You are an expert at summarizing long articles. Proceed to summarize the following text: 
it is currently admitted that processes of structure formation in the universe are mainly dominated at large scale by an invisible component called dark matter . although the microscopic origin of dark matter is still unclear , it is macroscopically described as a self - gravitating collisionless fluid following the collisionless boltzmann or vlasov equation in a cosmological background , @xmath0 f({\mbox{\boldmath$x$}},{\mbox{\boldmath$v$}},t ) = 0 , \label{eq : vp_3d}\end{aligned}\ ] ] supplemented with the poisson equation for the newton potential @xmath1 , @xmath2 , \label{eq : poisson_3d}\end{aligned}\ ] ] where @xmath3 is the phase - space density at comoving position @xmath4 , peculiar velocity @xmath5 and time @xmath6 , @xmath7 is the expansion factor of the universe , @xmath1 is the gravitational potential and @xmath8 is the average dark matter density .    in the standard picture of structure formation , dark matter was initially cold i.e. with a virtually null local velocity dispersion , so the six - dimensional phase - space distribution is effectively reduced to a three - dimensional hyper - surface and this remains true at all times thanks to the hamiltonian nature of the system . at early times , dark matter thus follows the single - stream flow regime , with a velocity field @xmath9 uniquely determined as a function of position , and its evolution is that of a pressure - less fluid with a phase - space distribution function given by @xmath10 , \label{eq : df_single - stream}\end{aligned}\ ] ] where @xmath11 is the density contrast of the dark matter distribution , initially of very small magnitude , as well as @xmath9 .    substituting eq . ( [ eq : df_single - stream ] ) into eqs . ( [ eq : vp_3d ] ) and ( [ eq : poisson_3d ] ) and taking the zeroth and first velocity moments yields the eulerian formulation of large - scale structure dynamics @xmath12=0 , \label{eq : continuity}\\ & \frac{\partial{\overline { \mbox{\boldmath$v$}}}}{\partial t}+\frac{1}{a}\,\bigl\{{\overline { \mbox{\boldmath$v$}}}\cdot\nabla\bigr\}\overline{\mbox{\boldmath$v$}}=-\frac{1}{a}\,\nabla\phi , \label{eq : euler}\\ & \frac{1}{a^2}\nabla^2\phi=4\pi\,g\,\overline{\rho}_{\rm m}\,\delta . \label{eq : poisson_v2}\end{aligned}\ ] ] although the single - stream flow is , strictly speaking , valid only during the early phase of structure formation , the above equations have been shown in practice to describe accurately nonlinear mode - coupling in the weakly nonlinear regime and provide a solid basis for perturbation theory ( pt ) calculations to predict statistical quantities of large - scale structure , such as the power spectrum or the two - point correlation function of the matter distribution ( see , e.g. * ? ? ? with the advent of cosmological observations aiming at precisely mapping galaxy clustering at large scales , techniques have been developed to improve on the slow convergence of the perturbative expansion as well as to accelerate higher - order calculations ( e.g. , * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? . improved calculations involving the next - to - next - to - leading order called two - loop give a remarkable agreement with cosmological @xmath13-body simulations at weakly nonlinear scales and have been applied to observations . however , this does not ensure the convergence of perturbative calculations including higher - order corrections , e.g. , three - loop . in fact , a direct calculation at three - loop order suggests a very large uv contribution to the large - scale modes through nonlinear mode - coupling @xcite and indicates a break down of higher - order perturbative expansions even at large scales . deficiency of pt calculations has been also highlighted in recent numerical analyses . @xcite directly measured the coupling between different scales in cosmological @xmath13-body simulations and found that the actual contribution from small scales to the large - scale modes is suppressed , as opposed to the prediction of pt based on the single - stream approximation . these facts imply that the validity of the single - stream treatment is questionable even at large scales , and higher - order perturbative corrections need to be remedied with a proper account of small - scale dynamics , where the multi - stream flow contribution is important .    one way to account for small - scale dynamics consists of using some ansatz to summarize the main physical effects in the multi - stream regime , such as burgers equation ( e.g. , * ? * ; * ? ? ? * ) and various alternatives ( e.g. , * ? ? ? * ) , which roughly amount to adding a source term in the right hand side of equation ( [ eq : euler ] ) . in particular , it was proposed recently to use effective fluid equations giving account of the non - vanishing stress tensor arising in eq . ( [ eq : euler ] ) when calculating the first velocity moment of vlasov equation in the multi - stream regime . this _ effective - field theory _ approach has attracted a lot of interest and has been studied in detail ( e.g. , * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? the drawback of this approach is however that the parameters in the stress tensor characterizing the small - scale dynamics need to be calibrated with @xmath13-body simulations to keep predictions of perturbative calculations under control . furthermore , these parameters generally vary with cosmology and with time and no prediction with pt is really possible in this framework independently of @xmath13-body simulations . an alternative approach that we consider in this paper consists in going back to a more fundamental description , i.e. the vlasov - poisson system , eqs . ( [ eq : vp_3d ] ) and ( [ eq : poisson_3d ] ) , and trying to follow accurately the phase - space distribution function in the multi - stream regime , which is essential for describing the formation of dark matter halos . for this purpose , it is useful to employ lagrangian pt , where the small parameter is the displacement field . we shall consider one - dimensional gravitational dynamics in standard ( three - dimensional ) cosmology . in this case , large - scale structure dynamics is described by the gravitational interaction of massive parallel infinite planes moving left and right along a fixed axis , while hubble expansion is taking place as usual in all the directions following standard friedman - lematre equations . despite its simplicity , one dimensional dynamics defined as such still displays a rich physical content which somewhat shares the same features as 3d clustering . this is partly the reason why the 1d model has recently attracted much attention ( e.g. , * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? in particular , the zeldovich solution provides an exact solution for the dynamics of massive sheets before shell - crossing , and thus , starting with zeldovich solution , a tractable perturbative treatment of multi - stream flows is made possible based on a lagrangian description . the analysis in the present paper is an extension of the method developed in @xcite to the cosmological setup . we will describe perturbatively post - collapse dynamics around the shell - crossing region and apply it to several cases including random initial conditions . in addition , we shall present a novel regularization scheme reducing the impact of small - scale clustering and improving greatly our post - collapse pt predictions of large - scale structure statistics . the idea is to apply adaptive smoothing to initial density peaks and to better capture the bulk properties of phase - space structures in the post - collapse regions , where the interaction or merger of halos is supposed to be significant . the idea is similar to the peak - patch treatment proposed by @xcite and subsequent works ( see , e.g. * ? ? ? * ) , but we here implement it in the pt prescription in order to better describe the late - time post - collapse dynamics . indeed , our perturbative approach will not allow us to follow post - collapse dynamics beyond next - crossing time , although an iterative prescription such as proposed by @xcite could make this possible but is out of the scope of the present work . here , we will show instead that the predictions with adaptive smoothing reproduce simulations remarkably well even at nonlinear scales , where zeldovich solution significantly deviates from simulations . note finally that an analytical study in 1d of course represents only a first step toward a proper description of 6d phase - space dynamics . apart from a few examples including self - similar solutions @xcite , little is indeed known analytically for the phase - space dynamics of vlasov - poisson systems , particularly in cosmology . therefore , further development of analytical treatment in this framework seems indispensable , complementary to simulations , and even helpful to cross check simulation codes . this is all the more motivated by the fact that simulations in 6d phase space have now become available thanks to recent efforts on the development of pure vlasov codes @xcite . this paper is organized as follows . in [ sec:1dcosmology ] , we begin by describing the basic setup of our calculations in one - dimensional cosmology . we then discuss in [ sec : postcollapse ] the analytic treatment beyond shell - crossing and develop post - collapse pt . [ sec : improvement ] introduces a novel filtering scheme to improve the perturbative description of post - collapse dynamics , suited for a system with interacting clusters . in   [ sec : comparison_nbody ] , analytic calculations in the post - collapse perturbative framework are tested against controlled @xmath13-body experiments . several cases are considered , including the single initial sine wave and random initial conditions . finally , [ sec : conclusion ] is devoted to discussion and conclusion . we consider an ensemble of massive parallel infinite planes moving along @xmath14 axis and interacting through gravitational force in the expanding universe . the lagrangian equations of motion of the planes are given by @xmath15 where @xmath16 and @xmath17 are respectively the comoving position and peculiar velocity of each plane , @xmath18 the gravitational potential , @xmath8 the average matter density , @xmath11 the density contrast and @xmath7 the expansion factor of the universe . to simplify the equations , it is useful to introduce the _ super - conformal time _ @xmath19 defined by ( e.g. , * ? ? ? * ; * ? ? ? * ) , @xmath20 the new velocity @xmath21 and potential @xmath1 : @xmath22 then , eqs . ( [ eq : eom0_x])([eq : eom0_delta ] ) are simplified : @xmath23 where @xmath24 is the matter density parameter of the universe and @xmath25 is the hubble constant .    to deal with this system , especially for the dynamics after shell - crossing , it is useful to introduce in the cold case considered here the lagrangian coordinate @xmath26 defining the initial position of the planes in the absence of perturbation @xmath27 and to express the subsequent position and velocity of each mass element as @xmath28 and @xmath29 . in particular , mass conservation implies @xmath30\,{\rm d}x \,\,\longrightarrow\,\ , \delta(x)=\left(\frac{\partial x}{\partial q}\right)^{-1}-1 . \label{eq : jacobian_delta}\end{aligned}\ ] ] this equation is valid in the single - stream regime , i.e. as long as @xmath31 remains monotonic as a function of @xmath26 .    with the new expressions above , the solution of the equations of motion can be formally written as : @xmath32 , \label{eq : sol1_p}\end{aligned}\ ] ] where @xmath33 and @xmath34 are the initial positions and velocities given at a starting time @xmath35 , which will be specified below .    in what follows , we consider the dynamics of the cosmological system in a finite - size box with periodic boundaries , @xmath36 . the solution of poisson equation , eq . ( [ eq : eom1_delta ] ) , can be expressed in an integral form as : @xmath37\,\delta(x ' ) . \label{eq : phi_poisson}\end{aligned}\ ] ] the derivation of this expression is presented in appendix [ sec : green_func ] . then , the force exerted on a mass element at position @xmath14 is given by : @xmath38 , \label{eq : force_poisson}\end{aligned}\ ] ] where function @xmath39 represents the heaviside step function . in the above , we used the fact that fluctuations averaged over space vanish , @xmath40 . in the limit @xmath41 , equation ( [ eq : force_poisson ] ) of course converges to the well - known result in the infinite space .      in one - dimension and in the cold case , zeldovich approximation is known to provide the exact solution for the dynamics of mass elements before shell - crossing . it can be explicitly written as @xmath42 here , function @xmath43 corresponds to the linear growth factor satisfying the following equation : @xmath44\,d_+(\tau)=0 . \label{eq : evolv_lin}\end{aligned}\ ] ] note that in terms of cosmic time @xmath6 , eq .  ( [ eq : evolv_lin ] ) reduces to the standard form of linear evolution equation : @xmath45\,d_+(t)=0.\end{aligned}\ ] ] the zeldovich solution in eq . ( [ eq : zeldovich_sol ] ) contains an arbitrary function @xmath46 that we call displacement field . it is related to the linear density field @xmath47 given at a very early time ( @xmath48 or @xmath49 ) through @xmath50 since zeldovich solution is exact until shell - crossing , we do not necessarily assume that the evolved density field @xmath51 is small . one may thus consider the situation where in a region around a lagrangian coordinate @xmath52 , the density field becomes large and the region will undergo shell - crossing at time @xmath53 . recalling the fact that in 1d , a shell - crossing point corresponds to an inflection point of the mapping from lagrangian to eulerian space , the actual conditions for shell - crossing are given by @xmath54 at the time of shell - crossing , @xmath53 , the solution ( [ eq : zeldovich_sol ] ) around the shell - crossing region can be expanded as follows : @xmath55 using eq . ( [ eq : def_deltal ] ) , the conditions ( [ eq : condition_shellcrossing ] ) for shell - crossing can be rewritten as @xmath56 in other words , shell - crossing takes place at local density peaks . hence , in the 1d case , the conditions for shell - crossing are equivalent to peak constraints and collapse take places exactly when the eulerian linear density contrast at the peak positions becomes equal to unity . the structure of the system around the density peak is symmetric with respect to the lagrangian coordinate @xmath57 , where @xmath52 is the lagrangian position of the shell - crossing point ( see also * ? ? ? [ fig : shellcrossing],width=321 ]    we are interested in the dynamics of mass elements after shell - crossing , when the zeldovich solution is no longer valid . in this section , extending the work of @xcite , we develop perturbative calculations to deal with the multi - stream motion around the shell - crossing region . the basic formalism to treat post - collapse dynamics can be described as follows . starting with the cold initial conditions in sec . [ sec : initial ] , we first follow pre - collapse dynamics with the exact zeldovich solution . then , in regions undergoing shell - crossing , we switch to a perturbative treatment and compute the backreaction to the zeldovich flow based on an explicit functional form of the displacement field in the shell - crossing region . to be precise , we compute the force from eq . ( [ eq : force_poisson ] ) exerted at each position using the extrapolation of the zeldovich flow . integrating the force over time , we obtain the correction for the velocity to the zeldovich motion from eq . ( [ eq : sol1_p ] ) . further integrating the corrected velocity over time , a correction for the position is obtained from eq . ( [ eq : sol1_x ] ) . throughout these calculations , we assume that the collapsing region is small , which allows us to express the phase - space configuration in terms of polynomial forms of low order of the lagrangian coordinate .    let us focus on a collapse point centered on lagrangian position @xmath52 and perform a perturbative description of the multi - stream flow around @xmath58 . as we mentioned earlier , the conditions for shell - crossing are given by eq . ( [ eq : condition_shellcrossing ] ) . shortly after shell - crossing time @xmath53 , the multi - stream flow has just started growing in a small lagrangian region around @xmath52 and the displacement field @xmath46 can be approximately described by a third - order polynomial in @xmath59 . since the deviation from zeldovich flow is small , the motion of a mass element around @xmath58 may be expanded at third - order in @xmath26 for @xmath60 as @xmath61 with the time - dependent coefficients @xmath62 , @xmath63 and @xmath64 defined by @xmath65 where the prime denotes the derivative with respect to @xmath26 . here , we used eqs . ( [ eq : def_deltal ] ) and ( [ eq : shellcrossing_peak ] ) . the above expressions imply that the shell - crossing point slightly moves from @xmath66 to @xmath67 and the shell - crossing structure develops as an `` @xmath68 '' shape , as shown in fig . [ fig : shellcrossing ] . we now introduce the new lagrangian coordinate , @xmath57 , for which the structure of the local density peak is symmetric . since we expand the displacement field at third - order , the equation @xmath69 inside the shell - crossing region has three solutions . an example is shown in fig . [ fig : shellcrossing ] for the equation @xmath70 with solutions @xmath71 . in general , the three ordered roots of the equation @xmath72 , when they exist , are related to each other through the following expressions : @xmath73 here , the lagrangian extent @xmath74 of the boundary of the multi - valued region , also related to the lagrangian position @xmath75 of the left caustic in fig . [ fig : shellcrossing ] through @xmath76 , is expressed as @xmath77 because of the time dependence of the coefficients ( see eqs . [ eq : def_b ] and [ eq : def_c ] ) , the position of the boundary gradually changes in time . the leading - order expression in powers of time @xmath78 reads @xmath79 with @xmath80 in other words , it also means that a given element of fluid of lagrangian position @xmath81 will enter the multi - valued region at some time @xmath82 defined by @xmath83 and then subsequently coincide exactly with the lagrangian position of the caustic at a time @xmath84 with @xmath85 . by inverting the relations @xmath86 and @xmath87 , we obtain @xmath88 and @xmath89 the quantities defined above thus play a crucial role to disentangle the single - stream regime from post - collapse dynamics , as well as the inner part of the `` @xmath90 '' shape of fig . [ fig : shellcrossing ] from its outer part , the lagrangian position of the caustic acting as a separator in the latter case . the expression of the force will indeed be different between the zeldovich single - flow regime @xmath91 and the multi - valued regime composed of an outer caustic region , @xmath92 , and an inner caustic region , @xmath93 . one important remark for the subsequent calculations is that the boundaries coordinates @xmath74 and @xmath94 are assumed to be small and that the backreaction to the position and velocity inside the multi - valued region can be described perturbatively in a polynomial form of @xmath95 [ partly allowing the fractional power of @xmath96 . thus , the resultant expressions are , rigorously speaking , only valid for a short period after collapse time , but we shall see in practice , as also shown by @xcite in the non - cosmological case , that they remain impressively accurate even up to next - crossing time ( see sec .  [ sec : comparison_nbody ] ) , which will allow us to set up the framework for a powerful self - adaptive scheme . to derive the corrections to the motion , we first compute the force exerted on a mass element inside the multi - valued region , @xmath97 ( see fig .  [ fig : shellcrossing ] ) . note again that the outer regions @xmath98 and @xmath99 are described by the zeldovich solution . the force in the multi - valued region is computed using eq . ( [ eq : force_poisson ] ) , dividing each integral of the right - hand - side into three contributions : @xmath100 assuming that the collapse region , @xmath101 , is small enough , the contributions to the integrals from each domain can be computed analytically , based on the geometrical setup in fig . [ fig : shellcrossing ] . the detailed calculations are presented in appendix [ sec : derivation ] . summing up all the contributions given in eqs . ( [ eq : integral_inner1 ] ) , ( [ eq : integral_inner2 ] ) , ( [ eq : integral_outer1 ] ) and ( [ eq : integral_outer2 ] ) , the force exerted on the mass element at @xmath102 inside the multi - valued region becomes @xmath103 \label{eq : force_inner}\end{aligned}\ ] ] with the functions @xmath104 and @xmath105 respectively defined by @xmath106 and @xmath107 where the quantities @xmath62 , @xmath63 , and @xmath64 are defined by eqs . ( [ eq : def_a])([eq : def_c ] ) . note that in deriving eq . ( [ eq : force_inner ] ) , we have assumed that the system follows zeldovich solution if @xmath108 . since the resultant expressions are written in terms of the local quantities characterizing the density peak at position @xmath52 and the shell - crossing time @xmath53 , eq . ( [ eq : force_inner ] ) is in fact still applicable to other shell - crossing regions possibly appearing elsewhere in the region @xmath108 , allowing one to generalize the result to a smooth random initial field with multiple peaks . one issue , discussed later , is then to treat mergers , i.e. the case when for instance two multi - stream regions start overlapping . given the explicit expression for the force in the multi - stream region and using the formal solution given by eqs . ( [ eq : sol1_x ] ) and ( [ eq : sol1_p ] ) , we now compute corrections to the zeldovich flow that we write as follows : @xmath109 we noticed in previous section that , depending on the lagrangian position of interest , the expression for the force is different and we have to divide the domain of the integrals in eqs . ( [ eq : delta_p_formal ] ) and ( [ eq : delta_x_formal ] ) into several pieces :    1 . @xmath110 : the position @xmath81 is located in the single - valued region ( i.e. , @xmath111 ) and the motion is still described by the zeldovich solution . we have @xmath112    1 . @xmath113 : the position @xmath81 is in the multi - valued region and satisfies @xmath114 , i.e. lies in the outer part of the caustic . in addition to the zeldovich flow , the corrections arising from the multi - stream flow need to be added and we note them as follows : @xmath115 in this region , which correspond to the tails of the @xmath90 shape in fig . [ fig : shellcrossing ] , the system is globally expanding in phase space and gaining energy , at variance with the central part , which on the contrary , contracts . this can be understood from computing the variation of energy of a typical test particle during a fraction of orbit . this phenomenon was studied in the non - cosmological case by @xcite and the results should not be fundamentally different in the cosmological case studied here .    1 . @xmath116 : this corresponds to @xmath117 , i.e. the position @xmath81 now lies in the inner part of the multi - stream region with respect to the caustic . similarly to the above case , the backreaction to zeldovich flow needs to be computed , including at present both the multi - stream dynamics of the inner part and the incoming flow from the outer part . we may write @xmath118 note that the corrections @xmath119 and @xmath120 partly come from outer part contributions , @xmath121 and @xmath122 . as mentioned in the previous point ( ii ) , in this region , which corresponds to the central part of the @xmath90 on fig . [ fig : shellcrossing ] , the system is globally contracting and losing energy in favor of the tails of the @xmath90 .    in what follows , we compute the backreaction to zeldovich flow and derive the expressions for @xmath123 and @xmath124 in each domain . the calculation for the corrected motion is rather straightforward but needs several steps . here , we give a brief sketch of the calculation , deferring details to appendix [ sec : pcpt_basic ] . the final results are given in eqs . ( [ eq : p_for_tau_c > tau > hattau_c ] ) and ( [ eq : x_for_tau_c > tau > hattau_c ] ) for the outer part , eqs .  ( [ eq : p_for_tau > tau_c ] ) and ( [ eq : x_for_tau > tau_c ] ) for the inner part , together with the coefficients in tables  [ tab : xv_pcpt1 ] and [ tab : xv_pcpt2 ] . note that although our expansion is rigorously valid only at third - order in @xmath81 , higher - order contributions will appear in the corrections to provide a continuous solution , up to fifth- and seventh - order for the velocity and the position , respectively . cll coefficients & velocity @xmath125 [ eq .  ( [ eq : p_for_tau_c > tau > hattau_c ] ) ] & position @xmath126 [ eq . ( [ eq : x_for_tau_c > tau > hattau_c ] ) ] +   + @xmath127 & @xmath128 & @xmath129 +   + @xmath130 & @xmath131 & @xmath132 +   + @xmath133 & @xmath134 & @xmath135 +   + @xmath136 & @xmath137 & @xmath138 +   + @xmath139 &    ' '' ''    & @xmath140 +   + @xmath141 & @xmath142^\tau_{\widehat{\tau}_{\rm c}(q)}}$ ] & @xmath143 +      in the outer part of the multi - valued region , ( ii ) , the correction to the velocity can be expressed as @xmath145 in the above , while the first integral is performed with the help of formulae in appendix [ sec : integral ] , the second integral is computed exactly . the whole derivation is provided in appendix [ subsec : delta_v_out ] . the resultant expression can be summarized as @xmath146 + \ , \widetilde{\epsilon}_1(\tau,\widehat{\tau}_{\rm c } ) , \label{eq : p_for_tau_c > tau > hattau_c}\end{aligned}\ ] ] with the time - dependent coefficients given in table  [ tab : xv_pcpt1 ] . note that the coefficient @xmath147 implicitly depends on the lagrangian position @xmath81 through @xmath148 but is not taylor expanded with respect to @xmath81 , for simplicity . a fully analytical theory , in particular to predict the power spectrum of the projected density field , would in principle require such a taylor expansion . performing it should not change significantly the performances of post - collapse pt as presented in this article . once @xmath125 is obtained , the expression for the correction @xmath126 is evaluated by further integrating eq . ( [ eq : p_for_tau_c > tau > hattau_c ] ) over time : @xmath149 the resultant expression becomes ( see appendix [ subsec : delta_x_out ] for derivation ) : @xmath150 \nonumber\\ & \quad\quad + \,\widetilde{\epsilon}_2(\tau,\widehat{\tau}_{\rm c } ) , \label{eq : x_for_tau_c > tau > hattau_c}\end{aligned}\ ] ] with the time - dependent coefficients presented in table  [ tab : xv_pcpt1 ] . note again the dependence on @xmath81 of parameter @xmath151 through @xmath82 and the fact that we did not taylor expand it in polynomials of @xmath81 , as it would be required for a fully analytical theory . let us consider now the inner part of the multi - valued region ( iii ) . for the correction to the velocity , @xmath153 , the expression to be evaluated is divided into three contributions : @xmath154 . \label{eq : delta_p_2 } \end{aligned}\ ] ] evaluating each contribution in appendix [ subsec : delta_v_in ] , we obtain @xmath155\,+\ , \widetilde{\epsilon}_1(\tau , \widehat{\tau}_{\rm c } ) , \label{eq : p_for_tau > tau_c}\end{aligned}\ ] ] with the time - dependent coefficients given in table  [ tab : xv_pcpt2 ] . note that the coefficient @xmath147 in the above equation is the same one as in eq . ( [ eq : p_for_tau_c > tau > hattau_c ] ) . cll coefficients & velocity @xmath153 [ eq .  ( [ eq : p_for_tau > tau_c ] ) ] & position @xmath119 [ eq .  ( [ eq : x_for_tau > tau_c ] ) ] +   + @xmath127 & @xmath156 & @xmath157 +   + @xmath130 & @xmath158 & @xmath159 +   + @xmath136 & @xmath160 & @xmath161 +   + @xmath139 &    ' '' ''    & @xmath162 +   + @xmath141 & @xmath142^\tau_{\widehat{\tau}_{\rm c}(q)}}$ ] & @xmath143 +    similarly , the correction to the position , @xmath119 , is divided into three pieces : @xmath163 where the function @xmath164 represents the correction of the velocity given in either eq . ( [ eq : p_for_tau_c > tau > hattau_c ] ) or eq . ( [ eq : p_for_tau > tau_c ] ) , but without the contribution corresponding to the zeldovich motion , i.e. , the last term in their expressions . the term @xmath165 is the same quantity as in eq . ( [ eq : x_for_tau_c > tau > hattau_c ] ) ( see table  [ tab : xv_pcpt1 ] ) .    based on the calculation presented in appendix [ subsec : delta_x_in ] , we obtain the expression for @xmath119 : @xmath166\ , + \,\widetilde{\epsilon}_2(\tau,\widehat{\tau}_{\rm c } ) , \label{eq : x_for_tau > tau_c}\end{aligned}\ ] ] with the time dependent coefficients given in table  [ tab : xv_pcpt2 ] . the perturbative description of post - collapse dynamics introduces two kinds of approximations around the singularities ( shell - crossing points ) . firstly , we assumed the position and velocity to be described by third - order polynomials in @xmath81 in the multi - stream regions . to preserve continuity between the multi - stream and the single - stream regime , we added higher - order terms in @xmath81 to the solutions for @xmath14 and @xmath167 , but strictly speaking , these solutions are valid locally only at third - order in @xmath81 around the singularities . secondly , these solutions are correct only very shortly after collapse time @xmath53 up to second - order in @xmath168 . our post - collapse dynamics can thus be improved at a twofold level :    ( a ) : :    _ improvement at the spatial level : _ we could try , starting from    initial conditions , to use a description of the @xmath90    shape of fig . [ fig : shellcrossing ] at higher - order in @xmath81 ,    e.g. fourth- or fifth- instead of third - order . in this latter case , deriving the solution of equation @xmath70 comes down    to find the roots of a third or fourth - order polynomial instead of a    second - order one , which makes , in addition to the symmetry breaking of    the @xmath90 shape , the calculations much more    cumbersome , even though a series expansion of the roots at fourth or    fifth - order in @xmath81 remains tractable and might bring    interesting improvements on the description of the tails of the    @xmath90 shape . ( b ) : :    _ higher - order in time : _ we can introduce higher - order corrections in    time to follow more accurately the evolution of the central part of    the singularity after collapse . note that after integration of the    equations of motion , such corrections bring out higher - order terms in    @xmath81 . strictly speaking , the right way to proceed would    consist of performing `` post - post - collapse '' dynamics , i.e. of    computing a counter - term to the zeldovich+leading - order post - collapse    solution . this procedure is cumbersome , and seems to make sense only    if performed simultaneously with a description at higher - order in    @xmath81 , in particular because the tails of the    @xmath90 structure in fig . [ fig : shellcrossing ] also    influence the dynamics of the central part of the system ( but this    remains to be proved ) . going beyond leading - order in post - collapse dynamics thus seems rather challenging even in the simple one - dimensional case considered in this article . however , we can still try to investigate some potential improvements of the theory while maintaining the complexity of the calculations at an acceptable level . for this purpose , we consider here mainly two approaches . the first one just consists in incorporating next - to - leading order terms into the expressions for @xmath74 and @xmath169 ( and @xmath75 and @xmath170 ) . these quantities , which determine the boundary between the single - stream and the multi - stream regime , are indeed crucial because they correspond to the main scale / timescale of the system . the second approach consists , in addition , of improving the evaluation of the time integrals in eqs . ( [ eq : delta_p_formal ] ) and ( [ eq : delta_x_formal ] ) . in sec . [ subsec : corrections_pcpt ] , when evaluating these time integrals , a part of the integrands is taylor - expanded and the results are presented at leading - order in time . it is thus possible to include higher - order corrections from the taylor - expanded integrands . again , one has to be aware that the two corrections we propose might not improve the results , since a correct higher - order time approach , as discussed in point ( b ) above , would actually require to compute counter terms to the standard post - collapse result , which we do not do here . in fact , our higher - order corrections might even _ worsen _ the results since the corrections brought by such additional terms might have opposite sign compared to the iterated counter - term approach . the derivation of the higher - order corrections we propose here and summarize just below is rather straightforward compared to its rigorous alternatives ( a ) and ( b ) above , but still cumbersome . the resultant expressions now include higher - order polynomials of @xmath81 as well as new contributions , which are given in appendix [ sec : pcpt_higher - order ] . to be more specific , we consider the following cases :    1 .   higher - order for the expressions of critical time @xmath171 and position @xmath172 ( _ hc _ ) : based on the basic results in sec . [ subsec : corrections_pcpt ] , we include the contributions arising from the next - to - leading order expressions for @xmath74 and @xmath169 ( as well as @xmath75 and @xmath173 ) . the calculation to derive higher - order terms is almost the same as presented in sec . [ subsec : corrections_pcpt ] , but we replace the leading - order expressions for @xmath74 and @xmath169 in eqs . ( [ eq : qchat_approx ] ) and ( [ eq : hattau_c ] ) with those including the next - to - leading order , given in eqs . ( [ eq : higher_order_qhat ] ) and ( [ eq : higher_order_hattau_c ] ) . the higher - order results , denoted by @xmath174 and @xmath175 , correspond to the expressions summarized in appendix [ sec : pcpt_higher - order ] , but where the terms involving the quantities @xmath176 or @xmath177 are all set to zero .    1 . higher - order in time ( _ ho _ ) : in addition to the above treatment , we also include higher - order corrections to the time - integrals , based on the formulas given in appendix [ sec : integral ] . the resultant expressions , @xmath178 and @xmath179 , are summarized in appendix [ sec : pcpt_higher - order ] .    finally , we also consider a semi - analytic approach in order to smoothly connect the inner part of the multi - stream region to the zeldovich solution :    1 . third - order spline ( _ spl _ ) : basic post - collapse pt is applied to the inner part , @xmath152 , but we take the higher - order expression for @xmath75 . in the outer part , @xmath180 , on the other hand , a third - order spline interpolation is used to smoothly connect the inner post - collapse region to the zeldovich solution , with again , the higher - order expression for @xmath74 . post - collapse pt as described in [ sec : postcollapse ] provides a way to describe the dynamics around an initial density peak shortly beyond collapse time . in practice , the quality of the description is expected to be good until next - crossing time , at least this is what is known in the non - cosmological case @xcite . using the leading - order expression given in eq . ( [ eq : x_for_tau > tau_c ] ) , the duration of time spent between collapse and next shell - crossing can be estimated from eq . ( [ eq : x_post - collapse_inner ] ) , @xmath181 which is derived by solving the equation @xmath182 .    beyond time @xmath183 , post - collapse pt is expected to fail , similarly as zeldovich approximation fails beyond @xmath184 . the only way to fix this problem at the fine grain level would be to use an iterative procedure similar to what was proposed by @xcite to reconstruct the spiral structure built up during the course of the dynamics . in addition , post - collapse pt assumes the initial density distribution to be sufficiently smooth . while this is certainly expected at very small scales in the standard cold dark matter scenario , the issue of mergers arises . indeed , post - collapse pt is valid only if the forming halo is isolated . if there are two adjacent collapsing regions , they might , at some point , overlap in lagrangian space , hence merge with each other in eulerian space . after the merger , a separate description of each halo becomes inadequate . in this case , the prediction of post - collapse pt , which is in principle able to handle multiple but non - overlapping structures in lagrangian space , will break down much earlier than expected . one simple way to fix both limitation in time and mergers issues , already applied in the past to the zeldovich approximation ( see , e.g. * ? ? ? * ; * ? ? ? * and references therein ) , consists in trying to describe the dynamics at the coarse - grained level by introducing a smoothing procedure in lagrangian space , so that post - collapse pt remains applicable , in particular that the condition @xmath185 remains valid for the system with coarse - grained initial conditions . this means that single or composite halos as well as mergers are locally summarized by a simplified @xmath90 shape ( fig . [ fig : shellcrossing ] ) corresponding to a single peak of the coarse - grained initial density field and such that @xmath185 , while retaining the correct large - scale dynamics . the hope is that this simplified structure still captures the most important features of the nonlinear system , e.g. its mass and its typical size .    of course , the coarsening scale should depend on environment , hence on lagrangian position @xmath81 , so for best results , _ adaptive smoothing _ of initial conditions has to be introduced . in this sense our method is very similar to the `` peak - patch '' treatment by @xcite ( see also * ? ? ? * ) , but we implemented it in the post - collapse pt framework . to be precise , the practical set up of post - collapse dynamics with adaptive smoothing is summarized as follows : the initial density field is smoothed at various scales , employing the sharp-@xmath186 filter function in fourier space : @xmath187 with @xmath39 being the heaviside step function . the cutoff wavenumber @xmath188 ranges from the fundamental mode @xmath189 to a certain value @xmath190 , where @xmath191 is an integer corresponding to the effective resolution of the simulation , as discussed more in details in [ subsec : simulation ] . in our investigations , we covered all the possible values of @xmath188 , i.e. @xmath192 with @xmath193 , but this is costly and probably not required , in practice .    1 . zeldovich solutions are constructed from the smoothed density fields . then we identify the critical points @xmath52 in lagrangian space that have undergone shell - crossing at time @xmath194 , where @xmath19 is the time of interest for the analyses . for each critical point , various quantities needed for perturbative calculations , such as @xmath195 and @xmath196 , are also computed .    1 . starting from the largest smoothing scale ( i.e. @xmath197 ) , we compute the next - crossing time @xmath198 associated to each critical point . for the points satisfying @xmath199 with @xmath200 a constant parameter , the post - collapse correction to the zeldovich flow is obtained in the interval @xmath201 $ ] , which is labeled .    1 .   the procedure ( iii ) is repeated , but with a smaller smoothing scale , i.e. , larger @xmath188 , and only for peaks lying in the unlabeled regions , until the cutoff wavenumber reaches @xmath202 . note that if two lagrangian regions overlap , e.g. @xmath203 and @xmath204 , corresponding to two different values of the smoothing scale , @xmath205 , the interval @xmath204 corresponding to the smallest smoothing scale , i.e. largest wavenumber @xmath206 , prevails . this is the way we solve the cloud - in - cloud problem , but not the only possible choice .    1 .   at the end of the iterative process , post - collapse dynamics is applied as well for the rest of the critical points corresponding to the highest resolution scale @xmath207 that still require treatment , i.e. for those with @xmath208 . the procedure is the same as described in point ( iv ) above , i.e. only critical points which do not fall in a labeled region are treated , and if it is the case , then their lagrangian interval @xmath209 prevails for applying post - collapse dynamics . for the rest of the unlabeled regions , standard zeldovich dynamics is applied and the mapping uses the initial density corresponding , of course , to the highest resolution scale . note that the approximation of the phase - space density distribution constructed this way is generally discontinuous , because we collect perturbative solutions for the displacement field in phase space coming from different coarse - graining scales without imposing smoothness at the transition between these solutions . in this respect , the procedure given above is not entirely consistent and remains to be improved in the future . nevertheless , the discontinuities only affect the results at small scales and in regions where the projected density is small . indeed , the discontinuities correspond to transitions between tails of the @xmath90 shape representing a halo and a non collapsed region or transitions between tails of two @xmath90 shapes in the case of a merger . as illustrated below through numerical examples , these regions of phase space have small density contrast and thus contribute weakly to second- or higher - order statistics , hence to the power spectrum which is the statistics under focus in this work . with an appropriate choice of @xmath210 , we shall show that the adaptive smoothing procedure provides a substantial improvement on the post - collapse prediction for the overall phase - space structure of the system and the power spectrum of the projected density . we are now in position to compare the predictions of post - collapse pt with simulation results . in   [ subsec : simulation ] , we briefly describe the one - dimensional cosmological @xmath13-body code we wrote to perform the numerical analyses of this article and discuss various important parameters settings . we then present the results of the comparison between theory and numerical experiments in various cases with cold initial conditions : formation of an isolated single halo (   [ subsec : single_cluster ] ) , merging of two neighboring halos (   [ subsec : merging_clusters ] ) , as well as random initial conditions with cold dark matter - like ( cdm ) and power - law power spectra (   [ subsec : cdm_initial ] and [ subsec : power - law_initial ] ) . in the later case , we shall discuss more in details highly nonlinear dynamics , in particular predictions from stable clustering , the expected outcome on the power spectrum from halos with power - law density profiles , as well as the statistical properties of the phase - space sheet resulting from a caustic distribution . to solve numerically vlasov - poisson equations in one dimension , we wrote a fortran90 particle - mesh code , vlafroid . in this code , publicly available through the following web page , www.vlasix.org , poisson equation is solved on a periodic grid of fixed resolution @xmath211 using fast fourier transform , after projecting the @xmath212 particles on the grid using simple cloud - in - cell interpolation ( see , e.g. * ? ? ? * ) . the particles positions and velocities are updated at each time step using a standard predictor - corrector integration scheme . after the predictor step , @xmath213 poisson equation is solved and the force is estimated on the grid by direct derivation in fourier space . the acceleration @xmath214 of each particle is computed by linear interpolation of the force . this is followed by the corrector step : @xmath215 the value of the time step @xmath216 is slowly varying with time and is bounded by 3 constraints :    1 .   a limit on the relative variations of the expansion factor during a time step to have a good description of the linear growing modes at early stages of the simulation , @xmath217 2 .   a dynamical condition which states that the time step should be a small fraction of the smallest harmonic oscillator time scale ( see , e.g. * ? ? ? * ) , @xmath218 where @xmath219 is the maximum value of the projected density ( normalized to be unity in average ) , 3 . the standard courant - friedrichs - lewy ( cfl ) condition @xmath220 with @xmath221 the maximum magnitude of the velocity , stating that a particle should not travel more than a fraction of a grid cell size during a time step .    for all the runs performed in this article , our practical and safe choice of the time step parameters was : @xmath222 the evolution of single or composite objects as studied in [ subsec : single_cluster ] and [ subsec : merging_clusters ] does not pose any problem from the numerical point of view . indeed , we consider here only a small number of dynamical times , so we do not expect significant defects due to the discrete nature of the particle distribution ( for instance , such as described in * ? ? ? the case of a random field , as studied in [ subsec : cdm_initial ] for cdm and in   [ subsec : power - law_initial ] for power - law power spectra , is less trivial . indeed , one needs accurate numerical modeling of the dynamics of each single structure , which has to be sampled with many particles and be resolved with a sufficient number of grid elements . in particular , to avoid contamination by the discrete nature of the particle distribution , we require , in average , many particles per grid element : @xmath223 to avoid the formation of small structures that would be insufficiently resolved by the computational grid , there must be a coarse - graining scale below which the initial conditions are smooth . this also facilitates , of course , the calculation of ( post - collapse ) pt predictions . here we impose a cut - off on the power spectrum of initial conditions using the sharp-@xmath186 filter introduced in eq . ( [ eq : sharp - k ] ) , with @xmath224 in practice , a choice which works well and adopted in this paper is @xmath225 which means that @xmath202 is five time smaller than the nyquist frequency of the grid . since the cut - off scale on the power spectrum can affect the dynamics compared to the true , un - smooth system , we shall perform some convergence studies to determine the available dynamic range in which the results are not influenced by @xmath202 . let us first present the simplest case of the formation of a single structure seeded with the following sinusoidal density perturbation : @xmath226 we adopt the einstein - de sitter cosmology ( @xmath227 , @xmath228 with @xmath229 ) and the box size of the simulation is arbitrarily set to @xmath230 . the simulation was started at expansion factor @xmath231 with initial amplitude @xmath232 and was run using @xmath233 . pt analyses were performed on the fly as part of routines of the code , when needed . [ fig : single - cluster_xv_basic ] shows snapshots of the phase - space structure ( upper inserts ) and projected density profiles ( lower inserts ) obtained from the simulation ( red ) and post - collapse pt ( blue ) based on the basic prescription described in sec . [ subsec : corrections_pcpt ] . also , for reference , the zeldovich solution is plotted as a green dashed line on each panel . although the three curves are identical at early time , differences show up , as expected , after collapse , which occurs at @xmath234 . a closer look at the central part of the system reveals however that post - collapse pt captures the trends of the @xmath13-body simulation quite well and succeeds to reproduce the structure of the inner part of the system until next - crossing time , which occurs at @xmath235 . this is in marked contrast with zeldovich solution , which , without the back reaction contribution present in the post - collapse correction , overshoots the displacement of the phase - space sheet after first crossing time , in particular of the caustics . after next - crossing time , @xmath236 , the phase - space structure given by post - collapse pt goes away from the @xmath13-body result . nevertheless , the predicted projected density profile still provides a reasonably good approximation of the evolved density structure , which is clearly not the case for the zeldovich solution . [ fig : single - cluster_xv_variant ] presents results of post - collapse pt when including the higher - order corrections discussed in [ subsec : higher - order ] . while the three prescriptions examined here , namely _ hc _ , _ ho _ and _ spl _ , account for higher - order terms in different manners , they are all based on the force expression in eq . ( [ eq : force_inner ] ) , which is valid only shortly after collapse . indeed , as discussed in [ subsec : higher - order ] , our implementation of higher - order corrections is not performed in a fully rigorous way , that would require to compute a back reaction to post - collapse pt , similarly as post - collapse pt introduces a back reaction correction to zeldovich solution . it also does not propose an iteration procedure that would allow us , similarly as in @xcite , to follow the system behind next - crossing time . hence , all of the treatments break down at the end of the simulation . nevertheless , taking the higher - order corrections to @xmath74 and @xmath169 into account improves the prediction of the outer boundary of the multi - valued region , and the phase - space description becomes visually better than that of the basic pt prescription in fig . [ fig : single - cluster_xv_basic ] . even after next - crossing time ( @xmath236 ) , the predicted phase - space structure still remains consistent , although the size of the halo is prone to be over - estimated , leading to an extended density profile , the best result being obtained with the _ spl _ prescription . still , one can notice that if only the projected density is under consideration , basic post - collapse pt does better than any of its higher - order counter - part , except for _ spl _ which , not surprisingly , provides nearly equivalent results , since it is exactly the same as post - collapse pt in the center of the system . as far as statistics depending on the projected density are concerned , such as the power spectrum that will be studied later , we thus expect from the analyses of this single halo that the simplest prescription of post - collapse pt is going to perform best . we examine now the second simplest case of two merging structures in einstein - de sitter cosmology ( @xmath227 , @xmath237 and @xmath229 ) . the initial density contrast is given by @xmath238 . \label{eq : delta_init_two_clusters}\end{aligned}\ ] ] the constant @xmath239 is determined so as to satisfy the condition @xmath240 . with @xmath241 , @xmath242 and @xmath243 , we have @xmath244 . the overall amplitude is set to @xmath245 and the simulation , of box size unity ( @xmath230 ) , was started at @xmath231 and was run using @xmath233 . the results are shown in fig . [ fig : merging - cluster_xv_basic ] for basic post - collapse pt and fig . [ fig : merging - cluster_xv_variant ] for the higher - order variants of post - collapse pt . here , we focus on the merger phase which takes place after collapse of each individual structure . figs . [ fig : merging - cluster_xv_basic ] and [ fig : merging - cluster_xv_variant ] demonstrate how the adaptive smoothing introduced in sec . [ sec : improvement ] improves pt prescriptions to describe the overall phase - space structure after the merger . in each figure , the upper and lower panels show the results without and with adaptive smoothing , respectively . the free parameter controlling adaptive smoothing , @xmath210 , is set here to @xmath246 for post - collapse pt and @xmath247 for the zeldovich solution .    without adaptive smoothing , as soon as the two clusters cross each other ( @xmath248 ) , both post - collapse pt and zeldovich solution start failing to describe the real dynamics in the @xmath13-body simulation . indeed , while post - collapse pt tries to provide a local perturbative correction for each single cluster , it is clearly unable to account for nonlinear couplings taking place at larger scales , where it simply follows zeldovich motion . as a result , the location of multi - valued regions predicted by post - collapse pt , as well as their shape , largely deviates from the correct one . on the other hand , implementing adaptive smoothing in the pt calculations significantly improves the predicted phase - space structure , starting already from the beginning of the merger phase ( @xmath248 ) . after merging , dynamics of the central part of the system is now described by the displacement field derived from the smoothed initial density field . by construction , it summarizes the composite structure made of two sub - halos into a single cluster . while this is a rough approximation of the real dynamics , the description of the outer part of the system is substantially improved as well as the prediction for projected density profiles . introducing both adaptive smoothing and higher - order corrections to post - collapse pt further improves the results from the visual point of view ( fig . [ fig : merging - cluster_xv_variant ] ) , although the predicted size of the halo tends to be slightly overestimated . as already noticed in [ subsec : single_cluster ] , the best results are obtained for _ spl _ , which also gives , after adaptive smoothing , comparable projected density profiles to basic post - collapse pt . the above results demonstrate that post - collapse pt with adaptive smoothing is effective in capturing the main features of phase - space structures . strictly speaking , it does not give an accurate prescription for the fine structure of high - density regions , but it provides a way of regularizing or mitigating the impact of highly nonlinear dynamics , keeping the location and size of halos reasonably accurate . as we will see next , adaptive smoothing brings a drastic improvement on the prediction of power spectra for systems with random initial conditions . furthermore , the introduction of adaptive smoothing makes the analytic calculations much less sensitive to the choice of the small - scale cutoff in the initial conditions so that pt predictions become more robust . in these respects , the determination of the extension of multi - streaming regions for a given smoothing scale of initial conditions , as performed in step ( iii ) of   [ sec : improvement ] , is essential , making the choice of @xmath210 crucial . concerning the particular merger experiment studied here , @xmath249 and @xmath247 seem to be the best choices for tuning post - collapse pt and zeldovich solutions , respectively , and we shall adopt these values in subsequent section . we shall however discuss again about the choice of @xmath210 later . let us now consider random cdm initial conditions . although there is no realistic setup in 1d , a way of mimicking as well as possible the 3d case may be to consider initial conditions given by a random gaussian field with the following power spectrum ( see , e.g. * ? ? ? * ) : @xmath250 with @xmath251 being the 3d matter power spectrum of the initial density fluctuations . here we use for @xmath251 the transfer function of @xcite who wrote the corresponding fortran module that we included in vlafroid . we set the cosmological parameters to those of the base @xmath252cdm model determined by planck @xcite : total matter density parameter @xmath253 , cosmological constant @xmath254 , baryonic matter density parameter @xmath255 , hubble constant @xmath256kms@xmath257mpc@xmath257 , effective slope @xmath258 and root mean square deviation @xmath259 . through @xmath260 with @xmath261mpc , @xmath262 and @xmath263 . [ foot : toto ] ] the simulations were performed with a box size @xmath264mpc and started at initial redshift @xmath265 . convergence of simulations results was tested by varying the number of particles @xmath266 , the resolution @xmath211 of the grid used to compute the force and the cutoff wavenumber of the initial power spectrum , @xmath202 , introduced in   [ subsec : simulation ] . , @xmath267 and @xmath268 . while the highest resolution runs with @xmath269 only allowed us to test the @xmath13-body results and zeldovich solution without adaptive smoothing , lower resolution runs included full pt predictions with and without adaptive smoothing . ] here , we present the results of the analyses for @xmath270 , @xmath271 , and @xmath272 mpc@xmath257 . our convergence study shows , unless specified otherwise below , that the results can be trusted for @xmath273 mpc@xmath257 , i.e. in the available dynamic range shown in left panels of figs . [ fig : pk_snap_cdm_midres ] and [ fig : pk_snap_cdm_midres2 ] . to have sufficient statistics for measuring power spectra , we ran @xmath274 simulations with different random seeds .                    figs . [ fig : pk_snap_cdm_midres ] and [ fig : pk_snap_cdm_midres2 ] show the measured power spectra at different times ( left panel of each figure ) as well as snapshots of phase space for a given realization ( right panel ) . right panel of fig . [ fig : pk_snap_cdm_midres ] is supplemented with fig . [ fig : snap_qxv_midres ] , which displays in approximately the same scale interval , eulerian positions and velocities as functions of lagrangian coordinates . one can first notice that , in contrast to the 3d case , the amplitude of the power spectrum at small scales is not strongly enhanced in 1d , but damped . for the cdm like initial conditions considered here , function @xmath275 becomes nearly flat in the highly nonlinear regime , which is probably a consequence of stable clustering as discussed further in details in next section , where the case of scale - free initial conditions is considered .    in fig . [ fig : pk_snap_cdm_midres ] , the simulations results are compared to basic pt predictions ( solid blue ) and zeldovich solution ( dotted green ) , while , in fig . [ fig : pk_snap_cdm_midres2 ] , they are compared to various higher - order variants of post - collapse pt . note , that except for the solid line that we explain below , all power spectra predicted by post - collapse pt or zeldovich solution are computed directly from the particle distribution evolved according to the displacement field predicted by each theoretical model under scrutiny . so it is important to mention that we did not actually compute analytically power spectra in these cases , although it would be in fact possible . the thin black solid lines , on the other hand , provide the analytic results for the zeldovich solution without adaptive smoothing , @xmath276 ( e.g. , * ? ? ? * ; * ? ? ? * ; * ? ? ? * ) : @xmath277\ , ; \nonumber\\ & \qquad\qquad\ i(q)=\int_0^\infty \frac{{\rm d}p}{\pi}\cos(p\,q)\frac{p_{\rm 1d}(p)}{p^2 } , \label{eq : pk_za_analytic}\end{aligned}\ ] ] where @xmath278 is the _ un - truncated _ power spectrum , i.e. without suppression of power at @xmath279 . the black solid curves therefore differ slightly from the green dotted ones on left panel of fig . [ fig : pk_snap_cdm_midres ] , because of the effect of the cut - off at @xmath202 ( and not because of discreteness effects or finite grid resolution ) . while the cut - off wavenumber @xmath202 has , as expected , a noticeable effect on the un - smoothed zeldovich solution as well as the un - smoothed post - collapse pt prediction and its variants ( a property that we do not show but that we verified ) , its influence becomes negligible in the range of values of @xmath186 under consideration , both for the simulations and the theoretical models with adaptive smoothing . as expected , without adaptive smoothing , both post - collapse pt prescriptions and zeldovich solution deviate from the simulations as soon as the nonlinear scale becomes larger than the cut - off scale @xmath280 , which is already the case in upper insert of left panel of figs . [ fig : pk_snap_cdm_midres ] and [ fig : pk_snap_cdm_midres2 ] , i.e. at @xmath281 . still , at this redshift , post - collapse pt ( or its variants ) without adaptive smoothing is able to capture the main features of phase space better than zeldovich solution , although this is not obvious to decipher on upper left insert of right panel of figs . [ fig : pk_snap_cdm_midres ] and [ fig : pk_snap_cdm_midres2 ] . at lower redshift , halos relax to a highly nonlinear state and can merge together , which can not be described well by zeldovich solution that gives elongated structures in phase space ( middle - left and lower - left inserts of right panels of fig . [ fig : pk_snap_cdm_midres ] . as already argued in previous section , post - collapse pt provides only local corrections to zeldovich flow and therefore can not account for nonlinear dynamics at large scales , hence presents the same defect . this elongation effect obviously leads to a strongly underestimated power spectrum . note however that this depends also on the choice of the cut - off scale @xmath207 : if this scale would be kept , for the theoretical predictions , variable so that it remains close to the typical transition scale to the nonlinear regime , the results would improve greatly both for zeldovich solution and for post - collapse pt , as it is the case for upper insert of left panels of figs . [ fig : pk_snap_cdm_midres ] and [ fig : pk_snap_cdm_midres2 ] . but here , instead of varying the effective smoothing scale @xmath207 for the pt predictions , we now focus on adaptive smoothing , which should be even better because it provides locally the `` optimal '' softening of initial conditions . when adaptive smoothing is employed , agreement between pt predictions and simulation improves strikingly both from the visual and statistical points of view , even for zeldovich solution which performs comparably well to post - collapse pt at @xmath282 , reproducing simulations results in the whole available dynamic range , including the highly nonlinear regime , with a slight underestimation of the power spectrum , though . of course , superiority of post - collapse pt shows up at @xmath283 but zeldovich solution remains quite good . keep in mind , however , that zeldovich solution with adaptive smoothing as presented here uses the calculation of next - crossing time to estimate the local smoothing scale , hence relies as well partly on post - collapse pt . these results confirm the conclusions of [ subsec : merging_clusters ] : adaptive smoothing allows one to summarize composite structures into a single s shape like halo with the right size so that two - point statistics matches the true one even in the highly nonlinear regime . in fact , pt predictions reproduce , at least partly , the plateau that can be observed at small scales in middle and lower insert of left panels of figs . [ fig : pk_snap_cdm_midres ] and [ fig : pk_snap_cdm_midres2 ] . in this regime that we will discuss more in details in next section , the system is most probably following stable clustering . of course , because of the way we implement our adaptive smoothing , continuity of the phase - space sheet is not preserved anymore , as illustrated by right panel of figs  [ fig : pk_snap_cdm_midres ] and [ fig : pk_snap_cdm_midres2 ] , but the overall description of phase - space structures is nevertheless improved tremendously . this is illustrated even more clearly by fig . [ fig : snap_qxv_midres ] , which also shows , as expected , that positions predicted by post - collapse pt or zeldovich solution match the simulation significantly better than predicted velocities . this is a trivial consequence of the fact that positions correspond to velocities integrated over time : by construction , velocities are better approximated by theoretical models at early than at late times , hence their integrated counterpart , which corresponds to some averaged behavior , will compare better to the simulations . note thus on right panel of fig . [ fig : snap_qxv_midres ] the excellent agreement with the simulation for the predicted position from basic post - collapse pt , even when significant mergers take place : one just needs to compare bottom right insert of right panel of fig . [ fig : pk_snap_cdm_midres ] to bottom insert of right panel of fig . [ fig : snap_qxv_midres ] to be convinced by this state of fact . as a consequence , we obtain an excellent match between basic post - collapse pt with adaptive smoothing and simulations measurements for the power spectrum , which is a two - point statistics not depending directly on the quality of representation of the velocity . these arguments also apply , to a lesser extent , obviously , to the zeldovich solution .    when focusing again on the power spectrum ( left panel of figs . [ fig : pk_snap_cdm_midres ] and [ fig : pk_snap_cdm_midres2 ] ) , a detailed comparison between various variant of post - collapse pt models reveals that the basic prescription remains the best among all of them , although _ spl _ performs nearly as well . this confirms the results of [ subsec : single_cluster ] and [ subsec : merging_clusters ] , where we noticed that higher - order prescriptions for post - collapse pt tended to over - predict the size of halos and mergers ( although this is not obvious at first sight when examining phase - space diagrams of right panels of figs . [ fig : pk_snap_cdm_midres ] and [ fig : pk_snap_cdm_midres2 ] ) , leading to a smaller power spectrum than the basic post - collapse pt prediction at small scales .    note that , in agreement with intuition , the results of pt with adaptive smoothing are rather insensitive to the choice of cutoff wavenumber @xmath202 ( as long as it is kept large enough compared to the wavenumber corresponding to transition towards nonlinearity ) , but the quality of the agreement with simulations depends , on the other hand , on the value of the parameter @xmath210 introduced in sec .  [ sec : improvement ] . here , we adopt the intuitive choices @xmath249 for post - collapse pt and @xmath247 for zeldovich solution , but , as illustrated in appendix [ sec : choice_of_f_cross ] , it is possible to optimize the value of @xmath210 . for instance at @xmath283 , we find for the analyses performed in this section that @xmath284 and @xmath285 would provide better agreement with simulation measurements . however we also noticed that the optimal value of @xmath210 depends slightly on redshift as well as on initial conditions , this is why we decided for the analyses presented in this article to stick to our fiducial , slightly suboptimal choice . we now turn to random initial conditions with a power - law power spectrum : @xmath286 this allows us to test post - collapse pt for different kinds of initial conditions and to understand more in details small - scale dynamics . assuming einstein - de sitter universe ( @xmath227 , @xmath237 and @xmath287 ) , we performed simulations for @xmath288 and @xmath246 , with the same resolution as our main cdm runs , i.e. @xmath270 , @xmath271 and @xmath289 mpc@xmath257 . although the initial conditions are scale free , we set an actual size for the simulation box , @xmath290 mpc , which allows us to normalize the power spectrum using @xmath291 ( see footnote [ foot : toto ] ) . with such a normalization , the level of nonlinearity of @xmath292 at the end of the simulations is roughly the same as in the cdm runs studied in previous section . to measure power spectra with sufficient statistics , we performed sets of 50 simulations with different initial random seed for each value of @xmath293 . [ fig : pk_snap_powerlaw ] displays the power spectra ( left panel ) and phase - space diagrams ( right panel ) measured at @xmath283 in the scale - free simulations . the numerical results are compared to basic post - collapse pt ( blue solid ) and to zeldovich solution ( green dotted ) with adaptive smoothing . for reference , linear theory ( black dashed ) and the analytic prediction given by eq . ( [ eq : pk_za_analytic ] ) for the zeldovich solution without smoothing ( black solid ) are also plotted in left panel . in eq . ( [ eq : pk_za_analytic ] ) , we introduced the actual cut - offs in the initial conditions of the simulations , i.e. below @xmath189 and above @xmath289mpc@xmath257 . ] in addition , predictions from stable clustering @xcite are displayed as a thin solid line and will be discussed more in details below . to complete visual inspection of phase space , fig . [ fig : snap_qxv_powerlaw ] plots velocities and positions as functions of lagrangian coordinate @xmath26 in the same region of phase space as in right panel of fig . [ fig : pk_snap_powerlaw ] . note that the simulations considered in right panel of fig . [ fig : pk_snap_powerlaw ] and in fig . [ fig : snap_qxv_powerlaw ] have the same random initial seed as for right panel of fig .  [ fig : pk_snap_cdm_midres ] and fig . [ fig : snap_qxv_midres ] , so one can compare directly phase - space features between scale - free and cdm - like simulations . examination of figs . [ fig : pk_snap_powerlaw ] and [ fig : snap_qxv_powerlaw ] demonstrates again the success of post - collapse pt with adaptive smoothing . while the results are nearly as good for @xmath288 as for the cdm - like case , post - collapse pt does slightly less well for @xmath294 . the same trend can be seen for zeldovich solution with adaptive smoothing . this result is not surprising , for two reasons . firstly , the level of nonlinearity in the @xmath294 simulations is slightly larger than for @xmath288 and the cdm - like initial conditions . indeed , deviation from linear theory or pure zeldovich solution happens at smaller values of @xmath186 for @xmath294 than for other cases . secondly and more importantly , adaptive smoothing using a sharp low - pass filter is less efficient , from the dynamical point of view , when considering a power spectrum with a high slope . the results would look even worse in the case @xmath295 that we simulated but do not show here . to explain this , one can just take the example of the single @xmath296 mode system , e.g. initial conditions with a single sine wave as shown in fig . [ fig : single - cluster_xv_basic ] . in this case , with a sharp-@xmath186 filter , the fluctuation at @xmath297 is either kept , hence post - collapse pt becomes invalid after the system reaches next - crossing time , or it is completely suppressed , which means that the procedure with adaptive smoothing is at the end equivalent to no smoothing at all . increasing the power spectrum slope @xmath293 basically makes us approaching the single mode case with @xmath298 , since we always introduce a cut - off at some value @xmath202 to have a sufficiently smooth representation of the phase - space sheet at very small scales . of course , this situation is unrealistic in practice , because no such problem occurs when considering cdm - like power spectra , but one has to stay aware of the fact that adaptive smoothing does not always work . a possible way to improve it could consist in changing the filtering window so that all the wavenumbers are affected by the softening procedure , even if it is performed predominantly at a given scale : to this respect , adaptive smoothing with a gaussian window might represent an interesting alternative to our sharp-@xmath186 filter . another interesting result of our measurements is the excellent agreement between the simulations and the stable clustering prediction from @xcite and @xcite in the nonlinear regime at very high-@xmath186 , irrespective of the initial power - law slope @xcite . ] , who predict that the slope of the power spectrum should be equal to @xmath299 in fact , even in the cdm - like case , which roughly corresponds , in one dimension , to an asymptotic value of @xmath293 equal to @xmath300 at large @xmath186 , stable clustering prediction works , since it predicts a plateau , @xmath301 , which is indeed what we observe on left panel of fig . [ fig : pk_snap_cdm_midres ] . note that what @xcite and @xcite call `` stable clustering '' is very specific . in the 1d case considered here , it is equivalent to assuming that relaxed objects become of constant size in the coordinate @xmath302 instead the physical coordinate @xmath303 as normally supposed in the three - dimensional case ( see , e.g. * ? ? ? * ; * ? ? ? it is important to notice here that post - collapse pt seems to be able to reproduce correctly the measured power spectrum in the beginning of the stable clustering regime . this suggests that stable clustering regime is established only in a few dynamical times , since post - collapse pt is only able to follow the dynamics of a single halo until next - crossing time .    on the other hand , if one considers a `` halo model '' approach ( see , e.g. * ? ? ? * ; * ? ? ? * ; * ? ? ? * ; * ? ? ? * and references therein ) and assumes that at large @xmath186 the power spectrum is dominated by the `` one - halo term '' , interesting properties are expected , at least when initial density fluctuations are sufficiently smooth at small scale . indeed , one dimensional simulations of single cold halos suggest that their central density profile nearly behaves like a power - law @xmath304 with @xmath305 . the actual value of @xmath306 is controversial . while @xcite postulated it to be equal to @xmath247 , detailed @xmath13-body simulations suggest instead @xmath307 @xcite . on the other hand , pure vlasov simulations results using the waterbag method suggest @xmath308 @xcite . other values to keep in mind are @xmath309 from the singularity occurring exactly at crossing times and @xmath310 from the caustics ( see , e.g. * ? ? ? * ) . these singularities imply , as discussed below , some expected behavior for @xmath311 at large @xmath186 @xcite . a slope @xmath312 for the density profile corresponds , after some trivial algebraic calculations , to a power - law slope @xmath313 for the power spectrum . for instance , one expects the zeldovich approximation to be dominated by the caustics at very small scales , i.e. @xmath314 constant at large @xmath186 , which is indeed the case on left panels of figs . [ fig : pk_snap_cdm_midres ] and [ fig : pk_snap_powerlaw ] in both analytic ( black ) and numerical ( green ) cases , irrespective of adaptive smoothing . turning to post - collapse pt with adaptive smoothing , since , by using @xmath249 , we push the ( smoothed ) system exactly to next - crossing time , we are dominated by the singularity with @xmath309 . hence , in this case one expects @xmath315 close to @xmath316 when @xmath186 becomes very large , as can be observed on all the figures . with a smaller value of @xmath210 , the system is dominated , similarly as for zeldovich solution , by the caustics at very small scales and then @xmath317 reaches a plateau at large @xmath186 as can be seen in fig . [ fig : pk_xv_fcross ] of appendix [ sec : choice_of_f_cross ] . without adaptive smoothing , one also intuitively expect post - collapse pt power spectrum to be dominated by caustics at very small scales , however the objects have a more complex structure than zeldovich solution as shown by lower right panel of fig . 2 , hence the plateau at large @xmath186 is not clearly visible on left panel of figs . [ fig : pk_snap_cdm_midres ] , [ fig : pk_snap_cdm_midres2 ] and [ fig : pk_snap_powerlaw ] probably because not reached as quickly as for the pure zeldovich solution when increasing @xmath186 .    turning to the actual nonlinear evolution of single objects , assuming as found in the literature @xmath308 and @xmath247 for the density profile of a halo would respectively give @xmath318 and @xmath319 for the `` one halo '' power spectrum . this of course supposes that statistical averaging over many halos of different masses does not affect the slope of the power spectrum predicted by the single halo . again , we have to be aware of the fact that at very large @xmath186 , one must be dominated by the caustics , since we introduce a cut - off at large @xmath186 . such a cut - off enforces some level of smoothness of the curve representing the system in phase space and thus actually imposes the existence of well - defined caustic structures in the nonlinear regime , hence a plateau for @xmath317 at very large @xmath186 , even for the simulations . however this plateau is not visible in the simulations , because it lies outside the range of @xmath186 shown in the figures . on the other hand , there must be , according to our `` halo model '' prescription , a regime where @xmath317 might decrease ( @xmath308 ) or present another plateau ( @xmath320 ) or at least some kind of inflection at some intermediate values of @xmath186 . this is what can indeed be observed e.g. on lower insert of left panel of fig . [ fig : pk_snap_cdm_midres ] for @xmath321 mpc@xmath257 for the cdm - like case and on lower insert of left panel of fig . [ fig : pk_snap_powerlaw ] for @xmath322 mpc@xmath257 for the @xmath294 simulation . there is no such signature obviously visible in the @xmath288 simulation , although the measurements are too noisy to conclude in this case . interestingly , when the inflection is visible , its extension in @xmath186 space is correctly predicted by basic post - collapse pt with adaptive smoothing . this just shows that basic post - collapse pt is able to predict approximately the correct size for the halos , which we already knew . note that here , we did not try to position our single halo term discussion in the context of stable clustering predictions . these latter can be considered as the outcome of some statistical averaging after weighting the one - halo term by the mass function . they are also expected to be valid only in a finite range of values of @xmath186 if one assumes a cut - off at @xmath202 : what we mean here is that our `` one - halo term '' discussion is in fact not incompatible with the concept of stable clustering although the link between the halo model and the stable clustering predictions remain to be performed in details in the present case . clearly , however , the small inflection discussed above marks a small deviation from stable clustering . in this paper , we have developed a lagrangian perturbation theory method for solving , beyond crossing time , vlasov - poisson equations in the expanding universe for initially cold systems in the infinite parallel planes geometric set up . the proposed approach captures post - collapse dynamics by computing at leading - order a counter term to the zeldovich solution just after first crossing time . it extends earlier work of @xcite to the cosmological case and , thanks to adaptive smoothing , to random initial conditions instead of a single halo . by performing a local taylor expansion of the velocity and position of each mass element as functions of lagrangian coordinate around each initial density peak , we are able to compute the force field in the multi - valued region . a correction to the single - stream flow motion , described by the zeldovich solution , is systematically obtained . the results are given as polynomial functions of lagrangian position and the coefficient of each power is expressed in terms of quantities associated to density peaks . a simple analysis of the formation of a single cluster shows that post - collapse pt reproduces well the phase - space structure up to next - crossing time and can describe the projected density profile quite well even at later time .    to deal with the bulk properties of phase - space structures , we propose an adaptive filtering scheme the idea is very close to the peak - patch treatment by @xcite . in particular , the local smoothing scale of initial conditions is determined by the calculation of a dynamical time , but in post - collapse pt this dynamical time is given by next - crossing time instead of collapse time . this allows us to produce a realistic coarse - grained version of the phase - space structure , even in the presence of mergers . after demonstrating that adaptive smoothing indeed works well in a simple setup with coalescing clusters , we applied it to more general , random initial conditions that include cdm - like as well as scale - free power spectra . the results of our analyses show that post - collapse pt predictions with adaptive smoothing reproduce the power spectrum of our numerical simulations remarkably well even at scales sufficiently nonlinear to probe the beginning of the stable clustering regime .    on important thing to keep in mind for further work is that the critical quantity in post - collapse pt is the time of next crossing , which can in fact be used to also create a spectacularly good pure zeldovich prescription with adaptive smoothing . while the calculations of post - collapse pt can be cumbersome , the expression ( [ eq : tau_cross ] ) of this time is relatively simple and this can be very interesting when trying to generalize results of post - collapse pt to the three - dimensional case . when analyzing equation ( [ eq : tau_cross ] ) , one indeed notes that @xmath198 depends both on collapse time @xmath53 and on the height @xmath323 of the initial density peak . since @xmath53 depends solely on @xmath323 , we can see that basically the dynamical time we consider for defining the local smoothing scale is a single function of @xmath323 , larger than @xmath53 . in fact , things are very simple : when performing adaptive smoothing one enforces @xmath324 to be equal to present time , or in other words , adaptive smoothing of initial conditions is performed so that the local peak height given by linear theory at present time , @xmath323 , is a number larger than unity ( instead of exactly unity for the collapse time ) . for instance , in the einstein - de sitter case , @xmath325 and @xmath237 , we have @xmath326 with @xmath210 chosen in this article to be equal to unity for post - collapse pt and one - half for the zeldovich solution . note that the optimal value of @xmath210 is smaller than these intuitive choices ( see appendix [ sec : choice_of_f_cross ] ) . the present post - collapse pt treatment together with the adaptive smoothing scheme seems a very promising tool . a successful extension to the 3d case should lead to a major breakthrough in the precision pt calculations of large - scale structure statistics beyond the single - stream approximation . our 1d study is however only a first step . toward a practical application to the 3d case , there are several issues to be addressed . one is the analytic calculation of statistical quantities such as the power spectrum . indeed , in the present work , we measured directly the power spectrum from a particle distribution following the dynamics prescribed by the various pt prescriptions under scrutiny . in fact , analytic calculations of the power spectrum are tractable , at least in the presence of fixed smoothing . adaptive smoothing itself can in principle be implemented in the analytic framework by inspiring for instance from the theory of excursion sets discussed in e.g. @xcite , although we can foresee that actual calculations will probably be very cumbersome . another issue , when generalizing to the 3d case , is that the zeldovich solution no longer provides a sufficiently accurate description of the dynamics before collapse . at least third - order lagrangian pt is required to estimate correctly collapse times as well as the power spectrum beyond tree - level approximation , which clearly adds some significant level of complexity . furthermore , the variety and complexity of topological configurations of singularities appearing in the 3d case might in fact represent a nearly impassable obstacle in three dimensions ( see , e.g. * ? ? ? we will tackle these issues in future work .    finally , while we have focused in this article on the impact of post - collapse dynamics in 1d on the phase - space structure and on the power spectrum of the matter distribution , it would be interesting to see quantitatively how the small - scale modes responsible for the multi - stream flows couple to large - scale fluctuations described by the single - stream flows . this can be presumably clarified by measuring the response of the evolved power spectrum to a small initial perturbation of high-@xmath186 modes , through the comparison between a pair of particle distributions modeled with post - collapse pt that slightly differ in their initial conditions @xcite . in the 3d case , a strong damping is found for the mode transfer from small to large scales , but the physical origin of it is still unclear . a detailed analysis in the 1d cosmological setup would provide important clues on the nature of mode transfers . at thanks patrick valageas for references and useful comments on the adhesion model . sc greatly acknowledges hospitality of yukawa institute for theoretical physics ( yitp ) , where this work was initiated . we also acknowledge the support of yitp in organizing the workshop `` vlasov - poisson : towards numerical methods without particles '' in kyoto , funded by grant yitp - t-15 - 02 , anr grant anr-13-monu-0003 and by institut lagrange de paris ( anr-10-labx-63 and anr-11-idex-0004 - 02 ) . this work was also supported in part by mext / jsps kakenhi grant numbers jp15h05889 and jp16h03977 ( at ) as well as anr grant anr-13-monu-0003 ( sc ) . in this appendix , we derive the expression in eq . ( [ eq : phi_poisson ] ) . to do this , we first consider the 1d green function which satisfies the following equation : @xmath327 with this green function , the solution of the poisson equation is formally expressed as @xmath328 to solve eq . ( [ eq : green_func ] ) under periodic boundaries condition , we move to harmonic space and express the green function as : @xmath329 recalling that the fourier transform of dirac s delta function reads @xmath330 $ ] , . ] substitution of eq . ( [ eq : g_fourier ] ) into eq . ( [ eq : green_func ] ) gives @xmath331 thus , we obtain @xmath332 note that there is a subtlety in the above equation . in deriving it , we ignore the contribution from @xmath288 mode in @xmath333 , while the @xmath288 mode does not vanish in dirac s delta function . the result remains however correct as long as the density field does not contain a @xmath288 mode , or equivalently , @xmath334 . keeping this point in mind , we can further simplify the expression . we use the following formula ( see , e.g. chap.1.44 of * ? ? ? * ; * ? ? ? * ) : @xmath335 the final expression for the green function is then given by @xmath336,\qquad |x - x'|\leq l.\end{aligned}\ ] ] substituting this into eq . ( [ eq : integral ] ) , we obtain the expression in eq . ( [ eq : phi_poisson ] ) . in this appendix , we derive the expression for the force in the multi - stream region , given in eqs . ( [ eq : force_inner ] ) , ( [ eq : force_inner_q ] ) and ( [ eq : force_inner_const ] ) .    to start with , we divide the integral in eq . ( [ eq : force_poisson ] ) into several pieces according to eq . ( [ eq : part_by_part_integrals ] ) : @xmath337 . \label{eq : force_pieces}\end{aligned}\ ] ] we shall below evaluate these integrals one by one .    consider the first integral in right - hand side of eq . ( [ eq : force_pieces ] ) . we evaluate it at position @xmath338 in fig . [ fig : shellcrossing ] . using eq . ( [ eq : jacobian_delta ] ) , the integral is rewritten @xmath339.\end{aligned}\ ] ] exploiting the geometric properties of the configuration in fig . [ fig : shellcrossing ] , substitution of the explicit form of @xmath340 ( eq .  [ eq : ballistic_shellcrossing ] ) into the above leads to @xmath341 here we used the fact that @xmath71 . recalling that @xmath342 and @xmath343 are related with each other through eqs . ( [ eq : relation_three - values_1 ] ) and ( [ eq : relation_three - values_2 ] ) , eq . ( [ eq : integral_inner1_mid ] ) can be recast as a single function of either @xmath342 or @xmath343 , which respectively provides the general expression valid for @xmath180 and @xmath344 . omitting the subscript @xmath345 or @xmath346 in either cases , we finally obtain @xmath347    next , we consider the second integral in eq . ( [ eq : force_pieces ] ) . rewriting the integrand in terms of the variable @xmath81 with the help of eqs . ( [ eq : jacobian_delta ] ) and ( [ eq : ballistic_shellcrossing ] ) , it is straightforward to obtain @xmath348    the remaining integrals to be computed , i.e. , the third and fourth terms in eq . ( [ eq : force_pieces ] ) , are performed over intervals outside the multi - valued region , which allows us to use zeldovich solution to evaluate the integrands . for the third integral in eq . ( [ eq : force_pieces ] ) , we obtain @xmath349 although this expression is exact , our perturbative treatment implicitly assumes that the multi - stream region is small and the quantities inside this region are written in polynomial forms of @xmath81 . we may expand the last two terms of eq . ( [ eq : integral_outer1_mid ] ) as @xmath350 where we used the fact that @xmath351 ( see eq .  [ eq : shellcrossing_peak ] ) . additionally , periodic boundaries imply @xmath352 . therefore , eq .  ( [ eq : integral_outer1_mid ] ) can be recast as @xmath353_{q_0}^{q(x = l ) } d_+(\tau ) . \label{eq : integral_outer1}\end{aligned}\ ] ]    in similar manner , we evaluate the fourth integral in eq . ( [ eq : force_pieces ] ) as follows : @xmath354 the first of these two integrals is performed exactly to give @xmath355_{q(x=0)}^{q(x = l ) } -\frac{1}{2 } \bigl[\psi(q')^2\bigr]_{q(x=0)}^{q(x = l)}\ , \nonumber\\ & \quad\quad + \int_{q(x=0)}^{q(x = l ) } { \rm d}q'\,\psi(q')\biggr\}\,d_+(\tau ) \nonumber\\ & \quad = -\psi(q(x=0))\,d_+(\tau ) ,   \label{eq : outer2_part1}\end{aligned}\ ] ] where we exploited properties of periodic boundaries , @xmath356 and @xmath357 , as well as momentum conservation , which implies @xmath358 . turning to the second integral in eq . ( [ eq : outer2_step1 ] ) , we use again the fact that the domain of integration is assumed to be narrow to taylor expand the integrand using the formula @xmath359 to obtain the following expression , valid at third - order in @xmath74 : @xmath360 here , eqs . ( [ eq : def_a])([eq : def_c ] ) are used to rewrite the expression in the last line . summing up the results in eqs . ( [ eq : outer2_part1 ] ) and ( [ eq : outer2_part2 ] ) , we obtain , at the end , @xmath361    plugging eqs . ( [ eq : integral_inner1 ] ) , ( [ eq : integral_inner2 ] ) , ( [ eq : integral_outer1 ] ) and ( [ eq : integral_outer2 ] ) into eq . ( [ eq : force_pieces ] ) , we finally obtain the expression of the force in the multi - valued region ( eqs . [ eq : force_inner ] , [ eq : force_inner_q ] and [ eq : force_inner_const ] ) . in this appendix , we detail the calculations of basic post - collapse pt predictions summarized in eqs . ( [ eq : p_for_tau_c > tau > hattau_c ] ) , ( [ eq : x_for_tau_c > tau > hattau_c ] ) , ( [ eq : p_for_tau > tau_c ] ) and ( [ eq : x_for_tau > tau_c ] ) , together with the coefficients in table  [ tab : xv_pcpt1 ] and [ tab : xv_pcpt2 ] . our starting point is eq . ( [ eq : v_ii ] ) , which involves two integrals . to evaluate the first integral , the approximate formula given by eq . ( [ eq : integ_formula_j ] ) is applied . recalling that we are considering the outer part of the multi - valued region , @xmath180 , we have @xmath363 , \label{eq : integ_j1}\end{aligned}\ ] ] with @xmath364 . here , we used the leading - order expression for critical time @xmath82 , @xmath365 , and ignore in eq . ( [ eq : integ_formula_j ] ) the higher - order terms involving the quantities @xmath176 or @xmath177 .    on the other hand , the second integral in eq . ( [ eq : v_ii ] ) is analytically performed . using the evolution equation in eq . ( [ eq : evolv_lin ] ) , we obtain @xmath366_{\widehat{\tau}_{\rm c}(q)}^\tau . \label{eq : integ_f}\end{aligned}\ ] ]    summing up the above two contributions , we obtain the expression for @xmath362 summarized in eq . ( [ eq : p_for_tau_c > tau > hattau_c ] ) .      to derive the expression for @xmath126 , we just need to integrate over time eq . ( [ eq : p_for_tau_c > tau > hattau_c ] ) , as written in eq . ( [ eq : x_ii ] ) : @xmath367 \nonumber \\ & ~\quad + \int_{\widehat{\tau}_{\rm c}(q)}^\tau { \rm d}\tau'\,\,\widehat{\epsilon}_1 \,(\tau',\widehat{\tau}_{\rm c } ) . \label{eq : x_ii_prime}\end{aligned}\ ] ] making use of the expressions summarized in table [ tab : xv_pcpt1 ] , we evaluate the first integral . a straightforward calculation leads to @xmath368 \nonumber\\ & \simeq   \frac{t^2}{2}\,q+ \bigl\{-\frac{\kappa}{8}\,t+\frac{\delta_{\rm l}''(q_0)}{6}d_+(\tau_0 ) \ , \bigr\}q^3 \nonumber\\ & \quad -\mbox{sgn}(q)\,\frac{\kappa^2}{80\sqrt{3}}\ , \left\{{\widehat{q}_{\rm c}}^2(\tau)-q^2\right\}^{5/2 } \nonumber\\ & \quad + \biggl[\frac{1}{2}\left(\frac{\kappa}{8}\right)^2-\left(\frac{\kappa}{8}\right)\,\frac{\delta_{\rm l}''(q_0)}{6}d_+(\tau_0)\,t\,\biggr]q^5 \nonumber\\ & \quad+ \frac{1}{2}\left(\frac{\kappa}{8}\right)^2 \frac{\delta_{\rm l}''(q_0)}{6}d_+(\tau_0)q^7 . \label{eq : tilde_x_i}\end{aligned}\ ] ] note again that we used the leading - order expression for @xmath82 ( see eq .  [ eq : hattau_c ] ) . the second integral in eq . ( [ eq : x_ii_prime ] ) is analytically performed to give @xmath369^{\tau'}_{{\widehat{\tau}_{\rm c}}(q ) } \nonumber\\ & \quad = \psi(q_0)\,\bigg\ { d_+(\tau)-d_+({\widehat{\tau}_{\rm c}}(q ) )   \nonumber\\ & \quad \quad \quad \quad -\left.\left.\frac{{\rm d}d_+}{{\rm d}\tau}\right|_{{\widehat{\tau}_{\rm c}}(q)}(\tau-{\widehat{\tau}_{\rm c}}(q))\right\}. \label{eq : integ_eps_1}\end{aligned}\ ] ] collecting eqs .  ( [ eq : tilde_x_i ] ) and ( [ eq : integ_eps_1 ] ) , we obtain eq . ( [ eq : x_for_tau_c > tau > hattau_c ] ) .      for @xmath153 , the expression to be evaluated is given by eq . ( [ eq : delta_p_2 ] ) , which involves three integrals . the first integral is obtained from eq . ( [ eq : integ_j1 ] ) by setting @xmath370 ( eq .  [ eq : tau_c ] ) . recalling the fact that @xmath371 ( eq .  [ eq : qchat_qc ] ) , we have @xmath372.\end{aligned}\ ] ] for the second integral in eq . ( [ eq : delta_p_2 ] ) , we use the approximate formula given by eq . ( [ eq : integ_formula_j ] ) in the regime @xmath373 . ignoring higher - order terms involving the quantities @xmath176 or @xmath177 , we obtain @xmath374.\end{aligned}\ ] ] finally , the last integral in eq . ( [ eq : delta_p_2 ] ) is evaluated using the formula given in eq . ( [ eq : integ_f ] ) . collecting the three contributions above , we obtain the expression for @xmath153 ( eq . [ eq : p_for_tau > tau_c ] ) . the starting expression for @xmath119 is given in eq . ( [ eq : x_iii ] ) , which involves two integrals . to evaluate the first integral in eq . ( [ eq : x_iii ] ) , we use eq . ( [ eq : tilde_x_i ] ) . setting @xmath375 , an expression at seventh - order in @xmath81 is obtained : @xmath376    to compute the second integral in eq . ( [ eq : x_iii ] ) , we substitute eq . ( [ eq : p_for_tau > tau_c ] ) into the integrand , but without the term @xmath377 . using the explicit expression for the coefficients in table [ tab : xv_pcpt2 ] , a straightforward calculation gives @xmath378 \nonumber\\ & ~=-\frac{3}{2}h_0 ^ 2\,\omega_{\rm m,0}\,a(\tau_0)\ , \nonumber\\ & \quad \times\biggl(\,-t^2\,q+ \left\{\frac{5\kappa}{8}\,t+\frac{\delta_{\rm l}''(q_0)}{12}\,d_+(\tau_0 ) t^2\right\}\,q^3 \nonumber\\ & \qquad\quad- \bigl\{\left(\frac{\kappa}{4}\right)^2\ , + \left(\frac{\kappa}{8}\right)\frac{\delta_{\rm d}''(q_0)}{6}\,d_+(\tau_0)\,t\ , \bigr\}\,q^5 \nonumber\\ & \qquad\quad -4\left(\frac{\kappa}{8}\right)^2\frac{\delta_{\rm l}''(q_0)}{6}\,d_+(\tau_0)\,q^7 \biggr ) . \label{eq : integ_x_ii_part2}\end{aligned}\ ] ]    combining the above two results , we obtain eq . ( [ eq : x_for_tau > tau_c ] ) , which gives the final expression for @xmath119 along with the coefficients in table [ tab : xv_pcpt2 ] . in this appendix , we detail calculations of the various higher - order corrections to post - collapse pt we propose in   [ subsec : higher - order ] .    one possible way to improve on basic post - collapse pt is to taylor expand at next - to - leading order the expression for @xmath74 ( as well as @xmath75 ) that determines the boundary of the post - collapse region ( eqs . [ eq : qchat_qc ] and [ eq : qchat_approx ] ) : @xmath379 with @xmath380 inverting the relation @xmath381 [ @xmath382 , we obtain the corresponding critical times at next - to - leading order , @xmath383 and @xmath384 using these expressions , we repeat below the same calculations as in sec .  [ subsec : corrections_pcpt ] . higher - order expressions for the correction to zeldovich flow are derived on top of basic post - collapse pt results and are denoted by @xmath385 and @xmath386 . in sec . [ subsec : higher - order ] , a semi - analytic treatment designed by _ spl _ is also discussed . in this case , we adopt the basic post - collapse pt results in the inner part of the multi - stream region , but the extension of this latter is calculated using the higher - order expression ( [ eq : higher_order_qhat ] ) for @xmath75 . then , a straightforward third - order spline interpolation , that we do not detail here , is used to connect the inner part to the zeldovich solution , assuming again that the lagrangian boundary @xmath74 of the multi - stream region is given by the higher - order expansion ( [ eq : higher_order_qhat ] ) . another improvement may come from the calculation of the integrals in eqs . ( [ eq : delta_p_formal ] ) and ( [ eq : delta_x_formal ] ) . in sec . [ subsec : corrections_pcpt ] , a part of the integrands is taylor - expanded in time and the integration is performed for the leading - order terms . the resultant expressions include time - dependent terms up to @xmath387 for the velocity and @xmath388 for the position . here , using the approximate integral formulae for the taylor - expanded integrands in appendix [ sec : integral ] , we keep terms up to @xmath388 for the velocity , @xmath389 for the position , and derive the corresponding higher - order corrections to the zeldovich flow , @xmath390 and @xmath391 . these corrections assume , of course , higher - order expansions for critical times @xmath82 and @xmath170 . below , we present the expressions for @xmath390 and @xmath391 . the expressions for @xmath385 and @xmath386 are obtained by simply setting @xmath176 and @xmath177 to zero . the expression for the higher - order correction to the velocity , @xmath393 , is @xmath394 , \label{eq : p_higher_tau_c > tau > hattau_c_higher}\end{aligned}\ ] ] with @xmath395 . \label{eq : def_delta_zeta_1_higher}\end{aligned}\ ] ] here , we have introduced the following function : @xmath396 where @xmath177 is the value of the hubble parameter at time @xmath53 . for the higher - order correction to the position , @xmath397 , we have @xmath398 , \label{eq : x_higher_tau_c > tau > hattau_c_higher}\end{aligned}\ ] ] with @xmath399 note the term proportional to @xmath400 , which is necessary to enforce continuity of the solution . the expression for the higher - order corrections to the velocity , @xmath402 , becomes @xmath403 , \label{eq : p_higher_tau > tau_c}\end{aligned}\ ] ] with @xmath404 , \\ & \widetilde{\beta}_3(\tau)= \frac{g(\tau_0)}{2}t^2\,\frac{\delta_{\rm l}''(q_0)}{6}\,d_+(\tau_0 ) , \\ &   \widetilde{\delta}_3(\tau_0)=\frac{47}{64}\kappa^2\,\eta   -\frac{\kappa^2}{128 }   \bigl[\,g(\tau_0)-48\{a(\tau_0)\}^2\,h(\tau_0)\bigr ] , \\ & \widetilde{\zeta}_3(\tau_0)=-\left(\frac{\kappa}{8}\right)^2\ , \left\{\,\frac{g(\tau_0)}{2}+\eta\,\right\ } \frac{\delta_{\rm l}''(q_0)}{6}\,d_+(\tau_0 ) \nonumber\\ & \qquad\qquad-\left(\frac{\kappa}{8}\right)^3\,\eta\ \bigl[\,g(\tau_0)-192\{a(\tau_0)\}^2\,h(\tau_0)\bigr].\end{aligned}\ ] ] finally , the higher - order corrections to the position , @xmath405 , are given by @xmath406 , \label{eq : x_higher - order_tau > tau_c}\end{aligned}\ ] ] with @xmath407 , \label{eq : def_delta_alpha_4 } \\ & \widetilde{\beta}_4(\tau)= \frac{\delta_{\rm l}''(q_0)}{6}\,d_+(\tau_0 ) \frac{g(\tau_0)}{6}\,t^3 ,   \label{eq : def_delta_beta_4 } \\ & \widetilde{\delta}_4(\tau)= -\left(\frac{\kappa}{4}\right)^2\,\left(\frac{g(\tau_0)}{8}\,t-6\{a(\tau_0)\}^2h(\tau_0)\,t\,\right ) \nonumber\\ & \qquad\quad+\,47\left(\frac{\kappa}{8}\right)^2\,\eta\,t , \label{eq : def_delta_delta_4 } \\ & \widetilde{\zeta}_4(\tau)= -\left(\frac{\kappa}{8}\right)^2 t \left(\frac{g(\tau_0)}{2}+\eta\right)\frac{\delta_{\rm l}''(q_0)}{6}\,d_+(\tau_0 ) \nonumber\\ & \qquad\quad + \left(\frac{\kappa}{8}\right)^3\biggl\{\,\frac{g(\tau_0)}{3}-64\,\{a(\tau_0)\}^2h(\tau_0 ) \nonumber\\ & \qquad\quad -\eta\,\bigl(95+g(\tau_0)\,t-192\,\{a(\tau_0)\}^2h(\tau_0)\,t\bigr)\,\biggr\ } , \label{eq : def_delta_zeta_4 } \\ & \widetilde{\mu}_4(\tau)=\kappa^3\,\delta_{\rm l}''(q_0)\,d_+(\tau_0 ) \left(\frac{g(\tau_0)}{9216}+\frac{1}{3072}\eta-\frac{g(\tau_0)\,\eta\,t}{3072}\right ) \nonumber\\ & \qquad\quad + \kappa^4\,\eta\,\biggl\{\frac{g(\tau_0)}{4096}-\frac{3}{16}\,\{a(\tau_0)\}^2h(\tau_0)-\frac{767\eta}{8092 } \nonumber\\ & \qquad\quad -\frac{g(\tau_0)\,\eta\,t}{8192}+\frac{3}{32}\eta\,\{a(\tau_0)\}^2h(\tau_0)\,t\biggr\ } , \nonumber\\ & \widetilde{\nu}_4(\tau)=\frac{1}{49152}\,\kappa^4\eta^2\delta_{\rm l}''(q_0)\,d_+(\tau_0 ) . \label{eq : def_delta_nu_3}\end{aligned}\ ] ] in this appendix , we provide useful formulae for primitive integrals over time , which are used to derive the post - collapse pt results in appendices [ sec : pcpt_basic ] and [ sec : pcpt_higher - order ] . our aim here is to evaluate the time - integral of the force in the multi - valued region , given in eq . ( [ eq : force_inner ] ) . the time - integral of the second term , involving function @xmath408 , is performed exactly ( see eq . [ eq : integ_f ] ) . on the other hand , the integral of the first term needs an approximate treatment which will be valid as long as the multi - valued region is sufficiently small , or equivalently , shortly after collapse time @xmath53 . below , we evaluate the integral by taylor - expanding the integrand . we first give the results term by term : @xmath409 where @xmath176 is given by eq . ( [ eq : def_t_g ] ) . for the term with fractional power of time , we obtain @xmath410 where we used the leading - order approximation for @xmath74 , i.e. , @xmath411 .    summing up the above formulas , we obtain the primitive of the time integrals of the force involving function @xmath412 . the result including polynomials of time up to @xmath413 is @xmath414_{\tau_{\rm i}}^{\tau_{\rm f}}\,q   \\ \\   \qquad+\frac{1}{6}\delta_{\rm l}''(q_0)\,d_+(\tau_0)\ , \\ \\ \qquad\times\ , \bigl[\ ,   ( \tau-\tau_0)+ \frac{1}{2}g(\tau_0)\,(\tau-\tau_0)^2 \bigr]_{\tau_{\rm i}}^{\tau_{\rm f}}\,q^3 \\ \qquad - \mbox{sgn}(q)\,{\displaystyle \frac{\kappa}{4\sqrt{3}}}\,\bigl[\ , \left({\widehat{q}_{\rm c}}(\tau)^2-q^2\right)^{3/2 } \bigr]_{\tau_{\rm i}}^{\tau_{\rm f } } \,\biggr\ } \\ \qquad\qquad\qquad\qquad\qquad\qquad;\quad q_{\rm c}<|q|<{\widehat{q}_{\rm c } } , \\ \\ a(\tau_0)\,\biggl\{\ , \bigl[\,-2 ( \tau-\tau_0)+ \frac{1}{2}\bigl\{-3\,h(\tau_0)\,a(\tau_0)^2 \\ \qquad + g(\tau_0)\bigr\}(\tau-\tau_0)^2\ , \bigr]_{\tau_{\rm i}}^{\tau_{\rm f}}\,q    +   \frac{1}{6}\delta_{\rm l}''(q_0)\,d_+(\tau_0 ) \\ \qquad \times \bigl[\,(\tau-\tau_0)+ \frac{1}{2}g(\tau_0)\,(\tau-\tau_0)^2\ , \bigr]_{\tau_{\rm i}}^{\tau_{\rm f}}\,q^3 \biggr\ } \\ \qquad\qquad\qquad\qquad\qquad\qquad ; \quad |q|<q_{\rm c}. \end{array } \right .   \label{eq : integ_formula_j}\end{aligned}\ ] ] in this appendix , we study the effect of varying the parameter @xmath210 on the performances of pt predictions with adaptive smoothing ( see sec .  [ sec : improvement ] for definition ) . in the main text , we adopt the intuitive setup @xmath249 for post - collapse pt and @xmath247 for zeldovich solution . however , pt predictions are prone , after shell - crossing , to get worse over time , so we can naively expect that choosing a smaller value of @xmath210 will provide a better agreement with simulations . [ fig : pk_xv_fcross ] shows the results obtained in the case of the cdm - like power spectrum studied in   [ subsec : cdm_initial ] when varying the value of @xmath210 . here , the calculations are performed using @xmath415 , @xmath416 and @xmath417mpc@xmath257 for a boxsize @xmath290 mpc , which is a relatively low resolution set up , but this will not have any consequence on the discussion that follows . left panel shows the power spectrum at @xmath283 and the resultant predictions for different values of @xmath210 are plotted as solid lines with different colors . as anticipated , decreasing @xmath210 results in an enhancement of small - scale power and a better agreement between post - collapse pt and the simulation is obtained . but a too small value of @xmath210 overshoots the simulation and an optimal choice of @xmath210 roughly corresponds to @xmath418 . the same trend can be also seen for the zeldovich solution , but the improvement of the power spectrum amplitude is rather mild : while @xmath419 seems to provide the best overall behavior for the zeldovich solution , the result is clearly still not as good as what post - collapse pt can provide , in particular in the regime @xmath420 mpc@xmath257 . right panel of fig . [ fig : pk_xv_fcross ] shows the phase - space portraits and compares the simulation results ( red ) to predictions obtained with our `` standard '' value of @xmath249 for post - collapse pt ( upper insert ) and @xmath421 for zeldovich solution ( lower insert ) and to predictions obtained with respective nearly optimal values of @xmath284 and @xmath285 . it seems difficult to judge by eye whether the prediction with a smaller value of @xmath210 really improves the description at small scales , although more structures seems to be captured for post - collapse pt with @xmath249 , which contradict the results obtained for the power spectrum . however the simple examination of this figure corresponding to a single realization of the random initial conditions is not conclusive , obviously : it does not preclude the fact that a smaller value of @xmath210 can give a better result for the power spectrum obtained from averaging over many realizations . finally , to conclude this section , although we do not show the results here , we noticed as well that the best choice of @xmath210 could vary slightly according to redshift or initial condition : for example , for the cdm cosmology considered here , a better choice of @xmath210 at redshift @xmath282 is @xmath422 and @xmath423 respectively for post - collapse pt and zeldovich solution . this finally explains why we decided to keep , for simplicity , the generic values @xmath249 and @xmath421 respectively for post - collapse pt and zeldovich solution , even if they are sub - optimal .
we develop a new perturbation theory ( pt ) treatment that can describe gravitational dynamics of large - scale structure after shell - crossing in the one - dimensional cosmological case . starting with cold initial conditions , 
 the motion of matter distribution follows at early stages the single - stream regime , which can , in one dimension , be described exactly by the first - order lagrangian perturbation , i.e. the zeldovich solution . 
 however , the single - stream flow no longer holds after shell - crossing and a proper account of the multi - stream flow is essential for post - collapse dynamics . in this paper , extending previous work by ( * ? ? ? 
 * mnras 446 , 2902 ) , we present a perturbative description for the multi - stream flow after shell - crossing in a cosmological setup . in addition , we introduce an adaptive smoothing scheme to deal with the bulk properties of phase - space structures . 
 the filtering scales in this scheme are linked to the next - crossing time in the post - collapse region , estimated from our pt calculations . 
 our pt treatment combined with adaptive smoothing is illustrated in several cases . 
 predictions are compared to simulations and we find that post - collapse pt with adaptive smoothing reproduces the power spectrum and phase - space structures remarkably well even at small scales , where zeldovich solution substantially deviates from simulations . 
 [ firstpage ]    large - scale structure  cosmology  vlasov - poisson equation