Saddle node bifurcation at a nonhyperbolic limit cycle in a phase locked loop

Saddle Node B$iucation at a Nonhyperbolic Limit Cycle in a Phase Locked Loop by JOHNSTENSBY Electrical and Computer Engineering Huntsville, AL 35899, U.S.A.

Department,

University of Alabama

in

ABSTRACT : New results are given on the phenomenon of false lock in second-order, Type I phase locked loops (PLLs) with a constant frequency reference of wi and a voltage controlled oscillator (VCO) quiescent frequency of wO. This loop has on[v one false lock state whose frequency error approaches Wi- co0 as Ioop gain 6 approaches zero, and this false lock state is stable. This false lock state corresponds to a stable, hyperbolic limit cycle X,(t ; 6) of the nonlinear equation describing the loop. As gain 6 is increasedfrom a value of zero, it is shown that a value 6, can be reached where X,Ybecomes semi-stable and nonhyperbolic. Furthermore, saddle node bifurcation occurs at 6,, and a second limit cycle X,(t; 6) branches from this bifurcation point. Limit cycle X, is unstable, and it corresponds to an unstable false lock state of the PLL. Furthermore, X, can be continued as a function of gain on an interval o2 < 6 < 6,, for some 6, > 0. Finally, ZXY and X, do not exist for 6 > 6,. Two numerical algorithms are given to anaIyse the faIse locked PLL under consideration. The first is useful for computing the above-mentioned limit cycles and the bifurcation point 6,. The second algorithm can calculate a Poincare map and its derivative which are useful in studying the saddle node bifurcation which occurs at 6 ,. Also given is a detailed decription of a laboratory experiment which was used to substantiate the theory and numerical techniques. The quahtative theory, numerical method and laboratory procedure are applied to a simple example. The numerical and empirical results are shown to be in close agreement. I. Introduction

Phase locked loops (PLLs) are an important part of modern electronic communication and control systems. They play a major role in systems where it is necessary to estimate the phase of a received signal. Also, they are used extensively in applications requiring the synthesis of highly stable sinusoidal signals. Many applications can be served by PLLs which can be modeled by the loop depicted by Fig. 1. The voltage controlled oscillator (VCO) in this second-order,

YCO

“CO center

Gain = Kv Fly. = w.

FIG. 1.Second-order Type I PLL.

The Franklin Institute 00164032/93

%6.00+0.00

775

J. Stenshy Type I loop produces a sinusoid whose instantaneous frequency is linearly related to control voltage e; let o. and K, denote the VCO’s quiescent frequency and gain, respectively. The loop’s analogue multiplier forms the product of the sinusoidal reference and VCO voltages; assume that this multiplier has unity gain. Filter F processes this product and supplies the result to the VCO. This PLL is considered to be phase locked when the loop’s phase error 4(t) = o,t-O,(t) is a constant, and the loop is in a stable equilibrium state. It is easy to show that phase lock is possible for

and that it cannot occur for 6 < 6, (1,2,3). This PLL can be in an anomalous state known the closed-loop phase error has the form #J(t ; 4 = dW+

as false lock (2,4).

44t; a

In this state

(2)

where $(t ; 6) is periodic with a fundamental frequency of wf, the apparent difference between the reference and VCO frequencies. When interpreted modulo 271, phase error 4 has a period of T(6) = 2n/or(6). Finally, phase error (2) corresponds to a limit cycle of the nonlinear equation which describes the PLL. It is known that this PLL has a stable false lock state corresponding to a hyperbolic limit cycle X,(t ; S) for 6 in a sufficiently small neighborhood of the origin. Furthermore, X,‘s fundamental frequency o,(6) approaches oos = o, - o. as 6 goes to zero (4). In Sections III and IV, it is shown that a value of gain 6, can be reached where X, becomes nonhyperbolic and semi-stable. Furthermore, saddle node bifurcation occurs at 6,. As is characteristic of this type of bifurcation, a second limit cycle X,, bifurcates from 6,. This cycle is unstable, and it can be continued on some interval a2 < 6 < 6,, where 6, is positive. Neither the stable or unstable limit cycles exist for 6 > 6,. This is shown by the use of a bifurcation diagram and a Poincare map, and algorithms are supplied for computing these graphical aids. These results on bifurcation in PLLs have not appeared in the literature, and they are the main contribution of this paper. Finally, Section VII describes a laboratory experiment which was used to validate the theory and numerical procedure. The experimental results are in close agreement with those obtaind by the numerical algorithms.

II. Dynamic Model The loop depicted

in Fig. 1 is described

d24

by

d4

~+(h+6Cos~)dl+a6sin~=bo,,, dtwhere 6 s AK,K,,

oOs 5 oi-CO”, and the pair K,, oO denotes

the VCO’s gain and Journal

776

c&the Franklin lnsulute Pergamon Press Ltd

Phase Locked Loops center frequency, respectively. phase error (2,3). Equation (3) can be written X(t,s)

Finally,

4(t) = oil--e,(t)

in state variable

= [x(t; S)

Y(t;6)lT

F(X;6)

=

the closed loop

form by defining

= [#(t;@

d4(t;6)/dtlr,

where bold face upper case letters denote vectors, transpose. The state equation for X(r; S), becomes dX = dt

denotes

and

(4)

superscript

T denotes

F(X;G)+G Y - [b + 6 cos (x)JY - a& sin (x)

The first variational equation with respect to a solution in what follows (5). This equation is given as dZ dr = F,(X;

1

X(t; S) of (5) is required

(6)

4%

where

F,(X;6) =

0

1

Y6 sin (x) - a6 cos (x)

-b - 6 cos (x)

1

is the Jacobian of F evaluated at solution X. Under false lock conditions (6) has periodic coefficients. In this case the characteristic exponents of this equation play a significant role in the structural stability of (5) and the ability to continue X,(t ; S) as a function of gain 6.

III. Properties

of the Stable False Lock State

For sufficiently

small 6, Eq. (5) has a stable solution,

i.e.

(7) where $ is periodic has the property

with a fundamental

frequency

X,(t ; 6) -+ X,(t; 0) =

of wf (4). Furthermore,

[-sr 1

X,(t ; S)

(8)

as 6 -+ 0. X,(t ; S) is periodic with fundamental frequency u_+(S) if the first entry of this vector is interpreted modulo 27~ [period T(6) = 271/o,-(6)]. This solution is a stable limit cycle of (5), and it corresponds to a stable false lock state of the PLL. Vol. 330, No. 5. pp. 77%7X6, 1993 Printed m Great Bntain

777

J. Stensby It is the only false lock state which exists for 6 in a sufficiently small neighborhood of the origin (4). First variation (6) with respect to X,(t ; 6) [in (6) replace X by X,] has periodic coefficients. This equation has a characteristic exponent 1, = 0 (5,6). Also, it has a second characteristic exponent : i.e. = -[b+6cos(o,t+$(t;13))],

Iz,(S) rTraceFY(X,;6)

(9)

where the overbar denotes that the average value is to be taken (6). Thus X, is a stable, hyperbolic limit cycle, and vector field F is structurally stable, as long as &(6) < 0 (7,8). Under these conditions, the false lock state corresponding to X, can be continued as a function of 6. Note that 1, given by (9) is proportional to Andronov’s J factor [see p. 103 of (7)].

IV. Saddle Node Bifurcation at 6, The PLL remains in a stable false lock state as long as i, remains in the left half plane. However, as 6 -+ 6;, characteristic exponent A,(6) + &(d ,) = 0 from the left, and X, becomes a nonhyperbolic limit cycle with period T(6,). Furthermore, vector field F becomes structurally unstable at 6, (7). The first entry of the nonhyperbolic limit cycle at 6, is the closed-loop phase error, and it has the form given by (2). Hence, the limit cycle passes through the y-axis in the phase plane ; let p0 = [0 y,lT denote the point. Let C denote a straight line which passes through p,, on the phase plane’s y-axis ; furthermore, we require that the transversal C be normal to the nonhyperbolic limit cycle at pO. For Eq. (5), let P(s;6) denote the Poincare map along the transversal Z [see Perko (8) for a discussion of this map]. The quantity s denotes the signed distance along C from the nonhyperbolic limit cycle. If the point p, = [x y]’ on C has y < y, (y > y,,), then distance s from p0 to p, is taken as negative (positive). Poincare map P can be used to establish the occurrence of saddle node bifurcation at 6,. First, note that

P(O;6,) = 0. Next, the derivative

(10)

of P at s = 0,6 = 6, is

1

T(6,)

=

exp

[S -

(b+h,cos4(t;h,))dt

,

0

where +(t ; S,) is the first entry in X,(t ; 6,), and T(6,) is the period cycle [see p. 332 of (S)]. Since n,(s,) = 0, Eqs (9) and (11) produce

&6,

6=6,

=

1.

of this limit

(12)

.s=0 Journal

778

(11)

of the Franklin Institute Pergamon Press Ltd

Phase Locked Loops Saddle node bifurcation cation that occurs under conditions,

of a hyperbolic limit cycle is the generic type of bifurconditions (10) and (12), (8, 9). In addition to these

dP ,,(s;6)

6=6, # 0

(13)

S=O

is the generic hypothesis

under which saddle node bifurcation

is certain

[see p. 64

of (9>1. As is characteristic of saddle node bifurcation at a nonhyperbolic limit cycle, an unstable periodic solution X, bifurcates from the point 6,. As a function of 6, it can be continued on some interval 6, < 6 d 6,, where d2 is positive [X, is the only limit cycle for sufficiently small 6 ; see Stensby (4)]. Equation (9) with X, substituted for X, yields a positive value of 2,. Both X, and X, become semi-stable, nonhyperbolic limit cycles at 6,, and neither cycle exists for 6 > 6 ,.

V. Computation

of X,, X, and the Bifurcation Diagram

In this section 6, X, and X, are considered to be functions of a user-supplied A2. Given Al, the algorithm described below computes the periodic solution with the specified exponent 1,. The solution is either X,, the nonhyperbolic limit cycle, or X, depending on whether the supplied 1, is negative, zero or positive, respectively. Also, the algorithm computes the fundamental frequency wf of the solution and the value of 6 at which the solution occurs. For example, if A2 = -b is specified, the algorithm produces the limit cycle given by (8) as well as 6 = 0 and wr = oos. If /2, = 0 is used as input, the result is the nonhyperbolic, semi-stable limit cycle that occurs at bifurcation ; also computed are the frequency wr of this limit cycle and the bifurcation point 6,. Finally, positive AZ results in the algorithm computing 6 and of for the unstable hyperbolic cycle XU. Let [0 y,lT denote the y-axis intercept point for the limit cycle. Note that using & = 0 results in y, = y, as was discussed in Section IV. For a specified I,, unknowns 6, or and y , must satisfy the nonlinear algebraic equations g,(&cUr,y,)

= x(2rL/cor;6)-2271=

gx(&of,v,)

= 1.,+(b+6cos(x(t;S)))

0

= 0,

(14)

where oos = w,--o. is assumed to be positive so that phase x increases by 271over one period T = 27c/or [the overbar in (14) denotes average value]. Note that the quantities x(2n/or; S) and y(27c/or; S) in (14) are obtained by numerically integrating (5) subject to the initial condition [0 y,lT. The nonlinear algebraic system given by (14) can be solved for 6, or and y1 by numerical integration coupled with standard routines for solving algebraic equations. Routine ZSYSTEM in the IMSL package (the more recent routine ZSPOW also works) was used to solve the algebraic equation. ZSYSTEM Vol. 330, No. 5, pp. 775~786, 1993 Printed m Great Bruin

779

J. Stensby

200 j

1’,“““~,,1,“‘rlllll 0 40

80

T

ll”‘l,“‘nl”l 120

160 Gain

FIG.

2. Saddle node bifurcation

diagram

6

with IL2appearing

parametrically.

accomplishes this by calling a user-supplied subroutine that produces x(27c/or; 6) and y(27c/or;6) by numerically integrating (5). The fourth-order Runge Kutta algorithm was used in the numerical integration. The algorithm was applied for the case, a = 100, b = 10 and woS = 20071. It was run for a spread of A>on the interval -b < iz < 41.5. The results for AZ = 0 yielded a computed value of 6, = 143.97 for the bifurcation point. Figure 2 depicts y , plotted as a function of 6 with i, appearing as a parametric parameter displayed on this bifurcation diagram. The pairs (y,, 6) were obtained as solutions of (14) with the user-specified values of AZ. Note that the solid line segment of the graph (with - 10 d AZ < 0 displayed parametrically) corresponds to the stable limit cycle X,. The dashed line segment (with 0 < A2 < 41.5) corresponds to the unstable limit cycle X,. Finally, the point 6, = 143.97, y, = 37 1.75 on the bifurcation diagram corresponds to the semi-stable limit cycle. Figure 3 depicts computed gain 6 as a function of A2, with or appearing as a parametric parameter displayed on the graph. With the plotted values of 6 to the left (to the right) of the vertical axis, (5) is satisfied by stable limit cycle X, (unstable limit cycle X,), and (6) has exponent AZ, - 10 < A2 < 0 (0 < A, < 41.5), as dis-

/I,”

Unstable 13 18

FIG. 3. Computed

gain 6 as a function

23 Limit 28 Cycle 33

-38

A2

of i, with wy appearing

parametrically. Journal

780

of the

Franklin lnstltule Pergamon Press Ltd

Phase Locked Loops played on the negative (positive) abscissa. The limit cycle has or (some values of which appear parametrically on the graph) as its fundamental frequency. Note that 6 = 0 at AZ = -b = - 10, and that bifurcation point 6, = 143.97 occurs at A2 = 0. On Figs ‘2 and 3 note that for i2 > 41.5 (corresponding to 6 < 121.5) the unstable limit cycle is not continued. The algorithm described in Section IV becomes very sensitive to initial guesses for the roots of (14), and for lb2 > 41.5 it would not converge. Hence, for 6 < 121.5 the exact limiting nature of unstable X, is unknown at this time. However, it is known that X, cannot be continued down to 6 = 0; that is, as a function of 6 it fails to exist for 0 < 6 < d2, where d2 is positive. This results from the fact that X, is the only limit cycle of (5) for sufficiently small 6 (4). At present, the quantity 6, is unknown. Also unknown is the limiting nature of X, as 6 + S:. However, it is conjectured that X, approaches a homoclinic loop (also known as a separatrix cycle) as 6 + 8:.

VI. Computation

of the Poincavt! Map and Derivative

Recall from Section IV that C represents a transversal which passes through the point p0 = [0 y,JT, where y, comes from the solution of (14) with EWz = 0. Furthermore, the nonhyperbolic limit cycle (which occurs at 6,) is normal to C at pO. Hence, the slope of Z is -I

1

_

(15)

which was obtained easily by dividing the components y-y,

of (5). Finally, the linear equation

= mx

(16)

represents C on the phase plane. Let the quantity s denote the signed distance along C from [0 y,jr as was discussed in Section IV. Also, let P(s ; S) denote the Poincare map along C for (5). An algorithm for computing P(s ; S) follows. At t = 0, start (5) at the point which lies s units along C; this initial condition which lies on C is

[x(O;S,s)

y(o,s,.~)]T

=

[g&

(

y”+,+l

ins r

)I T

(17)

Note that state vector notation (4) has been modified to indicate the dependence on initial conditions through s. Then at time t = z(s;6) the system returns to the point [x(z(s; S) ; S, s) -271 y(z(s; S) ; 6, s)]’ on C; this first return point is at distance P(s; 6). From Eq. (16) it is easy to calculate that P(s;6)

= Jm2 + I[x(t(s;

S) ; 6,s) -2711.

It is important to remember that s enters the right-hand side of (18) through the initial conditions (17) used to compute x(z(s ; S) ; 6, s), Vol. 330. No. 5, pp. 775-186, Pruned in Great Britam

(18) T and

1993

781

J. Stensby Note that first return

time r(s; 6) must satisfy

g(s) = {y(s ; 6, s) -y,}

- m(x(r ; 6, s) - 27r) = 0,

(19)

where it is assumed that the system was started at the given initial point (17) which lies s units along X. Equation (19) can be solved for r(s ; S) by using ZSYSTEM and integrating (5). Then the numerical value P(s ; 6) can be computed by using (18). From (18) the derivative

of P can be calculated

as

1

dP &==&FTT

(20) The partial derivative with respect to s reflects the dependence of the state on s through the initial condition (17). Development of the derivative dt/ds on the right-hand side of (20) follows. From (16) note that

which leads to

ax dt -_= ds

ay

mas-as dy z-my

ax ay dz ~=~ ds

mas-as - [b + 6 cos (x) + m] y - a6 sin (x) + bwOs’

(21)

where the dependence of (x, y) on (7, 6, s) has been suppressed. An equation follows which can be integrated to produce partial x, = ax/as and y, E ay/& for use in (20) and (21). Take the partial derivative of (5) to obtain

dx,_ - Ys dt dys = 6[ y sin (x) -a dt

xslr=o= ysl,=o =

782

cos (x)]x, - [b + 6 cos (x)] y,,

+ m2+1 J&i

(22) Journal of the Frankhn Institute Pergamon Press Ltd

Phase Locked Loops

FIG. 4. Poincark

map P(s ; S,) for nonhyperbolic

limit cycle.

The partial

derivatives x, and y, are obtained by integrating (22) simultaneously with (5) on the time interval [0, z(s ; 6 ,)I. They are used in (21) to compute dr/ds and in (20) to compute dP/ds. This completes the algorithm for calculating dP/ds. The procedure based on (18) and (19) for calculating P(s ; 6 ,) was applied to the case discussed in Section V where a = 100, b = 10 and ooS = 2001~ Also, the derivative dP(s ; 6 ,)/ds was computed using (20), (21) and (22). Figures 4 and 5 depict the results for these loop parameters. These figures contain a lot of information about the nonhyperbolic limit cycle at 6 ,. Note that P(s; 6,) has the fixed point s = 0 which corresponds to this limit cycle. Also, Fig. 5 implies that the Poincart map lies below the unity slope straight line which passes through the origin in Fig. 4. This implies that the nonhyperbolic cycle is semi-stable ; trajectories above (below) the cycle converge to (diverge from)

FIG. 5. Derivative Vol. 330, No. 5, pp. 775-786, 1993 Printed in Great Britain

of PoincarC map P(s; 6 ,) for the nonhyperbolic

limit cycle.

783

J. Stensby 0.3

-0.5

,P(O;6)

3

FIG. 6. Poincark map P(0 ; 6) as a function

of 6 in a neighborhood

of 6, = 143.97.

the cycle. Finally, note that dP(0 ; 6 ,)/ds = 1 ; under this condition, saddle node is the generic type of bifurcation which occurs at a nonhyperbolic limit cycle. Generic condition (13) is verified by inspection of Fig. 6, a plot of P(0 ; S) for 6 in a neighborhood of bifurcation point 6 ,. This plot was constructed by using (18) and (19) in a manner similar to that used to produce Fig. 4. However, this time s was set to zero, and 6 was varied. Finally, note from Figs 5 and 6 that d*P(O ; 6 J/d? and dP(0 ; 6 J/d6 have the same sign. This implies that X, and X, exist for some 6 < 6, and not for 6 > 6, [see Theorem 1, p. 331 of Perko (S)]. Also, this conclusion can be made from inspection of Fig. 2.

VII. Experimental

Results

The PLL discussed in Sections V and VI was set up in the laboratory. The reference and VCO were implemented using Exact Model 126 signal generators running in their sinusoidal mode. The output amplitude of each generator was set to 2 V RMS. The reference generator was set to 1.1 kHz, and the center frequency of the VCO was set to 1.O kHz so that oos = 2007r. The phase comparator consisted of an RCA3091D four-quadrant multiplier followed by a wide-band operational amplifier ; the gain of this phase comparator was set to unity. An active loop filter with F(s) = (s+ lOO)/(s+ 10) was used ; the filter’s break frequencies were set by hand with particular emphasis on establishing accurately the filter’s DC gain F(0). The ability to set accurately the overall DC loop gain was critical to the success of the experiment. Hence, a high-grade, IO-turn, linear (to within 0.25%) potentiometer and precision IO-turn counter knob (calibrated to read from 1 to 10 in steps of 0.01) were used to adjust loop gain. Loop gain calibration was established by (1) providing enough gain to achieve phase lock, and (2) slowly decreasing the gain until the loop dropped out of lock. At this point the IO-turn counter knob reading corresponded to 6, = 62.83 given 784

Phase Locked Loops by(2) ; this is the smallest value of loop gain for which phase lock is possible. Now, loop gain at any knob setting could be calculated easily by assuming a linear potentiometer so that gain is proportional to knob setting. This method was used to obtain the results outlined below. The first measurement yielded the knob setting at which loss of lock occurs. Once phase locked, the loop would fall out of lock when the precision lo-turn gain potentiometer was reduced slowly to a value of 2.3 1. Hence, the empirically-derived formula : 6 = (knob setting)

hw& 1=

T [.

27.20 (knob setting)

(23)

gives loop gain in terms of counter knob setting. The second measurement produced the knob setting corresponding to 6,. The knob was slowly increased from a starting point of zero until the stable false lock state was broken and the loop locked up. This occurred at 5.28 on the counter knob. Hence, the empirically determined value of 6, is 143.6. Note that this compares favorably with the computed value 143.97 discussed in Sections V and VI. The methods and results outlined here can serve as an independent check on an existing approximate technique for computing 6, [Section 10-2.4 and p. 469 of Lindsey (2)]. The need for such a check results from the fact that it is not always clear for what values of loop parameters the existing technique is applicable (due to the approximations made). Replace A0 and AK in the cited reference by woSand i?(a/b), respectively. Also, use z, = a- ’ and z2 = b- ‘. Then 6, can be approximated as the solution of (24) Use c()~>= 2007c, a = 100 and b = 10 in (24) to obtain to 6,.

140.5 as an approximation

VIII. Conclusions New results were given on the phenomenon of false lock in second-order, Type I PLLs. The results show that saddle node bifurcation can occur where two hyperbolic limit cycles bifurcate from a nonhyperbolic, semi-stable one. This bifurcation is a mechanism responsible for breaking false lock in the loop. The limit cycles which bifurcate have opposite stability, and they correspond to false lock states in the PLL. The stable limit cycle corresponding to the stable false lock state starts at 6 = 0 with a closed loop frequency error wr of was, where o,,~ E o, - w. is the difference between the frequency of the reference and the VCO’s center frequency. Error wi decreases with increasing 6, and it was shown that a value 6, could be reached where the limit cycle becomes nonhyperbolic and semi-stable. Furthermore, the limit cycle and its associated false lock state do not exist for 6 > 6,. Vol. 330, No. 5, pp. 775-786,

Printed

in Great

Britam

1993

785

J. Stensby The stable limit cycle and its corresponding false lock state were analysed numerically and experimentally. An algorithm was given for computing the limit cycle, its period and the bifurcation point 6 ,. An experimental method was discussed for determining 6,. The computed value of 6, was within 1% of the measured result. When in the stable false lock state, the PLL constructed for the laboratory experiment consistently phase locked when its gain was increased through the value 6 ,. Saddle node bifurcation was shown to occur at 6,. This was accomplished by using a Poincare map to analyse solutions near the nonhyperbolic limit cycle at 6 ,. As is characteristic of such bifurcations, an unstable limit cycle bifurcates at 6, in addition to the above-mentioned stable limit cycle. The unstable limit cycle does not exist for 6 > 6,. As a function of 6, the unstable limit cycle cannot be continued down to 6 = 0. That is, it fails to exist on an interval 0 < 6 < b2, where d2 is positive. Furthermore, it is believed that the unstable limit cycle bifurcates from a homoclinic loop (also known as a separatrix cycle) at 6,. More work is required to establish this conjecture and extend the results contained in this paper to higher-order PLLs.

References (1) F. M. Gardner, “Phaselock Techniques”, 2nd Ed., John Wiley, New York, 1979. (2) W. C. Lindsey, “Synchronization Systems in Communication and Control”, Prentice Hail, Englewood Cliffs, NJ, 1972. (3) A. J. Viterbi, “Principles of Coherent Communications”, McGraw-Hill, New York, 1966. (4) J. L. Stensby, “False lock and bifurcation in the phase locked loop”, SIAM J. A@. Math., Vol. 47, No. 6, pp. 1177-1184, 1987. (5) E. A. Coddington and N. Levinson, “Theory of Ordinary Differential Equations”, McGraw-Hill, New York, 1955. (6) B. D. Hassard and Y. H. Wan, “Theory and Application of Hopf Bifurcation”, Cambridge University Press, Cambridge, 198 1. (7) A. A. Andronov, E. A. Leontovich, I. I. Gordon and A. G. Maier (translated by D. Louvish), “Theory of Bifurcations of Dynamic Systems on a Plane”, John Wiley, New York, 1973. (8) L. Perko, “Differential Equations and Dynamical Systems”, Springer, New York, 1991. (9) D. Ruelle, “Elements of Differentiable Dynamics and Bifurcation Theory”, Academic Press, New York, 1989. Received : 15 February 1993 Accepted : 30 April 1993

Journal

786

OC the Frankim Institute Pergamon Press Ltd