C BGRMAS SOURCE CHAT 05/01/12 21:39:48 5004 C======================================================================= C= B G R M A S = C= ----------- = C= = C= Fonction : = C= ---------- = C= Calcul de la matrice de gradients BGR au point de Gauss iGau de = C= l'element fini de type iTEl. = C= Le jacobien est egalement evalue au point de Gauss pour verifier = C= ulterieurement si l'element fini n'est pas trop distordu. = C= = C= Parametres : (E)=Entree (S)=Sortie = C= ------------ = C= iGau (E) Numero du point de Gauss considere = C= iTEl (E) Type de l'element fini (cf. NOMTP dans bdata.eso) = C= NBNO (E) Nombre de noeuds de l'element fini = C= LRE (E) Nombre de DDL de l'element fini = C= IFOU (E) Mode de calcul utilise (cf. IFOUR dans CCOPTIO) = C= NGRA (E) Nombre de composantes du gradient = C= NN (E) Numero du mode de Fourier (si IFOU=1) = C= XEL (E) Coordonnees des noeuds de l'element fini etudie = C= Excen (E) Excentrement dans le cas d'un element coque = C= SHPTOT (E) Fonctions de forme et leurs derivees = C= SHP (S) Fonctions de forme et leurs derivees actuelles = C= BB (S) Matrice de travail = C= BGR (S) Matrice de gradients (B) calcule au point de Gauss = C= DJac (S) Jacobien au point de Gauss etudie = C= IIPDPG (E) Numero du noeud support en modes GENERALISEs = C= = C= Remarque : = C= ---------- = C= SHPTOT(2 a 4,*) contient les DERIVEES des fonctions de forme par = C= rapport aux coordonnees de REFERENCE Qsi,Eta,Dzeta. = C= En sortie du sousprogramme, SHPTOT(2 a 4,*) contient les DERIVEES = C= des fonctions de FORME par rapport aux coordonnees REELLES x,y,z. = C======================================================================= SUBROUTINE BGRMAS (iGau,iTEl,NBNO,LRE,IFOU,NGRA,NN, . XEL,Excen,SHPTOT,SHP,BB,BGR,DJAC,IIPDPG) IMPLICIT INTEGER(I-N) IMPLICIT REAL*8 (A-H,O-Z) -INC PPARAM -INC CCOPTIO -INC CCREEL -INC SMCOORD DIMENSION XEL(3,*),BGR(NGRA,*),SHP(6,*),SHPTOT(6,NBNO,*),BB(2,*) DIMENSION GEOM(20),XX(3),YY(3),BBF(4,9) DATA XX / .5, .0, .5 / DATA YY / .0, .5, .5 / C= Quelques constantes (2.Pi et 4.Pi) PARAMETER (X2Pi=6.283185307179586476925286766559D0) PARAMETER (X4Pi=12.566370614359172953850573533118D0) CALL ZERO(BGR,NGRA,LRE) IF (iTEl.EQ.28.OR.iTEl.EQ.45) GOTO 28 IF (iTEl.GE.57.AND.iTEl.LE.68) GOTO 57 C 1 - Elements MASSIFS en MECANIQUE C =================================== IFOR=IFOU+4 GOTO (100,110,110,120,130,140,150,150,150,150, . 150,150,150,150,150,160,160,160,170 ),IFOR RETURN C ===== C 1.1 - Elements MASSIFS 2D PLAN GENE C ===== 100 XXX=XZero YYY=XZero DO i=1,NBNO SHP(1,i)=SHPTOT(1,i,iGau) SHP(2,i)=SHPTOT(2,i,iGau) SHP(3,i)=SHPTOT(3,i,iGau) XXX=XXX+SHP(1,i)*XEL(1,i) YYY=YYY+SHP(1,i)*XEL(2,i) ENDDO CALL JACOBI(XEL,SHP,2,NBNO,DJAC) K=1 DO i=1,NBNO BGR(1,K)=SHP(2,i) BGR(2,K)=SHP(3,i) BGR(4,K+1)=SHP(2,i) BGR(5,K+1)=SHP(3,i) K=K+2 ENDDO IREF=(IIPDPG-1)*(IDIM+1) BGR(9,K)=1. BGR(9,K+1)=XCOOR(IREF+1)-XXX BGR(9,K+2)=YYY-XCOOR(IREF+2) RETURN C ===== C 1.2 - Elements MASSIFS 2D PLAN DEFO ou PLAN CONT C ===== 110 DO i=1,NBNO SHP(1,i)=SHPTOT(1,i,iGau) SHP(2,i)=SHPTOT(2,i,iGau) SHP(3,i)=SHPTOT(3,i,iGau) ENDDO CALL JACOBI(XEL,SHP,2,NBNO,DJAC) K=1 DO i=1,NBNO BGR(1,K)=SHP(2,i) BGR(2,K)=SHP(3,i) BGR(4,K+1)=SHP(2,i) BGR(5,K+1)=SHP(3,i) K=K+2 ENDDO RETURN C ===== C 1.3 - Elements MASSIFS 2D AXISymetrique C ===== 120 DO i=1,NBNO SHP(1,i)=SHPTOT(1,i,iGau) SHP(2,i)=SHPTOT(2,i,iGau) SHP(3,i)=SHPTOT(3,i,iGau) ENDDO CALL JACOBI(XEL,SHP,2,NBNO,DJAC) CALL DISTRR(XEL,SHP,NBNO,RR) C*OF DJAC=X2Pi*DJAC*RR DJAC=DJAC*RR K=1 DO i=1,NBNO BGR(1,K)=SHP(2,i) BGR(2,K)=SHP(3,i) BGR(4,K+1)=SHP(2,i) BGR(5,K+1)=SHP(3,i) BGR(9,K)=SHP(1,i)/RR K=K+2 ENDDO RETURN C ===== C 1.4 - Elements MASSIFS 2D FOURIER C ===== 130 DO i=1,NBNO SHP(1,i)=SHPTOT(1,i,iGau) SHP(2,i)=SHPTOT(2,i,iGau) SHP(3,i)=SHPTOT(3,i,iGau) ENDDO CALL JACOBI(XEL,SHP,2,NBNO,DJAC) CALL DISTRR(XEL,SHP,NBNO,RR) C*OF DJAC=X2Pi*DJAC*RR C*OF IF (NN.EQ.0) DJAC=0.5*DJAC DJAC=DJAC*RR XNSUR=DBLE(NN)/RR K=1 DO i=1,NBNO BGR(1,K)=SHP(2,i) BGR(2,K)=SHP(3,i) BGR(3,K)=-SHP(1,i)*XNSUR BGR(3,K+2)=-SHP(1,i)/RR BGR(4,K+1)=SHP(2,i) BGR(5,K+1)=SHP(3,i) BGR(6,K+1)=-SHP(1,i)*XNSUR BGR(7,K+2)=SHP(2,i) BGR(8,K+2)=SHP(3,i) BGR(9,K)=SHP(1,i)/RR BGR(9,K+2)=SHP(1,i)*XNSUR K=K+3 ENDDO RETURN C ===== C 1.6 - Elements MASSIFS 3D C ===== 140 DO i=1,NBNO SHP(1,i)=SHPTOT(1,i,iGau) SHP(2,i)=SHPTOT(2,i,iGau) SHP(3,i)=SHPTOT(3,i,iGau) SHP(4,i)=SHPTOT(4,i,iGau) ENDDO CALL JACOBI(XEL,SHP,3,NBNO,DJAC) K=1 DO i=1,NBNO BGR(1,K)=SHP(2,i) BGR(2,K)=SHP(3,i) BGR(3,K)=SHP(4,i) BGR(4,K+1)=SHP(2,i) BGR(5,K+1)=SHP(3,i) BGR(6,K+1)=SHP(4,i) BGR(7,K+2)=SHP(2,i) BGR(8,K+2)=SHP(3,i) BGR(9,K+2)=SHP(4,i) K=K+3 ENDDO RETURN C ===== C 1.7 - Elements MASSIFS 1D C ===== C= 1.7.1 - Modes 1D PLAN 150 DO i=1,NBNO SHP(1,i)=SHPTOT(1,i,iGau) SHP(2,i)=SHPTOT(2,i,iGau) ENDDO CALL JACOBI(XEL,SHP,1,NBNO,DJAC) DO i=1,NBNO BGR(1,i)=SHP(2,i) ENDDO IF (IFOU.GE.7) THEN BGR(9,i)=1. IF (IFOU.EQ.11) BGR(5,i+1)=1. ENDIF RETURN C= 1.7.2 - Modes 1D AXIS 160 DO i=1,NBNO SHP(1,i)=SHPTOT(1,i,iGau) SHP(2,i)=SHPTOT(2,i,iGau) ENDDO CALL JACOBI(XEL,SHP,1,NBNO,DJAC) CALL DISTRR(XEL,SHP,NBNO,RR) C*OF DJAC=X2Pi*DJAC*RR DJAC=DJAC*RR DO i=1,NBNO BGR(1,i)=SHP(2,i) BGR(9,i)=SHP(1,i)/RR ENDDO IF (IFOU.EQ.14) BGR(5,i)=1. RETURN C= 1.7.3 - Mode 1D SPHE 170 DO i=1,NBNO SHP(1,i)=SHPTOT(1,i,iGau) SHP(2,i)=SHPTOT(2,i,iGau) ENDDO CALL JACOBI(XEL,SHP,1,NBNO,DJAC) CALL DISTRR(XEL,SHP,NBNO,RR) C*OF DJAC=X4Pi*DJAC*RR*RR DJAC=DJAC*RR*RR DO i=1,NBNO BGR(1,i)=SHP(2,i) ZZ=SHP(1,i)/RR BGR(5,i)=ZZ BGR(9,i)=ZZ ENDDO RETURN C 2 - Elements MASSIFS en THERMIQUE C =================================== 57 IFOR=IFOU+4 GOTO (200,200,200,210,210,220,230,230,230,230, . 230,230,230,230,230,230,230,230,230 ),IFOR RETURN C ===== C 2.1 - Elements MASSIFS 2D C ===== 200 DO i=1,NBNO SHP(1,i)=SHPTOT(1,i,iGau) SHP(2,i)=SHPTOT(2,i,iGau) SHP(3,i)=SHPTOT(3,i,iGau) ENDDO CALL JACOBI(XEL,SHP,2,NBNO,DJAC) DO i=1,NBNO BGR(1,i)=SHP(2,i) BGR(2,i)=SHP(3,i) ENDDO RETURN C ===== C 2.2 - Elements MASSIFS 2D AXISymetrique C 2.3 - Elements MASSIFS 2D FOURIER C ===== 210 DO i=1,NBNO SHP(1,i)=SHPTOT(1,i,iGau) SHP(2,i)=SHPTOT(2,i,iGau) SHP(3,i)=SHPTOT(3,i,iGau) ENDDO CALL JACOBI(XEL,SHP,2,NBNO,DJAC) CALL DISTRR(XEL,SHP,NBNO,RR) C*OF DJAC=X2Pi*DJAC*RR DJAC=DJAC*RR DO i=1,NBNO BGR(1,i)=SHP(2,i) BGR(2,i)=SHP(3,i) ENDDO RETURN C ===== C 2.4 - Elements MASSIFS 3D C ===== 220 DO i=1,NBNO SHP(1,i)=SHPTOT(1,i,iGau) SHP(2,i)=SHPTOT(2,i,iGau) SHP(3,i)=SHPTOT(3,i,iGau) SHP(4,i)=SHPTOT(4,i,iGau) ENDDO CALL JACOBI(XEL,SHP,3,NBNO,DJAC) DO i=1,NBNO BGR(1,i)=SHP(2,i) BGR(2,i)=SHP(3,i) BGR(3,i)=SHP(4,i) ENDDO RETURN C ===== C 2.5 - Elements MASSIFS 1D C ===== 230 DO i=1,NBNO SHP(1,i)=SHPTOT(1,i,iGau) SHP(2,i)=SHPTOT(2,i,iGau) ENDDO CALL JACOBI(XEL,SHP,1,NBNO,DJAC) IF (IFOUR.GE.12.AND.IFOUR.LE.14) THEN CALL DISTRR(XEL,SHP,NBNO,RR) C*OF DJAC=X2Pi*DJAC*RR DJAC=DJAC*RR ELSE IF (IFOUR.EQ.15) THEN CALL DISTRR(XEL,SHP,NBNO,RR) C*OF DJAC=X4Pi*DJAC*RR*RR DJAC=DJAC*RR*RR ENDIF DO i=1,NBNO BGR(1,i)=SHP(2,i) ENDDO RETURN C 3 - Element COQUE DKT en MECANIQUE C ==================================== 28 DO i=1,NBNO SHP(1,i)=SHPTOT(1,i,iGau) SHP(2,i)=SHPTOT(2,i,iGau) SHP(3,i)=SHPTOT(3,i,iGau) ENDDO CALL JACOBI(XEL,SHP,2,NBNO,DJAC) K=1 DO i=1,NBNO BGR(1,K)=SHP(2,i) BGR(2,K)=SHP(3,i) BGR(4,K+1)=SHP(2,i) BGR(5,K+1)=SHP(3,i) K=K+6 ENDDO C= Prise en compte de l'excentrement des elements CALL GEOCST(XEL,GEOM) CALL BBGFDK(XX(iGau),YY(iGau),GEOM,BBF) K=2 KK=0 DO i=1,NBNO DO j=1,3 BGR(1,K+j)=Excen*BBF(1,j+KK) BGR(2,K+j)=Excen*BBF(2,j+KK) BGR(4,K+j)=Excen*BBF(3,j+KK) BGR(5,K+j)=Excen*BBF(4,j+KK) ENDDO K=K+6 KK=KK+3 ENDDO C= Fin du traitement de l'excentrement CALL GERADK(XEL,GEOM) CALL BBGRDK(XX(iGau),YY(iGau),GEOM,BB) K=2 KK=0 DO i=1,3 DO j=1,3 BGR(3,K+j)= BB(1,j+KK) BGR(6,K+j)= BB(2,j+KK) BGR(7,K+j)=-BB(1,j+KK) BGR(8,K+j)=-BB(2,j+KK) ENDDO K=K+6 KK=KK+3 ENDDO RETURN END