NDSolve "Coupled differential equation"

Question

I am trying to solve the given coupled differential equation using NDsolve.

$$ \begin{array}{c} n''(r)+\frac{n'(r) T'(r)}{T(r)}+\frac{n'(r)}{r}+\frac{n(r) T''(r)}{T(r)}-\frac{n(r) T'(r)^2}{T(r)^2}+\frac{n(r) T'(r)}{r T(r)}-n(r)+e^{-r^2}=0\\ T(r) n''(r)+2 n'(r) T'(r)+\frac{T(r) n'(r)}{r}+n(r) T''(r)+\frac{n(r) T'(r)}{r}-n(r)+e^{-r^2}=0 \end{array} $$

eq1 = Div[Grad[n[r] , {r, h, z}, "Cylindrical"], {r, h, z}, "Cylindrical"] + Div[(n[r]/T[r])*Grad[T[r], {r, h, z}, "Cylindrical"], {r, h, z}, "Cylindrical"] - n[r] + Exp[-(r^2)] == 0 ; 

eq2 = Div[T[r]*Grad[n[r], {r, h, z}, "Cylindrical"], {r, h, z}, "Cylindrical"] + Div[n[r]*Grad[T[r] , {r, h, z}, "Cylindrical"], {r, h, z}, "Cylindrical"] - n[r] + Exp[-(r^2)] == 0 ;

initialConditions = {n'[50] == 0, n'[0] == 0, T'[50] == 0, T'[0] == 0};

solution = NDSolve[{eq1, eq2, initialConditions}, {n, T}, {r, 0, 50}];

Plot[Evaluate[{n[r], T[r]} /. %], {r, 0, 50}, PlotLegends -> {"n[r]", "T[r]"}]

Mathematica returns an error - NDSolve::bvdae: Differential-algebraic equations must be given as initial value problems

I have two questions regarding it,
Q1: Can we solve such equations accurately with NDSolve?
Q2: How do I get initial values for the two quantities, n and T?

Answer 1

The pair of ODEs can be solved as follows. Begin by eliminating the second-derivative terms from eq2.

eq3 = Simplify[T[r] Subtract @@ eq1 - Subtract @@ eq2] == 0

(* E^-r^2 (-1 + T[r] - E^r^2 n'[r] T'[r]) + n[r] (1 - T[r] - T'[r]^2/T[r]) == 0 *)

Because n[r] and its derivatives enter the ODEs at most linearly, they can be eliminated to obtain an ODE for T[r] as follows.

Solve[{eq3, D[eq3, r], D[eq3, r, r]}, {n'''[r], n''[r], n'[r]}] // 
    Simplify // Flatten;
Simplify[{eq2, D[eq2, r]} /. %];
eqs = Equal @@@ (Simplify@Solve[%, {T'''[r], n[r]}] // Flatten)

(* {T’’’[r] == (-4 r^2 (-1 + r^2) T[r]^5 + r T[r]^4 (12 r (-1 + r^2) 
    + (-7 + 8 r^2) T’[r] + r (5 - 4 r^2) T’’[r]) + 
    r T’[r]^2 ((-2 + 4 r^2) T’[r] + 4 r T’[r]^2 + r (-1 + 2 T’’[r])) + 
    T[r] ((2 + 4 r^2 - 4 r^4) T’[r]^2 - 8 r^3 T’[r]^3 + 2 r^2 T’[r]^4 - 
    r^2 (-1 + T’’[r]) T’’[r] + T’[r] (r - 2 r^3 - (r + 2 r^3) T’’[r])) + 
    T[r]^3 (2 (1 - 5 r^2 + 2 r^4) T’[r]^2 + 
    r T’[r] (15 - 18 r^2 + (-1 + 10 r^2) T’’[r]) + 
    r^2 (-12 (-1 + r^2) + (-9 + 8 r^2) T’’[r] - T’’[r]^2)) + 
    T[r]^2 ((2 r - 8 r^3) T’[r]^3 + T’[r]^2 (-4 + 7 r^2 - 2 r^2 T’’[r]) + 
    r T’[r] (-9 + 12 r^2 + (2 - 8 r^2) T’’[r]) + 
    r^2 (4 (-1 + r^2) + (3 - 4 r^2) T’’[r] + 
    2 T’’[r]^2)))/(r^2 (-1 + T[r])^2 T[r] (2 r T[r] - T’[r])), 
    n[r] == (E^-r^2 T[r] (r + r T[r]^2 + (1 - 2 r^2) T’[r] - 
    2 r T’[r]^2 + T[r] ((-1 + 2 r^2) T’[r] + r (-2 + T’’[r])) - 
    r T’’[r]))/(r T[r]^3 - r T’[r]^2 + T[r]^2 (-T’[r] + r (-2 + T’’[r])) + 
    T[r] (r + T’[r] - r T’[r]^2 - r T’’[r]))} *)

An important point is that it is a third-order system, so only three boundary conditions can be specified, not four as in the question. In the absence of suitable boundary conditions, I offer the following solution with made-up initial conditions:

NDSolveValue[{eqs, T[1] == 2, T'[1] == 0, T''[1] == 0}, {T[r], n[r]}, {r, 1, 3}];
Plot[%, {r, 1, 3}, AxesLabel -> {r, "T,n"}, LabelStyle -> {12, Bold, Black}]

Solution without eqs

An even simpler approach, at least for the sample solution, is

NDSolveValue[{eq2, eq3, T[1] == 2, T'[1] == 0, n'[1] == -.7}, 
    {T[r], n[r]}, {r, 1, 3}];

which produces an answer very close to that plotted above.

Answer 2

This seems an overdetermined system.You get

     res1 = Solve[eq1 , T''[r] ][[1]]

$$\left\{T''(r)\to \frac{e^{-r^2} \left(-e^{r^2} r T(r)^2 n''(r)-e^{r^2} r T(r) n'(r) T'(r)-e^{r^2} T(r)^2 n'(r)-e^{r^2} n(r) T(r) T'(r)+e^{r^2} r n(r) T'(r)^2+e^{r^2} r n(r) T(r)^2-r T(r)^2\right)}{r n(r) T(r)}\right\}$$

But there is no solution for n''

    Simplify[eq2 /. res1 ] // FullSimplify

$$n'(r) T'(r)+n(r) \left(\frac{T'(r)^2}{T(r)}+T(r)-1\right)=e^{-r^2} (T(r)-1)$$

Using this equation as a replacement for n', there remains an undefined n in the equation for T''. Conclusion: There is an error in the algebraic ansatz.