1 : C TESTMA SOURCE BP208322 12/10/16 21:15:00 7530 2 : SUBROUTINE TESTMA(IPCHE,IMEL,COND,CONST,IRET,IMODI) 3 : *______________________________________________________________________ 4 : * 5 : * test sur les maillages 6 : * 7 : * entrees : 8 : * --------- 9 : * ipche chamelem dont il faut vérifier le maillage (type mchaml) 10 : * imel maillage du modèle: sert de référence (type meleme) 11 : * Maillage actif en Entree (?) et en Sortie (oui) 12 : * cond vrai si il faut tenir compte de const 13 : * const nom de constituant 14 : * 15 : * 16 : * sortie : 17 : * -------- 18 : * iret chamelem réordonné si nécessaire 19 : * = 0 si pb 20 : * imodi =1 il a fallu decouper le maillage imel 21 : * 22 : *______________________________________________________________________ 23 : * 24 : * declarations 25 : * 26 : IMPLICIT REAL*8(A-H,O-Z) 27 : IMPLICIT INTEGER(I-N) 28 : * 29 : -INC CCOPTIO 30 : * 31 : -INC SMCHAML 32 : -INC SMELEME 33 : -INC SMCOORD 34 : -INC SMINTE 35 : * 36 : SEGMENT ICPR(NPT) 37 : SEGMENT ICPEL(NPT) 38 : SEGMENT ITRAV(NNNT) 39 : SEGMENT NIMEL(NNNT) 40 : SEGMENT IPERM(NNNT) 41 : SEGMENT NZONE(NSOUS) 42 : SEGMENT ICOM(NBEL) 43 : * 44 : SEGMENT ITRA 45 : INTEGER ICC (IA),IRE(IA,IMA) 46 : ENDSEGMENT 47 : * 48 : SEGMENT WTRAV 49 : REAL*8 QSIREF(NBPGAU),QSICHAM(NBPGAU),ETAREF(NBPGAU) 50 : REAL*8 ETACHAM(NBPGAU),DZEREF(NBPGAU),DZECHAM(NBPGAU) 51 : REAL*8 XECHAM(3,NP1),XEREF(3,NP1) 52 : INTEGER TABOK(NBPGAU),TAB(NBPGAU) 53 : ENDSEGMENT 54 : LOGICAL INIT,COND 55 : CHARACTER*(*) CONST 56 : * 57 : * executable 58 : * 59 : NPT = XCOOR(/1)/(IDIM+1) 60 : ICPR = 0 61 : ICPEL = 0 62 : * 63 : IRET = 0 64 : IMODI = 0 65 : * 66 : IPT1 = IMEL 67 : SEGACT IPT1 68 : NBSOU1=IPT1.LISOUS(/1) 69 : * La routine ne travaille qu'avec un maillage elementaire 70 : IF ( NBSOU1 .GT. 1) THEN 71 : SEGDES,IPT1 72 : RETURN 73 : ENDIF 74 : MELEME = IPT1 75 : * 76 : * champ servant à stocker provisoirement le resultat 77 : * 78 : MCHELM = IPCHE 79 : SEGACT MCHELM 80 : SEGINI,MCHEL1=MCHELM 81 : * 82 : * segment servant a stocker l'ordre de stockage dans mchel1 83 : NS01=0 84 : NSOUS=ICHAML(/1) 85 : SEGINI NZONE 86 : * 87 : * creation de wtrav : on l'ajustera ensuite si necessaire 88 : * mais on sort la creation destruction de la boucle 89 : * la definition doit etre coherente avec celle dans tabgau 90 : nbpgau=0 91 : np1=0 92 : wtrav=0 93 : segini wtrav 94 : * 95 : * boucle conditionnelle sur les sous zones du maillage imel 96 : * 97 : IBOU1 = 0 98 : * 99 : 1 CONTINUE 100 : IBOU1 = IBOU1 + 1 101 : * 102 : * Initialisations de ICPR, ICPEL la premiere fois puis 103 : * remises a zero pour les fois suivantes 104 : IF (ICPR.EQ.0) THEN 105 : SEGINI ICPR,ICPEL 106 : ELSE 107 : DO i = 1, NPT 108 : ICPR(i) = 0 109 : ICPEL(i) = 0 110 : ENDDO 111 : ENDIF 112 : * 113 : IF (NBSOU1.NE.0) THEN 114 : MELEME = IPT1.LISOUS(IBOU1) 115 : ENDIF 116 : SEGACT MELEME 117 : * 118 : * nnnt :nb d'elements de imel pour la sous zone consideree 119 : * 120 : NNNT=NUM(/2) 121 : * 122 : IA=0 123 : IMA=0 124 : * 125 : * fabrication de icpel qui dira combien de fois un point apparait 126 : * de icpr qui donne une numerotation locale 127 : * 128 : DO 3 J=1,NNNT 129 : DO 3 K=1,NUM(/1) 130 : ID=NUM(K,J) 131 : ICPEL(ID)=ICPEL(ID)+1 132 : * on retient le max des occurences des points 133 : IMA=MAX(ICPEL(ID),IMA) 134 : IF(ICPR(ID).NE.0) GO TO 3 135 : * numerotation locale 136 : IA=IA+1 137 : ICPR(ID)=IA 138 : 3 CONTINUE 139 : * 140 : * fabrication de itra : 141 : * icc donne le nombre d'elements a partir du numero local 142 : * ire donne les numeros d'elements de imel 143 : * a partir du numero local et de l'occurence 144 : * 145 : SEGINI ITRA 146 : DO 8 J=1,NNNT 147 : DO 8 K=1,NUM(/1) 148 : ID=ICPR(NUM(K,J)) 149 : ICC(ID)=ICC(ID)+1 150 : IA=ICC(ID) 151 : IRE(ID,IA)=J 152 : 8 CONTINUE 153 : * 154 : * il faut maintenant regarder si dans un imache il existe 155 : * les elements de la sous zone du meleme imel 156 : * boucle sur les sous zone du chamelem 157 : * 158 : IFLAG=0 159 : * initialisons ipt2 160 : IPT2=-1 161 : DO 10 I=1,ICHAML(/1) 162 : IF (COND) THEN 163 : IF (CONCHE(I) .NE. CONST) GOTO 10 164 : ENDIF 165 : MCHAML=ICHAML(I) 166 : IPT2 =IMACHE(I) 167 : IF (IPT2 .EQ. MELEME) THEN 168 : IFLAG=1 169 : NS01=NS01+1 170 : NZONE(NS01)=I 171 : MCHEL1.ICHAML(NS01)=MCHEL1.ICHAML(I) 172 : MCHEL1.IMACHE(NS01)=MELEME 173 : MCHEL1.CONCHE(NS01)=MCHEL1.CONCHE(I) 174 : GOTO 10 175 : ENDIF 176 : SEGACT,MELEME,IPT2 177 : IF(ITYPEL.NE.IPT2.ITYPEL) GO TO 11 178 : * 179 : NBEL = NUM(/2) 180 : SEGINI ICOM 181 : * 182 : INIT=.FALSE. 183 : * 184 : ICO=0 185 : NP1=IPT2.NUM(/1) 186 : DO 12 K=1,IPT2.NUM(/2) 187 : DO 13 L=1,NP1 188 : ID=IPT2.NUM(L,K) 189 : IDD = ICPR(ID) 190 : IF(IDD.EQ.0) GO TO 12 191 : 13 CONTINUE 192 : * 193 : * ok l'element k possede ses noeuds dans imel 194 : * 195 : IF(ITYPEL.EQ.1) THEN 196 : ID=IPT2.NUM(1,K) 197 : IDD=ICPR(ID) 198 : IRE1=IRE(IDD,1) 199 : GOTO 20 200 : ENDIF 201 : * 202 : * ces noeuds correspondent-t-ils a un meme element ire1 203 : * 204 : ID = IPT2.NUM(1,K) 205 : IDD=ICPR(ID) 206 : IDE=ICC(IDD) 207 : c boucle sur les occurences de idd dans les elements 208 : DO 14 L=1,IDE 209 : IRE1=IRE(IDD,L) 210 : DO 15 M=2,NP1 211 : IDD2=ICPR(IPT2.NUM(M,K)) 212 : IF(IDD2.EQ.0) GO TO 12 213 : IDE2=ICC(IDD2) 214 : DO 16 N=1,IDE2 215 : IF(IRE(IDD2,N).EQ.IRE1) GO TO 15 216 : 16 CONTINUE 217 : GO TO 14 218 : 15 CONTINUE 219 : GOTO 20 220 : 14 CONTINUE 221 : * 222 : * l'element n'appartient pas a meleme 223 : GOTO 12 224 : * 225 : * itrav donne le numero de l'element dans le chamelem 226 : * nimel donne le numero de l'element dans la sous zone 227 : * de imel 228 : cbp IPERM = 0 si pas de permutation des noeuds 229 : * 230 : 20 CONTINUE 231 : c 232 : IF(.NOT.INIT) THEN 233 : INIT=.TRUE. 234 : SEGINI ITRAV 235 : SEGINI NIMEL 236 : SEGINI IPERM 237 : ENDIF 238 : ICO=ICO+1 239 : ITRAV(ICO)=K 240 : NIMEL(ICO)=IRE1 241 : ICOM(IRE1)=1 242 : IPERM(ICO)=0 243 : do iii=1,NP1 244 : if(IPT2.NUM(iii,K).ne.NUM(iii,IRE1)) then 245 : IPERM(ICO)=iii 246 : goto 12 247 : endif 248 : enddo 249 : 12 CONTINUE 250 : * 251 : IF(INIT) THEN 252 : * 253 : * le nb d'elements est-il egal pour la sous zone imel et 254 : * la sous zone mchaml 255 : * 256 : IF (ICO.EQ.NNNT) THEN 257 : IFLAG=1 258 : NS01=NS01+1 259 : NZONE(NS01)=I 260 : MINTE=INFCHE(I,4) 261 : IPMINT=MINTE 262 : IF (IPMINT.NE.0) SEGACT,MINTE 263 : MCHAML=ICHAML(I) 264 : SEGINI,MCHAM1=MCHAML 265 : MCHEL1.ICHAML(NS01)=MCHAM1 266 : MCHEL1.IMACHE(NS01)=MELEME 267 : MCHEL1.CONCHE(NS01)=MCHEL1.CONCHE(I) 268 : * 269 : DO 24 ICOMP=1,MCHAM1.IELVAL(/1) 270 : MELVAL=MCHAM1.IELVAL(ICOMP) 271 : SEGACT MELVAL 272 : IMEM=0 273 : IF (MCHAM1.TYPCHE(ICOMP).EQ.'REAL*8') THEN 274 : N1PTEL=VELCHE(/1) 275 : N1EL=ICO 276 : IF (VELCHE(/2).EQ.1) THEN 277 : N1EL=1 278 : IMEM=1 279 : ENDIF 280 : NEL=N1EL 281 : N2PTEL=0 282 : N2EL =0 283 : ELSE 284 : N2PTEL=IELCHE(/1) 285 : N2EL=ICO 286 : IF (IELCHE(/2).EQ.1) THEN 287 : N2EL=1 288 : IMEM=1 289 : ENDIF 290 : NEL=N2EL 291 : N1PTEL=0 292 : N1EL =0 293 : ENDIF 294 : * 295 : SEGINI MELVA1 296 : * 297 : DO 21 J=1,NEL 298 : K =ITRAV(J) 299 : IRE1=NIMEL(J) 300 : KTEST=K 301 : ITEST=IRE1 302 : IF (NEL.EQ.1) THEN 303 : IRE1=1 304 : IF (IMEM.EQ.1) THEN 305 : K=1 306 : ENDIF 307 : ENDIF 308 : * 309 : * dans le cas des champs qui ne sont pas définis aux noeuds 310 : * vérification de la bonne orientation des maillages 311 : * sinon création d'un tableau de correspondance 312 : * 313 : NBPGAU=MAX(N1PTEL,N2PTEL) 314 : cbp -cas ou l'on doit identifier la position des pts d'integration 315 : cbp (ne marche pas avec des coques epaisses MFR=5, 316 : cbp ni avec DKT dont nombre de points >1 dans l'épaisseur) 317 : IF((IPMINT.NE.0).and.(IPERM(J).ne.0)) THEN 318 : 319 : if (qsiref(/1).ne.nbpgau .or. 320 : & xeref(/2).ne.np1) then 321 : segadj,wtrav 322 : endif 322 : CALL TABGAU(IPMINT,ITEST,KTEST,NP1,MELEME, 323 : + IPT2,NBPGAU,IRET1,wtrav) 324 : IF (IRET1.EQ.0) THEN 324 : CALL ERREUR(757) 325 : SEGSUP MCHEL1,ITRA 326 : if (imel.ne.ipt2) SEGDES IPT2 327 : if (imel.ne.meleme) segdes meleme 328 : GOTO 900 329 : ENDIF 330 : IF (N1PTEL.EQ.0) THEN 331 : DO 22 IGAU=1,N2PTEL 332 : IGAU1=TAB(IGAU) 333 : MELVA1.IELCHE(IGAU,IRE1)=IELCHE(IGAU1,K) 334 : 22 CONTINUE 335 : ELSE 336 : DO 23 IGAU=1,N1PTEL 337 : IGAU1=TAB(IGAU) 338 : MELVA1.VELCHE(IGAU,IRE1)=VELCHE(IGAU1,K) 339 : 23 CONTINUE 340 : iimpi=0 341 : ENDIF 342 : 343 : cbp -cas ou on recopie directement sans chercher la permutation 344 : ELSE 345 : IF (N1PTEL.EQ.0) THEN 346 : DO 25 IGAU=1,N2PTEL 347 : MELVA1.IELCHE(IGAU,IRE1)=IELCHE(IGAU,K) 348 : 25 CONTINUE 349 : ELSE 350 : DO 26 IGAU=1,N1PTEL 351 : MELVA1.VELCHE(IGAU,IRE1)=VELCHE(IGAU,K) 352 : 26 CONTINUE 353 : iimpi=0 354 : ENDIF 355 : 356 : ENDIF 357 : * 358 : 21 CONTINUE 359 : MCHAM1.IELVAL(ICOMP)=MELVA1 360 : SEGDES MELVA1,MELVAL 361 : 24 CONTINUE 362 : SEGDES MCHAM1,MCHAML 363 : SEGSUP ITRAV,NIMEL,IPERM 364 : ELSE 365 : * 366 : * le maillage est a cheval sur plusieurs sous zones 367 : * du mchaml 368 : * il va falloir scinder le maillage en deux 369 : * et recommencer le processus 370 : * 371 : IMODI = 1 372 : * creation des deux maillages dont l'union est le maillage meleme 373 : NBNN = NUM(/1) 374 : NBELEM = ICO 375 : NBREF = 0 376 : NBSOUS = 0 377 : SEGINI IPT3 378 : NBELEM = NUM(/2) - ICO 379 : IF (NBELEM.LE.0) THEN 379 : CALL ERREUR(26) 380 : GOTO 900 381 : ENDIF 382 : SEGINI IPT4 383 : IPT3.ITYPEL = ITYPEL 384 : IPT4.ITYPEL = ITYPEL 385 : I3 = 0 386 : I4 = 0 387 : DO 45 II=1,NUM(/2) 388 : IF (ICOM(II) .EQ. 1) THEN 389 : * l'element est commun aux deux maillages 390 : I3=I3+1 391 : IPT3.ICOLOR(I3)=ICOLOR(II) 392 : DO 43 JJ=1,NUM(/1) 393 : IPT3.NUM(JJ,I3)=NUM(JJ,II) 394 : 43 CONTINUE 395 : ELSE 396 : * l'element appartient seulement a meleme 397 : I4=I4+1 398 : IPT4.ICOLOR(I4)=ICOLOR(II) 399 : DO 44 JJ=1,NUM(/1) 400 : IPT4.NUM(JJ,I4)=NUM(JJ,II) 401 : 44 CONTINUE 402 : ENDIF 403 : 45 CONTINUE 404 : * 405 : * modification de ipt1 406 : * 407 : IF (IPT1.LISOUS(/1) .LE. 1) THEN 408 : NBSOUS=2 409 : NBREF=0 410 : NBELEM=0 411 : NBNN=0 412 : SEGINI IPT5 413 : IPT5.LISOUS(1)=IPT3 414 : IPT5.LISOUS(2)=IPT4 415 : if (imel.ne.ipt1) SEGDES IPT1 416 : IPT1=IPT5 417 : ELSE 418 : NBSOUS=IPT1.LISOUS(/1)+1 419 : NBREF=0 420 : NBELEM=0 421 : NBNN=0 422 : SEGADJ IPT1 423 : IPT1.LISOUS(IBOU1)=IPT3 424 : IPT1.LISOUS(NBSOUS)=IPT4 425 : if (imel.ne.meleme) SEGDES MELEME 426 : ENDIF 427 : IBOU1=IBOU1-1 428 : NBSOU1=NBSOUS 429 : SEGSUP ITRAV,NIMEL,IPERM 430 : GOTO 1 431 : ENDIF 432 : ENDIF 433 : SEGSUP ICOM 434 : 11 CONTINUE 435 : if (imel.ne.ipt2) SEGDES IPT2 436 : 10 CONTINUE 437 : * 438 : * on n'a pas trouve de sous zone dans le mchaml qui 439 : * correspondent a une sous zone de imel 440 : * 441 : IF (IFLAG.EQ.0) THEN 442 : if (imel.ne.meleme) segdes meleme 443 : SEGSUP MCHEL1 444 : IF(IPT2.LE.0) THEN 445 : SEGSUP,ITRA 446 : ELSE 447 : if (imel.ne.ipt2) SEGDES IPT2 448 : MOTERR(1:8)='MAILLAGE' 449 : MOTERR(9:16)='MCHAML' 449 : CALL ERREUR(135) 450 : ENDIF 451 : GOTO 900 452 : ELSE 453 : SEGSUP,ITRA 454 : if (imel.ne.meleme) segdes meleme 455 : if (imel.ne.ipt2) SEGDES IPT2 456 : ENDIF 457 : * 458 : * fin de la boucle sur les sous zones du maillage 459 : * 460 : IF (IBOU1 .LT. NBSOU1 ) GOTO 1 461 : * 462 : IF(NS01.NE.NSOUS) THEN 463 : N1=NS01 464 : N3=MCHEL1.INFCHE(/2) 465 : L1=MCHEL1.TITCHE(/1) 466 : SEGINI MCHEL2 467 : MCHEL2.TITCHE=TITCHE 468 : MCHEL2.IFOCHE=IFOCHE 469 : DO 31 ISOUS=1,NS01 470 : MCHEL2.ICHAML(ISOUS) = MCHEL1.ICHAML(ISOUS) 471 : MCHEL2.IMACHE(ISOUS) = MCHEL1.IMACHE(ISOUS) 472 : MCHEL2.CONCHE(ISOUS) = MCHEL1.CONCHE(ISOUS) 473 : NSS=NZONE(ISOUS) 474 : DO 33 N33=1,N3 475 : MCHEL2.INFCHE(ISOUS,N33) = MCHEL1.INFCHE(NSS,N33) 476 : 33 CONTINUE 477 : 31 CONTINUE 478 : SEGSUP MCHEL1 479 : MCHEL1=MCHEL2 480 : ENDIF 481 : * SEGDES,MCHEL1 482 : IRET=MCHEL1 483 : 484 : 900 CONTINUE 485 : SEGDES MCHELM 486 : SEGSUP,ICPR,ICPEL,NZONE,WTRAV 487 : 488 : RETURN 489 : END 490 : 491 : 492 :
© Cast3M 2003 - All rights reserved.
Disclaimer