1 : C NOR4D3 SOURCE FANDEUR 13/01/29 21:16:06 7683 2 : SUBROUTINE NOR4D3( 3 : & MELEFA,MELEFL,MLECEN,MELEFP,MLESOM,MPONOR, 4 : & MPOSUR,MELTFA,MLEFA,MLEFA2,MPOTEN,MPOCHP,MLENCL, 5 : & MPOVCL,ICHTE,ICHCL,ICHCO,MPOVCO,IOP, 6 : & IPO2,SCMB,INDLI,VAL1,VAL2,IND22,IND2,IND, 7 : & IPO3,TAB,MPOGRA,MELVA1,MELVA2, 8 : & NSOMM,NBMAX,NBFAC,NBCOT,MCHEL2,MCHAM2) 9 : C 10 : C PROJET : CASTEM 2000 11 : C 12 : C NOM : NORV4 13 : C 14 : C DESCRIPTION : Appelle par NORV1 15 : C 16 : C LANGAGE : FORTRAN 77 + ESOPE 2000 (avec extensions CISI) 17 : C 18 : C AUTEUR : C. LE POTIER, DM2S/SFME/MTMS 19 : C 20 : C************************************************************************ 21 : C 22 : IMPLICIT INTEGER(I-N) 23 : IMPLICIT real*8 (a-h,o-z) 24 : -INC SMLENTI 25 : -INC SMELEME 26 : -INC SMCHPOI 27 : -INC CCOPTIO 28 : -INC SMCOORD 29 : -INC SMLREEL 30 : POINTEUR MELEFL.MELEME, MELEFP.MELEME, MELEFA.MELEME, 31 : & MELTFA.MELEME,MELEP2.MELEME 32 : POINTEUR MPOSUR.MPOVAL, MPONOR.MPOVAL, 33 : & MPOCHP.MPOVAL, MPOVCL.MPOVAL, MPGSOM.MPOVAL, MPVOSO.MPOVAL, 34 : & MPOGRA.MPOVAL,MPOTEN.MPOVAL,MPOVCO.MPOVAL 35 : POINTEUR MLENCL.MLENTI, MLECEN.MLENTI, MLESOM.MLENTI, 36 : & MLEFA.MLENTI,MLEFA2.MLENTI 37 : -INC SMCHAML 38 : INTEGER NBNN,NBREF,NBMAX 39 : 40 : 41 : 42 : C 43 : C**** Variables de COOPTIO 44 : C 45 : C 46 : C**** Variables de COOPTIO 47 : C 48 : INTEGER IPLLB, IERPER, IERMAX, IERR, INTERR 49 : & ,IOTER, IOLEC, IOIMP, IOCAR, IOACQ 50 : & ,IOPER, IOSGB, IOGRA, IOSAU, IORES 51 : & ,IECHO, IIMPI, IOSPI 52 : & ,IDIM 53 : C & ,MCOORD 54 : & ,IFOMOD, NIFOUR, IFOUR, NSDPGE, IONIVE 55 : & ,NGMAXY, IZROSF, ISOTYP, IOSCR,LTEXLU 56 : & ,NORINC,NORVAL,NORIND,NORVAD 57 : & ,NUCROU, IPSAUV 58 : 59 : C**** Variable de SMLENTI, SMCHPOI 60 : C 61 : INTEGER JG, N, NC, NSOUPO, NAT, NBSOUS, NBNO,NBELEM 62 : C 63 : C**** Les includes 64 : C 65 : INTEGER I1,ICOMP,ICOMGR,IGEOM 66 : & ,IOP1,ICEN,ISOMM,IFAC,IFACEL,IFACEP,INORM 67 : & ,ISURF,IMAIL,ICHPO,ICHCL,ICHGRA,ICOEFF 68 : & ,NTOT,NSOMM,NCOMP,NFAC,NCEN 69 : & ,NLCF,NGCF,NGCF1,NGCF2,NGCG,NGCD,NLCG,NLCD,NGS1,NGS2 70 : & ,NLS1,NLS2,NLFCL 71 : & ,ISOUS,IELEM,INOEUD,ICELL,NCON,NFIN,INDGA,INDDR 72 : INTEGER ICEN2 73 : REAL*8 SCNX,SCNY,SCNZ,SURF,VOL,VAL,VALX,VALY,VALZ,XG,XD,XF,XS1,XS2 74 : & ,YG,YD,YF,YS1,YS2,ZG,ZD,ZF,ZS1,ZS2, 75 : & PSCA,XNORM,VECX,VECY,PSCAGX,PSCAGY,PSCAGZ, 76 : & PSCADX,PSCADY,PSCADZ,K11G,K22G,K21G,K31G,K32G,K33G, 77 : & K11D,K22D,K21D,K31D,K32D,K33D,VXG1,VXG2, 78 : & VXAU,VYAU,VZAU,VXD1,VXD2,VYG1,VYG2,VYD1,VYD2,VZG1,VZG2,VZD1,VZD2, 79 : & TRG1,TRG2,TRG3, 80 : & TRD1,TRD2,TRD3,TRG,TRD,AUX,AUY,AUZ,AUXMA,THETA,COEFDD 81 : REAL*8 XLONG,AG1,AG2,AD1,AD2,PSCAG1,PSCAG2,PSCAD1,PSCAD2, 82 : & COEF1,COEF2,COEF3,COEF4,SCN1X,SCN1Y,SCN1Z,VX,VY,VZ,COEF1X,COEF2X, 83 : & COEF1Y,COEF2Y,COEF1Z,COEF2Z,CX,CY,CZ,ANCX,ANCY,ANCZ, 84 : & DIFFX,DIFFY,DIFFZ,XLONGG,XLONGD 85 : & VALD,VALG,COEF,GX,GY,GZ,UN,EXPR1,EXPR2 86 : 87 : c REAL*8 VECXG1(2),VECYG1(2) 88 : c REAL*8 VECXG2(2),VECYG2(2) 89 : c REAL*8 VECXD1(2),VECYD1(2) 90 : c REAL*8 VECXD2(2),VECYD2(2) 91 : 92 : REAL*8 VECXG(4,4),VECYG(4,4),VECZG(4,4) 93 : REAL*8 VECXD(4,4),VECYD(4,4),VECZD(4,4) 94 : REAL*8 VOLUG(4),SURFAGX(4),SURFAGY(4),SURFAGZ(4),COEFG(4) 95 : REAL*8 VOLUD(4),SURFADX(4),SURFADY(4),SURFADZ(4),COEFD(4) 96 : REAL*8 PX(3,4),PY(3,4),PZ(3,4),TRGAUX(4) 97 : REAL*8 XPRO,TRACE 98 : INTEGER NGS(4),NLS(4),XS(4),YS(4),ZS(4) 99 : INTEGER NLOCFG(4,4),NLOCFD(4,4) 100 : INTEGER IGNS(4) 101 : REAL*8 EPS 102 : INTEGER ICRIT 103 : CHARACTER*(4) NOMCOM(3) 104 : CHARACTER*8 TYPE 105 : C 106 : DATA NOMCOM /'P1DX','P1DY','P1DZ'/ 107 : 108 : INTEGER NDIM 109 : SEGMENT MMAT1 110 : REAL*8 PM(NDIM,NDIM),PM1(NDIM,NDIM),XSOL(NDIM) 111 : INTEGER IC(NDIM) 112 : ENDSEGMENT 113 : 114 : INTEGER K1,K2 115 : SEGMENT INDICE 116 : INTEGER NUME(K1,K2) 117 : ENDSEGMENT 118 : POINTEUR IND.INDICE,IND2.INDICE,IND22.INDICE 119 : 120 : SEGMENT MATRICE 121 : REAL*8 MAT(K1,K2) 122 : ENDSEGMENT 123 : POINTEUR VAL1.MATRICE,VAL2.MATRICE,SCMB.MATRICE 124 : 125 : INTEGER K3 126 : SEGMENT POINT2 127 : INTEGER POINT(K3) 128 : ENDSEGMENT 129 : POINTEUR IPO2.POINT2 130 : 131 : SEGMENT MATRICE2 132 : REAL*8 MAT2(K1,K2) 133 : ENDSEGMENT 134 : POINTEUR MATR1.MATRICE2,MATR2.MATRICE2 135 : 136 : SEGMENT POINT3 137 : INTEGER POINT33(K3) 138 : ENDSEGMENT 139 : POINTEUR IPO3.POINT3 140 : 141 : SEGMENT INDICE3 142 : INTEGER NU(K1,K2) 143 : ENDSEGMENT 144 : POINTEUR INDIC.INDICE3 145 : 146 : SEGMENT REP 147 : INTEGER ID(K3) 148 : ENDSEGMENT 149 : POINTEUR TAB.REP,INDLI.REP 150 : 151 : INTEGER K5 152 : SEGMENT NBFAC 153 : INTEGER NBFACEL(K5) 154 : INTEGER IMELEM(K5) 155 : ENDSEGMENT 156 : 157 : INTEGER K6 158 : SEGMENT NBCOT 159 : INTEGER NBCOTE(K6) 160 : INTEGER IMECOTE(K6) 161 : ENDSEGMENT 162 : 163 : INTEGER K7,K8 164 : SEGMENT MISZERO 165 : INTEGER TABL(K7) 166 : INTEGER TABL2(K7) 167 : INTEGER TABL1(K8),IPOS(K8),ICOURANT(K8) 168 : REAL*8 XMAX(K7) 169 : ENDSEGMENT 170 : POINTEUR ITAB.MISZERO 171 : 172 : 173 : c WRITE(6,*) 'COUCOU NORV4' 174 : 175 : 176 : NMAI1 = 0 177 : NMAI2 = 0 178 : NMAI3 = 0 179 : NMAI4 = 0 180 : NBSO = MAX(1,MELTFA.LISOUS(/1)) 181 : c WRITE(6,*) 'NBSO MAILLE= ',NBSO 182 : c WRITE(6,*) 'MELTFA= ',MELTFA 183 : IELTFA = MELTFA 184 : IF (NBSO.EQ.1) THEN 185 : K5 = MELTFA.NUM(/2) 186 : ELSEIF (NBSO.EQ.2) THEN 187 : IPT1 = MELTFA.LISOUS(1) 188 : SEGACT IPT1 189 : N1 = IPT1.NUM(/2) 190 : NMAI1 = N1 191 : SEGDES IPT1 192 : IPT2 = MELTFA.LISOUS(2) 193 : SEGACT IPT2 194 : N2 = IPT2.NUM(/2) 195 : NMAI2 = N2 196 : SEGDES IPT2 197 : K5 = N1 + N2 198 : ELSEIF (NBSO.EQ.3) THEN 199 : IPT1 = MELTFA.LISOUS(1) 200 : SEGACT IPT1 201 : N1 = IPT1.NUM(/2) 202 : NMAI1 = N1 203 : SEGDES IPT1 204 : c WRITE(6,*) 'N1= ',N1 205 : IPT2 = MELTFA.LISOUS(2) 206 : SEGACT IPT2 207 : N2 = IPT2.NUM(/2) 208 : NMAI2 = N2 209 : SEGDES IPT2 210 : c WRITE(6,*) 'N2= ',N2 211 : IPT3 = MELTFA.LISOUS(3) 212 : SEGACT IPT3 213 : N3 = IPT3.NUM(/2) 214 : NMAI3 = N3 215 : c WRITE(6,*) 'N3= ',N3 216 : SEGDES IPT3 217 : K5 = N1 + N2 + N3 218 : ELSEIF (NBSO.EQ.4) THEN 219 : IPT1 = MELTFA.LISOUS(1) 220 : SEGACT IPT1 221 : N1 = IPT1.NUM(/2) 222 : NMAI1 = N1 223 : SEGDES IPT1 224 : IPT2 = MELTFA.LISOUS(2) 225 : SEGACT IPT2 226 : N2 = IPT2.NUM(/2) 227 : NMAI2 = N2 228 : SEGDES IPT2 229 : IPT3 = MELTFA.LISOUS(3) 230 : SEGACT IPT3 231 : N3 = IPT3.NUM(/2) 232 : NMAI3 = N3 233 : SEGDES IPT3 234 : IPT4 = MELTFA.LISOUS(4) 235 : SEGACT IPT4 236 : N4 = IPT4.NUM(/2) 237 : NMAI4 = N4 238 : SEGDES IPT4 239 : K5 = N1 + N2 + N3 + N4 240 : ENDIF 241 : c WRITE(6,*) 'K5= ',K5 242 : 243 : 244 : 245 : IF (NBSO.EQ.1) THEN 246 : DO I = 1,K5 247 : NTYPE = MELTFA.ITYPEL 248 : c WRITE(6,*) 'NTYPE= ',NTYPE 249 : IF (NTYPE .EQ. 16) THEN 250 : NBFACEL(I) = 6 251 : IMELEM(I) = MELTFA 252 : ELSEIF (NTYPE .EQ. 25) THEN 253 : NBFACEL(I) = 5 254 : IMELEM(I) = MELTFA 255 : ELSEIF (NTYPE .EQ. 23) THEN 256 : NBFACEL(I) = 4 257 : IMELEM(I) = MELTFA 258 : ELSEIF (NTYPE .EQ. 9) THEN 259 : NBFACEL(I) = 5 260 : IMELEM(I) = MELTFA 261 : ENDIF 262 : c SEGDES MELTFA 263 : ENDDO 264 : ELSEIF (NBSO.EQ.2) THEN 265 : IPT1 = MELTFA.LISOUS(1) 266 : SEGACT IPT1 267 : IPT2 = MELTFA.LISOUS(2) 268 : SEGACT IPT2 269 : DO I = 1,K5 270 : N1 = IPT1.NUM(/2) 271 : IF (I.LE.N1) THEN 272 : NTYPE = IPT1.ITYPEL 273 : IF (NTYPE .EQ. 16) THEN 274 : NBFACEL(I) = 6 275 : IMELEM(I) = IPT1 276 : ELSEIF (NTYPE .EQ. 25) THEN 277 : NBFACEL(I) = 5 278 : IMELEM(I) = IPT1 279 : ELSEIF (NTYPE .EQ. 23) THEN 280 : NBFACEL(I) = 4 281 : IMELEM(I) = IPT1 282 : ELSEIF (NTYPE .EQ. 9) THEN 283 : NBFACEL(I) = 5 284 : IMELEM(I) = IPT1 285 : ENDIF 286 : ELSE 287 : NTYPE = IPT2.ITYPEL 288 : IF (NTYPE .EQ. 16) THEN 289 : NBFACEL(I) = 6 290 : IMELEM(I) = IPT2 291 : ELSEIF (NTYPE .EQ. 25) THEN 292 : NBFACEL(I) = 5 293 : IMELEM(I) = IPT2 294 : ELSEIF (NTYPE .EQ. 23) THEN 295 : NBFACEL(I) = 4 296 : IMELEM(I) = IPT2 297 : ELSEIF (NTYPE .EQ. 9) THEN 298 : NBFACEL(I) = 5 299 : IMELEM(I) = IPT2 300 : ENDIF 301 : ENDIF 302 : ENDDO 303 : ELSEIF (NBSO.EQ.3) THEN 304 : c WRITE(6,*) 'COUCOU' 305 : IPT1 = MELTFA.LISOUS(1) 306 : SEGACT IPT1 307 : NTYPE = IPT1.ITYPEL 308 : c WRITE(6,*) 'NTYPE= ',IPT1.ITYPEL 309 : IPT2 = MELTFA.LISOUS(2) 310 : SEGACT IPT2 311 : NTYPE = IPT2.ITYPEL 312 : c WRITE(6,*) 'NTYPE= ',IPT2.ITYPEL 313 : IPT3 = MELTFA.LISOUS(3) 314 : SEGACT IPT3 315 : NTYPE = IPT3.ITYPEL 316 : c WRITE(6,*) 'NTYPE= ',IPT3.ITYPEL 317 : N1 = IPT1.NUM(/2) 318 : N2 = IPT2.NUM(/2) 319 : N3 = IPT3.NUM(/2) 320 : DO I = 1,K5 321 : IF (I.LE.N1) THEN 322 : NTYPE = IPT1.ITYPEL 323 : IF (NTYPE .EQ. 16) THEN 324 : NBFACEL(I) = 6 325 : IMELEM(I) = IPT1 326 : ELSEIF (NTYPE .EQ. 25) THEN 327 : NBFACEL(I) = 5 328 : IMELEM(I) = IPT1 329 : ELSEIF (NTYPE .EQ. 23) THEN 330 : NBFACEL(I) = 4 331 : IMELEM(I) = IPT1 332 : ELSEIF (NTYPE .EQ. 9) THEN 333 : NBFACEL(I) = 5 334 : IMELEM(I) = IPT1 335 : ENDIF 336 : ELSEIF (I.LE.(N1+N2)) THEN 337 : NTYPE = IPT2.ITYPEL 338 : IF (NTYPE .EQ. 16) THEN 339 : NBFACEL(I) = 6 340 : IMELEM(I) = IPT2 341 : ELSEIF (NTYPE .EQ. 25) THEN 342 : NBFACEL(I) = 5 343 : IMELEM(I) = IPT2 344 : ELSEIF (NTYPE .EQ. 23) THEN 345 : NBFACEL(I) = 4 346 : IMELEM(I) = IPT2 347 : ELSEIF (NTYPE .EQ. 9) THEN 348 : NBFACEL(I) = 5 349 : IMELEM(I) = IPT2 350 : ENDIF 351 : ELSE 352 : NTYPE = IPT3.ITYPEL 353 : IF (NTYPE .EQ. 16) THEN 354 : NBFACEL(I) = 6 355 : IMELEM(I) = IPT3 356 : ELSEIF (NTYPE .EQ. 25) THEN 357 : NBFACEL(I) = 5 358 : IMELEM(I) = IPT3 359 : ELSEIF (NTYPE .EQ. 23) THEN 360 : NBFACEL(I) = 4 361 : IMELEM(I) = IPT3 362 : ELSEIF (NTYPE .EQ. 9) THEN 363 : NBFACEL(I) = 5 364 : IMELEM(I) = IPT3 365 : ENDIF 366 : ENDIF 367 : ENDDO 368 : ELSEIF (NBSO.EQ.4) THEN 369 : IPT1 = MELTFA.LISOUS(1) 370 : SEGACT IPT1 371 : NTYPE = IPT1.ITYPEL 372 : c WRITE(6,*) 'NTYPE= ',IPT1.ITYPEL 373 : IPT2 = MELTFA.LISOUS(2) 374 : SEGACT IPT2 375 : NTYPE = IPT2.ITYPEL 376 : c WRITE(6,*) 'NTYPE= ',IPT2.ITYPEL 377 : IPT3 = MELTFA.LISOUS(3) 378 : SEGACT IPT3 379 : NTYPE = IPT3.ITYPEL 380 : c WRITE(6,*) 'NTYPE= ',IPT3.ITYPEL 381 : IPT4 = MELTFA.LISOUS(4) 382 : SEGACT IPT4 383 : NTYPE = IPT4.ITYPEL 384 : c WRITE(6,*) 'NTYPE= ',IPT4.ITYPEL 385 : N1 = IPT1.NUM(/2) 386 : N2 = IPT2.NUM(/2) 387 : N3 = IPT3.NUM(/2) 388 : N4 = IPT4.NUM(/2) 389 : DO I = 1,K5 390 : IF (I.LE.N1) THEN 391 : NTYPE = IPT1.ITYPEL 392 : IF (NTYPE .EQ. 16) THEN 393 : NBFACEL(I) = 6 394 : IMELEM(I) = IPT1 395 : ELSEIF (NTYPE .EQ. 25) THEN 396 : NBFACEL(I) = 5 397 : IMELEM(I) = IPT1 398 : ELSEIF (NTYPE .EQ. 23) THEN 399 : NBFACEL(I) = 4 400 : IMELEM(I) = IPT1 401 : ELSEIF (NTYPE .EQ. 9) THEN 402 : NBFACEL(I) = 5 403 : IMELEM(I) = IPT1 404 : ENDIF 405 : ELSEIF (I.LE.(N1+N2)) THEN 406 : NTYPE = IPT2.ITYPEL 407 : IF (NTYPE .EQ. 16) THEN 408 : NBFACEL(I) = 6 409 : IMELEM(I) = IPT2 410 : ELSEIF (NTYPE .EQ. 25) THEN 411 : NBFACEL(I) = 5 412 : IMELEM(I) = IPT2 413 : ELSEIF (NTYPE .EQ. 23) THEN 414 : NBFACEL(I) = 4 415 : IMELEM(I) = IPT2 416 : ELSEIF (NTYPE .EQ. 9) THEN 417 : NBFACEL(I) = 5 418 : IMELEM(I) = IPT2 419 : ENDIF 420 : ELSEIF (I.LE.(N1+N2+N3)) THEN 421 : NTYPE = IPT3.ITYPEL 422 : IF (NTYPE .EQ. 16) THEN 423 : NBFACEL(I) = 6 424 : IMELEM(I) = IPT3 425 : ELSEIF (NTYPE .EQ. 25) THEN 426 : NBFACEL(I) = 5 427 : IMELEM(I) = IPT3 428 : ELSEIF (NTYPE .EQ. 23) THEN 429 : NBFACEL(I) = 4 430 : IMELEM(I) = IPT3 431 : ELSEIF (NTYPE .EQ. 9) THEN 432 : NBFACEL(I) = 5 433 : IMELEM(I) = IPT3 434 : ENDIF 435 : ELSE 436 : NTYPE = IPT4.ITYPEL 437 : IF (NTYPE .EQ. 16) THEN 438 : NBFACEL(I) = 6 439 : IMELEM(I) = IPT4 440 : ELSEIF (NTYPE .EQ. 25) THEN 441 : NBFACEL(I) = 5 442 : IMELEM(I) = IPT4 443 : ELSEIF (NTYPE .EQ. 23) THEN 444 : NBFACEL(I) = 4 445 : IMELEM(I) = IPT4 446 : ELSEIF (NTYPE .EQ. 9) THEN 447 : NBFACEL(I) = 5 448 : IMELEM(I) = IPT4 449 : ENDIF 450 : ENDIF 451 : ENDDO 452 : ENDIF 453 : 454 : c CAS OU LES FACES SONT DES TRIANGLES OU DES FACES 455 : MAUX = MELEFP 456 : NFAI1 = 0 457 : NBSOF = MAX(1,MELEFP.LISOUS(/1)) 458 : c WRITE(6,*) 'NBSOF FACE= ',NBSOF 459 : IELTFA = MELTFA 460 : IF (NBSOF.EQ.1) THEN 461 : K6 = MELEFP.NUM(/2) 462 : ELSEIF (NBSOF.EQ.2) THEN 463 : IPT5 = MELEFP.LISOUS(1) 464 : SEGACT IPT5 465 : N1 = IPT5.NUM(/2) 466 : NFAI1 = N1 467 : SEGDES IPT5 468 : IPT6 = MELEFP.LISOUS(2) 469 : SEGACT IPT6 470 : N2 = IPT6.NUM(/2) 471 : NFAI2 = N2 472 : SEGDES IPT6 473 : K6 = N1 + N2 474 : ENDIF 475 : 476 : C ON EST ICI 477 : IF (NBSOF.EQ.1) THEN 478 : DO I = 1,K6 479 : NTYPE = MELEFP.ITYPEL 480 : c WRITE(6,*) 'NTYPE= ',NTYPE 481 : IF (NTYPE .EQ. 5) THEN 482 : NBCOTE(I) = 3 483 : IMECOTE(I) = MELEFP 484 : ELSE 485 : NBCOTE(I) = 4 486 : IMECOTE(I) = MELEFP 487 : ENDIF 488 : c SEGDES MELTFA 489 : ENDDO 490 : ELSEIF (NBSOF.EQ.2) THEN 491 : c WRITE(6,*) 'POINT2' 492 : IPT5 = MELEFP.LISOUS(1) 493 : SEGACT IPT5 494 : IPT6 = MELEFP.LISOUS(2) 495 : SEGACT IPT6 496 : c WRITE(6,*) 'IPT5= ',IPT5.ITYPEL 497 : c WRITE(6,*) 'IPT6= ',IPT6.ITYPEL 498 : DO I = 1,K6 499 : N1 = IPT5.NUM(/2) 500 : C MISE A JOUR DE MLEFA.LECT 501 : IF (I.LE.N1) THEN 502 : N0 = IPT5.NUM(/1) 503 : NGFAUX = IPT5.NUM(N0,I) 504 : MLEFA2.LECT(NGFAUX) = I 505 : c WRITE(6,*) 'NGFAUX = ',NGFAUX, 506 : c & 'MLEFA2=',MLEFA2.LECT(NGFAUX) 507 : IF (IPT5.ITYPEL .EQ. 5) THEN 508 : NBCOTE(I) = 3 509 : IMECOTE(I) = IPT5 510 : ELSE 511 : NBCOTE(I) = 4 512 : IMECOTE(I) = IPT5 513 : ENDIF 514 : c SEGDES IPT5 515 : ELSE 516 : N0 = IPT6.NUM(/1) 517 : NGFAUX = IPT6.NUM(N0,I-NFAI1) 518 : MLEFA2.LECT(NGFAUX) = I 519 : IF (IPT6.ITYPEL .EQ. 5) THEN 520 : NBCOTE(I) = 3 521 : IMECOTE(I) = IPT6 522 : ELSE 523 : NBCOTE(I) = 4 524 : IMECOTE(I) = IPT6 525 : ENDIF 526 : c SEGDES IPT6 527 : ENDIF 528 : 529 : ENDDO 530 : ENDIF 531 : 532 : 533 : 534 : NFAC=MELEFL.NUM(/2) 535 : 536 : 537 : c ON MAJORE SUPERIEUREMENT NBNN : ON LE REAJUSTERA PAR LA SUITE 538 : NCON = ((2*NBMAX)) + 1 539 : c NCON = ((3*NBMAX)/2) + 1 540 : NBNN = NCON 541 : NESSAI = NCON 542 : NBNN = NESSAI 543 : c WRITE(6,*) 'NBMAX= ',NBMAX 544 : c WRITE(6,*) 'NBNN= ',NBNN 545 : c WRITE(6,*) 'NFAC= ',NFAC 546 : 547 : 548 : 549 : C DEFINITION DES PARAMETRES DU CHAMELEM DES COEFFICIENTS 550 : 551 : c INITIALISATION DU CHAMELEM 552 : N1=1 553 : N2=1 554 : c WRITE(6,*) 'N2= ',N2 555 : N3=6 556 : L1=8 557 : SEGINI MCHELM 558 : ICOEFF = MCHELM 559 : MCHELM.TITCHE='Gradient' 560 : MCHELM.IFOCHE=IFOUR 561 : C 562 : ISOUS=0 563 : NBSOUS=0 564 : NBREF=0 565 : NBELEM = NFAC 566 : ISOUS=ISOUS+1 567 : SEGINI MELEME 568 : C ITYPEL=32 -> 'POLY' 569 : ITYPEL=32 570 : MCHELM.IMACHE(ISOUS)=MELEME 571 : SEGINI MCHAML 572 : MCHELM.ICHAML(ISOUS)=MCHAML 573 : MCHAML.NOMCHE(1)='SCAL' 574 : MCHAML.TYPCHE(1)='REAL*8 ' 575 : N1PTEL=NESSAI 576 : N1EL=NBELEM 577 : N2PTEL=0 578 : N2EL=0 579 : SEGINI MELVA1 580 : MCHAML.IELVAL(1)=MELVA1 581 : 582 : 583 : 584 : 585 : CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 586 : C CALCUL DE LA VITESSE EN CHAQUE FACE 587 : AA = 0.0 588 : BB = 0.0 589 : C INDICE QUI COMPTE LES COEFFICIENTS POUR CHAQUE FACE 590 : NAUX2 = 0 591 : NMOY = 0 592 : DO NLCF= 1, NFAC, 1 593 : NCON = 1 594 : 595 : NGCF=MELEFL.NUM(2,NLCF) 596 : NGCG=MELEFL.NUM(1,NLCF) 597 : NGCD=MELEFL.NUM(3,NLCF) 598 : NLCG=MLECEN.LECT(NGCG) 599 : NLCD=MLECEN.LECT(NGCD) 600 : 601 : SCNX=MPONOR.VPOCHA(NLCF,1) 602 : SCNY=MPONOR.VPOCHA(NLCF,2) 603 : SCNZ=MPONOR.VPOCHA(NLCF,3) 604 : SCN1X = SCNX 605 : SCN1Y = SCNY 606 : SCN1Z = SCNZ 607 : SURF=0.5D0*MPOSUR.VPOCHA(NLCF,1) 608 : SCNX=SCNX*SURF 609 : SCNY=SCNY*SURF 610 : SCNZ=SCNZ*SURF 611 : 612 : C 4=IDIM+1 613 : ICELL=(4*(NGCG -1))+1 614 : XG=MCOORD.XCOOR(ICELL) 615 : YG=MCOORD.XCOOR(ICELL+1) 616 : ZG=MCOORD.XCOOR(ICELL+2) 617 : ICELL=(4*(NGCD -1))+1 618 : XD=MCOORD.XCOOR(ICELL) 619 : YD=MCOORD.XCOOR(ICELL+1) 620 : ZD=MCOORD.XCOOR(ICELL+2) 621 : ICELL=(4*(NGCF -1))+1 622 : XF=MCOORD.XCOOR(ICELL) 623 : YF=MCOORD.XCOOR(ICELL+1) 624 : ZF=MCOORD.XCOOR(ICELL+2) 625 : 626 : 627 : MELEME.NUM(1,NLCF)=NGCF 628 : MELVA1.VELCHE(1,NLCF)=0.0D0 629 : 630 : DO J= 1,NBNN 631 : MELVA1.VELCHE(J,NLCF) = 0.0D0 632 : ENDDO 633 : 634 : C MISE A ZERO DE NLOC 635 : DO JA=1,4 636 : DO IA=1,3 637 : NLOCFG(IA,JA) = 0 638 : NLOCFD(IA,JA) = 0 639 : ENDDO 640 : ENDDO 641 : 642 : 643 : MELTFA = IMELEM(NLCG) 644 : NBF = NBFACEL(NLCG) 645 : IF (NLCG.LE.NMAI1) THEN 646 : NGAUX = NLCG 647 : ELSEIF ((NLCG.GT.NMAI1).AND.(NLCG.LE.(NMAI1+NMAI2))) THEN 648 : NGAUX = NLCG - NMAI1 649 : ELSEIF ((NLCG.GT.(NMAI1+NMAI2)).AND. 650 : & (NLCG.LE.(NMAI1+NMAI2+NMAI3))) THEN 651 : NGAUX = NLCG - (NMAI1+NMAI2) 652 : ELSEIF (NLCG.GT.(NMAI1+NMAI2+NMAI3)) THEN 653 : NGAUX = NLCG - (NMAI1+NMAI2+NMAI3) 654 : ENDIF 655 : 656 : 657 : 658 : c ON TIENT COMPTE DU CHANGEMENT DE NUMEROTATION 659 : 660 : 661 : 662 : NLCF1 = MLEFA2.LECT(NGCF) 663 : NBNO = NBCOTE(NLCF1) 664 : MELEFP = IMECOTE(NLCF1) 665 : IF (NLCF1.GT.NFAI1) THEN 666 : NLCF1AUX = NLCF1 - NFAI1 667 : ELSE 668 : NLCF1AUX = NLCF1 669 : ENDIF 670 : 671 : 672 : DO IA=1,NBNO 673 : NGS(IA) = MELEFP.NUM(IA,NLCF1AUX) 674 : NLS(IA) = MLESOM.LECT(NGS(IA)) 675 : ICELL=(4*(NGS(IA) -1))+1 676 : XS(IA)=nint(MCOORD.XCOOR(ICELL)) 677 : YS(IA)=nint(MCOORD.XCOOR(ICELL+1)) 678 : ZS(IA)=nint(MCOORD.XCOOR(ICELL+2)) 679 : ENDDO 680 : 681 : 682 : 683 : 684 : C ON REPERE LES VECTEURS PRINCIPAUX DE LA BASE 685 : 686 : DO JA = 1,NBNO 687 : NGS(JA) = MELEFP.NUM(JA,NLCF1AUX) 688 : c WRITE(6,*) 'NGAUX= ',NGAUX,'JA= ',JA,'NGS= ',NGS(JA) 689 : c WRITE(6,*) 'NGCF= ',NGCF,'NLCF= ',NLCF 690 : 691 : ICOUR = 0 692 : DO J = 1,NBF 693 : N1 = MELTFA.NUM(J,NGAUX) 694 : NL1 = MLEFA2.LECT(N1) 695 : NBNO2 = NBCOTE(NL1) 696 : MELEP2 = IMECOTE(NL1) 697 : IF (NL1.GT.NFAI1) THEN 698 : NL1AUX = NL1 - NFAI1 699 : ELSE 700 : NL1AUX = NL1 701 : ENDIF 702 : c WRITE(6,*) 'N1= ',N1,'NL1= ',NL1,'NL1AUX= ',NL1AUX 703 : 704 : 705 : DO IA =1,NBNO2 706 : NSOM1 = MELEP2.NUM(IA,NL1AUX) 707 : c WRITE(6,*) 'NBNO2= ',NBNO2,'IA= ',IA,'NSOM1= ',NSOM1 708 : IF (NSOM1.EQ.NGS(JA)) THEN 709 : 710 : ICELL=(4*(N1 -1))+1 711 : XF=MCOORD.XCOOR(ICELL) 712 : YF=MCOORD.XCOOR(ICELL+1) 713 : ZF=MCOORD.XCOOR(ICELL+2) 714 : 715 : ICOUR = ICOUR + 1 716 : VECXG(ICOUR,JA) = (XF - XG) 717 : VECYG(ICOUR,JA) = (YF - YG) 718 : VECZG(ICOUR,JA) = (ZF - ZG) 719 : NLOCFG(ICOUR,JA) = N1 720 : C ON PERMUTE 721 : IF (N1.EQ.NGCF) THEN 722 : NAUX = NLOCFG(1,JA) 723 : VXAU = VECXG(1,JA) 724 : VYAU = VECYG(1,JA) 725 : VZAU = VECZG(1,JA) 726 : VECXG(1,JA) = (XF - XG) 727 : VECYG(1,JA) = (YF - YG) 728 : VECZG(1,JA) = (ZF - ZG) 729 : NLOCFG(1,JA) = N1 730 : VECXG(ICOUR,JA) = VXAU 731 : VECYG(ICOUR,JA) = VYAU 732 : VECZG(ICOUR,JA) = VZAU 733 : NLOCFG(ICOUR,JA) = NAUX 734 : ENDIF 735 : ENDIF 736 : ENDDO 737 : ENDDO 738 : ENDDO 739 : 740 : MELTFA = IMELEM(NLCD) 741 : NBF = NBFACEL(NLCD) 742 : IF (NLCD.LE.NMAI1) THEN 743 : NDAUX = NLCD 744 : ELSEIF ((NLCD.GT.NMAI1).AND.(NLCD.LE.(NMAI1+NMAI2))) THEN 745 : NDAUX = NLCD - NMAI1 746 : ELSEIF ((NLCD.GT.(NMAI1+NMAI2)).AND. 747 : & (NLCD.LE.(NMAI1+NMAI2+NMAI3))) THEN 748 : NDAUX = NLCD - (NMAI1+NMAI2) 749 : ELSEIF (NLCD.GT.(NMAI1+NMAI2+NMAI3)) THEN 750 : NDAUX = NLCD - (NMAI1+NMAI2+NMAI3) 751 : ENDIF 752 : 753 : 754 : C ON REPERE LES VECTEURS PRINCIPAUX DE LA BASE 755 : DO JA = 1,NBNO 756 : NGS(JA) = MELEFP.NUM(JA,NLCF1AUX) 757 : c WRITE(6,*) 'NDAUX= ',NDAUX,'JA= ',JA,'NGS= ',NGS(JA) 758 : c WRITE(6,*) 'NGCF= ',NGCF,'NLCF= ',NLCF 759 : ICOUR = 0 760 : DO J = 1,NBF 761 : N1 = MELTFA.NUM(J,NDAUX) 762 : NL1 = MLEFA2.LECT(N1) 763 : c WRITE(6,*) 'N1= ',N1,'NL1= ',NL1 764 : 765 : NBNO2 = NBCOTE(NL1) 766 : MELEP2 = IMECOTE(NL1) 767 : IF (NL1.GT.NFAI1) THEN 768 : NL1AUX = NL1 - NFAI1 769 : ELSE 770 : NL1AUX = NL1 771 : ENDIF 772 : 773 : 774 : DO IA =1,NBNO2 775 : NSOM1 = MELEP2.NUM(IA,NL1AUX) 776 : c WRITE(6,*) 'NBNO2= ',NBNO2,'IA= ',IA,'NSOM1= ',NSOM1 777 : IF (NSOM1.EQ.NGS(JA)) THEN 778 : 779 : ICELL=(4*(N1 -1))+1 780 : XF=MCOORD.XCOOR(ICELL) 781 : YF=MCOORD.XCOOR(ICELL+1) 782 : ZF=MCOORD.XCOOR(ICELL+2) 783 : 784 : ICOUR = ICOUR + 1 785 : VECXD(ICOUR,JA) = (XF - XD) 786 : VECYD(ICOUR,JA) = (YF - YD) 787 : VECZD(ICOUR,JA) = (ZF - ZD) 788 : NLOCFD(ICOUR,JA) = N1 789 : C ON PERMUTE 790 : IF (N1.EQ.NGCF) THEN 791 : NAUX = NLOCFD(1,JA) 792 : VXAU = VECXD(1,JA) 793 : VYAU = VECYD(1,JA) 794 : VZAU = VECZD(1,JA) 795 : VECXD(1,JA) = (XF - XD) 796 : VECYD(1,JA) = (YF - YD) 797 : VECZD(1,JA) = (ZF - ZD) 798 : NLOCFD(1,JA) = N1 799 : VECXD(ICOUR,JA) = VXAU 800 : VECYD(ICOUR,JA) = VYAU 801 : VECZD(ICOUR,JA) = VZAU 802 : NLOCFD(ICOUR,JA) = NAUX 803 : ENDIF 804 : ENDIF 805 : ENDDO 806 : ENDDO 807 : c WRITE(6,*) 'JA= ',JA 808 : c WRITE(6,*) 'ICOUR= ',ICOUR 809 : ENDDO 810 : 811 : MPOGRA.VPOCHA(NLCF,1) = 0.D0 812 : DO JA = 1,NBNO 813 : c WRITE(6,*) 'XPRO= ',XPRO 814 : 815 : 816 : DO KA = 1,ICOUR 817 : C PRODUIT MIXTES 818 : C PRODUIT VECTORIEL 819 : IF (KA.EQ.1) THEN 820 : PSCAGX = (VECYG(2,JA)*VECZG(3,JA)) - 821 : & (VECZG(2,JA)*VECYG(3,JA)) 822 : PSCAGY = (VECZG(2,JA)*VECXG(3,JA)) - 823 : & (VECXG(2,JA)*VECZG(3,JA)) 824 : PSCAGZ = (VECXG(2,JA)*VECYG(3,JA)) - 825 : & (VECYG(2,JA)*VECXG(3,JA)) 826 : VOLUG(JA) = (VECXG(1,JA)* PSCAGX) + 827 : & (VECYG(1,JA)* PSCAGY) + 828 : & (VECZG(1,JA)* PSCAGZ) 829 : SURFAGX(KA) = 0.5D0* PSCAGX 830 : SURFAGY(KA) = 0.5D0* PSCAGY 831 : SURFAGZ(KA) = 0.5D0* PSCAGZ 832 : IF ( VOLUG(JA).GT.0) THEN 833 : SURFAGX(KA) = - SURFAGX(KA) 834 : SURFAGY(KA) = - SURFAGY(KA) 835 : SURFAGZ(KA) = - SURFAGZ(KA) 836 : ENDIF 837 : VOLUG(JA) = 1.D0/6.D0*ABS(VOLUG(JA)) 838 : ENDIF 839 : 840 : IF (KA.EQ.2) THEN 841 : PSCAGX = (VECYG(3,JA)*VECZG(1,JA)) - 842 : & (VECZG(3,JA)*VECYG(1,JA)) 843 : PSCAGY = (VECZG(3,JA)*VECXG(1,JA)) - 844 : & (VECXG(3,JA)*VECZG(1,JA)) 845 : PSCAGZ = (VECXG(3,JA)*VECYG(1,JA)) - 846 : & (VECYG(3,JA)*VECXG(1,JA)) 847 : SURFAGX(KA) = 0.5D0* PSCAGX 848 : SURFAGY(KA) = 0.5D0* PSCAGY 849 : SURFAGZ(KA) = 0.5D0* PSCAGZ 850 : PSCA = (VECXG(2,JA)* PSCAGX) + (VECYG(2,JA)* PSCAGY) + 851 : & (VECZG(2,JA)* PSCAGZ) 852 : IF ( PSCA.GT.0) THEN 853 : SURFAGX(KA) = - SURFAGX(KA) 854 : SURFAGY(KA) = - SURFAGY(KA) 855 : SURFAGZ(KA) = - SURFAGZ(KA) 856 : ENDIF 857 : ENDIF 858 : 859 : 860 : IF (KA.EQ.3) THEN 861 : PSCAGX = (VECYG(1,JA)*VECZG(2,JA)) - 862 : & (VECZG(1,JA)*VECYG(2,JA)) 863 : PSCAGY = (VECZG(1,JA)*VECXG(2,JA)) - 864 : & (VECXG(1,JA)*VECZG(2,JA)) 865 : PSCAGZ = (VECXG(1,JA)*VECYG(2,JA)) - 866 : & (VECYG(1,JA)*VECXG(2,JA)) 867 : 868 : SURFAGX(KA) = 0.5D0* PSCAGX 869 : SURFAGY(KA) = 0.5D0* PSCAGY 870 : SURFAGZ(KA) = 0.5D0* PSCAGZ 871 : PSCA = (VECXG(3,JA)* PSCAGX) + (VECYG(3,JA)* PSCAGY) + 872 : & (VECZG(3,JA)* PSCAGZ) 873 : IF ( PSCA.GT.0) THEN 874 : SURFAGX(KA) = - SURFAGX(KA) 875 : SURFAGY(KA) = - SURFAGY(KA) 876 : SURFAGZ(KA) = - SURFAGZ(KA) 877 : ENDIF 878 : ENDIF 879 : 880 : c CALCUL DE MATRICE POUR LE NOEUD D INDICE NS1 881 : IF (ICHTE.EQ.0) THEN 882 : PX(KA,JA) = SURFAGX(KA)/(3.D0*VOLUG(JA)) 883 : PY(KA,JA) = SURFAGY(KA)/(3.D0*VOLUG(JA)) 884 : PZ(KA,JA) = SURFAGZ(KA)/(3.D0*VOLUG(JA)) 885 : 886 : ELSE 887 : IF (MPOTEN.VPOCHA(/2) .EQ.6) THEN 888 : C TENSEUR ANISOTROPE 889 : K11G = MPOTEN.VPOCHA(NLCG,1) 890 : K22G = MPOTEN.VPOCHA(NLCG,2) 891 : K33G = MPOTEN.VPOCHA(NLCG,3) 892 : K21G = MPOTEN.VPOCHA(NLCG,4) 893 : K31G = MPOTEN.VPOCHA(NLCG,5) 894 : K32G = MPOTEN.VPOCHA(NLCG,6) 895 : ELSEIF (MPOTEN.VPOCHA(/2) .EQ.1) THEN 896 : C TENSEUR DIAGONAL 897 : K11G = MPOTEN.VPOCHA(NLCG,1) 898 : K22G = K11G 899 : K33G = K11G 900 : K21G = 0.0D0 901 : K31G = 0.0D0 902 : K32G = 0.0D0 903 : ELSE 904 : WRITE(6,*) 'TENSEUR NON PREVU' 905 : STOP 906 : ENDIF 907 : 908 : PSCAGX = (K11G*SURFAGX(KA)) + (K21G*SURFAGY(KA)) + 909 : & (K31G*SURFAGZ(KA)) 910 : PSCAGY = (K21G*SURFAGX(KA)) + (K22G*SURFAGY(KA)) + 911 : & (K32G*SURFAGZ(KA)) 912 : PSCAGZ = (K31G*SURFAGX(KA)) + (K32G*SURFAGY(KA)) + 913 : & (K33G*SURFAGZ(KA)) 914 : PX(KA,JA) = PSCAGX/(3.D0*VOLUG(JA)) 915 : PY(KA,JA) = PSCAGY/(3.D0*VOLUG(JA)) 916 : PZ(KA,JA) = PSCAGZ/(3.D0*VOLUG(JA)) 917 : ENDIF 918 : ENDDO 919 : ENDDO 920 : 921 : XPRO = 1.D0/NBNO 922 : DO JA = 1,NBNO 923 : 924 : MARQ = 0 925 : DO I5 = 1,INDLI.ID(NLS(JA)) 926 : INDAUX = IND2.NUME(I5,NLS(JA)) 927 : IF (INDAUX.EQ.NLOCFG(2,JA)) THEN 928 : IAFF = I5 929 : IG2 = IAFF 930 : GOTO 62 931 : ENDIF 932 : ENDDO 933 : 62 CONTINUE 934 : 935 : TRG2 = SCMB.MAT(IAFF,NLS(JA)) 936 : c WRITE(6,*) 'TR ','IAFF= ',IAFF,TRG2 937 : 938 : DO I5 = 1,INDLI.ID(NLS(JA)) 939 : INDAUX = IND2.NUME(I5,NLS(JA)) 940 : IF (INDAUX.EQ.NLOCFG(3,JA)) THEN 941 : IAFF = I5 942 : IG3 = IAFF 943 : GOTO 629 944 : ENDIF 945 : ENDDO 946 : 629 CONTINUE 947 : TRG3 = SCMB.MAT(IAFF,NLS(JA)) 948 : c WRITE(6,*) 'TR ','IAFF= ',IAFF,TRG3 949 : 950 : MARQ = 0 951 : DO I5 = 1,INDLI.ID(NLS(JA)) 952 : INDAUX = IND2.NUME(I5,NLS(JA)) 953 : IF (INDAUX.EQ.NGCF) THEN 954 : IAFF = I5 955 : IG = IAFF 956 : GOTO 63 957 : ENDIF 958 : ENDDO 959 : 63 CONTINUE 960 : TRG = SCMB.MAT(IAFF,NLS(JA)) 961 : TRGAUX(JA) = TRG 962 : IGNS(JA) = IAFF 963 : 964 : c WRITE(6,*) 'TR ','IAFF= ',IAFF,TRG 965 : 966 : VAL = MPOCHP.VPOCHA(NLCG,1) 967 : VALD = MPOCHP.VPOCHA(NLCD,1) 968 : c ICI 969 : 970 : AUX = (XPRO*( 971 : & (PX(1,JA) * (VAL - TRG)) 972 : & + (PX(2,JA) * (VAL - TRG2)) 973 : & + (PX(3,JA) * (VAL - TRG3)))) 974 : AUY = (XPRO*( 975 : & (PY(1,JA) * (VAL - TRG)) 976 : & + (PY(2,JA) * (VAL - TRG2)) 977 : & + (PY(3,JA) * (VAL - TRG3)))) 978 : AUZ = (XPRO*( 979 : & (PZ(1,JA) * (VAL - TRG)) 980 : & + (PZ(2,JA) * (VAL - TRG2)) 981 : & + (PZ(3,JA) * (VAL - TRG3)))) 982 : MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) + 983 : & (AUX*SCN1X) + (AUY*SCN1Y) + 984 : & (AUZ*SCN1Z) 985 : 986 : c IF (NLCF.EQ.791) THEN 987 : c WRITE(6,*) 'NLCF= ',NLCF,'GR1= ',MPOGRA.VPOCHA(NLCF,1) 988 : c WRITE(6,*) 'NLCF= ',NLCF,'GR2= ',MPOGRA.VPOCHA(NLCF,2) 989 : c WRITE(6,*) 'NLCF= ',NLCF,'GR3= ',MPOGRA.VPOCHA(NLCF,3) 990 : c WRITE(6,*) 'PX= ',PX(1,JA),PX(2,JA),PX(3,JA) 991 : c WRITE(6,*) 'PY= ',PY(1,JA),PY(2,JA),PY(3,JA) 992 : c WRITE(6,*) 'PZ= ',PZ(1,JA),PZ(2,JA),PZ(3,JA) 993 : c WRITE(6,*) 'TR ',TRG,TRG2,TRG3,'VAL',VAL 994 : c WRITE(6,*) 'NLS= ',NLS(JA), 995 : c & 'IG= ',IG,'IG2= ',IG2,'IG3= ',IG3 996 : c WRITE(6,*) 'NGCF= ',NGCF 997 : c ENDIF 998 : 999 : 1000 : ITROUVE = 0 1001 : INDIC = IPO3.POINT33(NLS(JA)) 1002 : SEGACT INDIC *MOD 1003 : MATR1 = IPO2.POINT(NLS(JA)) 1004 : SEGACT MATR1 *MOD 1005 : c NLS1 = NLS(JA) 1006 : c DO IAUX = 1,INDLI.ID(NLS1) 1007 : c DO IAUX2 = 1,TAB.ID(NLS1) 1008 : c WRITE(6,*) 'NLS1= ',NLS1,'IAUX= ',IAUX ,'IAUX2= ', 1009 : c & IAUX2,'VAUX',MATR1.MAT2(IAUX,IAUX2) 1010 : c & ,'IND3= ',INDIC.NU(IAUX,IAUX2) 1011 : c ENDDO 1012 : c ENDDO 1013 : 1014 : DO ICOUR = 1,TAB.ID(NLS(JA)) 1015 : IA = ICOUR 1016 : J1 = INDIC.NU(IG,IA) 1017 : DO IAUX2 = 2,NCON 1018 : J2 = MELEME.NUM(IAUX2,NLCF) 1019 : IF (J1.EQ.J2) THEN 1020 : IAUX = IAUX2 1021 : ITROUVE = 1 1022 : GOTO 5119 1023 : ENDIF 1024 : ENDDO 1025 : C ON A RIEN TROUVE : ON INCREMENTE LE COMPTEUR 1026 : NCON = NCON + 1 1027 : IAUX = NCON 1028 : 5119 CONTINUE 1029 : 1030 : 1031 : CX = MATR1.MAT2(IG,IA) 1032 : CY = MATR1.MAT2(IG,IA) 1033 : CZ = MATR1.MAT2(IG,IA) 1034 : IF (J1.EQ.NGCG) THEN 1035 : CX = CX - 1.D0 1036 : CY = CY - 1.D0 1037 : CZ = CZ - 1.D0 1038 : INDGA = IAUX 1039 : ENDIF 1040 : IF (J1.EQ.NGCD) THEN 1041 : INDDR = IAUX 1042 : ENDIF 1043 : c MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - 1044 : c & (XPRO*PX(1,JA)*CX) 1045 : c MELVA2.VELCHE(IAUX,NLCF) = MELVA2.VELCHE(IAUX,NLCF) - 1046 : c & (XPRO*PY(1,JA)*CY) 1047 : c MELVA3.VELCHE(IAUX,NLCF) = MELVA3.VELCHE(IAUX,NLCF) - 1048 : c & (XPRO*PZ(1,JA)*CZ) 1049 : AUX = (XPRO*PX(1,JA)*CX*SCN1X) + 1050 : & (XPRO*PY(1,JA)*CY*SCN1Y) + 1051 : & (XPRO*PZ(1,JA)*CZ*SCN1Z) 1052 : MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - 1053 : & AUX 1054 : 1055 : INDAUX = INDIC.NU(IG,IA) 1056 : MELEME.NUM(IAUX,NLCF) = INDAUX 1057 : ENDDO 1058 : 1059 : 1060 : ITROUVE = 0 1061 : DO ICOUR = 1,TAB.ID(NLS(JA)) 1062 : IA = ICOUR 1063 : J1 = INDIC.NU(IG2,IA) 1064 : DO IAUX2 = 2,NCON 1065 : J2 = MELEME.NUM(IAUX2,NLCF) 1066 : IF (J1.EQ.J2) THEN 1067 : IAUX = IAUX2 1068 : ITROUVE = 1 1069 : GOTO 511 1070 : ENDIF 1071 : ENDDO 1072 : C ON A RIEN TROUVE : ON INCREMENTE LE COMPTEUR 1073 : NCON = NCON + 1 1074 : IAUX = NCON 1075 : 511 CONTINUE 1076 : 1077 : CX = MATR1.MAT2(IG2,IA) 1078 : CY = MATR1.MAT2(IG2,IA) 1079 : CZ = MATR1.MAT2(IG2,IA) 1080 : IF (J1.EQ.NGCG) THEN 1081 : CX = CX - 1.D0 1082 : CY = CY - 1.D0 1083 : CZ = CZ - 1.D0 1084 : ENDIF 1085 : c MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - 1086 : c & (XPRO*PX(2,JA)*CX) 1087 : c MELVA2.VELCHE(IAUX,NLCF) = MELVA2.VELCHE(IAUX,NLCF) - 1088 : c & (XPRO*PY(2,JA)*CY) 1089 : c MELVA3.VELCHE(IAUX,NLCF) = MELVA3.VELCHE(IAUX,NLCF) - 1090 : c & (XPRO*PZ(2,JA)*CZ) 1091 : AUX = (XPRO*PX(2,JA)*CX*SCN1X) + 1092 : & (XPRO*PY(2,JA)*CY*SCN1Y) + 1093 : & (XPRO*PZ(2,JA)*CZ*SCN1Z) 1094 : MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - 1095 : & AUX 1096 : INDAUX = INDIC.NU(IG2,IA) 1097 : MELEME.NUM(IAUX,NLCF) = INDAUX 1098 : ENDDO 1099 : 1100 : ITROUVE = 0 1101 : DO ICOUR = 1,TAB.ID(NLS(JA)) 1102 : IA = ICOUR 1103 : J1 = INDIC.NU(IG3,IA) 1104 : DO IAUX2 = 2,NCON 1105 : J2 = MELEME.NUM(IAUX2,NLCF) 1106 : IF (J1.EQ.J2) THEN 1107 : IAUX = IAUX2 1108 : ITROUVE = 1 1109 : GOTO 5118 1110 : ENDIF 1111 : ENDDO 1112 : C ON A RIEN TROUVE : ON INCREMENTE LE COMPTEUR 1113 : NCON = NCON + 1 1114 : IAUX = NCON 1115 : 5118 CONTINUE 1116 : 1117 : CX = MATR1.MAT2(IG3,IA) 1118 : CY = MATR1.MAT2(IG3,IA) 1119 : CZ = MATR1.MAT2(IG3,IA) 1120 : IF (J1.EQ.NGCG) THEN 1121 : CX = CX - 1.D0 1122 : CY = CY - 1.D0 1123 : CZ = CZ - 1.D0 1124 : ENDIF 1125 : c MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - 1126 : c & (XPRO*PX(3,JA)*CX) 1127 : c MELVA2.VELCHE(IAUX,NLCF) = MELVA2.VELCHE(IAUX,NLCF) - 1128 : c & (XPRO*PY(3,JA)*CY) 1129 : c MELVA3.VELCHE(IAUX,NLCF) = MELVA3.VELCHE(IAUX,NLCF) - 1130 : c & (XPRO*PZ(3,JA)*CZ) 1131 : AUX = (XPRO*PX(3,JA)*CX*SCN1X) + 1132 : & (XPRO*PY(3,JA)*CY*SCN1Y) + 1133 : & (XPRO*PZ(3,JA)*CZ*SCN1Z) 1134 : MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - 1135 : & AUX 1136 : INDAUX = INDIC.NU(IG3,IA) 1137 : MELEME.NUM(IAUX,NLCF) = INDAUX 1138 : ENDDO 1139 : 1140 : c WRITE(6,*) 'NLCF= ',NLCF,'JA= ',JA,'VOLUG(JA) = ',VOLUG(JA), 1141 : c & 'XPRO= ',XPRO 1142 : 1143 : c FIN DE LA BOUCLE SUR LES NOEUDS 1144 : SEGDES MATR1 1145 : SEGDES INDIC 1146 : ENDDO 1147 : 1148 : c MPOGRA.VPOCHA(NLCF,1) = 0.0D0 1149 : c ISUP = NCON 1150 : c DO J= ISUP+1,NBNN 1151 : c DO J= 2,NBNN 1152 : c MELVA1.VELCHE(J,NLCF) = 0.0D0 1153 : c MELEME.NUM(J,NLCF) = MELEME.NUM(ISUP,NLCF) 1154 : c ENDDO 1155 : ISUP = NCON 1156 : DO J= ISUP+1,NBNN 1157 : MELVA1.VELCHE(J,NLCF) = 0.0D0 1158 : MELEME.NUM(J,NLCF) = MELEME.NUM(ISUP,NLCF) 1159 : ENDDO 1160 : 1161 : CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 1162 : C ON RAJOUTE LE CONVECTIF 1163 : 1164 : C CALCUL PLUS PRECIS 1165 : 1166 : 1167 : c IF (NLCD.NE.NLCG) THEN 1168 : c MELVA1.VELCHE(INDDR,NLCF) = - MELVA1.VELCHE(INDGA,NLCF) 1169 : c ELSE 1170 : c DO J= 1,ISUP 1171 : c ICENT = MELEME.NUM(J,NLCF) 1172 : c ICEN = MLECEN.LECT(ICENT) 1173 : c DIFF = MELVA1.VELCHE(J,NLCF) + MELVA1.VELCHE(INDGA,NLCF) 1174 : c DIFF = ABS(DIFF) 1175 : c XX = ABS(MELVA1.VELCHE(INDGA,NLCF)) 1176 : c IF ((ICEN.EQ.0).AND.(DIFF.LT.(1e-5*XX))) THEN 1177 : c MELVA1.VELCHE(J,NLCF) = - MELVA1.VELCHE(INDGA,NLCF) 1178 : c ENDIF 1179 : c ENDDO 1180 : c ENDIF 1181 : 1182 : 1183 : IF (ICHCO.GT.0) THEN 1184 : C BOUCLE POUR CALCUER INDGA,INDDR 1185 : INDFR = 0 1186 : DO J= 1,ISUP 1187 : IF (MELEME.NUM(J,NLCF).EQ.NGCG) INDGA = J 1188 : IF (MELEME.NUM(J,NLCF).EQ.NGCD) INDDR = J 1189 : IF (MELEME.NUM(J,NLCF).EQ.NGCF) INDFR = J 1190 : ENDDO 1191 : 1192 : 1193 : 1194 : c WRITE(6,*) 'NGCF= ',NGCF 1195 : UN = MPOVCO.VPOCHA(NLCF,1) 1196 : c WRITE(6,*) 'UN= ',UN 1197 : C OPTION CENTRE 1198 : IF (IOP.EQ.2) THEN 1199 : IF (NLCD.NE.NLCG) THEN 1200 : VAL = 0.5D0*(MPOCHP.VPOCHA(NLCG,1) + 1201 : & MPOCHP.VPOCHA(NLCD,1))*UN 1202 : MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) - VAL 1203 : ELSE 1204 : C CONDITIONS AUX LIMITES : TRACE CALCULEE PAR LA DIFFUSION 1205 : 1206 : XPRO = 1.D0/NBNO 1207 : TRACE = 0.0D0 1208 : DO JA = 1,NBNO 1209 : TRACE =TRACE + (XPRO*TRGAUX(JA)) 1210 : ENDDO 1211 : c WRITE(6,*) 'NLCF= ',NLCF,'TRACE= ',TRACE 1212 : VAL = TRACE*UN 1213 : MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) - VAL 1214 : ENDIF 1215 : 1216 : C ON COMPLETE MELVA1 POUR LE CONVECTIF 1217 : IF (NLCD.NE.NLCG) THEN 1218 : MELVA1.VELCHE(INDGA,NLCF) = MELVA1.VELCHE(INDGA,NLCF) - 1219 : & (0.5D0*UN) 1220 : MELVA1.VELCHE(INDDR,NLCF) = MELVA1.VELCHE(INDDR,NLCF) - 1221 : & (0.5D0*UN) 1222 : C CONDITION AUX LIMITES : ON RAJOUTE LES DEPENDENCES DES TRACES 1223 : c POUR LES CONDITIONS MIXTES OU DE NEUMAN 1224 : 1225 : 1226 : ELSE 1227 : 1228 : NLFCL = MLENCL.LECT(NGCF) 1229 : 1230 : C ON RAJOUTE CECI POUR L OPTION GRADGEO 1231 : IF (NLFCL.NE.0) THEN 1232 : MELVA1.VELCHE(NCON+1,NLCF) = - UN 1233 : MELEME.NUM(NCON+1,NLCF) = NGCF 1234 : ENDIF 1235 : 1236 : IF (NLFCL.EQ.0) THEN 1237 : XPRO = 1.D0/NBNO 1238 : DO JA = 1,NBNO 1239 : INDIC = IPO3.POINT33(NLS(JA)) 1240 : SEGACT INDIC *MOD 1241 : MATR1 = IPO2.POINT(NLS(JA)) 1242 : SEGACT MATR1 *MOD 1243 : 1244 : DO ICOUR = 1,TAB.ID(NLS(JA)) 1245 : IA = ICOUR 1246 : J1 = INDIC.NU(IGNS(JA),IA) 1247 : DO IAUX2 = 2,NCON 1248 : J2 = MELEME.NUM(IAUX2,NLCF) 1249 : IF (J1.EQ.J2) THEN 1250 : IAUX = IAUX2 1251 : GOTO 5169 1252 : ENDIF 1253 : ENDDO 1254 : 5169 CONTINUE 1255 : 1256 : 1257 : CX = MATR1.MAT2(IGNS(JA),IA) 1258 : MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - 1259 : & (UN*CX*XPRO) 1260 : ENDDO 1261 : SEGDES INDIC *MOD 1262 : SEGDES MATR1 *MOD 1263 : ENDDO 1264 : ENDIF 1265 : ENDIF 1266 : 1267 : C OPTION UPWIND 1268 : ELSEIF (IOP.EQ.1) THEN 1269 : IF (NLCD.NE.NLCG) THEN 1270 : IF (UN.GT.0.0D0) THEN 1271 : VAL = (MPOCHP.VPOCHA(NLCG,1)*UN) 1272 : ELSE 1273 : VAL = (MPOCHP.VPOCHA(NLCD,1)*UN) 1274 : ENDIF 1275 : c WRITE(6,*) 'VAL= ',VAL 1276 : MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) - VAL 1277 : ELSE 1278 : C CONDITIONS AUX LIMITES : TRACE CALCULEE PAR LA DIFFUSION 1279 : c ANCIENNE VERRUE 1280 : IF (UN.GT.0.0D0) THEN 1281 : VAL = (MPOCHP.VPOCHA(NLCG,1)*UN) 1282 : ELSE 1283 : XPRO = 1.D0/NBNO 1284 : TRACE = 0.0D0 1285 : DO JA = 1,NBNO 1286 : TRACE =TRACE + (XPRO*TRGAUX(JA)) 1287 : ENDDO 1288 : 1289 : c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF 1290 : c WRITE(6,*) 'UN= ',UN,'TRACE=',TRACE 1291 : 1292 : 1293 : VAL = TRACE*UN 1294 : ENDIF 1295 : 1296 : 1297 : MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) - VAL 1298 : ENDIF 1299 : 1300 : C ON COMPLETE MELVA1 POUR LE CONVECTIF 1301 : IF (NLCD.NE.NLCG) THEN 1302 : c WRITE(6,*) 'UN= ',UN 1303 : IF (UN.GT.0.0D0) THEN 1304 : MELVA1.VELCHE(INDGA,NLCF) = MELVA1.VELCHE(INDGA,NLCF) - 1305 : & (UN) 1306 : ELSE 1307 : MELVA1.VELCHE(INDDR,NLCF) = MELVA1.VELCHE(INDDR,NLCF) - 1308 : & (UN) 1309 : ENDIF 1310 : C CONDITION AUX LIMITES : ON RAJOUTE LES DEPENDENCES DES TRACES 1311 : c POUR LES CONDITIONS MIXTES OU DE NEUMAN 1312 : 1313 : ELSE 1314 : c ANCIENNE VERRUE 1315 : IF (UN.GT.0.0D0) THEN 1316 : MELVA1.VELCHE(INDGA,NLCF) = MELVA1.VELCHE(INDGA,NLCF) - 1317 : & (UN) 1318 : ELSE 1319 : 1320 : NLFCL = MLENCL.LECT(NGCF) 1321 : 1322 : C ON RAJOUTE CECI POUR L OPTION GRADGEO 1323 : C IL Y A PROBABLEMENT REDONDANCE 1324 : IF (NLFCL.NE.0) THEN 1325 : MELVA1.VELCHE(NCON+1,NLCF) = - UN 1326 : MELEME.NUM(NCON+1,NLCF) = NGCF 1327 : ENDIF 1328 : * ENDIF 1329 : 1330 : IF (NLFCL.EQ.0) THEN 1331 : XPRO = 1.D0/NBNO 1332 : DO JA = 1,NBNO 1333 : INDIC = IPO3.POINT33(NLS(JA)) 1334 : SEGACT INDIC *MOD 1335 : MATR1 = IPO2.POINT(NLS(JA)) 1336 : SEGACT MATR1 *MOD 1337 : 1338 : DO ICOUR = 1,TAB.ID(NLS(JA)) 1339 : IA = ICOUR 1340 : J1 = INDIC.NU(IGNS(JA),IA) 1341 : DO IAUX2 = 2,NCON 1342 : J2 = MELEME.NUM(IAUX2,NLCF) 1343 : IF (J1.EQ.J2) THEN 1344 : IAUX = IAUX2 1345 : GOTO 5129 1346 : ENDIF 1347 : ENDDO 1348 : 5129 CONTINUE 1349 : 1350 : 1351 : CX = MATR1.MAT2(IGNS(JA),IA) 1352 : MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - 1353 : & (UN*CX*XPRO) 1354 : ENDDO 1355 : 1356 : SEGDES INDIC *MOD 1357 : SEGDES MATR1 *MOD 1358 : ENDDO 1359 : ENDIF 1360 : ENDIF 1361 : 1362 : 1363 : ENDIF 1364 : 1365 : 1366 : C OPTION UPWICENT 1367 : ELSEIF (IOP.EQ.3) THEN 1368 : IF (NLCD.NE.NLCG) THEN 1369 : c CALCUL DE THETA 1370 : c DANS MELVA1 : INDDR EST LA NUM de NLCD 1371 : C COEFD = MELVA1(INDDR,NLCF) 1372 : COEFDD = MELVA1.VELCHE(INDDR,NLCF) 1373 : 1374 : 1375 : THETA = 0.5D0 1376 : EXPR1 = - (THETA*UN) - COEFDD 1377 : EXPR2 = - ((1.D0 - THETA)*UN) + COEFDD 1378 : c WRITE(6,*) 'EXPR1',EXPR1 1379 : c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF, 1380 : c & 'COEFDD=',COEFDD,'THETA= ',THETA,'UN= ',UN 1381 : c WRITE(6,*) 'EXPR2',EXPR2 1382 : c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF, 1383 : c & 'COEFDD=',COEFDD,'THETA= ',THETA,'UN= ',UN 1384 : IF (EXPR1.GT.0.0D0) THEN 1385 : IF (ABS(UN) .GT. 1e-20) THEN 1386 : THETA = - (0.5D0*COEFDD / UN) 1387 : ENDIF 1388 : ENDIF 1389 : IF (EXPR2.LT.0.0D0) THEN 1390 : IF (ABS(UN) .GT. 1e-20) THEN 1391 : THETA = 1.D0 - (0.5D0*COEFDD / UN) 1392 : ENDIF 1393 : ENDIF 1394 : THETA = MIN(1.D0,THETA) 1395 : THETA = MAX(0.D0,THETA) 1396 : 1397 : EXPR1 = (THETA*UN) + COEFDD 1398 : EXPR2 = ((1.D0 - THETA)*UN) - COEFDD 1399 : c IF (EXPR1.LT.0.0D0) THEN 1400 : c WRITE(6,*) 'EXPR1',EXPR1 1401 : c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF, 1402 : c & 'COEFDD=',COEFDD,'THETA= ',THETA,'UN= ',UN 1403 : c ENDIF 1404 : c IF (EXPR2.GT.0.0D0) THEN 1405 : c WRITE(6,*) 'EXPR2',EXPR2 1406 : c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF, 1407 : c & 'COEFDD=',COEFDD,'THETA= ',THETA,'UN= ',UN 1408 : c ENDIF 1409 : c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF, 1410 : c & 'EXPR1= ',EXPR1,'EXPR2= ',EXPR2 1411 : 1412 : 1413 : VAL = (MPOCHP.VPOCHA(NLCG,1)*UN*THETA) + 1414 : & (MPOCHP.VPOCHA(NLCD,1)*(1.D0-THETA)*UN) 1415 : MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) - VAL 1416 : ELSE 1417 : C CONDITIONS AUX LIMITES : TRACE CALCULEE PAR LA DIFFUSION 1418 : IF (INDFR.NE.0) THEN 1419 : COEFDD = MELVA1.VELCHE(INDFR,NLCF) 1420 : ELSE 1421 : COEFDD = 0.0D0 1422 : ENDIF 1423 : C C'EST LE THETA LE MIEUX ADPATE POUR LES CONDITIONS AUX LIMITES 1424 : THETA = 0.0D0 1425 : EXPR1 = - (THETA*UN) - COEFDD 1426 : EXPR2 = - ((1.D0 - THETA)*UN) + COEFDD 1427 : 1428 : IF (EXPR1.GT.0.0D0) THEN 1429 : IF (ABS(UN) .GT. 1e-20) THEN 1430 : THETA = - (0.5D0*COEFDD / UN) 1431 : ENDIF 1432 : ENDIF 1433 : IF (EXPR2.LT.0.0D0) THEN 1434 : IF (ABS(UN) .GT. 1e-20) THEN 1435 : THETA = 1.D0 - (0.5D0*COEFDD / UN) 1436 : ENDIF 1437 : ENDIF 1438 : THETA = MIN(1.D0,THETA) 1439 : THETA = MAX(0.D0,THETA) 1440 : C ANCIENNE VERRUE 1441 : c THETA = 0.0D0 1442 : 1443 : XPRO = 1.D0/NBNO 1444 : TRACE = 0.0D0 1445 : DO JA = 1,NBNO 1446 : TRACE =TRACE + (XPRO*TRGAUX(JA)) 1447 : ENDDO 1448 : VAL = TRACE*UN 1449 : VAL = (MPOCHP.VPOCHA(NLCG,1)*UN*THETA) + 1450 : & (TRACE*(1.D0-THETA)*UN) 1451 : MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) - VAL 1452 : 1453 : c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF, 1454 : c & 'GR1= ',MPOGRA.VPOCHA(NLCF,1),'VAL = ',VAL, 1455 : c & 'TRACE = ',TRACE,'UN= ',UN 1456 : c MPOGRA.VPOCHA(NLCF,1) = MPOGRA.VPOCHA(NLCF,1) - VAL 1457 : c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF, 1458 : c & 'GR1= ',MPOGRA.VPOCHA(NLCF,1) 1459 : ENDIF 1460 : 1461 : C ON COMPLETE MELVA1 POUR LE CONVECTIF 1462 : IF (NLCD.NE.NLCG) THEN 1463 : MELVA1.VELCHE(INDGA,NLCF) = MELVA1.VELCHE(INDGA,NLCF) - 1464 : & (THETA*UN) 1465 : MELVA1.VELCHE(INDDR,NLCF) = MELVA1.VELCHE(INDDR,NLCF) - 1466 : & ((1.D0-THETA)*UN) 1467 : C CONDITION AUX LIMITES : ON RAJOUTE LES DEPENDANCES DES TRACES 1468 : c POUR LES CONDITIONS MIXTES OU DE NEUMAN 1469 : 1470 : ELSE 1471 : 1472 : MELVA1.VELCHE(INDGA,NLCF) = MELVA1.VELCHE(INDGA,NLCF) - 1473 : & (THETA*UN) 1474 : NLFCL = MLENCL.LECT(NGCF) 1475 : C ON EST ICI : CORRIGER 1476 : C ON RAJOUTE CECI POUR L OPTION GRADGEO 1477 : c IF (NLFCL.NE.0) THEN 1478 : c MELVA1.VELCHE(NCON+1,NLCF) = - UN 1479 : c MELEME.NUM(NCON+1,NLCF) = NGCF 1480 : c ENDIF 1481 : 1482 : c IF (NLFCL.EQ.0) THEN 1483 : XPRO = 1.D0/NBNO 1484 : DO JA = 1,NBNO 1485 : INDIC = IPO3.POINT33(NLS(JA)) 1486 : SEGACT INDIC *MOD 1487 : MATR1 = IPO2.POINT(NLS(JA)) 1488 : SEGACT MATR1 *MOD 1489 : 1490 : DO ICOUR = 1,TAB.ID(NLS(JA)) 1491 : IA = ICOUR 1492 : J1 = INDIC.NU(IGNS(JA),IA) 1493 : DO IAUX2 = 2,NCON 1494 : J2 = MELEME.NUM(IAUX2,NLCF) 1495 : IF (J1.EQ.J2) THEN 1496 : IAUX = IAUX2 1497 : GOTO 5159 1498 : ENDIF 1499 : ENDDO 1500 : 5159 CONTINUE 1501 : 1502 : 1503 : CX = MATR1.MAT2(IGNS(JA),IA) 1504 : MELVA1.VELCHE(IAUX,NLCF) = MELVA1.VELCHE(IAUX,NLCF) - 1505 : & ((1.D0-THETA)*UN*CX*XPRO) 1506 : ENDDO 1507 : 1508 : SEGDES INDIC *MOD 1509 : SEGDES MATR1 *MOD 1510 : ENDDO 1511 : c ENDIF 1512 : ENDIF 1513 : 1514 : ENDIF 1515 : ENDIF 1516 : 1517 : c DO J= 1,NBNN 1518 : c PSCA = (MPOGRA.VPOCHA(NLCF,1)*SCN1X) + 1519 : c & (MPOGRA.VPOCHA(NLCF,2)*SCN1Y) + 1520 : c & (MPOGRA.VPOCHA(NLCF,3)*SCN1Z) 1521 : c MPOGRA.VPOCHA(NLCF,1) = PSCA*SCN1X 1522 : c MPOGRA.VPOCHA(NLCF,2) = PSCA*SCN1Y 1523 : c MPOGRA.VPOCHA(NLCF,3) = PSCA*SCN1Z 1524 : c ENDDO 1525 : c IF (NLCF.EQ.791) THEN 1526 : c WRITE(6,*) 'NLCF= ',NLCF,'GR1= ',MPOGRA.VPOCHA(NLCF,1) 1527 : c ENDIF 1528 : c WRITE(6,*) 'NLCF= ',NLCF,'GR2= ',MPOGRA.VPOCHA(NLCF,2) 1529 : c WRITE(6,*) 'NLCF= ',NLCF,'GR3= ',MPOGRA.VPOCHA(NLCF,3) 1530 : 1531 : 1532 : c DO J= 1,NBNN 1533 : c PSCA = (MELVA1.VELCHE(J,NLCF)*SCN1X) + 1534 : c & (MELVA2.VELCHE(J,NLCF)*SCN1Y) + 1535 : c & (MELVA3.VELCHE(J,NLCF)*SCN1Z) 1536 : c MELVA1.VELCHE(J,NLCF) = PSCA*SCN1X 1537 : c MELVA2.VELCHE(J,NLCF) = PSCA*SCN1Y 1538 : c MELVA3.VELCHE(J,NLCF) = PSCA*SCN1Z 1539 : c ENDDO 1540 : c DO J= 1,NBNN 1541 : c WRITE(6,*) 'NLCF= ',NLCF,'J=',J,'MELVA1= ', 1542 : c & MELVA1.VELCHE(J,NLCF) 1543 : c WRITE(6,*) 'NLCF= ',NLCF,'J=',J,'MELVA2= ', 1544 : c & MELVA2.VELCHE(J,NLCF) 1545 : c WRITE(6,*) 'NLCF= ',NLCF,'J=',J,'MELVA3= ', 1546 : c & MELVA3.VELCHE(J,NLCF) 1547 : c WRITE(6,*) 'MELEME=',MELEME.NUM(J,NLCF) 1548 : c ENDDO 1549 : C IF (1.EQ.0) THEN 1550 : C WRITE(6,*) 'NGCF= ',NGCF 1551 : C WRITE(6,*) 'NLCF= ',NLCF,'GR1= ',MPOGRA.VPOCHA(NLCF,1) 1552 : C WRITE(6,*) 'NLCF= ',NLCF,'GR2= ',MPOGRA.VPOCHA(NLCF,2) 1553 : C WRITE(6,*) 'NLCF= ',NLCF,'GR3= ',MPOGRA.VPOCHA(NLCF,3) 1554 : C VALD = MPOCHP.VPOCHA(NLCD,1) 1555 : C WRITE(6,*) 'NLCF= ',NLCF,'VALD= ',VALD 1556 : C WRITE(6,*) 'NLCF= ',NLCF,'VAL= ',VAL 1557 : C WRITE(6,*) 'NLCF= ',NLCF,'TRG= ',TRG 1558 : C WRITE(6,*) 'NLCF= ',NLCF,'TRG2= ',TRG2 1559 : C WRITE(6,*) 'NLCF= ',NLCF,'TRG3= ',TRG3 1560 : C WRITE(6,*) 'NLCF= ',NLCF,'TRG= ',TRGAUX 1561 : C WRite(6,*) 'AG1=',AG1 1562 : C WRite(6,*) 'AG2=',AG2 1563 : C WRite(6,*) 'AD1=',AD1 1564 : C WRite(6,*) 'AD2=',AD2 1565 : C ENDIF 1566 : NAUX2 = MAX(NAUX2,NCON) 1567 : NMOY = NMOY + NCON 1568 : 1569 : c WRITE(6,*) 'INDGA= ',INDGA 1570 : c WRITE(6,*) 'INDDR= ',INDDR 1571 : c DO J= 1,NBNN 1572 : c WRITE(6,*) 'NLCF= ',NLCF,'J=',J,'MELVA1= ', 1573 : c & MELVA1.VELCHE(J,NLCF) 1574 : c ENDDO 1575 : 1576 : 1577 : ENDDO 1578 : NMOY = NMOY/(NFAC*1.D0) 1579 : 1580 : IF (1.EQ.0) THEN 1581 : DO NLCF = 1, NFAC, 1 1582 : MPOGRA.VPOCHA(NLCF,1) = 0.D0 1583 : SCNX=MPONOR.VPOCHA(NLCF,1) 1584 : SCNY=MPONOR.VPOCHA(NLCF,2) 1585 : SCNZ=MPONOR.VPOCHA(NLCF,3) 1586 : SCN1X = SCNX 1587 : SCN1Y = SCNY 1588 : SCN1Z = SCNZ 1589 : NGCF=MELEFL.NUM(2,NLCF) 1590 : DO IVOI=2,MELEME.NUM(/1) 1591 : ICENT = MELEME.NUM(IVOI,NLCF) 1592 : ICEN = MLECEN.LECT(ICENT) 1593 : VAL = 0.0D0 1594 : IF (ICEN.NE.0) THEN 1595 : c WRITE(6,*) 'INTERIEUR' 1596 : VAL = MPOCHP.VPOCHA(ICEN,1) 1597 : ELSE 1598 : ICENL = MLENCL.LECT(ICENT) 1599 : IF (ICENL.GT.0) THEN 1600 : c WRITE(6,*) 'DIRICHLET' 1601 : VAL = MPOVCL.VPOCHA(ICENL,1) 1602 : ENDIF 1603 : ENDIF 1604 : 1605 : COEF1 = MELVA1.VELCHE(IVOI,NLCF) 1606 : MPOGRA.VPOCHA(NLCF,1)= MPOGRA.VPOCHA(NLCF,1) + 1607 : & (COEF1 * VAL) 1608 : 1609 : c WRITE(6,*) 'NLCF= ',NLCF,'VAL= ',VAL 1610 : c WRITE(6,*) 'IVOI= ',IVOI,'MELEME= ', MELEME.NUM(IVOI,NLCF), 1611 : c & 'COEF1 = ',COEF1,'COEF2= ',COEF2,'COEF3= ',COEF3 1612 : ENDDO 1613 : c DO J= 1,NBNN 1614 : c PSCA = (MPOGRA.VPOCHA(NLCF,1)*SCN1X) + 1615 : c & (MPOGRA.VPOCHA(NLCF,2)*SCN1Y) + 1616 : c & (MPOGRA.VPOCHA(NLCF,3)*SCN1Z) 1617 : c MPOGRA.VPOCHA(NLCF,1) = PSCA*SCN1X 1618 : c MPOGRA.VPOCHA(NLCF,2) = PSCA*SCN1Y 1619 : c MPOGRA.VPOCHA(NLCF,3) = PSCA*SCN1Z 1620 : c ENDDO 1621 : c NGCF=MELEFL.NUM(2,NLCF) 1622 : c WRITE(6,*) 'NLCF= ',NLCF,'NGCF= ',NGCF 1623 : c WRITE(6,*) 'MPOGRA1= ', MPOGRA.VPOCHA(NLCF,1) 1624 : c WRITE(6,*) 'MPOGRA2= ', MPOGRA.VPOCHA(NLCF,2) 1625 : c WRITE(6,*) 'MPOGRA3= ', MPOGRA.VPOCHA(NLCF,3) 1626 : ENDDO 1627 : ENDIF 1628 : 1629 : IF (NBSO.EQ.2) THEN 1630 : SEGDES IPT1 1631 : SEGDES IPT2 1632 : ELSEIF (NBSO.EQ.3) THEN 1633 : SEGDES IPT1 1634 : SEGDES IPT2 1635 : SEGDES IPT3 1636 : ELSEIF (NBSO.EQ.4) THEN 1637 : SEGDES IPT1 1638 : SEGDES IPT2 1639 : SEGDES IPT3 1640 : SEGDES IPT4 1641 : ENDIF 1642 : IF (NBSOF.EQ.2) THEN 1643 : SEGDES IPT5 1644 : SEGDES IPT6 1645 : ENDIF 1646 : 1647 : C ON REJUSTE LE CHAMELEM 1648 : 1649 : 1650 : c WRITE(6,*) 'NAUX2= ',NAUX2 1651 : c WRITE(6,*) 'NMOY2= ',NMOY 1652 : c IF (NAUX2.GT.NESSAI) THEN 1653 : c WRITE(6,*) 'NAUX2 = ',NAUX2,'NESSAI = ',NESSAI 1654 : c WRITE(6,*) 'NESSAI TROP PETIT' 1655 : c STOP 1656 : c ENDIF 1657 : N1PTEL=NAUX2 1658 : N1EL=NBELEM 1659 : N2PTEL=0 1660 : N2EL=0 1661 : SEGADJ MELVA1 1662 : 1663 : NBNN = NAUX2 1664 : SEGADJ MELEME 1665 : 1666 : K7 = NFAC 1667 : K8 = NAUX2 1668 : SEGINI ITAB 1669 : 1670 : c ON SUPPRIME LES ZEROS INTERIEURS 1671 : IND = 2 1672 : DO NLCF = 1,NFAC 1673 : IND = 2 1674 : NFIN = NAUX2 1675 : C ON CALCULE D4ABORD LE MAXIMUM DE LA LIGNE 1676 : AUXMA = 0 1677 : DO J=2,NFIN 1678 : ICEN = MLECEN.LECT(MELEME.NUM(J,NLCF)) 1679 : 1680 : IF (ICEN.NE.0) THEN 1681 : AUXMA = MAX(ABS(MELVA1.VELCHE(J,NLCF)),AUXMA) 1682 : ENDIF 1683 : ENDDO 1684 : ITAB.XMAX(NLCF) = AUXMA 1685 : 1686 : DO J=IND, NFIN 1687 : ICEN = MLECEN.LECT(MELEME.NUM(J,NLCF)) 1688 : IF (ABS(MELVA1.VELCHE(J,NLCF)).Le.(1e-14*AUXMA) 1689 : & .AND.(ICEN.NE.0) ) THEN 1690 : 1691 : 1692 : DO K=1,NFIN-J 1693 : AUX = MELVA1.VELCHE(J+K,NLCF) 1694 : ICEN = MLECEN.LECT(MELEME.NUM(J+K,NLCF)) 1695 : C ON DECALE DE K CRAN 1696 : IF (ABS(AUX).gt.(1e-14*AUXMA).OR.(ICEN.EQ.0)) THEN 1697 : DO I=J,NFIN - K 1698 : MELVA1.VELCHE(I,NLCF) = MELVA1.VELCHE(I+K,NLCF) 1699 : MELEME.NUM(I,NLCF) = MELEME.NUM(I+K,NLCF) 1700 : ENDDO 1701 : C MISE A ZERO DES TERMES DU BOUT 1702 : DO I=NFIN-K+1,NFIN 1703 : MELVA1.VELCHE(I,NLCF) = 0.0 1704 : ENDDO 1705 : GOTO 2000 1706 : ENDIF 1707 : ENDDO 1708 : 2000 CONTINUE 1709 : ENDIF 1710 : ENDDO 1711 : ENDDO 1712 : 1713 : c IF (1.EQ.0) THEN 1714 : c ON RECONSTUIT UN CHAMELEM EN SUPPRIMANT LES 0 1715 : 1716 : 1717 : 1718 : 1719 : 1720 : c TABLEAU CALCULANT LE NOMBRE DE VOISIN NON NUL POUR CHAQUE FACE 1721 : NMOY = 0 1722 : INF = NAUX2 1723 : DO NLCF = 1,NFAC 1724 : IREC = 3 1725 : DO J=NAUX2,1,-1 1726 : IF (ABS(MELVA1.VELCHE(J,NLCF)).gt. 1727 : & (1e-14*ITAB.XMAX(NLCF))) THEN 1728 : IREC = J 1729 : GOTO 1111 1730 : ENDIF 1731 : ENDDO 1732 : 1111 CONTINUE 1733 : ITAB.TABL(NLCF) = IREC 1734 : NMOY = NMOY + IREC 1735 : ENDDO 1736 : NMOY = NMOY/(NFAC*1.D0) 1737 : c WRITE(6,*) 'NEWMOY2= ',NMOY 1738 : 1739 : 1740 : C TAB1(U) TABLEAU QUI CONTIENT LE NOMBRE DE FACE AYANT U VOISIN 1741 : NMAX = 0 1742 : DO ICOUR = 1,NAUX2 1743 : ITAB.TABL1(ICOUR) = 0 1744 : ENDDO 1745 : 1746 : DO NLCF = 1,NFAC 1747 : ICOUR = ITAB.TABL(NLCF) 1748 : ITAB.TABL1(ICOUR)=ITAB.TABL1(ICOUR) + 1 1749 : ENDDO 1750 : 1751 : C ON COMPTE LE NOMVRE DE SOUS DOMAINE 1752 : NTSOUS = 0 1753 : DO ICOUR = 1,NAUX2 1754 : IF (ITAB.TABL1(ICOUR) .NE.0) NTSOUS = NTSOUS +1 1755 : ENDDO 1756 : C IPOS INDICE DE LA PREMIERE FACE AYANT I VOISIN 1757 : C ICOUR INDICE COURANT INITIALISE A IPOS 1758 : 1759 : 1760 : ITAB.IPOS(1) = 1 1761 : DO I = 2,NAUX2 1762 : ITAB.IPOS(I) = ITAB.IPOS(I-1) + ITAB.TABL1(I-1) 1763 : ITAB.ICOURANT(I) = ITAB.IPOS(I) 1764 : ENDDO 1765 : 1766 : 1767 : c TABL2 TABLEAU QUI RANGE DANS L4ODRES DES SOUS DOMAINES LES FACES NLCF 1768 : 1769 : DO I =1,NFAC 1770 : IHELP = ITAB.TABL(I) 1771 : IAUX = ITAB.ICOURANT(IHELP) 1772 : ITAB.TABL2(IAUX) = I 1773 : IAUX2 = ITAB.TABL(I) 1774 : ITAB.ICOURANT(IAUX2) = ITAB.ICOURANT(IAUX2) + 1 1775 : ENDDO 1776 : 1777 : 1778 : C**** Initialisation du MCHELM 1779 : C 1780 : N1=NTSOUS 1781 : N2=1 1782 : N3=6 1783 : L1=8 1784 : SEGINI MCHEL2 1785 : MCHEL2.TITCHE='Gradient' 1786 : MCHEL2.IFOCHE=IFOUR 1787 : C 1788 : ISOUS=0 1789 : NBSOUS=0 1790 : NBREF=0 1791 : DO I1 = 1, NAUX2, 1 1792 : NBELEM=ITAB.TABL1(I1) 1793 : IF(NBELEM .GT. 0)THEN 1794 : ISOUS=ISOUS+1 1795 : NBNN=I1 1796 : SEGINI IPT8 1797 : C ITYPEL=32 -> 'POLY' 1798 : ITYPEL=32 1799 : MCHEL2.IMACHE(ISOUS)=IPT8 1800 : SEGINI MCHAM2 1801 : MCHEL2.ICHAML(ISOUS)=MCHAM2 1802 : MCHAM2.NOMCHE(1)='SCAL' 1803 : MCHAM2.TYPCHE(1)='REAL*8 ' 1804 : N1PTEL=NBNN 1805 : N1EL=NBELEM 1806 : N2PTEL=0 1807 : N2EL=0 1808 : SEGINI MELVA2 1809 : MCHAM2.IELVAL(1)=MELVA2 1810 : 1811 : DO I2=1,NBELEM,1 1812 : DO I3=1,NBNN,1 1813 : IAUX = ITAB.IPOS(I1) 1814 : IP= ITAB.TABL2(I2 + IAUX - 1) 1815 : c WRITE(6,*) 'IP= ',IP,'I1= ',I1,'I2=',I2,'I3=',I3 1816 : IPT8.NUM(I3,I2)= MELEME.NUM(I3,IP) 1817 : MELVA2.VELCHE(I3,I2)=MELVA1.VELCHE(I3,IP) 1818 : ENDDO 1819 : ENDDO 1820 : C 1821 : SEGDES MCHAM2 1822 : SEGDES IPT8 1823 : SEGDES MELVA2 1824 : ENDIF 1825 : ENDDO 1826 : SEGDES MCHEL2 1827 : SEGDES ITAB 1828 : c ENDIF 1829 : 1830 : 1831 : c SEGDES MCHAML 1832 : c SEGDES MELEME 1833 : c SEGDES MELVA1 1834 : c SEGDES MCHELM 1835 : 1836 : C 3009 1837 : SEGSUP MCHAML 1838 : SEGSUP MELEME 1839 : SEGSUP MELVA1 1840 : SEGSUP MCHELM 1841 : 1842 : 1843 : C SUPRESSION DES SEGMENTS 1844 : K3 = NSOMM 1845 : DO I = 1,K3 1846 : MATR1 = IPO2.POINT(I) 1847 : SEGACT MATR1 1848 : SEGSUP MATR1 1849 : ENDDO 1850 : 1851 : SEGSUP IPO3 1852 : SEGSUP IPO2 1853 : SEGSUP INDLI 1854 : SEGSUP TAB 1855 : SEGSUP IND2 1856 : c SEGSUP IND 1857 : c SEGSUP IND22 1858 : c SEGSUP VAL1 1859 : c SEGSUP VAL2 1860 : SEGSUP SCMB 1861 : C 3009 1862 : SEGSUP NBCOT 1863 : 1864 : 1865 : 1866 : 9999 CONTINUE 1867 : RETURN 1868 : END 1869 : 1870 : 1871 : 1872 : 1873 : 1874 : 1875 : 1876 : 1877 : 1878 : 1879 : 1880 : 1881 : 1882 : 1883 :
© Cast3M 2003 - All rights reserved.
Disclaimer