src/tool/ES_Xfig/pppma3.F

   1       subroutine pppma3 ( nbtrvi, nbquvi,
   2      >                    nntrvi, nnquvi,
   3      >                    coopro,
   4      >                    posini, xyzfac, tabaux, nivsup,
   5      >                    ulsort, langue, codret )
   6 c ______________________________________________________________________
   7 c
   8 c                             H O M A R D
   9 c
  10 c Outil de Maillage Adaptatif par Raffinement et Deraffinement d'EDF R&D
  11 c
  12 c Version originale enregistree le 18 juin 1996 sous le numero 96036
  13 c aupres des huissiers de justice Simart et Lavoir a Clamart
  14 c Version 11.2 enregistree le 13 fevrier 2015 sous le numero 2015/014
  15 c aupres des huissiers de justice
  16 c Lavoir, Silinski & Cherqui-Abrahmi a Clamart
  17 c
  18 c    HOMARD est une marque deposee d'Electricite de France
  19 c
  20 c Copyright EDF 1996
  21 c Copyright EDF 1998
  22 c Copyright EDF 2002
  23 c Copyright EDF 2020
  24 c ______________________________________________________________________
  25 c
  26 c     Post-Processeur - Preparation du MAillage - phase 3
  27 c     -    -            -              --               -
  28 c ______________________________________________________________________
  29 c
  30 c but : recherche de l'ordre d'affichage des faces pour les cacher
  31 c
  32 c       par defaut, on affiche les objets vus par l'observateur avec
  33 c       l'axe (oz+) dans l'oeil, donc regardant de z>0 vers z<0.
  34 c       on utilise l'algorithme dit du peintre, ou encore du z-buffer
  35 c       utilise par les affichages graphiques standards, consistant a
  36 c       imprimer les objets dans l'ordre inverse de leur eloignement,
  37 c       i.e. a imprimer en dernier, et donc par dessus le reste, les
  38 c       objets les plus proches. (donc avec les z les plus grands)
  39 c ______________________________________________________________________
  40 c .        .     .        .                                            .
  41 c .  nom   . e/s . taille .           description                      .
  42 c .____________________________________________________________________.
  43 c . nbtrvi . e   .   1    . nombre triangles visualisables             .
  44 c . nbquvi . e   .   1    . nombre de quadrangles visualisables        .
  45 c . nntrvi . e   .10nbtrvi. 1 : niveau du triangle a afficher          .
  46 c .        .     .        . 2 : numero HOMARD du triangle              .
  47 c .        .     .        . 3, 4, 5 : numeros des noeuds p1            .
  48 c .        .     .        . 6 : famille du triangle                    .
  49 c .        .     .        . 7, 8, 9 : numeros des noeuds p2            .
  50 c .        .     .        . 10 : numero du noeud interne               .
  51 c . nnquvi . e   .12nbquvi. 1 : niveau du quadrangle a afficher        .
  52 c .        .     .        . 2 : numero HOMARD du quadrangle            .
  53 c .        .     .        . 3, 4, 5, 6 : numeros des noeuds p1         .
  54 c .        .     .        . 7 : famille du quadrangle                  .
  55 c .        .     .        . 8, 9, 10, 11 : numeros des noeuds p2       .
  56 c .        .     .        . 12 : numero du noeud interne               .
  57 c . coopro . e   .   3*   . coordonnees projetees de :                 .
  58 c .        .     .nbnot+12. le triedre : -8:O ; -9:I ; -10:J ; -11:K   .
  59 c .        .     .        . la fenetre de zoom : de -7 a 0 en 3D ou    .
  60 c .        .     .        .                      de -3 a 0 en 2D       .
  61 c .        .     .        . les noeuds de 1 a nbnoto                   .
  62 c . posini . aux . nbquvi . tableau auxiliaire de renumerotation des   .
  63 c .        .     .+nbtrvi . faces en fonction de l'affichage           .
  64 c . xyzfac .  /  .nbtrvi+ . tableau de travail reel                    .
  65 c .        .     .nbtrvi,9.                                            .
  66 c . tabaux . aux . nbquvi . tableau auxiliaire                         .
  67 c .        .     .+nbtrvi .                                            .
  68 c . nivsup . e   .   1    . niveau superieur present dans le maillage  .
  69 c . ulsort . e   .   1    . unite logique de la sortie generale        .
  70 c . langue . e   .    1   . langue des messages                        .
  71 c .        .     .        . 1 : francais, 2 : anglais                  .
  72 c . codret .  s  .    1   . code de retour des modules                 .
  73 c .        .     .        . 0 : pas de probleme                        .
  74 c ._____________________________________________________________________
  75 c
  76 c====
  77 c 0. declarations et dimensionnement
  78 c====
  79 c
  80 c 0.1. ==> generalites
  81 c
  82       implicit none
  83       save
  84 c
  85       character*6 nompro
  86       parameter ( nompro = 'PPPMA3' )
  87 c
  88 #include "nblang.h"
  89 c
  90 c 0.2. ==> communs
  91 c
  92 #include "envex1.h"
  93 #include "envca1.h"
  94 #include "nombno.h"
  95 c
  96 c 0.3. ==> arguments
  97 c
  98       integer nivsup
  99       integer nbtrvi, nbquvi
 100       integer nntrvi(10,nbtrvi)
 101       integer nnquvi(12,nbquvi)
 102       integer posini(nbtrvi+nbquvi)
 103       integer tabaux(nbtrvi+nbquvi)
 104 c
 105       double precision coopro(3,-11:nbnoto)
 106       double precision xyzfac(nbtrvi+nbquvi,9)
 107 c
 108       integer ulsort, langue, codret
 109 c
 110 c 0.4. ==> variables locales
 111 c
 112       integer nutrvi, nuquvi
 113       integer nbfavi
 114       integer rangfa, rangqu
 115       integer nbfast, nbfas0
 116       integer rangfd
 117       integer iaux, jaux
 118 c
 119       double precision daux
 120       double precision borne
 121       double precision xminfa, xmaxfa
 122       double precision yminfa, ymaxfa
 123       double precision zminfa, zmaxfa
 124       double precision coface(2,3)
 125       double precision v1(3), v2(3), v3(3)
 126 c
 127       integer nbmess
 128       parameter ( nbmess = 10 )
 129       character*80 texte(nblang,nbmess)
 130 c
 131 c 0.5. ==> initialisations
 132 c
 133       data borne / 1.d-8 /
 134 c
 135 c 0.5. ==> initialisations
 136 c_______________________________________________________________________
 137 c
 138 c====
 139 c 1. prealables
 140 c====
 141 c
 142 #include "impr01.h"
 143 c
 144 #ifdef _DEBUG_HOMARD_
 145       write (ulsort,texte(langue,1)) 'Entree', nompro
 146       call dmflsh (iaux)
 147 #endif
 148 c
 149 #include "impr03.h"
 150 c
 151 #ifdef _DEBUG_HOMARD_
 152       write (ulsort,90002) 'nbtrvi', nbtrvi
 153       write (ulsort,90002) 'nbquvi', nbquvi
 154       write (ulsort,90002) 'sdim  ', sdim
 155       write (ulsort,90002) 'nivsup', nivsup
 156 #endif
 157 c
 158       codret = 0
 159 c
 160 c====
 161 c 2. En dimension 2, les faces sont reclassses par niveau croissant
 162 c    A la fin de cette etape, posini(1) contient l'indice de la face a
 163 c    tracer en premier, parce que la plus loin, posini(2) contient
 164 c    l'indice de la suivante a tracer et ainsi de suite, jusqu'a
 165 c    posini(nbfavi) qui contient l'indice de la derniere face
 166 c    a tracer parce que la plus proche de l'observateur.
 167 c    Attention : posini ne contient pas les numeros des faces mais une
 168 c    indirection dans la liste des faces a traiter.
 169 c====
 170 cgn      do 34 , jaux = 1 , nbnoto
 171 cgn        write(*,1797)jaux,(coopro(iaux,jaux),iaux=1,sdim)
 172 cgn   34 continue
 173 c
 174       if ( sdim.le.2 ) then
 175 c
 176         jaux = 0
 177 c
 178         do 21 , iaux = 0 , nivsup+1
 179 c
 180           do 211 , nutrvi = 1 , nbtrvi
 181             if ( nntrvi(1,nutrvi).eq.iaux ) then
 182               jaux = jaux + 1
 183               posini(jaux) = nutrvi
 184             endif
 185   211     continue
 186 c
 187           do 212 , nuquvi = 1 , nbquvi
 188             if ( nnquvi(1,nuquvi).eq.iaux ) then
 189               jaux = jaux + 1
 190               posini(jaux) = nbtrvi+nuquvi
 191             endif
 192   212     continue
 193 c
 194    21   continue
 195 c
 196 c====
 197 c 2. En dimension 3, on applique l'algorithme du z-buffer
 198 c====
 199 c
 200       else
 201 c
 202 c====
 203 c 3. Calcul des dimensions extremes des faces et des vecteurs normaux
 204 c====
 205 c 3.1. ==> Extrema pour un triangle
 206 c
 207       do 31 , nutrvi = 1 , nbtrvi
 208 c
 209         xminfa = min(coopro(1,nntrvi(3,nutrvi)),
 210      >           coopro(1,nntrvi(4,nutrvi)),coopro(1,nntrvi(5,nutrvi)))
 211         xmaxfa = max(coopro(1,nntrvi(3,nutrvi)),
 212      >           coopro(1,nntrvi(4,nutrvi)),coopro(1,nntrvi(5,nutrvi)))
 213         xyzfac(nutrvi,1) = xminfa + borne * (xmaxfa - xminfa)
 214         xyzfac(nutrvi,2) = xmaxfa - borne * (xmaxfa - xminfa)
 215 c
 216         yminfa = min(coopro(2,nntrvi(3,nutrvi)),
 217      >           coopro(2,nntrvi(4,nutrvi)),coopro(2,nntrvi(5,nutrvi)))
 218         ymaxfa = max(coopro(2,nntrvi(3,nutrvi)),
 219      >           coopro(2,nntrvi(4,nutrvi)),coopro(2,nntrvi(5,nutrvi)))
 220         xyzfac(nutrvi,3) = yminfa + borne * (ymaxfa - yminfa)
 221         xyzfac(nutrvi,4) = ymaxfa - borne * (ymaxfa - yminfa)
 222 c
 223         zminfa = min(coopro(3,nntrvi(3,nutrvi)),
 224      >           coopro(3,nntrvi(4,nutrvi)),coopro(3,nntrvi(5,nutrvi)))
 225         zmaxfa = max(coopro(3,nntrvi(3,nutrvi)),
 226      >           coopro(3,nntrvi(4,nutrvi)),coopro(3,nntrvi(5,nutrvi)))
 227         xyzfac(nutrvi,5) = zminfa + borne * (zmaxfa - zminfa)
 228         xyzfac(nutrvi,6) = zmaxfa - borne * (zmaxfa - zminfa)
 229 cgn      write(*,1798)nutrvi,(nntrvi(iaux,nutrvi),iaux=2,5),
 230 cgn     > (xyzfac(nutrvi,iaux),iaux=1,6)
 231 c
 232    31 continue
 233 c
 234 c 3.2. ==> Extrema pour un quadrangle
 235 c
 236       do 32 , nuquvi = 1 , nbquvi
 237 c
 238         xminfa = coopro(1,nnquvi(3,nuquvi))
 239         xmaxfa = xminfa
 240         yminfa = coopro(2,nnquvi(3,nuquvi))
 241         ymaxfa = yminfa
 242         zminfa = coopro(3,nnquvi(3,nuquvi))
 243         zmaxfa = zminfa
 244         do 321 , iaux = 4 , 6
 245           daux = coopro(1,nnquvi(iaux,nuquvi))
 246           xminfa = min ( daux , xminfa )
 247           xmaxfa = max ( daux , xmaxfa )
 248           daux = coopro(2,nnquvi(iaux,nuquvi))
 249           yminfa = min ( daux , yminfa )
 250           ymaxfa = max ( daux , ymaxfa )
 251           daux = coopro(3,nnquvi(iaux,nuquvi))
 252           zminfa = min ( daux , zminfa )
 253           zmaxfa = max ( daux , zmaxfa )
 254   321   continue
 255         rangfa = nbtrvi+nuquvi
 256         xyzfac(rangfa,1) = xminfa + borne * (xmaxfa - xminfa)
 257         xyzfac(rangfa,2) = xmaxfa - borne * (xmaxfa - xminfa)
 258         xyzfac(rangfa,3) = yminfa + borne * (ymaxfa - yminfa)
 259         xyzfac(rangfa,4) = ymaxfa - borne * (ymaxfa - yminfa)
 260         xyzfac(rangfa,5) = zminfa + borne * (zmaxfa - zminfa)
 261         xyzfac(rangfa,6) = zmaxfa - borne * (zmaxfa - zminfa)
 262 cgn      write(*,1799)nuquvi,(nnquvi(iaux,nuquvi),iaux=2,6),
 263 cgn     > (xyzfac(nuquvi,iaux),iaux=1,6)
 264 c
 265    32 continue
 266 c
 267 c 3.3. ==> Vecteurs normaux
 268 c
 269       nbfavi = nbtrvi + nbquvi
 270 cgn      write (ulsort,1795) 'nbfavi', nbfavi
 271       do 33 , rangfa = 1 , nbfavi
 272 c
 273         if ( rangfa.le.nbtrvi ) then
 274           v1(1) = coopro(1,nntrvi(3,rangfa))
 275           v1(2) = coopro(2,nntrvi(3,rangfa))
 276           v1(3) = coopro(3,nntrvi(3,rangfa))
 277           v2(1) = coopro(1,nntrvi(4,rangfa))
 278           v2(2) = coopro(2,nntrvi(4,rangfa))
 279           v2(3) = coopro(3,nntrvi(4,rangfa))
 280           v3(1) = coopro(1,nntrvi(5,rangfa))
 281           v3(2) = coopro(2,nntrvi(5,rangfa))
 282           v3(3) = coopro(3,nntrvi(5,rangfa))
 283         else
 284           rangqu = rangfa - nbtrvi
 285           v1(1) = coopro(1,nnquvi(3,rangqu))
 286           v1(2) = coopro(2,nnquvi(3,rangqu))
 287           v1(3) = coopro(3,nnquvi(3,rangqu))
 288           v2(1) = coopro(1,nnquvi(4,rangqu))
 289           v2(2) = coopro(2,nnquvi(4,rangqu))
 290           v2(3) = coopro(3,nnquvi(4,rangqu))
 291           v3(1) = coopro(1,nnquvi(5,rangqu))
 292           v3(2) = coopro(2,nnquvi(5,rangqu))
 293           v3(3) = coopro(3,nnquvi(5,rangqu))
 294         endif
 295 cgn          write(*,1796)'v1',v1
 296 cgn          write(*,1796)'v2',v2
 297 cgn          write(*,1796)'v3',v3
 298 c
 299         coface(1,1) = v2(1) - v1(1)
 300         coface(1,2) = v2(2) - v1(2)
 301         coface(1,3) = v2(3) - v1(3)
 302         coface(2,1) = v3(1) - v1(1)
 303         coface(2,2) = v3(2) - v1(2)
 304         coface(2,3) = v3(3) - v1(3)
 305 cgn          write(*,1796)'v1',v1
 306 cgn          write(*,1796)'v2',v2
 307 cgn          write(*,1796)'v3',v3
 308 c
 309         xyzfac(rangfa,7) = coface(1,2)*coface(2,3)
 310      >                   - coface(1,3)*coface(2,2)
 311         xyzfac(rangfa,8) = coface(1,3)*coface(2,1)
 312      >                   - coface(1,1)*coface(2,3)
 313         xyzfac(rangfa,9) = coface(1,1)*coface(2,2)
 314      >                   - coface(1,2)*coface(2,1)
 315 cgn        write(*,1796)'normale',(xyzfac(rangfa,iaux),iaux=7,9)
 316 c
 317    33 continue
 318 c
 319  1795 format(3(a,' =',i5,' , '))
 320  1796 format(a,6f12.5)
 321  1797 format(i5,' *',6f12.5)
 322  1798 format(i4,' :',i5,' *',3i4,' *',6f12.5)
 323  1799 format(i4,' :',i5,' *',4i4,' *',6f12.5)
 324 c
 325 c====
 326 c 4. On parcourt toutes les faces a afficher et on les range de la
 327 c    plus eloignee a la plus proche du point de vue.
 328 c    Groso modo, les plus proches sont avec les z les plus grands.
 329 c    A la fin de cette etape, posini(1) contient l'indice de la face a
 330 c    tracer en premier, parce que la plus loin, posini(2) contient
 331 c    l'indice de la suivante a tracer et ainsi de suite, jusqu'a
 332 c    posini(nbfavi) qui contient l'indice de la derniere face
 333 c    a tracer parce que la plus proche de l'observateur.
 334 c    Attention : posini ne contient pas les numeros des faces mais une
 335 c    indirection dans la liste des faces a traiter.
 336 c====
 337 c
 338 c 4.1. ==> A priori, les faces ne sont en contact avec aucune autre
 339 c
 340       do 41 , rangfa = 1 , nbfavi
 341         tabaux(rangfa) = 0
 342    41 continue
 343 c
 344 c 4.2. ==> A priori, la 1ere face est en premiere position
 345 c
 346       nbfast = 1
 347       posini(nbfast) = 1
 348       tabaux(nbfast) = 1
 349 c
 350 c 4.3. ==> On examine toutes les faces en ne classant que celles qui
 351 c          ont une partie commune avec celles deja classees
 352 c          On n'examine que les faces non encore stockees
 353 c          On boucle jusqu'a ce qu'il n'y en ait plus
 354 c          Role de tabaux :
 355 c          . tabaux(i) = 1 : la face i est stockee
 356 c          . tabaux(i) = 0 : la face i n'est pas stockee
 357 c
 358       rangfd = nbfast + 1
 359 c
 360   430 continue
 361 c
 362 cgn      write (ulsort,*) '---------- NOUVELLE SERIE ---------'
 363 cgn      write (ulsort,1795) 'Depart avec rangfd', rangfd,'nbfast',nbfast
 364 c
 365       nbfas0 = nbfast
 366 c
 367       do 43 , rangfa = rangfd , nbfavi
 368 c
 369         if ( codret.eq.0 ) then
 370 c
 371         if ( tabaux(rangfa).eq.0 ) then
 372 c
 373           iaux = rangfa
 374 #ifdef _DEBUG_HOMARD_
 375       write (ulsort,texte(langue,3)) 'PPPMA4', nompro
 376 #endif
 377           call pppma4 ( iaux, nbfast,
 378      >                  nbtrvi, nbquvi,
 379      >                  nntrvi, nnquvi,
 380      >                  coopro,
 381      >                  posini, xyzfac, tabaux,
 382      >                  ulsort, langue, codret )
 383 c
 384             if ( rangfa.eq.-70)  stop
 385         endif
 386 c
 387         endif
 388 c
 389    43 continue
 390 c
 391 c 4.4. ==> Bilan
 392 c
 393         if ( codret.eq.0 ) then
 394 c
 395 cgn        write (ulsort,*) 'nbfas0, nbfast =', nbfas0, nbfast
 396 cgn        if ( nbfast.eq.222 ) then
 397 cgn         nbquvi=nbfast
 398 cgn          goto 45
 399 cgn        endif
 400 c
 401 c       . Si le nombre de faces stockees est egal au nombre de faces a
 402 c         visualiser, c'est fini
 403 c
 404 cgn        write (ulsort,*) 'Bilan'
 405         if ( nbfast.eq.nbfavi ) then
 406 cgn          write (ulsort,1795) 'Fin avec nbfast', nbfast
 407           goto 45
 408 c
 409 c          . Si le nombre de faces stockees a change, c'est que l'on
 410 c            se trouve encore dans une suite de faces qui se superposent
 411 c          . Sinon, c'est qu'il ne reste plus de faces se superposant
 412 c            a l'une de la serie en cours. Il faut commencer une
 413 c            autre serie.
 414 c          Dans les deux cas, on va refaire un tour en partant de la
 415 c          premiere face non superposee.
 416 c
 417         else
 418 cgn      do 998 , iaux = 1 , nbfavi
 419 cgn        if ( iaux.le.nbtrvi ) then
 420 cgn          jaux = nntrvi(2,iaux)
 421 cgn        else
 422 cgn          jaux = -nnquvi(2,iaux-nbtrvi)
 423 cgn        endif
 424 cgn        if ( tabaux(iaux).ne.0 ) then
 425 cgn          write (*,1999) iaux,tabaux(iaux),jaux
 426 cgn 1999 format('tabaux(',i5,') =',i5,', face',i5)
 427 cgn        endif
 428 cgn 998  continue
 429 c
 430           do 44 , rangfa = rangfd , nbfavi
 431             if ( tabaux(rangfa).eq.0 ) then
 432               iaux = rangfa
 433               goto 441
 434             endif
 435    44     continue
 436   441     continue
 437           rangfd = iaux
 438 c
 439 cgn      write (ulsort,1795) 'nbfas0', nbfas0, 'nbfast', nbfast,
 440 cgn     > 'rangfd', rangfd
 441           if ( nbfast.eq.nbfas0 ) then
 442 cgn         write (ulsort,*) '---------- NOUVELLE SERIE ---------'
 443 cgn         write (ulsort,*) '---------- FACE', -nnquvi(2,rangfd-nbtrvi)
 444             nbfast = nbfast + 1
 445 cgn        write (ulsort,888)
 446 cgn     > (rangfa,posini(rangfa),rangfa =1,nbtrvi+nbquvi)
 447 cgn  888 format(5(i3,i4,' * '))
 448             posini(nbfast) = rangfd
 449             tabaux(rangfd) = 1
 450             rangfd = rangfd + 1
 451           endif
 452 c
 453 cgn            write (ulsort,*) 'goto 430   goto 430   goto 430'
 454           goto 430
 455 c
 456         endif
 457 c
 458         endif
 459 c
 460 c 4.5. ==> C'est fini
 461 c
 462    45     continue
 463 c
 464 cgn              nbfavi=nbtrvi+nbquvi
 465 cgn      do 999 , iaux = 1 , nbfast
 466 cgn        if ( posini(iaux).le.nbtrvi ) then
 467 cgn          jaux = nntrvi(2,posini(iaux))
 468 cgn        else
 469 cgn          jaux = -nnquvi(2,posini(iaux)-nbtrvi)
 470 cgn        endif
 471 cgn        write (*,*) iaux,posini(iaux),'-eme face :',jaux
 472 cgn 999  continue
 473 c
 474       endif
 475 c
 476 c====
 477 c 5. la fin
 478 c====
 479 c
 480       if ( codret.ne.0 ) then
 481 c
 482 #include "envex2.h"
 483 c
 484       write (ulsort,texte(langue,1)) 'Sortie', nompro
 485       write (ulsort,texte(langue,2)) codret
 486 c
 487       endif
 488 c
 489 #ifdef _DEBUG_HOMARD_
 490       write (ulsort,texte(langue,1)) 'Sortie', nompro
 491       call dmflsh (iaux)
 492 #endif
 493 c
 494       end