1

I am numerically solving this ODE, and I do get a solution. I then need to use this solution to find derivatives and integrals of it, so I define functions using this solution. But when I try to obtain a parametric plot, I do not get a result. I am sure that the problem is with the part that I am Nintegrating but I do not understand the reason because I am able to plot that separately.

sol2 = NDSolve[{1 - Derivative[1][r][t]^2 - r[t]*Derivative[2][r][t] - r[t]*(1 - Derivative[1][r][t]^2)^(1/2) == 0, Derivative[1][r][0] == 1, Derivative[1][r][1] == 0}, r, {t, 0, 1}]
func[t_] := r[t] /. First[sol2]
func1[u_] := (D[r[t] /. sol2, t] /. t -> u)^2
ParametricPlot[{func[t] /. t -> u, NIntegrate[(1 - func1[y])^(1/2), {y, 0, u}]}, {u, 0, 1}]

Here's the result that I get when I plot that part separately:

Plot[NIntegrate[(1 - func1[y])^(1/2), {y, 0, u}], {u, 0, 1}]

enter image description here

juv95
  • 351
  • 1
  • 7

1 Answers1

2

We can use NumericQ to evaluate integral as follows

sol2 = NDSolve[{1 - Derivative[1][r][t]^2 - r[t]*Derivative[2][r][t] -
       r[t]*(1 - Derivative[1][r][t]^2)^(1/2) == 0, 
    Derivative[1][r][0] == 1, Derivative[1][r][1] == 0}, r, {t, 0, 1}];
func[t_] := r[t] /. First[sol2];
func1[u_] := (D[r[t], t] /. sol2[[1]] /. t -> u)^2;
g[u_?NumericQ] := NIntegrate[(1 - func1[y])^(1/2), {y, 0, u}];

ParametricPlot[{func[t] /. t -> u, g[u]}, {u, 0, 1}, Frame -> True, 
 TicksStyle -> Directive["Label", 8], AspectRatio -> 1/2]

Figure 1

Note, that BVP solution with Neuman boundary condition on two ends computed with NDSolve up to arbitrary constant of about $4.70033\times 10^9$ in this case. If you need func[t] of about 1, then BVP could be regularized like here.

Alex Trounev
  • 44,369
  • 3
  • 48
  • 106
  • I had worked on this problem and noticed the solution does not satisfy the boundary conditions: When taking the derivative of the solution via funcD[t_]=D[func[t],t], then funcD[0]=0.873 and funcD[1]=0.061 and increasing WorkingPrecision did not improve the results. – josh Sep 13 '22 at 10:02
  • @josh Yes, it is why constant $4.70033\times 10^9$ so big. – Alex Trounev Sep 13 '22 at 10:44