Thread Subject: high stiff problem

Hi,

I am trying to simulate an almost spontaneous reaccion.

this means that the differential equation is very very stiff.

dCdt = k*Ca*Cb

Cb =~ 0 and k = 1e7.

because the concentration of Cb is so near to 0 and, due to numerical error becommes negative sometimes, the solvers make to small steps that need many iterations to solve the problem.

I have tried nonnegative and other ways to avoid Cb being negative, but the solvers doesnt seem to be able to solve it fast.

I am also prety shure I tryed all the solvers available.

Is there a way to overcome this problem without having to change my equation?

Is there I solver I might have overseen?

Have you tried any of our solvers meant specifically for stiff systems,
such as ode15s, ode23s, ode23t, ode23tb?

http://www.mathworks.com/access/helpdesk/help/techdoc/ref/ode23.html

Also, I'm not sure if I understand the equation exactly. Is it:

dCa/dt = k*Ca*Cb

where k and Cb are constants? Is Ca the only state variable? If so, the
sign of k*Cb must be negative for the system to be stable. Can I assume
from what you say that Cb < 10e-7? If this is the case, and k and Cb are
constants, then abs(k*Cb) < 1, and the system has slow dynamics and is
not stiff (but perhaps unstable?).

-James

nico cruz wrote:
> Hi,
>
> I am trying to simulate an almost spontaneous reaccion.
>
> this means that the differential equation is very very stiff.
>
> dCdt = k*Ca*Cb
>
> Cb =~ 0 and k = 1e7.
>
> because the concentration of Cb is so near to 0 and, due to numerical
> error becommes negative sometimes, the solvers make to small steps that
> need many iterations to solve the problem.
>
> I have tried nonnegative and other ways to avoid Cb being negative, but
> the solvers doesnt seem to be able to solve it fast.
>
> I am also prety shure I tryed all the solvers available.
>
> Is there a way to overcome this problem without having to change my
> equation?
>
> Is there I solver I might have overseen?

Hi,
thanx for answering.

I have tryed all these solvers.
k is a costant, but Cb is the concentration of the component B

the whole equation looks like this

 dCd/dt = k*Ca*Cb
 dCa/dt = -k*Ca*Cb+D*Caext -> (a difussion term and the external concentration of component A)
 dCb/dt = -k*Ca*Cb

The problem is that I get a very ugly curve near cero for Ca. it is very stiff and sometimes negative.

Non of the solvers on the matlab help can do it without reducing the step drastically.

do you know any other solution?

James Allison <james.allison@mathworks.com> wrote in message <hlme2a$ak2$1@fred.mathworks.com>...
> Have you tried any of our solvers meant specifically for stiff systems,
> such as ode15s, ode23s, ode23t, ode23tb?
>
> http://www.mathworks.com/access/helpdesk/help/techdoc/ref/ode23.html
>
> Also, I'm not sure if I understand the equation exactly. Is it:
>
> dCa/dt = k*Ca*Cb
>
> where k and Cb are constants? Is Ca the only state variable? If so, the
> sign of k*Cb must be negative for the system to be stable. Can I assume
> from what you say that Cb < 10e-7? If this is the case, and k and Cb are
> constants, then abs(k*Cb) < 1, and the system has slow dynamics and is
> not stiff (but perhaps unstable?).
>
> -James
>
> nico cruz wrote:
> > Hi,
> >
> > I am trying to simulate an almost spontaneous reaccion.
> >
> > this means that the differential equation is very very stiff.
> >
> > dCdt = k*Ca*Cb
> >
> > Cb =~ 0 and k = 1e7.
> >
> > because the concentration of Cb is so near to 0 and, due to numerical
> > error becommes negative sometimes, the solvers make to small steps that
> > need many iterations to solve the problem.
> >
> > I have tried nonnegative and other ways to avoid Cb being negative, but
> > the solvers doesnt seem to be able to solve it fast.
> >
> > I am also prety shure I tryed all the solvers available.
> >
> > Is there a way to overcome this problem without having to change my
> > equation?
> >
> > Is there I solver I might have overseen?

With k so large the dynamics may be very fast, and you need a step size
small enough to capture the dynamics of the system. There may not be a
way to avoid a very small step size. How is Caext calculated, and what
is D? I can have a deeper look if you provide enough information to
reproduce this.

-James

