Matrix analysis of linear conservative systems

MATRIX ANALYSIS OF LINEAR CONSERVATIVE SYSTEMS BY

L O U I S A. P I P E S t ABSTRACT

The loop differential equations of any lossless, reciprocal electric circuit are integrated by the use of functions of matrices. It is shown how the transient analysis of such circuits is a generalization of the analysis of the simple single-loop circuit composed of a linear inductance in series with a linear elastance. The scaler trigonometric functions that appear in the analysis of the single-loop circuit are replaced by matric trigonometric functions in the analysis of the general circuit of n loops. The matric functions are then expanded by a form of Sylvester's theorem. T h e case in which the circuit has multiple eigenvalues is particularly easy to analyze by the use of functions of matrices. The analysis is extended to the forced oscillations of the n-loop linear conservative circuit. I. THE CLASSICAL THEORY OF LINEAR CONSERVATIVE CIRCUITS

Before discussing the free and forced oscillations of the general n-loop linear lossless circuit, let the simple single-loop circuit of Fig. 1 be considered. The circuit depicted in Fig. 1 consists of a linear inductance, L, connected in series with a linear elastance S. The parameters L and S are in series with a driving electromotive force E ( t ) . If the circulating charge in this circuit is q, the circulating current i is given by dq i = - = 4.

(1)

dt

The differential equation that determines the circulating charge q ( t ) of the circuit is Lii + Sq = E(t) .

(2)

The magnetic energy T of the circuit is T = a_ L 2.

2

(3)

The electric energy of the circuit is V = -1- S q z. 2

(4)

! Professor of Engineering, University' of California, Los Angeles, Calif.; Consultant, Aerospace Corporation. 3OO

Oet. 1962. ]

M A T R I X ANALYSIS OF LINEAR CONSERVATIVE SYSTEMS

3Ol

If there is no electromotive force applied to the circuit, the E ( t ) = 0 and Eq. 2 reduces to L~ + Sq = 0. This is the differential equation for the free oscillations of the circuit. the following notation is introduced S

w ~ = -L

(5) If

(6)

then the differential Eq. 5 for the free oscillations of the circuit takes the form,

q + W2q = O.

(7)

L

D J-, T Fro. 1. Single-loop L-C circuit

If q0 is the initial charge at t = 0 on the elastance and i0 is the initial current that flows in the circuit at t = 0, then the solution of Eq. 7 subject to these initial values of charge and current is i0 q = q0cos(wt) + - - sin(wt). w

(8)

After this brief review of well-known results for the simple single-loop circuit, let us now consider a general linear conservative circuit. Given any such circuit, there are well-established rules for choosing a set of independent loops and for writing the loop equations (1).2 Let the n circulating loop charges of the general n-loop circuit be q~, q2, q3, ' ' ' , q~" These loop charges m a y be taken as the elements of an n th order column matrix (q) of the form, ql

q2 (q) =

q3

"'

2 The boldface numbers in parentheses refer to the references appended to this paper.

(9)

3o2

Louis A. PIPES

[,I.F.I.

T h e column matrix (q) is the charge matrix of the circuit; the time derivative of the matrix (q) is the current matrix of the circuit (i), so that, (i) = (4).

(10)

T h e elements of the current matrix (i) are the various circulating currents in the n loops of the general lossless circuit. Let (E) be a column matrix whose elements E~ (t), E2(t), E3(t), . . . , Eo(t), are the n electromotive forces applied to the various loops of the circuit; then in the above matrix notation (2), the loop equations of the circuit may be written as a single matrix differential equation of the form, ILl (q) + IS] (q) = (E).

(11)

T h e matrix [L] is a symmetric matrix whose elements are the self and mutual inductances of the n-loop circuit. This matrix is the inductance matrix of the circuits. T h e matrix [S] is a symmetric matrix whose elements are the self and mutual elastances of the circuit. [S] is the elastance matrix of the circuit. Equation 11 m a y be regarded as a generalization of the simple scalar Eq. 2 for the simple single-loop circuit. It is easy to show (2) that the magnetic energy T of the n-loop circuit may be expressed in the compact form, T=

1 ~- (q)' [L] (q)

(12)

where (4)' is the transpose of the matrix (4). Similarly, the electric energy V of the general circuit m a y be expressed in the form V = _1 (q)'[S] (q). 2

(13)

