First I will give the relevant information for my question, and then I'll ask the question.$\large{\textrm{Background}}$For evolving the nonlinear Schrodinger equation (NLS), one typically uses [a variant of] the Split-Operator, in which one split the nonlinear and linear components of the differential equation and treats them seperately, one after another. The version of the NLS I'm concerned with is$$\frac{\partial \psi}{\partial t}=\frac{i}{2}\frac{\partial^2 \psi}{\partial x^2}+i|\psi|^2\psi=[\hat{D}+\hat{N}(x,t)]\psi$$where we have defined $$\psi=\psi(x,t),\,\,\hat{D}=\frac{i}{2}\frac{\partial^2}{\partial x^2},\,\,\hat{N}(x,t)=i|\psi|^2$$As outlined quite succinctly in the previously cited wiki page, the general evolution $\psi(x,t)\to\psi(x,t+\Delta t)$ for sufficiently small step size $\Delta t$ is given by$$\psi(x,t+\Delta t)=e^{\Delta t(\hat{D}+\hat{N})}\psi(x,t)$$To actually evolve this with a computer, we can split/factor the exponential operator in the line above to a obtain an algorithm that is correct to order $(\Delta t)^2$.$$\psi(x,t+\Delta t)\approx e^{\Delta t(\hat{D})}e^{\Delta t(\hat{N})}\psi(x,t)$$Realizing that the second derivative in the $\hat{D}$ operator turns into a multiplication by $-k^2$ in Fourier-space, we can perform this evolution on a computer by calculating, in a discretized fashion of course, $$\psi(x,t+\Delta t)\approx\mathcal{F}^{-1}\left(e^{-\frac{i}{2}\Delta t (k^2)}\mathcal{F}\left(e^{i\Delta t |\psi^2|}\psi(x,t)\right)[k,t]\right)[x,t]$$This can be interpreted as follows: we first evolve the nonlinear part of the equation, and then we take a time step $\Delta t$ by the evolution of the linear part of the equation (in Fourier space). With this scheme, one can start with an initial wave-form/function $\psi(x,t_0)$ and evolve forward in time.$\large{\textrm{Question}}$How does one use higher order algorithms on the NLS? By using the Baker-Campell-Hausdorff formula for factorizing exponentials of operators, one can come up with a simple third order algorithm,$$\psi(x,t+\Delta t)\approx e^{\frac{\Delta t}{2}\hat{D}}e^{\Delta t(\hat{N})}e^{\frac{\Delta t}{2}\hat{D}}\psi(x,t)$$However, if we try to attack this in the same way we did before, we now have to evaluate $\hat{N}(x,t+\frac{\Delta t}{2})$, which involves evaluating $\psi(x,t+\frac{\Delta t}{2})$, which is the function we are trying to evolve! More precisely, if we only start off with an initial wave-form/function $\psi(x,t)$, then we don't have access to $\psi(x,t+\frac{\Delta t}{2})$. How can I resolve this problem.Edit: After writing this I realized that by evaluating $e^{\frac{\Delta t}{2}\hat{D}}\psi(x,t)$, it's possible we are obtaining an approximation for $\psi(x,t+\frac{\Delta t}{2})$, which would answer my question affirmatively. However I'm going to keep this question up in case I'm incorrect.
Evolving nonlinear Schrodinger equation with higher-order algorithms?
nonlinear equations;operator splitting
Your added edit to the question points exactly the way this needs to be done: You first evolve your field using a half-step with the $\hat D$ operator, then a full step with the $\hat N$ operator, and then another half step with $\hat D$. These are done one after the other, i.e., the result of one operation is the input for the next.If you're interested in finding more about this, the scheme you describe is commonly referred to as the Strang splitting.