A perturbation method for the solution of linear matrix differential equations

A perturbation method for the solution of linear matrix differential equations

Journal of The Franklin Institute DEVOTED TO SCIENCE AND THE MECHANIC ARTS Volume 283, Numhw 5 May 1967 A Perturbation Method for the Solution...

737KB Sizes 90 Downloads 103 Views

Journal of The Franklin Institute DEVOTED

TO SCIENCE

AND

THE

MECHANIC

ARTS

Volume 283, Numhw 5

May 1967

A Perturbation Method for the Solution of Linear Matrix DiJkrential Equations by

LOUIS A. PIl?ES

Department of Engineering, University of California, The Aerospace Corporation

Los Angeles, and Consultant,

ABSTRACT: The mathematicalformulation of problems in electric circuit theory, structura engineering, and dynamics by the “state-vector” approach leads to a matrix di$eTential equation of the form dB/dx = A(x)S(x) + C+(X)>.This paper presents a useful method for the approximate

solution

of the matrix

di,erential

equation

boundary conditions. The general procedure is illustrated uniform electrical transmission lines of a general type.

I.

subject

to initial

by applying

and two-point

it to the studg of nm-

Introd7.4ction,

The matrix differential equation @ = A(x)X(s) da

+ 4(x)

I

X(0) = &I

(1)

where X(z) is a coIumn matrix or vector of the nth order, A (2) is a square matrix of the nth order and $(x) is a column matrix or vector of the nth order with the initial condition, S( 0) = X0 arises in the mathematica1 study of a diversified group of physical problems. For example, equations of the type of Eq. 1 arise in the study of nonuniform structural problems, dynamic systems with variable elements, time-varying electrical networks, non-uniform transmission lines and in many other problems in physics and engineering (1-5). When the elements of the matrix A (x) are functions of z of a general type, the solution of Eq. 1 cannot be effected by analytical procedures and the investigator must then resort to methods by which approximate solutions of Eq. 1

357

Louis A. Pipes

may be obtained. In this discussion we present certain procedures of the perturbation type which are useful in obtaining approximate solutions of equations of the form Eq. 1. In the physical problems enumerated above, the elements of the matrices A (5) and Q (2) are known functions of x and, it is desired to determine the vector S(z) by solving Eq. 1 subject to the given indiuE vector 80. II,

P(x)

The Transfer Matrix

The solution of Eq. I may be effected by assuming that the vector S(z) may be represented as the sum of two vectors in the form,

(2)

S(s) = SC(X)+ S,(x) with the initial conditions, S,(O) = X0,

and

S,(O)

= 0.

(3)

The vector X,(x) represents the complementary function of the differential Eq. I and the vector S,(z) may be regarded as the particular integral of Eq. 1. If Eq. 2 is now substituted into Eq, 1, we obtain the following:

dS zc

=

A(z)&(s),

(4)

and

ZP=A(~)+%(x) + 4bJ). In order to solve Eq. 4 we introduce the transfer matrix P(z),

S,(x)

= P(x)&,

P(0)

= I,

P(X)

and write,

= matrix of order (n, n)

(6)

with the stipulation that P(0) = I, the nth order unit matrix. If Eq. 6 is substituted into Eq. 4, the result may be written in the form,

(7) Since the elements of the initial vector SOare arbitrary, Eq. 7 can only be satisfied if dP = A(ci~)P(x), dx The square nth order tram&r matrix P(z) condition P(0) = I.

358

P(0)

= I.

thus satisfies Eq. 8 with the initial

_Toumalof The Franklin Institute

If it is possible to solve Eq. 8 for the transfer matrti P(z), by the following procedure: Let

we obtain X,(z)

S,(x) = P(X)&(~)

(9)

where Q(X) is a column matrix to be determined+ If Eq. 9 is now substituted into Eq. 5, the result is,

OF

The first term of the left member of Eq. II vanishes as a consequence of Eq. 8, and we have

P(x)

2 = &cc).

(12)

By now integrating Eq. 12, we obtain

