1 : C NUMOP2 SOURCE PV 09/03/30 21:15:00 6349 2 : C RACINE DE LA NUMEROTATION POUR LA SORTIE SUR FAC 3 : C 4 : C methode utilisee: NESTED DISSECTION. 5 : 6 : C ce programme refait l'indicage de la matrice afin de minimiser 7 : C le profil. 8 : 9 : 10 : SUBROUTINE NUMOP2(MELEME,ICPR,NODES) 11 : IMPLICIT INTEGER(I-N) 12 : -INC SMELEME 13 : -INC SMCOORD 14 : -INC CCOPTIO 15 : -INC CCASSIS 16 : 17 : SEGMENT JMEM(NODES+1),JMEMN(NODES+1) 18 : C JMEM et JMEMN contiennent le nombre d'element auquel appartient un noeud 19 : 20 : SEGMENT JNT(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) 28 : SEGMENT JADJC(0) 29 : C IADJ(i) pointe sur JADJC 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 ete deja mentionne dans la liste 41 : C des voisins JADJC. 42 : 43 : SEGMENT IMEMOIR(NBV),LMEMOIR(NBV) 44 : C contient les elements appartenant a chaque noeud,sous forme de liste. 45 : 46 : INTEGER ELEM 47 : C nom d'un element 48 : 49 : INTEGER N 50 : 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*3) 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 : C 68 : C segments utilisés dans sepa2 69 : C 70 : SEGMENT NRELONG(NODES*nbthr) 71 : C NRELONG contient pour chaque noeud sa profondeur. 72 : 73 : SEGMENT NOELON(NODES*nbthr) 74 : SEGMENT NOEL2(NODES) 75 : SEGMENT LONDIM(NODES*nbthr) 76 : C NOELON contient les noeuds de profondeur LONG. 77 : C DIMLON= dimension de NOELON. 78 : 79 : 80 : C********************************** 81 : 82 : C debut du program 83 : 84 : C********************************** 85 : 86 : 87 : 88 : C initialisation 89 : C******************************* 90 : 91 : 92 : C norme d'erreur 93 : SEGACT ICPR*MOD 94 : NODES=ICPR(/1) 95 : SEGACT MELEME 96 : C icpr: numero des noeuds. 97 : C meleme: objet de maillage (cf assem2.eso) 98 : 99 : DO 10 I=1,ICPR(/1) 100 : ICPR(I)=0 101 : 10 CONTINUE 102 : 103 : IPT1=MELEME 104 : IKOU=0 105 : NBV=0 106 : NB1=0 107 : NB2=0 108 : 109 : DO 100 IO=1,MAX(1,LISOUS(/1)) 110 : IF (LISOUS(/1).GT.0) THEN 111 : IPT1=LISOUS(IO) 112 : SEGACT IPT1 113 : ENDIF 114 : C on cree la numerotation des noeuds. 115 : C 'nb noeuds/element'=IPT1.NUM(/1) 116 : C 'nb element'=IPT1.NUM(/2) 117 : IF(abs(IPT1.ITYPEL).EQ.22) THEN 118 : NB1=NB1+IPT1.NUM(/2) 119 : NB2=MAX(NB2,IPT1.NUM(/1)) 120 : C NB1= nbre d'éléments de type 22. 121 : C NB2=nbre de noeuds/élément maximum parmi 122 : C les éléments de type 22. 123 : ENDIF 124 : DO 150 J=1,IPT1.NUM(/2) 125 : DO 150 I=1,IPT1.NUM(/1) 126 : IJ=IPT1.NUM(I,J) 127 : C IJ est le Ième noeud du Jème élément. 128 : IF (ICPR(IJ).EQ.0) THEN 129 : C s'il est déjà numéroté, on ne fait rien. 130 : IKOU=IKOU+1 131 : ICPR(IJ)=IKOU 132 : ENDIF 133 : 150 CONTINUE 134 : 100 CONTINUE 135 : 136 : NODES=IKOU 137 : NB=NB2*NB1 138 : 139 : C***** initalisation des segments********* 140 : 141 : SEGINI IADJ,JADJC,JMEM,JMEMN,LAGRAN 142 : SEGINI BOOLEEN,JNT 143 : 144 : DO 20 I=1,NODES+1 145 : IADJ(I)=0 146 : JMEM(I)=0 147 : JMEMN(I)=0 148 : 20 CONTINUE 149 : 150 : C****************************************** 151 : 152 : IPT1=MELEME 153 : IADJ(1)=1 154 : INC=0 155 : DO 200 IO=1,MAX(1,LISOUS(/1)) 156 : 157 : IF (LISOUS(/1).NE.0) THEN 158 : IPT1=LISOUS(IO) 159 : ENDIF 160 : 161 : DO 210 J=1,IPT1.NUM(/2) 162 : IF(abs(IPT1.ITYPEL).EQ.22) THEN 163 : is=sign(1,ipt1.itypel) 164 : DO 220 I=1,IPT1.NUM(/1) 165 : C chaque element de type 22 a au plus NB2 noeuds. 166 : LAGRAN(INC*NB2+I)=ICPR(IPT1.NUM(I,J))*is 167 : C les noeuds de l'elements de type 22 168 : C sont ranges dans le vecteur LAGRAN. 169 : 220 CONTINUE 170 : DO 225 I=IPT1.NUM(/1)+1,NB2 171 : LAGRAN(INC*NB2+I)=0 172 : C comme on a alloue la place memoire maximale, 173 : C on remplit les cases restantes avec des 0. 174 : 225 CONTINUE 175 : INC=INC+1 176 : C INC=le nbre d'elements de type 22. 177 : ENDIF 178 : 179 : DO 230 I=1,IPT1.NUM(/1) 180 : IJ=ICPR(IPT1.NUM(I,J))+1 181 : JMEM(IJ)=JMEM(IJ)+1 182 : C JMEM(I+1): nb elements auquel le noeud I appartient 183 : 230 CONTINUE 184 : 210 CONTINUE 185 : 186 : 200 CONTINUE 187 : 188 : 189 : JMEM(1)=1 190 : DO 30 I=1,NODES 191 : JMEM(I+1)=JMEM(I)+JMEM(I+1) 192 : C JMEM(I+1)=indice de depart des elements 193 : C auxquels le noeud I appartient. 194 : 30 CONTINUE 195 : NBV=JMEM(NODES+1) 196 : C NBV= dimension de IMEMOIR. 197 : SEGINI IMEMOIR,LMEMOIR 198 : 199 : 200 : 201 : IPT1=MELEME 202 : 203 : DO 300 IO=1,MAX(1,LISOUS(/1)) 204 : IF (LISOUS(/1).NE.0) THEN 205 : IPT1=LISOUS(IO) 206 : ENDIF 207 : DO 350 J=1,IPT1.NUM(/2) 208 : DO 350 I=1,IPT1.NUM(/1) 209 : IJ=ICPR(IPT1.NUM(I,J)) 210 : JMEMN(IJ+1)=JMEMN(IJ+1)+1 211 : IMEMOIR(JMEM(IJ)+JMEMN(IJ+1)-1)=J 212 : LMEMOIR(JMEM(IJ)+JMEMN(IJ+1)-1)=IO 213 : C on range dans IMEMOIR tous les elements des sous-objets 214 : C IO auxquels appartient le noeud ICPR(IPT1.NUM(I,J)). 215 : C On connait pour chaque element, le sous-objet auquel 216 : C il appartient grace a LMEMOIR 217 : 350 CONTINUE 218 : 300 CONTINUE 219 : 220 : DO 410 J=1,NODES 221 : BOOL(J)=.FALSE. 222 : 410 CONTINUE 223 : DO 400 I=1,NODES 224 : IADJ(I+1)=IADJ(I) 225 : DO 420 J=JMEM(I),JMEM(I+1)-1 226 : ELEM=IMEMOIR(J) 227 : C ELEM=element auquel appartient le noeud I. 228 : 229 : IPT1=MELEME 230 : IF (LISOUS(/1).NE.0) IPT1=LISOUS(LMEMOIR(J)) 231 : DO 430 K=1,IPT1.NUM(/1) 232 : C k representatif du nb de noeuds par elements. 233 : IK=ICPR(IPT1.NUM(K,ELEM)) 234 : IF ((I.NE.IK).AND. 235 : & (.NOT.(BOOL(IK)))) THEN 236 : C si i n'est pas egal a un des nouveaux numeros des noeuds 237 : C de l'element ELEM et si il n'appartient pas deja a l'ens des 238 : C voisins du noeud i(jadjc(i)),alors on le rajoute. 239 : C JADJC(IADJ(I+1))=IK 240 : JADJC(**)=IK 241 : IADJ(I+1)=IADJ(I+1)+1 242 : BOOL(IK)=.TRUE. 243 : ENDIF 244 : 430 CONTINUE 245 : 420 CONTINUE 246 : * remise a faux de bool 247 : DO 412 J=IADJ(I),IADJ(I+1)-1 248 : IK=JADJC(J) 249 : BOOL(IK)=.FALSE. 250 : 412 CONTINUE 251 : 400 CONTINUE 252 : 253 : SEGSUP JMEM,JMEMN,IMEMOIR,LMEMOIR,BOOLEEN 254 : 255 : 256 : 257 : C************************************************************************** 258 : 259 : 260 : C affectation 261 : C************************ 262 : 263 : 264 : SEGINI IPOS,MASQUE 265 : IPOSMAX=0 266 : 267 : DO 50 I=1,NODES 268 : MASQ(I)=.TRUE. 269 : IPOS(I)=0 270 : IPOS(NODES+I)=0 271 : 50 CONTINUE 272 : C initialement, les noeuds ne sont pas masques,ont donc 273 : C une position nulle. 274 : 275 : DIM=0 276 : C le nombre de noeuds renumerotes DIM est initialement egal a zero. 277 : 278 : 279 : C **************************************** 280 : C boucle principale 281 : NS=NODES 282 : 283 : nbthr=min(64,nbesc) 284 : nbthr=max(1,nbesc) 284 : if (nbthr.gt.1) call threadii 285 : SEGINI NRELONG,NOELON,noel2,londim 286 : ** write (6,*) ' avant appel sepa2 ' 287 : DO 500 I=1,NODES 288 : 550 IF(.NOT.MASQ(I)) GOTO 500 289 : C si le noeud est masque alors ne rien faire: il est deja 290 : C renumerote. On passe au noeud suivant. 291 : 292 : PIVOT=I 293 : 293 : CALL SEPA2(IADJ,JADJC,PIVOT,MASQUE,DIMSEP,NS, 294 : > IPOS,NODES,IPOSMAX,nrelong,noelon,noel2, 295 : > londim,nbthr) 296 : C separe le domaine d'etude en 2 parties. 297 : C on decrit le domaine d'etude a partir du pivot et on cherche la 298 : C longueur maximale en decrivant les voisins de pivot, et leurs 299 : C voisins... jusqu'a rencontrer un voisin masque. On cree alors 300 : C une nouvelle separation. 301 : C les noeuds masques delimitent la separation du domaine. 302 : 303 : 304 : DIM=DIM+DIMSEP 305 : NS=NS-DIMSEP 306 : C la dimension de noeuds renumerotes est augmente de DIMSEP. 307 : C Celle de noeuds a renumeroter est diminue de DIMSEP. 308 : 309 : * IF (DIM.GE.NODES) GOTO 600 310 : C si tous les noeuds ont ete renumerotes, on arrete. 311 : 312 : GOTO 550 313 : 314 : 500 CONTINUE 315 : ** write (6,*) ' apres appel sepa2 ' 316 : 317 : SEGSUP NRELONG,NOELON,noel2,londim 317 : if (nbesc.gt.1) call threadis 318 : 319 : 600 CONTINUE 320 : * 321 : * tri dans chaque zone 322 : do 610 lpoint=1,nodes 323 : mdomn=ipos(lpoint+nodes) 324 : iposv=ipos(mdomn+1) 325 : iposi=3*nodes 326 : do 620 kk=iadj(lpoint),iadj(lpoint+1)-1 327 : k=jadjc(kk) 328 : iposk=ipos(ipos(k+nodes)+1) 329 : if (iposk.ne.iposv) then 330 : iposi=min(iposi,iposk) 331 : endif 332 : 620 continue 333 : ipos(lpoint+2*nodes)=iposi 334 : 610 continue 335 : 336 : * 337 : * mise à la bonne place des multiplicateurs de Lagrange 338 : do 700 J=0,NB1-1 339 : iposvs=2**30 340 : iposvr=-2**30 341 : ipr=0 342 : ips=0 343 : mdomnr=0 344 : mdomns=0 345 : if (nb2.eq.2) goto 700 346 : * write (6,*) 'numop2 ',(lagran(J*NB2+il),il=1,nb2) 347 : do 800 il=3,nb2 348 : ip = abs(LAGRAN(J*NB2+il)) 349 : * if (ip.eq.0) write (6,*) ' prob numop2 ' 350 : if (ip.eq.0) goto 800 351 : mdomn=ipos(ip+nodes) 352 : iposv=ipos(mdomn+1) 353 : if (iposv.gt.iposvr) then 354 : ipr=ip 355 : iposvr=iposv 356 : mdomnr=mdomn 357 : endif 358 : if (iposv.lt.iposvs) then 359 : ips=ip 360 : iposvs=iposv 361 : mdomns=mdomn 362 : endif 363 : 800 continue 364 : *pv if (ipr.eq.0) write (6,*) ' probleme ipr numop2' 365 : *pv if (ips.eq.0) write (6,*) ' probleme ips numop2' 366 : ip=abs(LAGRAN(J*NB2+2)) 367 : IPOS(IP+2*NODES)= 3*nodes+2 368 : * numéroter le frottement apres le contact 369 : if (LAGRAN(J*NB2+1).lt.0) IPOS(IP+2*NODES)=IPOS(IP+2*NODES)+1 370 : IPOS(IP+nodes)=mdomnr 371 : ip=abs(LAGRAN(J*NB2+1)) 372 : IPOS(IP+2*NODES)= -2 373 : if (LAGRAN(J*NB2+2).lt.0) IPOS(IP+2*NODES)=IPOS(IP+2*NODES)-1 374 : IPOS(IP+nodes)=mdomns 375 : 700 continue 376 : SEGSUP IADJ,JADJC,MASQUE 377 : * ok maintenant on trie 377 : CALL SORTI2(IPOS,JNT,NODES) 378 : 379 : C*************************************************************************** 380 : 381 : DO 860 I=1,ICPR(/1) 382 : IF(ICPR(I).NE.0) ICPR(I)=JNT(ICPR(I)) 383 : C numerotation finale. 384 : 860 CONTINUE 385 : 386 : 387 : SEGSUP JNT,IPOS,LAGRAN 388 : 389 : 390 : RETURN 391 : END 392 : 393 : 394 : 395 : 396 : 397 : 398 : 399 : 400 : 401 : 402 : 403 : 404 : 405 : 406 : 407 : 408 : 409 :
© Cast3M 2003 - All rights reserved.
Disclaimer