A Method for Solving Systems of First Order Linear Homogeneous Differential Equations when the Elements of the Forcing Vector are Modelled as Step Functions

International Finance Discussion Papers Number 275

February 1986

A METHOD FOR SOLVING SYSTEMS OF FIRST ORDER LINEAR HOMOGENEOUS DIFFERENTIAL EQUATIONS WHEN THE ELEMENTS OF THE FORCING VECTOR ARE MODELLED AS STEP FUNCTIONS

by

Robert A. Johnson

NOTE: International Finance Discussion Papers are preliminary materials circulated to stimulate discussion and critical comment. References in publications to International Finance Discussion Papers (other than an acknowledgment by a writer that he has had access to unpublished material) should be cleared with the author or authors.

Abstract

This paper presents a method for solving a system of first order linear differential equations with constant coefficients when the elements of the forcing vector are step functions. The analysis presented in the text has been programmed for use in the computer simulation of linear continous time rational expectations models using any combination of anticipated and unanticipated, permanent or temporary shocks. The program

entitled "JAB" is available from the author upon request.

A Method for Solving Systems of First Order Linear Homogeneous Differential Equations When the Elements of the Forcing Vector are Modelled as Step Functions by

Robert A. Johnson*

This paper presents analytical expressions for the solution of a system of first order linear differential equations with constant coefficients when the elements of the vector of exogenous forcing variables are restricted to being modelled as step functions. The solution of such a system can then be programmed quite easily for use on the computer. As a result, simulation of deterministic continuous time rational expectations models can be performed under a wide variety of circumstances including any combination of anticipated or unanticipated shocks, permanent or transitory shocks and future or present shocks.

The first section below develops the notation and presents the discussion of the basic solution of first order linear differential equation systems with constant coefficients. This section constitutes the background for the focus of this paper and is essentially a review of the material presented in Buiter (1984). Section three begins with the general solution form developed in the previous section and derives the specific solution when the elements of the foreing vector are restricted to step functions. A brief conclusion is followed by an appendix exhibiting the computer program which incorporates the analysis

contained in the body of the paper.

* The author is grateful to seminar participants at the Board of Governors of the Federal Reserve and to William Buiter for helpful comments. I would also like to thank Ruby Brooks for her patience and skillful typing of this paper. This paper represents the views of the author and should not be interpreted as reflecting the views

of the Board of Governors of the Federal Reserve System or other members of its staff.

-2-

SECTION II. Background In this section we will develop the notation and the basic solution method that will be used below.

We begin with the linear first order differential system,

x(t) x(t)

where: x(t) is an Ny vector of predetermined variables y(t) is an N-N, vector of nonpredeterminted variables z(t) is a k vector of exogenous variables A is a constant matrix of dimension N¥*N B is a constant matrix of dimension N*K E denotes the expectations operator

wu denotes the derivative with respect to time d/dt I will assume that the A matrix can be diagonalized as follows.

(2) A-z=v7lav or vav7l =a.

where: V is the matrix of left eigenvectors of A of dimension N¥N. A is a diagonal matrix whose elements are the

eigenvalues of matrix A. I will further assume that there are N, eigenvalues with nonpositive real components and N-N, eigenvalues with positive real components. Following the exposition of Buiter (1984) I will partition the A,E,V,

v-' and A matrices to conform to the dimensions of x(t) and y(t). [ a ae | o-[ ] ;

A21 A22 Bo

Yiq Vi2 vile{ War Way yl] A, 0

Vo1 Vo Wo1 Wo O A °

(3) A

where Aij» Vy1, and Wy, are Ny * N; matrices. Ai2, Vis and Wi2 are Ny * (N-N,) matrices. Aa» V21 and Wo, are (N-N,) * Ny, matrices. A22, Va2 and Woo are (N-N,) * (N-N,) matrices. Aj, is an Ny * Ny matrix containing the stable roots of A.

Ag is an (N-N1) * (N-Ny) matrix containing the unstable roots

of A. By is an Nj * K constant matrix.

B2 is an (N-N,) * K constant matrix.

Let variables p and q be defined as follows. 4 Py =v (*) or (%) = v-1 (P), (4) (F) (Y? CY a

p is an Ny vector and q and N-Nj vector.

Substituting expression (2) for A, and taking the expectation

of (1) yields,

Ey x(t) _y-] Ez, x(t) (5) E xe) | =V lav Ex y(t) | * BE, z(t). Premultiplying both sides of equation (5) by V using the

definitions presented in (4) and partitioning the A, V and B matrices as

in (3) one obtains.

(6) at Pity | = 1 4 nt nea + yi yl2 a! B,2(t). t Wt) t 4 21 ‘22 2

Let (7) De V21B, + Vo2Bo.

Then CS) Fy g(t) = Ap Eeq(t) + DEg z(t).

The forward looking solution for equations (8) is = Ans _ 5? QA (s-T) > (9) E,a(s) = Ke 2 Sy e*, DE,2(t)dt, s2t.

where Ko is an N-N, vector of arbitrary constants. Provided that he expectation of z(t) is bounded over the interval [t, + »), the second term of equation (9) will exist. But the expression in (9) will not converge unless Ko=0. Assuming that the elements of Ko are all set

equal to zero (9) becomes.

(10) Epa(s) = - sg e825") « pee, a( rar.

From the definition of q(t) given in (4) above. (11) q(t) = Voz, * x(t) + Voo*y(t).

or

vit ® g(t) - Vv

(12) y(t) 29

-1 ¥ % 22" Vo, * x(t).

Evaluating (10) at s=t using the weak consistency assumption that

Erq(t) = q(t) and substituting into (12) one obtains

(143) y(t) =- vi!

-1 @ Ap(t-t) oo * Var¥x(t)- V4, S, e°? *D¥E, 2(t)dt.

22

which is equation (14) from Buiter (1984). Recall that by definition (3),

Vaw = I, Thus,

V = 0.

21 M41 * Yoo Wo,

Premultiplying by Vo5 and postmultiplying by wi yields,

-1 -1 (14) Wo, Wy = Vo5 Vo4°

Equation (13) can be expressed as

(13") y(t) = Way *W"tx(t) - VS prehee Danas, z(t)dt.

From equation (1) x(t) = Ayq*x(t) + Aypoty(t) + Byz(t).

Substitute (13) in for y(t) and one obtains .

° -1 (15) x(t )=(A) 4-Ay VooV5, x(t) + Byz(t)

-1 i ela (t-1)

“Ay oVo50 j DE,z(t)dt.

The compound term which multiplies x(t) can be simplified as follows.

From (2), A = WAV. (where W=V71) Aqq=WyqAqVqq + Wy2haVay.

Aya=W11AqVi2 + Wy2AoV20

-1 -1 Avi~Ay2aV22 Vaq=Wy 4 AqV 41 *Wy ahaV2q ~ WyyAqVyaVoaVay - Wy 5A,V,,-

-4 (16) = Wi AV, ~ VyoVo0Voy]

Recall again that

WV = I.

tt tet .

W11Vq14Wy2Vo4

= Vat evor = Mate

and Wi1V12 + W12Vo0 = 0.

Premultiplying by wi!

"1 and post multiplying by V55 yields,

~ -4 Way Wyo = “Yao Voo-

Putting this back into (17) creates

viv. ew!

a4 7 V2 Yoo Voy = Way

and (16) becomes

1

