Sources Esope | Cast3M

bgfcq2
C BGFCQ2    SOURCE    AM        10/09/20    21:15:00     6756
      SUBROUTINE BGFCQ2(B,NGRA,DJAC,IGAU,IFOU,XEL,NN,T,P,EXCEN,DIM3,
     .     IARR,XDPGE,YDPGE)
C==================================================================
C           CALCUL DE LA MATRICE B DES GRADIENTS DE FLEXION
C                  DES COQUES @ 2 NOEUDS
C==================================================================
C ENTREES
C     IGAU=NUMERO DU POINT DE GAUSS
C     IFOU=IFOUR DE CCOPTIO
C     XEL=COORDONNEES LOCALE DE L'ELEMENT
C     NN=NUMERO DU MODE DE FOURIER
C     T(IGAU)=POSITION DU POINT DE GAUSS
C     P(IGAU)=POIDS DU POINT DE GAUSS
C     EXCEN  = EXCENTREMENT
C     DIM3   = EPAISSEUR DANS L'AUTRE DIMENSION
C     XDPGE,YDPGE : COORDONNEE DU POINT AUTOUR DUQUEL
C          FAIT LE MOUVEMENT EN DEFO PLAN GENE
C SORTIE
C     B(2,*)=MATRICE B AU POINT DE GAUSS
C     DJAC=JACOBIEN AU POINT DE GAUSS=POIGAU*LONG/2 (*R(IGAU), SI
C          IFOU EST SUPERIEUR OU EGAL A ZERO)
C     IARR=0 SI OK 1 SI LONGUEUR ELEMENT NULLE
C                  2 SI R / D INFERIEUR A 10-3
C==================================================================
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8(A-H,O-Z)
C    Include contenant quelques constantes dont XPI, XZERO :
-INC CCREEL
      PARAMETER(UNDE=.5D0,UN=1.D0,DEUX=2.D0,TRS=3.D0)
      PARAMETER(QUTR=4.D0,SIX=6.D0,DOUZ=12.D0)
      DIMENSION B(NGRA,*),T(*),XEL(3,*),P(*)
C
C     ---------------------------------INITIALISATION
      IARR=0
      DJAC=XZERO
      IF(IFOU.GT.0)  THEN
         LRE=8
         CALL ZERO(B,2,8)
      ELSE IF(IFOU.LE.0) THEN
         LRE=6
         CALL ZERO(B,2,6)
      ENDIF
C
      D=SQRT((XEL(1,2)-XEL(1,1))**2+(XEL(2,2)-XEL(2,1))**2)
      IF(D.EQ.0) THEN
         IARR=1
         GOTO 4
      ENDIF
      DD=UN/D
      RO=(XEL(1,1)+XEL(1,2))*UNDE
      R1=(XEL(2,1)+XEL(2,2))*UNDE
      SP=(XEL(1,2)-XEL(1,1))/D
      CP=(XEL(2,2)-XEL(2,1))/D
      SP2=SP*SP
      CP2=CP*CP
      SPCP=SP*CP
C           X FONCTION FORME NOEUD 2     RRRR RAYON
C           Y FONCTION FORME NOEUD 1     D LONGUEUR DD INVERSE LONGUEUR
      X=UNDE+UNDE*T(IGAU)
      Y=UNDE-UNDE*T(IGAU)
      RRRR=RO+UNDE*D*SP*T(IGAU)
C     ---------------------------------
C
C     TEST D'ERREUR
C
      IF(IFOU.GE.0) THEN
         IF(ABS(RRRR/D).LE.1.D-03) THEN
            IARR=2
            GOTO 4
         ENDIF
      ENDIF
C
C     ---------------------------------CALCULS
C
      IF(IFOU.LT.0) RRRR =UN
      U=X/RRRR
      V=Y/RRRR
C
C          AXISYMETRIQUE DEFORMATIONS PLANES CONTRAINTES PLANES
C
      IF(IFOU.LE.0) THEN
C
C         RT,S
C
         AUX = SIX*T(IGAU)*DD*DD
         B(1,1)= CP*AUX
         B(1,2)=-SP*AUX
         B(1,3)=(QUTR-SIX*X)*DD
         B(1,4)=-B(1,1)
         B(1,5)=-B(1,2)
         B(1,6)=(DEUX-SIX*X)*DD
C
C          FOURIER
C
      ELSE IF(IFOU.GT.0) THEN
         AN=DBLE(NN)
C
C           RT,S
C
         AUX= SIX*T(IGAU)*DD*DD
         B(1,1)= CP*AUX
         B(1,2)=-SP*AUX
         B(1,4)=(QUTR-SIX*X)*DD
         B(1,5)=-B(1,1)
         B(1,6)=-B(1,2)
         B(1,8)=(DEUX-SIX*X)*DD
C
C          RT,T
C
         AUX1=DOUZ*U*Y*AN*DD
         AUX =-AUX1
         B(3,1)=-CP*AUX
         B(3,2)= SP*AUX
         B(3,3)= 0.
         B(3,4)= AN*V*(DEUX-SIX*X)
         AUX= AUX1
         B(3,5)=-CP*AUX
         B(3,6)= SP*AUX
         B(3,7)= 0.
         B(3,8)= AN*U*(SIX*X-QUTR)
      ENDIF
      IF(IFOU.EQ.0.OR.(IFOU.EQ.1.AND.NN.EQ.0)) THEN
         DJAC=D*UNDE*P(IGAU)*RRRR*2*XPI
      ELSEIF(IFOU.EQ.1.AND.NN.NE.0) THEN
         DJAC=D*UNDE*P(IGAU)*RRRR*XPI
      ELSE
         DJAC=D*UNDE*P(IGAU)*RRRR*DIM3
      ENDIF
*
 4    CONTINUE
      RETURN
      END