# DSolve - incorrect solution to PDE

Bug introduced in 12.2, persisting through 12.3.1.

I'm trying to solve a simple Schrödinger equation in external field:

DSolve[{-Derivative[2, 0][psi][x,t]/2 + F Sin[π t]x  psi[x,t]==
I Derivative[0, 1][psi][x,t],psi[x,0]==Exp[-x^2]}, psi, {x,t},
Assumptions-> t ∈ Reals && F ∈ Reals && x ∈ Reals]


And the solution that Mathematica spits out appears to be incorrect:

E^(I x (x/(-I + 2 t) - F t Sin[π t]))/Sqrt[1 + 2 I t]


When I plug it to the original equation and FullSimplify I get:

(E^(I x (x/(-I+2 t)-F t Sin[π t])) F t (2 π (I-2 t) x Cos[π t]+Sin[π t]
(-4 x+F t (-I+2 t) Sin[π t])))/(Sqrt[1+2 I t] (-I+2 t))==0


which obviously doesn't appear to be zero for all F values (although it is for F=0).

Here's a screenshot from my notebook:

I'm using MMA 12.3.0 on MacOS.

• what's psi^(2,0) (etc.)? I think something might have gone wrong while copying... May 25, 2021 at 20:22
• ah, Derivative[2,0][psi] and Derivative[0,1][psi], probably, i see. May 25, 2021 at 20:24
• This is partial derivative alternative notion. Basically the same as D[psi[x,t],{x,2}], also psi^(0,1) would be D[psi[x,t],{t,1}] May 25, 2021 at 20:25
• right, it just doesn't parse as mathematica input syntax, so I couldn't copy-paste your code! I'll edit it so others can (and I checked to make sure it gives the same answer) May 25, 2021 at 20:29
• @thorimur thanks! BTW, Is there a way to easily copy compatible code from MMA? May 25, 2021 at 20:34

Since DSolve cannot solve the problem up to now i.e. v12.3.0, I'd like to add a (imperfect) work-around. The idea is similar to that in the linked paper. We first transform the problem to a initial value problem of 1st order PDE with the new-in-12.3 BilateralLaplaceTransform. (Notice bilateral Laplace transform is essentially a Fourier transform with special coefficient. )

With[{psi = psi[x, t]},
eq = -D[psi, x, x]/2 + F Sin[π t] x psi == I D[psi, t];
ic = psi == Exp[-x^2] /. t -> 0];

tsys = BilateralLaplaceTransform[{eq, ic}, x, s] /.
HoldPattern@BilateralLaplaceTransform[h_[x, t_], __] :> h[s, t] /. psi -> Ψ
(*
{(-(1/2)) s^2 Ψ[s, t] - F Sin[Pi t] Derivative[1, 0][Ψ][s, t] ==
I Derivative[0, 1][Ψ][s, t], Ψ[s, 0] == E^(s^2/4) Sqrt[Pi]}
*)


Then solve it with DSolve. DSolve spits out ifun warning, which isn't too surprising, because you've chosen the periodic F Sin[π t] as $$E(t)$$:

tsol = Ψ[s, t] /. DSolve[tsys, Ψ, {s, t}]
(*
{E^((π (I F + π s - I F Cos[π t])^2 -
I ArcSin[Cos[π t]] (F^2 - 2 (π s - I F Cos[π t])^2) +
1/2 I π (2 F^2 - 2 π^2 s^2 + 4 I F π s Cos[π t] +
F^2 Cos[2 π t]) + F (4 π s - 3 I F Cos[π t]) Sqrt[Sin[π t]^2])/(
4 π^3)) Sqrt[π],
E^((π (I F + π s - I F Cos[π t])^2 +
I ArcSin[Cos[π t]] (F^2 - 2 (π s - I F Cos[π t])^2) -
1/2 I π (2 F^2 - 2 π^2 s^2 + 4 I F π s Cos[π t] +
F^2 Cos[2 π t]) + I F (4 I π s + 3 F Cos[π t]) Sqrt[Sin[π t]^2])/(
4 π^3)) Sqrt[π]}
*)


The tsol is probably valid only for certain interval of $$t$$ (maybe $$0, if I have to guess), but let's proceed anyway. The last step in principle is to transform back with

InverseBilateralLaplaceTransform[tsol, s, x]


But sadly, InverseBilateralLaplaceTransform cannot handle tsol. Given the paper doesn't include solution in the time domain either, I think this is acceptable.

"OK, but how do you know the method is correct? " A rigorous validation isn't easy, but experimental validation is. Since definition of inverse bilateral Laplace transform of $$F(s)$$ is $$\frac{1}{2\pi\mathbb{i}} \int_{\gamma-\mathbb{i}\infty}^{\gamma+\mathbb{i}\infty}F(s)e^{st}ds$$, we can validate for certain $$(x,t)$$ as follows:

integrand2[x_, t_] = 1/(2 Pi I) I tsol[[2]] Exp[I w x] /. s -> I w;

(*validate the i.c.: *)
Integrate[integrand2[x, 0], {w, -Infinity, Infinity}]
(* Exp[-x^2] *)

(* validate for x == 1, t == 1/2 *)
eq /. psi -> integrand2 /. x -> 1 /. t -> 1/2 /. F -> 1 // Simplify;
Integrate[Subtract @@ %, {w, -Infinity, Infinity}]
(* 0 *)


BTW, if a simpler $$E(t)$$ is chosen, the symbolic solution without integral can be found. For $$E(t)=F t$$:

(* Still, DSolve cannot correctly handle this: *)
With[{psi = psi[x, t]},
eq = -D[psi, x, x]/2 + F t x psi == I D[psi, t];
ic = psi == Exp[-x^2] /. t -> 0];

tsys = BilateralLaplaceTransform[{eq, ic}, x, s] /.
HoldPattern@BilateralLaplaceTransform[h_[x, t_], __] :> h[s, t] /. psi -> Ψ

tsol = Ψ[s, t] /. DSolve[tsys, Ψ, {s, t}]

{f1[x_, t_], f2[x_, t_]} = InverseBilateralLaplaceTransform[tsol, s, x]
(*
{1/(E^((8 I F^(5/2) t^6 + 9 (F t^2)^(5/2) + 180 F^(3/2) t^2 x + 240 I (F t^2)^(3/2) x -
360 I Sqrt[F] x^2)/(360 (-I Sqrt[F] + 2 Sqrt[F t^2]))) Sqrt[
1 + (2 I Sqrt[F t^2])/Sqrt[F]]), E^((
8 I F^(5/2) t^6 - 9 (F t^2)^(5/2) + 180 F^(3/2) t^2 x - 240 I (F t^2)^(3/2) x -
360 I Sqrt[F] x^2)/(360 (I Sqrt[F] + 2 Sqrt[F t^2])))/Sqrt[
1 - (2 I Sqrt[F t^2])/Sqrt[F]]}
*)


This isn't the end. Substituting the solution back to the PDE, only f1 turns out to be correct for $$t>0$$:

eq /. psi -> f1 // Simplify[#, t > 0]&
(* True *)
eq /. psi -> f2 // Simplify[#, t < 0]&
(* True *)

• Nice solution, thank you! I'm accepting this as a workaround. May 26, 2021 at 15:57