Sources Esope | Cast3M

bdbs1
C BDBS1     SOURCE    CHAT      05/01/12    21:37:51     5004
 
C=======================================================================
C=                             B D B S 1                               =
C=                             ---------                               =
C=                                                                     =
C=  Fonction :                                                         =
C=  ----------                                                         =
C=   Calcul de la contribution du point de Gauss iGau a la matrice de  =
C=   RIGIDITE elastique REL d'un element fini. Traitement de cas par-  =
C=   ticuliers suivant l'element fini et le type de materiau utilise.  =
C=   Contribution : Jacobien.transpose(BGENE)*DHOOK*BGENE              =
C=                                                                     =
C=  Parametres :  (E)=Entree  (S)=Sortie                               =
C=  ------------                                                       =
C=   BGENE    (E)   Matrice B (B-Barre dans le cas incompressible)     =
C=   DJac     (E)   "Volume" associe au point de Gauss                 =
C=   DHOOK    (E)   Matrice de Hooke au point de Gauss                 =
C=   LRE      (E)   Nombre de DDL de l'element fini                    =
C=   NSTRS    (E)   Nombre de composantes de contraintes/deformations  =
C=   REL     (E/S)  Matrice de rigidite de l'element fini              =
C=   MFR      (E)   Formulation associee a l'element fini              =
C=   IFOU     (E)   Mode de calcul utilise (cf. IFOUR dans CCOPTIO)    =
C=   MATE     (E)   Type du materiau utilise                           =
C=   iGau     (E)   Numero du point de Gauss considere                 =
C=   Excen    (E)   Excentrement (coque avec couplage membrane-flexion)=
C=   IMAT     (E)   =1 si                                              =
C=                                                                     =
C=  Remarque : Actuellement, NSTRS doit etre inferieur ou egal a 10 !  =
C=======================================================================
 
      SUBROUTINE BDBS1 (BGENE,DJac,DHOOK,LRE,NSTRS,REL,MFR,IFOU,MATE,
     .                  iGau,IMAT,Excen)
 
      IMPLICIT INTEGER(I-N)
      IMPLICIT REAL*8 (A-H,O-Z)
 
-INC CCREEL
 
      DIMENSION DHOOK(NSTRS,*),REL(LRE,*),BGENE(NSTRS,*)
      DIMENSION WORK1(10),WORK2(10),TRD1(3,3),TRD2(3,3),TRG(3,3)
 
C  1 - Determination de la methode de multiplication matricielle la
C      plus adaptee a l'element fini et au type de materiau utilise
C ==================================================================
      IMULT=0
      IF (IMAT.EQ.1) THEN
        IF (MFR.EQ.1.AND.MATE.EQ.1) THEN
          IF (IFOU.EQ.-1.OR.IFOU.EQ.-2) THEN
            IMULT=1
          ELSE IF (IFOU.EQ. 2) THEN
            IMULT=2
          ENDIF
        ELSE IF (MFR.EQ.3) THEN
          IF (MATE.EQ.1.AND.Excen.EQ.XZero) IMULT=3
C*      ELSE IF (MFR.EQ.31.AND.MATE.EQ.1) THEN
C*        IF (IFOU.EQ.-1.OR.IFOU.EQ.-2) THEN
C*          IMULT=1
C*        ELSE IF (IFOU.EQ. 2) THEN
C*          IMULT=2
C*        ENDIF
        ENDIF
      ENDIF
 
C  2 - Calcul de la contribution DJac*transpose(BGENE)*DHOOK*BGENE
C =================================================================
C  2.1 - Multiplication matricielle normale
C =====
      IF (IMULT.EQ.0) THEN
        DO IX=1,LRE
          DO j=1,NSTRS
            WORK1(j)=DJac*BGENE(j,IX)
          ENDDO
          DO i=1,NSTRS
            CC=XZero
            DO j=1,NSTRS
              CC=CC+WORK1(j)*DHOOK(i,j)
            ENDDO
            WORK2(i)=CC
          ENDDO
          DO IY=1,IX
            CC=XZero
            DO i=1,NSTRS
              CC=CC+BGENE(i,IY)*WORK2(i)
            ENDDO
            REL(IX,IY)=REL(IX,IY)+CC
            REL(IY,IX)=REL(IX,IY)
          ENDDO
        ENDDO
