Fwd: Help-octave Digest, Vol 32, Issue 51

[Francis Poulin]

Hello Chris,

The problem you are having is not with lsode but it's with the approach.  If you are using a random number you should not use a standard package to solve this equation.  The random number you are using is white noise and is not continuous so the solution is more complicated.  For example, if you are solving the following ode,

dx/dt = xi(t)

where xi(t) is a random number and you wanted to solve this numerically, there is a special Euler method for Stochastic DEs.  What you do is,

x^{n+1} = x^{n} + \sqrt{\Delta t} xi^n.

Notice that there is a square root because xi is stochastic and not a \Delta t as you would get for the regular Euler method.   As far as I know there are not any packages to solve SDEs.

I suggest you learn a bit more about Stochastic DEs.  A great introduction that I used is "An Introduction to Stochastic Processes" by Lemons.  

Good luck,

Francis

Message: 1

Date: Sun, 23 Nov 2008 08:57:28 +0100

From: Matthias Brennwald <address@hidden>

Subject: Re: lsode and random term

To: address@hidden

Message-ID: <address@hidden>

Content-Type: text/plain; charset=US-ASCII; delsp=yes; format=flowed

On 23.11.2008, at 06:32, address@hidden wrote:

Message: 1

Date: Sat, 22 Nov 2008 13:10:42 -0500

From: "William Christopher Carleton" <address@hidden>

Subject: lsode and random term

To: address@hidden

Message-ID: <address@hidden>

Content-Type: text/plain; charset="UTF-8"

Hi all,

I have a function (below), which includes a random term (I wrote a  

function 'random(x)' that just uses rand() to produce a value  

between -x and x b/c I couldn't find such an existing function in  

Octave's manual) and when I use lsode to solve it occasionally  

works, but often produces as error about convergence;

LSODE--  AT T (=R1) AND STEP SIZE H (=R2), THE

       CORRECTOR CONVERGENCE FAILED REPEATEDLY

       OR WITH ABS(H) = HMIN

      In above,  R1 =  0.3220694060927D+03   R2 =  0.1207113904987D-05

error: lsode: repeated convergence failures (t = 322.069perhaps bad  

jacobian supplied or wrong choice of integration method or tolerances)

error: evaluating assignment _expression_ near line 10, column 2

Can anyone tell me why the random term has caused this problem and  

what I can do about it? I want the random term, or I suppose an  

unpredictable function with magnitude of fluctuations determined by  

a variable selected before solving, to simulate unpredictable  

environmental fluctuations. Thanks,

Chris

Octave Code:

function xdot=f(x,t);

k=0.01;

c=.9;

p=0.025;

l=0.5;

C=1000;

r=0.001;

s=1.1;

d=0.002;

v=1;

E=random(1) * v;

xdot(1)=x(1) * (k * (x(3) - x(1)));

xdot(2)=(x(2) * -d * x(1) * E) + (r * x(2) * ((C - x(2)) / x(2)));

xdot(3)=x(2) * (x(1) * l) * p * (((x(2) * s) - x(3)) / (x(2) * s))  

- (c * x(1));

endfunction

I've tried something similar, and it did not work either. With the  

random term, the function value will never be the same for a given  

argument. So, if the solver evaluates the function twice for (almost)  

the same argument, the function value will not be (almost) the same.  

In other words, there solver cannot converge.

I haven't looked closely into your function with the random term, but  

if the random part is here to introduce some noise, then a solution  

would be to define the noise-term values before starting the solver.  

With this approach, the noise term will always be the same for a  

given argument, so the solver will be able to converge.

IHTH

Matthias