How to use NDSolve with limited RAM?

Question

How can I use NDSolve when running out of RAM?

The code I am implementing:

sol = NDSolve[ pdes, {A, B, C},
    {x, -L, L}, {t, Tf, Tf} ]

in order to obtain the end solution of NDSolve to be exported to the next run (due to limited RAM).

In previous attempts the solution of NDSolve was exported as an interpolating function, such as:

Export["sol_Tf.wdx",sol]

sol_Tf.wdx is the export file of the NDSolve run.

Also see my previous question.

I then imported it as:

Import["sol_Tf.wdx"]

and used the boundary conditions:

A[x, Ti] == A[x, Tf]
B[x, Ti] == B[x, Tf] 
C[x, Ti] == C[x, Tf]

But the results were different from running it all in one go. What went wrong here? Is there a better way of doing this?

For example, in the heat equation:

sol = NDSolve[{D[u[t, x], t] == D[u[t, x], x, x], u[0, x] == 0,
u[t, 0] == Sin[t], u[t, 5] == 0}, u, {t, 0, 10}, {x, 0, 5}]

I would want to continue it with the end results of the previous run, as the new boundary conditions of the next run; RAM free, and by running it again from t = 10 to t = 20.

I think it would be useful if you added a simple PDE exmaple such that people can experiment with this. As the question stands now, I'd first need to come up with an example to try this. — user21, Sep 16 '19 at 05:56

Edmund · Answer 1 · 2021-07-29T13:47:14.707

You may use the NDSolve Components and Data Structures tutorial to control the memory usage and save down the state of intermediate runs.

Initialise the NDSolve`StateData for the complete range you need to create solutions. Below is done for 0 <= t <= 30 and I will iterate in chunks of 10.

ndsStateData = First@NDSolve`ProcessEquations[
    {
     D[u[t, x], t] == D[u[t, x], x, x],
     u[0, x] == 0,
     u[t, 0] == Sin[t],
     u[t, 5] == 0
     },
    u,
    {t, 0, 30}, {x, 0, 5}
    ];

Next I iterate in 3 chunks of 10. ndsStateData can be saved down after the current chunk's solution is extracted with NDSolve`ProcessSolutions. Below I reassign ndsStateData with the reinitialised NDSolve`StateData instead.

sols = {};
Module[{step = #},
    NDSolve`Iterate[ndsStateData, {(step - 1)*10, step 10}];
    AppendTo[sols, u /. NDSolve`ProcessSolutions[ndsStateData]];
    ndsStateData =
     First@NDSolve`Reinitialize[ndsStateData, {u[step 10, x] == sols[[step]][step 10, x]}];
    ] & /@ Range@3;

sols contains the 3 solutions over each chunk. Notice the difference in the domains.

sols

Plotting the solutions.

Plot3D[sols[[#]][t, x], {t, (# - 1) 10, # 10}, {x, 0, 5},
    PlotRange -> Full,
    PlotStyle -> ColorData[109][#]
    ] & /@ Range@3 // Show[#, PlotRange -> All] &

Hope this helps.

How did you make a gif of you rotating the figure? – Derek H Jul 29 '21 at 14:10 — Derek H, Jul 29 '21 at 14:10

Tim Laska · Answer 2 · 2019-09-18T02:29:22.860

The following throws some inconsistency warnings, but appears to work (i.e., no obvious discontinuities at the seam). Also, I did not change the boundary condition, but I fed the final condition of the first solution to the "initial" condition of the second solution.

uif1 = NDSolveValue[{D[u[t, x], t] == D[u[t, x], x, x], u[0, x] == 0, 
    u[t, 0] == Sin[t], u[t, 5] == 0}, u, {t, 0, 10}, {x, 0, 5}];
uif2 = NDSolveValue[{D[u[t, x], t] == D[u[t, x], x, x], 
    u[10, x] == uif1[10, x], u[t, 0] == Sin[t], u[t, 5] == 0}, 
   u, {t, 10, 20}, {x, 0, 5}];
plt1 = Plot3D[uif1[t, x], {t, 0., 10.}, {x, 0.`, 5.`}, 
   PlotRange -> {{0, 20}, {0, 5}, {-1, 1}}];
plt2 = Plot3D[uif2[t, x], {t, 10.`, 20.`}, {x, 0.`, 5.`}, 
   PlotRange -> {{0, 20}, {0, 5}, {-1, 1}}, 
   PlotStyle -> RGBColor[0.22`, 1.`, 0.44`]];
Show[{plt1, plt2}]

How to use NDSolve with limited RAM?

2 Answers2