1
$\begingroup$

I have a second order ODE of the form :

u''[x]+ p[x]u'[x]+ q[x] u[x]==0 ;

with,

p[x]= -(((1.15+ 4.3 x - 6 x^2))/(-0.03333333333333333 + 1.15 x + 2.15 x^2 - 2 x^3))
q[x]= ((6 -(0.06666666666666667/x^2) +(1.15/x) + 2 x))/(-0.03333333333333333+ 1.15 x + 2.15 x^2 - 2 x^3)

With Boundary condition 

u[0]==0,u'[0]==0

When I solved it using NDSolve, I got the following errors :

Power::infy: Infinite expression 1/0.^2 encountered. >>
Power::infy: Infinite expression 1/0. encountered. >>
Infinity::indet: Indeterminate expression 6.+ComplexInfinity+ComplexInfinity
encountered. >>
NDSolve::ndnum: Encountered non-numerical value for a derivative at x == 0.`. >>

Could anyone please help me in solving this or similar equation ? Any suggestion in this regard is appreciated. Thanks in advance

I have used the following code:

Clear[u, t, x, eq, sol, v];
L = 6; m = -0.1; n = 1.15; k = 1.15; 
p = (-2 x^3 + (x^2) (1 + k) + n x + m/3); 
eq = -((L + (2 x + (n/x) + (2 m)/(3 x^2)) )/p) u[x] 
-D[p, x]/p D[u[x], x] +   D[D[u[x], x], x]; 
eqn = {eq == 0, u[0] == 0, D[u[x],x]/.x->0 == 0}; 
sol = NDSolve[eqn, u, {x, 0, 1}]
Improve this question
$\endgroup$
6
  • $\begingroup$ Welcome to Mathematica.SE! I suggest the following: 1) As you receive help, try to give it too, by answering questions in your area of expertise. 2) Take the tour! 3) When you see good questions and answers, vote them up by clicking the gray triangles, because the credibility of the system is based on the reputation gained by users sharing their knowledge. Also, please remember to accept the answer, if any, that solves your problem, by clicking the checkmark sign! $\endgroup$
    – Michael E2
    Oct 17, 2015 at 15:39
  • $\begingroup$ You can format inline code and code blocks by selecting the code and clicking the {} button above the edit window. The edit window help button ? is also useful for learning how to format your questions and answers. You may also find this this meta Q&A helpful $\endgroup$
    – Michael E2
    Oct 17, 2015 at 15:39
  • 1
    $\begingroup$ Isn't the solution $u(t) \equiv 0$? $\endgroup$
    – Michael E2
    Oct 17, 2015 at 15:47
  • $\begingroup$ No, actually the solution may come in the range of the order of 10^-3 or so. Here is the code i have used : Clear[u, t, x, eq, sol, v] L = 6; m = -0.1; n = 1.15; k = 1.15; p = (-2 x^3 + (x^2) (1 + k) + n x + m/3); eq = -((L + (2 x + (n/x) + (2 m)/(3 x^2)) )/p) u[x] - D[p, x]/p D[u[x], x] + D[D[u[x], x], x]; eqn = {eq == 0, u[0] == 0, u'[0] == 0}; sol = NDSolve[eqn, u, {x, 0, 1}] $\endgroup$
    – Pisces
    Oct 17, 2015 at 15:50
  • $\begingroup$ Please edit your post with the code, properly formatted, rather than put it in a comment. Comments are by design transitory on stack exchange, so it's better to have all the info in the post itself. $\endgroup$
    – march
    Oct 17, 2015 at 16:02

1 Answer 1

Reset to default
1
$\begingroup$

Using the original code

eq = u''[x] + p[x] u'[x] + q[x] u[x];

p[x] = -(((1.15 + 4.3 x - 6 x^2))/(-0.03333333333333333 + 1.15 x + 
        2.15 x^2 - 2 x^3)) // Rationalize;
q[x] = ((6 - (0.06666666666666667/x^2) + (1.15/x) + 
       2 x))/(-0.03333333333333333 + 1.15 x + 2.15 x^2 - 2 x^3) // Rationalize;

we can use frobeniusNDSolve, which applies Frobenius' method, from my answer to Attempting to use NDSolve to plot harmonic oscillator solutions, to find two independent solutions.

usol = frobeniusNDSolve[eq, u, {x, 0, 1}]

Mathematica graphics

Note that the indicial equation has two conjugate complex roots ${1\over2}\big(1\pm i\sqrt{7}\big)$. We can get two independent real solutions by taking the real and imaginary parts of the Frobenius solution corresponding to one of the roots.

Plot[(Through[{Re, Im}@u[x]] /. First@usol) // Evaluate, {x, 0, 1}]
Plot[(Through[{Re, Im}@u[x]] /. First@usol) // Evaluate, {x, 0, 0.01}]

Mathematica graphics Mathematica graphics

While the value of u[0] is 0 in both cases, the derivative u'[x] approaches infinity (and is roughly asymptotic to 1/Sqrt[x]).

{Through[{Re, Im}@u[0]], Through[{Re, Im}[u'[1*^-12]]]} /. First@usol
(*  {{0., 0.}, {-999941., -1.00006*10^6}}  *)

Consequently, the only solution satisfying the initial condition u[0] == u'[0] == 0 is the trivial solution.

Improve this answer
$\endgroup$
7
  • $\begingroup$ @Pisces You're welcome. $\endgroup$
    – Michael E2
    Oct 21, 2015 at 12:32
  • $\begingroup$ I've checked the above code with the help of frobeniusNDSolve. But it shows error. $\endgroup$
    – Pisces
    Oct 21, 2015 at 15:27
  • $\begingroup$ @Pisces I'll take a look when I get a chance. I updated the code earlier, but maybe I forgot to include something. Sorry about the inconvenience. $\endgroup$
    – Michael E2
    Oct 21, 2015 at 15:47
  • $\begingroup$ Yea, Sure and no inconvenience at all. You have already helped me a lot. Thanks so much for spending your time $\endgroup$
    – Pisces
    Oct 21, 2015 at 15:51
  • $\begingroup$ @Pisces I checked with a fresh kernel and the code worked without trouble. Some possibilities: Maybe you have a lurking definition that interferes with my code? Perhaps you made a change or tried a different example? What version of Mathematica are you using? What errors do you get? $\endgroup$
    – Michael E2
    Oct 22, 2015 at 10:58

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.