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$
7
  • $\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$ Commented Aug 27, 2015 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$ Commented Aug 27, 2015 at 11:01
  • $\begingroup$ A pretty good solution is obtained with initial conditions v[1] == 5, v'[1] == -5.740226373077693`50. $\endgroup$ Commented Aug 27, 2015 at 11:08
  • $\begingroup$ to Guess who it is, the warning message is :Power::infy: Infinite expression 1/0.^3 encountered. >> $\endgroup$
    – zeo
    Commented Aug 27, 2015 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
    Commented Aug 27, 2015 at 12:07

2 Answers 2

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$
1
  • $\begingroup$ Mariusz, So many thanks for your wonderful answer. $\endgroup$
    – zeo
    Commented Aug 28, 2015 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$
1
  • $\begingroup$ Splendidly and Simple :) $\endgroup$ Commented Aug 28, 2015 at 16:48

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

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