# NDSolve for Complex Algebraic-Differential Equation

Let consider the following complex equations: $$\frac{i x(s)}{2\pi} - \frac{\log(1+ e^{-y(s)}) - \log(1 + e^{y(s)})}{\dot x(s)}$$ which I will cal eq1 and $$\frac{2 \left(\text{Li}_2\left(-e^{-y(s)}\right)+\text{Li}_2\left(-e^{y(s)}\right)-\frac{3 \pi ^2}{2}\right) x''(s)}{x'(s)^3}-\frac{y'(s) \log \left(e^{-y(s)}+1\right)-y'(s) \log \left(e^{y(s)}+1\right)}{x'(s)^2}+\frac{i y(s)}{2 \pi }$$ which I will call eq2. (btw, notice that since we are in the complex domain I didn't use on purpose sum/subraction formula for polylogs, since in principle I do not know in which branch I'm in) More specifically

eq1 = -((I*x[s])/(2*Pi)) - (2*(Log[1 - E^(-y[s])] - Log[1 + E^y[s]]))/Derivative[1][x][s]
eq2 = -((I*y[s])/(2*Pi)) + (2*(-Log[1 - E^(-y[s])] + Log[1 + E^y[s]])*Derivative[1][y][s])/Derivative[1][x][s]^2 + (4*(-(Pi^2/3) + PolyLog[2,E^(-y[s])] +PolyLog[2, -E^y[s]])*Derivative[2][x][s])/Derivative[1][x][s]^3


I would like to solve numerically them in the interval $$s\in (-1/2,1/2)$$ with some initial condition, say $$x(0)=0$$, $$y(0)=0$$ and $$\dot x(0) = 1$$. I try to input them in

NDSolve[{eq1 == 0, eq2 == 0, x[0] == 0, y[0] == 0, x'[0] == 1}, {x, y}, {s, -1/2, 1/2}]


I get the following warnings:

NDSolve::ntdvdae: Cannot solve to find an explicit formula for the derivatives. NDSolve will try solving the system as differential-algebraic equations.
NDSolve::mconly: For the method IDA, only machine real code is available. Unable to continue with complex values or beyond floating-point exceptions.
NDSolve::mconly: For the method IDA, only machine real code is available. Unable to continue with complex values or beyond floating-point exceptions.


I think that the first warns me about the method it will use to try to solve the equations, and the second and the third tell me that in order to use IDA method I have to use real valued functions. In order to deal with this I defined

eq1r = (eq1 /. {x -> (xr[#] + I xi[#] &), y -> (yr[#] + I yi[#] &)}) //Re // ComplexExpand;
eq1i = (eq1 /. {x -> (xr[#] + I xi[#] &), y -> (yr[#] + I yi[#] &)}) //Im // ComplexExpand;
eq2r = (eq2 /. {x -> (xr[#] + I xi[#] &), y -> (yr[#] + I yi[#] &)}) //Re // ComplexExpand;
eq2i = (eq2 /. {x -> (xr[#] + I xi[#] &), y -> (yr[#] + I yi[#] &)}) //Im // ComplexExpand;


and then, when I try to execute

NDSolveValue[{eq1r == 0, eq1i == 0, eq2r == 0, eq2i == 0, xr[0] == 0, xi[0] == 0, yr[0] == 0, yi[0] == 0, xr'[0] == 1, xi'[0] == 0}, {xr, xi, yr, yi}, {s, -1/2, 1/2}]


apart for the usual warning I get some other messages like

NDSolveValue::nrnum1: The function value 2. -3.10849*10^-8 I is not a real number when the arguments are {-1.95313*10^-7,-7.51414*10^-33+0. I,3.84724*10^-26+0. I,-1.95313*10^-7+0. I,1. +0. I,3.10849*10^-8+0. I,-1.19591*10^-33+0. I,1,1}.


which I really cannot understand since I explicitly use Re and Im. I think that maybe I have been too naive.

## QUESTION

Do you have any suggestion to numerically solve these equations with Mathematica?

• Add code for eq1,eq2. Oct 25, 2019 at 9:39
• @AlexTrounev. I added the code.
– MaPo
Oct 25, 2019 at 9:48
• This is a misprint PolyLog[2,E^-Y[s]]? Should it be PolyLog[2,E^-y[s]]? Oct 25, 2019 at 10:32
• Your code for eq1, eq2 does not work: it needs to be posted it in InputForm. See this meta Q&A for help. Also you have both y and Y -- You might try differentiating the algebraic equation. Oct 25, 2019 at 10:36
• Sorry for the misprint. I edited the code for eq1 and eq2
– MaPo
Oct 25, 2019 at 10:48

## 1 Answer

The first equation can be solved with respect to y[s]

Solve[eq1 == 0, y[s]]
(*{{y[s] ->
Log[1/2 (-1 + E^((I x[s] Derivative[1][x][s])/(4 \[Pi])) - Sqrt[
1 - 6 E^((I x[s] Derivative[1][x][s])/(4 \[Pi])) + E^((
I x[s] Derivative[1][x][s])/(2 \[Pi]))])]}, {y[s] ->
Log[1/2 (-1 + E^((I x[s] Derivative[1][x][s])/(4 \[Pi])) + Sqrt[
1 - 6 E^((I x[s] Derivative[1][x][s])/(4 \[Pi])) + E^((
I x[s] Derivative[1][x][s])/(2 \[Pi]))])]}}*)


We define two equations on two branches

eq3 = eq2 /. {y[s] ->
Log[1/2 (-1 + E^((I x[s] Derivative[1][x][s])/(4 \[Pi])) - Sqrt[
1 - 6 E^((I x[s] Derivative[1][x][s])/(4 \[Pi])) + E^((
I x[s] Derivative[1][x][s])/(2 \[Pi]))])],
y'[s] ->
D[Log[1/2 (-1 + E^((I x[s] Derivative[1][x][s])/(4 \[Pi])) - Sqrt[
1 - 6 E^((I x[s] Derivative[1][x][s])/(4 \[Pi])) + E^((
I x[s] Derivative[1][x][s])/(2 \[Pi]))]), s]]};
eq4 = eq2 /. {y[s] ->
Log[1/2 (-1 + E^((I x[s] Derivative[1][x][s])/(4 \[Pi])) + Sqrt[
1 - 6 E^((I x[s] Derivative[1][x][s])/(4 \[Pi])) + E^((
I x[s] Derivative[1][x][s])/(2 \[Pi]))])],
y'[s] ->
D[Log[1/2 (-1 + E^((I x[s] Derivative[1][x][s])/(4 \[Pi])) + Sqrt[
1 - 6 E^((I x[s] Derivative[1][x][s])/(4 \[Pi])) + E^((
I x[s] Derivative[1][x][s])/(2 \[Pi]))]), s]]};


We exclude the point s = 0 from the solution. We solve the equations eq,e4 on the segments {s,-1/2,-x0}, {s,x0,1/2}, where x0 =10^-10. For example

x0 = 10^-10; X2 =
NDSolve[{eq3 == 0, x[x0] == 0 + I 0, x'[x0] == 1 + I 0},
x, {s, x0, 1/2}]

{Plot[Evaluate[ReIm[x[s] /. X2]], {s, x0, 1/2}],
Plot[Evaluate[
ReIm[Log[
1/2 (-1 + E^((I x[s] Derivative[1][x][s])/(4 \[Pi])) + Sqrt[
1 - 6 E^((I x[s] Derivative[1][x][s])/(4 \[Pi])) + E^((
I x[s] Derivative[1][x][s])/(2 \[Pi]))])]] /. X2], {s, x0,
1/2}]}