1

I've got a mathematica notebook that uses RK4 to successively approximate solutions to a vector valued DE. The output is a nested list similar to:

solutions := {{a[t0],b[t0],c[t0],d[t0]},{a[t0+dt],b[t0+dt],c[t0+dt],d[t0+dt]},...}

The list that I want to use as the x-axis is simply {t0,t0+dt,t0+2dt,...}, and I want to plot it against the a column of my solutions matrix. Perhaps I'm missing it, but I couldn't get ParametricPlot to work with these lists.

Here's a link to my notebook: http://pastebin.com/ZK0pUBSz

Sample data:

t:={0,20,40,60,80,100,120,140,160,180,200}

y:={200.,177.847,162.025,150.607,142.306,136.239,131.788,128.518,126.111,124.341,123.04}

Desired output is a plot with y as the dependent variable and t as the independent variable.


(edit: Solved! I had previously defined solutions as a matrix without any argument, and when I tried to turn it into something that took an argument, the kernel threw up errors. Solved by running a simple Remove@solutions.)

An an extra aside, my end goal is to see the effect of step size on solution resolution. I tried to accomplish this using the following code:

Tlist[step_] := 
Transpose[{Range[0, tmax, step], Part[solutions[step], All, 1]}];
ListLinePlot[{TList[.1], Tlist[.5], Tlist[1]}]

I'm trying to create a plot with each step size using Tlist, but the evaluation fails with write protection. Maybe this isn't the right way to tackle this problem in Mathematica, and I'm taking too much of a procedural approach. Any guidance on this would be welcome!

Here's the code for solutions:

RungeKutta[func_List, yinit_List, y_List, step_] := 
 Module[{k1, k2, k3, k4}, 
  k1 = step N[func /. MapThread[Rule, {y, yinit}]];
  k2 = step N[func /. MapThread[Rule, {y, k1/2 + yinit}]];
  k3 = step N[func /. MapThread[Rule, {y, k2/2 + yinit}]];
  k4 = step N[func /. MapThread[Rule, {y, k3 + yinit}]];
  yinit + Total[{k1, 2 k2, 2 k3, k4}]/6]
solutions[step_] := 
 NestList[RungeKutta[func, #, y, step] &, N[yinit], Round[tmax/step]]

RungeKutta code from here: Solving a system of ODEs with the Runge-Kutta method

  • Welcome! It would be easier to help you if you included some actual data to plot with, but I suspect Thread or MapThread or Transpose may come in handy here. – Yves Klett Oct 06 '15 at 17:05
  • I added a link to the notebook, sorry about that! I'll look at those functions! – ijustlovemath Oct 06 '15 at 17:08
  • No problem! If you post a small sample (or just bogus data with the same structure) here, you will get more answers. This will be quicker to browse, more convenient and secure for users to access, and permanent (your link may vanish). – Yves Klett Oct 06 '15 at 17:12
  • Okay, I'll add some in. – ijustlovemath Oct 06 '15 at 17:12
  • Welcome to Mathematica.SE! I hope you will become a regular contributor. To get started, 1) take the introductory [tour] now, 2) 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, 3) remember to accept the answer, if any, that solves your problem, by clicking the checkmark sign, and 4) give help too, by answering questions in your areas of expertise. – bbgodfrey Oct 06 '15 at 17:17
  • Thanks, bbgodfrey. @YvesKlett, the information you asked for has been included. Although, I'm not exactly sure why the second value is repeated over and over in the yvalues list. – ijustlovemath Oct 06 '15 at 17:20
  • For your sample data and alsoalso the pastebin code: consider using Set (=) instead of using SetDelayed (:=) for your sample lists for instance. Not so important in this case but that will avoid unnecessary evaluations of expressions and can tremensously reduce runtime of code – Lukas Oct 06 '15 at 17:35
  • Thanks for the tip @Lukas. I've added a minor edit to my question about functions in mathematica if you have time to look. – ijustlovemath Oct 06 '15 at 17:50
  • I do not see anything wrong. Of course, you first have to properly define your solutions function so that it actualy is a fubction of step. But the plotting should be okay as far as i can see. Edit: just realized that you tried to run thr code... You have to define the a function solutions[step_] before you use it so that Mathematica knowd what to do – Lukas Oct 06 '15 at 18:33
  • @Lukas, added the solutions[] code that throws the errors. Maybe I'm doing it differently than I should be. – ijustlovemath Oct 06 '15 at 18:35
  • Unfortunately I cannot run the code on a computer right now to help ypu. Using a mobile all the time. You'll have to wait for another one or until the day after tomorrow for me to actually run the code – Lukas Oct 06 '15 at 18:45
  • 1
    @ijustlovemath It would be very nice to give a link to the author of the RungeKutta code. solving-a-system-of-odes-with-the-runge-kutta-method RunnyKine –  Oct 06 '15 at 19:19
  • @Lukas, I figured it out. I had to do a Remove@solutions before it would allow me to pass it a parameter. – ijustlovemath Oct 06 '15 at 22:57

1 Answers1

3

Suppose xColumn is the column of x-values and yColumn the y-values. Then - as pointed to by YvesKlett - you could use ListPlot[Transpose[{xColumn,yColumn}]], for instance.

In your case:

xColumn=Table[t0+i*dt,{i,0,n}];
yColumn=solutions[[All,1]];

Here, n clearly has to be adjusted to your code, depending on how long your time interval is. Then use for instance the ListPlot command from above with option Joined->True or use ListLinePlot instead.

Lukas
  • 2,702
  • 1
  • 13
  • 20
  • I see. So here, we transpose the nested list to make it an appropriate list for listplot? Thanks so much for the clarification – ijustlovemath Oct 06 '15 at 17:28
  • Well, since you want to plot only the first column of your solution matrix you first extract it via solutions[[All,1]]. Then you essentially form the x-y tuples using Transpose so that it is in the required form for the plot commands. Glad I could help – Lukas Oct 06 '15 at 17:31
  • Would it be possible to make that Transpose a function of step size, for plotting multiple plots over varying step sizes? – ijustlovemath Oct 06 '15 at 17:52