Ways to speed up user implemented RK4Speed up Numerical IntegrationSpeed of convergence for NIntegrateTough Calculation, novice mathematica userNumerical integration's speedNumerical integral speedImprove the speed of Gaussian quadrature integrationSolving an unstable BVP numerically, accurately and efficientlyHow to speed up integral of results of PDE modelSolve BVP involving user defined functionUser defined ArcTan function

Is a roofing delivery truck likely to crack my driveway slab?

How was Earth single-handedly capable of creating 3 of the 4 gods of chaos?

What would happen if the UK refused to take part in EU Parliamentary elections?

Have I saved too much for retirement so far?

Ways to speed up user implemented RK4

Applicability of Single Responsibility Principle

How can I get through very long and very dry, but also very useful technical documents when learning a new tool?

Is there an Impartial Brexit Deal comparison site?

Mapping a list into a phase plot

How do I keep an essay about "feeling flat" from feeling flat?

Transcription Beats per minute

Failed to fetch jessie backports repository

What would be the benefits of having both a state and local currencies?

Is there a problem with hiding "forgot password" until it's needed?

Your magic is very sketchy

How to be diplomatic in refusing to write code that breaches the privacy of our users

How do I rename a LINUX host without needing to reboot for the rename to take effect?

Where in the Bible does the greeting ("Dominus Vobiscum") used at Mass come from?

HashMap containsKey() returns false although hashCode() and equals() are true

Can somebody explain Brexit in a few child-proof sentences?

Increase performance creating Mandelbrot set in python

Valid Badminton Score?

Displaying the order of the columns of a table

Lay out the Carpet



Ways to speed up user implemented RK4


Speed up Numerical IntegrationSpeed of convergence for NIntegrateTough Calculation, novice mathematica userNumerical integration's speedNumerical integral speedImprove the speed of Gaussian quadrature integrationSolving an unstable BVP numerically, accurately and efficientlyHow to speed up integral of results of PDE modelSolve BVP involving user defined functionUser defined ArcTan function













3












$begingroup$


So, I've implemented RK4, and I'm wondering what I can do to make it more efficient? What I've got so far is below. I wish to still record all steps. I think AppendTo is doing the most damage to the time, is there a faster alternative?



rk4[f_, variables_, valtinit_, tinit_, tfinal_, nsteps_] := 
Module[table, xlist, ylist, step, k1, k2, k3, k4,
xlist = tinit;
step = N[(tfinal - tinit)/(nsteps)];
ylist = valtinit;
table = xlist, ylist;
Table[
k1 = step* f /. MapThread[Rule, variables, ylist]; (*
Equivalent to step* f/.Thread[Rule[variables,ylist]]*)
k2 = step*f /. MapThread[Rule, variables, k1/2 + ylist];
k3 = step*f /. MapThread[Rule, variables, k2/2 + ylist];
k4 = step*f /. MapThread[Rule, variables, k3 + ylist];
ylist += 1/6 (k1 + 2 (k2 + k3) + k4);
xlist += step;
AppendTo[table, xlist, ylist];
xlist, ylist, nsteps];
table
];


Example Input:



funclist = -x + y, x - y;
initials = 1, 2;
variables = x, y;
init = 0;
final = 200;
nstep = 20000;
approx = rk4[funclist, variables, initials, init, final, nstep]//AbsoluteTiming;



3.59932,...




I'd love some suggestions!










share|improve this question









