1 : C JACONO SOURCE BP208322 13/02/05 21:15:00 7692 2 : SUBROUTINE JACONO(IPMODL,INORM,IPCHE,IRET) 3 : C======================================================================= 4 : C ENTREES : 5 : C --------- 6 : C IPMODL= pointeur sur un MMODEL 7 : C INORM = 1 si les vecteurs doivent etre normes 0 sinon 8 : C 9 : C SORTIES : 10 : C -------- 11 : C 12 : C IPCHE = CHAMELEM contenant les JACOBIENS 13 : C ( = normale aux faces des elements dans le cas des coques) 14 : C ( = tangente a la fibre neutre dans le cas des poutres) 15 : C IRET = 1 si succes 0 sinon 16 : C 17 : C Passage au nouveau Chamelem PAR S.RAMAHANDRY le 11/09/90 18 : C 19 : C 2013-01-02 (BP) : ajout zones cohesives (ZCO2,3 et 4 => coque mince) 20 : C + calcul de la tangente pour les poutres 21 : C 22 : C 23 : C===================================================================== 24 : IMPLICIT INTEGER(I-N) 25 : IMPLICIT REAL*8(A-H,O-Z) 26 : PARAMETER(UN=1.D0,XZER=0.D0) 27 : -INC CCOPTIO 28 : -INC CCHAMP 29 : -INC SMCHAML 30 : -INC SMMODEL 31 : -INC SMELEME 32 : -INC SMCOORD 33 : -INC SMINTE 34 : C 35 : SEGMENT TRA 36 : REAL*8 XEL(3,NBNN) ,SHP(6,NBNN) ,XE(3,NBNN) 37 : ENDSEGMENT 38 : C 39 : SEGMENT TR1 40 : REAL*8 TH(NBN1) ,TXR(3,3,NBN1) ,XJ(3,3) 41 : ENDSEGMENT 42 : C 43 : SEGMENT INFO 44 : INTEGER INFELL(JG) 45 : ENDSEGMENT 46 : C 47 : SEGMENT MPTVAL 48 : INTEGER IPOS(NS) , NSOF(NS) 49 : INTEGER IVAL(NCOSOU) 50 : CHARACTER*16 TYVAL(NCOSOU) 51 : ENDSEGMENT 52 : C 53 : DIMENSION BPSS(3,3) 54 : CHARACTER*8 CMATE 55 : C 56 : REAL*8 XU(3) 57 : REAL*8 XV(3) 58 : REAL*8 XW(3) 59 : C 60 : NHRM=NIFOUR 61 : IRET=1 62 : C 63 : C ACTIVATION DU MODELE 64 : C 65 : MMODEL= IPMODL 66 : SEGACT MMODEL 67 : NSOUS=KMODEL(/1) 68 : C 69 : C CREATION DU MCHELM 70 : C 71 : N1=NSOUS 72 : L1=8 73 : N3=6 74 : SEGINI MCHELM 75 : IPCHE=MCHELM 76 : IF (INORM .EQ. 1) THEN 77 : TITCHE='NORMALES' 78 : ELSE 79 : TITCHE='VECTEURS SURFACE' 80 : ENDIF 81 : IFOCHE=IFOUR 82 : IRT8=1 83 : C____________________________________________________________________ 84 : C 85 : C DEBUT DE LA BOUCLE SUR LES DIFFERENTES ZONES 86 : C____________________________________________________________________ 87 : C 88 : DO 500 ISOUS=1,NSOUS 89 : C 90 : C ON RECUPERE L INFORMATION GENERALE 91 : C 92 : IMODEL=KMODEL(ISOUS) 93 : SEGACT IMODEL 94 : IPMAIL=IMAMOD 95 : IMACHE(ISOUS)=IPMAIL 96 : CONCHE(ISOUS)=CONMOD 97 : C 98 : C TRAITEMENT DU MODELE 99 : C 100 : MELE=NEFMOD 101 : MELEME=IMAMOD 102 : NFOR=FORMOD(/2) 103 : NMAT=MATMOD(/2) 104 : C____________________________________________________________________ 105 : C 106 : C INFORMATION SUR L'ELEMENT FINI 107 : C____________________________________________________________________ 108 : C 109 : if(infmod(/1).lt.7) then 109 : CALL ELQUOI(MELE,0,5,IPINF,IMODEL) 110 : IF (IERR.NE.0) THEN 111 : SEGDES IMODEL,MMODEL 112 : SEGSUP MCHELM 113 : IRET=0 114 : RETURN 115 : ENDIF 116 : INFO=IPINF 117 : MELE =INFELL(1) 118 : MFR =INFELL(13) 119 : MINTE=INFELL(11) 120 : MINTE1=INFELL(12) 121 : segsup info 122 : else 123 : MELE =INFELE(1) 124 : MFR =INFELE(13) 125 : MINTE=INFMOD(7) 126 : MINTE1=INFMOD(8) 127 : endif 128 : C 129 : INFCHE(ISOUS,1)=0 130 : INFCHE(ISOUS,2)=0 131 : INFCHE(ISOUS,3)=NHRM 132 : INFCHE(ISOUS,4)=MINTE 133 : INFCHE(ISOUS,5)=0 134 : INFCHE(ISOUS,6)=5 135 : C 136 : C INITIALISATION DE MINTE 137 : C 138 : SEGACT MINTE 139 : NBPGAU=POIGAU(/1) 140 : C 141 : C ACTIVATION DU MELEME 142 : C 143 : SEGACT MELEME 144 : NBNN =NUM(/1) 145 : NBELEM=NUM(/2) 146 : C 147 : C RECHERCHE DE LA TAILLE DES MELVAL A ALLOUER 148 : C 149 : IF (MFR.EQ.3.OR.MFR.EQ.5.OR.MFR.EQ.9 .OR. MFR.EQ.77) THEN 150 : N1PTEL=NBPGAU 151 : N1EL=NBELEM 152 : ELSEIF(MFR.EQ.7 .OR. MFR.EQ.13) THEN 153 : N1PTEL=NBPGAU 154 : N1EL=NBELEM 155 : ELSE 156 : N1PTEL = 0 157 : N1EL = 0 158 : ENDIF 159 : C 160 : C CREATION DU MCHAML DE LA SOUS ZONE 161 : C 162 : NJAC=IDIM 163 : N2 = IDIM 164 : SEGINI MCHAML 165 : ICHAML(ISOUS)=MCHAML 166 : NS=1 167 : NCOSOU=NJAC 168 : SEGINI MPTVAL 169 : IVAJAC=MPTVAL 170 : C 171 : C 2 OU 3 COMPOSANTES 172 : C 173 : ICOMP=1 174 : IF (IFOMOD .EQ. 0) THEN 175 : NOMCHE(ICOMP)='VR ' 176 : ELSE 177 : NOMCHE(ICOMP)='VX ' 178 : ENDIF 179 : TYPCHE(ICOMP)='REAL*8' 180 : N2PTEL=0 181 : N2EL=0 182 : SEGINI MELVA1 183 : IELVAL(ICOMP)=MELVA1 184 : IVAL(ICOMP)=MELVA1 185 : C 186 : ICOMP=2 187 : IF (IFOMOD .EQ. 0) THEN 188 : NOMCHE(ICOMP)='VZ ' 189 : ELSE 190 : NOMCHE(ICOMP)='VY ' 191 : ENDIF 192 : TYPCHE(ICOMP)='REAL*8' 193 : N2PTEL=0 194 : N2EL=0 195 : SEGINI MELVA2 196 : IELVAL(ICOMP)=MELVA2 197 : IVAL(ICOMP)=MELVA2 198 : C 199 : IF (IDIM .EQ. 3) THEN 200 : ICOMP=3 201 : NOMCHE(ICOMP)='VZ ' 202 : TYPCHE(ICOMP)='REAL*8' 203 : N2PTEL=0 204 : N2EL=0 205 : SEGINI MELVA3 206 : IELVAL(ICOMP)=MELVA3 207 : IVAL(ICOMP)=MELVA3 208 : ENDIF 209 : C 210 : 44 CONTINUE 211 : C 212 : SEGINI TRA 213 : C 214 : C ================ FORMULATION MASSIVE ======================= 215 : C 216 : C 217 : IF(MFR.EQ.1.OR.MFR.EQ.33) THEN 218 : GOTO 520 219 : C 220 : C ================ FORMULATION COQUE MINCE ===================== 221 : C 222 : C 223 : ELSE IF(MFR.EQ.3.OR.MFR.EQ.9 .OR. MFR.EQ.77) THEN 224 : IDI2=IDIM-1 225 : DO 3000 IB=1,NBELEM 226 : *--------------Calcul de la normale a l'élément 227 : IREF1 = (IDIM + 1)*(MELEME.NUM(1, IB) - 1) 228 : IREF2 = (IDIM + 1)*(MELEME.NUM(2, IB) - 1) 229 : DO 28 IC = 1, IDIM 230 : 28 XU(IC) = XCOOR(IREF2+IC)-XCOOR(IREF1+IC) 231 : XNORU = 0. 232 : DO 31 IC = 1, IDIM 233 : 31 XNORU = XNORU + XU(IC)*XU(IC) 234 : XNORU = SQRT(XNORU) 235 : DO 32 IC = 1, IDIM 236 : 32 XU(IC) = XU(IC)/XNORU 237 : IF (IDIM .EQ. 2) THEN 238 : XW(1) = -XU(2) 239 : XW(2) = XU(1) 240 : ELSE 241 : IN = 3 242 : 33 IREF3 = (IDIM + 1)*(MELEME.NUM(IN, IB) - 1) 243 : DO 34 IC = 1, IDIM 244 : 34 XV(IC) = XCOOR(IREF3+IC)-XCOOR(IREF1+IC) 245 : XNORV = 0. 246 : DO 35 IC = 1, IDIM 247 : 35 XNORV = XNORV + XV(IC)*XV(IC) 248 : XNORV = SQRT(XNORV) 249 : DO 36 IC = 1, IDIM 250 : 36 XV(IC) = XV(IC)/XNORV 251 : XW(1) = XU(2)*XV(3) - XU(3)*XV(2) 252 : XW(2) = XU(3)*XV(1) - XU(1)*XV(3) 253 : XW(3) = XU(1)*XV(2) - XU(2)*XV(1) 254 : XNORW = 0. 255 : DO 37 IC = 1, IDIM 256 : 37 XNORW = XNORW + XW(IC)*XW(IC) 257 : IF (XNORW .LT. 1.E-4) THEN 258 : IN = IN + 1 259 : if(IN.le.NBNN) GOTO 33 260 : ENDIF 261 : XNORW = SQRT(XNORW) 261 : IF(XNORW .LT.1.E-4) call erreur(345) 262 : DO 38 IC = 1, IDIM 263 : 38 XW(IC) = XW(IC)/XNORW 264 : ENDIF 265 : *--------------Fin du calcul de la normale a l'élément 265 : CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE) 266 : C 267 : IF(IDIM.EQ.2)THEN 267 : CALL VPAST2(XE,BPSS) 268 : ELSE IF(IDIM.EQ.3) THEN 268 : CALL VPAST(XE,BPSS) 269 : ENDIF 269 : CALL VCORL1(XE,XEL,BPSS,NBNN) 270 : DO 3002 IC=1,NBPGAU 271 : IF (INORM .EQ. 1) THEN 272 : DJAC = 1. 273 : ELSE 274 : DO 3003 ID=1,6 275 : DO 3003 IE=1,NBNN 276 : SHP(ID,IE)=SHPTOT(ID,IE,IC) 277 : 3003 CONTINUE 277 : CALL JACOBI(XEL,SHP,IDI2,NBNN,DJAC) 278 : ENDIF 279 : MPTVAL=IVAJAC 280 : MELVAL = IVAL(1) 281 : IBMN=MIN(IB, MELVAL.VELCHE(/2)) 282 : IGMN=MIN(IC,MELVAL.VELCHE(/1)) 283 : MELVAL.VELCHE(IGMN,IBMN)=DJAC*XW(1) 284 : MELVAL = IVAL(2) 285 : IBMN=MIN(IB, MELVAL.VELCHE(/2)) 286 : IGMN=MIN(IC,MELVAL.VELCHE(/1)) 287 : MELVAL.VELCHE(IGMN,IBMN)=DJAC*XW(2) 288 : IF (IDIM .EQ. 3) THEN 289 : MELVAL = IVAL(3) 290 : IBMN=MIN(IB, MELVAL.VELCHE(/2)) 291 : IGMN=MIN(IC,MELVAL.VELCHE(/1)) 292 : MELVAL.VELCHE(IGMN,IBMN)=DJAC*XW(3) 293 : ENDIF 294 : 3002 CONTINUE 295 : 3000 CONTINUE 296 : GOTO 520 297 : C 298 : C ================ FORMULATION POUTRE ET TUYAU ==================== 299 : C 300 : C 301 : ELSE IF(MFR.EQ.7.OR.MFR.EQ.13) THEN 302 : IDI2=IDIM-1 303 : DO 4000 IB=1,NBELEM 304 : *-----------Calcul de la tangente a l'élément 305 : IREF1 = (IDIM + 1)*(MELEME.NUM(1, IB) - 1) 306 : IREF2 = (IDIM + 1)*(MELEME.NUM(2, IB) - 1) 307 : DO 41 IC = 1, IDIM 308 : 41 XU(IC) = XCOOR(IREF2+IC)-XCOOR(IREF1+IC) 309 : XNORU = 0.D0 310 : DO 42 IC = 1, IDIM 311 : 42 XNORU = XNORU + XU(IC)*XU(IC) 312 : XNORU = SQRT(XNORU) 313 : DO 43 IC = 1, IDIM 314 : 43 XU(IC) = XU(IC)/XNORU 315 : *-----------Fin du calcul de la tangente a l'élément 316 : *-----------Calcul de la tangente a l'élément en chaque point de Gauss 317 : c BP : On suppose le jacobien constant dans l'element (idem POUJAC.eso) 318 : C => on sort le calcul du jacobien de la boucle sur les points de Gauss. 319 : C Cela implique que la POUTre de Bernoulli n'est pas isoparamétrique... 320 : IF (INORM .EQ. 1) THEN 321 : DJAC = 1.D0 322 : ELSE 323 : DJAC = 1.D0/DBLE(NBPGAU) 324 : ENDIF 325 : DO 4002 IC=1,NBPGAU 326 : MPTVAL=IVAJAC 327 : MELVAL = IVAL(1) 328 : IBMN=MIN(IB, MELVAL.VELCHE(/2)) 329 : IGMN=MIN(IC,MELVAL.VELCHE(/1)) 330 : MELVAL.VELCHE(IGMN,IBMN)=DJAC*XU(1) 331 : MELVAL = IVAL(2) 332 : IBMN=MIN(IB, MELVAL.VELCHE(/2)) 333 : IGMN=MIN(IC,MELVAL.VELCHE(/1)) 334 : MELVAL.VELCHE(IGMN,IBMN)=DJAC*XU(2) 335 : IF (IDIM .EQ. 3) THEN 336 : MELVAL = IVAL(3) 337 : IBMN=MIN(IB, MELVAL.VELCHE(/2)) 338 : IGMN=MIN(IC,MELVAL.VELCHE(/1)) 339 : MELVAL.VELCHE(IGMN,IBMN)=DJAC*XU(3) 340 : ENDIF 341 : 4002 CONTINUE 342 : 4000 CONTINUE 343 : 344 : GOTO 520 345 : C 346 : C ================ FORMULATION COQUE EPAISSE ==================== 347 : C 348 : C 349 : ELSE IF(MFR.EQ.5) THEN 350 : SEGACT MINTE1 351 : NBPGA1=MINTE1.POIGAU(/1) 352 : NBN1 =MINTE1.SHPTOT(/2) 353 : SEGINI TR1 354 : C 355 : C UNE PETITE HORREUR ON CONSIDERE LES EPAISSEURS CONSTANTES 356 : C 357 : DO 5010 IC=1,NBNN 358 : TH(IC)=UN 359 : 5010 CONTINUE 360 : DO 5000 IB=1,NBELEM 361 : *--------------Calcul de la normale a l'élément 362 : IREF1 = (IDIM + 1)*(MELEME.NUM(1, IB) - 1) 363 : * IREF2 = (IDIM + 1)*(MELEME.NUM(2, IB) - 1) 364 : * bp : les EF de coque epaisse etant quadratiques (coq6 et coq8), on 365 : * prend les noeuds "coins" pour eviter pb avec les noeuds 1,2,3 si courbures 366 : IREF2 = (IDIM + 1)*(MELEME.NUM(3, IB) - 1) 367 : DO 8 IC = 1, IDIM 368 : 8 XU(IC) = XCOOR(IREF2+IC)-XCOOR(IREF1+IC) 369 : XNORU = 0. 370 : DO 11 IC = 1, IDIM 371 : 11 XNORU = XNORU + XU(IC)*XU(IC) 372 : XNORU = SQRT(XNORU) 373 : DO 12 IC = 1, IDIM 374 : 12 XU(IC) = XU(IC)/XNORU 375 : IF (IDIM .EQ. 2) THEN 376 : XW(1) = -XU(2) 377 : XW(2) = XU(1) 378 : ELSE 379 : * IN = 3 380 : IN = 5 381 : 13 IREF3 = (IDIM + 1)*(MELEME.NUM(IN, IB) - 1) 382 : DO 14 IC = 1, IDIM 383 : 14 XV(IC) = XCOOR(IREF3 + IC) - XCOOR(IREF1 + IC) 384 : XNORV = 0. 385 : DO 15 IC = 1, IDIM 386 : 15 XNORV = XNORV + XV(IC)*XV(IC) 387 : XNORV = SQRT(XNORV) 388 : DO 16 IC = 1, IDIM 389 : 16 XV(IC) = XV(IC)/XNORV 390 : XW(1) = XU(2)*XV(3) - XU(3)*XV(2) 391 : XW(2) = XU(3)*XV(1) - XU(1)*XV(3) 392 : XW(3) = XU(1)*XV(2) - XU(2)*XV(1) 393 : XNORW = 0. 394 : DO 17 IC = 1, IDIM 395 : 17 XNORW = XNORW + XW(IC)*XW(IC) 396 : XNORW = SQRT(XNORW) 397 : IF (XNORW .LT. 1.E-4) THEN 398 : if(IN.LT.NBNN) then 399 : IN = IN + 1 400 : GOTO 13 401 : else 402 : write(6,*) 'Difficultes pour etablir la normale de' 403 : write(6,*) 'l element',IB 404 : write(6,*) 'Verifiez votre maillage' 405 : GOTO 9990 406 : endif 407 : ENDIF 408 : DO 18 IC = 1, IDIM 409 : 18 XW(IC) = XW(IC)/XNORW 410 : ENDIF 411 : *--------------Fin du calcul de la normale a l'élément 412 : IF (INORM .EQ. 0) THEN 412 : CALL DOXE(XCOOR,IDIM,NBNN,NUM,IB,XE) 413 : C 413 : CALL CQ8LOC(XE,NBN1,MINTE1.SHPTOT,TXR,IRR) 414 : ENDIF 415 : C 416 : DO 5002 IC=1,NBPGAU 417 : IF (INORM .EQ. 1) THEN 418 : DJAC = 1. 419 : ELSE 420 : E=DZEGAU(IC) 420 : CALL COQ8JC(IC,NBN1,E,XE,TH,TXR,SHPTOT,XJ,DJAC,IRR) 421 : ENDIF 422 : MPTVAL=IVAJAC 423 : MELVAL = IVAL(1) 424 : IBMN=MIN(IB, MELVAL.VELCHE(/2)) 425 : IGMN=MIN(IC, MELVAL.VELCHE(/1)) 426 : MELVAL.VELCHE(IGMN,IBMN)=DJAC*XW(1) 427 : MELVAL = IVAL(2) 428 : IBMN=MIN(IB, MELVAL.VELCHE(/2)) 429 : IGMN=MIN(IC, MELVAL.VELCHE(/1)) 430 : MELVAL.VELCHE(IGMN,IBMN)=DJAC*XW(2) 431 : IF (IDIM .EQ. 3) THEN 432 : MELVAL = IVAL(3) 433 : IBMN=MIN(IB, MELVAL.VELCHE(/2)) 434 : IGMN=MIN(IC, MELVAL.VELCHE(/1)) 435 : MELVAL.VELCHE(IGMN,IBMN)=DJAC*XW(3) 436 : ENDIF 437 : 5002 CONTINUE 438 : 5000 CONTINUE 439 : GOTO 520 440 : ENDIF 441 : C--------------------------------------------------------------------- 442 : C DESACTIVATION DES SEGMENTS PROPRES A LA ZONE GEOMETRIQUE ISOUS 443 : C--------------------------------------------------------------------- 444 : C 445 : 520 CONTINUE 446 : MPTVAL=IVAJAC 447 : DO 515 IO=1,NJAC 448 : IF(IVAL(IO).NE.0) THEN 449 : MELVAL=IVAL(IO) 450 : SEGDES MELVAL 451 : ENDIF 452 : 515 CONTINUE 453 : SEGSUP MPTVAL 454 : C 455 : SEGDES MINTE 456 : * SEGSUP INFO 457 : C 458 : SEGDES IMODEL 459 : SEGDES MELEME 460 : SEGDES MCHAML 461 : C 462 : IF (MFR.EQ.5) THEN 463 : SEGDES MINTE1,TR1 464 : ENDIF 465 : C 466 : SEGSUP TRA 467 : C 468 : 500 CONTINUE 469 : SEGDES MMODEL 470 : SEGDES MCHELM 471 : RETURN 472 : C 473 : C 474 : 9990 CONTINUE 475 : * 476 : C------------------------------------------------------------------- 477 : C ERREUR DANS UNE ZONE , DESACTIVATION ET RETOUR 478 : C------------------------------------------------------------------- 479 : * 480 : * 481 : * 482 : MPTVAL=IVAJAC 483 : DO 9993 IO=1,NJAC 484 : IF (IVAL(IO).NE.0) THEN 485 : MELVAL=IVAL(IO) 486 : SEGDES MELVAL 487 : ENDIF 488 : 9993 CONTINUE 489 : SEGSUP MPTVAL 490 : * 491 : * 492 : * 493 : SEGDES MELEME 494 : SEGDES IMODEL 495 : SEGDES MCHAML 496 : C 497 : SEGDES MMODEL 498 : SEGDES MCHELM 499 : SEGDES MINTE 500 : * SEGSUP INFO 501 : RETURN 502 : END 503 : 504 : 505 : 506 : 507 : 508 : 509 : 510 : 511 : 512 : 513 : 514 : 515 : 516 :
© Cast3M 2003 - All rights reserved.
Disclaimer