C GAJOME    SOURCE    CHAT      05/01/13    00:15:19     5004
      SUBROUTINE GAJOME(A,M,N,B,IRET,IOIMP)
*
*     Cette subroutine associe a la matrice rectangulaire A de
*     rang M et de dimension MxN ou M<N la matrice carree B de
*     dimension NxN telle que A.B=[I:0] ou I est la matrice
*     identite de dimension MxM et 0 la matrice nulle de
*     dimension Mx(N-M).
*
*     P.PEGON (ISPRA) AOUT 1996
*
*************************************
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8 (A-H,O-Z)
      DIMENSION A(M,N),B(N,N)
*
      IRET=0
      IF (M.GT.N)THEN
        WRITE(IOIMP,*)'GAJOME: reduction impossible'
        IRET=1
        RETURN
      ENDIF
*
*     Initialisation de B a la matrice identite
*
      DO J=1,N
        DO I=1,N
          B(I,J)=0.D0
        ENDDO
        B(J,J)=1.D0
      ENDDO
*
*     Boucle sur les lignes du systeme original
*
      DO I=1,M
*
*     recherche du pivot dnas la ligne I
*
        BIG=1.D-6
        ICOL=0
        DO J=I,N
          IF (ABS(A(I,J)).GE.BIG)THEN
            BIG=ABS(A(I,J))
            ICOL=J
          ENDIF
        ENDDO
        IF(ICOL.EQ.0)THEN
          WRITE(IOIMP,*)'GAJOME: Il y a un mecanisme global'
          IRET=1
          RETURN
        ENDIF
*
*     on place la colone pivot en position I et on fait de meme dans B
*
        IF (ICOL.NE.I)THEN
          DO II=1,M
            DUM=A(II,ICOL)
            A(II,ICOL)=A(II,I)
            A(II,I)=DUM
          ENDDO
          DO II=1,N
            DUM=B(II,ICOL)
            B(II,ICOL)=B(II,I)
            B(II,I)=DUM
          ENDDO
        ENDIF
*
*     on pivote, on normalise la colonne et on reduit la ligne
*
        PIVINV=1.D0/A(I,I)
        DO II=1,M
          A(II,I)=A(II,I)*PIVINV
        ENDDO
        A(I,I)=1.D0
        DO II=1,N
          B(II,I)=B(II,I)*PIVINV
        ENDDO
        DO J=1,N
          IF(J.NE.I)THEN
            DUM=A(I,J)
            IF(DUM.NE.0.D0)THEN
              A(I,J)=0.D0
              DO II=1,M
                A(II,J)=A(II,J)-A(II,I)*DUM
              ENDDO
              DO II=1,N
                B(II,J)=B(II,J)-B(II,I)*DUM
              ENDDO
            ENDIF
          ENDIF
        ENDDO
*
      ENDDO
*
      RETURN
      END


