# Solve an ODE involving an Interpolating function with Mathematica

Using Mathematica, I want to solve an ODE of an Interpolating Function. First I solve the ODE obtaining the function y[s] (or Y[s])

aa = .4; sol =  NDSolve[{-4 + 4/y[s]^2 - y[s]^2 + y[s]^4 + y'[s]^2 ==  0, y == 1.207}, y[s], {s, 0, aa}];
ParametricPlot [Evaluate[{s, y[s]} /. sol[]], {s, 0, aa},  PlotRange -> All]
Y[s_] = y[s] /. sol[];


(with the function Y[s] I can plot the solution y[s])

My problem is solving the ODE w'[s]= Sqrt[y[s]^2 + 2] with initial condition w = 0. My first attempt has been

w1 = NDSolve[{w'[s] == Sqrt[Y[s]^2 + 2], w == 0}, w[s], {s, 0, aa}]


Next, I put W[s_] = w[s] /. w1[] and I plot the function W[s] by means of Plot[W[s], {s, 0, aa}]. However, I think that it is wrong because the plot of W[s] seems to be a straight-line. In fact, I hope that Mathematica solves the ODE for w[s] by taking only the first value for y[s], namely, y=1.207 and then solves the ODE w'[s]=Sqrt[1.207^2+2], w=0 (whose solution is a linear function). I have computed the last ODE:

p1 = NDSolve[{p'[s] == Sqrt[1.207^2 + 2], p == 0},   p[s], {s, 0, aa}]


and the plots of both functions p[s] and w[s] are 'very similar'. Please, help.

## migrated from math.stackexchange.comJul 25 '18 at 11:07

This question came from our site for people studying math at any level and professionals in related fields.

The first equation has two branches of the solution, so you need to select a branch, define a function, and then solve the second equation. The solutions of the second equation can be compared for each branch with a linear function. It is obvious that the solutions differ among themselves and differ from the linear function.

aa = .417527; sol =
NDSolve[{-4 + 4/y[s]^2 - y[s]^2 + y[s]^4 + y'[s]^2 == 0,
y == 1.207}, y, {s, 0, aa}];
Y1[s_] := y[s] /. First[sol]
Y2[s_] := y[s] /. Last[sol]
sol1 = NDSolveValue[{w'[s] == Sqrt[2 + Y1[s]^2], w == 0},
w, {s, 0, aa}];
sol2 = NDSolveValue[{w'[s] == Sqrt[2 + Y2[s]^2], w == 0},
w, {s, 0, aa}];
sol3 = NDSolveValue[{w'[s] == Sqrt[2 + 1.207^2], w == 0},
w, {s, 0, aa}];
{Plot[{Y1[x], Y2[x]}, {x, 0, aa}, PlotLegends -> Automatic,
PlotRange -> All],
Plot[{sol1[x], sol2[x], sol3[x]}, {x, 0, aa},
PlotLegends -> Automatic, PlotRange -> All]} 