C =====
C  2.2 - Cas des elements massifs 2D incompressibles et isotropes
C =====
      ELSE IF (IMULT.EQ.1) THEN
        LRE2=LRE-1
        D1=DHOOK(1,1)*DJac
        D2=DHOOK(1,2)*DJac
        G=DHOOK(4,4)*DJac
        DO j=1,LRE2,2
          j1=j+1
          DO i=1,LRE2,2
            i1=i+1
            Tr11=BGENE(1,i )*BGENE(1,j )
            Tr21=BGENE(2,i1)*BGENE(1,j )
            Tr12=BGENE(1,i )*BGENE(2,j1)
            Tr22=BGENE(2,i1)*BGENE(2,j1)
            REL(i ,j )=REL(i ,j )+D1*Tr11+G*Tr22
            REL(i1,j )=REL(i1,j )+D2*Tr21+G*Tr12
            REL(i ,j1)=REL(i ,j1)+D2*Tr12+G*Tr21
            REL(i1,j1)=REL(i1,j1)+D1*Tr22+G*Tr11
          ENDDO
        ENDDO
C =====
C  2.3 - Cas des elements massifs 3D incompressibles et isotropes
C =====
      ELSE IF (IMULT.EQ.2) THEN
        LRE3=LRE-2
        D1=DHOOK(1,1)*DJac
        D2=DHOOK(1,2)*DJac
        G=DHOOK(4,4)*DJac
        DO j=1,LRE3,3
          j1=j+1
          j2=j+2
          DO i=1,LRE3,3
            i1=i+1
            i2=i+2
            TR11=BGENE(1,i )*BGENE(1,j )
            TR12=BGENE(1,i )*BGENE(2,j1)
            TR13=BGENE(1,i )*BGENE(3,j2)
            TR21=BGENE(2,i1)*BGENE(1,j )
            TR22=BGENE(2,i1)*BGENE(2,j1)
            TR23=BGENE(2,i1)*BGENE(3,j2)
            TR31=BGENE(3,i2)*BGENE(1,j )
            TR32=BGENE(3,i2)*BGENE(2,j1)
            Tr33=BGENE(3,i2)*BGENE(3,j2)
            REL(i ,j )=REL(i ,j )+D1*Tr11+G*(Tr22+Tr33)
            REL(i1,j )=REL(i1,j )+D2*Tr21+G*Tr12
            REL(i2,j )=REL(i2,j )+D2*Tr31+G*Tr13
            REL(i ,j1)=REL(i ,j1)+D2*Tr12+G*Tr21
            REL(i1,j1)=REL(i1,j1)+D1*Tr22+G*(Tr11+Tr33)
            REL(i2,j1)=REL(i2,j1)+D2*Tr32+G*Tr23
            REL(i ,j2)=REL(i ,j2)+D2*Tr13+G*Tr31
            REL(i1,j2)=REL(i1,j2)+D2*Tr23+G*Tr32
            REL(i2,j2)=REL(i2,j2)+D1*Tr33+G*(Tr11+Tr22)
          ENDDO
        ENDDO
C =====
C  2.4 - Cas des coques minces 3D isotropes
C =====
      ELSE IF (IMULT.EQ.3) THEN
        LRE6=LRE-5
        D1=DHOOK(1,1)*DJac
        D2=DHOOK(1,2)*DJac
        G=DHOOK(3,3)*DJac
        DH1=DHOOK(4,4)*DJac
        DH2=DHOOK(5,4)*DJac
        GH=DHOOK(6,6)*DJac
        DO j=1,LRE6,6
          j1=j+1
          DO i=1,LRE6,6
            i1=i+1
            Tr11=BGENE(1,i )*BGENE(1,j )
            Tr21=BGENE(2,i1)*BGENE(1,j )
            Tr12=BGENE(1,i )*BGENE(2,j1)
            Tr22=BGENE(2,i1)*BGENE(2,j1)
            REL(i ,j )=REL(i ,j )+D1*Tr11+G*Tr22
            REL(i1,j )=REL(i1,j )+D2*Tr21+G*Tr12
            REL(i ,j1)=REL(i ,j1)+D2*Tr12+G*Tr21
            REL(i1,j1)=REL(i1,j1)+D1*Tr22+G*Tr11
            DO l=1,3
              jl=j1+l
              DO k=1,3
                ik=i1+k
                TrD1(k,l)=( BGENE(4,ik)*BGENE(4,jl)
     .                     +BGENE(5,ik)*BGENE(5,jl))*DH1
                TrD2(k,l)=( BGENE(4,ik)*BGENE(5,jl)
     .                     +BGENE(5,ik)*BGENE(4,jl))*DH2
                TrG(k,l)=BGENE(6,ik)*BGENE(6,jl)*GH
              ENDDO
            ENDDO
            DO l=1,3
              jl=j1+l
              DO k=1,3
                ik=i1+k
                REL(ik,jl)=REL(ik,jl)+TrD1(k,l)+TrD2(k,l)+TrG(k,l)
              ENDDO
            ENDDO
          ENDDO
        ENDDO
      ENDIF
 
      RETURN
      END