0
$\begingroup$

Here is the code:

a=4;b=10;
s = NDSolve[{y^2 v[y]^3 v''[y] + v[y]^3 y v'[y] + 2 y^2 v[y]^2 v'[y]^2 - a v[y]^4 +b== 0, 
v[1.0] == 5.0, v'[50.0] == 0.0}, v, {y, 1.0, 50.0}]

v1[y_] = Evaluate[v[y] /. s];

LogLinearPlot[v1[y], {y, 1, 50}]

For a=4; b=10; mathematica can't solve this equation. However, if I change the values of a and b (say a=1; b=1;), I can obtain the results. For b=0, mathematic can solve this equation with any given positive value of a. What's the problem here? Thanks!

$\endgroup$
  • $\begingroup$ "can't solve" - what errors are you getting? Anyway, try replacing the inexact numbers (e.g. 50.0) with exact ones (e.g. 50) and report back. $\endgroup$ – J. M. will be back soon Aug 27 '15 at 10:37
  • 1
    $\begingroup$ Your equation has some sort of singularity in the region of v[1] == 5, v'[1] == -5.740225. For much more negative values of $v'[1]$, the equations become non-amenable to solution. You seem to be trying to solve the DE at a singularity. $\endgroup$ – Patrick Stevens Aug 27 '15 at 11:01
  • $\begingroup$ A pretty good solution is obtained with initial conditions v[1] == 5, v'[1] == -5.740226373077693`50. $\endgroup$ – Patrick Stevens Aug 27 '15 at 11:08
  • $\begingroup$ to Guess who it is, the warning message is :Power::infy: Infinite expression 1/0.^3 encountered. >> $\endgroup$ – zeo Aug 27 '15 at 12:04
  • $\begingroup$ to Patrick Stevens, the problem is that I don't the value of v'[y] for y=1, but for larger y, say y>10, v'[y]=0. $\endgroup$ – zeo Aug 27 '15 at 12:07
3
$\begingroup$

NDSolve It has problems with equation.Why? I don't No.Solution with $MAPLE.$

enter image description here

Let's try analytically.

$$-a v(y)^4+b+y^2 v(y)^3 v''(y)+2 y^2 v(y)^2 v'(y)^2+y v(y)^3 v'(y)=0 \tag{1}$$

By making the substitution $y=exp(t)$. When subbing for the independent variable, we need to make careful use of the chain rule to express $v(y),v′(y),v″(y)$ in terms of $v(t),v′(t),v″(t)$.

 ode = y^2*v[y]^3*v''[y] + v[y]^3 *y*v'[y] + 2 y^2*v[y]^2*v'[y]^2 - 
 a*v[y]^4 + b;
 subs = Exp[t]; newode = 
 ode /. {v[y] -> v[t], v'[y] -> v'[t]/D[subs, t], 
 v''[y] -> (v''[t] - D[subs, {t, 2}] v'[t]/D[subs, t])/D[subs, t]^2,
 y -> subs};
 newode2 = Expand[Collect[FullSimplify[newode], v''[t]]]

$-a v(t)^4+b+v(t)^3 v''(t)+2 v(t)^2 v'(t)^2=0 \tag{2}$

Using the substitution $v'(t)=k(v)$,where $v$ plays the role of the independent variable, and taking into account the relations $v''(t)=k'(t)=k'(v)*v'(t)=k'(v)*k(v)$, can reduce to a first-order equation:

newode3 = newode2 /. {v'[t] -> k[v], v''[t] -> k[v]*k'[v], v[t] -> v}

$-a v^4+b+v^3 k(v) k'(v)+2 v^2 k(v)^2=0 \tag{3}$

ODE = DSolve[newode3 == 0, k[v], v]

$\left\{k(v)\to -\frac{\sqrt{a v^6-3 b v^2+3 c_1}}{\sqrt{3} v^2},k(v)\to \frac{\sqrt{a v^6-3 b v^2+3 c_1}}{\sqrt{3} v^2}\right\} \tag{4}$

then $v'(t)=k(v)$ and $t=Log(y)$ :

{Integrate[1/k[v] /. ODE[[1]], v] == Log[y] + C[2], 
Integrate[1/k[v] /. ODE[[2]], v] == Log[y] + C[2]}

We have analyticall solution:

$\left\{-\sqrt{3} \int \frac{v^2}{\sqrt{a v^6-3 b v^2+3 c_1}} \, dv=c_2+\log (y),\sqrt{3} \int \frac{v^2}{\sqrt{a v^6-3 b v^2+3 c_1}} \, dv=c_2+\log (y)\right\}$

$c_1$ and $c_2$ integration constants. Mathematica can't find this Integrated's $:($

Try's with Series.

We need a equation $(2)$.By making the substitution $y=exp(t)$ and solve with parametr $t$ $t=Log(y)$ we have a new boundary conditions:

$$\{v(0)=5,v(\log (50))=0\}$$

a = 4;
b = 10;
n = 13;
v = Sum[c[i]*t^i, {i, 0, n}] + O[t]^(n + 1);
ODE = v^3*D[v, {t, 2}] + 2*D[v, t]^2*v^2 - a*v^4 + b == 0;
V = FullSimplify@
Normal[v /. Solve[LogicalExpand[ODE], Table[c[i], {i, 1, n}]]] // 
Quiet;
Sol = First@
Solve[{V == 5 /. t -> 0, D[V, t] == 0 /. t -> Log[50]}, {c[0], 
c[1]}];
VV = V /. Sol /. t -> Log[y] // N

LogLinearPlot[VV, {y, 1, 5}, PlotRange -> All, 
PlotLegends -> {"Series"}, PlotStyle -> {Thick, Dashed}]

Solution is not accurates ,you must increase $n=20$,mayby $50$.

enter image description here

$\endgroup$
  • $\begingroup$ Mariusz, So many thanks for your wonderful answer. $\endgroup$ – zeo Aug 28 '15 at 14:20
1
$\begingroup$

This equation also can be solved using the substitution v[y] -> u[y]^(1/4).

Clear[a, b];
eq = y^2 v[y]^3 v''[y] + v[y]^3 y v'[y] + 2 y^2 v[y]^2 v'[y]^2 - a v[y]^4 + b;
First@Solve[{u''[y] == D[v[y]^4, {y, 2}], u'[y] == D[v[y]^4, {y, 1}]}, {v''[y], v'[y]}];
equ = eq /. % /. v[y] -> u[y]^(1/4)
(* 10 - 4*u[y] + (y*Derivative[1][u][y])/4 + (y^2*Derivative[1][u][y]^2)/(8*u[y]) + 
   (y^2*(-3*Derivative[1][u][y]^2 + 4*u[y]*Derivative[2][u][y]))/(16*u[y]) *)

a = 4; b = 10;
su = NDSolveValue[{equ == 0, u[1] == 5^4, u'[50] == 0}, u, {y, 1, 50}, 
    WorkingPrecision -> 30];
Plot[su[y]^(1/4), {y, 1, 50}, PlotRange -> All, AxesLabel -> {v, y}]

enter image description here

Addendum: Method -> "Shooting"

Method -> "Shooting" also works, but here only with an excellent initial guess.

sv = NDSolveValue[{eq == 0, v[1] == 5, v'[50] == 0}, v, {y, 1, 50}, 
  Method -> {"Shooting", "StartingInitialConditions" -> {v[50] == 1.25745, v'[50] == 0}}]
$\endgroup$
  • $\begingroup$ Splendidly and Simple :) $\endgroup$ – Mariusz Iwaniuk Aug 28 '15 at 16:48

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.