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
ThreadorMapThreadorTransposemay come in handy here. – Yves Klett Oct 06 '15 at 17:05Set(=) instead of usingSetDelayed(:=) 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:35step. 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 functionsolutions[step_]before you use it so that Mathematica knowd what to do – Lukas Oct 06 '15 at 18:33Remove@solutionsbefore it would allow me to pass it a parameter. – ijustlovemath Oct 06 '15 at 22:57