! GAUSSIAN ELIMINATION DEMONSTRATION - SIMPLE VERSION ! PROGRAM P141 IMPLICIT NONE REAL :: M1(5,6),SOL(5) INTEGER :: N INTERFACE SUBROUTINE MATSIN(M1,N) IMPLICIT NONE REAL ,INTENT(IN OUT) :: M1(:,:) INTEGER ,INTENT(IN OUT) :: N END SUBROUTINE MATSIN SUBROUTINE GAUSS(M1,SOL,N) IMPLICIT NONE REAL ,INTENT(IN OUT) :: M1(:,:),SOL(:) INTEGER ,INTENT(IN OUT) :: N SUBROUTINE PRNMAT(M3,N) IMPLICIT NONE REAL ,INTENT(IN OUT) :: M3(:,:) INTEGER ,INTENT(IN OUT) :: N END SUBROUTINE PRNMAT END SUBROUTINE GAUSS END INTERFACE ! N=5 ! NUMBER OF EQUATIONS ! PRINT *,' This is Program P141 - Gaussian elimination' PRINT *,'PROGRAM IS READING DATA INTO ARRAYS' CALL MATSIN(M1,N) PRINT *,'SOLVING SYSTEM OF EQUATIONS' CALL GAUSS(M1,SOL,N) PRINT *,'SOLUTION:' PRINT 90,SOL 90 FORMAT(' | ',F8.3,' |') STOP END PROGRAM P141 SUBROUTINE MATSIN(M1,N) IMPLICIT NONE REAL ,INTENT(IN OUT) :: M1(:,:) INTEGER ,INTENT(IN OUT) :: N INTEGER :: I,J ! ! Tell program where data for READ is coming from OPEN(UNIT=5, FILE='P141.DAT') ! UNIT=5 is the default input ! ! READ IN M1 ! ONE ROW PER CARD ! L1: DO I=1,N READ 27,(M1(I,J),J=1,N+1) END DO L1 27 FORMAT(10(F5.2)) RETURN END SUBROUTINE MATSIN ! SUBROUTINE PRNMAT(M3,N) IMPLICIT NONE REAL ,INTENT(IN OUT) :: M3(:,:) INTEGER ,INTENT(IN OUT) :: N INTEGER :: I,J ! ! L3: DO I=1,N PRINT 202, (M3(I,J),J=1,N+1) END DO L3 202 FORMAT(10(' ',F8.3)) ! RETURN END SUBROUTINE PRNMAT ! SUBROUTINE GAUSS(M1,SOL,N) IMPLICIT NONE REAL ,INTENT(IN OUT) :: M1(:,:),SOL(:) INTEGER ,INTENT(IN OUT) :: N ! ! THIS ROUTINE PERFORMS GAUSSIAN ELIMINATION AND BACKSUBSTITUTION. ! WE HAVE SIMPLIFIED THE PROBLEM BY NOT WORRYING ABOUT MATRICES ! WITHOUT A SOLUTION. IF WE WANTED TO CONSIDER THAT CASE, ALL THE ! CODE THAT FOLLOWS WOULD STILL BE VALID, BUT WE WOULD HAVE TO ADD ! MORE TESTS. ! REAL :: M2(:,:),TEMP INTEGER :: I,J,K ! ! ! INSTEAD OF MODIFYING THE ORIGINAL ARRAY, WE WILL PRODUCE A WORKING COPY ! OF IT ! M2 = M1 ! L1: DO I=1,N ! L2: DO J=I+1,N TEMP=M2(J,I)/M2(I,I) L3: DO K=I,N+1 M2(J,K)= M2(J,K) - TEMP*M2(I,K) END DO L3 END DO L2 END DO L1 ! PRINT *,'TRIANGULARIZED MATRIX' CALL PRNMAT(M2,5) ! ! MATRIX IS NOW TRIANGULAR. USE BACKSUBSTITUTION TO SOLVE ! SOL(N)=M2(N,N+1)/M2(N,N) L4: DO I=N-1,1,-1 TEMP=0.0 L5: DO K=N,I+1,-1 TEMP=TEMP+M2(I,K)*SOL(K) END DO L5 SOL(I)=(M2(I,N+1) - TEMP)/M2(I,I) END DO L4 RETURN END SUBROUTINE GAUSS