| Author |
Comment/Response |
John Barber
|
 09/13/01 3:47pm
Hi Folks:
I have a small system of stochastic ODE's to solve. Here is the system:
q'[t]==p[t]
p'[t]==-(1+Q[t]^2)q[t]
here Q[t] is a random variable, i.e., it takes on a different random value at every instant of time. What makes this difficult is that the statistics of Q depend on the value of q[t] at that moment. Specifically, at each instant in time, Q is sampled from a normal distribution with mean zero and variance 1/(1+q[t]^2).
I wish to be able to numerically solve this ODE, and thus generate random trajectories in the (q,p) plane. (Naturally, each one will be different.) I can't seem to get it to work. I've tried the following, to no avail. (Yes, I remembered to load the Statistics`ContinuousDistributions` package):
NDSolve[{q'[t] == p[t],
p'[t] == -(1 + (Random[NormalDistribution[0, Sqrt[1/(1 + q[t]^2)]]])^2)q[
t], q[0] == 1, p[0] == 0}, {q, p}, {t, 0, 20}]
This generates the warning that the second argument of NormalDistribution[] is not a machine sized number.
Then, I tried creating a function and doing the following:
Q[b_] := Random[NormalDistribution[0, Sqrt[1/(1 + b^2)]]];
NDSolve[{q'[t] == p[t], p'[t] == -(1 + Q[q[t]]^2)q[t], q[0] == 1,
p[0] == 0}, {q, p}, {t, 0, 20}]
This generated the same error message.
Is it possible to do what I am trying to do in Mathematica, or do I have to write my own code?
(Just for reference: a correct solution to the above ODE should trace out a wiggly, Brownian-looking trajectory in the (q,p) plane that goes around in a very rough circle in the clockwise direction.)
Thanks,
John Barber
URL: , |
|
