NotebookEvaluate[NotebookDirectory[]<>"/split_step.nb"];

HarmonicTrap1D[];
sizes={256};
timeStep=0.001;
k=-0.5;
contactInteractionFactor=4 Pi k;
max=5ho;
min=-max;
initialize[];
gridReport

moments:=Join[Array[mean,dim],Array[sigma,dim],Table[mom[n,4],{n,dim}],Table[mom[n,6],{n,dim}],{norm}];
momentslegend:=Join[Table["mean"<>dir[i],{i,dim}],Table["sigma"<>dir[i],{i,dim}],
Table["mom4"<>dir[i],{i,dim}],Table["mom6"<>dir[i],{i,dim}],{"norm"}];

calcAllTable
p1=plotProjectionsOpts[PlotStyle->Orange];

myPrint=Null&;
AbsoluteTiming[evolve["ite",5,10]]

calcAllTable
p2=plotProjectionsOpts[PlotStyle->Black]

Show[p1,p2,PlotRange->{{-3,3},All}]

reportEnergies

threeBodyLosses=True;
lossRate=0.05;
threeBodyLossesFactor=(-I/2)lossRate;
tmax=10;
reset;
myPrint=Null&;
AbsoluteTiming[evolve["rte",tmax,100,100]]

plotM["sigmaX"]

plotM["mom4X"]

plotM["mom6X"]

plotM["A0"]

plotM["norm"]

Show[lplots[[1;;;;10,2]],PlotRange->{{-1,1},All}]