1 : C QUAHM2 SOURCE CHAT 05/01/13 02:40:56 5004 2 : C QUAHM2 SOURCE AM1 95/11/24 22:54:38 1918 3 : SUBROUTINE QUAHM2(IGAU,ITEL,MFR,NBNO,XEL,SHPTOT,SHP,IFOU, 4 : # NHARM,VKL12,VKL23,VKL33,POIGAU,ISDJC,LRE,REL,IRET) 5 : C======================================================================= 6 : C 7 : C CALCULE LES TERMES EN P * PI ,PI * (UR,RT) ,(UR,RT) *(UR,RT) 8 : C (UT,RR) * (UT,RR) , PI * (UT,RR) DE LA MATRICE 9 : C MASSE DANS LE CAS AXISYMETRIQUE OU FOURIER POUR 10 : C LA FORMULATION (37) HOMOGENE 11 : C======================================================================= 12 : C INPUT 13 : C IGAU=NUMERO DU POINT DE GAUSS 14 : C ITEL=NUMERO DE L ELEMENT DANS NOMTP 15 : C MFR =NUMERO DE LA FORMULATION 16 : C NBNO=NOMBRE DE NOEUDS 17 : C XEL =COORDONNEES DE L ELEMENT 18 : C IFOU=IFOUR DE CCOPTIO 19 : C NHARM=NUMERO DU MODE DE FOURIER 20 : C VKL12=-((COEFPI*COEFPR)/(RHOF*C**2))*SFLU/SCEL 21 : C VKL23=(BET11+BET22)*COEFPI/(2.*SCEL) 22 : C VKL33=(RHOS*2.+RHOF*(BET11+BET22))/SCEL 23 : C POIGAU=MINTE.POIGAU(IGAU) 24 : C LRE =NOMBRE DE D.D.L DE LA MATRICE DE RIGIDITE 25 : C SHPTOT(6,NBNO,NBGAU)=FONCTIONS DE FORMES ET DERIVEES 26 : C ISDJC = INDICATEUR SUR LE SIGNE DU JACOBIEN 27 : C ZONE DE TRAVAIL 28 : C SHP(5,NBNO)=TABLEAU DE TRAVAIL 29 : C OUTPUT 30 : C ISDJC = INDICATEUR SUR LE SIGNE DU JACOBIEN 31 : C REL=MATRICE DE MASSE 32 : C IRET:INDICATEUR = 1 : SUCCES 33 : C 0 : ECHEC (ELEMENT MELE INCOMPATIBLE ) 34 : C 2 : ECHEC (JACOBIEN NUL ) 35 : C 3 :ECHEC (ROUTINE N EST VALABLE QU 36 : C EN AXISYMETRIQUE OU FOURIER ) 37 : C======================================================================= 38 : IMPLICIT INTEGER(I-N) 39 : IMPLICIT REAL*8(A-H,O-Z) 40 : DIMENSION XEL(3,*),SHP(6,*),SHPTOT(6,NBNO,*),REL(LRE,*) 41 : IF (ITEL.EQ.126) GOTO 10 42 : C 43 : C ERREUR : TYPE D' ELEMENT INCOMPATIBLE AVEC LA FORMULATION 44 : C 45 : IRET = 0 46 : GOTO 666 47 : 10 CONTINUE 48 : IF (IFOU.EQ.0.OR.IFOU.EQ.1) GOTO 11 49 : C 50 : C MESSAGE D ERREUR : ROUTINE N EST VALABLE QU EN FOURIER 51 : C OU EN AXISYMETRIQUE 52 : C 53 : IRET = 3 54 : GOTO 666 55 : 11 CONTINUE 56 : C 57 : C ELEMENTS HOMOGENEISES QUAH EN AXISYMETRIE OU EN FOURIER 58 : C NBDL = LRE/NBNO NOMBRE DE D.D.L PAR NOEUD 59 : C 60 : NBDL = LRE /NBNO 61 : C 62 : C SHP(1,I) : FONCTION DE FORME 63 : C SHP(2,I) : DERIVEE % R DE LA FONCTION DE FORME 64 : C SHP(3,I) : DERIVEE % Z DE LA FONCTION DE FORME 65 : C 66 : DO 101 NP=1,NBNO 67 : SHP(1,NP)=SHPTOT(1,NP,IGAU) 68 : SHP(2,NP)=SHPTOT(2,NP,IGAU) 69 : SHP(3,NP)=SHPTOT(3,NP,IGAU) 70 : 101 CONTINUE 71 : C 71 : CALL DEVOLU(XEL,SHP,MFR,NBNO,IFOU,NHARM,2,1.D0,RR,DJAC) 72 : IF (DJAC.EQ.0.) GOTO 667 73 : IF (DJAC.LT.0.) ISDJC = ISDJC + 1 74 : C 75 : C FONCTIONS DE FORME 76 : C 77 : SHP(4,1)=SHP(1,1) + SHP(1,4) 78 : SHP(4,2)=SHP(1,2) + SHP(1,3) 79 : SHP(4,3)=0.D0 80 : SHP(4,4)=0.D0 81 : C 82 : XZH=SHP(1,3) + SHP(1,4) 83 : C 84 : A1=MIN(XEL(2,1),XEL(2,2)) 85 : A2=MIN(XEL(2,1),XEL(2,3)) 86 : A3=MIN(XEL(2,1),XEL(2,4)) 87 : A4=MIN(A1,A2) 88 : A5=MIN(A1,A3) 89 : Z1=MIN(A4,A5) 90 : C 91 : B1=MAX(XEL(2,1),XEL(2,2)) 92 : B2=MAX(XEL(2,1),XEL(2,3)) 93 : B3=MAX(XEL(2,1),XEL(2,4)) 94 : B4=MAX(B1,B2) 95 : B5=MAX(B1,B3) 96 : Z2=MAX(B4,B5) 97 : C 98 : RH=Z2-Z1 99 : C 100 : S1=SHP(1,1) + SHP(1,4) 101 : S2=SHP(1,2) + SHP(1,3) 102 : C 103 : C FONCTIONS DE FORME EN Z 104 : C 105 : XZH= SHP(1,3) + SHP(1,4) 106 : C 107 : SHP(5,1)=1.D0-(3.D0*XZH*XZH)+(2.D0*XZH*XZH*XZH) 108 : SHP(5,2)=3.D0*(XZH**2)-(2.D0*XZH*XZH*XZH) 109 : SHP(5,3)=(RH*XZH)*(1.D0-2.D0*XZH+XZH*XZH) 110 : SHP(5,4)=(RH*XZH)*(XZH*XZH-XZH) 111 : C 112 : C FONCTIONS DE FORME POUR LA FLEXION 113 : C 114 : SHP(6,1)=S1*SHP(5,1) 115 : SHP(6,2)=S2*SHP(5,1) 116 : SHP(6,3)=S2*SHP(5,2) 117 : SHP(6,4)=S1*SHP(5,2) 118 : C 119 : SHP(6,5)=S1*SHP(5,3) 120 : SHP(6,6)=S2*SHP(5,3) 121 : SHP(6,7)=S2*SHP(5,4) 122 : SHP(6,8)=S1*SHP(5,4) 123 : C 124 : C 125 : C1=MIN(XEL(1,1),XEL(1,2)) 126 : C2=MIN(XEL(1,1),XEL(1,3)) 127 : C3=MIN(XEL(1,1),XEL(1,4)) 128 : C4=MIN(C1,C2) 129 : C5=MIN(C1,C3) 130 : R1=MIN(C4,C5) 131 : C 132 : D1=MAX(XEL(1,1),XEL(1,2)) 133 : D2=MAX(XEL(1,1),XEL(1,3)) 134 : D3=MAX(XEL(1,1),XEL(1,4)) 135 : D4=MAX(D1,D2) 136 : D5=MAX(D1,D3) 137 : R2=MAX(D4,D5) 138 : C 139 : DELTAR = R2 - R1 140 : DV = (DELTAR*RH)/4.d0 141 : C 142 : DJAC = ABS(DJAC) 143 : C 144 : C 145 : C TERMES EN P * PI 146 : C 147 : DJAC1 = ABS(DJAC)*POIGAU 148 : IX1=0 149 : IY1=0 150 : DO 102 IX=2,LRE ,NBDL 151 : IX1=IX1 + 1 152 : DO 103 IY=1,IX ,NBDL 153 : IY1=IY1 + 1 154 : REL(IY,IX) = REL(IY,IX) + VKL12*DJAC1*SHP(1,IX1)*SHP(1,IY1) 155 : REL(IX,IY) = REL(IY,IX) 156 : 103 CONTINUE 157 : IY1=0 158 : 102 CONTINUE 159 : DO 104 IX=2+NBDL,LRE ,NBDL 160 : IX2=IX - NBDL 161 : DO 105 IY=1,IX2 ,NBDL 162 : REL(IY+1,IX-1) = REL(IY,IX) 163 : REL(IX-1,IY+1) = REL(IY+1,IX-1) 164 : 105 CONTINUE 165 : 104 CONTINUE 166 : C 167 : C TERMES EN PI * (UR , RT ) 168 : C 169 : C 170 : IX1=0 171 : IY1=0 172 : DO 106 IX=3,LRE ,NBDL 173 : IX1=IX1 + 1 174 : DO 107 IY=2,IX ,NBDL 175 : IY1=IY1 + 1 176 : REL(IY,IX) = REL(IY,IX) + VKL23*DJAC1*SHP(2,IY1) 177 : #*SHP(6,IX1) 178 : REL(IY,IX+1) = REL(IY,IX+1) + VKL23*DJAC1*SHP(2,IY1) 179 : #*SHP(6,IX1+4) 180 : REL(IX,IY) = REL(IY,IX) 181 : REL(IX+1,IY) = REL(IY,IX+1) 182 : 107 CONTINUE 183 : IY1=0 184 : 106 CONTINUE 185 : IX1=1 186 : IY1=0 187 : DO 108 IX=2+NBDL,LRE ,NBDL 188 : IX1=IX1 + 1 189 : DO 109 IY=3,IX ,NBDL 190 : IY1=IY1 + 1 191 : REL(IY,IX) = REL(IY,IX) + VKL23*DJAC1*SHP(2,IX1) 192 : #*SHP(6,IY1) 193 : REL(IY+1,IX) = REL(IY+1,IX) + VKL23*DJAC1*SHP(2,IX1) 194 : #*SHP(6,IY1+4) 195 : REL(IX,IY) = REL(IY,IX) 196 : REL(IX,IY+1) = REL(IY+1,IX) 197 : 109 CONTINUE 198 : IY1=0 199 : 108 CONTINUE 200 : C 201 : IF ( IFOU.EQ.1) THEN 202 : C 203 : C TERMES EN PI * (UT , RR ) 204 : C NON NULS QU EN FOURIER 205 : C 206 : DJAC2 = ABS(DJAC)*POIGAU 207 : C 208 : VKL25 = -1.D0* VKL23*NHARM / RR 209 : IX1=0 210 : IY1=0 211 : DO 110 IX=5,LRE ,NBDL 212 : IX1=IX1 + 1 213 : DO 111 IY=2,IX ,NBDL 214 : IY1=IY1 + 1 215 : REL(IY,IX) = REL(IY,IX) + VKL25*DJAC2*SHP(1,IY1)*SHP(6,IX1) 216 : REL(IY,IX+1) = REL(IY,IX+1) + VKL25*DJAC2*SHP(1,IY1) 217 : #*SHP(6,IX1+4)*(-1.D0) 218 : REL(IX,IY) = REL(IY,IX) 219 : REL(IX+1,IY) = REL(IY,IX+1) 220 : 111 CONTINUE 221 : IY1=0 222 : 110 CONTINUE 223 : IX1=1 224 : IY1=0 225 : DO 112 IX=2+NBDL,LRE ,NBDL 226 : IX1=IX1 + 1 227 : DO 113 IY=5,IX ,NBDL 228 : IY1=IY1 + 1 229 : REL(IY,IX) = REL(IY,IX) + VKL25*DJAC2*SHP(1,IX1)*SHP(6,IY1) 230 : REL(IY+1,IX) = REL(IY+1,IX) + VKL25*DJAC2*SHP(1,IX1) 231 : #*SHP(6,IY1+4)*(-1.D0) 232 : REL(IX,IY) = REL(IY,IX) 233 : REL(IX,IY+1) = REL(IY+1,IX) 234 : 113 CONTINUE 235 : IY1=0 236 : 112 CONTINUE 237 : ENDIF 238 : C 239 : C TERMES EN (UR,RT ) * (UR , RT ) 240 : C 241 : C 242 : IX1=0 243 : IY1=0 244 : DO 114 IX=3,LRE ,NBDL 245 : IX1=IX1 + 1 246 : DO 115 IY=3,IX ,NBDL 247 : IY1=IY1 + 1 248 : C 249 : REL(IY,IX) = REL(IY,IX) + VKL33*DJAC1*SHP(6,IY1)*SHP(6,IX1) 250 : REL(IY,IX+1) = REL(IY,IX+1) + VKL33*DJAC1*SHP(6,IY1)*SHP(6,IX1+4) 251 : REL(IY+1,IX) = REL(IY+1,IX) + VKL33*DJAC1*SHP(6,IY1+4)*SHP(6,IX1) 252 : REL(IY+1,IX+1) = REL(IY+1,IX+1)+VKL33*DJAC1*SHP(6,IY1+4) 253 : #*SHP(6,IX1+4) 254 : REL(IX,IY) = REL(IY,IX) 255 : REL(IX+1,IY) = REL(IY,IX+1) 256 : REL(IX,IY+1) = REL(IY+1,IX) 257 : REL(IX+1,IY+1) = REL(IY+1,IX+1) 258 : 115 CONTINUE 259 : IY1=0 260 : 114 CONTINUE 261 : C 262 : IF ( IFOU.EQ.1) THEN 263 : C 264 : C TERMES EN (UT,RR ) * (UT , RR ) 265 : C NON NULS QU EN FOURIER 266 : C 267 : IX1=0 268 : IY1=0 269 : DO 303 IX=5,LRE ,NBDL 270 : IX1=IX1 + 1 271 : DO 403 IY=5,IX ,NBDL 272 : IY1=IY1 + 1 273 : REL(IY,IX) = REL(IY,IX) + VKL33*DJAC1*SHP(6,IY1)*SHP(6,IX1) 274 : REL(IY,IX+1) = REL(IY,IX+1) + VKL33*DJAC1*SHP(6,IY1)*SHP(6,IX1+4) 275 : #*(-1.D0) 276 : REL(IY+1,IX) = REL(IY+1,IX) + VKL33*DJAC1*SHP(6,IY1+4)*SHP(6,IX1) 277 : #*(-1.D0) 278 : REL(IY+1,IX+1) = REL(IY+1,IX+1)+VKL33*DJAC1*SHP(6,IY1+4) 279 : #*SHP(6,IX1+4) 280 : REL(IX,IY) = REL(IY,IX) 281 : REL(IX+1,IY) = REL(IY,IX+1) 282 : REL(IX,IY+1) = REL(IY+1,IX) 283 : REL(IX+1,IY+1) = REL(IY+1,IX+1) 284 : 403 CONTINUE 285 : IY1=0 286 : 303 CONTINUE 287 : C 288 : ENDIF 289 : IRET = 1 290 : GOTO 666 291 : C 292 : C MESSAGE D ERREUR : ELEMENT A SURFACE NULLE 293 : C 294 : 667 CONTINUE 295 : IRET = 2 296 : GOTO 666 297 : C 298 : 666 CONTINUE 299 : RETURN 300 : END 301 : 302 : 303 : 304 :
© Cast3M 2003 - All rights reserved.
Disclaimer