Approximation of evolutional system using singular forcing

Computers and Chemical Engineering 26 (2002) 1013– 1021 www.elsevier.com/locate/compchemeng

Approximation of evolutional system using singular forcing Mikhail Skliar *, Prashant Tathireddy Department of Chemical and Fuels Engineering, Uni6ersity of Utah, 84112 Salt Lake City, UT, USA Accepted 22 February 2002

Abstract A method of finding an asymptotic approximation of the solution of a general class of linear evolutional time-varying systems is developed and applied to illustrative problems, including ordinary and partial differential evolutional systems. The evolutional problem is first forced into the singularly perturbed form by coordinate re-scaling. Singular perturbation methods are then used to find an approximate solution of the forced evolutional problem. The developed approach can be used to: find an approximate solution, often in analytical form, for a broad class of ordinary and partial differential equations; analyze steady-state and transient behavior of the ODE and PDE models; develop hybrid asymptotic-numerical methods for solving time-varying evolutional problems. © 2002 Elsevier Science Ltd. All rights reserved. Keywords: Evolutional system; Singular forcing; Perturbation methods

1. Introduction Perturbation methods (Kevorkian & Cole, 1981; Nayfeh, 1973; Hinch, 1991) are concerned with the study of the behavior of algebraic, integral, or differential problems, which depend on a small perturbation parameter m (or, inversely, large parameter 1/m), and their approximate solution. The problem is said to be regularly perturbed if solutions of the perturbed (finite m) and unperturbed (m “0) problems are close in some sense. Otherwise, the problem is perturbed singularly. An example of the singularly perturbed problem is the case when m multiplies the highest order derivative with respect to one of the independent variables: With m= 0 all initial or boundary conditions of the original problem can no longer be satisfied, and the difference between solutions of the perturbed and unperturbed problems may be arbitrarily large at least in some region of the problem domain. From the physical perspective, singularly perturbed ordinary and partial differential problems (O’Malley, 1991; Eckhause, 1979; Kokotovic, Khalil, & O’Reilly, 1986) describe systems which exhibit either vastly different gradients (different spatial scales) in different * Corresponding author. E-mail address: [email protected] (M. Skliar).

regions of spatial domain (Il’in, 1992), or vastly different speed of response (time scales) during their time evolution (Kevorkian & Cole, 1996; Holems, 1999). In either case, to apply the singular perturbation methods, the perturbation parameter must be explicitly introduced into the mathematical description of the applied problem. For a single-equation model, the perturbation parameter is typically introduced by casting the original model in a re-scaled and non-dimensional form. The model re-scaling usually involves the introduction of ‘characteristic’ time(s), length(s) or other selectable, non-unique parameters. For the system of equations exhibiting multi-scale behavior, the model in the standard singular perturbation form is typically obtained using coordinate transformation (Marino & Kokotovich, 1988; Kumar, Christofides, & Daoutidis, 1998) and such a transformation is in general not unique. Perturbation methods give the systematic way of approximating a singularly perturbed system by a set of simpler unperturbed subsystems. The approximation obtained using singular perturbation methods is asymptotic in nature, which means that the expected accuracy is achieved only when m“ 0. For small but finite m, an approximate solution in the form of the series of terms of progressively higher power of m may actually be divergent, and only a truncated m-dependent series may provide a useful approximation. The dependence of the

0098-1354/02/$ - see front matter © 2002 Elsevier Science Ltd. All rights reserved. PII: S 0 0 9 8 - 1 3 5 4 ( 0 2 ) 0 0 0 2 3 - 6

M. Skliar, P. Tathireddy / Computers and Chemical Engineering 26 (2002) 1013–1021

1014