Q(x) = [= P”(u)+(u)

du.

(13)

If Eq. 13 is now substituted into Eq. 9 the result is

x,(z)

= P(x)

/EWt+$44 du, 0

S,(O)

= 0.

(14)

The vector S(z) is now obtained by substituting Eq. 6 and Eq. 14 into Eq. 2, and thus obtaining X(z)

= P(X) [Ho + [=F(u)~~(u) 0

du].

(15)

This is the formal solution of Eq. 1 with X(O) = X0. Note that once the truns@‘ermat& P(z) is obtained, S(z) is found by the use of Eq. 15. Thus it is apparent that the general solution of Eq. 1 for the vector S(z) depends entirely on being able to obtain the lrunsfer mat& P(z) by solving Eq. 8 with the initial condition P(0) = I. III.

Problems

Involving

Two-Point

Boundary

Conditions

In certain problems of practical importance, the initial vector 80 is not known but information regarding X(z) is given at x = 0 and at 17;= L of the form:

MX(0)

Vol. 283,No. 5, May 1967

+ ivS(I.4) = Q

(16) 359

where M and N are square matrices of order n and G is a column matrix or vector of order n. If we place x = L in Eq. 15 we obtain the vector = P(L)[So

X(L)

w =

+ W] where W is the vector,

(17)

JLP-yu)rb(u) ah.

(18)

0

Therefore, if X(z) as given by Eq. 15 is now substituted boundary condition Eq. 16, the result is MSo + NP(L)[So Equation result :

into the two-point

+ W] = G.

(19)

19 may now be solved for the unknown

SO with the following

So = [M + NP(L)]-l[G

- &Y=(L)W].

(20)

With this for SO, Eq. 15 now gives the solution in the range 0 5 2 2 L. IV.

Computation

of the Transfer Matrix

by a Perturbation

Method

As discussed in Section II, in order to obtain the general solution of the matrix differential Eq. 1, first it is necessary to solve the equation, dP = A(x)P(z)

(211

dX

subject to the initial condition P(O) = I to obtain the &an& mat&x P(x). In the usual application to physical problems, Eq. 21 must be solved in the range 0 5 z $ L. An effective perturbation procedure for solving Eq. 21 consists of first obtaining an nth order square matrix of constants Ao, whose elements are the averages of the elements of the matrix A (x) in the range 0 5 x 6 I;, ao that

Having obtained AO by the use of Eq. 22, we now write A (x) as follows: A(s)

= AO + [A(s)

- AD] = AO + a(x)

(23)

where a(z) is the square matrix of the nth order given by a(x) = [A(z)

360

- AoJ.

(24)

Journ;tl of The Franklin Institute

The matrix a(x) is, therefore, a sort of perturbutiolz matrix about the constant matrix Ao. If the elements of the original matrix A (z) are nearly constant or vary slowly in the range 0 $ z 6 L, then the elements of the matrix a(x) are small relative to the corresponding elements of Ao. If Eq. 23 is substituted into Eq. 21, the equation satisfied by the transfer matrix P(X) is,

g

=

P(0) = I.

A0Pb) + ab)P(x),

(25)

The matrix differential Eq. 25 with the initial condition P(0) transformed (2) into the following integral equation:

[I + Is exp 0

P(x) = exp (AOX) The function exp (&x)

is the m&k exp (AOX) =

( -Aou)a(u)P(u)

= I may be

du],

(26)

exponential junction defined by

2

(AG)~/~c!

(27)

k=O

In order to solve the integral Eq. 26 we try a soIution of the form,

P(x)

=

2 Pk(2).

(28)

b=O

If Eq. 28 is now substituted into Eq. 26 we obtain,

2 pk(x> = exp(A&) p

k-0

+

J:exp ( -Aou)a(u)

g IS(U)

du]. (29)

Equation 29 may be satisfied if we let, (30) (31)

P,(x)

= PO(X)

J’Po( -Ub(u)Pl(u) a

Pk+l(u) = PO(X)

du

k = ~‘ P*(-u)a(u)Pk(u) du, cl

(32)

1, 2, 3,

?? ?? ??

(33)

The sequence of matrix functions P,(s) may be shown (6) to form the uniformly convergent series Eq. 28 under suitable restrictions and is, therefore, the solution of the integral Eq. 26. Vol. 283, No. 5, May 1967

361

Louis A. Pipes

FIG. 1. Non-uniform

V. Application

of the Perturbation

to ccNon-Uniform

Transmission

line.

Method Line

As an interesting and important exampIe of the method described in Section IV, consider the uniform transmission line depicted by Fig. I. The folIowing notation is used to describe the line: (3) x = a space coordinate measured along the line

V(z) 1(x) Z(z) Y(X)

= = = =

S(x)

=

the the the the

complex potential across the line conductors at z current flow at the point o of the line impedance per unit length of the line admittance per unit length of the Iine = the state vector of the line

[II 0 A(x)

-Y(z) P(x)

--Z(x)

=[I

I/ = the “A” matrix of the line

0

= the transfer matrix of the line

The state vector X (5) of the non-uniform line satisfies the following differential equation: (3)

2 = A(zi~)X(z).

(34)

The initial state vector X(0) = SO is now the vector whose elements are the potential Y0 and current 10 at the end z = 0 of the line. The state vector of the line X(o) is related to the trunsfer matrk~ of the line P(z) by the equation, S(x) The transfer matrix P(z)

(35)

= P(x)Xo.

of the line, satisfies the equation dP = A(x)P(x)

(36)

dX

362

Journal

of The Franklin

Institute

Equation 36 will now be solved by the method described in Section XV. We first compute the matrix A0 by the equation,

Ao = i /‘A(s) L 0

0

-zJ

[ -Y*

0

d;z: =

I

(37)

where 20 and YO are the merage values of Z(z) and Y(s) over the length L of the line. We now determine the a(z) matrix by the following equation:

a(z)

= CA(z)

x(x)

= Z(a)

-AI,

= [_;CXi

(38)

-yX’]

where - ZO,

&J(z) = Y(x)

by Eq. 26, the integral equation satisfied by P(a)

-

Yo

(39)

is

P(s) = exp L&a) I + /= exp ( -Aou)a(u)P(u)

1

du .

0

(40)

In this special case we have (3), cash (OX) PO(z) = exp (AG)

-&

i - sinh ( f?z)/Ko

where &

sinh (OX) (41)

= cash (8~)

0 = (ZOYO)II2 = a propagation function = (ZO/YO)~‘~= a characteristic impedance.

1 (42) (43)

The matrix PO(Z) is the transfer matrix of a uniform line whose impedance and admittance per unit length are ZO and Yo, the average values of Z(z) and Y(z) over the length of the line. The matrix P,(r) given by Eq. 31 is the first-order correction matrix which must be added to PO(z) to take into account the variability of the line parameters. If the parameters Z(z) and Y Cx) vary slow@ throughout the length of the line, it us usually sufficient to take

P(s) G PO(Z) t- PI(S) as a good approximation to the true transfer matrix P(s) the first-order correction matrix is P,(x)

= PO(X) I’&(

-u)a(u)Po(u)

(44) of the line. By Eq. 31,

du.

(45)

0

The matrix u(u) is obtained by replacing the letter 2 by u in Eq. 38 and is,

Vot.283.No.5,May 1967

363

Louis A. Pipes therefore, u(u)

-a’j@

= c_YOu,

(46)

I’,( -u) is obtained by replacing x by --u in Eq. 41. Accordingly, multiplication and some aIgebraic reductions, we have 2 /0

a,( -zc)a(u)Po(u)

B(x)

du = [

L-E”(~)

-KID%

C(~IlIKo

-

by direct

+ C(z)1 (47)

--B(x)

where F(a)

= 5 I= cz(u)/& 0

- &y(u)]

cash (20~) du

(48)

B(X)

= i 1% [X(U)/& 0

-

sinh (28u) du

(49)

du.

(50)

C(x) = f

P$Yzc)]

J’C4u)/& + &y(u)] 0

The first-order correction matrix PI(S) by Eq, 45 is cash (0x) PI(S)

as given

1 -&C~b) B(2) +Cb)l 1