1 - (18) A,.-A W AWW,

11 AyaVoaVay = May As a result (15) can be reexpressed as

° 71 -1,@ A5(t-1) . (19) x(t )=Wy 1 AqWy 4 x(t )#Byz(t)-AyaV oS, e DE,2(1)dt.

Given an Ny vector of initial conditions for each element of the vector

x(t) at time t=t, one can express the solution to (19) as

_ + yl t ~e\ynl (20) x(t) = WizeAq(t to)W) x(t. )+ Feoli1edt(t SW, , B,2z(s)ds

t © - 7 Fel nett (t-soWs 1a, v5! # sete’ pe oc ryaras. s

The expression W554, oVoo is then transformed using (2),

Ay2 = WyqAqVq2 + WrohoVoo.

Premultiplication by Ww and postmultiplication by V povields,

“1 lav vil 4 wily

(21) Wis AioVS5 = 112"22 * W4,W) 545.

and the solution for x(t) becomes

(20)' x(t) = welt (tt OW Te(t, past oP, en fe Ou tB B.z(s)ds

~ ty Ay (t-s), “Ty x *7”,42(s-T) Seo a 1e CAy 1V12V55 + WWW 1245] fe DE,2(t)dtds.

These two equations, (20) and (20') are the expressions developed by Buiter (1984).

The solution to x(t), the predetermined variable, reflects three influences. First, the initial condition X(to), the relative influence of which diminishes as t increases. The second influence, which is expressed in the Second term of (20) or (20') above, is that exerted by the elements of the forcing vector from the initial time, to until time t. Finally the third influence can be expressed as, what was believed about the future course of the forcing elements z(t) at each point in the past from to to t. Let us turn now to the evolution of the Solution components when the elements of the forcing vector z are |

modelled as step functions.

Section III

We can now examine the components of the solution for x(t) and y(t). First the solution for x(t) is developed and then the sclution to y(t) is analyzed.

To develop a solution that is analytically tractable and that can be implemented on computer I will limit the elements of the forcing vector to change only in steps.

Z=Zo for t § to,

Z=Z1 for to < ts ty,

Z=Zo for ty < ts to,

Z=Z; for ty-, < ts ti,

Z=Z™ for th < t S$ @,

Where n= the number of columns in the Z matrix minus one.

The solution to x(t) has three components. The first, representing the n; vector of initial conditions, can be evaluated without manipulation. The second term represents the influence of the forcing vector from the initial time to, to time to the present.

t Ay(t-s),, -1 Leo W44e Way B,2(s)ds.

If the forcing vector is constant between sucessive points as assumed above then this integral becomes

J (21) J “Why eAt (Erba) Ay (ttya4)

jay + i=1

-y ph taelt (ty),

1 Jtl:

~9-

Where J is the number of transitions of the forcing vector z from time to to time t. One then has to compute the third term of the solution,

eM (t-S ) exp? h2(s- 1)

t SJ, oW to 1 s

% 1 D E2(t)dtds.

al -4 where M=A1V15V55 + Way Woh,

First of all one analyses the inside integral,

(22) £(s)= sreh2(S-D ape a(n)dt.

The relationship between s and the transition times tt; is crucial in the development of this expression. Let us call tt, the first

transition time after s. Thus our integral becomes.

f(s) = peste’ xg a(t)de

*DE | z(t)dt.

r ttstit+l Ao(s-t) +) Sissi e 1

izo When evaluated this becomes IMAX

tts i=o

(23) f(s) = a,!

2 Zetstit1 “tts+i?*

-10-

where IMAX = number of transitions of the z vector after point s.

Note that for s greater than the last transition time one obtains .

f(s) = s2 ehats-™) pe (z(t) )dt.

Now we must evaluate the entire expression

Ay (t-s)

t (24) SeoW 1° * M¥f(s)ds.

1

The difficulty arises because the function f(s) is continuous across transition times but it is nondifferentiable at the transition times.

For instance at the first transition time, tt; for e€ small

_ a7 -1 A, (s-tt,), _ f(s) ttt, A, DZ,+ A, e2 1°* D(Z,-Z,) “1 Aj (s-tt,) yore _ + A, e2 2 D(z, Z,) + oe. -1 “1 Aj (s-tt,) yoo _ (3) | ett yeg= he Zot Ap e 2 2° *D[Z,-2,] + «.. As one takes the limit as «70, f(s)_ = f(s),.

The function f(s) is continuous. But differentiating the two

expressions and evaluating as e>0,

of(s) _ of(s), _ ., _ 38 3s TD 22) tatty.

Which is only zero (i.e. the expression is differentiable) if zo = 24,

in other words, if this is not a transition point.

-ll-

Thus, the expression (24) must be broken up into a number of

integrals.

t tie] Ay(t-s), ¥ i+] (25) Br fitt wie M* £(s),"'ds,

Where f(s); is the appropriate function between te and t

i+]’

titel

and nt = the number of transitions of the z vector up to the current

time t.

Consider a representative time t which is between tt; and

ttj+1- One obtains,

sely ety (ts) ro eM‘ Sue(s)2as

1 toni Mf(s) (ds + J

wc eeeee + st W eh (tS) e(gI*1G, ty 11 J

Referring back to equation (23), we can recognize that a given

point t between ty and Uj+1 all terms of the form,

-1, Ao(s-ty.) y A, De K+1 DZ oe Kad

for JS KS » will be common to all subfunctions f(s);"! between

t) and t. One can then recombine terms and evaluate the expression t A, (t-S)y yyy (A5(s-t, nr.

D (26) Jens 12 1 M A, e2 jrit+1’Dlz

gear *y+14i

~12- Thus

peeky ef (ts)

-1_A,(s-tk) * * ty M11 M A, e2

* . dl2,, 4724]

becomes Wi1x DEZK+1 - 2x], where yx is an N; *(N-N1) matrix whose element 6, o equals,

6 Oy 6 _ Y. =@ xCehy ft ty) - eh, ft t? x eholto tK)y, 1,6,0 6,0

Now for every K from 1 to J there will be an expression of this type.

Thus the contribution of these past terms can be expressed as

qQ

(28) re Wi, DIZ... - z].

1 K

The final set of terms are derived from the elements of

f(s);7! Which are of the form AL!

2 Day and they evaluate to expressions of the fora, J -1 A,(t-t.)__ A, (t-t, ) (29) ) - Wy, bAoJ * M¥Le 1 i’-e “1 i i=o

¥ bz, ,J-

J+1 and finally between tj and t one must evaluate the component of" f(s)y

of the form A3'Dz,, ,which integrates to become

-1 -1 _ oA (t-t3) (30) W,,CA,] * CA] Mm * Ci-el NS’) # Dee,

-~13-

To solve this integral one forms a matrix @ whose element in row 6

20lumn o is equal to

_ GO _ ,8471y, x raoq7l 8 5,0" CA. Av d Meg [A.J where 136 8Ny, and 1 S$ o S N-Ny-

The Ny * (N-N,) matrix © is then used to formulate another Ny * (N-N1) matrix By, where the subscript i denotes the number of forcing element transitions beyond the first transition occuring after the current time

.

0 ) 0 B, . = 6, * [ehatt tt i ed -ehy Ct tos ehalty tt oeied]. 1,6,0 6,6 “he solution to (26) then becomes NZ-nt-1 ¥ - (27) Hy eZ od 7 Byer ei! L=0 4 similar technique is employed for all terms f(s)! from 0 to K

