1 : C MHYBR4 SOURCE CHAT 11/03/16 21:27:58 6902 2 : SUBROUTINE MHYBR4(IMAIL,IPMODE,IPRIGI,IPGEOS) 3 : C----------------------------------------------------------------------- 4 : C Calcul de la matrice masse hybride utilisé pour prendre en compte les forces 5 : C de volume. 6 : C Traitement du cas des elements finis massifs a integration numerique 7 : C pour un maillage elementaire et une formulation hybride. 8 : C----------------------------------------------------------------------- 9 : C 10 : C--------------------------- 11 : C Parametres Entree/Sortie : 12 : C--------------------------- 13 : C 14 : C E/ IMAIL : Numero du maillage elementaire considere, 15 : C dans l'objet modele. 16 : C E/ IPMODE : Pointeur sur un segment IMODEL. 17 : C E/ IPGEOS : Pointeur sur le maillage sommet 18 : C E/S IPRIGI : Pointeur sur l'objet rigidite resultat. 19 : C 20 : C---------------------- 21 : C Variables en COMMON : 22 : C---------------------- 23 : C 24 : C E/ XCOOR : VOIR SMCOORD 25 : C E/ IERR : VOIR CCOPTIO 26 : C E/ IDIM : VOIR CCOPTIO 27 : C E/ INTERR : VOIR CCOPTIO 28 : C E/ IFOMOD : VOIR CCOPTIO 29 : C E/ XPETIT : VOIR CCREEL 30 : C 31 : C---------------------- 32 : C Tableaux de travail : 33 : C---------------------- 34 : C 35 : C NBNN : Nombre de noeuds dans l'element considere 36 : C NEFHYB : Numero de l'element fini dans NOMTP. 37 : C NEF : Numero de l'element fini support geometrique 38 : C dans NOMTP (voir CCHAMP) 39 : C NBELEM : Nombre d'element dans le maillage elementaire 40 : C NBPGAU : Nombre de points de gauss pour l'element fini NEF 41 : C** CEL : Matrice de conductivite elementaire 42 : C CEL : Matrice des masses elementaires 43 : C XE : Coordonnees des noeuds dans le repere global 44 : C** CMAT : Matrice de permeabilite dans le repere global 45 : C SHP : Tableau de travail contenant les fonctions de forme au 46 : C point de gauss utilise + derivees 47 : C SHY : Contient les fonctions de base hybride en un point 48 : C mais pas les derivees de la fonction de base. 49 : C** VALMAT : Valeurs des coeff. de la matrice CMAT et des 50 : C** cosinus directeurs des axes d'ortho. / repere local 51 : C** XGLOB : Cosinus directeurs des axes d'ortho. / repere global 52 : C** XLOC : Cosinus directeurs des axes d'ortho. / repere local 53 : C** TXR : Cosinus directeurs des axes locaux / repere global 54 : C 55 : C 56 : C----------------------------------------------------------------------- 57 : C 58 : C Langage : ESOPE + FORTRAN77 59 : C 60 : C 02/96 L.V.BENET : fonction propre à l'option 'MASSE' 61 : C 62 : C----------------------------------------------------------------------- 63 : IMPLICIT INTEGER(I-N) 64 : IMPLICIT REAL*8(A-H,O-Z) 65 : * 66 : -INC CCHAMP 67 : -INC CCOPTIO 68 : -INC CCREEL 69 : -INC SMCOORD 70 : -INC SMINTE 71 : -INC SMMODEL 72 : -INC SMRIGID 73 : -INC SMELEME 74 : ****-INC SMCHAML 75 : * 76 : * 77 : SEGMENT MMAT1 78 : REAL*8 CEL(NBDDL,NBDDL),CEL1(NBDDL,NBDDL),XE(3,NBNN) 79 : REAL*8 SHP(6,NBNN),SHY(IDIM,NBDDL),ZJAC(IDIM,IDIM) 80 : REAL*8 CMAT(IDIM,IDIM),CMAT1(IDIM,IDIM),CMAT2(IDIM,IDIM) 81 : INTEGER ICSTO(NBDDL) 82 : ENDSEGMENT 83 : * 84 : SEGMENT NOTYPE 85 : CHARACTER*16 TYPE(NBTYPE) 86 : ENDSEGMENT 87 : * 88 : SEGMENT HYBSTO 89 : REAL*8 HYBASE(NDIM,NBDDL,NBPP) 90 : ENDSEGMENT 91 : * 92 : CHARACTER*8 CNM 93 : PARAMETER(NINF=3) 94 : INTEGER INFOS(NINF) 95 : * 96 : * Recup. des caracteristiques geometriques du maillage elementaire 97 : * et du maillage hybride dual 98 : * 99 : IMODEL = IPMODE 100 : SEGACT IMODEL 101 : IPMAIL = IMAMOD 102 : MELEME = IPGEOS 103 : SEGACT MELEME 104 : NBNN = NUM(/1) 105 : NBELEM = NUM(/2) 106 : NEFHYB = NEFMOD 107 : NEF = NUMGEO(NEFHYB) 108 : MFR = NUMMFR(NEFHYB) 109 : * 110 : MRIGID = IPRIGI 111 : SEGACT MRIGID 112 : IPT1 = IRIGEL(1,IMAIL) 113 : SEGACT IPT1 114 : NBDDL = IPT1.NUM(/1) 115 : SEGDES IPT1 116 : * 117 : * Recup. des caracteristiques d'integration de l'EF support geometrique 118 : * de l'EF hybride 119 : * 119 : CALL TSHAPE(NEF,'GAUSS',IPINTE) 120 : IF (IERR.NE.0) THEN 121 : SEGDES IMODEL , MELEME 122 : RETURN 123 : ENDIF 124 : * 125 : * Recup. des fonctions de bases hybrides 126 : * 126 : CALL HYSHPT(NEFHYB,NBDDL,IPINTE,IPTHYB) 127 : IF (IERR.NE.0) THEN 128 : SEGDES IMODEL , MELEME 129 : RETURN 130 : ENDIF 131 : * 132 : * Activation des segments "d'integration" 133 : * 134 : MINTE = IPINTE 135 : SEGACT MINTE 136 : NBPGAU = POIGAU(/1) 137 : HYBSTO = IPTHYB 138 : SEGACT HYBSTO 139 : * 140 : * Recup. des caracteristiques d'integration au centre de l'EF 141 : * 141 : CALL RESHPT(1,NBNN,NEF,NEF,0,IPT1,IRT1) 142 : MINTE1 = IPT1 143 : SEGACT MINTE1 144 : * 145 : * Initialisation des chapeaux de l'objet rigidité 146 : * 147 : xMATRI = IRIGEL(4,IMAIL) 148 : SEGACT xMATRI*MOD 149 : NLIGRP = NBDDL 150 : NLIGRD = NBDDL 151 : * 152 : * Remplissage du tableau INFOS (informations sur element). 153 : * 154 : INFOS(1) = 0 155 : INFOS(2) = 0 156 : INFOS(3) = NIFOUR 157 : * 158 : * Initialisation des tableaux de travail 159 : * 160 : NDIM = IDIM * (IDIM+1) 161 : SEGINI MMAT1 161 : CALL ZERO(CMAT2,IDIM,IDIM) 162 : DO 5 I=1,IDIM 163 : CMAT2(I,I)=1.D0 164 : 5 CONTINUE 165 : * 166 : *------------------------------------------------------- 167 : * BOUCLE SUR LES ELEMENTS DU MAILLAGE ELEMENTAIRE IMAIL 168 : *------------------------------------------------------- 169 : * 170 : DO 30 IEL=1,NBELEM 171 : * 172 : * Initialisations 173 : * 174 : IFOIS = 0 175 : IFOI2 = 0 175 : CALL ZERO(CEL,NBDDL,NBDDL) 176 : * 177 : * Recuperation des coordonnees globales des noeuds de l'element IEL 178 : * 178 : CALL DOXE(XCOOR,IDIM,NBNN,NUM,IEL,XE) 179 : *-------------------------------- 180 : * BOUCLE SUR LES POINTS DE GAUSS 181 : *-------------------------------- 182 : DO 20 IGAU=1,NBPGAU 183 : * 184 : * Calcul de la matrice jacobienne de la fonction de passage du 185 : * repere local au repere global, de son determinant ET recup. 186 : * des fonctions de base hybride au point de gauss IGAU. 187 : * 187 : CALL MHYBR3(IGAU,NBNN,NBDDL,NDIM,IDIM,IDIM,XE,HYBASE, 188 : S SHPTOT,SHY,SHP,ZJAC,DJAC) 189 : * 190 : * Controle du maillage 191 : * 192 : IF (DJAC.LT.0.D0) IFOIS = IFOIS + 1 193 : IF (ABS(DJAC).LT.XPETIT) THEN 194 : IFOI2 = IFOI2 + 1 195 : DJAC = XPETIT 196 : ENDIF 197 : * 198 : * Calcul du poids d'integration global affecte dans DJAC. 199 : * 200 : DJAC = POIGAU(IGAU) / ABS(DJAC) 201 : * 202 : *- Calcul de la contribution du point de gauss IGAU a la matrice 203 : *- elementaire CEL de l'element IEL : 204 : *- POIGAU/DJAC * transpose( ZJAC*SHY ) *inv(CMAT)* ( ZJAC*SHY ) 205 : *- On calcule CMAT2=inv(CMAT) avec INVRS puis 206 : *- on calcule CMAT1=transpose(ZJAC)*CMAT2*ZJAC avec PRODT puis 207 : *- on somme POIGAU/DJAC * transp.(SHY)*CMAT1*SHY avec BDBST. 208 : * 208 : CALL PRODT(CMAT1,CMAT2,ZJAC,IDIM,IDIM) 208 : CALL BDBST(SHY,DJAC,CMAT1,NBDDL,IDIM,CEL) 209 : 20 CONTINUE 210 : * 211 : * Le jacobien est negatif --> MAILLAGE INCORRECT 212 : * 213 : IF (IFOIS.NE.0.AND.IFOIS.NE.NBPGAU) THEN 214 : INTERR(1) = IEL 214 : CALL ERREUR(195) 215 : SEGSUP xMATRI , MRIGID 216 : GOTO 40 217 : ENDIF 218 : * 219 : * Le jacobien est tres petit --> MAILLAGE INCORRECT 220 : * 221 : IF (IFOI2.EQ.NBPGAU) THEN 222 : INTERR(1) = IEL 222 : CALL ERREUR(259) 223 : SEGSUP xMATRI , MRIGID 224 : GOTO 40 225 : ENDIF 226 : * 227 : * Remplissage de XMATRI 228 : * 229 : * SEGINI XMATRI 230 : * IMATTT(IEL) = XMATRI 230 : CALL REMPMT(CEL,NBDDL,RE(1,1,iel)) 231 : * SEGDES XMATRI 232 : 30 CONTINUE 233 : * 234 : * Desactivation des segments 235 : * 236 : SEGDES xMATRI , MRIGID 237 : 40 CONTINUE 238 : SEGSUP MMAT1 , HYBSTO 239 : SEGDES MELEME 240 : SEGDES IMODEL 241 : SEGDES MINTE , MINTE1 242 : RETURN 243 : END 244 : 245 : 246 :
© Cast3M 2003 - All rights reserved.
Disclaimer