If there are no electromotive forces acting on the circuit, then (E) is a zero matrix and Eq. 11 reduces to [L](q) + IS] (q) = (0).

(14)

This matrix differential equation governs the behavior of the free oscillations of the general n-loop lossless linear circuit. Equation 14 is a generalization of the simple scalar equation (5) for the free oscillation of the single-loop case. II. THE FREE OSCILLATIONS OF THE N-LOOP CIRCUIT

T h e study of the free oscillations of the general n-loop lossless linear circuit subject to aft initial distribution of charges and currents reduces to the solution of the matrix differential equation (14) subject to the following initial conditions: (q) = (q)o} att = 0. (4) (i)o

(15)

Oct. ,962.] ]VIATRIX ANALYSIS OF LINEAR CONSERVATIVE SYSTEMS

30:3

In (15), (q),, is a c o l u m n matrix of the initial loop charges and (i)~, is a colu m n matrix of the initial loop currents. In general, the i n d u c t a n c e matrix [L] is non-singular a n d therefore, its inverse exists. If Eq. 14 is pre-multiplied by [L] ~, the resulting equation is (16)

+ [L]-, IS] (9) = (0). Now let [L] -I [S] = [o2]2.

(17)

With this notation, Eq. 16 takes the form (q) + [02]2(q) = (0).

(18)

