Re: How to solve a stiff ODE system?

Orlando wrote:
Narasimham wrote:
Orlando wrote:
d Y
---- = f(t, Y), Y(t_0) = Y_0
d t
Warning: Failure at t=1.533227e-003. Unable to meet integration
tolerances without reducing the step size below the smallest value
allowed (5.447117e-018) at time t.

Suspect stiffness is not a problem,cannot say further as you do not
indicate f(t,Y).For eps =10^-8, if initial condition Y(t_eps) = Y_0 is
used and ode behaves well, then problem was with usage of first order
equation.To fix such cases switch over to a second order differential
equation providing two initial conditions.

Thank you very much for your answer!

The ode system of mine is as follows:

dy(1) / dt = y(2)
dy(2) / dt = - (y(3) - y(4)) / epsilon0 / E^2
dy(3) / dt = 1 / epsilon0 / y(2) / E^2 * (y(3)^2 - y(3) * y(4) * (1 -
epsilon0 * Ri / mu_+ / e))
dy(4) / dt = 1 / epsilon0 / y(2) / E^2 * (y(3) * y(4) * (1 - epsilon0 *
Ri / mu_- / e) - y(4)^2)

Where,
epsilon0: permittivity of vacuum 8.85e-12
E = E(t): the electric field. it is a function of t. After integrating
one step of the ODE, the electric field should be recalculated from
y(1), y(2), y(3), y(4), and this formula doesn't belong to the ODE
system, so it is not given here.
Ri: the recombination coefficient of bipolar charge. 2.2e-12
mu_+, mu_-: the mobility of positive and negative charge, 1.6e-4
e: elementary charge 1.6e-19

The above system is used to model the space charge flow of the high
voltage direct current power line.

I still don't quite follow you: how should I switch to a second order
differential equation? What is the reason of doing this? And how to
apply this procedure to my ODE system?

Thank you again!!

Tried to solve the ODE numerically. It crashes at start of any chosen
time interval, even as an initial value problem,found no errors.
Perhaps you need to at first check the physics modeling even as an
IVP.In a simulation program start with the simplest electrostatic
situations like equipotentials or force lines and equations depicting
them.When it works,gradually increase its complexity in steps.That
traps the error at a particular stage.

Narasimham

.

Relevant Pages

  • Re: How to solve a stiff ODE system?
    ... Unable to meet integration ... The ode system of mine is as follows: ... I still don't quite follow you: how should I switch to a second order ... apply this procedure to my ODE system? ...
    (sci.math)
  • Re: How to solve a stiff ODE system?
    ... Unable to meet integration ... the recombination coefficient of bipolar charge. ... I still don't quite follow you: how should I switch to a second order ... apply this procedure to my ODE system? ...
    (sci.math)
  • Re: How to deal with such ODE with fortran?
    ... I have to an integrate different ODE for each kp. ... without causing confusion in the recursive procedure. ... RHscalls P, but the subroutine doing integration only ... I am not familiar with this kind of equation, so I am now writing codes ...
    (comp.lang.fortran)
  • Re: Runge-Kutta question
    ... point of phase transitions) and, therefore, a numerical solution goes crazy. ... Try regularizing the heat capacity by replacing the jump with a steep, ... Your ODE will then have an analytical ... the integration step to really small values. ...
    (sci.math.num-analysis)
  • Re: How to deal with such ODE with fortran?
    ... I have to an integrate different ODE for each kp. ... without causing confusion in the recursive procedure. ... RHscalls P, but the subroutine doing integration only ... available literature for suitable numerical methods? ...
    (comp.lang.fortran)