nico cruz wrote:
> Hi, thanx for answering.
>
> I have tryed all these solvers.
> k is a costant, but Cb is the concentration of the component B
>
> the whole equation looks like this
>
> dCd/dt = k*Ca*Cb
> dCa/dt = -k*Ca*Cb+D*Caext -> (a difussion term and the external
> concentration of component A)
> dCb/dt = -k*Ca*Cb
>
> The problem is that I get a very ugly curve near cero for Ca. it is very
> stiff and sometimes negative.
>
> Non of the solvers on the matlab help can do it without reducing the
> step drastically.
>
> do you know any other solution?
>
> James Allison <james.allison@mathworks.com> wrote in message
> <hlme2a$ak2$1@fred.mathworks.com>...
>> Have you tried any of our solvers meant specifically for stiff
>> systems, such as ode15s, ode23s, ode23t, ode23tb?
>>
>> http://www.mathworks.com/access/helpdesk/help/techdoc/ref/ode23.html
>>
>> Also, I'm not sure if I understand the equation exactly. Is it:
>>
>> dCa/dt = k*Ca*Cb
>>
>> where k and Cb are constants? Is Ca the only state variable? If so,
>> the sign of k*Cb must be negative for the system to be stable. Can I
>> assume from what you say that Cb < 10e-7? If this is the case, and k
>> and Cb are constants, then abs(k*Cb) < 1, and the system has slow
>> dynamics and is not stiff (but perhaps unstable?).
>>
>> -James
>>
>> nico cruz wrote:
>> > Hi,
>> > > I am trying to simulate an almost spontaneous reaccion.
>> > > this means that the differential equation is very very stiff.
>> > > dCdt = k*Ca*Cb
>> > > Cb =~ 0 and k = 1e7.
>> > > because the concentration of Cb is so near to 0 and, due to
>> numerical > error becommes negative sometimes, the solvers make to
>> small steps that > need many iterations to solve the problem.
>> > > I have tried nonnegative and other ways to avoid Cb being
>> negative, but > the solvers doesnt seem to be able to solve it fast.
>> > > I am also prety shure I tryed all the solvers available.
>> > > Is there a way to overcome this problem without having to change
>> my > equation?
>> > > Is there I solver I might have overseen?

One more suggestion: it looks like with slight modification you could
put your system into linear state-space form. You then could take
advantage of the efficiency of lsim. ODEXX algorithms are for general
nonlinear problems; lsim exploits the linearity of a system, and has
impressive speed and accuracy.

If you can linearize your system then you could also diagonalize it,
which enables you to simulate each component of the system
independently. You will still have at least one fast component, and it
will require many times steps to simulate, but the slower components
could be simulated much faster. Afterwards you can transform the results
back into the original state space.

-James

nico cruz wrote:
> Hi, thanx for answering.
>
> I have tryed all these solvers.
> k is a costant, but Cb is the concentration of the component B
>
> the whole equation looks like this
>
> dCd/dt = k*Ca*Cb
> dCa/dt = -k*Ca*Cb+D*Caext -> (a difussion term and the external
> concentration of component A)
> dCb/dt = -k*Ca*Cb
>
> The problem is that I get a very ugly curve near cero for Ca. it is very
> stiff and sometimes negative.
>
> Non of the solvers on the matlab help can do it without reducing the
> step drastically.
>
> do you know any other solution?
>
> James Allison <james.allison@mathworks.com> wrote in message
> <hlme2a$ak2$1@fred.mathworks.com>...
>> Have you tried any of our solvers meant specifically for stiff
>> systems, such as ode15s, ode23s, ode23t, ode23tb?
>>
>> http://www.mathworks.com/access/helpdesk/help/techdoc/ref/ode23.html
>>
>> Also, I'm not sure if I understand the equation exactly. Is it:
>>
>> dCa/dt = k*Ca*Cb
>>
>> where k and Cb are constants? Is Ca the only state variable? If so,
>> the sign of k*Cb must be negative for the system to be stable. Can I
>> assume from what you say that Cb < 10e-7? If this is the case, and k
>> and Cb are constants, then abs(k*Cb) < 1, and the system has slow
>> dynamics and is not stiff (but perhaps unstable?).
>>
>> -James
>>
>> nico cruz wrote:
>> > Hi,
>> > > I am trying to simulate an almost spontaneous reaccion.
>> > > this means that the differential equation is very very stiff.
>> > > dCdt = k*Ca*Cb
>> > > Cb =~ 0 and k = 1e7.
>> > > because the concentration of Cb is so near to 0 and, due to
>> numerical > error becommes negative sometimes, the solvers make to
>> small steps that > need many iterations to solve the problem.
>> > > I have tried nonnegative and other ways to avoid Cb being
>> negative, but > the solvers doesnt seem to be able to solve it fast.
>> > > I am also prety shure I tryed all the solvers available.
>> > > Is there a way to overcome this problem without having to change
>> my > equation?
>> > > Is there I solver I might have overseen?