i where K in turn runs from 0 to J, the number of the transition time

preceding the current time t

For all f(s);"", i=0,K the term

-1 JA (s-t e

) _ Ay. > k D{z,,, Zz

K!

will be a common element.

~14-- The solution for x(t) in total is

(31) x(t) = W, eM EEO WT CE) ~ £(27) + (28) + (29) +

1

+ (30)] + (21).

One can then use x(t) in equation (13') and that leaves only the second term in that expression to evaluate

' _ -1 cul p@4 (13") y(t) = W,*W,, * x(t)-V,,* Se 2

(t-t)y D¥E, 2(t)dt.

The second term reflects the expectation of future forcing vectors as of time t presuming once again that E,z(t) = z(t), (i.e. correct expectations), one can evaluate

i ehalt-™) « peg ea(r)dt.

At a time t between ty and ty,, one obtains,

nt 1

“Vyns ~ xfe4 (32) \, D 25447 7 A, [e 2

(t-te ebalt-ty,. dys

D¥lZysi+1 ~ Zy+il

Where nt is the number of transitions of the forcing vector between t and » . Equation (32) and the first term of equation (13'), which can be computed using the results obtained in equation (31), are

combined to solve for y(t).

-15-

Conclusion:

The analysis in the preceding sections can be readily implemented on the computer to permit one to conduct policy simulation experiments which illustrate the effects of any combination of anticipated or unanticipated, permanent or transitory changes in exogenous variables. The program called "JAB" written by Johnson, Austin and Buiter, which is the offspring of "Saddlepoint", has recently been completed and is available from the author. A copy of the program

follows as an Appendix.

-16-

References Austin, G. P., and Willem H. Buiter: "Saddlepoint" A Programme for Solving Continuous Time Linear Rational Expectations Models", University of Bristol Discussion Papaer No. 82/132, November, 1982. Buiter, Willem H., "Saddlepoint Problems in Continuous Time Rational Expectations Models: A General Method and Some Macroeconomic

Examples," Econometrica, Vol. 52 (May 84), pp. 665-80.

-17- Appendix Presented here is the program "JAB" which implements the

results of the paper. This program is an outgrowth of "Saddlepoint";

see Austin and Buiter (1982).

~18-

SSET LINEINFO

SRESET FREE LIST

SSET AUTOBIND

SBIND=FROM XNAGF/=;

FILE 19(KIND=DISK, TITLE="INFILE', FILETYPE=7)

FILE 18CKIND=DISK, TITLE="OUTFILE', PROTECTION=SAVE) FILE 99(KIND=REMOTE,MYUSE=O0UT)

FILE 98¢CKIND=REMOTE,MYUSE=IN)

THIS IS THE "JAB" PROGRAM WRITTEN BY G.P. AUSTIN, W. BUITER AND ROBERT JOHNSON. ONE MUST HAVE THE N.A.G. LIBRARY TO BE ABLE TO IMPLEMENT THIS PROGRAM.

ANNAAAANAANA|NN|NANO

IMPLICIT REAL®8(CA~-H,0-Z)

REALS E1(024,24),E2(24,24),E3(24,24),E4(24,24),E5( 24,24) 1,£6(24,24),E7(24,24),E8(24,24),F1024,24) ,F2024, 24) ,F3(24, 24), FG(24 1,24), BC24,24),2(24, 24), WKSPCE( 24), FOIDAF,X,Y,ZZ,1T,CR»CI,XEND

COMPLEX*16 C1(24,24),C2(24, 24) 1 ,C3€24,24) ,C4(24,24) ,C5(24,24) ,C6 (24,24) ,C7( 24,24) ,C8( 124,24),C9(24, 24) ,C10(24,24),C11(24, 24) ,C12(24, 24) ,C13(24 1,24),C14(24, 24) ,C15( 24,24) ,C16( 24,24) ,CC(24, 24), CX

INTEGER INTGER(24),N1,N2,N3,NZ,NSR,NUSR,NIC,IA

CHARACTER¥128 INFILE,OUTFILE

PRINT, "WHAT IS YOUR DATA FILE CALLED?

READC98,X®)INFILE

OPENC19, FILE=INFILE, FORM="FORMATTED!® )

PRINT, "WHERE SHOULD I WRITE THE OUTPUT?

READ( 98, )OUTFILE

OPEN (18,FILE=QUTFILE, FORM="FORMATTED! ) READC19,X)N1,N2,N3,NZ

TA=MAX0(N1,N2,N3,NZ)

READ SIZE OF MODEL AND PASS LARGEST DIMENSION TO ALL SUBROUTINES

Cats Meo KE}, E2,E5,E4,E5,£6,67,E8,F1,F2,F3,F4,B,Z,WKSPCE,C1,C2,C3, 1€4,C5,C6,C7, 1 C8,C9,010,C11,C12,C13,C14,C15,C16,CC,CX, INTGER,N1,N2,N3,NZ, IA) RETURN E SUBROUTINE SUBSCE1,E2,E£3,£4,E£5,£6,£7,£8,F1,F2,F3,F4,B,Z,WKSPCE,Cl, 1C2,C3,C4,C5,C6, 1 C7,C8,C9,C10,C11,C12,C13,C14,C15,C16,CC,CX, INTGER, N1,N2,N3,NZ,IA) IMPLICIT REAL®8(A-H,0-Z) REAL%8 E1(IA,IA),E2CIA,1IA),E3CIA,IA),EGCIA, IA), E5(IA,IA) 1,E6CIA,IA),E7CIA,IA), 1 E8CIA, IA), FICIA, IA), F2CIA,1A),F3CIA,IA),F4CIA,IA),BCIA,IA),Z(IA,1 . 1A),WKSPCECIA) LEME LEX#16 C1CIA,IA),C2CIA,IA),C3CIA,IA),C4CIA,IA),C5CIA, IA) ,C6C A,IA), 1 C7CIA,IA),C8CIA,IA),C9CIA,IA),C1OCIA,IA),CLICIA,IA), 1 C12(TA,1TA),C13CIA,IA),C14CIA,IA),C15(IA,IA),C16CIA, IA), 1 CCCIA,IA),CX REALX8 BLANK 7! 4 INTEGER IA,N1,N2,N3,NZ,NUSR,NIC,IFLAG,KSTR, INTGER . DIMENSION INTGERCIA) , 99995 FORMAT(1X,A6) 200 FORMAT(/* SOLN. TO MATRIX A',7) 210 FORMAT(/* SOLN. TO MATRIX B',7)

ANNANAANN

aan

220 230 240 250 260

270

FORMAT(/# FORMAT(/* FORMAT(/# FORMAT(/# FORMATC*

Lt)

