# Set of ODEs with Integrate inside NDsolve?

I have set of ODEs with integration and to be solved with NDsolve. But because of existance of Exp[-2 qr[z] r^2], I can not find the solution. By the way, if Exp[-2 qr[z] r^2] is replaced by Exp[-2 r^2], or (TeV*Exp[-2 qr[z] r^2])^(-3/2) is replaced by (TeV*Exp[-2 qr[z] r^2])^(-2), the equations is solved. Why? Could anyone know how to solve it?

 ClearAll["Global*"];
Z = 5/2; λ = 1; a0 = 0.01; σ = 1; n = 100; B = 10; TeV = 100;
νe[z_, r_] := 1.72*Z*n*λ^-1*(TeV*Exp[-2  r^2])^(-3/2)
ϵr[z_, r_] := 1 - (n (1 - σ*B))/((1 - σ*B)^2 + νe[z, r]^2)
ϵi[z_, r_] := (n νe[z, r])/((1 - σ*B)^2 + νe[z,r]^2)

K1[z_] := Integrate[((1 - ϵr[z, r])*Cos[qi[z] r^2] -ϵi[z, r]*Sin[qi[z] r^2]) Exp[-qr[z] r^2] r, {r, 0, 100}];
K2[z_] := Integrate[((1 - ϵr[z, r])*Sin[qi[z] r^2] + ϵi[z, r]*Cos[qi[z] r^2]) Exp[-qr[z] r^2] r, {r, 0, 100}];
K0[z_] :=Integrate[ϵi[z, r]*Exp[-2 qr[z]*r^2] r, {r, 0, 100}];

equ = {
-ν[z] D[qr[z], z] + 2 qr[z] D[ν[z], z] + 4 K0[z] qr[z]^2ν[z] == 0,
ν[z] D[qi[z], z] - qi[z] D[ν[z], z] + ν[z] (K1[z] (qi[z]^2 - qr[z]^2) - 2 K2[z] qi[z] qr[z]) == 0,
ν[z] D[qr[z], z] - qr[z] D[ν[z], z] + ν[z] (K2[z] (qi[z]^2 - qr[z]^2) + 2 K1[z] qi[z] qr[z]) == 0
};
ic = {ν[0] == 1, qr[0] == 1, qi[0] == 0};

{νsol, qrsol, qisol} =NDSolveValue[{equ, ic}, {ν, qr, qi}, {z, 0, 100}]
Plot[{νsol[z], qrsol[z], qisol[z]}, {z, 0, 1},PlotTheme -> "Detailed"]

• It looks like you're splitting the equations into their real and imaginary parts. Have you tired writing the whole thing in terms of complex-valued functions for $\nu$, $K$, $q$, etc.? Mathematica can handle complex-valued functions in NDSolve quite easily, and leaving the quantities in terms of arbitrary complex-valued functions might help Mathematica perform the integrations you're asking it to perform. Jun 15, 2021 at 14:54
• Yes, \[Epsilon]r and \[Epsilon]i are real and imaginary part of \[Epsilon], but they do stay alone in the ODEs. So I will not use complex-valued function. Thanks. Jun 15, 2021 at 15:10
• You try to solve system of integrodifferential equations. It is not solvable problem in the current version of NDSolve. Actually we need to develop some algorithm based on FDM and collocation method. Jun 15, 2021 at 18:31
• What a sad information! FDM and collocation method are too complicated and require lot of coding. Thanks. Jun 16, 2021 at 0:27

## 1 Answer

This system of integrodifferential equations can be solved with using Gauss quadrature rule for integral approximation as follows (we made some simplification)

ClearAll["Global*"]; Needs["DifferentialEquationsNDSolveProblems"];
Needs["DifferentialEquationsNDSolveUtilities"]; \
Get["NumericalDifferentialEquationAnalysis"];

np = 200; g = GaussianQuadratureWeights[np, 0, 100]; points =
g[[All, 1]]; weights = g[[All, 2]];
Z = 5/2; \[Lambda] = 1; a0 = 0.01; \[Sigma] = 1; n = 100; B = 10; TeV \
= 100;
\[Nu]e = 1.72*Z*n*\[Lambda]^-1*(TeV^(-3/2)*Exp[3 qr[z] r^2]);
\[Epsilon]r = 1 - (n (1 - \[Sigma]*B))/((1 - \[Sigma]*B)^2 + \[Nu]e^2);
\[Epsilon]i = (n \[Nu]e)/((1 - \[Sigma]*B)^2 + \[Nu]e^2);
K1 = Sum[((1 - \[Epsilon]r)*Cos[qi[z] r^2] - \[Epsilon]i*
Sin[qi[z] r^2]) Exp[-qr[z] r^2] r weights[[i]] /.
r -> points[[i]], {i, Length[points]}];
K2 = Sum[((1 - \[Epsilon]r)*Sin[qi[z] r^2] + \[Epsilon]i*
Cos[qi[z] r^2]) Exp[-qr[z] r^2] r weights[[i]] /.
r -> points[[i]], {i, Length[points]}];
K0 = Sum[\[Epsilon]i*Exp[-2 qr[z]*r^2] r weights[[i]] /.
r -> points[[i]], {i, Length[points]}];

equ = {-\[Nu][z] D[qr[z], z] + 2 qr[z] D[\[Nu][z], z] +
4 K0 qr[z]^2 \[Nu][z] ==
0, \[Nu][z] D[qi[z], z] -
qi[z] D[\[Nu][z], z] + \[Nu][
z] (K1 (qi[z]^2 - qr[z]^2) - 2 K2 qi[z] qr[z]) ==
0, \[Nu][z] D[qr[z], z] -
qr[z] D[\[Nu][z], z] + \[Nu][
z] (K2 (qi[z]^2 - qr[z]^2) + 2 K1 qi[z] qr[z]) == 0};
ic = {\[Nu][0] == 1, qr[0] == 1, qi[0] == 0};

{\[Nu]sol, qrsol, qisol} =
NDSolveValue[{equ, ic}, {\[Nu], qr, qi}, {z, 0, 100},
Method -> {"EquationSimplification" -> "Residual"}]

Plot[{\[Nu]sol[z], qrsol[z], qisol[z]}, {z, 0, 1},
PlotTheme -> "Detailed", PlotRange -> All]
`

• Thanks a lot!!! It works even after a revision of the equations. By the way, is it the universal skill to numerically solve the complicated ODE and PDE in mma? For further learning, do you have some suggestions? Jun 16, 2021 at 5:46
• @sixpenny This forum is a very good resource to study numerical approach with Mathematica. Jun 16, 2021 at 10:35