'***************************************************** '* Calculate the determinant of a real square matrix * '* by the LU decomposition method * '* ------------------------------------------------- * '* SAMPLE RUN: * '* (Calculate the determinant of real square matrix: * '* | 1 0.5 0 0 0 0 0 0 0 | * '* |-1 0.5 1 0 0 0 0 0 0 | * '* | 0 -1 0 2.5 0 0 0 0 0 | * '* | 0 0 -1 -1.5 4 0 0 0 0 | * '* | 0 0 0 -1 -3 5.5 0 0 0 | * '* | 0 0 0 0 -1 -4.5 7 0 0 | * '* | 0 0 0 0 0 -1 -6 8.5 0 | * '* | 0 0 0 0 0 0 -1 -7.5 10 | * '* | 0 0 0 0 0 0 0 -1 -9 | ) * '* * '* * '* Determinant = 1 * '* * '* * '* Basic version by J-P Moreau, Paris. * '* (www.jpmoreau.fr) * '***************************************************** ' Exact value is: 1 '----------------------------------------------------- DEFINT I-N DEFDBL A-H, O-Z CLS : PRINT n = 9 DIM A(n, n) DIM INDX(n) 'matrix A(i,j) column by column DATA 1,-1,0,0,0,0,0,0,0 DATA 0.5,0.5,-1,0,0,0,0,0,0 DATA 0,1,0,-1,0,0,0,0,0 DATA 0,0,2.5,-1.5,-1,0,0,0,0 DATA 0,0,0,4,-3,-1,0,0,0 DATA 0,0,0,0,5.5,-4.5,-1,0,0 DATA 0,0,0,0,0,7,-6,-1,0 DATA 0,0,0,0,0,0,8.5,-7.5,-1 DATA 0,0,0,0,0,0,0,10,-9 'read matrix A(i,j) FOR j = 1 TO n FOR i = 1 TO n READ A(i, j) NEXT i NEXT j GOSUB 1000 'call LUDCMP(A,N,INDX,ID,CODE) FOR j = 1 TO n D = D * A(j, j) NEXT j PRINT PRINT " Determinant = "; D PRINT END '*************************************************************** '* Given an N x N matrix A, this routine replaces it by the LU * '* decomposition of a rowwise permutation of itself. A and N * '* are input. INDX is an output vector which records the row * '* permutation effected by the partial pivoting; D is output * '* as -1 or 1, depending on whether the number of row inter- * '* changes was even or odd, respectively. This routine is used * '* in combination with LUBKSB to solve linear equations or to * '* invert a matrix. Return code is 1, if matrix is singular. * '*************************************************************** 1000 'PROCEDURE LUDCMP; NMAX = 100: TINY = 1.5E-16 DIM VV(n) D = 1: CODE = 0 FOR i = 1 TO n AMAX = 0# FOR j = 1 TO n IF ABS(A(i, j)) > AMAX THEN AMAX = ABS(A(i, j)) NEXT j IF (AMAX < TINY) THEN CODE = 1 RETURN END IF VV(i) = 1# / AMAX NEXT i 'i loop FOR j = 1 TO n FOR i = 1 TO j - 1 SUM = A(i, j) FOR K = 1 TO i - 1 SUM = SUM - A(i, K) * A(K, j) NEXT K A(i, j) = SUM NEXT i 'i loop AMAX = 0# FOR i = j TO n SUM = A(i, j) FOR K = 1 TO j - 1 SUM = SUM - A(i, K) * A(K, j) NEXT K A(i, j) = SUM DUM = VV(i) * ABS(SUM) IF DUM >= AMAX THEN IMAX = i AMAX = DUM END IF NEXT i 'i loop IF j <> IMAX THEN FOR K = 1 TO n DUM = A(IMAX, K) A(IMAX, K) = A(j, K) A(j, K) = DUM NEXT K 'k loop D = -D VV(IMAX) = VV(j) END IF INDX(j) = IMAX IF ABS(A(j, j)) < TINY THEN A(j, j) = TINY END IF IF j <> n THEN DUM = 1# / A(j, j) FOR i = j + 1 TO n A(i, j) = A(i, j) * DUM NEXT i END IF NEXT j 'j loop RETURN 'subroutine LUDCMP 'End of file deter1.bas