DECLARE SUB CMUL (Z1!(), Z2!(), Z!()) DECLARE SUB RSHCGT (eps!, N%, R0 AS INTEGER, M0%) DECLARE FUNCTION CABS! (Z!()) DECLARE SUB CDIV (Z1!(), Z2!(), Z!()) '************************************************************** '* Solving a complex homogeneous linear system by Gauss * '* Method with full pivoting. * '* ---------------------------------------------------------- * '* SAMPLE RUN: * '* | -1 i -i | * '* Solve AX = 0 with A = | -i -1 1 | * '* | i 1 -1 | * '* * '* Basic Family of system solutions: * '* Solution 1 * '* X 1 = ( 0 , 1 ) * '* X 2 = ( 1 , 0 ) * '* X 3 = ( 0 , 0 ) * '* * '* Solution 2 * '* X 1 = ( 0 ,-1 ) * '* X 2 = ( 0 , 0 ) * '* X 3 = ( 1 , 0 ) * '* * '* ---------------------------------------------------------- * '* Ref.: "Algèbre, Algorithmes et programmes en Pascal * '* By Jean-Louis Jardrin, DUNOD Paris, 1988". * '* * '* Quick Basic Release By J-P Moreau, Paris. * '* (www.jpmoreau.fr) * '*************************************************************} 'Program TRSHCGT DEFINT I-N OPTION BASE 1 DIM SHARED NMAX NMAX = 20 DIM R0 AS INTEGER 'Complex matriczs qhared with main subroutine RSHCGT DIM SHARED A(NMAX, NMAX, 2), VX(NMAX, NMAX, 2) N = 3 'size of complex matrix A A(1, 1, 1) = -1!: A(1, 1, 2) = 0! A(1, 2, 1) = 0!: A(1, 2, 2) = 1! A(1, 3, 1) = 0!: A(1, 3, 2) = -1! A(2, 1, 1) = 0!: A(2, 1, 2) = -1! A(2, 2, 1) = -1!: A(2, 2, 2) = 0! A(2, 3, 1) = 1!: A(2, 3, 2) = 0! A(3, 1, 1) = 0!: A(3, 1, 2) = 1! A(3, 2, 1) = 1!: A(3, 2, 2) = 0! A(3, 3, 1) = -1!: A(3, 3, 2) = 0! eps = 1E-10 RSHCGT eps, N, R0, M0 CLS PRINT IF R0 = N THEN PRINT " System solution:" FOR i = 1 TO N PRINT " X"; i; " = ("; VX(i, 1, 1); ","; VX(i, 1, 2); ")" NEXT i ELSE PRINT " Basic Family of system solutions:" FOR l = 1 TO M0 PRINT " Solution "; l FOR i = 1 TO N PRINT " X"; i; " = ("; VX(i, l, 1); ","; VX(i, l, 2); ")" NEXT i NEXT l END IF PRINT END 'of main program FUNCTION CABS (Z()) CABS = SQR(Z(1) ^ 2 + Z(2) ^ 2) END FUNCTION SUB CADD (Z1(), Z2(), Z()) Z(0) = Z1(0) + Z2(0) Z(1) = Z1(1) + Z2(1) END SUB SUB CDIF (Z1(), Z2(), Z()) Z(0) = Z1(0) - Z2(0) Z(1) = Z1(1) - Z2(1) END SUB SUB CDIV (Z1(), Z2(), Z()) DIM C1(2) D = Z2(1) ^ 2 + Z2(2) ^ 2 IF D > 2E-12 THEN C1(1) = Z2(1): C1(2) = -Z2(2) CMUL Z1(), C1(), Z() Z(1) = Z(1) / D: Z(2) = Z(2) / D END IF END SUB SUB CMUL (Z1(), Z2(), Z()) Z(1) = Z1(1) * Z2(1) - Z1(2) * Z2(2) Z(2) = Z1(1) * Z2(2) + Z1(2) * Z2(1) END SUB SUB CSwap (Z1(), Z2()) DIM T(2) T(1) = Z2(1): T(2) = Z2(2) Z2(1) = Z1(1): Z2(2) = Z1(2) Z1(1) = T(1): Z1(2) = T(2) END SUB '********************************************************** '* Procedure RSHCGT solves the complex homogeneous linear * '* system AX = 0 by Gauss method with full pivoting. * '* ------------------------------------------------------ * '* INPUTS: * '* eps: required precision * '* N : size of complex matrix A (NxN) * '* A : Complex Matrix shared with main program * '* OUTPUTS: * '* R0: rank of A in output * '* M0: dimension of system solutions space * '* (if R0 <> N) * '* VX: contains in output the M0 solution * '* vectors (stored in columns 1..M0) when * '* R0<>N, or the unique solution vector * '* (stored in column 1) when R0=N. * ''* Shared with mainprogram. * '********************************************************** SUB RSHCGT (eps, N, R0 AS INTEGER, M0) DIM C0(2), C1(2), P0(2), Q0(2), S(2), TMP(2), TMP1(2)'Complex numbers DIM Z0(2), Z1(2) DIM I0(NMAX) Z0(1) = 0!: Z0(2) = 0! Z1(1) = 1!: Z1(2) = 0! FOR K = 1 TO N I0(K) = K NEXT K R0 = 0: K = 1 WHILE R0 = 0 AND K <= N P0(1) = A(K, K, 1): P0(2) = A(K, K, 2) l0 = K: k0 = K FOR i = K TO N FOR j = K TO N TMP(1) = A(i, j, 1): TMP(2) = A(i, j, 2) IF CABS(TMP()) > CABS(P0()) THEN P0(1) = A(i, j, 1): P0(2) = A(i, j, 2) l0 = i: k0 = j END IF NEXT j NEXT i IF CABS(P0()) < eps THEN R0 = K - 1 M0 = N - R0 ELSE IF K = N THEN R0 = N: M0 = 0 ELSE IF l0 <> K THEN FOR j = K TO N 'CSwap(A(k,j),A(l0,j)) TMP(1) = A(l0, j, 1): TMP(2) = A(l0, j, 2) A(l0, j, 1) = A(K, j, 1): A(l0, j, 2) = A(K, j, 2) A(K, j, 1) = TMP(1): A(K, j, 2) = TMP(2) NEXT j END IF IF k0 <> K THEN FOR i = 1 TO N 'CSwap(A(i,k),A(i,k0)) TMP(1) = A(i, k0, 1): TMP(2) = A(i, k0, 2) A(i, k0, 1) = A(i, K, 1): A(i, k0, 2) = A(i, K, 2) A(i, K, 1) = TMP(1): A(i, K, 2) = TMP(2) NEXT i END IF N0 = I0(K): I0(K) = I0(k0): I0(k0) = N0 END IF END IF FOR i = K + 1 TO N 'Q0 = A(i,k)/P0 TMP(1) = A(i, K, 1): TMP(2) = A(i, K, 2) CDIV TMP(), P0(), Q0() FOR j = K + 1 TO N C0(1) = A(i, j, 1): C0(2) = A(i, j, 2) TMP1(1) = A(K, j, 1): TMP1(2) = A(K, j, 2) 'C1 = Q0 * A(k,j) CMUL Q0(), TMP1(), C1() A(i, j, 1) = C0(1) - C1(1) A(i, j, 2) = C0(2) - C1(2) NEXT j NEXT i K = K + 1 WEND IF R0 = N THEN FOR i = 1 TO N VX(i, 1, 1) = Z0(1): VX(i, 1, 2) = Z0(2) NEXT i ELSE FOR l = 1 TO M0 FOR i = N TO R0 + 1 STEP -1 IF i = R0 + l THEN VX(i, l, 1) = Z1(1): VX(i, l, 2) = Z1(2) ELSE VX(i, l, 1) = Z0(1): VX(i, l, 2) = Z0(2) END IF NEXT i FOR i = R0 TO 1 STEP -1 S(1) = Z0(1): S(2) = Z0(2) FOR j = i + 1 TO N C0(1) = S(1): C0(2) = S(2) 'C1 = A(i,j) * VX(j,l) TMP(1) = A(i, j, 1): TMP(2) = A(i, j, 2) TMP1(1) = VX(j, l, 1): TMP1(2) = VX(j, l, 2) CMUL TMP(), TMP1(), C1() S(1) = C0(1) + C1(1) S(2) = C0(2) + C1(2) NEXT j 'C0 = S/A(i,i) TMP(1) = A(i, i, 1): TMP(2) = A(i, i, 2) CDIV S(), TMP(), C0() VX(i, l, 1) = -C0(1) VX(i, l, 2) = -C0(2) NEXT i NEXT l FOR i = 1 TO N - 1 j = i + 1 WHILE j <= N IF I0(j) = i THEN FOR l = 1 TO M0 'CSwap(VX(j,l),VX(i,l)) TMP(1) = VX(i, l, 1): TMP(2) = VX(i, l, 2) VX(i, l, 1) = VX(j, l, 1): VX(i, l, 2) = VX(j, l, 2) VX(j, l, 1) = TMP(1): VX(j, l, 2) = TMP(2) NEXT l I0(j) = I0(i): I0(i) = i j = N + 1 ELSE j = j + 1 END IF WEND NEXT i END IF END SUB 'RSHCGT