Esope source files

   1 : C CHMSMQ    SOURCE    CHAT      05/01/12    22:00:15     5004           
   2 :       SUBROUTINE CHMSMQ(Z,Y,N,NXDIM,IER)
   3 : C=======================================================================        
   4 : C        ISSU DE TRIOEF ( TRSIMQ)
   5 : C                                                                               
   6 : C OBJET: TROUVE LA SOLUTION D'UN ENSEMBLE DE N<=NXDIM EQUATIONS LINEAIRES       
   7 : C NON HOMOGENES Y=Z.X: TROUVE X                                                 
   8 : C                                                                               
   9 : C ARGUMENTS: Z(NXDIM,NXDIM): MATRICE DES COEFFICIENTS                           
  10 : C            Y(NYDIM): TABLEAU DES RESIDUS                                      
  11 : C            N :DIMENSION ACTUELLE DU TABLEAU Z                                 
  12 : C            NXDIM: DIMENSION MAXIMUM DU TABLEAU Z                              
  13 : C                                                                               
  14 : C EST APPELLEE PAR CHMSLV
  15 : C IER=8  SI LA MATRICE EST SINGULIERE
  16 : C                                                                               
  17 : C    FORTRAN 77
  18 : C                                                                               
  19 : C=======================================================================        
  20 :       IMPLICIT INTEGER(I-N)
  21 :       IMPLICIT REAL*8(A-H,O-Z)                                                  
  22 :       DIMENSION Z(NXDIM,NXDIM),Y(NXDIM)                                         
  23 :                                                                                 
  24 : C     PROVISION FOR N=1                                                         
  25 :       IF(N.EQ.1) THEN
  26 :       Y(1)=Y(1)/Z(1,1)                                                          
  27 :       RETURN                                                                    
  28 :       ENDIF
  29 :                                                                                 
  30 : C     ELEMENT OF ELIMINATION                                                    
  31 :                                                                                 
  32 :       L1=N-1
  33 :       DO 10 M=1,L1
  34 :       ZMAX=0.D0
  35 :       IMAX=0                                                                    
  36 :                                                                                 
  37 : C     FIND MAX OF COLUMN                                                        
  38 :                                                                                 
  39 :       DO 20 I=M,N                                                               
  40 :         IF(ABS(Z(I,M)).LE.ZMAX) GO TO 20                                        
  41 :         IMAX=I                                                                  
  42 :         ZMAX=ABS(Z(I,M))                                                        
  43 : 20    CONTINUE                                                                  
  44 :                                                                                 
  45 : C     ERROR RETURN                                                              
  46 :                                                                                 
  47 :       IF(IMAX.EQ.0) THEN
  48 :       IER=8
  49 :       RETURN
  50 :       ENDIF
  51 :                                                                                 
  52 : C     ROW INTERCHANGE                                                           
  53 :                                                                                 
  54 :       IF(IMAX.NE.M) THEN
  55 :                                                                                 
  56 :             V=Y(M)                                                                    
  57 :             Y(M)=Y(IMAX)                                                              
  58 :             Y(IMAX)=V                                                                 
  59 :                                                                                       
  60 :             DO 40 J=M,N                                                               
  61 :                V=Z(M,J)                                                               
  62 :                Z(M,J)=Z(IMAX,J)                                                       
  63 :                Z(IMAX,J)=V                                                            
  64 : 40          CONTINUE
  65 :       ENDIF
  66 :                                                                                 
  67 :                                                                                 
  68 : C     DIAGONALIZE                                                               
  69 :                                                                                 
  70 :       M1=M+1                                                                    
  71 :       DO 70 I=M1,N                                                              
  72 :          V=Z(I,M)/Z(M,M)                                                        
  73 :          Y(I)=Y(I)-V*Y(M)                                                       
  74 :          DO 65 J=M,N                                                            
  75 :             Z(I,J)=Z(I,J)-V*Z(M,J)                                              
  76 : 65       CONTINUE                                                               
  77 : 70    CONTINUE                                                                  
  78 :                                                                                 
  79 : 10    CONTINUE                                                                  
  80 :                                                                                 
  81 : C     BACK SUBSTITUTE                                                           
  82 :                                                                                 
  83 :       Y(N)=Y(N)/Z(N,N)                                                          
  84 :       L1=N-1
  85 :       DO 100 K=1,L1
  86 :          I=N-K                                                                  
  87 :          I1=I+1                                                                 
  88 :          DO 90 J=I1,N                                                           
  89 :             Y(I)=Y(I)-Y(J)*Z(I,J)                                               
  90 :  90      CONTINUE                                                               
  91 :          Y(I)=Y(I)/Z(I,I)                                                       
  92 : 100   CONTINUE                                                                  
  93 :                                                                                 
  94 :       RETURN                                                                    
  95 :       END                                                                       
  96 :  
  97 :