1 : C ELLP21 SOURCE CHAT 05/01/12 23:35:21 5004 2 : SUBROUTINE ELLP21 (CARACT,INP,ZS,COOR,GAMA,ICORR,ZC1) 3 : C 4 : IMPLICIT INTEGER(I-N) 5 : IMPLICIT REAL*8 (A-H,O-Y) 6 : IMPLICIT COMPLEX*16 (Z) 7 : C 8 : C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 9 : C 10 : C OPERATEUR ELFE LAPLACE POUTRE 11 : C 12 : C CALCUL POUR LA POUTRE N› INP DE LA MATRICE CI TEL QUE CI * U = SM 13 : C 14 : C PARAMETRES : 15 : C 16 : C CARACT : TABLEAU DES CARACTERISTIQUE DES POUTRES (9,NP) 17 : C INP : NUMERO DE LA POUTRE 18 : C ZS : VALEUR DE S = S0 + I W 19 : C COOR : TABLEAU DES COORDONNEES DES NOEUDS (3,2*NP) 20 : C GAMA : TABLEAU DES VECTEUR DEFINISSANT LE PLAN LOCALE OXY (3,2*NP) 21 : C ICORR : TABLEAU DE CORRESP. DE NOUEDS REELS <--> NOEUDS FICTIFS 22 : C (2*NP) 23 : C SORTIES : 24 : C 25 : C ZC1 : MATRICE 12*24 CALCULEE EN DIMENSIONNEL DANS LE REPERE 26 : C GLOBLAL 27 : C 28 : C 29 : C AUTEUR : SAINT-DIZIER 30 : C DATE : 30 NOVEMBRE 1989 (VERSION DU 03 AVRIL 1990) 31 : C 32 : C CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 33 : C 34 : DIMENSION CARACT(12,*),COOR(3,*),GAMA(3,*),ICORR(*) 35 : DIMENSION ZC1(12,24) 36 : DIMENSION ICHANG(24) 37 : DIMENSION ZR(3,3),ZD(12),ZE(12,12),ZA(12,24),ZA1(12,24) 38 : DIMENSION ZER(24,24) 39 : C 40 : C -------------- ACQUISITION DES CARACTERISTIQUES DE LA POUTRE 41 : C --------------------------------------------- 42 : CE = CARACT( 1,INP) 43 : CNU = CARACT( 2,INP) 44 : CRO = CARACT( 3,INP) 45 : CSE = CARACT( 4,INP) 46 : CC = CARACT( 5,INP) 47 : CIP = CARACT( 6,INP) 48 : CIY = CARACT( 7,INP) 49 : CIZ = CARACT( 8,INP) 50 : CKCY = CARACT( 9,INP) 51 : CKCZ = CARACT(10,INP) 52 : CAM = CARACT(11,INP) 53 : CETA = CARACT(12,INP) 54 : C 55 : ZI = (0.D0 , 1.D0) 56 : C 57 : ZCE = CE*(1.D0+CETA*ZI) 58 : ZCP = SQRT(ZCE/CRO) 59 : ZCMU = ZCE/(2.D0*(1.D0+CNU)) 60 : ZCR = SQRT(ZCMU/CRO) 61 : C 62 : ZALFA2= ZS*ZS*CRO*CSE + CAM*ZS 63 : ZALFA2= ZALFA2 / (ZCE*CSE) 64 : ZALFAX= SQRT(ZALFA2) 65 : C 66 : ZALFA2= ZS*ZS*CRO*CIP + CAM*ZS 67 : ZALFA2= ZALFA2 / (ZCMU*CC) 68 : ZALFAR= SQRT(ZALFA2) 69 : C 70 : ZALFA4= CRO*CRO*CIZ*ZS*ZS*ZS*ZS/ZCMU/CKCY 71 : * + CAM*CRO*CIZ*ZS*ZS*ZS/ZCMU/CKCY/CSE 72 : * + CRO*CSE*ZS*ZS 73 : * + CAM*ZS 74 : ZALFA4= ZALFA4/(ZCE*CIZ) 75 : ZALFA2= SQRT(ZALFA4) 76 : ZALFAY= SQRT(ZALFA2) 77 : C 78 : ZZ3 = CRO*CIZ*(1.D0+ZCE/(ZCMU*CKCY))*ZS*ZS 79 : * + CAM*ZCE*CIZ*ZS/(ZCMU*CKCY*CSE) 80 : ZZ3 = ZZ3 / (ZCE*CIZ*ZALFA2) 81 : C 82 : ZALFA4= CRO*CRO*CIY*ZS*ZS*ZS*ZS/ZCMU/CKCZ 83 : * + CAM*CRO*CIY*ZS*ZS*ZS/ZCMU/CKCZ/CSE 84 : * + CRO*CSE*ZS*ZS 85 : * + CAM*ZS 86 : ZALFA4= ZALFA4/(ZCE*CIY) 87 : ZALFA2= SQRT(ZALFA4) 88 : ZALFAZ= SQRT(ZALFA2) 89 : C 90 : ZZ4 = CRO*CIY*(1.D0+ZCE/(ZCMU*CKCZ))*ZS*ZS 91 : * + CAM*ZCE*CIY*ZS/(ZCMU*CKCZ*CSE) 92 : ZZ4 = ZZ4 / (ZCE*CIY*ZALFA2) 93 : C 94 : N1 = 2*INP-1 95 : N2 = 2*INP 96 : X1 = COOR (1,N1) 97 : Y1 = COOR (2,N1) 98 : Z1 = COOR (3,N1) 99 : X2 = COOR (1,N2) 100 : Y2 = COOR (2,N2) 101 : Z2 = COOR (3,N2) 102 : C 103 : XL = SQRT((X2-X1)**2 + (Y2-Y1)**2 + (Z2-Z1)**2) 104 : C 105 : C ------------------------------ VECTEUR UNITAIRE OX REPERE LOCALE 106 : C --------------------------------- 107 : XI1 = (X2-X1)/XL 108 : XI2 = (Y2-Y1)/XL 109 : XI3 = (Z2-Z1)/XL 110 : C 111 : C ------------------------- VECTEUR UNITAIRE OY REPERE LOCALE 112 : C --------------------------------- 113 : GX = GAMA(1,INP) 114 : GY = GAMA(2,INP) 115 : GZ = GAMA(3,INP) 116 : GG = SQRT(GX*GX + GY*GY + GZ*GZ) 117 : GX = GX/GG 118 : GY = GY/GG 119 : GZ = GZ/GG 120 : C 121 : DELTA = SQRT (1.D0 - (XI1*GX + XI2*GY + XI3*GZ)**2) 122 : C 123 : DET=-(GX*XI3-GZ*XI1)**2-(GY*XI1-GX*XI2)**2-(GY*XI3-GZ*XI2)**2 124 : C 125 : IF (ABS(DET).LT.1E-12) THEN 126 : XJ1 = -XI2 127 : XJ2 = 0. 128 : XJ3 = 0. 129 : ELSE 130 : XJ1 = (XI2*(GY*XI1-GX*XI2)-XI3*(GX*XI3-GZ*XI1))*DELTA/DET 131 : XJ2 = (XI3*(GZ*XI2-GY*XI3)-XI1*(GY*XI1-GX*XI2))*DELTA/DET 132 : XJ3 = (XI1*(GX*XI3-GZ*XI1)-XI2*(GZ*XI2-GY*XI3))*DELTA/DET 133 : ENDIF 134 : C 135 : C ---------------------------- VECTEUR UNITAIRE OZ REPERE LOCALE 136 : C --------------------------------- 137 : XK1 = XI2*XJ3 - XI3*XJ2 138 : XK2 = XI3*XJ1 - XI1*XJ3 139 : XK3 = XI1*XJ2 - XI2*XJ1 140 : C 141 : ZR(1,1) = XJ2*XK3 - XJ3*XK2 142 : ZR(1,2) = XJ3*XK1 - XJ1*XK3 143 : ZR(1,3) = XJ1*XK2 - XJ2*XK1 144 : ZR(2,1) = XI3*XK2 - XI2*XK3 145 : ZR(2,2) = XI1*XK3 - XI3*XK1 146 : ZR(2,3) = XI2*XK1 - XI1*XK2 147 : ZR(3,1) = XI2*XJ3 - XI3*XJ2 148 : ZR(3,2) = XI3*XJ1 - XI1*XJ3 149 : ZR(3,3) = XI1*XJ2 - XI2*XJ1 150 : C 151 : C 152 : C ----------------------------- MATRICE DIMENSIONNEMENT ZD 153 : C -------------------------- 154 : ZD( 1) = XL**2/ZCP 155 : ZD( 2) = XL**2/ZCP 156 : ZD( 3) = XL**2/ZCP 157 : ZD( 4) = XL/ZCR 158 : ZD( 5) = XL**2/ZCP*ZALFAZ 159 : ZD( 6) = XL**2/ZCP*ZALFAY 160 : ZD( 7) = XL**2/ZCP*ZALFAX 161 : ZD( 8) = XL**2/ZCP*(ZALFAY**3) 162 : ZD( 9) = XL**2/ZCP*(ZALFAZ**3) 163 : ZD(10) = XL/ZCR*ZALFAR 164 : ZD(11) = XL**2/ZCP*(ZALFAZ**2) 165 : ZD(12) = XL**2/ZCP*(ZALFAY**2) 166 : C 167 : C ------------------------------- MATRICE DES EFFORTS ZE 168 : C ---------------------- 169 : DO 20 I = 1 , 12 170 : DO 21 J = 1 , 12 171 : ZE(I,J) = 0.D0 172 : 21 CONTINUE 173 : 20 CONTINUE 174 : C 175 : ZEPSY = CRO*CIZ*ZS*ZS/(ZCMU*CKCY*CSE) 176 : ZEPSZ = CRO*CIY*ZS*ZS/(ZCMU*CKCZ*CSE) 177 : C 178 : ZA5 = 1.D0 179 : ZA6 = 1.D0 180 : ZB6 = 1.D0 / (ZCMU*CKCY*CSE) 181 : ZB5 = 1.D0 / (ZCMU*CKCZ*CSE) 182 : C 183 : ZA8 = (CRO*CSE*ZS*ZS+CAM*ZS)/(ZCMU*ZCMU*CKCY*CKCY*CSE*CSE) 184 : * - 1.D0/(ZCE*CIZ) 185 : ZA9 = (CRO*CSE*ZS*ZS+CAM*ZS)/(ZCMU*ZCMU*CKCZ*CKCZ*CSE*CSE) 186 : * - 1.D0/(ZCE*CIY) 187 : ZB8 = CRO*ZS*ZS*(1.D0/(ZCMU*CKCY)+1.D0/ZCE) 188 : * + CAM*ZS/(ZCMU*CKCY*CSE) 189 : ZB9 = CRO*ZS*ZS*(1.D0/(ZCMU*CKCZ)+1.D0/ZCE) 190 : * + CAM*ZS/(ZCMU*CKCZ*CSE) 191 : C 192 : ZA12= 1.D0/(ZCE*CIZ) 193 : ZA11= 1.D0/(ZCE*CIY) 194 : ZB12= (CRO*CSE*ZS*ZS+CAM*ZS)/(ZCMU*CKCY*CSE) 195 : ZB11= (CRO*CSE*ZS*ZS+CAM*ZS)/(ZCMU*CKCZ*CSE) 196 : C 197 : C 198 : ZE( 1, 1) = 1.D0 199 : ZE( 2, 2) = 1.D0 200 : ZE( 3, 3) = 1.D0 201 : ZE( 4, 4) = 1.D0 202 : ZE( 5, 5) = ZA5 203 : ZE( 5, 9) = ZB5 204 : ZE( 6, 6) = ZA6 205 : ZE( 6, 8) = ZB6 206 : ZE( 7, 7) = 1.D0/(ZCE*CSE) 207 : ZE( 8, 6) = ZB8 208 : ZE( 8, 8) = ZA8 209 : ZE( 9, 5) = ZB9 210 : ZE( 9, 9) = ZA9 211 : ZE(10,10) =-1.D0/(ZCMU*CC) 212 : ZE(11, 3) = ZB11 213 : ZE(11,11) = ZA11 214 : ZE(12, 2) = ZB12 215 : ZE(12,12) = ZA12 216 : C 217 : C -- CALCUL DES FONCTIONS DE GREEN 218 : C 218 : CALL ELLP31(0.D0,INP,ZS,ZALFAX,ZALFAR,ZALFAY,ZALFAZ, 219 : * ZG010,ZG110,ZG020,ZG120, 220 : * ZG030,ZG130,ZG230,ZG330,ZG430,ZZ3, 221 : * ZG040,ZG140,ZG240,ZG340,ZG440,ZZ4) 222 : C 222 : CALL ELLP31(XL,INP,ZS,ZALFAX,ZALFAR,ZALFAY,ZALFAZ, 223 : * ZG011,ZG111,ZG021,ZG121, 224 : * ZG031,ZG131,ZG231,ZG331,ZG431,ZZ3, 225 : * ZG041,ZG141,ZG241,ZG341,ZG441,ZZ4) 226 : C 227 : C -- REMPLISSAGE DE LA MATRICE A LOCAL DANS L'ORDRE INITIAL 228 : C 229 : DO 50 J = 1,24 230 : DO 51 I = 1,12 231 : ZA(I,J)= 0. 232 : 51 CONTINUE 233 : 50 CONTINUE 234 : C 235 : ZA( 1, 1) = ZG110+1.D0 236 : ZA( 1, 2) = ZG010 237 : ZA( 1, 3) = ZG111 238 : ZA( 1, 4) =-ZG011 239 : ZA( 2, 1) = ZG111 240 : ZA( 2, 2) = ZG011 241 : ZA( 2, 3) = ZG110+1.D0 242 : ZA( 2, 4) =-ZG010 243 : C 244 : ZA( 3, 5) = ZG120+1.D0 245 : ZA( 3, 6) = ZG020 246 : ZA( 3, 7) = ZG121 247 : ZA( 3, 8) =-ZG021 248 : ZA( 4, 5) = ZG121 249 : ZA( 4, 6) = ZG021 250 : ZA( 4, 7) = ZG120+1.D0 251 : ZA( 4, 8) =-ZG020 252 : C 253 : ZA( 5, 9) =-ZG330+ZZ3*ZG130+1.D0 254 : ZA( 5,10) =-ZG230+ZZ3*ZG030 255 : ZA( 5,11) =-ZG130 256 : ZA( 5,12) =-ZG030 257 : ZA( 5,13) =-ZG331+ZZ3*ZG131 258 : ZA( 5,14) = ZG231-ZZ3*ZG031 259 : ZA( 5,15) =-ZG131 260 : ZA( 5,16) = ZG031 261 : ZA( 6, 9) =-ZG331+ZZ3*ZG131 262 : ZA( 6,10) =-ZG231+ZZ3*ZG031 263 : ZA( 6,11) =-ZG131 264 : ZA( 6,12) =-ZG031 265 : ZA( 6,13) =-ZG330+ZZ3*ZG130+1.D0 266 : ZA( 6,14) = ZG230-ZZ3*ZG030 267 : ZA( 6,15) =-ZG130 268 : ZA( 6,16) = ZG030 269 : ZA( 7, 9) =-ZG430+ZZ3*ZG230 270 : ZA( 7,10) =-ZG330+ZZ3*ZG130+1.D0 271 : ZA( 7,11) =-ZG230 272 : ZA( 7,12) =-ZG130 273 : ZA( 7,13) = ZG431-ZZ3*ZG231 274 : ZA( 7,14) =-ZG331+ZZ3*ZG131 275 : ZA( 7,15) = ZG231 276 : ZA( 7,16) =-ZG131 277 : ZA( 8, 9) =-ZG431+ZZ3*ZG231 278 : ZA( 8,10) =-ZG331+ZZ3*ZG131 279 : ZA( 8,11) =-ZG231 280 : ZA( 8,12) =-ZG131 281 : ZA( 8,13) = ZG430-ZZ3*ZG230 282 : ZA( 8,14) =-ZG330+ZZ3*ZG130+1.D0 283 : ZA( 8,15) = ZG230 284 : ZA( 8,16) =-ZG130 285 : C 286 : ZA( 9,17) =-ZG340+ZZ4*ZG140+1.D0 287 : ZA( 9,18) =-ZG240+ZZ4*ZG040 288 : ZA( 9,19) =-ZG140 289 : ZA( 9,20) =-ZG040 290 : ZA( 9,21) =-ZG341+ZZ4*ZG141 291 : ZA( 9,22) = ZG241-ZZ4*ZG041 292 : ZA( 9,23) =-ZG141 293 : ZA( 9,24) = ZG041 294 : ZA(10,17) =-ZG341+ZZ4*ZG141 295 : ZA(10,18) =-ZG241+ZZ4*ZG041 296 : ZA(10,19) =-ZG141 297 : ZA(10,20) =-ZG041 298 : ZA(10,21) =-ZG340+ZZ4*ZG140+1.D0 299 : ZA(10,22) = ZG240-ZZ4*ZG040 300 : ZA(10,23) =-ZG140 301 : ZA(10,24) = ZG040 302 : ZA(11,17) =-ZG440+ZZ4*ZG240 303 : ZA(11,18) =-ZG340+ZZ4*ZG140+1.D0 304 : ZA(11,19) =-ZG240 305 : ZA(11,20) =-ZG140 306 : ZA(11,21) = ZG441-ZZ4*ZG241 307 : ZA(11,22) =-ZG341+ZZ4*ZG141 308 : ZA(11,23) = ZG241 309 : ZA(11,24) =-ZG141 310 : ZA(12,17) =-ZG441+ZZ4*ZG241 311 : ZA(12,18) =-ZG341+ZZ4*ZG141 312 : ZA(12,19) =-ZG241 313 : ZA(12,20) =-ZG141 314 : ZA(12,21) = ZG440-ZZ4*ZG240 315 : ZA(12,22) =-ZG340+ZZ4*ZG140+1.D0 316 : ZA(12,23) = ZG240 317 : ZA(12,24) =-ZG140 318 : C 319 : C -- TRANSFORMATION DE CETTE MATRICE DANS L'ORDRE U R F M 320 : C 321 : C - ICHANG : MATRICE DE TRANSFORMATION DES COLONNES 322 : C 323 : ICHANG( 1) = 1 324 : ICHANG( 2) = 7 325 : ICHANG( 3) = 13 326 : ICHANG( 4) = 19 327 : ICHANG( 5) = 4 328 : ICHANG( 6) = 10 329 : ICHANG( 7) = 16 330 : ICHANG( 8) = 22 331 : ICHANG( 9) = 2 332 : ICHANG(10) = 6 333 : ICHANG(11) = 12 334 : ICHANG(12) = 8 335 : ICHANG(13) = 14 336 : ICHANG(14) = 18 337 : ICHANG(15) = 24 338 : ICHANG(16) = 20 339 : ICHANG(17) = 3 340 : ICHANG(18) = 5 341 : ICHANG(19) = 11 342 : ICHANG(20) = 9 343 : ICHANG(21) = 15 344 : ICHANG(22) = 17 345 : ICHANG(23) = 23 346 : ICHANG(24) = 21 347 : C 348 : DO 10 I = 1,24 349 : K = ICHANG(I) 350 : DO 11 J = 1,12 351 : ZA1(J,K) = ZA(J,I) 352 : 11 CONTINUE 353 : 10 CONTINUE 354 : C 355 : C 356 : C -- CACUL DE ZCI = ZA1*ZR (DANS CET ORDRE CAR LE PRODUIT N'EST 357 : C -- PAS COMMUTATIF 358 : C 359 : C 360 : C - CALCUL DE ZR : 361 : C 362 : DO 26 I = 1 , 24 363 : DO 25 J = 1 , 24 364 : ZER(I,J)= 0.D0 365 : 25 CONTINUE 366 : 26 CONTINUE 367 : C 368 : DO 30 L = 1,4 369 : DO 31 I = 1,12 370 : DO 32 J = 1,3 371 : DO 33 K = 1,3 372 : ZER(I,(L-1)*3+J) = ZER(I,(L-1)*3+J) + ZE(I,(L-1)*3+K)*ZR(K,J) 373 : 33 CONTINUE 374 : 32 CONTINUE 375 : 31 CONTINUE 376 : 30 CONTINUE 377 : C 378 : C 379 : DO 34 L = 1,2 380 : DO 35 I = 1,12 381 : DO 36 J = 1,3 382 : DO 37 K = 1,3 383 : ZER(I+12,(L-1)*3+J+12) = ZER(I+12,(L-1)*3+J+12) 384 : * + ZE(I,(L-1)*3+K)*ZR(K,J) 385 : 37 CONTINUE 386 : 36 CONTINUE 387 : 35 CONTINUE 388 : 34 CONTINUE 389 : C 390 : C 391 : DO 54 L = 3,4 392 : DO 55 I = 1,12 393 : DO 56 J = 1,3 394 : DO 57 K = 1,3 395 : ZER(I+12,(L-1)*3+J+12) = ZER(I+12,(L-1)*3+J+12) 396 : * - ZE(I,(L-1)*3+K)*ZR(K,J) 397 : 57 CONTINUE 398 : 56 CONTINUE 399 : 55 CONTINUE 400 : 54 CONTINUE 401 : C 402 : C 403 : DO 38 I = 1 , 12 404 : DO 39 J = 1 , 24 405 : ZER(I ,J) = ZER(I ,J)/ZD(I ) 406 : ZER(I+12,J) = ZER(I+12,J)/ZD(I ) 407 : 39 CONTINUE 408 : 38 CONTINUE 409 : C 410 : C --------------------------------- CALCUL DE ZA1*ZER 411 : C ----------------- 412 : DO 60 J = 1,24 413 : DO 61 I = 1,12 414 : ZC1(I,J)= (0.D0,0.D0) 415 : 61 CONTINUE 416 : 60 CONTINUE 417 : C 418 : DO 40 I = 1,12 419 : DO 41 J = 1,24 420 : DO 42 K = 1,24 421 : ZC1(I,J) = ZC1(I,J) + ZA1(I,K)*ZER(K,J) 422 : 42 CONTINUE 423 : 41 CONTINUE 424 : 40 CONTINUE 425 : C 426 : END 427 :
© Cast3M 2003 - All rights reserved.
Disclaimer