13
$\begingroup$

After some googling, i've found similar problems around, but didn't find a 100% satisfactory answer, so let me ask here:

I'd like to solve a 1+1 problem using the method of lines. In spherical symmetry, I discretize the radial direction and end up with a system of ODEs that I'd like to solve with NDSolve. Problem is, the original 1+1 problem is a system of ~10 first-order PDEs and eventually I need to solve ~10*N first-order ODEs, where N are the radial grid points. To achieve enough resolution, N ranges from O(10^3) to O(10^4).

My problem is that NDSolve requires a lot of memory to solve the system. Also, I need the solution at any time, not only at the final time. What I would like to do is having NDSolve storing intermediate values at each time step in a data file, thus leaving the memory free. Is that possible? Is there any other solution to this problem? Frankly, I find this problem quite restrictive, given that a similar odeint in c would require essentially no memory to run, unless i'm doing some silly mistake of course.

Here is an example:

$HistoryLength = 0;
Timing[
Nv = 10^4;
tF = 40;
variables = Table[Subscript[ϕ, i], {i, 1, Nv}];
eqs = Table[
Subscript[ϕ, i]''[t] + i Subscript[ϕ, i]'[t]^3 + 
1/i Subscript[ϕ, i][t] == 0, {i, 1, Nv}];
initc = Table[{Subscript[ϕ, i][0] == 1, 
Subscript[ϕ, i]'[0] == 1}, {i, 1, Nv}] // Flatten;

ruleP = {MaxStepSize -> 0.01, MaxSteps -> 10^5};
sol = NDSolve[Union[eqs, initc], variables, {t, 0, tF}, ruleP];]
MemoryInUse[]/1024^2 // N

In my machine, this takes ~30 seconds and approx. 1 GByte to run! I'm not concerned about the speed, but the memory rapidly increases when the parameters Nv and tF increase and the notebook quickly sucks all machine's memory...

Improve this question
$\endgroup$
5
  • 1
    $\begingroup$ I think you should check MaxMemoryUsed to get an idea about how much memory is really used in intermediate steps (that's even worse :-). While I've seen cases where creating the final result for PDEs took a lot of time and memory that doesn't seem to be a problem for ODEs where MaxMemoryUsed isn't much different whether you create the full final result or not (the differece in MemoryInUse seems to result from the 10^4 interpolating functions in sol). One thing you could try is to use a matrix for for your system, but I have no idea how that effects memory usage... $\endgroup$
    – Albert Retey
    Sep 20, 2013 at 9:05
  • 1
    $\begingroup$ I just tried my suggestion to formulate the same system in matrix form and it turns out that it is using about 60% of the intermediate memory (MaxMemoryInUse), about the same final memory (MemoryInUse) but is about a factor 8 slower, which I find hard to understand. Another thing you could try is to stop and restart in ceratin intervals either as described here or manually with a fresh call to NDSolve with the final results as initial conditions. $\endgroup$
    – Albert Retey
    Sep 20, 2013 at 9:30
  • $\begingroup$ Thanks for your reply Albert. Indeed, I had tried to split NDSolve as in the link you sent. In that case, it saves a lot of memory but at expenses of being way slower. Doing it manually would be even worse, because each time the full NDSolve is called it will also reprocess the system of equations (whereas doing a proper splitting one can avoid to repeat that operation for each time step). However, what I'd like to do is really to store all data in a file, thus leaving the memory potentially free. Again, i guess one might to that manually, but I'm afraid it will be veeeery slow. $\endgroup$
    – Paolo
    Sep 20, 2013 at 18:25
  • $\begingroup$ My guess is that saving to file will be so slow that the overhead for stop/restart might be less relevant if not neglectable. The art will be to find the right balance between speed and memory consumption by adjusting the intervals between saving/restarting. Honestly it is out of my scope to assess whether the memory consumption of NDSolve for these cases is within the expected range or indicates that something isn't right. The difference to an optimized (less general?) solver written in C might well be the price you pay for all the extras that NDSolve provides... $\endgroup$
    – Albert Retey
    Sep 21, 2013 at 15:21
  • $\begingroup$ Hi, I am also facing similar problem. Do you guys able to solve this problem. $\endgroup$
    – santosh
    Jan 29, 2014 at 0:41

1 Answer 1

Reset to default
1
$\begingroup$

After digging around for a bit I found the tutorials in the standard documentation: FEMDocumentation/tutorial/SolvingPDEwithFEM and FEMDocumentation/tutorial/FiniteElementProgramming. Both have sections on methods to trade memory/time/accuracy, but the first tutorial is way easier for people new to Mathematica!

The most promissing proposals seem to be:

Method -> {"PDEDiscretization" -> {"FiniteElement",
"MeshOptions" -> {"MaxCellMeasure" -> 0.00025}}}

Method -> {"PDEDiscretization" -> {"FiniteElement", 
"MeshOptions" -> {"MeshOrder" -> 1}}}

Method -> {"PDEDiscretization" -> {"FiniteElement", 
"LinearSolveMethod" -> {"Pardiso"}}}

One can add these as options in NDSolve. The last one, changing LinearSolveMethod seems to have many options, and I am still trying to figureout how to use the example shown by Mr.Wizard on how to find those options

Improve this answer
$\endgroup$

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.