- KO sinh (0~)

= -sinh

of the transfer matrix P(z)

cash (0~)

(0X)/& X

[F(x)

-

-B(x)

C(x>]/~~o



(51)

If the two matrices of the right member of Eq. 51 are multiplied, the correction matrix PI(Z) is seen to be (52)

Plbl = [;;;;;; ;;::::I where p,,(z)

= B(x)

cash (&IT) + [C(z) - P(X)] sinh (03) sinh (h) - [F(s) + C(z)] co& (0~) JIG ~12b) = Pt C(X)] cash (0~) - B(x) sinh (0~) }&-I P&> = f V(z) ~~~(2) = [F(X) + C(x)]sinh (k> - B(s) cash (es). 2 )

If the variation of the elements of the matrix A(z) is small, then we have I 4~)

364

I <<

20

and

I y(x)

I <<

20

(53) (54) (55) (56)

in the range 0 5 5 S L

for 0 5 x 5 L.

Jaurnal of The Franklin

(57)

Institute

Linear Matrix Di$erential Equations

In such a case the transfer makix P(z) a sufficient degree of accuracy.

is given by P(z)

= PO(X) + PI(S) to

VI. Use of the Transfer Parameter Next we discuss an alternative method which is useful in the solution of equations of the type: (58)

Equation 58 describes the distribution of potential and current in a non-uniform electrical transmission line (3). We now define B(z), the variable propagation function of the line by the equation,

