Geology Reference
In-Depth Information
N2=2*N+1
NP1=N+1
NM1=N-1
C Insert the matrix A.
DO 10 I=1,N
DO 10 J=1,N
C(I,J)=A(I,J)
10 CONTINUE
C Insert the unknown vector.
DO 11 J=1,N
C(J,NP1)=B(J)
11 CONTINUE
C Augment with the identity matrix.
DO 12 I=1,N
DO 12 J=1,N
JJ=J+NP1
C(I,JJ)=0.D0
IF(I.EQ.J)C(I,JJ)=1.D0
12 CONTINUE
C Begin Gauss-Jordan elimination.
DET=1.D0
DO 13 I=1,N
IM1=I-1
C Search for largest pivot.
PIVOT=0.D0
DO 14 J=1,N
DO 14 K=1,N
C Scan IR and IC vectors for invalid pivot subscripts.
IF(I.EQ.1)GO TO 15
DO 16 ISC=1,IM1
DO 16 JSC=1,IM1
IF(J.EQ.IR(ISC))GO TO 14
IF(K.EQ.IC(JSC))GO TO 14
16 CONTINUE
15 IF(DABS(C(J,K)).LE.DABS(PIVOT))GO TO 14
PIVOT=C(J,K)
IR(I)=J
IC(I)=K
14 CONTINUE
DET=DET*PIVOT
C Label next row.
IRP1=IR(I)+1
C Label previous row.
IRM1=IR(I)-1
C Divide row by pivot.
DO 17 J=1,N2
C(IR(I),J)=C(IR(I),J)/PIVOT
17 CONTINUE
C Reduce elements below pivot to zero.
IF(IRP1.GT.N)GO TO 18
DO 19 K=IRP1,N
ELL=C(K,IC(I))
DO 20 J=1,N2
C(K,J)=C(K,J)-C(IR(I),J)*ELL
20 CONTINUE
19 CONTINUE
18 CONTINUE
C Reduce elements above pivot to zero.
IF(IRM1.LT.1)GO TO 21
DO 22 K=1,IRM1
ELL=C(K,IC(I))
Search WWH ::




Custom Search