This equation is analogous to the simple scalar equation (8). T h e solution of Eq. 18 subject to the initial conditions (15) m a y be written in terms of functions of matrices in the form, (q) = cos([0211) (q)0 + [02]-1 sin([02]t) (i)0-

(19)

This result m a y be regarded as a generalization of the solution (Eq. 8) for the single-loop case. T h e scalar trigonometric functions in the solution of the single-loop case are replaced by matric trigonometric functions in the solution (Eq. 19) for the general n-loop case. T h e matric functions cos ([02]t) and sin ([02]t) are defined (2) by the following convergent matric power series: ([02]/)2

cos ([02]0 = u sin ([02]t) -

2!

([02]/)4 + . . . . . 4!

([021t)

([021t)

1[

3!

+

([021tp

,

(20) (21)

5!

If Eqs. 20 and 21 are differentiated with respect to the scalar variable t, it is easy to obtain the following results, d -

-

cos ([02]0 = -[w] sin ([0211)

(22)

sin([02]t) = [02] cos ([021t).

(23)

dt d -

-

dt

If Eqs. 22 a n d 23 are differentiated again with respect to the scalar variable t, the following results are obtained: d2 - - cos ([02]t) = _[02]2 cos ([wlt),

(24)

dt 2 d2

- - sin[02]t) = --[wp sin ([021t).

dt z

(25)

304

Louis A. PIPES

[J.V.l.

As a consequence of Eqs. 23 and 24, it is evident that Eq. 19 is a solution of Eq. 18 subject to the initial conditions (Eq. 15). III. EXPANSION OF THE MATRIC FUNCTION SOLUTION

Although the solution (Eq. 19) is an elegant and formally complete solution, it is not well adapted for practical computations. If the solution is desired for small values of the time t, then the matric trigonometric functions may be c o m p u t e d by means of the series (20) and (21). However, this procedure is impractical for large values of time. T h e practical procedure which is valid for all values of the time t > 0, and that gives an insight into the physical behavior of the circuit to expand the solution (Eq. 19) in terms of the eigenvectors and eigenvalues of the n-loop circuit. In order to do this, we first turn our attention to the expressions for the magnetic and electric energies T a n d V given by Eqs. 12 and 13. These energies are given by the two quadratic forms, T : _1 (~)'[L] (q) 2

(26)

V = _1 (q), [S] (q). 2

(27)

T h e inductance matrix [L] and the elastance matrix [S] are symmetric matrices. T is a positive-definite quadratic form and V is non-negative. It is a well-known result (2), that if [L] is definite, an n'h order square matrix [HI m a y be found such that simultaneously [HI'[L] [H] = u

(28)

[HI' [S] [/4] : [d].

(29)

In Eq. 28, u is the n'h order unit matrix and the matrix [d] is a diagonal matrix that has the form . , 0

[d]

I

0

=

.

#z

0

0

0

0

0

0

0

0

0/

. . . . . . . . . . . . . . . . . . . . . . . . .

L 0

(30)

! 0

0

0

0

0

.oA

#1, #2, #3, -.., gn are the eigenvalues of the matrix ([L] -1 [S]) = [co]2. If Eq. 28 is post-multiplied by [H] -1 and the resulting equation by [L] 1, the following result is obtained: [H]' = [H1-1 [L] -1.

(31)

If now Eq. 31 is substituted into Eq. 29, the following equation is obtained: [H] -1 ]L] -1 [S] [/4] = [H] -1 [co]z [H] = [d].

(32)

Oct. 1962.1 ~IATRIX ANALYSIS OF LINEAR CONSERVATIVE SYSTEMS

305

The columns of the square matrix [H] are the normalized eigenvectors of the matrix [,0]2. It is easy to show (2) that the eigenvalues #i in Eq. 30 are the squares of the natural angular frequencies of the n-loop electric circuit, .02, that is, it, = *oz. Let [V] = [L]-'[S] (33) so that

[V]'/2.

[*0] =

(34)

With this notation Eq. 32 may be written in the form .

L 0

0

.

.

0

.

.

.

...

.0~

Since the matrix [H] reduces [ V] to the diagonal form (35) by the transformation [H]-~[V]]H], this same transformation will also reduce any function of[ V], F([ V]) to the following diagonal form [H]-' [F(*0)] [H] = [/4] ' F([ V)v2] [H] = IF(,0,)

0

...

0

-]

0

F(,0:,)

...

0

t

0

0

1

... F(,0,)

(36) J

[H] 'cos ([,0)) H] =

co,;t> ' - - ' " ; ............ ; ......... [[[ .... co'si;[llii

= [alt.]

(37)

= [B(t)].

(38)

[H]-' [.0]-i sin ([,0]t) [/4] = s i n (.01t)

0

...

0

sin (*0Q)

...

0

0.)1

0

/

*02 0 __

0

...

sin (*0,,t)_ COn

__J

306

Louis A. PIPES

[.I.F.I.

The matrix functions cos ([w]t) and [60] 1 sin ([oat]) may therefore be expressed in the form,

[H] [A(t)] [H]-'

(39)

[00] 1sin ([co]t) = [H] [B(t)] [H]-'.

(40)

cos ([oa]t) =

IfEqs. 39 and 40 are now substituted into the solution (Eq. 19) for the free oscillation of the n-loop circuit, we see that it may be expressed in the form (q) = [HI [A(t)] [H]-' (x)o + [H] [B(t)] [H]

1

(i)0 '

(41)

This is the classical solution. The eigenvalues ~ of Eq. 30 are the roots of the characteristic equation of the matrix [V] = [L] -1 [S], that is, they are the roots of the equation det (~u - [V] = 0.

(42)

To each eigenvalue #,, there corresponds the square of one of the natural frequencies of the circuit, oa, so that 0) 2, = t{.fi, i =

1, 2, 3 , . . . , n.

(43)

The columns of the transformation matrix [HI are the normalized eigenvectors corresponding to the ratios in the columns (H), that satisfy the equations, iV] (H), = #i(H)i, i = 1, 2, 3, . . . , n

(44)

and the normalizing condition, (H)l(H)i

=

1,i = 1 , 2 , . . . , n .

(45)

IV. THE DIRECT EVALUATION OF THE MATRIC FUNCTIONS

The matric functions cos ([ [q'/21) and [ V]1/2sin ([ V]'/2t) m a y be evaluated directly. If the matrix [V] = [L] '[S] has n distinct eigenvalues#~, #2,#3, --., #., a concise expression for the expansion of a matrix function, F ([ V] may be formulated in a symbolic determinant form (4). If, for example, F([ V]) is a matrix of the fourth order whose eigenvalues are #1, #2, tt3, and #4, supposed all distinct, the function F([ V]) may be obtained by expanding the following determinant, 1

1

1

1

u

#1

~-L2

].L3

~.L4

[ V]

#2

#23

tt2

[ V]2

F(/*2)

F(t*3)

u~

F(u,)

F(ttn)F([V])

= 0

(46)

Oct. ,962.]

M A T R I X ANALYSIS OF LINEAR CONSERVATIVE SYSTEMS

307

In this d e t e r m i n a n t , u is the fourth order unit m a t r i x a n d the desired expression for F([ V]) follows by e x p a n d i n g the d e t e r m i n a n t (46) in terms of the last row. It is easily seen t h a t {(IV]) - #~u} {([v] - u2u} . . . {([v] /~4u} is a factor of the above d e t e r m i n a n t a n d hence by the C a y l e y - H a m i l ton theorem, this is the zero matrix. If, for example, the m a t r i x [V] = [L] ~[S] is of the third order with eigenvalues #~

=

w~,#~

w2,#~

=

=

co 23

(47)

t h e n by e x p a n d i n g a d e t e r m i n a n t similar to (46), it is easy to show that

cos([V],/2t) = ([V] - o02u)([V] - 6o2~u) cos

(w,t)

(04 - 0 4 ) ( ~ 4 - oo~) +

+

([ v]

-

o~2u) ([

([ V]

-

L

v]

~4") ([ v]

- o)~)

cos ( % 1 )

- ,4u)

L

cos (w3t).

(48)

L

,7 , 7 ,T Fro. 2.

Three-loop circuit.

A similar expression is o b t a i n e d for the function IV] 1/2 sin ([1/]1/2t); the m a t r i x coefficients for this function are the same as those for cos ([ V] 1/21), but the functions cos (oait) are replaced by sin (w,t)/oo,, i = 1, 2, 3. As an example of the general procedure, consider the oscillations of the three loop circuit of Fig. 2. T h e i n d u c t a n c e m a t r i x [L] of this circuit has the form

[L] = L

1 __

(49)

= Lu

0

where u is the third order unit matrix. T h e elastance m a t r i x of the circuit of Fig. 2 has the form

[S] =

2Io 2S

-S

-S

2S_J

, = S

[i, i] 2

-

1

-

(50)

308

L o u I s A. PIPES

[J.F.I.

In this case, the matrix [ V] has the form _2 S [V] = [L]-'[S] = ~-

-1

1

2

0 = 002

I_ 21

- 12

0

-1

-1

- il

(51)

~000= ~ '

The eigenvalues ofg~,/*2, ~3 of[V] are u, = ( 2 -

v~4

#2 = 2°a2 = °a2

(52)

u3 = (2 + v~5o~2 = , 4 .

The expression for cos ([ V]l/2t) may be obtained by the use of Eq. 48 and written in the form

cos([V]'/2t)

= [M]lcos(oa,/) + [M]2cos(60it) -4- [M]3cis(oa3t)

(53)

where the matrices [M],, [M]2 , and [M]3 have the following form

[M]I = ~-

x/T

2

1

x/g-

1

1

0

-1-~

[M] 2 = ~-

0

0

0

0

1

-1

[M]3 The function

1 ~-

1

-~/Y-

vT-

2

1

-x/-2

(54)

-

[V]-l/2t) is given by the equation

[V] -1/2 sin

([V]l/Zt) = [M], -sin- (Oalt) 601

+ [M],

+

[M12

sin (w2t)

002

(55)

sin (o~3t)

603

If the initial loop charges and loop currents of the circuit at t = 0 are given by the vectors (q)0 and (i)o in the form

Oct.

,962.] MATRIX ANALYSIS OF LINEAR CONSERVATIVE SYSTEMS

/;,\

/q0\

(q)o= tqo~)

and

(i),~ = \,i.,;~)

3o9

(56)

\qo=/ then the subsequent oscillations of the circuit loop charges are given by (q) = cos ([ V]t/2t) (q)o + { V} -'/2 sin ([ V]t) (i) o

(57)

where the matric trigonometric functions have the expanded forms (53) and (54). V. THE CASE OF MULTIPLE EIGENVALUES OF THE MATRIX [V]

The case of coincident roots of the characteristic equation of the matrix [ V] so that the equation, det (#u - IV]) = 0

(58)

has multiple roots, of historic interest. In this case, some of the natural frequencies of the circuit now coincide. If the equations governing the oscillations of'the loop charges are of the form (~/) + {Vl (q) = (0),[Vl = [LI-'[S1

(59)

each loop charge q~ will satisfy an equation obtained by eliminating the remaining variables, which is simply, (D 2 + w~)(D 2 + w~) ... (D 2 -[- (.1)2) qi

=

0

(60)

where D 2 = dZ/dt 2 and where some of the w's are now assumed to coincide. The general solution m a y be expected to contain " s e c u l a r " terms of the form h(t) cos (~ot) corresponding to a repeated frequency w, h(t) and k(t) being polynomials of degree one less than the multiplicity of the root in addition to the usual terms in cos (oat) and sin (oat) for the single roots. This is indeed what Lagrange and others who followed him supposed. However, it turns out, that on substituting the assumed solution into the original differential equation (59), in accordance with the usual procedure for determining the relations between the constants, all the coefficients in h(t) and k(t) vanish identically with the exception of their constant terms. T h e reason for this is not directly obvious, although of course the matter has been well understood, since Weierstrass discussed it in 1865 (4). The matrix solution concisely demonstrates the disappearance of the " s e c u l a r " terms in the dynamical problem, as will now be shown. W h e n repeated roots occur, (q) = cos ([V]UZt) (q)o + [V] -'/2 sin ([ V]~12t) (i)o

(61)

is still formally the solution of the differential equation (59) with the appropriate initial conditions. However, since some of the eigenvalues/*i now coincide, symbolic determinants of the type (46) cannot be used to inter-

3 IO

LouIs

A . PIPES

[J.F.I.

pret the matrix functions, but a confluent form of Sylvester's identity must now be employed. Let it be assumed that the matrix [ V] is of the fourth order and that its eigenvalues are ]"gl

=

a, #2 = b, ]~3

=

b, ~4 = b.

(62)

T h e peculiarity of the dynamical case is that if the matrix [ V] has repeated eigenvalues, it satisfies an identical equation, that is of lower degree t h a n the order of [ V]. This equation is called the reduced Cayley-Hamilton equation, and it is obtained by counting each linear factor in the ordinary Cayley-Hamilton equation once only. In the above example, the ordinary Cayley-Hamilton equation has the form ([V] - au)([V] - bu) 3 where u is the fourth order unit matrix. equation is, in this case,

=

[0],

T h e reduced Cayley-Hamilton

(IV] - au)([V] - bu) = [0].

FIa. 3.

(63)

(64)

Symmetric circuit.

T h e fact that [V] satisfies the reduced equation follows at once from the fact that [V] can be reduced to the diagonal form whether it contains repeated eigenvalues or not. F r o m a knowledge of the reduced CayleyH a m i l t o n equation, an analytic function F([ V]) of the matrix [ V] m a y be expanded in the form F([V] = Cou + C~[V],

(65)

where u is the fourth order unit matrix, Co and C~ are constants that m a y be obtained from the equations F(a) = Co + Cla F(b) = Co + C~b.

(66)

Oct. i962.1

M A T R I X ANALYSIS OF LINEAR CONSERVATIVE SYSTEMS

31 I

IfEqs. 66 are solved for the constants Co and C1 and the results substituted into Eq. 65, then the matric function F([V]) may be written in the form F([V] = [ b F ( a ) - a F ( b ) ] u + [F(b) - F(a)] [V]. (b - a) (b - a)

(67)

VI. THE OSCILLATIONS OF A SYMMETRIC CIRCUIT

An example that illustrates the occurrence of a vanishing eigenvalue, a pair of repeated eigenvalues, and a single eigenvalue of the matrix IV] will now be given in order to clarify the ideas discussed in Section V. Consider the symmetric electric circuit of Fig. 3. The circuit of Fig. 3 consists of four equal inductances L and four equal elastances in the symmetric arrangement shown in the figure. The generalized coordinates of the system may be taken as the four circulating charges ql, q2, q3, and q4 of the internal loops of the circuit as shown. In this case, the inductance matrix [L] has the form [L] = L u ,

(68)

where u is the fourth order unit matrix. The elastance matrix [S] has the symmetric form,

2

[s]=s

-1

-1 0

2 -

, =S[B] -~_~

(69)

1

where [B] is the fourth order numerical matrix in (69). The matrix [ V] in this case is given by S IV] = [L]-'[S] = £ [B].

(7o)

The characteristic equation of the matrix [V] is det(#u-

[V]) = 0.

(71)

If we introduce the notation S L,

002 = b

(72)

[v] = ,4o [B].

(73)

then

The characteristic equation (71) may be written in the form

312

LouIs A. PIPES =

[J.F.I.

(74)

0

If we let

¢ = L

(75)

then Eq. 74 m a y be written in the form det (~bu - [B]) --- 0.

(76)

This is the characteristic e q u a t i o n of the numerical matrix [B]. If 76 is expanded, it will be found that the resulting polynomial in q~ can be factored and the eigenvalues of [B] are found to be ~b, = 0, ~b2 = 2, 4)3 = 2, 4)4 = 4.

(77)

As a c o n s e q u e n c e of the relation (75), these eigenvalues of [B] corlespond to the following eigenvalues of [ V]. ~1

:

0 , #2

=

2oa02,# 3

:

2~0~, # 4

:

4'o 2.

(78)

T h e differential e q u a t i o n for the free oscillations of the circuit has the form (/j) + [V] (q) = 0.

(79)

Since [ V] = oa2 [B], this e q u a t i o n m a y be written in the form,

if) + ~ [B] (q) = (0).

(80)

T h e solution of this matrix differential e q u a t i o n m a y now be written in the form,

(q) = cos(oaot[B] '/2) (q)o + ~o' [B]-'/2sin(~%t[B] '/2) (i)o.

(81)

T h e vectors (q)0 and (i)o have as elements the initial mesh charges and mesh currents in the circuit at t = 0. Since the matrix [B] has the eigenvalues (77), it has a reduced C a y l e y - H a m i l t o n e q u a t i o n of the form [B] ([B] - 2u)([B] - 4u) = [0],

(82)

where u is the fourth order unit matrix. If the same technique described in Section V is now used to c o m p u t e a matric function of [B], F([B]), the result is 1

F([B] = -~ ([B] - 2u) ([B] - 4u) F ( 0 ) -

1 -- [B] ([B]

4

-

4u) F(2)

+ - - 1 [B] ([B] - 2u) F ( 4 ) . 8

(83)

Oct. 1962.] MATRIX ANALYSIS OF LINEAR CONSERVATIVE SYSTEMS

313

This result may now be used to expand the two matrix functions of the right member of Eq. 81 to obtain (q)

=

1 [C1 +

+ 1__ [Clt +

4

2C2 cos

(x/2coot) + C~cos (2coot)] (q)o

-2C2 - Sm(v~ot ) + C~sin (2co0t)] (i),, 2o)(} 2co(}

(84)

where C1, ~ and C3are

=

1

1

1

1

1

1

1

1

1

1

1

1

1

1

1/ 1_.[__J

I I

C2=

1 0 -1

0

0

I

-1

C~ =

0 -1 1 0 -1

1 -1

1

1 -1 -1

1

-1

0t -1

1

0

0

1

(85)

1 -1 1 1

1 -1 1

1

The ratios of the coordinates in the various vibrational modes can be read off from the elements of the matrices C1, C2, C~, to give

(H) I =

1

1

' co = 0, (zero frequency mode)

(86)

\1

01or/0 (10/\°1 )

, degenerate modes of frequency, co = 2co (87)

314

(H)4 =

f 1) 1 1

Louis A. PIPES

[J.F.I.

f r e q u e n c y o~ -- 29.

(88)

1 VII. THE FORCED OSCILLATIONS OF THE CIRCUIT

It is possible to write the solution of the equation ILl (~) + IS] (q) = {E(t)}

(89)

by the use of functions of matrices. In order to do this, we pre-multiply Eq. 89 by the inverse of the inductance matrix and obtain the equation (~/) + IV] (q) = IF(t)}

(90)

[V] = [L]-' [S] and{F(t)} = [L] -I {E(t)}.

(91)

where

T h e solution of Eq. 90, subject to the initial conditions, (q) = (q)0 a t t = 0, and (~) = (i)0att = 0,

(92)

m a y be written in the form (q) = cos ([ V]llZt) (q)o + ([ V] -t sin ([ V]l/2t) (i)o + ([V])-'

£

(93)

sinl[V] '/z (t - u)} { f ( u ) l d u REFERENCES

(1) E. A. GUILLEmN,"Introductory Circuit Theory," New York, John Wiley & Sons, Inc., 1953. (2) E. A. GUILLEMIN, "The Mathematics of Circuit Analysis," New York, John Wiley & Sons, Inc., 1949. (3) L. A. PIPES, "Applied Mathematics for Engineers and Physicists," New York, McGraw-Hill Book Co., Inc., 1958, Chapter VII. (4) H. W. TURNBULLAND A. C. AITKEN,"An Introduction to the Theory of Canonical Matrices," London, Blackie and Son, Limited, 1932, pp. 76-79.