C BBAR SOURCE FANDEUR 10/08/31 21:15:13 6735 SUBROUTINE BBAR(IGAU,NBPGAU,POIGAU,QSIGAU,ETAGAU,DZEGAU, 1 MELE,NBNO,LRE,IFOU,NST,XE,DJAC,A,BB,BGENE) C======================================================================= C= B B A R = C= ----------- = C= Fonction : = C= ---------- = C= CALCUL DE LA MATRICE B-BARRE - ELEMENTS INCOMPRESSIBLES 'IC--' = C= - EN 2D : ELEMENTS ICT3, ICT6, ICQ4, ICQ8 = C= - EN 3D : ELEMENTS ICC8, IC20 = C= ELEMENTS ICT4, IC10, ICP6, IC15 -> A voir = C= = C= Calcul de la matrice B reliant les deformations en un point d'un = C= element fini aux ddls de deplacement aux noeuds de cet element = 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= NBPGAU (E) Nombre de points de Gauss de l'element fini = C= POIGAU,QSIGAU (E) | Poids et coordonnees de reference des = C= ETAGAU,DZEGAU (E) | differents points de Gauss de l'element = C= MELE (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= NST (E) Nombre de composantes de deformations = C= XE (E) Coordonnees des noeuds de l'element fini etudie = C= DJAC (S) Jacobien au point de Gauss etudie = C= BGENE (E/S) Matrice de gradients (B) calcule au point de Gauss = C= A,BB (S) Matrices utiles pour la methode BBar = C======================================================================= C LRE = NOMBRE DE COLONNES DE LA MATRICE B C A = COEFFICIENTS DE MODIFICATION C BB(3,LRE) = MATRICE B-BARRE DILATATION C BGENE(6,LRE) = MATRICE B C BGENE(6,LRE) = MATRICE B-BARRE C======================================================================= IMPLICIT INTEGER(I-N) IMPLICIT REAL*8(A-H,O-Z) PARAMETER (XZero=0., XUn=1., X1s2=0.5, & X1s3=0.333333333333333333333333333333333333333333) DIMENSION XE(3,*),BGENE(NST,*), & POIGAU(*),QSIGAU(*),ETAGAU(*),DZEGAU(*), & A(4,*),BB(3,*) *W WRITE(*,*) 'Entree dans le sousprogramme BBAR' *W WRITE(*,*) 'Element',MELE *W WRITE(*,*) 'Point Gauss: xsi=',QSIGAU(IGAU),'eta=',ETAGAU(IGAU), *W & 'dze=',DZEGAU(IGAU) IF ( MELE.LT.69 .OR. MELE.GT.78 ) GOTO 666 C----------------------------------------------------------------------- C---- Elements incompressibles (MFR = 31) ------------------------------ C----------------------------------------------------------------------- C NOM : ICT3, ICQ4, ICT6, ICQ8, ICC8, ICT4, ICP6, IC20, IC10, IC15 C MELE : 69 , 70 , 71 , 72 , 73 , 74 , 75 , 76 , 77 , 78 C----------------------------------------------------------------------- GOTO ( 69, 70, 71, 72, 73, 74, 75, 76, 77, 78 ) , (MELE-68) GOTO 999 C----------------------------------------------------------------------- C----- Element massif bidimensionnel ICT3 ------------------------------ C----------------------------------------------------------------------- 69 CONTINUE *W WRITE(*,*) 'Element ICT3',IFOU IFR = IFOU+4 GOTO (666,691,691,693,666),IFR GOTO 999 *----- Contraintes planes ou deformations planes (IFOU=-2,-1) 691 CONTINUE K = 1 DO i = 1,NBNO BGENE(2,K) = X1s2 * (BB(1,K)-BGENE(1,K)) BGENE(1,K) = X1s2 * (BB(1,K)+BGENE(1,K)) K = K+1 BGENE(1,K) = X1s2 * (BB(2,K)-BGENE(2,K)) BGENE(2,K) = X1s2 * (BB(2,K)+BGENE(2,K)) K = K+1 ENDDO GOTO 666 *----- Axisymetrie (IFOU=0) 693 CONTINUE K = 1 DO i = 1,NBNO r_z = X1s3 * (BB(1,K)-BGENE(1,K)+BB(3,K)-BGENE(3,K)) BGENE(1,K) = BGENE(1,K) + r_z BGENE(2,K) = r_z BGENE(3,K) = BGENE(3,K) + r_z K = K + 1 r_z = X1s3 * (BB(2,K)-BGENE(2,K)) BGENE(1,K) = r_z BGENE(2,K) = BGENE(2,K) + r_z BGENE(3,K) = r_z K = K + 1 ENDDO GOTO 666 C----------------------------------------------------------------------- C----- Element massif bidimensionnel ICQ4 ------------------------------ C----------------------------------------------------------------------- 70 CONTINUE *W WRITE(*,*) 'Element ICQ4',IFOU IFR = IFOU+4 GOTO (666,701,701,703,666),IFR GOTO 999 *----- Contraintes planes ou deformations planes (IFOU=-2,-1) 701 CONTINUE K = 1 DO i = 1,NBNO BGENE(2,K) = X1s2 * (BB(1,K)-BGENE(1,K)) BGENE(1,K) = X1s2 * (BB(1,K)+BGENE(1,K)) K = K + 1 BGENE(1,K) = X1s2 * (BB(2,K)-BGENE(2,K)) BGENE(2,K) = X1s2 * (BB(2,K)+BGENE(2,K)) K = K + 1 ENDDO GOTO 666 *----- Axisymetrie (IFOU=0) 703 CONTINUE K = 1 DO i = 1,NBNO r_z = X1s3 * (BB(1,K)-BGENE(1,K)+BB(3,K)-BGENE(3,K)) BGENE(1,K) = BGENE(1,K) + r_z BGENE(2,K) = r_z BGENE(3,K) = BGENE(3,K) + r_z K = K + 1 r_z = X1s3 * (BB(2,K)-BGENE(2,K)) BGENE(1,K) = r_z BGENE(2,K) = BGENE(2,K) + r_z BGENE(3,K) = r_z K = K + 1 ENDDO GOTO 666 C----------------------------------------------------------------------- C----- Element massif bidimensionnel ICT6 ------------------------------ C----------------------------------------------------------------------- 71 CONTINUE *W WRITE(*,*) 'Element ICT6',IFOU IFR = IFOU+4 GOTO (666,711,711,713,666),IFR GOTO 999 *----- Contraintes planes ou deformations planes (IFOU=-2,-1) 711 CONTINUE K = 1 DO i = 1,NBNO BGENE(2,K) = X1s2 * (BB(1,K)-BGENE(1,K)) BGENE(1,K) = X1s2 * (BB(1,K)+BGENE(1,K)) K = K + 1 BGENE(1,K) = X1s2 * (BB(2,K)-BGENE(2,K)) BGENE(2,K) = X1s2 * (BB(2,K)+BGENE(2,K)) K = K + 1 ENDDO GOTO 666 *----- Axisymetrie (IFOU=0) 713 CONTINUE K = 1 DO i = 1,NBNO r_z = X1s3 * (BB(1,K)-BGENE(1,K)+BB(3,K)-BGENE(3,K)) BGENE(1,K) = BGENE(1,K) + r_z BGENE(2,K) = r_z BGENE(3,K) = BGENE(3,K) + r_z K = K + 1 r_z = X1s3 * (BB(2,K)-BGENE(2,K)) BGENE(1,K) = r_z BGENE(2,K) = BGENE(2,K) + r_z BGENE(3,K) = r_z K = K + 1 ENDDO GOTO 666 C----------------------------------------------------------------------- C----- Element massif bidimensionnel ICQ8 ------------------------------ C----------------------------------------------------------------------- 72 CONTINUE *W WRITE(*,*) 'Element ICQ8',IFOU QSI_z = QSIGAU(IGAU) ETA_z = ETAGAU(IGAU) IFR = IFOU+4 GOTO (666,721,721,723,666),IFR GOTO 999 *----- Contraintes planes ou deformations planes (IFOU=-2,-1) 721 CONTINUE *----- Calculs des composantes de BB-DILATATION BB(1, 1) = A(1, 1) + A(2, 1)*QSI_z + A(3, 1)*ETA_z BB(2, 2) = A(1, 2) + A(2, 2)*QSI_z + A(3, 2)*ETA_z BB(1, 3) = A(1, 3) + A(2, 3)*QSI_z + A(3, 3)*ETA_z BB(2, 4) = A(1, 4) + A(2, 4)*QSI_z + A(3, 4)*ETA_z BB(1, 5) = A(1, 5) + A(2, 5)*QSI_z + A(3, 5)*ETA_z BB(2, 6) = A(1, 6) + A(2, 6)*QSI_z + A(3, 6)*ETA_z BB(1, 7) = A(1, 7) + A(2, 7)*QSI_z + A(3, 7)*ETA_z BB(2, 8) = A(1, 8) + A(2, 8)*QSI_z + A(3, 8)*ETA_z BB(1, 9) = A(1, 9) + A(2, 9)*QSI_z + A(3, 9)*ETA_z BB(2,10) = A(1,10) + A(2,10)*QSI_z + A(3,10)*ETA_z BB(1,11) = A(1,11) + A(2,11)*QSI_z + A(3,11)*ETA_z BB(2,12) = A(1,12) + A(2,12)*QSI_z + A(3,12)*ETA_z BB(1,13) = A(1,13) + A(2,13)*QSI_z + A(3,13)*ETA_z BB(2,14) = A(1,14) + A(2,14)*QSI_z + A(3,14)*ETA_z BB(1,15) = A(1,15) + A(2,15)*QSI_z + A(3,15)*ETA_z BB(2,16) = A(1,16) + A(2,16)*QSI_z + A(3,16)*ETA_z *W WRITE (*,*) 'Ecriture de la matrice BB[ ]' *W WRITE (*,4) (('BB(',I,',',J,')=',BB(I,J),I=1,3),J=1,16) *W 4 FORMAT (3(A,I1,A,I1,A,F8.4,2X)) *----- Calcul des composantes de B-BARRE K = 1 DO i = 1,NBNO BGENE(2,K) = X1s2 * (BB(1,K)-BGENE(1,K)) BGENE(1,K) = X1s2 * (BB(1,K)+BGENE(1,K)) K = K + 1 BGENE(1,K) = X1s2 * (BB(2,K)-BGENE(2,K)) BGENE(2,K) = X1s2 * (BB(2,K)+BGENE(2,K)) K = K + 1 ENDDO GOTO 666 *----- Axisymetrie (IFOU=0) 723 CONTINUE *----- Calculs des composantes de BB-DILATATION BB(1, 1) = A(1, 1) + A(2, 1)*QSI_z + A(3, 1)*ETA_z BB(2, 2) = A(1, 2) + A(2, 2)*QSI_z + A(3, 2)*ETA_z BB(3, 1) = A(1, 3) + A(2, 3)*QSI_z + A(3, 3)*ETA_z BB(1, 3) = A(1, 4) + A(2, 4)*QSI_z + A(3, 4)*ETA_z BB(2, 4) = A(1, 5) + A(2, 5)*QSI_z + A(3, 5)*ETA_z BB(3, 3) = A(1, 6) + A(2, 6)*QSI_z + A(3, 6)*ETA_z BB(1, 5) = A(1, 7) + A(2, 7)*QSI_z + A(3, 7)*ETA_z BB(2, 6) = A(1, 8) + A(2, 8)*QSI_z + A(3, 8)*ETA_z BB(3, 5) = A(1, 9) + A(2, 9)*QSI_z + A(3, 9)*ETA_z BB(1, 7) = A(1,10) + A(2,10)*QSI_z + A(3,10)*ETA_z BB(2, 8) = A(1,11) + A(2,11)*QSI_z + A(3,11)*ETA_z BB(3, 7) = A(1,12) + A(2,12)*QSI_z + A(3,12)*ETA_z BB(1, 9) = A(1,13) + A(2,13)*QSI_z + A(3,13)*ETA_z BB(2,10) = A(1,14) + A(2,14)*QSI_z + A(3,14)*ETA_z BB(3, 9) = A(1,15) + A(2,15)*QSI_z + A(3,15)*ETA_z BB(1,11) = A(1,16) + A(2,16)*QSI_z + A(3,16)*ETA_z BB(2,12) = A(1,17) + A(2,17)*QSI_z + A(3,17)*ETA_z BB(3,11) = A(1,18) + A(2,18)*QSI_z + A(3,18)*ETA_z BB(1,13) = A(1,19) + A(2,19)*QSI_z + A(3,19)*ETA_z BB(2,14) = A(1,20) + A(2,20)*QSI_z + A(3,20)*ETA_z BB(3,13) = A(1,21) + A(2,21)*QSI_z + A(3,21)*ETA_z BB(1,15) = A(1,22) + A(2,22)*QSI_z + A(3,22)*ETA_z BB(2,16) = A(1,23) + A(2,23)*QSI_z + A(3,23)*ETA_z BB(3,15) = A(1,24) + A(2,24)*QSI_z + A(3,24)*ETA_z *W WRITE (*,*) 'Ecriture de la matrice BB[ ]' *W WRITE (*,4) (('BB(',I,',',J,')=',BB(I,J),I=1,3),J=1,16) *----- Calcul des composantes de B-BARRE K = 1 DO i = 1,NBNO r_z = X1s3 * (BB(1,K)-BGENE(1,K)+BB(3,K)-BGENE(3,K)) BGENE(1,K) = BGENE(1,K) + r_z BGENE(2,K) = r_z BGENE(3,K) = BGENE(3,K) + r_z K = K + 1 r_z = X1s3 * (BB(2,K)-BGENE(2,K)) BGENE(1,K) = r_z BGENE(2,K) = BGENE(2,K) + r_z BGENE(3,K) = r_z K = K + 1 ENDDO GOTO 666 C----------------------------------------------------------------------- C----- Element massif tridimensionnel ICC8 ----------------------------- C----------------------------------------------------------------------- 73 CONTINUE *W WRITE(*,*) 'Element ICC8',IFOU IF (IFOU.NE.2) GOTO 999 QSI_z = QSIGAU(IGAU) ETA_z = ETAGAU(IGAU) DZE_z = DZEGAU(IGAU) GOTO 666 C----------------------------------------------------------------------- C----- Element massif tridimensionnel ICT4 ----------------------------- C----------------------------------------------------------------------- 74 CONTINUE *W WRITE(*,*) 'Element ICT4',IFOU IF (IFOU.NE.2) GOTO 999 QSI_z = QSIGAU(IGAU) ETA_z = ETAGAU(IGAU) DZE_z = DZEGAU(IGAU) GOTO 666 C----------------------------------------------------------------------- C----- Element massif tridimensionnel ICP6 ----------------------------- C----------------------------------------------------------------------- 75 CONTINUE *W WRITE(*,*) 'Element ICP6',IFOU IF (IFOU.NE.2) GOTO 999 QSI_z = QSIGAU(IGAU) ETA_z = ETAGAU(IGAU) DZE_z = DZEGAU(IGAU) GOTO 666 C----------------------------------------------------------------------- C----- Element massif tridimensionnel IC20 ----------------------------- C----------------------------------------------------------------------- 76 CONTINUE *W WRITE(*,*) 'Element IC20',IFOU IF (IFOU.NE.2) GOTO 999 QSI_z = QSIGAU(IGAU) ETA_z = ETAGAU(IGAU) DZE_z = DZEGAU(IGAU) GOTO 666 C----------------------------------------------------------------------- C----- Element massif tridimensionnel IC10 ----------------------------- C----------------------------------------------------------------------- 77 CONTINUE *W WRITE(*,*) 'Element IC10',IFOU IF (IFOU.NE.2) GOTO 999 QSI_z = QSIGAU(IGAU) ETA_z = ETAGAU(IGAU) DZE_z = DZEGAU(IGAU) GOTO 666 C----------------------------------------------------------------------- C----- Element massif tridimensionnel IC15 ----------------------------- C----------------------------------------------------------------------- 78 CONTINUE *W WRITE(*,*) 'Element IC15',IFOU IF (IFOU.NE.2) GOTO 999 QSI_z = QSIGAU(IGAU) ETA_z = ETAGAU(IGAU) DZE_z = DZEGAU(IGAU) GOTO 666 C----------------------------------------------------------------------- C----- ERREUR : Donnees incompatibles ---------------------------------- C----------------------------------------------------------------------- 999 CONTINUE *W WRITE(*,*) 'Mode de calcul ',IFOU,' incompatible avec ', *W & 'element ',MELE CALL ERREUR(5) C----------------------------------------------------------------------- C----- FIN du sousprogramme -------------------------------------------- C----------------------------------------------------------------------- 666 CONTINUE *W WRITE (*,*) 'Sortie du sousprogramme BBAR' RETURN END