1 : C NUMOP1 SOURCE CHAT 05/01/13 02:03:25 5004 2 : C RACINE DE LA NUMEROTATION POUR LA SORTIE SUR FAC 3 : C 4 : C methode utilisee: NESTED DISSECTION & MULTIFRONTALE 5 : 6 : C ce programme decompose un domaine d'etude par elements finis 7 : C et utilise ensuite la methode multifrontale. 8 : 9 : 10 : SUBROUTINE NUMOP1(MELEME,ICPR,NODES) 11 : IMPLICIT INTEGER(I-N) 12 : -INC SMELEME 13 : -INC SMCOORD 14 : -INC CCOPTIO 15 : 16 : SEGMENT JMEM(NODES+1),JMEMN(NODES+1) 17 : C JMEM et JMEMN contiennent le nombre d'elements auquel appartient 18 : C un noeud 19 : 20 : SEGMENT JNT(NODES),NUMERO(NODES) 21 : C JNT contient la nouvelle numerotation 22 : 23 : SEGMENT ICPR(XCOOR(/1)/(IDIM+1)) 24 : C ICPR au debut contient l'ancienne numerotation, 25 : C a la fin la nouvelle. 26 : 27 : SEGMENT IADJ(NODES+1),IVOIS(0) 28 : SEGMENT NPVOIS(2*NBENS),NVOIS(0) 29 : C IADJ(i) pointe sur IVOIS qui contient les voisins de i entre 30 : C IADJ(i) et IADJ(i+1)-1 31 : 32 : SEGMENT LAGRAN(NB) 33 : C contient les noeud de lagrange et les noeuds les suivant directement 34 : C cf element de type 22 35 : SEGMENT INVLAG(NB1) 36 : 37 : SEGMENT BOOLEEN 38 : LOGICAL BOOL(NODES) 39 : ENDSEGMENT 40 : C BOOL(i) = true si le noeud i a etait deja memtioner dans la liste 41 : C des voisin IVOIS 42 : 43 : SEGMENT IMEMOIR(NBV) 44 : C contient les element appartenant a chaque noeud,sous forme de liste, 45 : C au premier, puis au second, etc... 46 : 47 : INTEGER ELEM 48 : C nom d'un element 49 : 50 : INTEGER N,NBENS 51 : 52 : SEGMENT MASQUE 53 : LOGICAL MASQ(NODES) 54 : ENDSEGMENT 55 : C MASQ(X)=.TRUE. si le noeud X n'a pas ete renumerote; 56 : C .FALSE. si il l'a ete. 57 : 58 : INTEGER DIM,DIMSEP 59 : C DIM= nombre de noeuds renumerotes. 60 : 61 : INTEGER PIVOT 62 : C PIVOT est le noeud utile a la division du domaine. 63 : 64 : SEGMENT IPOS(NODES*2) 65 : C est le vecteur contenant la numerotation dans les deux sens,de 1 a NODES 66 : C puis de NODES+1 a 2* NODES, cf la subroutine SEPAR 67 : 68 : SEGMENT ICOMPT(NODES),IORDRE(2*NBENS) 69 : 70 : 71 : 72 : C initialisation 73 : C******************************* 74 : 75 : 76 : 77 : IENORM=2000000000 78 : C norme d'erreur. 79 : 80 : NODES=XCOOR(/1)/(IDIM+1) 81 : SEGACT ICPR*MOD 82 : SEGACT MELEME 83 : C icpr: numero des noeuds. 84 : C meleme: objet de maillage (cf assem2.eso) 85 : 86 : DO 10 I=1,ICPR(/1) 87 : 10 ICPR(I)=0 88 : 89 : IPT1=MELEME 90 : IKOU=0 91 : NBV=0 92 : NB1=0 93 : NB2=0 94 : 95 : DO 100 IO=1,MAX(1,LISOUS(/1)) 96 : IF (LISOUS(/1).NE.0) THEN 97 : IPT1=LISOUS(IO) 98 : SEGACT IPT1 99 : ENDIF 100 : C on cree la numerotation des noeuds. 101 : C 'nb noeuds/element'=IPT1.NUM(/1) 102 : C 'nb element'=IPT1.NUM(/2) 103 : IF(IPT1.ITYPEL.EQ.22) THEN 104 : NB1=NB1+IPT1.NUM(/2) 105 : NB2=MAX(NB2,IPT1.NUM(/1)) 106 : C NB1= nbre d'elements de type 22. 107 : C NB2=nbre de noeuds/element maximum parmi 108 : C les elements de type 22. 109 : ENDIF 110 : DO 150 I=1,IPT1.NUM(/1) 111 : DO 150 J=1,IPT1.NUM(/2) 112 : IJ=IPT1.NUM(I,J) 113 : C IJ est le Ieme noeud du Jeme element. 114 : IF (ICPR(IJ).NE.0) GOTO 150 115 : C s'il est deja numerote, on ne fait rien. 116 : IKOU=IKOU+1 117 : ICPR(IJ)=IKOU 118 : 150 CONTINUE 119 : NODES=IKOU 120 : C NODES= nbre de noeuds. 121 : 100 CONTINUE 122 : 123 : NB=NB2*NB1 124 : C NB= dimension(LAGRAN). 125 : 126 : C***** initalisation des segments********* 127 : 128 : SEGINI IADJ,IVOIS,JMEM,JMEMN,LAGRAN 129 : SEGINI BOOLEEN,JNT 130 : 131 : DO 20 I=1,NODES+1 132 : IADJ(I)=0 133 : JMEM(I)=0 134 : JMEMN(I)=0 135 : 20 CONTINUE 136 : C****************************************** 137 : 138 : 139 : 140 : IPT1=MELEME 141 : NGRAND=0 142 : IADJ(1)=1 143 : INC=0 144 : 145 : 146 : DO 200 IO=1,MAX(1,LISOUS(/1)) 147 : 148 : IF (LISOUS(/1).NE.0) THEN 149 : IPT1=LISOUS(IO) 150 : ENDIF 151 : 152 : DO 210 J=1,IPT1.NUM(/2) 153 : 154 : IF(IPT1.ITYPEL.EQ.22) THEN 155 : DO 220 I=1,IPT1.NUM(/1) 156 : C chaque element de type 22 a au plus NB2 noeuds. 157 : LAGRAN(INC*NB2+I)=ICPR(IPT1.NUM(I,J)) 158 : C les noeuds de l'elements de type 22 159 : C sont ranges dans le vecteur LAGRAN. 160 : 220 CONTINUE 161 : DO 225 I=IPT1.NUM(/1)+1,NB2 162 : LAGRAN(INC*NB2+I)=0 163 : C comme on a alloue la place memeoire maximale, 164 : C on remplit les cases restantes avec des 0. 165 : 225 CONTINUE 166 : INC=INC+1 167 : C INC=le nbre d'elements de type 22. 168 : ENDIF 169 : 170 : DO 230 I=1,IPT1.NUM(/1) 171 : IJ=ICPR(IPT1.NUM(I,J))+1 172 : JMEM(IJ)=JMEM(IJ)+1 173 : C JMEM(I+1): nb elements auquel le noeud I appartient 174 : 230 CONTINUE 175 : 176 : 210 CONTINUE 177 : 178 : NGRAND=MAX(NGRAND,IPT1.NUM(/2)) 179 : 180 : 200 CONTINUE 181 : 182 : 183 : 184 : NGRAND=NGRAND+1 185 : JMEM(1)=1 186 : DO 30 I=1,NODES 187 : JMEM(I+1)=JMEM(I)+JMEM(I+1) 188 : C JMEM(I+1)=indice de depart des elements 189 : C auxquels le noeud I appartient. 190 : 30 CONTINUE 191 : NBV=JMEM(NODES+1) 192 : C NBV= dimension de IMEMOIR. 193 : SEGINI IMEMOIR 194 : 195 : 196 : IPT1=MELEME 197 : 198 : IMATOT=LISOUS(/1)*(NGRAND+1) 199 : IF(IMATOT.GT.IENORM) THEN 199 : CALL ERREUR(382) 200 : RETURN 201 : ENDIF 202 : 203 : 204 : 205 : DO 300 IO=1,MAX(1,LISOUS(/1)) 206 : IF (LISOUS(/1).NE.0) THEN 207 : IPT1=LISOUS(IO) 208 : ENDIF 209 : DO 350 I=1,IPT1.NUM(/1) 210 : DO 350 J=1,IPT1.NUM(/2) 211 : JMEMN(ICPR(IPT1.NUM(I,J))+1)=JMEMN(ICPR(IPT1.NUM(I,J))+1)+1 212 : IMEMOIR(JMEM(ICPR(IPT1.NUM(I,J)))+ 213 : & JMEMN(ICPR(IPT1.NUM(I,J))+1)-1)=J+IO*NGRAND 214 : C on range dans IMEMOIR tous les elements des sous-objets 215 : C IO auxquels appartient le noeud ICPR(IPT1.NUM(I,J)). 216 : C On connait pour chaque element, le sous-objet auquel 217 : C il appartient. 218 : 350 CONTINUE 219 : 300 CONTINUE 220 : 221 : 222 : DO 400 I=1,NODES 223 : IADJ(I+1)=IADJ(I) 224 : DO 410 J=1,NODES 225 : BOOL(J)=.FALSE. 226 : 410 CONTINUE 227 : DO 420 J=JMEM(I),JMEM(I+1)-1 228 : ELEM=IMEMOIR(J) 229 : C ELEM=element auquel appartient le noeud I. 230 : 231 : IPT1=MELEME 232 : IF (LISOUS(/1).NE.0) IPT1=LISOUS(ELEM/NGRAND) 233 : ELEM=MOD(ELEM,NGRAND) 234 : DO 430 K=1,IPT1.NUM(/1) 235 : C k representatif du nb de noeuds par elements. 236 : IK=ICPR(IPT1.NUM(K,ELEM)) 237 : IF ((I.NE.IK).AND. 238 : & (.NOT.(BOOL(IK)))) THEN 239 : C si i n'est pas egal a un des nouveaux numeros des noeuds 240 : C de l'element ELEM et si il n'appartient pas deja a l'ens des 241 : C voisins du noeud i(jadjc(i)),alors on le rajoute. 242 : C IVOIS(IADJ(I+1))=IK 243 : IVOIS(**)=IK 244 : IADJ(I+1)=IADJ(I+1)+1 245 : BOOL(IK)=.TRUE. 246 : ENDIF 247 : 430 CONTINUE 248 : 249 : 250 : 420 CONTINUE 251 : 400 CONTINUE 252 : 253 : SEGSUP JMEM,JMEMN,BOOLEEN 254 : SEGSUP IMEMOIR 255 : 256 : 257 : C affectation 258 : C************************ 259 : 260 : 261 : N=NODES 262 : SEGINI IPOS,MASQUE 263 : 264 : DO 50 I=1,N 265 : MASQ(I)=.TRUE. 266 : IPOS(I)=0 267 : IPOS(NODES+I)=0 268 : 50 CONTINUE 269 : C initialement, les noeuds ne sont pas masques,ont donc 270 : C une position nulle. 271 : 272 : DIM=0 273 : C le nombre de noeuds renumerotes DIM est initialement egal a zero. 274 : 275 : 276 : C **************************************** 277 : C boucle principale 278 : MPOS=0 279 : DO 500 I=1,N 280 : 281 : 550 IF(.NOT.MASQ(I)) GOTO 500 282 : C si le noeud a ete renumerote, on passe au suivant. 283 : 284 : PIVOT=I 285 : C sinon, il devient PIVOT. 286 : 286 : CALL SEPAR(IADJ,IVOIS,PIVOT,MASQUE,DIMSEP,N,IPOS,NODES) 287 : C separe le domaine d'etude en 2 parties. 288 : C on decrit le domaine d'etude a partir du pivot et on cherche la 289 : C longueur maximale en decrivant les voisins de pivot, et leurs 290 : C voisins... jusqu'a rencontrer un voisin masque. On cree alors 291 : C une nouvelle separation. 292 : C les noeuds masques delimitent la separation du domaine. 293 : 294 : DIM=DIM+DIMSEP 295 : 296 : N=N-DIMSEP 297 : 298 : IF (DIM.GE.NODES) GOTO 600 299 : C si tous les noeuds ont une position non nulle,on arrete. 300 : 301 : GOTO 550 302 : 303 : 500 CONTINUE 304 : 305 : 600 SEGSUP MASQUE 306 : SEGINI NUMERO 307 : 307 : CALL PREPA(IPOS,NPVOIS,NVOIS,NBENS,NUMERO,NODES,IADJ,IVOIS) 308 : C procedure donnant une nouvelle gestion de donnees. 309 : C NBENS= nbre d'ensembles. 310 : C un ensemble est compose de noeuds ayant meme NUMERO. 311 : C NPVOIS(I,2)-NPVOIS(I,1)-1=nbre d'ensemble voisins a l'ensemble I. 312 : C NVOIS(J), pour J compris entre NPVOIS(I,1) et NPVOIS(I,2)-1, 313 : C contient les ensembles voisins de I. 314 : 314 : CALL MINDEG(NPVOIS,NVOIS,NBENS,IORDRE) 315 : C on recherche l'ensemble de minimum degre. 316 : C on le numerote et on recommence. 317 : 317 : CALL SORTID(IADJ,IVOIS,IORDRE,NUMERO,NODES,NBENS,JNT) 318 : C JNT(I) est le nouveau numero du NOEUD I. 319 : 320 : SEGSUP IORDRE 321 : 322 : C****************************************************** 323 : C insertion des lagrangiens. 324 : C****************************************************** 325 : C principe: on cherche le premier noeud auquel s'applique les 2 326 : C coefficients de lagrange, on les insere de part et d'autre de 327 : C ce noeud dans la numerotation. 328 : 329 : DO 1 I=1,NODES 330 : IPOS(JNT(I))=I 331 : 1 CONTINUE 332 : 333 : IF(NB2.NE.0) THEN 334 : N=NB/NB2 335 : ELSE 336 : N=0 337 : ENDIF 338 : 339 : NORME=0 340 : DO 700 I=1,NODES 341 : NOEUD=IPOS(I) 342 : IF (NOEUD.EQ.0) GOTO 700 343 : DO 750 J=0,N-1 344 : IF(NOEUD.EQ.LAGRAN(J*NB2+1)) THEN 345 : NORME=NORME+1 346 : IPOS(NODES+NOEUD)=0 347 : GOTO 700 348 : ENDIF 349 : 350 : IF(NOEUD.EQ.LAGRAN(J*NB2+2)) THEN 351 : NORME=NORME+1 352 : IPOS(NODES+NOEUD)=0 353 : GOTO 700 354 : ENDIF 355 : 750 CONTINUE 356 : IPOS(NODES+NOEUD)=I-NORME 357 : C IPOS(NODES+NOEUD)= numero d'elimination obtenu en tenant 358 : C compte des lagrangiens. 359 : 360 : 700 CONTINUE 361 : 362 : NCOMPT=0 363 : NORM=0 364 : SEGINI ICOMPT,INVLAG 365 : 366 : DO 800 I=1,NODES 367 : 368 : NOEUD=IPOS(I) 369 : C NOEUD =Ieme noeud elimine. 370 : 371 : IF(IPOS(NODES+NOEUD).EQ.0) GOTO 800 372 : 373 : ICOMPT(NOEUD)=0 374 : 375 : DO 850 J=0,NB1-1 376 : IF(LAGRAN(J*NB2+1).NE.0) THEN 377 : K=3 378 : 855 IF((K.LE.NB2).AND.(LAGRAN(J*NB2+K).NE.0)) THEN 379 : 380 : IF(NOEUD.EQ.LAGRAN(J*NB2+K)) THEN 381 : 382 : INVLAG(ICOMPT(NOEUD)+1)=J 383 : C INVLAG mémorise l'élément de type 22. 384 : 385 : C on change le numero d'élimination de noeud, 386 : C en insérant de part et d'autre les coef de lagrange. 387 : 388 : JNT(LAGRAN(J*NB2+1))= 389 : & IPOS(NOEUD+NODES)+NORM-ICOMPT(NOEUD) 390 : 391 : LAGRAN(J*NB2+1)=0 392 : C permet de savoir si le lagrangien a ete insere. 393 : 394 : NORM=NORM+1 395 : JNT(NOEUD)= 396 : & IPOS(NODES+NOEUD)+NORM-ICOMPT(NOEUD) 397 : C renumerotation de noeud. 398 : 399 : DO 852 IJ=1,ICOMPT(NOEUD) 400 : IL=INVLAG(IJ) 401 : JNT(LAGRAN(IL*NB2+2))=JNT(LAGRAN(IL*NB2+2))+1 402 : 852 CONTINUE 403 : 404 : NORM=NORM+1 405 : JNT(LAGRAN(J*NB2+2))=IPOS(NODES+NOEUD)+NORM 406 : ICOMPT(NOEUD)=ICOMPT(NOEUD)+1 407 : GOTO 850 408 : ELSE 409 : K=K+1 410 : GOTO 855 411 : ENDIF 412 : 413 : ENDIF 414 : 415 : ENDIF 416 : 417 : 850 CONTINUE 418 : 419 : IF(ICOMPT(NOEUD).EQ.0) THEN 420 : JNT(NOEUD)=IPOS(NODES+NOEUD)+NORM 421 : ENDIF 422 : 423 : 800 CONTINUE 424 : 425 : SEGSUP INVLAG 426 : C****************************************************************** 427 : DO 860 I=1,XCOOR(/1)/(IDIM+1) 428 : IF(ICPR(I).NE.0) THEN 429 : ICPR(I)=JNT(ICPR(I)) 430 : C numerotation finale. 431 : ENDIF 432 : 860 CONTINUE 433 : 434 : SEGSUP IADJ,IVOIS,JNT,IPOS,ICOMPT,NUMERO 435 : SEGSUP LAGRAN,NPVOIS,NVOIS 436 : 437 : RETURN 438 : END 439 : 440 :
© Cast3M 2003 - All rights reserved.
Disclaimer