approximate solution on non-unique m, a consequence of non-unique re-scaling or coordinate transformation, should not be used to cast a doubt on the usefulness of the approximations obtained using singular perturbation methods. On the contrary, a well-established theory of singular perturbation methods and their successful practical application despite non-unique and sometimes poorly motivated introduction of perturbation parameters leads us to question if the method of singular perturbations can be used to obtain approximate solutions for a broad class of unperturbed problems. For the case of linear evolutional systems, this paper gives a positive answer to the above question. An arbitrarily small perturbation parameter is introduced into the general unperturbed evolutional problem using a re-scaling procedure (referred to as singular forcing), and without assuming that the underlying process evolves in multiple time scales. The method of boundary functions (Vishik & Lyusternik, 1957; Trenogin, 1970; Vasil´eva, Butuzov, & Kalachev, 1995) is then used to find an asymptotic approximation of the singularly perturbed problem obtained by singular forcing. To insure that the asymptotic solution is independent from an arbitrary parameter m, the time-varying system operator is replaced by its Taylor series approximation. Finally, several examples are presented to illustrate the application of the developed method to ordinary and partial differential evolutional systems.

x(t,s;z)= U(t,s)x0(z)

(4)

Operator U satisfies the integral Volterra equation U(t,s)=I+

&

t

(5)

A(~)U(~,s)d~

s

