Group Abstract Group Abstract

Message Boards Message Boards

Index an arbitrary array in NDSolve to simulate white noise?

Sandip Roy
Posted 7 years ago

I'm trying to introduce a randomized element in a simulation, using NDSolve to do the integration. The random element is supposed to simulate a random thermal force. The way I've tried to implement this is by creating an array of random values and then trying to select values from the array during each time step. The following is a snippet of my current code.

(*Initialising some time values that will allow me to select particular values from my random array*)
Tfinal = 10;
SampleFreq = 2500;
NoT = Tfinal*SampleFreq + 1;
dt = Tfinal/(NoT - 1);

(*Array of random values and the function I'll be using in NDSolve*)
WN = RandomVariate[NormalDistribution[0, 1], NoT];
WNoise[t_] := Indexed[WN, Floor[Floor[t/dt]] + 1]

(*NDSolve function using white noise function*)
sol = x[t] /. 
  NDSolve[{D[x[t], t, t] == -(w^2) x[t] + WNoise[t], x[0] == 0, 
    x'[0] == 2 }, x, {t, 0, 10}]

The function above (without white noise) is the classic F = -kx equation. I've just tried to introduce a white noise driving force force to it as well. Unfortunately, it doesn't work when the WNoise function is added. The error message I typically get is the following.

Indexed::ind: The index 1/2 is not a nonzero integer.
NDSolve::ndnum: Encountered non-numerical value for a derivative at t == 0.`.

I'm just wondering if anyone has faced a similar problem and if this method of introducing white noise is possible or if I'm just being silly with this implementation.

Thank you for your help.
POSTED BY: Sandip Roy
7 Replies
Sort By:

Other responses are probably a better way to go but I want to mention that this can be made to work with very simple modifications.

(1) Give w a value.

(2) Define WNoise so that it only evaluates for explicit numeric values of t.

WNoise[t_?NumberQ] := Indexed[WN, Floor[t/dt] + 1]

(3) Also it would be good to maintain consistency with time parameters, so run NDSolve from 0 to Tfinal.

A possible advantage to this approach is that it does not try to create a differentiable function from the random impulses. Which is no guarantee NDSolve won't treat it as differentiable using numeric finite differences (though I do not think it will do that).
POSTED BY: Daniel Lichtblau
Anonymous User
Anonymous User
Posted 7 years ago

I forgot the most important principle (always a bad thing to do).

Observe F sin(wt+B) as the (counter/added) force function to -kx in time, and that it must have some ENERGY (work) (or F is 0, nothing).

If the force is limited, the effect of it will be limited if there is any DAMPING (imagine a spring oscillating in water, the water will stop the motion of the oscillating spring. If there were a "added force" to the oscillation it would need to BE MORE THAN the water damping to have any OVERALL IMPACT).

The overall analysis should be done first, which is that the energy of the random "addition" matters, and more so if the damping is greater than randomness matters more than that.
POSTED BY: Anonymous User
Anonymous User
Anonymous User
Posted 7 years ago
POSTED BY: Anonymous User

Try searching http://reference.wolfram.com/language/ref/ItoProcess.html (the doc. page for ItoProcess) for "white" (or "white noise"). You'll should find an example very similar to what you're trying to do.
POSTED BY: Michael Rogers
POSTED BY: Neil Singer
Sandip Roy
Sandip Roy
Posted 7 years ago

Hi, Neil.

Thank you so much! Your solution is exactly what I was looking for.

However, I am quite new to Mathematica and I was wondering if you could explain to me what exactly the function Interpolation[Transpose[{Range[0, 10, 1/SampleFreq], WN}]] is doing. Naively, to me, it looks like it's taking the WN array and interpolating between the different values, with a frequency of SampleFreq from t=0 to t=10. However, I don't quite understand the purpose of the Transpose function.

Thank you again for your help, Neil.

Best Regards,

Sandip

POSTED BY: Sandip Roy

Sandip,

You can't mix a digitized noise signal with a continuous integration that way. The derivatives are evaluated along continuous time points and are not on fixed intervals. You have several options. You can use WhenEvent to trigger an event at each sample point and add in the noise or you can use Interpolation to make your noise continuous. I'll show the later. You also need to decide the frequency of the noise (SampleFreq) because a very high rate will significantly slow down the integration.

Tfinal = 10;
SampleFreq = 2500;
NoT = Tfinal*SampleFreq + 1;
dt = Tfinal/(NoT - 1);

(*NDSolve function using white noise function*)

WN = RandomVariate[NormalDistribution[0, 1], NoT];
noise = Interpolation[Transpose[{Range[0, 10, 1/SampleFreq], WN}]];
sol = x[t] /. 
   NDSolve[{D[x[t], t, t] == -(w^2) x[t] + 100*noise[t], x[0] == 0, 
     x'[0] == 2}, x, {t, 0, 10}];
Plot[sol[[1]], {t, 0, 10}]

Regards,

Neil
POSTED BY: Neil Singer
Reply to this discussion
Community posts can be styled and formatted using the Markdown syntax.
Tag limit exceeded
Note: Only the first five people you tag will receive an email notification; the other tagged names will appear as links to their profiles.
Reply Preview
Attachments
Remove
or Discard