How to take derivative of the argument of an interpolating function

Question

I am trying to plot the derivative of the argument of an interpolating function u[t, x] with respect to $x$. Here, u[t, x] is the solution of Non-linear Schrodinger equation. It seems that there is a problem while taking the derivative of Arg[u[t, x]]. In fact, the derivative of Arg[u[t, x]] calculated for arbitrary $x,t$ gives a complex number, but it is supposed to be real.

Can anyone please help?

Here is the code I was trying.

sol = NDSolve[
       Evaluate[{I D[u[t, x], t] + D[u[t, x], x, x] + 
         1/2 Abs[u[t, x]]^2 u[t, x] == 0, u[0, x] == 2 Sech[2 x]}], 
       u, {t, 0, 10}, {x, -10, 10}]

g[t_, x_] = Evaluate[Arg[u[t, x]] /. sol]
h[t_, x_] = D[g[t, x], x]

DensityPlot[h[t, x], {t, 0, 10}, {x, -10, 10}]

For example, h[2, 4] returns {0.28502 + 0.459296 I}.

Answer 1

The main problem here seems to be that Mathematica is giving strange results when you ask it to take the derivative of the function Arg. As an example of this, suppose you take $\arg(e^{ix})$, where $x$ is a real number. You would expect that $\arg(e^{ix}) = x \; (\mathrm{mod} \; 2\pi)$, and therefore the derivative of this function should be a constant for most values of $x$. But this isn't what Mathematica does:

Plot[Evaluate[ReIm[D[Arg[Exp[I x]], x]]], {x, 0, 2}] 

You do get the expected results if you insert a ComplexExpand, and if your expression is something that simplifies when you use ComplexExpand:

Plot[Evaluate[ReIm[D[ComplexExpand[Arg[Exp[I x]]], x]]], {x, 0, 2}]

However, ComplexExpand can't do much with an arbitrary InterpolatingFunction, so this method isn't suitable here. There may be some well-justified complex analysis reason why the first snippet above returns the given plot (it appears that the result is $i e^{2ix}$), but I'm not sure what it is just now. Alternately, it may be a bug.

To work around this, you can rewrite your equations so that they solve for the real and imaginary parts of $u$ separately, and then use the usual sorts of formulas to find the phase of $u$ and its derivatives:

eqns = {I D[u[t, x], t] + D[u[t, x], x, x] + 
 1/2 Abs[u[t, x]]^2 u[t, x] == 0, u[0, x] == 2 Sech[2 x]};
expandedeqns = Flatten[ComplexExpand[
    ReIm[eqns /. {u[x_, t_] -> uR[x, t] + I uI[x, t], 
      Derivative[y__][u][x_, t_] -> 
       Derivative[y][uR][x, t] + I Derivative[y][uI][x, t]}]]]
sol = NDSolve[Evaluate[expandedeqns], {uR, uI}, {t, 0, 10}, {x, -10, 10}]

g[t_, x_] = ArcTan[uR[t, x], uI[t, x]] /. First[sol]
h[t_, x_] = D[g[t, x], x]
DensityPlot[h[t, x], {t, 0, 10}, {x, -10, 10}]

Finally, as noted by xzczd in the comments, you should not expect this answer to be correct unless you have fully specified the boundary conditions.

Answer 2

D[Arg[x + y*I], x]
(*Derivative[1][Arg][x + I y]*)

D[ComplexExpand[Arg[x + y*I], TargetFunctions -> {Re, Im}], x]
(*-(y/(x^2 + y^2))*)

But in case it is coming zero

D[ComplexExpand[Arg[u[t, x]], TargetFunctions -> {Re, Im}], x]
(*0*)

Answer 3

Another possible solution is, extract the data from InterpolatingFunction, modify it and build a new interpolating function:

solfunc = u /. sol[[1]];

solfuncarg = ListInterpolation[solfunc["ValuesOnGrid"] // Arg, solfunc["Coordinates"]]

h2 = Derivative[0, 1][solfuncarg]

The structure of InterpolatingFunction is discussed here.

Finally I'd like to emphasize once again that the missing of boundary condition (b.c.) is a serious problem and you really need to add proper b.c. to obtain the desired answer.