the solution of which can be found using the method of successive approximations. We now consider the following singularly perturbed problem: (6(t%,z) m =A(t%)6(t%,z) (t% i(t%)6(t%,y)= hb(t%,y), 6(0,z)=x0(z),

(6) (7)

y(d

x0 D(A)

(8)

where 0B m 1. The solution of Eqs. (6)–(8) is searched in the same set of function as the solution of the problem (1)–(3). We are interested in the following question: When does the solution of the above problem (6)–(8) coincide with the solution of problem (1)–(3)? Operator 1/mA(t%) gives rise to the evolutional operator Um (t%,s), which is the semigroup operator in parameter t%. The following identity: Um (t%,s)= I+

& &

s

=I+

t%

1 A(~)U(~,s)d~ m

t%/m

A(my)U(my,s)dy

s/m

shows that 2. Singular forcing Consider problem:

the

Um (t%,s)= U(t/m,s) following

homogeneous

Cauchy

(1)

where operator A(t) is generally unbounded, and assume that for all t  [0,T] A(t) is defined in a dense domain D(A(t))=D(A) of Banach space E and is closed in it. Function x(t, z) is the solution of problem (1) if: (a) it has its value in D(A); (b) x(t, z) is defined for each t in the spatial domain d R n with the boundary (d, where it satisfies the boundary condition i(t)x(t,y)=hb(t,y),

y (d

(2)

(c) it has strong derivative in time; (d) it satisfies Eq. (1) on the interval [0, t]; and (e) at the initial time x(0,z)= x0(z),

For some t%j = mtj, and assuming zero initial time, 6(t%j,z)=Um (t%j,0)x0(z)= U(mtj /m,0)x0(z)= x(tj,z)

(x(t,z) = A(t)x(t,z) (t

x0  D(A)

(3)

For the uniformly correct Cauchy problem (Krein, 1967) there exists an evolutional operator U(t, s) (where s is the initial time), which for every t maps each x0  D(A) into the solution of the Cauchy problem (1) –(3) on the interval [s, T]:

(9)

(10)

Therefore, the solution of the Cauchy problem (1)–(3) on the interval [0, T] can be found from the solution of the problem (6)–(8) on the interval [0,mT]: x(t,z)= 6(mt,z)

(11)

In essence, casting the general evolutional problem in the singularly perturbed form of Eqs. (6)–(8) involves a simple re-scaling of the time variable. It is instructive to interpret the described re-scaling procedure from the perspective of different time scales. The multiplication of the derivative by a small parameter m without changing the right-hand part of Eq. (1) means that to maintain the equality in Eq. (6), the value of the derivative with respect to the ‘fast’ time t’ is ‘forced’ to be 1/m times larger than the corresponding derivative in natural time t. For an observer in natural time, it appears that the system (6) is forced to achieve its final state 1/m times faster. We refer to the described process of converting the evolutional problem into the singularly perturbed form as singular forcing.

M. Skliar, P. Tathireddy / Computers and Chemical Engineering 26 (2002) 1013–1021

Turning our attention to the case of inhomogeneous evolutional systems, it is easy to see that if f(z) is a smooth function with values in D(A), the solution of the following problem on the interval [0, T] (x(t,z) = A(t)x(t,z)+ f(z) (t

(12)

can be found from the solution of the following singularly forced problem with the same initial and boundary conditions, and defined on the interval [0,mT]: (6(t%,z) m =A(t%)6(t%,z) +f(z) (t%

(13)

and that x(t,z)=6(mt,z). For the case of the time varying smooth function f, which for each t[0,T] takes its value from D(A), the solution of the inhomogeneous problem (x(t,z) =A(t)x(t,z)+ f(t,z) (t

(14)

can be expressed in the Volterra integral form using the evolutional operator U: x(t,z)=U(t,0)x0(z)+

&

and ~=t%/m. When m“ 0, xf tends to zero everywhere within the interval [l,mT], where l is the small neighborhood of the initial time, known as temporal boundary layer. Therefore, we can interpret xf as the boundary correction to the regular solution xs(t%, z). The approximate solution xs(t%,z)+ xf(~,z) must satisfy the initial condition (8). This can be achieved by setting [0] x [0] s (0,z)+ x f (0,z)= x0(z)

(21)

x [is ](0,z)+ x [if ](0,z)= 0,

(22)

i= 1,2,…

[0] If initial conditions cannot be satisfied with x [0] s and x f terms, then require that i−1

x [is ](0,z)+ x [if ](0,z)= x0(z)− % [x [s j ](0,z)+ x [f j ](0,z)] j=0

(23) Boundary conditions (7) can be satisfied by requiring that i(t%)6(t%,y)=i(t%)[xs(t%,y)+xf(t%/m,y)]= hb(t%,y), (24)

y(d

t

U(t,s)f(s,z)ds

(15)

0

The last expression easily yields that solution of the problem (14) can be expressed as x(t,z) = 6(mt,z), where 6 is the solution of the following singularly forced problem defined on the interval [0,mT]: (6(t%,z) m =A(t%)6(t%,z) +f(t%/m,z) (t%

(16)

which can be achieved if [0] i(t%)x [0] s (t%,y)+i(t%)x f (t%/m,y)= hb(t%,y),

x [is ](t%,y)+x [if ](t%/m,y)= 0,

i =1,2,…,

3. Asymptotic solution of the forced problem

m

6(t%,z)= xs(t%,z)+ xf(~,z) +RN (t%,z)

(17)

where N

xs(t%,z)= % m x (t%,z) [i ] s

(18)

i=0 N

xf(~,z)= % m ix [if ](~,z)

(19)

i=0

RN =O(m N + 1)

when m “ + 0

(20)

(25)

y(d

(26)

If boundary conditions cannot be satisfied with x [0] s and x [0] terms, we require that for y(d f i(t%)x [is ](t%,y)+i(t%)x [if ](t%/m,y)

According to Vishik– Lyusternik’s method of boundary functions (Vishik & Lyusternik, 1957; Trenogin, 1970; Vasil´eva et al., 1995), we seek the solution of the problem (16), (7), (8) in the following form:

y(d

and

with the previously specified initial and boundary conditions. Note that since the state x(t, z) of the system (14) is achieved 1/m times faster by the system (16), the source term f must be synchronized accordingly during the forcing procedure.

i

1015

i−1

= hb(t%,y)− % [i(t%)x [s j ](t%,y)+i(t%)x [f j ](t%/m,y)]

(27)

j=0

Formal substitution of Eq. (17) in Eq. (16) yields: ( [0] 2 [2] [x (t%,z)+mx [1] s (t%,z)+m x s (t%,z)+…] (t% s

+

( [0] 2 [2] [x (~,z)+ mx [1] f (~,z)+ m x f (~,z)+ …] (~ f

2 [2] [1] = A(t%)[x [0] s (t%,z)+mx s (t%,z)+ m x s (t%,z)+ …] 2 [2] [1] +A(m~)[x [0] f (~,z)+ mx f (~,z)+ m x f (~,z)+ …] + f(t%/m,z) (28)

A meaningful approximate solution given by equations Eqs. (18) and (19) should not depend on an artificially introduced parameter m. To obtain m-independent solution in the case of time-varying operator A, we introduce its approximation by a series of stationary operators. Assuming that A(t) is N+ 1 times differentiable, obtain that

M. Skliar, P. Tathireddy / Computers and Chemical Engineering 26 (2002) 1013–1021

1016

A(t) =A(0)+ tA1 +t 2A2 +… +t NAN +O(t N + 1) (29) where each An is defined as



where the inverse operator A − 1(t) is formally introduced to explicitly express the solution for the second term in Eq. (18). For the kth terms (k=1,N),

Á à [k] ( [k − 1] −1 s[k] (t%,z) , i(t%)x [k] y(d Ãx s (t%,z)=A (t%) (t%x s s (t%,y)= h b (t%,y), à k Ã(x [k] f (~,z) [k] −i ] =A(0)x (~,z) + % ~ iAi x [k (~,z) f f à (~ i=1 Í k−1 Ãx [k](0,z)= x (z)− % [x [ j ](0,z) +x [ j ](0,z)] −x [k](0,z) f 0 s f s à j=0 à k−1 [j] [j] s[k] Ãi(~)x [k] f (~,y) =hb(~,y) − % [i(t%)x s (t%,y) + i(t%)x f (t%/m,y)]− h b (~,y), à j=0 Ä Anx(·) =



n

1 ( nA(t) n! (t n

n

x(·),

n =1,N

(30)

t=0

Therefore, for ~-dependent terms in Eq. (28), obtain that N (xf(~,z) =A(0)xf(~,z) + % m n~ nAnxf(~,z) +O(m N + 1) (~ n=1 (31)

Using the last expression in Eq. (28), and equating terms in the same power of m, separately for functions of t% and ~, we arrive at a sequence of problems, which determine the terms in the exponential expansion (17)– (19). For the first terms of xs and xf, obtain: [0] s ÁA(t%)x [0] s (t%,z) + f(t%/m,z)=0, i(t%)x s (t%,y)=h b(t%,y), y (d Ã(x [0](~,z) f [0] x [0] = A(0)x [0] Í f (~,z), f (0,z) =x0(z) −x s (0,z) Ã (~ s y(d Äi(~)x [0] f (~,y) =hb(~,y)−h b(~,y),

(32) where h sb is selected so that the corresponding equilibrium problem is well posed. For example, if the original Eq. (14) is parabolic with Neumann or mixed boundary conditions, then solution of the elliptic problem (first equation of system (32)) will exist only if an integral of fluxes over the entire boundary is equal to zero, which may not be true with boundary conditions (2). For the second terms of the asymptotic solution (N=1),



(34)

y(d

The preceding development of the singular forcing method, followed by the application of the boundary functions method to obtain the asymptotic approximation of the forced evolutional problem are summarized in the form of the following theorem. Theorem 1. Assume that the following conditions hold: 1. Linear operator A for each t[0,T] is defined and closed in the dense, time-independent domain D(A) of Banach space E. 2. In domain D(A), operator A(t) is N +1 times strongly differentiable with respect to t. 3. Operator A(t) generates a uniformly correct Cauchy problem (1) – (3). 4. Function f(t, z) is smooth on the inter6al [0, T], and takes its 6alues from D(A). Then the solution of the Cauchy problem (14), (2), (3) on the inter6al t[0,T] can be found as x(t,z)= 6(t%,z),

t%= mt

(35)

where 6(t%, z) is the solution of the singularly perturbed Cauchy problem (16) on the inter6al t%= [0,mT]. Furthermore, the solution x(t, z) can be expressed as x(t,z)= xs(t%,z)+xf(~,z)+ RN (t%,z)

(36)

t% = mt,

(37)

~= t%/m =t

where xs and xf are defined by the formal decompositions (18) and (19), members of which are found as the solu-

n

( [0] Á [1] s[1] Ã x s (t%,z)= A − 1(t%) x s (t%,z) , i(t%)x [1] y (V s (t%,y)=h b (t%,y), (t% Ã [1] (x (~,z) [0] Ã f = A(0)x [1] f (~,z) +~A1x f (~,z) (t Í 0 [ j] [ j] [1] Ã x [1] f (0,z)=x0(z)− % [x s (0,z) +x f (0,z)] −x s (0,z) j = 0 Ã 0 Ã i(~)x [1] (~,y)=h (~,y) − % [i(t%)x [s j](t%,y) + i(t%)x [f j](t%/m,y)]− h s[1] f b b (~,y), Ä j=0

y(V

(33)

M. Skliar, P. Tathireddy / Computers and Chemical Engineering 26 (2002) 1013–1021

tions of problems (32) and (34), and the norm of the remainder RN is of the order O(m N + 1) as m “0. Remark 1. Note that for time independent system operator A the approximation (29)– (30) is not needed, and −i ] Ai 0, i= 1,2,…, A(0) =A, ki = 1~ iAi x [k (~,z)=0, f [k] and all x f terms are found from the solution of the identical homogeneous evolutional problems, but with different initial and boundary conditions. Remark 2. The approximation of the time varying operator A given by Eqs. (29) and (30) is used to insure that the asymptotic approximation of the original evolutional problem does not depend on an arbitrary (and indeterminate) parameter m. For the finite m, the error of the truncated approximation (29)– (30) will grow without bound as t“ . Therefore, for the propagation problems with time varying operator A, the approximation (17)– (19) must be obtained sequentially on subintervals t[tm,tm +Dt], m =0,1,2,…, where Dt is selected to provide the desired accuracy in approximation of the operator A(t) according to the Taylor series expansion (29)– (30) (Zgurovskiy & Skliar, 1991a,b).

1017

y; = − 6y− 25x +sint For x(0)= 1 and y(0)= 0, the exact solution is given by

 

n

Á 305 103 2 1 − 3t + sint − cost Ãx= 408sin4t + 102cos4t e 51 102 Í Ãy= − 2563sin4t − 2 cos4t e − 3t + 1 sint + 2 cost Ä 408 51 102 51

n

and is plotted in Fig. 1(a). According to the Theorem 1, the first-order asymptotic solution is sought in the form [0] xas = x [0] s +xf

where

xas =

 n xas yas

and

 n

 n

Remark 3. It is straightforward to extend the singular forcing method for the case of systems of operator equations:

x [0] s =

(x1(t,z) = A1(t)x1(t,z) +A2(t)x2(t,z) + f1(t,z) (t

(38)

[0] Terms x [0] s and y s , found from the algebraic equation

(x2(t,z) = A3(t)x1(t,z) +A4(t)x2(t,z) + f2(t,z) (t

(39)

x [0] s , y [0] s

x [0] f =

x [0] f y [0] f

Ax [0] s + F=0

where both x1 and x2 can be vectors, and operators Ai (t), i= 1,4 and functions f1, f2 possess properties described in Theorem 1. Singular forcing in this case will produce an equivalent system in the standard singularly perturbed form, asymptotic approximation of which can be obtained by applying Vishik– Lyusternik’s theory of boundary functions.

with

4. Examples

x [0] s =

4.1. Approximation of the ODE model

The fast component of the asymptotic solution, x [0] f , is found from the homogeneous evolutional problem

First, we apply the described approach to the approximation of the following linear pendulum model: x¨ + 2hx; + … 2x =sint



0 − 25

n

1 , −6

F=

 n 0 sint

are equal sint , 25

y [0] s =0

[0] x; [0] f = Ax f

(40)

which for h =3 and … =5 can be equivalently written as x; = y

A=

with initial conditions

xf(0)= x(0)− x [0] s (0)=

n 1 0

M. Skliar, P. Tathireddy / Computers and Chemical Engineering 26 (2002) 1013–1021

1018

leading to the following overall first-order asymptotic approximation of the exact solution: xas =



n

1 3 sint + sin4t + cos4t e − 3t 25 4



yas = −

n

25 sin4t e − 3t 4

Using the result of Theorem 1, the higher order asymptotic approximations can be similarly obtained. Fig. 1(b) shows that the integral (t  [0,10]) of the absolute error of the approximate solution decays exponentially with the increased number of terms in the asymptotic solution. Fig. 1(c) and (d) illustrate the rapid decay of time-dependent approximation errors for x and y as the number of terms in the asymptotic solution is increased.

4.2. Approximation of the PDE model

(T(z,t) = AT(z,t)+F(z) (t (2 (z 2

A=h 2

T(z,0) = sin(yz),

T(0,t)=T(1,t)=0,

05t5 tf

where F(z)= sin(3yz) describes the heat generation. The analytical solution of this problem is given by: 2

T(z,t)=e − (yh) tsin(yz)+

1 2 [1− e − (3yh) t]sin(3yz) (3yh)2

The truncated asymptotic solution of the singularly perturbed problem is given as: [0] [1] [1] (z,t%)= x [0] s (z,t%)+ x f (z,~)+ mx s (z,t%)+ mx f (z,~)

+ O(m 2)

4.2.1. Time in6ariant operator A Consider the following inhomogeneous heat transfer problem:

05 z5 1

(41)

The first-order approximation is found from the following simplified problems:

Fig. 1. (a) Exact solution; (b) integral of the absolute error vs. the order of approximation, N + 1; (c) error in approximation solution of x as a function of time and the order of approximation, N+ 1; (d) error in approximation solution of y as a function of time and the order of approximation, N +1.

M. Skliar, P. Tathireddy / Computers and Chemical Engineering 26 (2002) 1013–1021

1019

[0] Fig. 2. (a) Exact solution, T(z, t); (b) errors in the first order approximation, T(z, t)−[x [0] s (z, t) +x f (z, t)].

(T(t,z) Q(t,z) = A(t)T(t,z)+ (t zc

Á 2( 2x [0] s [0] [0] Ãh (z 2 +sin(3yz)=0, x s (0,t%)=x s (1,t%)=0 Ã [0] 2 [0] Í(x f (z,~) =h 2( x f (z,~), x [0](z,0)= sin(yz)− x [0](z,0) f s (z 2 Ã (~ Ãx [0](0,~)=x [0](1,~)=0 f Ä f

A(t)= \1 (42)

\1(t)=

(2 (2 + \ 2 (z 21 (z 22

u1(t) , zc

\2(t)=

(47)

(48) u2(t) zc

(49)

Á 2( x (z,t%) (x (z,t%) [1] = , x [1] Ãh s (0,t%) = x s (1,t%) = 0 2 (z (t% Ã [1] ( 2x [1] Ã(x f (z,~) f (z,~) 2 =h Í (~ (z 2 Ã [1] [1] [0] [0] Ãx f (z,0) = sin(yz)− x s (z,0) −x s (z,0) −x f (z,0) [1] Ãx [1] f (0,~)=x f (1,~) = 0 Ä (43)

defined on the interval t[0,tf], with boundary conditions at infinity. Here Q(t, z) describes the heat generation, c is the heat capacity and z is the density. Singular forcing converts this problem to the following singularly perturbed problem defined on the interval t% [0,mtf]:

Fig. 2 shows an excellent approximation with the 1st [0] order terms, x [0] given by the following s +x f , expressions:

(6(t%,z) m =A(t%)6(t%,z)+ Q(t%/m,z) (t%

2

[1] s

[0] s

[0] x [0] s (z,t)=x s (z)=

1 sin(3yz) (3y)2 2

− (nyh) ~ x [0] sin(nyz), f (z,~)= % Cn e

~= t

&

1

Cn = 2

6(t%,z)= 0

as z1,z2 “ 9

(50)

as z1,z2 “ 9

6(0,z)=0,

(45)

where m is the perturbation parameter. According to our results, the relationship between the original and singularly perturbed problems is given by equation T(t,z)= 6(t%,z)

ƒ(z)sin(nyz)dx,

0

with ƒ(z)= T(z,0) =sin(yz)

T(t,z)= 0

(44)

n=1

where

T(0,z)=0,

(46)

Note that analytical solutions of Eqs. (42) and (43) can be found easily for any source function F(z, t).

The truncated asymptotic solution of the singularly perturbed problem is given as [0] [1] [1] 6(t%,z)= x [0] s (t%,z)+x f (~,z)+ mx s (t%,z)+mx f (~,z)

+ O(m 2)

4.2.2. Time-6arying operator A Consider the following heat transfer problem with time varying coefficients:

when t% =mt

(51)

where ~= t%/m = t and terms in the approximation (51) are found from the following simplified problems:

1020

M. Skliar, P. Tathireddy / Computers and Chemical Engineering 26 (2002) 1013–1021

accuracy. For time varying operators, for which analytical solution is usually unavailable, we developed the procedure for piece-wise approximation using singular forcing, which may yield an analytical solution. The developed method can be used to find the hybrid asymptotic-numerical solution of evolutional problems. A numerical approach, which starts with the conversion of the original problem into simplified subproblems and the subsequent numerical solution of subproblems can be classified as a splitting method (Marchuk, 1988).

Á ( 2x [0] ( 2x [0] Q(t%/m,z) s (t%,z) s (t%,z) + \2(t%) + =0 Ã\1(t%) 2 2 (z 1 (z 2 zc à Ãb.c.:x [0] as z1,z2 “ 9 s (t%,z)= 0 à [0] 2 [0] (x (~,z) ( x ( 2x [0] f (~,z) f (~,z) Í f =\1(0) +\ (0) 2 2 à (~ (z 1 (z 22 à [0] [0] Ãi.c.:x f (0,z)=x0(z) − x s (0,z) Ãb.c.:x [0] as z1,z2 “ 9 f (~,z)=0 Ä (52)

Á ( 2x [1] ( 2x [1] (x [0] s (t%,z) s (t%,z) s (t%,z) + \2(t%) = Ã\1(t%) 2 2 (z (z (t% 1 2 à [1] Ãb.c.:x s (t%,z)=0 as z1,z2 “ 9 à [1] 2 [1] 2 [1] 2 [0] 2 [0] Í(x f (~,z) = \ (0)( x f (~,z) +\ (0)( x f (~,z) +~ ( \ (0)( x f (~,z) + ~ ( \ (0)( x f (~,z) 1 2 1 2 2 2 2 à (~ (z 1 (z 2 (z 1 (z 22 (~ (~ à [1] [0] [0] [1] Ãi.c.:x f (0,z)=x0(z) − x s (0,z) −x f (0,z) −x s (0,z) Ãb.c.:x [1] as z1,z2 “ 9 f (~,z)=0 Ä



The solution of the original problem is now approximated as 2 [0] [1] [1] T(t,z)=T [0] s + T f + mT s +mT f +O(m )

(54)

[0] [1] [1] [0] [0] where T [0] s = x s (mt,z), T s =x s (mt,z), T f =x f (t, z) [1] [1] and T f = x f (t, z). As previously indicated, the utilization of the Taylor series approximation of the operator A(t) implies that the accuracy of the approximate solution (51) will deteriorate with time. To mitigate this problem, the interval of the definition can be discretized with the step Dt, small enough to give the desired accuracy in approximation (29) which, with a two-term asymptotic solution, will be of order O(Dt 2). The overall asymptotic solution on the interval [0, tf] can now be found as a combination of asymptotic solutions sequentially found for each subinterval [tm,tm +Dt], m= 0,M, tf = MDt.

5. Conclusions For linear evolutional systems, we have shown that singular perturbation methods can be used to approximate the solution of unperturbed systems, and have developed the systematic method of singular forcing for obtaining such an approximation. The approximate solution of the original inhomogeneous evolutional problem is found from a sequence of simpler subproblems. These subproblems are either of equilibrium or evolutional type, analytical solution for which may be found more easily than for the original problem. Examples with ODE and PDE models show that even first order singular forcing approximation may capture the behavior of the original problem with a high degree of

n



n

(53)

Since subproblems defining terms of an approximate solution are only mildly coupled, the hybrid asymptotic-numerical approach to the solution of evolutional problems is well-suited for computational parallelization. The singular forcing method can also be used to convert a complex, legacy equilibrium solver (i.e. a code that solves problems described by the first equation of system (32)) into a solver of unsteady-state problems. To this end, the solution of the second equation in (32) provides a boundary (in space and time) correction of the equilibrium solution. Equilibrium solution plus correction now forms an approximation of the unsteadystate solution. The accuracy of the approximation can be improved by finding additional equilibrium and correction terms described by system (34). The developed method has potential implications for controller design. For instance, the singular forcing method suggests that one approach to design a controller for PDE systems is to first obtain a feedback controller, which modifies operator A to achieve the desired behavior of the homogeneous dynamic system governed by the second equation of system (32), followed by the design of the feedforward controller to provide the desired equilibrium behavior described by of the first equation in the same system.

References Eckhause, W. (1979). Asymptotic analysis of singular perturbations. Amsterdam: North-Holland. Hinch, E. J. (1991). Perturbation methods. Cambridge: Cambridge University Press.

M. Skliar, P. Tathireddy / Computers and Chemical Engineering 26 (2002) 1013–1021 Holems, M. H. (1999). The method of multiple scales. Proceeding of Symposia in Applied Mathematics, 56, 23 –46. Il’in, A. M. (1992). Matching of assymptotic expensions of solutions of boundary 6alue problems. Providence: American Mathematical Society. Kevorkian, J., & Cole, J. D. (1981). Perturbation methods in applied mathematics. New York: Springer –Verlag. Kevorkian, J., & Cole, J. D. (1996). Multiple scale and singular perturbation methods. New York: Springer –Verlag. Kokotovic, P. V., Khalil, H. K., & O’Reilly, J. (1986). Singular perturbations in control: analysis and design. London: Academic Press. Krein, S. G. (1967). Linear differential equation in Banach space. Moscow: Nauka (in Russian). Kumar, A., Christofides, P. D., & Daoutidis, P. (1998). Singular perturbation modeling of nonlinear processes with nonexplicit time-scale multiplicity. Chemical Engineering Science, 53, 1491 – 1504. Marchuk, G. I. (1988). Splitting methods. Moscow: Nauka (in Russian). Marino, R., & Kokotovich, P. V. (1988). A geometric approach to

1021

nonlinear singularly perturbed control systems. Automatica, 24, 31 – 41. Nayfeh, A. H. (1973). Perturbation methods. New York: Wiley. O’Malley, R. E. Jr (1991). Singular perturbation methods for ordinary differential equations. New York: Springer – Verlag. Trenogin, V. A. (1970). Development and application of the Lyusternik – Vishik asymptotic method. Uspekhi Matematicheskikh Nauk, 25 (4:154), 123 – 156. Vasil´eva, A.B., Butuzov, V.F., & Kalachev, L.V. (1995). The boundary function method for singular perturbation problems. SIAM Studies in Applied Mathematics, vol. 14. Vishik, M. I., & Lyusternik, L. A. (1957). Regular degeneration and the boundary layer for linear differential equations with a small parameter. Uspekhi Matematicheskikh Nauk, 12 (5:77), 3 –122. Zgurovskiy, M. Z., & Skliar, M. V. (1991a). Multiscale expansion of the solution of a singularly perturbed systems of evolutional operator equation. Part 1. Journal of Automation and Information Sciences, 24 (5), 61 – 71. Zgurovskiy, M. Z., & Skliar, M. V. (1991b). Multiscale expansion of the solution of a singularly perturbed systems of evolutional operator equation. Part 2. Journal of Automation and Information Sciences, 24 (6), 18 – 29.