$endgroup$







  • 1




    $begingroup$
    AppendTo is quadratic time complexity. Might be better to preallocate and set by index. Also it'll be much faster to not use Rule and instead code stuff up a little bit more explicitly. As a general rule, too, use vectorized operators. Those can be very fast. And if everything can be totally functional over "packed arrays" (look them up here) it'll be very quick too.
    $endgroup$
    – b3m2a1
    2 hours ago










  • $begingroup$
    I'll work on implementing it more explicity, this is what came to find first though. It'll require some changes to the inputs, I'll have to ponder this. And preallocating the list is a quick change that won't be an issue to do, I can't believe I forgot that's faster :(. Thanks though!
    $endgroup$
    – Shinaolord
    2 hours ago










  • $begingroup$
    Shinaoloard, using Join[ xlist, ylist, Table[ k1 = step*f /. MapThread[Rule, variables, ylist]; k2 = step*f /. MapThread[Rule, variables, k1/2 + ylist]; k3 = step*f /. MapThread[Rule, variables, k2/2 + ylist]; k4 = step*f /. MapThread[Rule, variables, k3 + ylist]; ylist += 1/6 (k1 + 2 (k2 + k3) + k4); xlist += step; xlist, ylist, nsteps ] ] as return value is already a first step. There is no point in appending if you use a Table anyways.
    $endgroup$
    – Henrik Schumacher
    2 hours ago










  • $begingroup$
    @HenrikSchumacher do you think it would be faster to Pre-allocate a list of length nsteps, and append the values, or to join the values using table? I can obviously change Table to Do to remove the time it takes to make the table list, going by b3m2a1's method, or I could use Join as you have suggested. I'm thinking your method may be faster, though. I've already removed the MapThread part, I am testing the speed increase granted by that at the moment. Just curious which path you think will be faster.
    $endgroup$
    – Shinaolord
    2 hours ago










  • $begingroup$
    I am currently testing the speed difference between the one in the post and rk4t2[f_, valtinit_, tinit_, tfinal_, nsteps_] := Module[test, table, xlist, ylist, step, k1, k2, k3, k4, xlist = tinit; step = N[(tfinal - tinit)/(nsteps)]; ylist = valtinit; table = xlist, ylist; test = Table[ k1 = step* f[ylist] ; k2 = step*f[k1/2 + ylist]; k3 = step*f[k2/2 + ylist]; k4 = step*f[k3 + ylist]; ylist += 1/6 (k1 + 2 (k2 + k3) + k4); xlist += step; xlist, ylist, nsteps]; Join[table, test] ];
    $endgroup$
    – Shinaolord
    2 hours ago















3












$begingroup$


So, I've implemented RK4, and I'm wondering what I can do to make it more efficient? What I've got so far is below. I wish to still record all steps. I think AppendTo is doing the most damage to the time, is there a faster alternative?



rk4[f_, variables_, valtinit_, tinit_, tfinal_, nsteps_] := 
Module[table, xlist, ylist, step, k1, k2, k3, k4,
xlist = tinit;
step = N[(tfinal - tinit)/(nsteps)];
ylist = valtinit;
table = xlist, ylist;
Table[
k1 = step* f /. MapThread[Rule, variables, ylist]; (*
Equivalent to step* f/.Thread[Rule[variables,ylist]]*)
k2 = step*f /. MapThread[Rule, variables, k1/2 + ylist];
k3 = step*f /. MapThread[Rule, variables, k2/2 + ylist];
k4 = step*f /. MapThread[Rule, variables, k3 + ylist];
ylist += 1/6 (k1 + 2 (k2 + k3) + k4);
xlist += step;
AppendTo[table, xlist, ylist];
xlist, ylist, nsteps];
table
];


Example Input:



funclist = -x + y, x - y;
initials = 1, 2;
variables = x, y;
init = 0;
final = 200;
nstep = 20000;
approx = rk4[funclist, variables, initials, init, final, nstep]//AbsoluteTiming;



3.59932,...




I'd love some suggestions!










share|improve this question









$endgroup$







  • 1




    $begingroup$
    AppendTo is quadratic time complexity. Might be better to preallocate and set by index. Also it'll be much faster to not use Rule and instead code stuff up a little bit more explicitly. As a general rule, too, use vectorized operators. Those can be very fast. And if everything can be totally functional over "packed arrays" (look them up here) it'll be very quick too.
    $endgroup$
    – b3m2a1
    2 hours ago










  • $begingroup$
    I'll work on implementing it more explicity, this is what came to find first though. It'll require some changes to the inputs, I'll have to ponder this. And preallocating the list is a quick change that won't be an issue to do, I can't believe I forgot that's faster :(. Thanks though!
    $endgroup$
    – Shinaolord
    2 hours ago










  • $begingroup$
    Shinaoloard, using Join[ xlist, ylist, Table[ k1 = step*f /. MapThread[Rule, variables, ylist]; k2 = step*f /. MapThread[Rule, variables, k1/2 + ylist]; k3 = step*f /. MapThread[Rule, variables, k2/2 + ylist]; k4 = step*f /. MapThread[Rule, variables, k3 + ylist]; ylist += 1/6 (k1 + 2 (k2 + k3) + k4); xlist += step; xlist, ylist, nsteps ] ] as return value is already a first step. There is no point in appending if you use a Table anyways.
    $endgroup$
    – Henrik Schumacher
    2 hours ago










  • $begingroup$
    @HenrikSchumacher do you think it would be faster to Pre-allocate a list of length nsteps, and append the values, or to join the values using table? I can obviously change Table to Do to remove the time it takes to make the table list, going by b3m2a1's method, or I could use Join as you have suggested. I'm thinking your method may be faster, though. I've already removed the MapThread part, I am testing the speed increase granted by that at the moment. Just curious which path you think will be faster.
    $endgroup$
    – Shinaolord
    2 hours ago










  • $begingroup$
    I am currently testing the speed difference between the one in the post and rk4t2[f_, valtinit_, tinit_, tfinal_, nsteps_] := Module[test, table, xlist, ylist, step, k1, k2, k3, k4, xlist = tinit; step = N[(tfinal - tinit)/(nsteps)]; ylist = valtinit; table = xlist, ylist; test = Table[ k1 = step* f[ylist] ; k2 = step*f[k1/2 + ylist]; k3 = step*f[k2/2 + ylist]; k4 = step*f[k3 + ylist]; ylist += 1/6 (k1 + 2 (k2 + k3) + k4); xlist += step; xlist, ylist, nsteps]; Join[table, test] ];
    $endgroup$
    – Shinaolord
    2 hours ago













3












3








3


1



$begingroup$


So, I've implemented RK4, and I'm wondering what I can do to make it more efficient? What I've got so far is below. I wish to still record all steps. I think AppendTo is doing the most damage to the time, is there a faster alternative?



rk4[f_, variables_, valtinit_, tinit_, tfinal_, nsteps_] := 
Module[table, xlist, ylist, step, k1, k2, k3, k4,
xlist = tinit;
step = N[(tfinal - tinit)/(nsteps)];
ylist = valtinit;
table = xlist, ylist;
Table[
k1 = step* f /. MapThread[Rule, variables, ylist]; (*
Equivalent to step* f/.Thread[Rule[variables,ylist]]*)
k2 = step*f /. MapThread[Rule, variables, k1/2 + ylist];
k3 = step*f /. MapThread[Rule, variables, k2/2 + ylist];
k4 = step*f /. MapThread[Rule, variables, k3 + ylist];
ylist += 1/6 (k1 + 2 (k2 + k3) + k4);
xlist += step;
AppendTo[table, xlist, ylist];
xlist, ylist, nsteps];
table
];


Example Input:



funclist = -x + y, x - y;
initials = 1, 2;
variables = x, y;
init = 0;
final = 200;
nstep = 20000;
approx = rk4[funclist, variables, initials, init, final, nstep]//AbsoluteTiming;



3.59932,...




I'd love some suggestions!










share|improve this question









$endgroup$




So, I've implemented RK4, and I'm wondering what I can do to make it more efficient? What I've got so far is below. I wish to still record all steps. I think AppendTo is doing the most damage to the time, is there a faster alternative?



rk4[f_, variables_, valtinit_, tinit_, tfinal_, nsteps_] := 
Module[table, xlist, ylist, step, k1, k2, k3, k4,
xlist = tinit;
step = N[(tfinal - tinit)/(nsteps)];
ylist = valtinit;
table = xlist, ylist;
Table[
k1 = step* f /. MapThread[Rule, variables, ylist]; (*
Equivalent to step* f/.Thread[Rule[variables,ylist]]*)
k2 = step*f /. MapThread[Rule, variables, k1/2 + ylist];
k3 = step*f /. MapThread[Rule, variables, k2/2 + ylist];
k4 = step*f /. MapThread[Rule, variables, k3 + ylist];
ylist += 1/6 (k1 + 2 (k2 + k3) + k4);
xlist += step;
AppendTo[table, xlist, ylist];
xlist, ylist, nsteps];
table
];


Example Input:



funclist = -x + y, x - y;
initials = 1, 2;
variables = x, y;
init = 0;
final = 200;
nstep = 20000;
approx = rk4[funclist, variables, initials, init, final, nstep]//AbsoluteTiming;



3.59932,...




I'd love some suggestions!







differential-equations numerical-integration






share|improve this question













share|improve this question











share|improve this question




share|improve this question










asked 2 hours ago









ShinaolordShinaolord

808




808







  • 1




    $begingroup$
    AppendTo is quadratic time complexity. Might be better to preallocate and set by index. Also it'll be much faster to not use Rule and instead code stuff up a little bit more explicitly. As a general rule, too, use vectorized operators. Those can be very fast. And if everything can be totally functional over "packed arrays" (look them up here) it'll be very quick too.
    $endgroup$
    – b3m2a1
    2 hours ago










  • $begingroup$
    I'll work on implementing it more explicity, this is what came to find first though. It'll require some changes to the inputs, I'll have to ponder this. And preallocating the list is a quick change that won't be an issue to do, I can't believe I forgot that's faster :(. Thanks though!
    $endgroup$
    – Shinaolord
    2 hours ago










  • $begingroup$
    Shinaoloard, using Join[ xlist, ylist, Table[ k1 = step*f /. MapThread[Rule, variables, ylist]; k2 = step*f /. MapThread[Rule, variables, k1/2 + ylist]; k3 = step*f /. MapThread[Rule, variables, k2/2 + ylist]; k4 = step*f /. MapThread[Rule, variables, k3 + ylist]; ylist += 1/6 (k1 + 2 (k2 + k3) + k4); xlist += step; xlist, ylist, nsteps ] ] as return value is already a first step. There is no point in appending if you use a Table anyways.
    $endgroup$
    – Henrik Schumacher
    2 hours ago










  • $begingroup$
    @HenrikSchumacher do you think it would be faster to Pre-allocate a list of length nsteps, and append the values, or to join the values using table? I can obviously change Table to Do to remove the time it takes to make the table list, going by b3m2a1's method, or I could use Join as you have suggested. I'm thinking your method may be faster, though. I've already removed the MapThread part, I am testing the speed increase granted by that at the moment. Just curious which path you think will be faster.
    $endgroup$
    – Shinaolord
    2 hours ago










  • $begingroup$
    I am currently testing the speed difference between the one in the post and rk4t2[f_, valtinit_, tinit_, tfinal_, nsteps_] := Module[test, table, xlist, ylist, step, k1, k2, k3, k4, xlist = tinit; step = N[(tfinal - tinit)/(nsteps)]; ylist = valtinit; table = xlist, ylist; test = Table[ k1 = step* f[ylist] ; k2 = step*f[k1/2 + ylist]; k3 = step*f[k2/2 + ylist]; k4 = step*f[k3 + ylist]; ylist += 1/6 (k1 + 2 (k2 + k3) + k4); xlist += step; xlist, ylist, nsteps]; Join[table, test] ];
    $endgroup$
    – Shinaolord
    2 hours ago












  • 1




    $begingroup$
    AppendTo is quadratic time complexity. Might be better to preallocate and set by index. Also it'll be much faster to not use Rule and instead code stuff up a little bit more explicitly. As a general rule, too, use vectorized operators. Those can be very fast. And if everything can be totally functional over "packed arrays" (look them up here) it'll be very quick too.
    $endgroup$
    – b3m2a1
    2 hours ago










  • $begingroup$
    I'll work on implementing it more explicity, this is what came to find first though. It'll require some changes to the inputs, I'll have to ponder this. And preallocating the list is a quick change that won't be an issue to do, I can't believe I forgot that's faster :(. Thanks though!
    $endgroup$
    – Shinaolord
    2 hours ago










  • $begingroup$
    Shinaoloard, using Join[ xlist, ylist, Table[ k1 = step*f /. MapThread[Rule, variables, ylist]; k2 = step*f /. MapThread[Rule, variables, k1/2 + ylist]; k3 = step*f /. MapThread[Rule, variables, k2/2 + ylist]; k4 = step*f /. MapThread[Rule, variables, k3 + ylist]; ylist += 1/6 (k1 + 2 (k2 + k3) + k4); xlist += step; xlist, ylist, nsteps ] ] as return value is already a first step. There is no point in appending if you use a Table anyways.
    $endgroup$
    – Henrik Schumacher
    2 hours ago










  • $begingroup$
    @HenrikSchumacher do you think it would be faster to Pre-allocate a list of length nsteps, and append the values, or to join the values using table? I can obviously change Table to Do to remove the time it takes to make the table list, going by b3m2a1's method, or I could use Join as you have suggested. I'm thinking your method may be faster, though. I've already removed the MapThread part, I am testing the speed increase granted by that at the moment. Just curious which path you think will be faster.
    $endgroup$
    – Shinaolord
    2 hours ago










  • $begingroup$
    I am currently testing the speed difference between the one in the post and rk4t2[f_, valtinit_, tinit_, tfinal_, nsteps_] := Module[test, table, xlist, ylist, step, k1, k2, k3, k4, xlist = tinit; step = N[(tfinal - tinit)/(nsteps)]; ylist = valtinit; table = xlist, ylist; test = Table[ k1 = step* f[ylist] ; k2 = step*f[k1/2 + ylist]; k3 = step*f[k2/2 + ylist]; k4 = step*f[k3 + ylist]; ylist += 1/6 (k1 + 2 (k2 + k3) + k4); xlist += step; xlist, ylist, nsteps]; Join[table, test] ];
    $endgroup$
    – Shinaolord
    2 hours ago







1




1




$begingroup$
AppendTo is quadratic time complexity. Might be better to preallocate and set by index. Also it'll be much faster to not use Rule and instead code stuff up a little bit more explicitly. As a general rule, too, use vectorized operators. Those can be very fast. And if everything can be totally functional over "packed arrays" (look them up here) it'll be very quick too.
$endgroup$
– b3m2a1
2 hours ago




$begingroup$
AppendTo is quadratic time complexity. Might be better to preallocate and set by index. Also it'll be much faster to not use Rule and instead code stuff up a little bit more explicitly. As a general rule, too, use vectorized operators. Those can be very fast. And if everything can be totally functional over "packed arrays" (look them up here) it'll be very quick too.
$endgroup$
– b3m2a1
2 hours ago












$begingroup$
I'll work on implementing it more explicity, this is what came to find first though. It'll require some changes to the inputs, I'll have to ponder this. And preallocating the list is a quick change that won't be an issue to do, I can't believe I forgot that's faster :(. Thanks though!
$endgroup$
– Shinaolord
2 hours ago




$begingroup$
I'll work on implementing it more explicity, this is what came to find first though. It'll require some changes to the inputs, I'll have to ponder this. And preallocating the list is a quick change that won't be an issue to do, I can't believe I forgot that's faster :(. Thanks though!
$endgroup$
– Shinaolord
2 hours ago












$begingroup$
Shinaoloard, using Join[ xlist, ylist, Table[ k1 = step*f /. MapThread[Rule, variables, ylist]; k2 = step*f /. MapThread[Rule, variables, k1/2 + ylist]; k3 = step*f /. MapThread[Rule, variables, k2/2 + ylist]; k4 = step*f /. MapThread[Rule, variables, k3 + ylist]; ylist += 1/6 (k1 + 2 (k2 + k3) + k4); xlist += step; xlist, ylist, nsteps ] ] as return value is already a first step. There is no point in appending if you use a Table anyways.
$endgroup$
– Henrik Schumacher
2 hours ago




$begingroup$
Shinaoloard, using Join[ xlist, ylist, Table[ k1 = step*f /. MapThread[Rule, variables, ylist]; k2 = step*f /. MapThread[Rule, variables, k1/2 + ylist]; k3 = step*f /. MapThread[Rule, variables, k2/2 + ylist]; k4 = step*f /. MapThread[Rule, variables, k3 + ylist]; ylist += 1/6 (k1 + 2 (k2 + k3) + k4); xlist += step; xlist, ylist, nsteps ] ] as return value is already a first step. There is no point in appending if you use a Table anyways.
$endgroup$
– Henrik Schumacher
2 hours ago












$begingroup$
@HenrikSchumacher do you think it would be faster to Pre-allocate a list of length nsteps, and append the values, or to join the values using table? I can obviously change Table to Do to remove the time it takes to make the table list, going by b3m2a1's method, or I could use Join as you have suggested. I'm thinking your method may be faster, though. I've already removed the MapThread part, I am testing the speed increase granted by that at the moment. Just curious which path you think will be faster.
$endgroup$
– Shinaolord
2 hours ago




$begingroup$
@HenrikSchumacher do you think it would be faster to Pre-allocate a list of length nsteps, and append the values, or to join the values using table? I can obviously change Table to Do to remove the time it takes to make the table list, going by b3m2a1's method, or I could use Join as you have suggested. I'm thinking your method may be faster, though. I've already removed the MapThread part, I am testing the speed increase granted by that at the moment. Just curious which path you think will be faster.
$endgroup$
– Shinaolord
2 hours ago












$begingroup$
I am currently testing the speed difference between the one in the post and rk4t2[f_, valtinit_, tinit_, tfinal_, nsteps_] := Module[test, table, xlist, ylist, step, k1, k2, k3, k4, xlist = tinit; step = N[(tfinal - tinit)/(nsteps)]; ylist = valtinit; table = xlist, ylist; test = Table[ k1 = step* f[ylist] ; k2 = step*f[k1/2 + ylist]; k3 = step*f[k2/2 + ylist]; k4 = step*f[k3 + ylist]; ylist += 1/6 (k1 + 2 (k2 + k3) + k4); xlist += step; xlist, ylist, nsteps]; Join[table, test] ];
$endgroup$
– Shinaolord
2 hours ago




$begingroup$
I am currently testing the speed difference between the one in the post and rk4t2[f_, valtinit_, tinit_, tfinal_, nsteps_] := Module[test, table, xlist, ylist, step, k1, k2, k3, k4, xlist = tinit; step = N[(tfinal - tinit)/(nsteps)]; ylist = valtinit; table = xlist, ylist; test = Table[ k1 = step* f[ylist] ; k2 = step*f[k1/2 + ylist]; k3 = step*f[k2/2 + ylist]; k4 = step*f[k3 + ylist]; ylist += 1/6 (k1 + 2 (k2 + k3) + k4); xlist += step; xlist, ylist, nsteps]; Join[table, test] ];
$endgroup$
– Shinaolord
2 hours ago










1 Answer
1






active

oldest

votes


















4












$begingroup$

Just to give you an impression how fast things may get when you use the right tools.



For given stepsize τ and given vectorfield F, this creates a CompiledFunction cStep that computes a single Runge-Kutta step



F = X [Function] -Indexed[X, 2], Indexed[X, 1];

τ = 0.01;
Block[YY, Y, k1, k2, k3, k4,

YY = Table[Compile`GetElement[Y, i], i, 1, 2];
k1 = τ F[YY];
k2 = τ F[0.5 k1 + YY];
k3 = τ F[0.5 k2 + YY];
k4 = τ F[k3 + YY];

cStep = With[code = YY + (k1 + 2. (k2 + k3) + k4)/6. ,
Compile[Y, _Real, 1,
code,
CompilationTarget -> "C",
RuntimeOptions -> "Speed"
]
]
];


Now we can apply it 2 million times with NestList and still need only 2 seconds.



nsteps = 20000000;
xlist = Range[0., step nsteps, step];
Ylist = NestList[cStep, initials, nsteps]; // AbsoluteTiming // First



2.08678







share|improve this answer









$endgroup$












  • $begingroup$
    Damn, you definitely know how to use Mathematica A LOT more efficiently than I do. Thanks!
    $endgroup$
    – Shinaolord
    2 hours ago










  • $begingroup$
    You're welcome.
    $endgroup$
    – Henrik Schumacher
    2 hours ago










  • $begingroup$
    I'll have to play around with Compile, it definitely seems like a massive speed up if used correctly.
    $endgroup$
    – Shinaolord
    2 hours ago










Your Answer





StackExchange.ifUsing("editor", function ()
return StackExchange.using("mathjaxEditing", function ()
StackExchange.MarkdownEditor.creationCallbacks.add(function (editor, postfix)
StackExchange.mathjaxEditing.prepareWmdForMathJax(editor, postfix, [["$", "$"], ["\\(","\\)"]]);
);
);
, "mathjax-editing");

StackExchange.ready(function()
var channelOptions =
tags: "".split(" "),
id: "387"
;
initTagRenderer("".split(" "), "".split(" "), channelOptions);

StackExchange.using("externalEditor", function()
// Have to fire editor after snippets, if snippets enabled
if (StackExchange.settings.snippets.snippetsEnabled)
StackExchange.using("snippets", function()
createEditor();
);

else
createEditor();

);

function createEditor()
StackExchange.prepareEditor(
heartbeatType: 'answer',
autoActivateHeartbeat: false,
convertImagesToLinks: false,
noModals: true,
showLowRepImageUploadWarning: true,
reputationToPostImages: null,
bindNavPrevention: true,
postfix: "",
imageUploader:
brandingHtml: "Powered by u003ca class="icon-imgur-white" href="https://imgur.com/"u003eu003c/au003e",
contentPolicyHtml: "User contributions licensed under u003ca href="https://creativecommons.org/licenses/by-sa/3.0/"u003ecc by-sa 3.0 with attribution requiredu003c/au003e u003ca href="https://stackoverflow.com/legal/content-policy"u003e(content policy)u003c/au003e",
allowUrls: true
,
onDemand: true,
discardSelector: ".discard-answer"
,immediatelyShowMarkdownHelp:true
);



);













draft saved

draft discarded


















StackExchange.ready(
function ()
StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fmathematica.stackexchange.com%2fquestions%2f194002%2fways-to-speed-up-user-implemented-rk4%23new-answer', 'question_page');

);

Post as a guest















Required, but never shown

























1 Answer
1






active

oldest

votes








1 Answer
1






active

oldest

votes









active

oldest

votes






active

oldest

votes









4












$begingroup$

Just to give you an impression how fast things may get when you use the right tools.



For given stepsize τ and given vectorfield F, this creates a CompiledFunction cStep that computes a single Runge-Kutta step



F = X [Function] -Indexed[X, 2], Indexed[X, 1];

τ = 0.01;
Block[YY, Y, k1, k2, k3, k4,

YY = Table[Compile`GetElement[Y, i], i, 1, 2];
k1 = τ F[YY];
k2 = τ F[0.5 k1 + YY];
k3 = τ F[0.5 k2 + YY];
k4 = τ F[k3 + YY];

cStep = With[code = YY + (k1 + 2. (k2 + k3) + k4)/6. ,
Compile[Y, _Real, 1,
code,
CompilationTarget -> "C",
RuntimeOptions -> "Speed"
]
]
];


Now we can apply it 2 million times with NestList and still need only 2 seconds.



nsteps = 20000000;
xlist = Range[0., step nsteps, step];
Ylist = NestList[cStep, initials, nsteps]; // AbsoluteTiming // First



2.08678







share|improve this answer









$endgroup$












  • $begingroup$
    Damn, you definitely know how to use Mathematica A LOT more efficiently than I do. Thanks!
    $endgroup$
    – Shinaolord
    2 hours ago










  • $begingroup$
    You're welcome.
    $endgroup$
    – Henrik Schumacher
    2 hours ago










  • $begingroup$
    I'll have to play around with Compile, it definitely seems like a massive speed up if used correctly.
    $endgroup$
    – Shinaolord
    2 hours ago















4












$begingroup$

Just to give you an impression how fast things may get when you use the right tools.



For given stepsize τ and given vectorfield F, this creates a CompiledFunction cStep that computes a single Runge-Kutta step



F = X [Function] -Indexed[X, 2], Indexed[X, 1];

τ = 0.01;
Block[YY, Y, k1, k2, k3, k4,

YY = Table[Compile`GetElement[Y, i], i, 1, 2];
k1 = τ F[YY];
k2 = τ F[0.5 k1 + YY];
k3 = τ F[0.5 k2 + YY];
k4 = τ F[k3 + YY];

cStep = With[code = YY + (k1 + 2. (k2 + k3) + k4)/6. ,
Compile[Y, _Real, 1,
code,
CompilationTarget -> "C",
RuntimeOptions -> "Speed"
]
]
];


Now we can apply it 2 million times with NestList and still need only 2 seconds.



nsteps = 20000000;
xlist = Range[0., step nsteps, step];
Ylist = NestList[cStep, initials, nsteps]; // AbsoluteTiming // First



2.08678







share|improve this answer









$endgroup$












  • $begingroup$
    Damn, you definitely know how to use Mathematica A LOT more efficiently than I do. Thanks!
    $endgroup$
    – Shinaolord
    2 hours ago










  • $begingroup$
    You're welcome.
    $endgroup$
    – Henrik Schumacher
    2 hours ago










  • $begingroup$
    I'll have to play around with Compile, it definitely seems like a massive speed up if used correctly.
    $endgroup$
    – Shinaolord
    2 hours ago













4












4








4





$begingroup$

Just to give you an impression how fast things may get when you use the right tools.



For given stepsize τ and given vectorfield F, this creates a CompiledFunction cStep that computes a single Runge-Kutta step



F = X [Function] -Indexed[X, 2], Indexed[X, 1];

τ = 0.01;
Block[YY, Y, k1, k2, k3, k4,

YY = Table[Compile`GetElement[Y, i], i, 1, 2];
k1 = τ F[YY];
k2 = τ F[0.5 k1 + YY];
k3 = τ F[0.5 k2 + YY];
k4 = τ F[k3 + YY];

cStep = With[code = YY + (k1 + 2. (k2 + k3) + k4)/6. ,
Compile[Y, _Real, 1,
code,
CompilationTarget -> "C",
RuntimeOptions -> "Speed"
]
]
];


Now we can apply it 2 million times with NestList and still need only 2 seconds.



nsteps = 20000000;
xlist = Range[0., step nsteps, step];
Ylist = NestList[cStep, initials, nsteps]; // AbsoluteTiming // First



2.08678







share|improve this answer









$endgroup$



Just to give you an impression how fast things may get when you use the right tools.



For given stepsize τ and given vectorfield F, this creates a CompiledFunction cStep that computes a single Runge-Kutta step



F = X [Function] -Indexed[X, 2], Indexed[X, 1];

τ = 0.01;
Block[YY, Y, k1, k2, k3, k4,

YY = Table[Compile`GetElement[Y, i], i, 1, 2];
k1 = τ F[YY];
k2 = τ F[0.5 k1 + YY];
k3 = τ F[0.5 k2 + YY];
k4 = τ F[k3 + YY];

cStep = With[code = YY + (k1 + 2. (k2 + k3) + k4)/6. ,
Compile[Y, _Real, 1,
code,
CompilationTarget -> "C",
RuntimeOptions -> "Speed"
]
]
];


Now we can apply it 2 million times with NestList and still need only 2 seconds.



nsteps = 20000000;
xlist = Range[0., step nsteps, step];
Ylist = NestList[cStep, initials, nsteps]; // AbsoluteTiming // First



2.08678








share|improve this answer












share|improve this answer



share|improve this answer










answered 2 hours ago









Henrik SchumacherHenrik Schumacher

57.9k579159




57.9k579159











  • $begingroup$
    Damn, you definitely know how to use Mathematica A LOT more efficiently than I do. Thanks!
    $endgroup$
    – Shinaolord
    2 hours ago










  • $begingroup$
    You're welcome.
    $endgroup$
    – Henrik Schumacher
    2 hours ago










  • $begingroup$
    I'll have to play around with Compile, it definitely seems like a massive speed up if used correctly.
    $endgroup$
    – Shinaolord
    2 hours ago
















  • $begingroup$
    Damn, you definitely know how to use Mathematica A LOT more efficiently than I do. Thanks!
    $endgroup$
    – Shinaolord
    2 hours ago










  • $begingroup$
    You're welcome.
    $endgroup$
    – Henrik Schumacher
    2 hours ago










  • $begingroup$
    I'll have to play around with Compile, it definitely seems like a massive speed up if used correctly.
    $endgroup$
    – Shinaolord
    2 hours ago















$begingroup$
Damn, you definitely know how to use Mathematica A LOT more efficiently than I do. Thanks!
$endgroup$
– Shinaolord
2 hours ago




$begingroup$
Damn, you definitely know how to use Mathematica A LOT more efficiently than I do. Thanks!
$endgroup$
– Shinaolord
2 hours ago












$begingroup$
You're welcome.
$endgroup$
– Henrik Schumacher
2 hours ago




$begingroup$
You're welcome.
$endgroup$
– Henrik Schumacher
2 hours ago












$begingroup$
I'll have to play around with Compile, it definitely seems like a massive speed up if used correctly.
$endgroup$
– Shinaolord
2 hours ago




$begingroup$
I'll have to play around with Compile, it definitely seems like a massive speed up if used correctly.
$endgroup$
– Shinaolord
2 hours ago

















draft saved

draft discarded
















































Thanks for contributing an answer to Mathematica Stack Exchange!


  • Please be sure to answer the question. Provide details and share your research!

But avoid


  • Asking for help, clarification, or responding to other answers.

  • Making statements based on opinion; back them up with references or personal experience.

Use MathJax to format equations. MathJax reference.


To learn more, see our tips on writing great answers.




draft saved


draft discarded














StackExchange.ready(
function ()
StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fmathematica.stackexchange.com%2fquestions%2f194002%2fways-to-speed-up-user-implemented-rk4%23new-answer', 'question_page');

);

Post as a guest















Required, but never shown





















































Required, but never shown














Required, but never shown












Required, but never shown







Required, but never shown

































Required, but never shown














Required, but never shown












Required, but never shown







Required, but never shown







Popular posts from this blog

Can not update quote_id field of “quote_item” table magento 2Magento 2.1 - We can't remove the item. (Shopping Cart doesnt allow us to remove items before becomes empty)Add value for custom quote item attribute using REST apiREST API endpoint v1/carts/cartId/items always returns error messageCorrect way to save entries to databaseHow to remove all associated quote objects of a customer completelyMagento 2 - Save value from custom input field to quote_itemGet quote_item data using quote id and product id filter in Magento 2How to set additional data to quote_item table from controller in Magento 2?What is the purpose of additional_data column in quote_item table in magento2Set Custom Price to Quote item magento2 from controller

How to solve knockout JS error in Magento 2 Planned maintenance scheduled April 23, 2019 at 23:30 UTC (7:30pm US/Eastern) Announcing the arrival of Valued Associate #679: Cesar Manara Unicorn Meta Zoo #1: Why another podcast?(Magento2) knockout.js:3012 Uncaught ReferenceError: Unable to process bindingUnable to process binding Knockout.js magento 2Cannot read property `scopeLabel` of undefined on Product Detail PageCan't get Customer Data on frontend in Magento 2Magento2 Order Summary - unable to process bindingKO templates are not loading in Magento 2.1 applicationgetting knockout js error magento 2Product grid not load -— Unable to process binding Knockout.js magento 2Product form not loaded in magento2Uncaught ReferenceError: Unable to process binding “if: function()return (isShowLegend()) ” magento 2

Nissan Patrol Зміст Перше покоління — 4W60 (1951-1960) | Друге покоління — 60 series (1960-1980) | Третє покоління (1980–2002) | Четверте покоління — Y60 (1987–1998) | П'яте покоління — Y61 (1997–2013) | Шосте покоління — Y62 (2010- ) | Посилання | Зноски | Навігаційне менюОфіційний український сайтТест-драйв Nissan Patrol 2010 7-го поколінняNissan PatrolКак мы тестировали Nissan Patrol 2016рвиправивши або дописавши її