5
$\begingroup$

I am using NDSolve for a Langevin dynamics problem. I want to the know long term behaviour of my system ($t>1$) but it has to be simulated with very small time steps ($dt\sim 10^{-9}$). An example code:

R[t_Real]:= RandomVariate[NormalDistribution[0,1]] NDSolve[ {(-x''[t] - k1*x'[t] + k2*R[t])==0, x[0]==0, x'[0]==0} //.values, x, {t, 0, 10}, StartingStepSize-> 10^-9, Method->{"FixedStep",Method->"ExplicitEuler"}, MaxSteps->\[Infinity] ]

The Problem: My computer runs out memory when trying to store $10^{10}$ data points necessary for this computation. Is there a way to sample and store only a small subset of all the integration points?

Improve this question
$\endgroup$
8
  • $\begingroup$ Have you tried [Mma's RandomVariate[...]](reference.wolfram.com/mathematica/ref/RandomVariate.html)? What does values do? Is there a reason you are using ExplicitEuler? You would benefit by letting Mathematica choose the solver automatically or use perhaps the magic-wand Method->"LSODA". Also, why not just let t range from 0 to a very large number (for long term behavior)? $\endgroup$ Commented Apr 24, 2014 at 13:46
  • $\begingroup$ Hi, the 'values' variable is used to store parameters k1 and k2. I do not trust the automatic solvers of Mathematica because 'R[t]' is not a continuos function. For me $t \sim 1$ is a long time behaviour. The problem presented here has a damping of momentum (second term) that happens on the timescale of $10^{-9}$. The momentum has to relax from thermal kicks obtained from 'R[t]' term. For this reason I am using constant time sampling. $\endgroup$ Commented Apr 24, 2014 at 22:32
  • $\begingroup$ May I ask if this molecular dynamics or some such? What is the physical phenomena that you are trying to model? $\endgroup$ Commented Apr 24, 2014 at 23:23
  • $\begingroup$ It is a simple 1D Molecular Dynamics model for free diffusion $\endgroup$ Commented Apr 24, 2014 at 23:40
  • $\begingroup$ I performed simulations on MD for my Master's degree and I am curious, if you don't mind. Are you simulating Argon? With periodic BCs? Is the eqn you are solving the convection-diffusion equation? In 2012 I had a similar question on SE and I am wondering if it would help you. Good luck! $\endgroup$ Commented Apr 25, 2014 at 11:56

1 Answer 1

Reset to default
4
$\begingroup$

You can have NDSolve do the integration but not save the results (set second argument to {}). Then use EvaluationMonitor to save points at whatever interval dt you please. This uses much less memory; in fact for dt = 10^-2, the memory use is negligible for the settings below.

Since parameters were not given by the OP, I used the simplest choices. Also, waiting for ten billion steps seemed silly for a proof-of-concept trial, so I lengthened the step size.

MaxMemoryUsed[] SeedRandom[1]; R[t_Real] := RandomVariate[NormalDistribution[0, 1]]; values = {k1 -> 1, k2 -> 1}; lastt = -1; dt = 10^-2; {foo, {pts}} = Reap@NDSolve[{(-x''[t] - k1*x'[t] + k2*R[t]) == 0, x[0] == 0, x'[0] == 0} /. values, {}, {t, 0, 10}, StartingStepSize -> 10^-5, Method -> {"FixedStep", Method -> "ExplicitEuler"}, MaxSteps -> ∞, EvaluationMonitor :> If[t >= lastt + dt, lastt = t; Sow[{t, x[t]}]] ]; sol = Interpolation[pts]; MaxMemoryUsed[] Plot[sol[t], {t, 0, 10}] (* 37243360 *) (* 37243360 *)

Mathematica graphics

When the integral x is requested, the memory use by NDSolve is significantly higher:

MaxMemoryUsed[] SeedRandom[1]; R[t_Real] := RandomVariate[NormalDistribution[0, 1]]; values = {k1 -> 1, k2 -> 1}; lastt = -1; dt = 10^-2; {foo, {pts}} = Reap@NDSolve[{(-x''[t] - k1*x'[t] + k2*R[t]) == 0, x[0] == 0, x'[0] == 0} /. values, x, {t, 0, 10}, StartingStepSize -> 10^-5, Method -> {"FixedStep", Method -> "ExplicitEuler"}, MaxSteps -> ∞, EvaluationMonitor :> If[t >= lastt + dt, lastt = t; Sow[{t, x[t]}]] ]; sol = Interpolation[pts]; MaxMemoryUsed[] Plot[x[t] /. First[foo] // Evaluate, {t, 0, 10}] (* 37243360 *) (* 127874624 *)

Mathematica graphics

Improve this answer
$\endgroup$

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.