B(x) = [Z(x) Y (x)-p.

(59)

The transfer parameter of the line #(CC) is defined by

95(x) =

/‘S(u) ch.

(60)

D

If Eq. SO is differentiated with respect to X, we obtain

d+ - = B(x). da

Using the relation

d & -d = _*_ - e(x) dx

d$ dx

5 &

(62)

we write Eq. 58 in the form,

-Zdr)] [;I or 0 -Y(x) We now introduce the notation:

-Z(X)/O(X) /O(x)

0

I1v

I[

(63)

. (64)

(65) (66) K ($) is the variable character&k

Vol.283, No. 5, May

1967

impedance of the non-uniform transmission

365

Louis A. Pipes

line. With this notation Eq, 64 may be written as (67) If we Iet and

X=

0

-K

i -l/R

0

B-

then Eq. 67 takes the compact form

1

(68)

(69)

-b=BS.

&

In order to simplify the form of Eq. 69, we introduce the following Iinear transformation :

If we now let

[I

j-p/Z 1 0

V

2:=

(70)

and

M=

i

[

0

K--l/a

(71)

the linear transformation Eq. 70 may be expressed in the form

s =

Alz.

(72)

Now by differentiating Eq. 72 with respect to the transfer parameter 4, and indicating differentiation with respect to 4 by a prime, we obtain S’ =

M'z +

Md.

(73)

Equation 69 may now be written in the form (74)

,.S = ES = BMz. If X’ as given by Eq. 73 is now substituted into Eq. 74, there results M’z + Mz’ = BMz.

(75)

We now solve Eq. 75 for .z’, and obtain 2’ =

366

[i,k13M

-

(76)

JP&-‘-J.

Joumal of The

Franklin

Institute

0

1t can easily be shown that

-1

= [ -1

M-IBM

and

0

1

(77)

(78)

If we now substitute Eqs. 77 and 78 into Eq. 76, we have -K’/2K

-1

er

In many cases of practical importance, it happens that the characteristic impedance K varies slowly throughout the range 0 5 x 6 L. In such cases, we have K’/2K 1 < 1. - (30) If we now define Co and C(9) by the equations,

/I 1 0

co=

-Kf/2K

-1

-1

0

ct+>

0

= 0

K’/2K

1 (811

we may then write Eq. 79 in the form,

$

= coz + C(l$)z.

(82)

Equation 82 is similar to Eq. 25 in form, therefore, it can be transformed into the integral equation, (33) where z(0) is the initial vector z at up= 0. In this case, we have (3)

co& (8) exp (Cd)

-sinh where u

Vol. 283, No. S, May

1967

-sinh

(4)

=

and

(4)

cash (#)

1

(34)

0(u) = [Z(U)Y(U)]~‘~.

367

Louis A. Pipes

The initial vector x(0) may be obtained by inverting Eq. 72 and placing 4 = 0, so that K-l’Z(O) 0 x(0) = .X(O) = M-l(O) S(0). (85) 0 Kl@ (0) i 1 The first term of the right member of integral Eq. 83 is

a(d

= exp (C04)2(0) = exp (GdJ-F(O)X(O).

(86)

Let the state vector S(4) that corresponds to the approximate solution x0($) be denoted by SO(4). We then have, as a consequence of Eqs. 72 and 86,

X0($> = IV($) exp (CO+>M-~(O)X(O). If the indicated matrix multiplications may be written as follows:

in Eq. 87 are performed,

(87) the result

L-[K(z)K(O) J-“z L

[K(z)/K(O)

J cash 9

-[K(x)K(O)]~/~

sinh 4

[K (0) /K(z)

sinh 4

-Jcash 0

(@I

where e(2) = C-z(z) U(z)-J1’2,

K(z)

= Z(z)/a(s).

The result, Eq. 88, is the W.K.B.J.L. (Wentsel-Kramers-Brillouin-JeffriesLiouville) approximation (7) . Note that the W.K.B.J.L. approximation reduces to cash (es) (89) -sinh

(6$)/K

in the simple cam of a uniform line in which the functions Z(X) and Y(z) become constants. The W.K.B.J.L. approximation is obtained by using onIy the first term of the right member of integrd Eq. 83. This integral equation has a solution of the form (90)

a(4) = &l(4) T&=0 where

(91)

.%i(+) = exp (cbQ)40) &I($) = exp (Cd) J’ exp (-Cd~)C(u)x~~(u) 0 368

du

n = 1, 2,3, ‘*a.

(92)

Journal of The Franklin Institute

Linear &la&ix DifferentialEquations The state-vector of the line rS, is given by

x = M(4)z($).

(93)

It can therefore be seen that the term ZO($) in Eq. 90 leads to the W.K.B.J.L. approximation. The terms zn($), n = 1, 2, 3, 4, lead to more accurate results than those obtained by the use of the W.K.B.J.L. method. The vector ZI($) is given by Eq. 92, in the form ?? ??

Q(4)

=

exp (C&) /“” exp (-C~U)C(U)ZO(U) 0

clu.

(94)

In cases where the chtiacteristic impedance K(X) varies slowly throughout the range 0 5 z 5 L, suflkient accuracy is usually obtained by the improved W.K.B.J.L. approximation, given by 8 = M(d)ko(&) VII.

Transmission

Non-Uai$wm

a(4>1.

+

Line with Constant

w Propagation

Function

Many non-uniform transmission lines of importance in practice have the property that their propagation function O(X) is a constant, so that 13(z) = cZ(z) Y(z)]~‘~

= 0

(a constant).

(96)

In these special cases, the lrccnsferparameter 4 discussed in Section VI is given by l#l = Jku 0

= ex.

(97)

In this case, it is convenient to write Eq. 58 in the following form: 0

-Z(z)/0

-Y(z)/0 where K(z)

0

0

2 0 I

/I -l/K(s)

--K(z) 0

I[1(981 Fi

I

is the variable characteristic impedance of the line given by K(x)

=

Ei$ = A-_ Y(x) *

In the notation of Eq. 68, Eq. 98 may be rewritten as

Vol. 283,No. 5, May 1967

369

Louis A. Pipes

of we again introduce the transformation Eq. 71, S = Mz, it is easy to show that &‘q. 100 is transformed into the equation,

By direct matrix multiplication,

0 (Ill-lSM>

= [ “,

we obtain

(ilk1%)

;I,

= i”‘?

__KU12K].

(102)

Therefore, Eq. 101 may be written in the form,

(103)

If we now let Nj”,

;],

C(X) =[-K;

K,;2J.

(104)

With this notation, Eq. 101 may be written as dx - = Nz + C(X)& dx

(105)

The integral equation form of Eq. 105 is, z=

exp (N%I /z(O) + J’

exp

(-Nu)c(u)z(u)

du

We have (5) cash (0X) exp(ivcc)

The solution of Eq. 106 is,

-sinh

(kc)

=

L -sinh

(106)

. I

0

(8~)

cash (0~)

I-

(107)

z(2) = 2 G.(x) 9-o

(108)

with .20(z) = exp (iv~)x(O)

(109)

and

znb) = exp (Nz) J

r* ~XP (-Nu)C(u)z,_~(u)

du,

n = 1, 2, 3, ...

(110)

0

370

~ournsl of The

Franklin

Institute

~~(2)

is given by the equation,

a(x) = exp (Nx) [‘exp 0

(-Nu)C(u)

(111)

exp (Nzc) &.X(O).

If we let P(s)

= K’(z)/2K(s)

exp (-Nu)C(u)

F(zc). 1 (112)

we have, - cash (28~)

sinh (20~)

- sinh (2824)

cash (28~)

(113)

exp (iVU) =

We now let,

(114) (1151 We may then write Eq. 111 in the expanded

form,

-= x(O)

(&)J(z)

- cash (~x)H(z)]

[cash (&)J(z)

[sinh

&)H(x)

- cash (Bz)d(x)]

[H(X) cash (8~) - sinh @XV(X)]

If 1 P(X) j = 1 K’(x)/~K(x) solution of Eq. 106 is

/ << 1, in the range 8(X) = z,(2)

The solution

+ G(X).

- sinh (&)A(x)]

* 1 (116)

[sinh a(x)

0 2 x $ L, a satisfactory (117)

of Eq. 98 is then given by,

where x(0)

= M-LAY(O).

(119)

Bibliography (1) L. A. Zadeh and C. A. Desoer, ‘
The State Vector. Approach,” New York, McGraw-Hill Book Co., 1983. P. M. DeRusso, R. J. Roy, and C. M. Close, “St& Variables for Engineers,” New York, John Wiley and Sons, Inc., 1965. L. A. Pipes, “Direct Computation of Transmission Matrices of Electrical Transmiesioo Lines,” Parts I and II, Jou~. Frank. Inst., Vol. 281, Nos. 4 and 5, April and May 1966. B. I<. Kinariwals, I.R.E. Convention Record, Pt. 4, pp. 268-276, 1961. L. A. Pipes, “Matrix Methods for Engineers,” Englewood Cliffs, N. J., Prentice-Hall, Inc., 1963. F. G. Trieomi, “Integral Equations,” New York, Interscience Publishers, 1857. J. C. Water, “Microwave Transmission,” New York, McGraw-Hill Book Co., 1942.

Vol. 283,No. 5, May 1967

371