Numerical solution of retarded functional differential equations as partial differential equations

IFAC

Copyright 'C IFAC Linear Time Delay Systems, Ancona , Italy , 2000

c:

0

C>

Publications www.elsevier.com..1ocate ·ifac

NUMERICAL SOLUTION OF RETARDED FUNCTIONAL DIFFERENTIAL EQUATIONS AS PARTIAL DIFFERENTIAL EQUATIONS Stefano Maset'

• Department 01 Mathematical Sciences, University 01 Trieste , Italy

Abstract : The initial value problem for Delay Differential Equations (DDEs) y'(t) = I(t, y(t), y(t - Td,·· · , y(t - Tn)), { y(t) = (t), t::; to .

(1)

is more complicated than the classical Cauchy problem from both theoretical and numerical point of view. In particular, the numerical integration of (1) requires the use of a continuous approximation such as, for instance, continuous Runge-Kutta methods. The need of continuous methods , rather than discrete methods, leads to new and difficult issues for the accuracy analysis and for the preservation of stability properties as well. Copyright :r 2000 IFAC Keywords : Delay Differential Equations , Functional Differential Equations, Runge-Kutta Methods , Methods of Lines , Asymptotic Stability

The classic approach for solving a Delay Differential Equation (DDE) y'(t)=/(t,y(t),y(t-T)), { y(t)=ip(t) -T::;t::;O

t2:0

and with an initial condition which depends on the initial data U(O,B)=ip(B), BE [- T,O] .

(2)

(5)

The solution of the DDE is then given by is the step-by-step integration where the delayed part y (t - T) is viewed as a known forcing term in an Ordinary Differential Equation (ODE) . A totally different approach is based on restate the DDE as a Partial Differential Equation (PDE) which consists of the advection equation

au Ft (t,B)

y (t

au

~~ (t, 0) = 1 (t , U (t , O) , U (t, -T)) ,

= U (t ,B) ,

t 2: 0, BE [-T ,O] . (6)

Indeed we can use this approach also on the more general class of the Retarded Functional Differential Equations (RFDEs) y'(t)=F(t , y(t) , ytl , t2:0 { y(t)=ip(t) -T::;t::;O

= aB (t , B) , t 2: 0, BE [-T , O] (3)

supplied with a boundary condition which depends on the DDE

+ B)

(7)

where Yt (B) = y (t + B) , B E [- T, O] . :\ow. the boundary condition (4) becomes

au

t 2: 0(4)

aB (t. 0)

133

= 1 (t , U (t , 0) • U (t , .)) ,

t

2: 0 (8)

By solving the PDE (3) with the method of lines we obtain a numerical solution for the DDE (2) or the RFDE (7). Schemes of discretization for (3) are characterized by triplets (d, w , W) , where d = (dl, ...,ds), 0 < d l < dz < ... < ds = 1, w E R S and W E RS xs , as follows For a given integer N 2: 1 we consider the mesh on the interval [-7, OJ N-I

nN

U {- (n + d

:= {O} U

i)

It I i = 1, ... , s},

n=O

U(t , _(n+ddh)) 1

for every

let us consider the partition

u = (u;r, [u()

'u~ (t) =

[U/J (t) {

(10)

U(O)='{!N' where F is a suitable finite-dimensional approximation of F, WI Wll . Wll

Ws WsI

Wll .

• 0

< Cl < ... < Cs

We will call such a RK method compatible. By a compatible method we obtain the triplet d := c,

. 'WI s

belongs to E RsNx (l+sN) and 'PN = ('P(8))OErl ".' The linear part of the ODE (10) is obtained froIll the PDE (3) by replacing the exact values

(12)

This approach for the resolution of DDEs could be superior to the classic step-by-step integration for the preservation of asymptotic stability. In fact , it is known (see (Maset, 2000a)) that no RK step-bystep integration preserves the asymptotic stability when applied to complex linear systems of DDE

TT!

t 2: 0 (13)

2: 2.

Unfortunately in (Maset, 2000b) it is proved that the systems of ODEs (10) relevant to complex scalar equations (13) does not preserve the asymptotic stability. In order to recover stability, we consider the full discretization of the method of lines. By using an A-stable RK method for timediscretization of (10) with constant stepsize tlt we conjecture (supported by numerical evidence) that in order to obtain preservation of asymptotic stability for systems (13) we have to choose tlt 2: K (It) where K(h) is bounded.

+ dIl h) 1. REFEREi\CES

(11) 88 (t, - (11

W:= _A-I .

y' (t) = Ly (t) + My (t - 1) , { y(t)=
Wj Wll

8U

w:= A-lIs,

A possible class of compatible RK methods is given by Radau HA methods. The case s = 1 (Implicit Euler method) is analyzed in (Bellen and Maset, 2000) for a constant coefficient linear system of DDE. Convergence analysis in the general case is studied in the forthcoming paper (Maset, 2000b).

of dimension

88 (t, - (n

= 1,

• A is invertible, • b= (UsI, .. ·,a ss )'

WI s

B=

8U

- 1.

Triplets (d, w, W) for discretization of (3) can be obtained by Runge-Kutta (RK) methods as follows. Let us consider an s-stage RK method (A , b, c) such that:

Wss WI

= 0, 1, .. . , N

U(t, -(~+ds)h)

T

where Uo E cm and [uJ E (Cm)flN \ {O} By discretizing the boundary condition (8) , the PDE (3) and the initial condition (5) we get the following (1 + sN) m-dimensional ordinary initial value problem

F (t, Uo (t) , [u J (t)) = (~B 0 Im) u (t)

Tt

.

(

It

(9) and the space X N := (Cm)flN ~ C (I+sN)m where m is the dimension of the RFD E (7) . For u E X N ,

.

+- (W 0 Im)

+ d s ) It)

with

~(W 0 Im).U(t , -nh)

Bellen , A. and S. Maset (2000) . Kumerical solution of constant coefficient linear delay differential equations as abstract Cauchy problems. Numer. Mathem. 84, 351- 374. ~laset , S. (2000a). Instability of runge-kutta methods when they are applied to linear systems of ddes. Numer. Mathem.

134

Maset, S. (2000b). Numerical solutions ofretarded functional differential equations as partial differential equations by higher order schemes .

135