Re: NDSolve -- n-body indexing ([[]]) problem
- To: mathgroup at smc.vnet.net
- Subject: [mg112233] Re: NDSolve -- n-body indexing ([[]]) problem
- From: Albert Retey <awnl at gmx-topmail.de>
- Date: Mon, 6 Sep 2010 04:13:54 -0400 (EDT)
- References: <i5vnri$2e1$1@smc.vnet.net>
Am 05.09.2010 11:28, schrieb Themis Matsoukas:
> The problem is that your differential equations are stochastic. I
> asked wolfram about it (see email below) using the following as a
> minimal example (which doesn't work):
>
> Remove@g1 g1[t_?NumberQ] := (RandomReal[] - 0.5)*t; sln = NDSolve[{
> y'[t] == g1[t], y[0] == 0 }, {y}, {t, 0, 1}] Plot[Evaluate[y /.
> sln[[1]]][t], {t, 0, 1}, PlotRange -> All, AxesOrigin -> {0, 0}]
I'm not an expert in the area of stochastic differential equation, but
it seems rather obvious that this isn't a meaningful formulation of one:
how often would you expect RandomReal to be updated, or g1 to be called?
I think what you sent can be interpreted as an approximation to a
stochastic differential equation _only if_ assuming a discretization
with fixed step size which will then fix the update rate and the actual
measure of randomness...
> Wolfram's answer is basically that you cannot do this NDSolve. Tech
> support suggested to use RecurrenceTable:
>
> Remove@g1 ; Clear [a, pp]; g1[t_?NumberQ] := (RandomReal[] - 0.5)*t;
> dt = 0.0001;
>
> pp = RecurrenceTable[{a[n + 1] - a[n] == dt*g1[n], a[0] == 0}, a[n],
> {n, 1, Round[1/dt]}] ; ListPlot[pp, PlotRange -> All]
>
> This, however, forces you to apply a simple Euler integration and
> that's unfortunate because you cannot take advantage of advanced
> integration algorithms (stiff systems, adjustable time steps, etc)
> unless you code them yourself. I hope that Mathematica will come up
> with a module similar to NDSolve that is truly numerical, i.e. it
> would work with tabulated data, not just functions that can be
> expressed in analytic form.
... but that would of course exclude variable step sizes, so I don't
think the suggestion of WRI isn't a limitation at all. Instead of
RecurrenceTable, you can of course force NDSolve to use a fixed step
solver, e.g.:
sln = NDSolve[{y'[t] == g1[t], y[0] == 0}, {y}, {t, 0, 1},
Method -> {"FixedStep", Method -> "ExplicitEuler"},
StartingStepSize -> 0.0001];
and then the result will look as you probably expected...
I think if you expect a variable step solver to make sense, you need to
rethink how to formulate your problem. I can't even imagine how coming
versions of Mathematica (or other systems) could make sense of the given
formulation if not assuming fixed step sizes.
You should also note that the same question has been discussed more than
once in this group, so when searching the archives (e.g. NDSolve Random)
you could find some more useful suggestions.
hth,
albert