I'm seeking to find solutions to a numerical integration with a large set of parameter combinations (basically, I'm doing a brute parameter sampling). Yet, I believe the memory of the computer is storing unnecessary information, but I can't understand why: the memory stores near 16Gb of memory coming from Mathematica Kernel, while the output is only ~25mb. This sometimes makes the computer run out of memory and Kernel quitting. I am also afraid this may slow my computations.
I will abbreviate a lot the code because I just want to know if I'm doing something fundamentally wrong with my rationale, specially on point 3 below. Basically, I am numerically integrating some equations for many parameter values. The output of the integration is then used to compute values of the dependent variables at specific time points. These values are my output ~25 mb.
Note that I simplified my code a lot, and also left some of it blank because it is unnecessary to explain what I am doing. What I'm doing is the following:
1- Considering a very large table
Table4 = Sort[ Flatten[ Table[ {10^b, 10^g, 10^v, 10^k}, {b, -5, 5, 0.5}, {g, -7, 3, 0.5}, {v, -11, -1, 0.5}, {k, -8, 2, 0.5} ], 3] ]; 2- I have a function deqns4 that is simply a system of ODEs to be numerically integrated and is called several times in point 3 (This is a random system of equations, my real one has 13 dependent variables).
deqns4[b_, g_, v_, k_] := {A'[t] == ki H[t] + k X[t] - b A[t], X'[t] == g A[t] - v X[t] - kc X[t], A[0] == 2, X[0] == 1 }; 3- Then, the following function calls deqns4 several times, to **numerically integrate the equations and output the values of the dependent variables at some time points, for the given parameter values in Table4:
lsim4[{b_, g_, v_, k_}, ba_] := Block[{ep, initconds, sssol0, tsim4, tsim3, tsimc50, tsimc100, c1, c2, c3, c4, c5, c6, c7}, (*This finds initial conditions*) sssol0 = Cases[ NSolve[ Rationalize[ Evaluate[{ ki H + k X - b A == 0, g A - v X - kc X == 0} /. ki -> 2 /. kc -> 0.1 ], 0 ], {A, X} ], {Rule[_, _?NonNegative] ..} ] /. {A -> A0, X -> X0}; (*This and the following commands use the initial conditions above to numerically integrate the system of equations "deqns4"*) tsim4 = Flatten[ NDSolve[ Rationalize[ Evaluate[ Join[ {H[t] == 0.009 Exp[-0.1 t] + ba }, deqns4[b, g, v, k] ] /. kc -> 0.1 /. sssol0], 0 ], {A, X, H}, {t, 0, 9720}, MaxSteps -> 1*^12 ], 1 ]; (*tsim3, tsimc50, and tsimc100 are also numerical integration of deqns4, similarly to tsim4*) tsim3 = Flatten[NDSolve[(*...left blank to simplify...*)], 1]; tsimc50 = Flatten[NDSolve[(*...left blank to simplify...*)], 1]; tsimc100 = Flatten[NDSolve[(*...left blank to simplify...*)], 1]; (*then I want the output to be a set of values computed from the \ numerical integration above*) c1 = (A[324] + X[324]) 100 /. tsim4; c2 = (A[9720]) 100/(A[0]) /. tsim4; c3 = (X[3240]) 100 /. tsim4; c4 = (X[9720]) 100 /. tsim4; c5 = (A[3240] + X[3240]) 100 /. tsim3; c6 = (X[3240]) 100 /. tsimc100; c7 = (X[3240]) 100 /. tsimc50; {c1, c2, c3, c4, c5, c6, c7} ]; lsimsols = ParallelMap[lsim4[#, 1.2] &, Table4]; (*output is a list with the form lsimsols={{0.1,0.99,100,0.52,76,100,100},{0.01,0.99,90,53,75,100},etc}*) My questions are the following:
- Is
NDSolvestoring information each time it runs, though it is within aBlockfunction? - The function
lsim4is only outputingc1toc7for each element ofTable4. This is what is being stored bylsimsols, and is the only information I am interested in. Yet, a lot more data is being stored by the kernel. Is there any way for me to avoid storage of so much information? - Even without specifying
MaxStepSize, takes me ~4 hours to run all the values of parameters. Is there any more efficient way of doing this parameter sampling other thanParallelMap?
Thanks in advance!
Edit:
- as suggested in the comments, I set
$HistoryLength=0which did not change all the memory allocated to Mathematica - I still get about 16 Gbs of stored information. - Also,
ParametricNDSolveseems to produce faster results thanNDSolve.
Edit #2:
- Someone suggested to use
Mapinstead ofParallelMap, as this would possibly allow less memory usage. This didn't make any improvements. - The only way I was able to attain a partial solution to the problem above is by splitting
Table4in multiple parts and rundeqns4for each part separately and then joining everything.
$HistoryLengthto something other thanInfinity(the default)? $\endgroup$ParametricNDSolve. As pointed out above try$HistoryLength=0$\endgroup$ParametricNDSolve. I think I've used it before. $\endgroup$$HistoryLength=0but Mathematica still steals all memory. UsingParametricNDSolvemade the computations faster though, but I then have to export the output, quit the kernel, import, and only then I can actually do something in mathematica $\endgroup$