C GEOLI4    SOURCE    GOUNAND   26/01/09    21:15:24     12441          
      SUBROUTINE GEOLI4(IESREL,IESREF,NDPOGO,NDELEM,JMAJAC,
     $     JJTJ,JJTJM1,
     $     JMIJAC,JDTJAC,
     $     IMPR,IRET)
      IMPLICIT REAL*8 (A-H,O-Z)
      IMPLICIT INTEGER (I-N)
C***********************************************************************
C NOM         : GEOLI4
C PROJET      : Noyau linéaire NLIN
C DESCRIPTION : Calcul du jacobien dans le cas où la matrice
C               jacobienne A n'est pas carrée. A(n,m)
C               On calcule sqrt (det (transpose(A) * A))
C               Ceci est effectué pour chaque point de Gauss d'un
C               élément.
C               On calcule aussi le pseudo-inverse de J
C
C LANGAGE     : Fortran 77 (sauf E/S)
C AUTEUR      : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/LTMF)
C               mél : gounand@semt2.smts.cea.fr
C***********************************************************************
C APPELES          : -
C APPELE PAR       : GEOLIN
C***********************************************************************
C ENTREES            : * IESREL (type entier) : dimension de l'espace
C                        réel (i.e. géométrique).
C                      * IESREF (type entier) : dimension de l'espace de
C                        référence.
C                      * NDNOEU (type entier) : nombre de ddl (par
C                        élément) de la transformation géométrique.
C                      * NDPOGO (type entier) : nombre de points
C                        d'intégration.
C                      * NDELEM (type entier) : nombre d'éléments du
C                        maillage élémentaire courant.
C                      * JMAJAC (type MCHEVA) : valeurs de la matrice
C                        jacobienne aux points de Gauss sur le maillage
C                        élémentaire courant.
C ENTREES/SORTIES    : * JDTJAC (type MCHEVA) : valeurs de
C                        sqrt(det(trans(J).J)) aux points de Gauss sur
C                        le maillage élémentaire courant.
C SORTIES            : -
C CODE RETOUR (IRET) : = 0 si tout s'est bien passé
C***********************************************************************
C VERSION    : v1, 26/10/06, version initiale
C HISTORIQUE : v1, 26/10/06, création
C HISTORIQUE :
C HISTORIQUE :
C***********************************************************************
C Prière de PRENDRE LE TEMPS de compléter les commentaires
C en cas de modification de ce sous-programme afin de faciliter
C la maintenance !
C***********************************************************************

-INC PPARAM
-INC CCOPTIO
      INTEGER IESREF,IESREL,NDPOGO,NDELEM
      REAL*8 JMAJAC(IESREL,IESREF,NDPOGO,NDELEM)
      REAL*8 JMIJAC(IESREF,IESREL,NDPOGO,NDELEM)
      REAL*8 JJTJ(IESREF,IESREF,NDPOGO,NDELEM)
      REAL*8 JJTJM1(IESREF,IESREF,NDPOGO,NDELEM)
      REAL*8 JDTJAC(NDPOGO,NDELEM)
*
      INTEGER IMPR,IRET
*
      INTEGER IELEM,IPOGO,IKF,IIF,IJL,IKL,IJF
*
* Executable statements
*
      IF (IMPR.GT.1) WRITE(IOIMP,*) 'Entrée dans geoli4'
*
* Calcul de JtJ
*
      DO IELEM=1,NDELEM
         DO IPOGO=1,NDPOGO
            DO IKF=1,IESREF
               DO IIF=1,IESREF
                  DO IJL=1,IESREL
                     JJTJ(IIF,IKF,IPOGO,IELEM)=
     $                    JJTJ(IIF,IKF,IPOGO,IELEM)+
     $                    (JMAJAC(IJL,IIF,IPOGO,IELEM)
     $                    *JMAJAC(IJL,IKF,IPOGO,IELEM))
                  ENDDO
               ENDDO
            ENDDO
         ENDDO
      ENDDO
*
* Inverse et déterminant
*
      CALL GEOLI2(IESREF,NDPOGO,NDELEM,JJTJ,
     $     JJTJM1,JDTJAC,
     $     IMPR,IRET)
      IF (IRET.GT.0) THEN
         GOTO 9999
      ELSEIF (IRET.EQ.-666.OR.IRET.EQ.-667) THEN
*         WRITE(IOIMP,*)
*     $        'Determinant change de signe avec matrice positive ??'
*     GOTO 9999
      ELSE
*
* Pseudo-inverse
*
         DO IELEM=1,NDELEM
            DO IPOGO=1,NDPOGO
               DO IKL=1,IESREL
                  DO IIF=1,IESREF
                     DO IJF=1,IESREF
                        JMIJAC(IIF,IKL,IPOGO,IELEM)=
     $                       JMIJAC(IIF,IKL,IPOGO,IELEM)+
     $                       (JJTJM1(IIF,IJF,IPOGO,IELEM)
     $                       *JMAJAC(IKL,IJF,IPOGO,IELEM))
                     ENDDO
                  ENDDO
               ENDDO
            ENDDO
         ENDDO
      ENDIF
*
* Racine du déterminant
*
      DO IELEM=1,NDELEM
         DO IPOGO=1,NDPOGO
            JDTJAC(IPOGO,IELEM)=SQRT(ABS(JDTJAC(IPOGO,IELEM)))
         ENDDO
      ENDDO
*
* Normal termination
*
      RETURN
*
* Format handling
*
*
* Error handling
*
* 9998 CONTINUE
*      WRITE(IOIMP,*) 'Déterminant de la matrice tAA nul ou négatif'
*      WRITE(IOIMP,*) 'IELEM=',IELEM,' IPOGO=',IPOGO
*      GOTO 9999
 9999 CONTINUE
      IRET=1
      WRITE(IOIMP,*) 'An error was detected in subroutine geoli4'
      RETURN
*
* End of subroutine GEOLI4
*
      END