Sources Esope | Cast3M

jcbien
C JCBIEN    SOURCE    CHAT      05/01/13    00:49:21     5004
      SUBROUTINE JCBIEN(XREF,ITYP,NDIM,NOEL1,XJAC,DET,IZSH)
*******************************************************************
*** SP 'JCBIEN' : calcul du jacobien pour le passage entre l'element
*** de reference et l'element reel.
***
*** APPELES 1 = 'SHAPE'
*** APPELES 2 = aucun
***
*** E = 'XREF' coordonnees de reference du pt de calcul
***     'ITYP' entier caracterisant le type de l'element
***     'NDIM' dimension de l'espace
***     'NOEL1' nbre de noeuds de l'element
***
*** S = 'XJAC' matrice Jacobienne au pt de calcul
***     'DET' valeur du Jacobien au pt de calcul
***
*** E/S = 'IZSH' segment content parametres utilises suivants :
***     'XYZL' coord reelles des noeuds de l'element (E)
***     'SHP(1,I)' fcts forme au pt calcul de reference (S)
***     'SHP(2,I)' derivee fcts forme / xsi au pt calcul reference (S)
***     'SHP(3,I)' derivee fcts forme / eta au pt calcul reference (S)
***     'SHP(4,I)' derivee fcts forme / dzeta au pt calcul reference (S)
***
***    Auteur  Patrick Meyniel
********************************************************************
 
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8 (A-H,O-Z)
      SEGMENT IZSH
          REAL*8 SHP(6,MNO9),SHY(12,MNO9),XYZL(3,MNO9)
      ENDSEGMENT
      DIMENSION XREF(*),XJAC(3,*)
***   calcul des fonctions de forme et de leur derivees
      CALL SHAPE(XREF(1),XREF(2),XREF(3),ITYP,SHP,IRET)
***   initialisation des coefficients de la matrice jacobienne
      DET=0.D0
      DO 10 I=1,3
         DO 20 J=1,3
            XJAC(I,J)=0.D0
 20      CONTINUE
 10   CONTINUE
***   calcul des coefficients de derivation de la matrice jacobienne
      IF (NDIM.EQ.2) THEN
***      cas 2D
         DO 30 I=1,NOEL1
            XJAC(1,1)=XJAC(1,1)+SHP(2,I)*XYZL(1,I)
            XJAC(1,2)=XJAC(1,2)+SHP(3,I)*XYZL(1,I)
            XJAC(2,1)=XJAC(2,1)+SHP(2,I)*XYZL(2,I)
            XJAC(2,2)=XJAC(2,2)+SHP(3,I)*XYZL(2,I)
 30      CONTINUE
         DET=XJAC(1,1)*XJAC(2,2)-XJAC(1,2)*XJAC(2,1)
      ELSE
***      cas 3D
         DO 40 I=1,NOEL1
            XJAC(1,1)=XJAC(1,1)+SHP(2,I)*XYZL(1,I)
            XJAC(1,2)=XJAC(1,2)+SHP(3,I)*XYZL(1,I)
            XJAC(1,3)=XJAC(1,3)+SHP(4,I)*XYZL(1,I)
            XJAC(2,1)=XJAC(2,1)+SHP(2,I)*XYZL(2,I)
            XJAC(2,2)=XJAC(2,2)+SHP(3,I)*XYZL(2,I)
            XJAC(2,3)=XJAC(2,3)+SHP(4,I)*XYZL(2,I)
            XJAC(3,1)=XJAC(3,1)+SHP(2,I)*XYZL(3,I)
            XJAC(3,2)=XJAC(3,2)+SHP(3,I)*XYZL(3,I)
            XJAC(3,3)=XJAC(3,3)+SHP(4,I)*XYZL(3,I)
 40     CONTINUE
         DET=XJAC(1,1)*(XJAC(2,2)*XJAC(3,3)-XJAC(3,2)*XJAC(2,3))
     $       +XJAC(1,2)*(XJAC(3,1)*XJAC(2,3)-XJAC(2,1)*XJAC(3,3))
     $       +XJAC(1,3)*(XJAC(2,1)*XJAC(3,2)-XJAC(2,2)*XJAC(3,1))
      ENDIF
      RETURN
      END