C GEOLIF    SOURCE    GOUNAND   06/08/04    21:15:57     5520           
      SUBROUTINE GEOLIF(DFFPG,JCOOR,SSFACT,NBELEF,
     $     JMAJAC,JMIJAC,JDTJAC,LERJ,
     $     IMPR,IRET)
      IMPLICIT REAL*8 (A-H,O-Z)
      IMPLICIT INTEGER (I-N)
C***********************************************************************
C NOM         : GEOLIF
C PROJET      : Noyau linéaire NLIN
C DESCRIPTION : Calcul de la matrice jacobienne d'une transformation
C               géométrique aux points de Gauss d'un élément de 
C               référence pour chaque élément réel.
C               On a :
C    Fonction f :      R^m      ->    R^n
C                 a=(a1...am)  |->  f(a)=(f1(a)...fn(a)   
C      =>         matjac(i,j)(a)=Dj fi (a) = dfi / dxj (a)
C
C  Par exemple, pour une surface en 3D, la matrice jacobienne
C  s'exprime comme :
C        
C                [  <xn>  ]    
C         [J] =  [  <yn>  ]   .  [ { dNg/d(ksi) } { dNg/d(eta) } ]
C                [  <zn>  ]
C
C        (3x2)   (3 x Nnoeuds)      (Nnoeuds x 2)
C
C    Ici, le nb de ddl est égal aux nbs de noeuds car l'interpolation
C    pour la géométrie est de type nodale.
C
C               Si la matrice jacobienne J est carrée, on calcule 
C               également l'inverse de la matrice jacobienne j ainsi que
C               son déterminant : det J (appelé Jacobien).
C               Si la matrice jacobienne n'est pas carrée, le Jacobien
C               est calculé comme :
C                                 sqrt(det(transpose(J)*J))
C               Le but de ce sous-programme est d'effectuer toutes les 
C               opérations de calcul dépendant uniquement de la
C               géométrie (ie de la variable : coordonnées des points).
C
C LANGAGE     : ESOPE
C AUTEUR      : Stéphane GOUNAND (CEA/DRN/DMT/SEMT/LTMF)
C               mél : gounand@semt2.smts.cea.fr
C***********************************************************************
C APPELES          : RSET8  (copie de tableaux de CHARACTER*8)
C                    PRCHVA (impression d'un segment de type MCHEVA)
C                    GEOLI1 (calcul de la matrice jacobienne (fortran 77
C                            ))
C                    GEOLI2 (calcul de l'inverse et du déterminant de la
C                            matrice jacobienne (fortran 77))
C                    GEOLI3 (calcul de sqrt(det(transpose(J)*J)), J
C                            étant la matrice jacobienne (fortran 77))
C APPELE PAR       : NLIN
C***********************************************************************
C ENTREES            : * DFFPG (type MCHEVA) : valeurs des dérivées
C                        premières des fonctions d'interpolation de la
C                        transformation géométrique aux points de gauss
C                        sur l'élément de référence. 
C                        Structure (cf.include SMCHAEL) : 
C                        (1, nb. ddl, 1, dim.esp.réf, nb. poi. gauss, 1) 
C                      * JCOOR  (type MCHEVA) : valeurs des ddl de la
C                        transformation géométrique sur le maillage
C                        élémentaire courant.
C                        Structure (cf.include SMCHAEL) : 
C                        (1, nb. ddl, 1, dim. esp. réel, 
C                         1, nb. éléments) 
C                      * NBELEV (type entier) : nombre d'éléments du
C                        maillage élémentaire courant.
C ENTREES/SORTIES    : -
C SORTIES            : * JMAJAC (type MCHEVA) : valeurs de la matrice
C                        jacobienne aux points de Gauss sur le maillage
C                        élémentaire courant. 
C                        Structure (cf.include SMCHAEL) : 
C                        (1, 1, dim. esp. réel, dim. esp. référence,
C                         nb. poi. gauss, nb. éléments) 
C                      * JMIJAC (type MCHEVA) : valeurs de l'inverse de
C                        la matrice jacobienne aux points de Gauss sur
C                        le maillage élémentaire courant. 
C                        Structure (cf.include SMCHAEL) : 
C                        (1, 1, dim. esp. référence, dim. esp. réel,
C                         nb. poi. gauss, nb. éléments) 
C                        JMIJAC n'existe que si dim.esp.réf=dim.esp.réel 
C                      * JDTJAC (type MCHEVA) : valeurs du déterminant
C                        de la matrice jacobienne aux points de Gauss
C                        sur le maillage élémentaire courant. 
C                        Structure (cf.include SMCHAEL) : 
C                        (1, 1, 1, 1, nb. poi. gauss, nb. éléments) 
C                        Si la matrice jacobienne n'est pas carrée, on
C                        a calculé sqrt(det(trans(J).J))
C CODE RETOUR (IRET) : = 0 si tout s'est bien passé
C***********************************************************************
C VERSION    : v2, 03/10/03, refonte complète (modif SMPOUET)
C VERSION    : v1, 09/08/99, version initiale
C HISTORIQUE : v1, 09/08/99, 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 CCOPTIO
CBEGININCLUDE SFACTIV
      SEGMENT FACTIV
      POINTEUR IFACTI(NBSOUV).SFACTI
      ENDSEGMENT
      SEGMENT SFACTI
      POINTEUR ISFACT(NBSOFV).SSFACT
      ENDSEGMENT
      SEGMENT SSFACT
      LOGICAL LFACTI(NBELFV,NBELEV)
      ENDSEGMENT
CENDINCLUDE SFACTIV
CBEGININCLUDE SMCHAEL
      SEGMENT MCHAEL
      POINTEUR IMACHE(N1).MELEME
      POINTEUR ICHEVA(N1).MCHEVA
      ENDSEGMENT
      SEGMENT MCHEVA
      REAL*8 VELCHE(NBLIG,NBCOL,N2LIG,N2COL,NBPOI,NBELM)
      ENDSEGMENT
      SEGMENT LCHEVA
      POINTEUR LISCHE(NBCHE).MCHEVA
      ENDSEGMENT
CENDINCLUDE SMCHAEL
      INTEGER NBLIG,NBCOL,N2LIG,N2COL,NBPOI,NBELM
      POINTEUR JCOOR.MCHEVA
      POINTEUR JMAJAC.MCHEVA
      POINTEUR JMIJAC.MCHEVA
      POINTEUR JDTJAC.MCHEVA
* Valeur des fns d'interpolation de l'espace géométrique aux pts de Gauss
      POINTEUR DFFPG.MCHEVA
*
      INTEGER NBELEV,NBELEF,NBELFV
      INTEGER IMPR,IRET
*
      LOGICAL LCARRE,LERJ
      INTEGER IESREF,IESREL,NBPOGO
*     
* Executable statements
*
      IF (IMPR.GT.1) WRITE(IOIMP,*) 'Entrée dans geolif'
      SEGACT SSFACT
      NBELFV=SSFACT.LFACTI(/1)
      NBELEV=SSFACT.LFACTI(/2)
      SEGACT JCOOR
      NBLIG=JCOOR.VELCHE(/1)
      NBCOL=JCOOR.VELCHE(/2)
      N2LIG=JCOOR.VELCHE(/3)
      N2COL=JCOOR.VELCHE(/4)
      NBPOI=JCOOR.VELCHE(/5)
      NBELM=JCOOR.VELCHE(/6)
* On pourrait envisager un cas où tous les éléments seraient 
* identiques (ex. carré) NBELM=1
      IF (NBLIG.NE.1.OR.N2LIG.NE.1
     $     .OR.NBPOI.NE.1.OR.NBELM.NE.NBELEV) THEN
         WRITE(IOIMP,*) 'Erreur dims JCOOR'
         GOTO 9999
      ENDIF
      NDDL=NBCOL
      IESREL=N2COL
      SEGACT DFFPG
      NBLIG=DFFPG.VELCHE(/1)
      NBCOL=DFFPG.VELCHE(/2)
      N2LIG=DFFPG.VELCHE(/3)
      N2COL=DFFPG.VELCHE(/4)
      NBPOI=DFFPG.VELCHE(/5)
      NBELM=DFFPG.VELCHE(/6)
      IF (NBLIG.NE.1.OR.NBCOL.NE.NDDL.OR.N2LIG.NE.1
     $     .OR.(NBELM.NE.NBELFV.AND.NBELM.NE.1)) THEN
         WRITE(IOIMP,*) 'Erreur dims DFFPG'
         GOTO 9999
      ENDIF
      IESREF=N2COL
      NBPOGO=NBPOI
      NLFVDF=NBELM
      IF (IESREF.LT.IESREL) THEN
         LCARRE=.FALSE.
      ELSEIF (IESREF.EQ.IESREL) THEN
         LCARRE=.TRUE.
      ELSE
         WRITE(IOIMP,*) 'Dim. esp. ref. > dim. esp. reel ???'
         WRITE(IOIMP,*) 'IESREF=',IESREF,' IESREL=',IESREL
         GOTO 9999
      ENDIF
*
* Initialisations...
*
      NBLIG=1
      NBCOL=1
      N2LIG=IESREL
      N2COL=IESREF
      NBPOI=NBPOGO
      NBELM=NBELEF
      SEGINI JMAJAC
*      SEGPRT,DFFPG
*
* On effectue le calcul de la matrice jacobienne aux points de Gauss
*
      CALL GEOLF1(IESREL,IESREF,NDDL,NBPOGO,NBELEV,NBELFV,NBELEF,
     $     NLFVDF,
     $     DFFPG.VELCHE,JCOOR.VELCHE,SSFACT.LFACTI,
     $     JMAJAC.VELCHE,
     $     IMPR,IRET)
      IF (IRET.NE.0) GOTO 9999
*      SEGPRT,JMAJAC
*
* Si la matrice est carrée, on tente de calculer
* son inverse et son déterminant...
*
      IF (LCARRE) THEN
         NBLIG=1
         NBCOL=1
         N2LIG=IESREF
         N2COL=IESREL
         NBPOI=NBPOGO
         NBELM=NBELEF
         SEGINI JMIJAC
         NBLIG=1
         NBCOL=1
         N2LIG=1
         N2COL=1
         NBPOI=NBPOGO
         NBELM=NBELEF
         SEGINI JDTJAC
         CALL GEOLI2(IESREF,NBPOGO,NBELEF,JMAJAC.VELCHE,
     $        JMIJAC.VELCHE,JDTJAC.VELCHE,LERJ,
     $        IMPR,IRET)
*         SEGPRT,JDTJAC
         IF (IRET.NE.0) THEN
            IF (LERJ) GOTO 9666
            GOTO 9999
         ENDIF
         SEGDES JDTJAC
         SEGDES JMIJAC
      ELSE
*
* Sinon on calcule : SQRT ( det (transpose(A) * A))
*
* Pas de matrice inverse
         JMIJAC=0
*
         NBLIG=1
         NBCOL=1
         N2LIG=1
         N2COL=1
         NBPOI=NBPOGO
         NBELM=NBELEF
         SEGINI JDTJAC
         CALL GEOLI3(IESREL,IESREF,NBPOGO,NBELEF,JMAJAC.VELCHE,
     $        JDTJAC.VELCHE,
     $        IMPR,IRET)
         IF (IRET.NE.0) GOTO 9999
         SEGDES JDTJAC
      ENDIF      
      SEGDES JMAJAC
      SEGDES SSFACT
      SEGDES JCOOR
      SEGDES DFFPG
*
* Normal termination
*
      IRET=0
      RETURN
*
* Format handling
*
*
* Error handling
*
 9666 CONTINUE
      IRET=666
      RETURN
 9999 CONTINUE
      IRET=1
      WRITE(IOIMP,*) 'An error was detected in subroutine geolif'
      RETURN
*
* End of subroutine GEOLIF
*
      END