Message Boards Message Boards

0
|
5590 Views
|
3 Replies
|
1 Total Likes
View groups...
Share
Share this post:

Integrating second order differential equation

Posted 10 years ago

Hi everyone, I'm a new user of Mathematica, and i'm having difficulties in solving a particular problem. My aim is to integrate this differential equation: $$ \frac{d^{2 }U}{dx^{2}} + \frac{2}{x} \frac{dU}{dx} = - \frac{I(U)}{I(U_{0})} - \lambda \frac{I(\mu U)}{I(\mu U_{0})}$$ with boundary conditions: $$ U(0)=U_{0} \; \frac{dU}{dx}=0 \; \mathrm{at} \;x=0$$ where $I(U)$ is defined as: $$I(U)= \frac{3}{4} \sqrt{\pi} -\frac{2}{3} e^{U} \Gamma(\frac{5}{2}, U) $$ where $\Gamma(a,U) $ is the incomplete gamma function. In Mathematica I tried this code, with given values for the parameters above:

f[x_]:= 3/4 Sqrt[Pi] - 2/3 Exp[x]* Gamma[5/2, x];
g[x_, p_, l_, m_]:= (f[x]/f[p]) + l* (f[m*x]/f[m*p]);

DSolve[{y''[x] + (2/x) y'[x] + g[y[x], 2, 2, 3]==0, y[0]==2, y'[0]==0}, y, {x,0, 1}]

but i get the following error: the function y appears with no arguments. Could you please give me any help? Please, also consider that my final aim is then to plot $g(U(x))$ as defined at different values of $U_{0}$, $\lambda$ and $\mu$. Thank you.

POSTED BY: Filippo Cova
3 Replies

Hello, I believe that you may have previous definitions of f,g, or y. When I try your code from a fresh kernel, I get the DSolve back, indicating that Mathematica failed to find a closed form solution.

Try:

Clear[f,g,y];
DSolve[{y''[x] + (2/x) y'[x] + g[y[x], 2, 2, 3]==0, y[0]==2, y'[0]==0}, y, {x,0, 1}]

I believe you may have to resort to NSolve:

Clear[y, f, g];
f[x_] := 3/4 Sqrt[Pi] - 2/3 Exp[x]*Gamma[5/2, x];
g[x_, p_?NumericQ, l_?NumericQ, m_?NumericQ] := (f[x]/f[p]) + 
   l*(f[m*x]/f[m*p]);
nsol = y /. 
  NDSolve[{y''[x] + (2/x) y'[x] + g[y[x], 2, 2, 3] == 0, y[1] == 2, 
     y'[1] == 0}, y, {x, 1, 2}][[1]]
Plot[nsol[t], {t, 1, 2}]

I changed the domain of x because there is a singularity at x=0. You could generalize nsol to take arguments for numerical values of p,l,and m.

Note the use of the ?NumericQ pattern restriction is because NSolve only plays nice when it calls functions that are defined for numerical values.

POSTED BY: W. Craig Carter
Posted 10 years ago

Thank you very much for you help, I tried your suggestions. However, my specific problem requests conditions for x=0, is there a way to do this avoiding the singularity?

POSTED BY: Filippo Cova

Hello, There is probably a clever way to do this by looking at the limits as $t_{init}\rightarrow 0$, but this will give you a visual clue about the behavior:

Clear[nsol]

nsol[init_ ] := 
 y /. NDSolve[{y''[x] + (2/x) y'[x] + g[y[x], 2, 2, 3] == 0, 
     y[init] == 2, y'[init] == 0}, y, {x, init, 2}][[1]]

Manipulate[
 Plot[nsol[init][t], {t, init, 2}, PlotRange -> {{0, 2}, All}], {init,
   0.0001, .1}]

or,

nsol2[eps_ ] := 
 y /. NDSolve[{ y''[x] + (2/(x + eps)) y'[x] + g[y[x], 2, 2, 3] == 0, 
     y[0] == 2, y'[0] == 0}, y, {x, 0, 2}][[1]]

Manipulate[
 Plot[nsol2[eps][t], {t, 0, 2}, 
  PlotRange -> {{0, 2}, All}], {{eps, 0.01}, .0001, .1}]
POSTED BY: W. Craig Carter
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Reply Preview
Attachments
Remove
or Discard

Group Abstract Group Abstract