'********************************************************* '* Calculate the coefficients of the characteristic * '* polynomial P(l)=A-l*I of a real matrix A(i,j) * '* by Souriau's method * '* ----------------------------------------------------- * '* Ref.: "Algebre, Algorithmes et programmes en Pascal * '* By Jean-Louis Jardrin, DUNOD Paris, 1988" [10].* '* ----------------------------------------------------- * '* SAMPLE RUN: * '* (Find the caracteristic polynomial of matrix: * '* 1 0 0 0 0 1 * '* 1 1 0 0 0 -1 * '* A = -1 1 1 0 0 1 * '* 1 -1 1 1 0 -1 * '* -1 1 -1 1 1 1 * '* 1 -1 1 -1 1 -1 ) * '* * '* Input matrix size (max 100): 6 * '* Line 1: * '* Element 1: 1 * '* Element 2: 0 * '* Element 3: 0 * '* Element 4: 0 * '* Element 5: 0 * '* Element 6: 1 * '* Line 2: * '* Element 1: 1 * '* Element 2: 1 * '* Element 3: 0 * '* Element 4: 0 * '* Element 5: 0 * '* Element 6: -1 * '* Line 3: * '* Element 1: -1 * '* Element 2: 1 * '* Element 3: 1 * '* Element 4: 0 * '* Element 5: 0 * '* Element 6: 1 * '* Line 4: * '* Element 1: 1 * '* Element 2: -1 * '* Element 3: 1 * '* Element 4: 1 * '* Element 5: 0 * '* Element 6: -1 * '* Line 5: * '* Element 1: -1 * '* Element 2: 1 * '* Element 3: -1 * '* Element 4: 1 * '* Element 5: 1 * '* Element 6: 1 * '* Line 6: * '* Element 1: 1 * '* Element 2: -1 * '* Element 3: 1 * '* Element 4: -1 * '* Element 5: 1 * '* Element 6: -1 * '* * '* Chracteristic polynomial P(l): * '* Coefficient for degree 6: 1.000000 * '* Coefficient for degree 5: -4.000000 * '* Coefficient for degree 4: 0.000000 * '* Coefficient for degree 3: 30.000000 * '* Coefficient for degree 2: -75.000000 * '* Coefficient for degree 1: 79.000000 * '* Coefficient for degree 0: -32.000000 * '* * '* (The characteristic polynomial of matrix A is: * '* P(l)=l^6-4*l^5+30*l^3-75*l^2+79*l-32 ) * '* * '* Basic version by J-P Moreau, Paris * '* (www.jpmoreau.fr) * '********************************************************* 'Program CARPOL2 defdbl a-h,o-z defint i-n NMAX = 20 NMAXP = 21 ' k,n: Integer dim A(NMAX,NMAX) dim P(NMAXP) cls gosub 1000 'Read data ' call routine to calculate coefficients gosub 2000 'PCMS(n,A,P) ' print results print print " Characteristic polynomial P(l):" For k=1 to n+1 print " Coefficient for degree "; n-k+1; ": "; print using "####.######"; P(k) next k print END 'read matrix A from screen 1000 print print " Input matrix size (max ";NMAX;"): "; : input "", n for i=1 to n print : print " Line "; i for j=1 to n print " Element "; j; ": "; : input "", A(i,j) next j next i return 'calculate the trace of real matrix C 'Function TRM(n:integer; C:matrix):Double 1500 t0=0.0 for i=1 to n t0=t0+C(i,i) next i return '***************************************************** '* The procedure PCMS calculates the coefficients of * '* the characteristic polynomial P(l)=A-l*I of a real* '* square matrix A(i,j). * '* * '* Note: the roots of P(l) are the eigenvalues of * '* matrix A(i,j). * '* ------------------------------------------------- * '* INPUTS: * '* n: zize of matrix (integer) * '* A: real matrix (n,n) * '* (see main program). * '* OUTPUT: * '* P: vector of coefficients of P(l) * '* * '***************************************************** 2000 dim B(NMAX,NMAX), C(NMAX,NMAX) if INT(n/2) <> n/2 then P(1)=-1# else P(1)= 1# end if for l=1 to n if (l=1) then for i=1 to n for j=1 to n C(i,j)=A(i,j) next j next i else for i=1 to n for j=1 to n s=0# for k=1 to n s=s+B(i,k)*A(k,j) next k C(i,j)=s next j next i end if gosub 1500 't0=TRM(n,C)/l t0=t0/L P(l+1)=-t0*P(1) if l