FORMATC/*

280 FORMATC/*

290 300 310 320

330

FORMAT(/* FORMAT(/! FORMAT(/* FORMATC!

-19-

SOLN. TO MATRIX C',7) SOLN. TO MATRIX D',7) REAL PART OF THE EIGENVALUES OF MATRIX A',/) IMAGINARY PARTS OF THE EIGENVALUES OF MATRIX A',/) REAL PARTS OF THE EIGENVECTORS OF MATRIX A (BY COLUMN)

IMAGINARY PARTS OF THE EIGENVECTORS OF MATRIX A'/)

SOLN. TO XBAR',7/)

SOLN. TO XBARBAR',7/)

SOLN. TO YBAR',/)

SOLN. TO YBARBAR',7) , A_PROGRAM FOR SOLVING FOR THE BEHAVIOUR OF LINEAR ',/?

Ll RATIONAL EXPECTATION MODELS ON THE STABLE MANIFOLD',//)

~FORMAT(!

DATA SUPPLIED ....E1,E2,E3,£4,E5,E6,E7,E8, Ll VECTORS',///) FORMAT(1X,16(F13.6)) FORMATC1X,16(€F10.6)) DO 54 I=1,I1A CC(I,I)= DCMPLX(1.0D0,0.0D0) 4 BCI,I)=1.0

SET ELEMENTS 0

READC19,*®) (CE READ(19,*)¢ READ(19, *)¢ READ(19,®)( READ(19, ¥)¢ READ(C19, x) ¢ ¢ ¢ ¢

~ .

= ~ ww

READ(19,x) READ(19,x) READC19, x) M=1 IFAIL=0

SET ALL NAG ROUTINES TO TAKE IFAIL PARAMETER AS ZERO

WRITEC18,320) WRITEC18,99995) BLANK WRITEC18,330) WRITEC18,99995) BLANK DO 81 I=1,N1 WRITEC18,100)(E1(I,J),J=1,N1) 81 CONTINUE WRITEC18,99995) BLANK DO 82 I=1,N1 WRITEC18,100)CE2C1I,J),J=1,N1) 82 CONTINUE ™ WRITEC18,99995) BLANK DO 83 I=1,N1 WRITEC18,100)CE3(1I,J),J=1,N2) CONTINUE WRITEC18,99995) BLANK DO 84 I=1,N1 WRITEC18,100)CE4CI,J),J=1,N3) CONTINUE WRITEC18,99995) BLANK DO 85 I=1,N2 WRITEC18,100)CE5(1I,J),J=1,N1) CONTINUE

fk ne ee Se Se >

HH nun aw wn tt

¢ ¢ ¢ ¢ ¢ ¢ ¢ 4

~weweyvuvvuves qewe we ve ewe vw

. BRE eee ee

zee seh eh YY ool Le |

Sd ed te

lokeke)

83

84

85

EXOGENOUS Z

86

87

88

ANOS No}

65

-20~

WRITEC18,99995) BLANK

DO 86 I=1,N2 WRITEC18,100)CE6CI,J),J=1,N1) CONTINUE

WRITEC18,99995) BLANK

DO 87 I=1,N2 WRITEC18,100)CE7(1I,J),J=1,N2) CONTINUE

WRITEC18,99995) BLANK

DO 88 I=1,N2

WRITEC18,100) CE8(1,J),J=1,N3) CONTINUE

WRITEC13,99995) BLANK

DO 89 I=1,N3 WRITEC18,100)(Z(1,J),J=1,NZ) CONTINUE

ECHO WHAT WE READ IN AS DATA

WRITEC18,99995) BLANK

CALL FOGAEFCE7,IA,B,IA,N2,N2,F1,IA,WKSPCE, F2,1A,F3,I1A,IFAIL) IOPT=

cape FOLCKFCE1,E3,F1,IA,IA,IA,WKSPCE,IA,IOPT,IFAIL) OPT=1

CALL FOICKFC(F2,F1,E6,I1A,IA,IA,WKSPCE,IA,IOPT,IFAIL)

CALL FOICKFCE7,E3,£6,1A,1IA,IA,WKSPCE,IA, IOPT,IFAIL)

CALL FOICEFCE2,E2,E7,1A,1A,IFAIL)

CALL FOGAEFCE2,IA,B,IA,N1,N1,F3,IA,WKSPCE,F4¢,1A,E6,IA,IFAIL)

CALL FOICKFCF4,E3,E5,1A,1IA,IA,WKSPCE,IA,IOPT,IFAIL)

CALL FOICEFCE1,E1,F4¢,1A,1A,IFAIL)

CALL FOICKFCE6,E3,E8,1A,IA,IA,WKSPCE,IA,IOPT,IFAIL)

CALL FOICEFCE4,E4,E6,I1A,IA,IFAIL)

IOPT=3

CALL FOICKFC(E1,F1,E8,1A,1IA,IA,WKSPCE,IA,IOPT,IFAIL)

cat FOLCKFCE1,F1,£5,1A,1A,IA,WKSPCE, IA, 10PT, IFAIL) 0 =

CALL FOICKFCE6,F2,F3,IA,1IA,IA,WKSPCE,IA,IOPT,IFAIL)

CALL FOICKFCE7,E6,E4,1A,1A,IA,WKSPCE,IA,IOPT,IFAIL)

CALL FOICEFCE8,E7,E8,IA,IA,IFAIL)

CALL FOICKFCE2,E6,£1,1A,IA,1IA,WKSPCE,IA,IOPT,IFAIL)

CALL FOICEFCE7,E2,E5,1A,1A,IFAIL)

CALL FOICKFCES5,F3,E1,1A,1IA,IA,WKSPCE,IA,IOPT,IFAIL)

CALL FOICKFCE6,F3,E4,1A,1IA,IA,WKSPCE,IA,IOPT,IFAIL) DO 65 I=1,1A

DO 65 J=1,I1A

E5(1I,J)=-E5(1I,J)

E6(I,J)=-E6CI,J)

CONTINUE

CALL FOGAEFCE5,IA,B,IA,N1,N1,F1,IA,WKSPCE,F2,1A,F3,I1A,IFAIL) CALL FOICKFCF2,F1,E6,1A,IA,IA,WKSPCE,IA,IOPT,IFAIL)

CALL FOICKFCF3,E7,F2,1A,1IA,1IA,WKSPCE,IA,IOPT,IFAIL)

CALL FOICEFCF4,E8,F3,1A,1IA,IFAIL)

CALL FOICKFCE1,F2,Z,1A,IA,IA,WKSPCE,IA,IOPT,IFAIL) CALL FOICKFCE2,F4,Z,I1A,IA,IA,WKSPCE,IA,IOPT,IFAIL)

WRITEC18,99995) BLANK DO 536 I=1,IA

DO 636 J=1,I1A Elci,J)=-E1CI,J)

636

501

502

637

an

55

56

57

loleked®,) oo

~21-

CONTINUE WRITEC18,99995) BLANK

MANIPULATE THE E MATRICES TO GIVE SOLN. TO A,B,C,D WRITEC18,99995) BLANK WRITEC99, 501)

FORMAT(7" SOLN TO XBAR') VIRITEC99,100)CE1(CJ,1),J=1,N1) MRITEC18,99995) BLANK HIRITEC99, 502)

FORMAT(“"® SOLN TO XBARBAR') WRITEC99,100)CE1CJ,NZ),J=1,N1) WIRITEC18,280) VIRITEC18,100)CE1(J,1),J=1,N1) VIRITEC18,99995) BLANK WRITEC18,290) WIRITEC18,100)CE1CJ,NZ),J=1,N1) WRITEC18,99995) BLANK WRITEC18,300) WRITEC18,100)CE2€J,1),J=1,N2) WRITEC18,99995) BLANK WRITEC18,310) WRITEC18,100)CE2CJ,NZ),J=1,N2) WRITEC18,99995) BLANK

DO 637 I=1,I1A DO 637 J=1,I1A F1CI,J)=-E1CI,J)

CONTINUE NRITE SOLN. TO XBAR,XBARBAR, YBAR, YBARBAR WRITEC18,200)

00 55 I=1,N1 WRITEC18,100)CE5(1I,J),J=1,N1) CONTINUE

WRITEC18,99995) BLANK WRITEC18,210)

DO 56 I=1,N1 WRITEC18,100)CE6CI,J),J=1,N3) CONTINUE

WRITEC18,99995) BLANK WRITEC18,220)

NO 57 I=1,N2 WRITEC18,100)CE7(1I,J),J=1,N1) CONTINUE

WRITEC18,99995) BLANK WRITEC18,230)

DO 58 I=1,N2 WRITEC18,100)CE8(1I,J),J=1,N3) CONTINUE

WRITE A,B,C,D TO OUTPUT FILE.

DO 10 IT=1,N1

DO 10 J=1,N1

(1201, J)=DCMPLXCE5(1I,J),0.0D0)

WRITEC18,99995) BLANK

(CALL FOZAGFCE5,IA,N1,E2,E3,E4,1A,F1,1IA,INTGER,IFAIL)

TO FIND THE E-VALUES AND VECTORS OF MATRIX A

NW

ANNO He

59

8838

15

-22-

ZZ=0.0

N11=N1-1

DO 12 I=1,N11 IFCE2(I,1).LE.E2¢I+1,1)) GO To 12 2Z=2Z+1.

X=E2(1I,1)

Y=E3¢€1I,1)

E2CI,1)=E2cI+1,1) E3(I,1)=E3(I+1,1)

CONTINUE CONTINUE IF (Z2Z.GT.0)GO TO 11

SORT THE E-VALUES INTO ASCENDING ORDER AND ARRANGE THE E-VECTORS CONFORMALLY.

WRITEC18,99995) BLANK WRITEC18,240) WRITEC18,100)(E2€I,1),1=1,N1) WRITEC18,99995) BLANK WRITE(18,250) WRITEC18,100)(E3€I,1),I=1,N1) WRITEC18,99995) BLANK WRITEC18,260)

DO 59 I=1,N1 WRITE(18,100)CE4(I,J),J=1,N1) WRITEC18,99995) BLANK WRITEC18,270)

DO 888 I=1,N1 WRITEC18,100)(F1(I,J),J=1,N1)

WRITE THE E-VECTORS AND E-VALUES TO THE OUTPUT FILE. WRITEC18,99995) BLANK

NSR=0

DO 14 I=1,N1

IF CE2€I,1).LE. O)NSR=NSRt1

CONTINUE

FIND THE NUMBER OF STABLE ROOTS.

NUSR=N1-NSR

DO 15 I=1,Nl DO 15 J=1,N1

X=E4(1I,J)

Y=F1CI,J)

C1CI,J)= DCMPLX(X,Y) C201I,J)= DCMPLX(X,Y) C3(I,J)= DCMPLXC(X,Y) DO 16 I=1,NSR

16 17

19

20

21

22

23

24

25

886

887

28

26

27

-23-

C4(1I,1)}= DCMPLXCE2(1I,1),E3(1,1))

DO 17 I=1,NUSR ‘

C5(I,1)= DCMPLXCE2CI+NSR,1),E3CI+NSR,1))

CALL FOGADF(C1,1A,CC,IA,NSR,NSR,C6,1A,WKSPCE,IFAIL) CALL FOGADF(C2,1A,CC,IA,N1,N1,C1,IA,WKSPCE,IFAIL) exe DCMPLXC0.0D0,0.0D0)

N11=NSR+1

DO 19 I=N11,N1

K=I-NSR

DO 19 J=1,NSR

CALL FOIDDFC(L,NSR,CX,1I,J,C3,IA,C6,1A,N,CR,CI) E3(K,J)=-CR

DO 20 I=1,N1

DO 20 J=1,N3

C7(I,J)= DCMPLXCE6(I,J),0.0D0)

DO 21 I=1,NSR

DO 21 J=1,N3 . CALL FOIDDF(L,NSR,CX,1,J,C6,IA,C7,1A,N,CR,CI) C8(I,J)= DCMPLX(-CR,-CI)

DO 22 I=N11,N1

K=I-NSR

DO 22 J=1,N3

CALL FOLDDF(L,N1,CX,1,J,C1,IA,C7,1A,N,CR,CI) C11(K,J)= DCMPLX(-CR,-CI)

DO 23 I=N11,N1

K=I-NSR

DO 23 J=N11,N1

N=J-NSR

C7(K,N)=C1CI,J)

CALL FO4GADF(C7,1A,CC,IA,NUSR,NUSR,C9,IA,WKSPCE, IFAIL) DO 24 I=1,NSR

DO 24 J=N11,N1

K=J-NSR

CALL FOIDDFC(L,NSR,CX,1,J,C6,IA,C12,1A,N,CR,CI) C1CI,K)= DCMPLX(-CR,-CI)

DO 25 I=1,NSR

DO 25 J=1,NUSR

CALL FOIDDF(L,NUSR,CX,1,J,C1,IA,C9,IA,N,CR,CI) C7(I,J)= DCMPLX(-CR,-CI)

DO 26 I=1,NSR

DO 886 J=1,NSR

C1(I,J)= DCMPLX(0.0D0,0.0D0)

X=DREAL(C4(1,1))

Y=DIMAG(C4(I,1))

ZZ=X¥*X+YVRY

C1C1I,1I)= DCMPLX((X/ZZ),-CY/ZZ))

CONTINUE

DO 27 I=1,NUSR

DO 887 J=1,NUSR

C2(1,J)=(0.,0.)

X=DREAL(C5(I,1))

Y=DIMAG(C5(I,1))

ZZ=X¥X+YRY

C2(1,1)= DCMPLX((X/ZZ),-CY/ZZ))

CONTINUE

DO 28 I=1,NSR

DO 28 J=1,NUSR

CALL FOLDDF(L,NUSR,CX,I,J,C7,1A,C2,IA,N,CR,CI) C10€I,J)= DCMPLX(-CR,-CI)

ioleokeke)

776

399 999

998

29 99991

777 41

778 781

779 99992

53 51

885

884 997

996

994

"88881

x!

~24-

MANIPULATE THE E-VECTOR MATRICES AND THE E-VALUE MATRICES TG GIVE THE COEFFICIENTS IN EQUATIONS 12B AND 13.

DO 776 I=1,N3

DO 776 J=1,NZ

CCCI, J)= DCMPLX(ZC1I,J),0.0D0)

X=0.0

WRITEC99,399)NSR

FORMATC/* THE NUMBER OF STABLE ROOTS IS ',12) WRITEC99,999)

FORMAT(/* HOW MANY INITIAL CONDITIONS?') READ(C98,X)NIC

M=NSR-NIC

IFCM .GE. 0) GO TO 29

WRITEC99,998)

FORMATC/" MORE INITIAL CONDITIONS THAN STABLE ROOTS") GO TO 50

WRITEC99,99991) FORMATC/* ENTER -1 TO MAKE XBAR-BAR THE INITIAL CONDITIONS '7 OR 0 TO USE XBAR, ELSE 1 TO INPUT YOUR OWN : =

READ(938,%*) Y

IFCY) 777,778,779

DO 41 I=1,NIC E1(I,1)=E1C1I,NZ)-E1(I,1) GO TO 53

DO 781 I=1,NIC El(1I,1)=-E1(1I,1)

GO TO 53 WRITEC6,99992) FORMAT(/" WHAT ARE YOUR INITIAL CONDITIONS := ?'7)

READC98,%) CE1CI,1),1=1,NIC IFCM .EQ. 0) GO TO 52

DO 885 I=1,1A

DO 885 J=1,IA

N11,NSR ElcI,1 . WRITEC99,997) FORMAT(/* WHAT ARE THE VALUES OF FI1C(I,J) (BY ROW)?") READ(98,%)(CF1CI,J),J=1,M),1=1,M) WRITE(99, 996) FORMAT(/* WHAT ARE THE VALUES OF F2C(I,J) (BY ROW)?") READC98,%)(CF2C1,J),J=1,NIC),I=1,M) WRITE(99, 995) FORMAT(/* WHAT ARE THE VALUES OF Pacl, J) CBY ROW)?7") READ(98, ¥)(CF3(1I,J),J=1,NUSR),I=1,M WRITE(99, 994) FORMAT(/* WHAT ARE THE VALUES OF FC(I,1) 7") READ(98,%)(F4(1;1),I=1,M) WRITE(18, 88881) FORMATC* COEFFICIENT MATRIX F1') WRITE(18,100)(CF1CI,J),I=1,M),J=1,M)

oltrmooooceoeo o-r oooooo

Web ton mato ott

~ OQwwuwvwwy

WRITEC18,99995) BLANK

oy

88882

88885

888384

30 31

32 33

aou © N

993

99993

35

38 500

WG =

780

42 989

988

-~25-

WRITE(18,88882)

FORMATC * COEFFICIENT MATRIX F2"')

WRITEC18,100)(CF2CI,J),1=1,M),J=1,NIC)

WRITEC18,99995) BLANK .

WRITEC18, 88883)

FORMAT (C * COEFFICIENT MATRIX F3" )

WRITEC18,100)(CF3(1,J),1=1,M),J=1,NUSR)

WRITEC18,99995) BLANK

WRITEC18,88884)

FORMAT C " R.H.S. VECTOR F* )

WRITEC18,100)(F4(1I,1),1=1,M)

WRITE(18,99995) BLANK ,

CALL FOGAEFCF1,IA,B,IA,M,M,E2,1IA,WKSPCE,E4,1IA,E5,IA,IFAIL)

CALL FOICKFCE6,F3,E3,IA,1IA,IA,WKSPCE,IA, IOPT, IFAIL)

DO 30 I=1,M

NIC] = NIC#+1

DO 30 J=NIC1,NSR

K=J-NIC

F1CI,K)=FOIDAFCL,M,X,1,J,E2,1A,E6,1A,N)

DO 31 I=1,M

FICI,I)=F1CI,1)+1.0

DO 33 I=1,M

Y=0.0

DO 32 J=1,NIC

Y=Y+CF201,J)+E6C1I,J))XE1CJ,1)

F2(1I,1)=Y

cape FOGAEFCF1,1A,B,1A,M,M,E6,1A,WKSPCE, E4,1A,E5,1A,IFAIL) OPT=2

CALL FOICKFCE1,E6,E2,1A,1A,IA,WKSPCE, IA, IOPT,IFAIL) E2€1,1)=0.0

SET T(0) TO BE ZERO

WRITEC99,993) ;

FORMAT(7* WHAT ARE THE TRANSITION TIMES FOR THE Z(I)?7"*) READ(98,*)CE2(1I,1),I=2,NZ)

WRITEC18,99993) CE2CI,1),1=2,NZ)

FORMAT(“* TRANSITION TIMES FOR THE Z(I) GIVEN AS := '7 * C1OF13.6)7)

E2(NZ+1,1)=1000.0

IFCM.EQ.0)GOTO 780

CALL A2(C2,C11,CC,E2,C5,C13,C14,C15,C16,E5,C9,CX,X,NZ, 1 NUSR,N3,IA,L)

DO 35 I=1,M

F3(I,1)=FOLDAFC(L,NUSR,X,1,L,F3,1A,E5,1A,N)

DO 38 I=1,M

ELCI+NIC,1)=FOIDAF(L,M,X,1,L,E6,1A,F3,1A,N)

DO 34 I=1,M E1LCI+NIC,1)=E1LCI+NIC,1)-FOIDAF(L,M,X,1,L,E6,1A,F2,IA,N) ELCI+NIC,1)=E1CI+NIC,1)+FO1DAFCL,M,X,1I,L,E6,1A,F4,IA,N) DO 42 I=1,NSR

C7(1I,1)= DCMPLXC0.0D0,0.0D0)

DO 42 J=1,NSR

C7(1,1)= DCMPLXCE1(CJ,1),0.0D0)*C6(1,J)+C7(1,1) WRITEC99,989)

FORMAT(/' HOW MANY TIME PERIODS TO INTEGRATE OVER?") READ( 98, )NN

WRITEC99, 988)

FORMAT(/" HOW BIG A TIME PERIOD?'*)

READ(98, *)XEND

116 44

45 47

444 48

110 111

112

113 119

115 117

~26~-

T=0.0

DO 44 I=1,NUSR E1CI+NSR,1)=FOIDAFCL,NSR,X,1,L,E3,1A,E1,1A,N) CALL A2(C2,C11,CC,E2,C5,C6,C13,C14,C15,E5,C9,CX,1T,NZ,NUSR,N3,IA,L) DO 45 I=1,NUSR

K=I+NSR

DO 47 I=1,N2 ElCT,2)=FOIDAF(L,N1,X,1,L,E£7,1A,E1,IA,N) IFLAG=0

DO 48 I=1,NZ

IFCT.LT.E2CI+1,1))GOTO 444

GOTO 48

IFCIFLAG.EQ.0)KSTR=I

IFLAG=1

CONTINUE

DO 110 J=1,N2

ElCJ,2)=E1(J,2)+FO1DAF(L,N3,X,J,KSTR,E8,IA,Z,IA,N) WRITEC18,199)T, CE1(I,1),1=1,N1),CELCI,2),1=1,N2) T=T+XEND

NN=NN-1

IFCNN.LT.0)GOTO117

CALL A1(C4,T,E2,C6,L,NSR,IA)

DO 112 I=1,NSR

CALL FOIDDF(L,NSR,CX,I,L,C6,1A,C7,1A,N,CR,CI) E4(I,1)=-CR

E4(I,2)=-clI

CALL A4(C4,C5,CC,C10,C11,E2,C6,C13,C14,C15,C16,C1,E5,CX, 1 T,NUSR,NSR,N3,NZ,IA,L)

CALL A6(C4,E2,CX,T,C6,C13,C8,CC,C14,C15,C1,NSR,IA,N3,L) DO 113 I=1,NSR C15€T,1)=DCMPLXCEG4(1,1),E4(1,2))4+C14(1,1)+DCMPLX(E5(1,1),E5(L,2)) DO 115 I=1,NSR

CALL FOIDDF(L,NSR,CX,1,L,C3,1A,C15,IA,N,CR,CI) E1(I,1)=-CR

GOTO 116

REWIND 18

REWIND 19

RETURN

END

SUBROUTINE A1(A,T,B,C,I,N,IA)

IMPLICIT REAL¥®8(A-H,0-Z)

COMPLEX¥16 ACIA,IA),CCIA,IA)

REAL%8 BCIA),X,Y,T,1T

TT=T-BCI)

DO 2 K=1, CCJ,K)= DCMPLX(0.0D0,0.0D0)

X=DEXPC (CDREALCACJ,J))) TT)

Y=TTXDIMAGCACJ,J))

CCJ,J)= DCMPLX(CXXDCOSCY)), (XXDSINCY)))

SUBROUTINE A2(A,B,C,D,F,G,H,P,Q,R,S,CX,T,NZ,NUSR,N3,IA,L) IMPLICIT REAL*®8(A-H,0-Z)

COMPLEX%16 ACIA,IA),BCIA,IA),CCIA,IA),FCIA,IA) 1,GCIA,IA),

1 HCTA,IA),PCIA,IA),QCIA,IA),SCIA, IA) ,CX

RG AL xs DCIA),RCIA,IA),T,CR,CI

-27-

CALL A5(B,C,CX,Q,NUSR,N3,I,IA,L) IFCDCI)-T)2,2,3

DO 1 J=1,NUSR

PCJ,1)=Q0J,1)

GOTO 14

CALL A1CF,T,D,H,I,NUSR,IA)

DO 4 J=1,NUSR

CALL FOIDDFCL,NUSR,CX,J,L,H,IA,Q,IA,N,CR,CI) PC = DCMPLX(¢ -CR,-C1)

I -

CALL A5(B,C,CX,Q,NUSR,N3,1,1A,L)

IF (DCI) 135,526

CALL A1CF,T, »I,NUSR,IA)

DO 7 J=1,N UR’

DO 8 J=1,NUSR

CALL FOLDDF(L,NUSR,CX,J,L,H,IA,Q,IA,N,CR,CI) PCJ,1)=PCJ,1)- DCMPLXCCR,CI)

DO 9 J=1,NUSR

HCJ,J)=GOJ,J)

I=I-1

GOTO 10

DO 11 J=1,NUSR H(CJ,J)=DCMPLX(1.0D0,0.0D0)-HCJ,J)

DO 13 J=1,NUSR

CALL FOILDDFCL,NUSR,CX,J,L,H,IA,Q,IA,N,CR,CI) PCJ,1)=DCMPLX(-CR,-CI)+P(J,1)

DO 15 J=1,NUSR

CALL FOIDDFCL, NUSR,CX,J,L,A,IA,P,IA,N,CR,CI) G(J,1)= DCMPLX¢~ CR,-CI)

DO 16 J=1,NUSR

CALL FOLDDFC(L,NUSR,CX,J,L,S,1A,G,1A,N,CR,CI) RCJ,1)=-CR

RETURN

END

SUBROUTINE ASCA,B»C,D,E,F»G/H,T,1,NUSR,NSR,L,IA,KFLAG) IMPLICIT REAL¥8(CA-H,0-Z)

COMPLEXX16 ACIA,IA),BCIA,IA),CCIA,IA),DCIA,IA),ECIA,IA), 1 FCIA,IA),GCIA,IA)

REAL x8 HCIA), X,XX,TT,»T,Y,ZR,Z1

TT=0.0

aire tia ean 4,5

CALL A1(A,T,H,C,1,NSR,IA) GOTO 6

CALL A1(B,T,H,C,1I,NUSR,IA) CALL A1(A,T,H,D,L,NSR, IA) CALL A1CB,TT,H,E,1I,NUSR,IA) IFCI-KFLAG)14,14,15

DO 1 J=1,NSR

DO 1 N=1,NUSR

XX=1.0 X=DREALCBCN,N)-ACJ,J)) Y=DIMAGCB(N,N)-ACJ,J))

CALL AO2ACFCXX,TT,X,Y,ZR,ZI)

X=DREAL(CCJ,J))-CDREAL (DCJ, J) (DREALCECN,N)) )-CDIMAGCD(J,J)))% * (DIMAGCECN,N)))

Y=DIMAG(CCJ,J))-CDREAL (DCJ, J) ) )*CDIMAGCECN,N)))-CDIMAG(D(J,J)))% % CDREALCECN,N)))

XX=DREAL (GC J,N) )®X-DIMAGCGCJ,N) )*Y

Y=DIMAGC(GCJ,N))*X+DREAL (GC J,N) )¥Y

15

ll

12 10

1

11 23

-28-

FCJ, NS 5 DOMPLXC CXX¥ZR- ~YXZ1) , CXXX¥ZI+YXZR) )

GOTO

DO 11 ae 1,NSR

DO 11 N=1,NUSR

XX=1.0

X=DREALCBCN,N)-ACJ,J))

Y=DIMAGCBCN,N)-ACJ,J))

CALL AQ2ACFC(XX,TT,X,Y,ZR,ZI) X=DREAL(CCON,N))-CDREALCDCJ,J)) XC DREALCECN,N) ) -CDIMAGCD(J,J)))% * CDIMAGCECN,N)))

Y=DIMAGCC(N,N) )-CDREALC DCJ, J) XC DIMAGCECN,N) ) )-CDIMAGCD(J,J)))% * CDREALCECN,N)))

XX=DREAL (GC J,N))*X-DIMAG(GCJ,N) XY Y=DIMAG(GCJ,N))*X+DREAL(GCJ,N) )*Y

FCJ,N)=DCMPLXC CXXX¥ZR-YXZI ) , CXXXZI+YXZR) )

CONTINUE

RETURN

END

SUBROUTINE AGCA,B,C,D,E,F,G,H,P,Q,R,S,U,CX,T,NUSR,NSR,N3,NZ,IA,L) IMPLICIT REAL*8(CA-H,0-Z)

COMPLEX¥16 ACIA,IA),BCIA,IA),CCIA,IA),DCIA,1IA),ECIA,IA), 1 GCIA, IA), HCIA,IA),PCIA, IA), QCIA, IA), RCIA, IA),SCIA,IA),CX REALX8 FCIA),UCIA,1IA),T,CR,CI

DO 11 I=1,NSR

UCI,1)=0.

UCI,2)=0.

KFLAG=1200

DO 12 I=1,NZ

IFCT.GE. FCI). AND. T.LT.FCI+1) )KFLAG=I

CONTINUE

I=NZ

CALL A3(A, B, P,Q,R,G,D,F,1,T,NUSR,NSR,L,IA,KFLAG)

DO 3 J=1,NSR

DO 3 K=1,N3

CALL FOIDDFCL,NUSR,CX,J,K,G,IA,E,IA,N,CR,CI) HCJ,K)=DCMPLX(-CR,-CI)

K=I-1

DO @ J=1,N3 GCJ,1)=DCMPLXCDREAL(CCJ,K))-DREAL(C(J,1)),0.0D0)

DO 5 J=1,NSR

CALL FOIDDF(L,N3,CX,J,L,H,IA,G,IA,N,CR,CI) UCJ,1)=-CRtHUCJ,1)

UCJ,2)=-CI+UCJ,2)

I=I-1

IFCI-2)2,10,10

IFCKFLAG.EQ.1)GOTO 1

NFLAG=KFLAG-1

DO 113 IM=1,NFLAG

JF=IM+1 CALL A1CA,T, ee A, IM,NSR,IA) CALL A1CA,T »JF,NSR,IA)

DO 112 J= 1, ASR” RASS ,J)= DCMPLX(DREAL(P(J, J)I-DREALCHCJ,J)), DIMAGCPCJ, JII-LIMAGCHCJ, CONTINUE

CALL A5CE,C,CX,G,NUSR,N3,IM,IA,L)

DO 111 J= WSR

DO 111 K=1,N5R

CALL FOIDDFCL NSR,CX, 3° K,S,IA,P,IA,N,CR,CI)

QC, K= DCMPLX(-CR, -CI)

115

119 113

g

%)

12

mo WwW Ne

DN

-29-

DO 115 J=1,NSR

DO 115 K=1,NUSR

CALL FOIDDFC(L,NSR,CX,J,K,Q,IA,D,1IA,N,CR,CI) PCJ,K)=DCMPLX(-CR,-CI)

DO 119 J=1,NSR

CALL FOIDDF(L,NUSR,CX,J,L,P,IA,G,IA,N,CR,CI) UCJ,1)=-CR+UCJ,1)

UCJ,2)=-CI+UCJ,2)

CONTINUE

CONTINUE

CALL A5(E,C,CX,G,NUSR,N3,KFLAG,IA,L)

CALL A1(CA,T,F,P,KFLAG,NSR, IA)

DO 6 J=1,NSR PCJ,J)=DCMPLX(1.0D0,0.0D0)-P(J,J)

DO 7 J=1,NSR

DO 7 K=1,NSR

CALL FOIDDF(L,NSR,CX,J,K,S,IA,P,IA,N,CR,CI) QCJ,K)=DCMPLX(-CR,-CI)

DO 8 J=1,NSR

DO 8 K=1,NUSR

CALL FOIDDFC(L,NSR,CX,J,K,Q,IA,D,IA,N,CR,CI) PCJ,K)=DCMPLX(-CR,-CI)

DO 9 J=1,NSR

CALL FOIDDF(L,NUSR,CX,J,L,P,IA,G,IA,N,CR,CI) UCJ,1)=-CRt+UCJ,1)

UCJ,2)=-CI+UCJ,2)

RETURN

END

SUBROUTINE A5(A,B,CX,C,NUSR,N3,1,IA,L) IMPLICIT REAL¥8(A-H,0-Z)

COMPLEXX16 ACIA,IA),BCIA,IA),CX,CCIA,IA) REAL%8& CR,CI

DO 1 J=1,NUSR

CALL FOIDDF(L,N3,CX,J,1,A,IA,B,IA,N,CR,CI) CCJ,1)=DCMPLX(-CR,-CI)

RETURN

END

SUBROUTINE A6(A,B,CX,T,D,E,F,G,H,P,Q,NSR,IA,N3,L) IMPLICIT REAL¥8(A~H,0-Z)

COMPLEXX16 ACTA, IA),CX,DCIA,IA),ECIA,IA),FCIA,IA),

1 GCIA, IA), HCIA,IA),PCIA,IA),QCIA,IA)

REAL®8 BCIA),T,CR,CI DO 12 J=1,NSR FOS 1)=C0.,0.)

DO 6 J=1,NSR DCJ,J)=DCJ,J5)-DCMPLX(1.0D0,0.0D0)

DO 7 J=1,NSR

DO .7.K=1,N3

CALL FOIDDFC(L,NSR,CX,J,K,D,IA,F,IA,N,CR,CI)

-30-

HC J,K)=DCMPLXC€-CR, -CI) K=I-]

DO 8 J=1,NSR

CALL FOIDDF(L,N3,CX,J,K,H,IA,G,IA,N,CR,CI) PCJ,1)=PCJ, 1)+DCMPLX¢ “CR, -CI)

IFCBCI). LT.T)GOTO 9

DO 11 J=1,NSR

CALL FOIDDFCL, NSR,CX,J,L,Q,IA,P,IA,N,CR,CI) HCJ,1)= DCMPLX(- CR, ~CI)

RETURN

END

31

International Finance Discussion Papers i cussion Papers

IF DP NUMBER TITLES 1986

276 Post~simulation Analysis of Monte Carlo Experiments: Interpreting Pesaran's (1974) Study of Non*nested Hypothesis Test Statistics

275 A Method for Solving Systems of First Order Linear Homogeneous Differential Equations When the Elements of the Forcing Vector are Modelled as Step Functions

274 International Comparisons of Fiscal Policy: The OECD and the IMF Measures of Fiscal Impulse

273 An Analysis of the Welfare Implications of Alternative Exchange Rate Regimes: An Intertemporal Model with an Application

1985

272 Expected Fiscal Policy and the Recession of 1982

271 Elections and Macroeconomic Policy Cycles

270 Assertion Without Empirical Basis: An Econometric Appraisal of Monetary Trends in... the United Kingdom by Milton Friedman and Anna J. Schwartz

269 Canadian Financial Markets: The Government's Proposal for Reform

268 Was It Real? The Exchange Rate Interest Differential Relation, 1973-1984

267 The U.K. Sector of the Federal Reserve's Multicountry Model: The Effects of Monetary and Fiscal Policies

266 Optimal Currency Basket in a World of

Generalized Floating: the Nordic Countries

An Application to

AUTHOR(s)

Neil R. Ericsson

Robert A. Johnson

Garry Schinasi

Andrew Feltenstein David Lebow Anne Sibert

William H. Branson Arminio Fraga Robert A. Johnson

Kenneth Rogoff Anne Sibert

David F. Hendry Neil R. Eriesson

Garry J. Schinasi

Richard Meese Kenneth Rogoff

Hali J. Edison

Hali J. Edison Erling Vardal

Please address requests for copies to International Finance Discussion Papers, Division of International Finance, Stop 24, Board of Governors of

the Federal Reserve Board, Washington, D.C. 20551.

IFDP NUMBER 265

264

263

262

261

260

259

258

257

256

255

254

253

252

251

250

32

International Finance Discussion Papers a et sston rapers

TITLES

Money Demand in Open Economies: Substitution Model for Venezuela

A Currency

Comparing Costs of Note Issuance Facilities and Eurocredits

Some Implications of the President's Tax Proposals for U.S. Banks with Claims on Developing Countries ,

Monetary Stabilization Policy in an Open Economy

Anticipatory Capital Flows and the Behaviour of the Dollar

Simulating Exchange Rate Shocks in the MPS and MCM Models: An Evaluation

Trade Policy for the Multiple Product Declining Industry

Long Memory Models of Interest Rates, the Term Structure, and Variance Bounds Tests

Currency Substitution and the New Divisia Monetary Aggregates: The U.S. Case

The International Transmission of Oil Price Effects and OPEC's Pricing Policy

U.S. Banks' Lending to Developing Countries: A Longer~Term View

Conditional Econometric Modelling: An Application to New House Prices in the United Kingdom

Low Pushing: Doctrine and Theory

1984 (partial listing) Postwar Financial Policies in Taiwan, China

Foreign Exchange Constraints and Growth Possibilities in LDCs

The Determination of Front-end Fees on Syndicated Eurocurrency Credits

AUTHOR(s)

Jaime Mar quez Rodney H. Mills

Allen B. Frankel

Marcus H. Miller Arnold Kling Arnold K_ing Catherine Mann Gary S. Shea Jaime Mar quez Jaime Mar quez Henry S. Terrell Rod Mills

Neil R. Ericsson

David F. Hendry

William Darity, Jr.

Robert F. Emery

Jaime Mar quez

Rodney H. Milis Henry S. Terreli