Mathematica Stack Exchange is a question and answer site for users of Mathematica. Join them; it only takes a minute:

Sign up
Here's how it works:
  1. Anybody can ask a question
  2. Anybody can answer
  3. The best answers are voted up and rise to the top

I am very new to Mathematica and have a question regarding plotting a ODE boundary value problem.

Let say I have an ODE

alpha = 1/2;
beta = -5;
eq = beta*f'''[x] + alpha*f[x]*f'[x] == 0
bc = {f[0] == 0, f'[0] == 0, f'[10] == 1}; 
NDSolve[{eq, bc}, f[x], {x, 0, 10}]

I want to plot alpha*f''[0] + (f''[0])^3 vs beta. In the plot, beta should take the range [-5, 5]. NDSolve should use shooting method for the numerical solution.

Also, how I can extract data from that plot as a .txt file.

share|improve this question
Could you try to write down some Mathematica code? – Dr. belisarius Jul 22 '13 at 4:27
No idea how to write mathematica code here. Sorry mate – MMM Jul 22 '13 at 4:34
you can take infinity to be some thing like 10 or 20 – MMM Jul 22 '13 at 4:40
Thanks mate. Now it looks very clear. Cheers – MMM Jul 22 '13 at 4:46
up vote 2 down vote accepted

This is not a complete solution as some issues remain. but may be will help you get started, as I have to go.

There is some issues when beta is takes some values (between -1 and 1). NDSolve complains for some values of beta.

Failed to converge to the requested accuracy or precision \
within 100 iterations.

This generates the plot you wanted, but you might have to play with NDSolve options or make sure your parameters are ok.

y[beta_] := 
  Module[{eq, sol, x, f, bc, max = 10, alpha = 0.5, r}, 
   eq = beta*f'''[x] + alpha*f[x]*f'[x] == 0;
   bc = {f[0] == 0, f'[0] == 0, f'[max] == 1};
   sol = First@NDSolve[{eq, bc}, {f[x], f'[x], f''[x]}, {x, 0, max}];
   r = f''[x] /. sol /. x -> 0;
   alpha*r + beta*r^3];
data = Table[{beta, y[beta]}, {beta, -5, 5, 1}];
ListPlot[data, Joined -> True, Mesh -> True, Frame -> True, 
 FrameLabel -> {{"y(beta)", None}, {beta, "y[beta] vs. beta"}}]

enter image description here

share|improve this answer
data = Table[{lambda, y[lambda]}, {lambda, -5, 0, 1}]; – MMM Jul 22 '13 at 5:48
Trying to plot alphar + betar^3 vs lambda here. Can we use shooting method? – MMM Jul 22 '13 at 5:49
I wanted the answer you already give me but I am unable to execute the code. STill I got an error. Could you please upload ur .nb file? As for as the change is concern, I was hope to give an ode with which NDsolve may have no issue. – MMM Jul 22 '13 at 6:03
let us continue this discussion in chat – MMM Jul 22 '13 at 6:09

Your Answer


By posting your answer, you agree to the privacy policy and terms of service.

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