  The large scale structure (LSS) of the universe is generated by the linear
density gaussian modes, which are evolved into the observed nonlinear LSS. The
posterior surface of the modes is convex in the linear regime, leading to a
unique global maximum (MAP), but this is no longer guaranteed in the nonlinear
regime. In this paper we investigate the nature of posterior surface using the
recently developed MAP reconstruction method, with a simplified but realistic
N-body simulation as the forward model. The reconstruction method uses
optimization with analytic gradients from back-propagation through the
simulation. For low noise cases we recover the initial conditions well into the
nonlinear regime ($k\sim 1$ h/Mpc) nearly perfectly. We show that the large
scale modes can be recovered more precisely than the linear expectation, which
we argue is a consequence of nonlinear mode coupling. For noise levels
achievable with current and planned LSS surveys the reconstruction cannot
recover very small scales due to noise. We see some evidence of non-convexity,
specially for smaller scales where the non-injective nature of the mappings:
several very different initial conditions leading to the same near perfect
final data reconstruction. We investigate the nature of these phenomena further
using a 1-d toy gravity model, where many well separated local maximas are
found to have identical data likelihood but differ in the prior. We also show
that in 1-d the prior favors some solutions over the true solution, though no
clear evidence of these in 3-d. Our main conclusion is that on very small
scales and for a very low noise the posterior surface is multi-modal and the
global maximum may be unreachable with standard methods, while for realistic
noise levels in the context of the current and next generation LSS surveys MAP
optimization method is likely to be nearly optimal.
