src/tool/AP_Conversion/pcset4.F

   1       subroutine pcset4 ( etan, etanp1, tehn, tehnp1, typint,
   2      >                    prfcan, prfcap,
   3      >                    hettet, filtet, nbante, anfite,
   4      >                    nteeca, ntesca,
   5      >                    nbfonc, vafoen, vafott,
   6      >                    ulsort, langue, codret )
   7 c ______________________________________________________________________
   8 c
   9 c                             H O M A R D
  10 c
  11 c Outil de Maillage Adaptatif par Raffinement et Deraffinement d'EDF R&D
  12 c
  13 c Version originale enregistree le 18 juin 1996 sous le numero 96036
  14 c aupres des huissiers de justice Simart et Lavoir a Clamart
  15 c Version 11.2 enregistree le 13 fevrier 2015 sous le numero 2015/014
  16 c aupres des huissiers de justice
  17 c Lavoir, Silinski & Cherqui-Abrahmi a Clamart
  18 c
  19 c    HOMARD est une marque deposee d'Electricite de France
  20 c
  21 c Copyright EDF 1996
  22 c Copyright EDF 1998
  23 c Copyright EDF 2002
  24 c Copyright EDF 2020
  25 c ______________________________________________________________________
  26 c
  27 c    aPres adaptation - Conversion de Solution Elements de volume -
  28 c     -                 -             -        -
  29 c                       Tetraedres d'etat anterieur 4
  30 c                       -                           -
  31 c ______________________________________________________________________
  32 c .        .     .        .                                            .
  33 c .  nom   . e/s . taille .           description                      .
  34 c .____________________________________________________________________.
  35 c . etan   . e   .    1   . ETAt du tetraedre a l'iteration N          .
  36 c . etanp1 . e   .    1   . ETAt du tetraedre a l'iteration N+1        .
  37 c . tehn   . e   .    1   . TEtraedre courant en numerotation Homard   .
  38 c .        .     .        . a l'iteration N                            .
  39 c . tehnp1 . e   .    1   . TEtraedre courant en numerotation Homard   .
  40 c .        .     .        . a l'iteration N+1                          .
  41 c . typint . e   .   1    . type d'interpolation                       .
  42 c .        .     .        .  0, si automatique                         .
  43 c .        .     .        .  elements : 0 si intensif, sans orientation.
  44 c .        .     .        .             1 si extensif, sans orientation.
  45 c .        .     .        .             2 si intensif, avec orientation.
  46 c .        .     .        .             3 si extensif, avec orientation.
  47 c .        .     .        .  noeuds : 1 si degre 1                     .
  48 c .        .     .        .           2 si degre 2                     .
  49 c .        .     .        .           3 si iso-P2                      .
  50 c . prfcan . e   .   *    . En numero du calcul a l'iteration n :      .
  51 c .        .     .        . 0 : l'entite est absente du profil         .
  52 c .        .     .        . i : l'entite est au rang i dans le profil  .
  53 c . prfcap . es  .   *    . En numero du calcul a l'iteration n+1 :    .
  54 c .        .     .        . 0 : l'entite est absente du profil         .
  55 c .        .     .        . 1 : l'entite est presente dans le profil   .
  56 c . hettet . e   . nbteto . historique de l'etat des tetraedres        .
  57 c . filtet . e   . nbteto . premier fils des tetraedres                .
  58 c . nbante . e   .   1    . nombre de tetraedres decoupes par   .
  59 c .        .     .        . conformite sur le maillage avant adaptation.
  60 c . anfite . e   . nbante . tableau filtet du maillage de l'iteration n
  61 c . nteeca . e   . reteto . numero des tetraedres dans le calcul entree.
  62 c . ntesca . e   . rsteto . numero des tetraedres dans le calcul sortie.
  63 c . nbfonc . e   .    1   . nombre de fonctions elements de volume     .
  64 c . vafoen . e   . nbfonc*. variables en entree de l'adaptation        .
  65 c .        .     . nbeven .                                            .
  66 c . vafott . es  . nbfonc*. variables en sortie de l'adaptation        .
  67 c .        .     . nbevso .                                            .
  68 c . ulsort . e   .   1    . numero d'unite logique de la liste standard.
  69 c . langue . e   .    1   . langue des messages                        .
  70 c .        .     .        . 1 : francais, 2 : anglais                  .
  71 c . codret . es  .    1   . code de retour des modules                 .
  72 c .        .     .        . 0 : pas de probleme                        .
  73 c .        .     .        . 1 : 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 = 'PCSET4' )
  87 c
  88 #include "nblang.h"
  89 #include "fracta.h"
  90 #include "fractc.h"
  91 #include "fractf.h"
  92 #include "fractg.h"
  93 #include "fracth.h"
  94 c
  95 c 0.2. ==> communs
  96 c
  97 #include "nombte.h"
  98 #include "nombsr.h"
  99 #include "nomber.h"
 100 c
 101 c 0.3. ==> arguments
 102 c
 103       integer etan, etanp1, tehn, tehnp1
 104       integer typint
 105       integer nbfonc
 106       integer prfcan(*), prfcap(*)
 107       integer hettet(nbteto), filtet(nbteto)
 108       integer nbante
 109       integer anfite(nbante)
 110       integer nteeca(reteto), ntesca(rsteto)
 111 c
 112       double precision vafoen(nbfonc,*)
 113       double precision vafott(nbfonc,*)
 114 c
 115       integer ulsort, langue, codret
 116 c
 117 c 0.4. ==> variables locales
 118 c
 119 c     tecnp1 = TEtraedre courant en numerotation du Calcul a l'it. N+1
 120 c
 121       integer tecnp1
 122 c
 123 c     f1hp = Fils 1er du tetraedre en numerotation Homard a l'it. N+1
 124 c     fihp = Fils ieme u tetraedre en numerotation Homard a l'it. N+1
 125 c     f1cp = Fils 1er du tetraedre en numerota. du Calcul a l'it. N+1
 126 c     f2cp = Fils 2eme du tetraedre en numerota. du Calcul a l'it. N+1
 127 c     f3cp = Fils 3eme du tetraedre en numerota. du Calcul a l'it. N+1
 128 c     f4cp = Fils 4eme du tetraedre en numerota. du Calcul a l'it. N+1
 129 c
 130       integer f1hp, fihp
 131       integer f1cp, f2cp, f3cp, f4cp
 132 c
 133 c     f1hn = Fils 1er du tetraedre en numerotation Homard a l'it. N
 134 c     f1cn = Fils 1er du tetraedre en numerotation du Calcul a l'it. N
 135 c     f2cn = Fils 2eme du tetraedre en numerotation du Calcul a l'it. N
 136 c     f3cn = Fils 3eme du tetraedre en numerotation du Calcul a l'it. N
 137 c     f4cn = Fils 4eme du tetraedre en numerotation du Calcul a l'it. N
 138 c
 139       integer f1hn
 140       integer f1cn, f2cn, f3cn, f4cn
 141 c
 142 c     f1fcp = Fils 1er du Fils en numerotation Calcul a l'it. N+1
 143 c     f2fcp = Fils 2eme du Fils en numerotation Calcul a l'it. N+1
 144 c     f3fcp = Fils 3eme du Fils en numerotation Calcul a l'it. N+1
 145 c     f4fcp = Fils 4eme du Fils en numerotation Calcul a l'it. N+1
 146 c
 147       integer f1fcp, f2fcp, f3fcp, f4fcp
 148 c
 149       integer nrofon
 150       integer coderr
 151 c
 152       integer iaux
 153       integer lglist, nrlist
 154       integer list(30)
 155 c
 156       double precision daux
 157       double precision daux1
 158 c
 159       integer nbmess
 160       parameter ( nbmess = 10 )
 161       character*80 texte(nblang,nbmess)
 162 c
 163 c 0.5. ==> initialisations
 164 c ______________________________________________________________________
 165 c
 166 c====
 167 c 1. initialisations
 168 c====
 169 c
 170 #include "impr01.h"
 171 c
 172 #ifdef _DEBUG_HOMARD_
 173       write (ulsort,texte(langue,1)) 'Entree', nompro
 174       call dmflsh (iaux)
 175 #endif
 176 c
 177       texte(1,4) =
 178      >'(/,''Tetr. en cours : numero a l''''iteration '',a3,'' : '',i10)'
 179       texte(1,5) =
 180      >'(  ''                  etat a l''''iteration '',a3,''   : '',i4)'
 181 c
 182       texte(2,4) =
 183      >'(/,''Current tetrahedron : # at iteration '',a3,''     : '',i10)'
 184       texte(2,5) =
 185      > '(  ''                     status at iteration '',a3,'' : '',i4)'
 186 c
 187       coderr = 0
 188 c
 189 c 1.2. ==> on repere les numeros dans le calcul pour ses quatre fils
 190 c          a l'iteration n
 191 c
 192       f1hn = anfite(tehn)
 193       f1cn = nteeca(f1hn)
 194       f2cn = nteeca(f1hn+1)
 195       f3cn = nteeca(f1hn+2)
 196       f4cn = nteeca(f1hn+3)
 197 c
 198       if ( prfcan(f1cn).gt.0 .and. prfcan(f2cn).gt.0 .and.
 199      >     prfcan(f3cn).gt.0 .and. prfcan(f4cn).gt.0 ) then
 200 c
 201 c====
 202 c 2. etan = 41, ..., 44 : le tetraedre etait coupe en 4
 203 c                         selon la face 1, 2, 3, 4
 204 c    etan = 45, 46, 47 : le tetraedre etait coupe en 4
 205 c                        selon une diagonale
 206 c====
 207 c
 208 c 2.1. ==> etanp1 = 0 : le tetraedre est reactive.
 209 c
 210       if ( etanp1.eq.0 ) then
 211 c
 212         tecnp1 = ntesca(tehnp1)
 213         prfcap(tecnp1) = 1
 214 c
 215         if ( typint.eq.0 ) then
 216           daux1 = unsqu
 217         else
 218           daux1 = 1.d0
 219         endif
 220         do 21 , nrofon = 1, nbfonc
 221           daux = daux1 * ( vafoen(nrofon,prfcan(f1cn))
 222      >                   + vafoen(nrofon,prfcan(f2cn))
 223      >                   + vafoen(nrofon,prfcan(f3cn))
 224      >                   + vafoen(nrofon,prfcan(f4cn)) )
 225           vafott(nrofon,tecnp1) = daux
 226    21   continue
 227 cgn        write (41,7777) tecnp1
 228 cgn        write (ulsort,7777) f1cn,f2cn,f3cn,f4cn,-1,tecnp1
 229 cgn7777   format(I3)
 230 c
 231 c 2.2. ==> etanp1 = 21, ...,26 : le tetraedre est decoupe en deux
 232 c             c'est ce qui se passe quand un decoupage de conformite
 233 c             est supprime au debut des algorithmes d'adaptation,
 234 c             puis reproduit a la creation du maillage car les faces
 235 c             autour n'ont pas change entre les deux iterations.
 236 c             on donne la valeur moyenne de la fonction sur les quatre
 237 c             anciens fils a chaque nouveau fils.
 238 c             remarque : on pourrait certainement faire mieux, avec des
 239 c                        moyennes ponderees en fonction du recouvrement
 240 c                        des anciennes et nouvelles filles. c'est trop
 241 c                        complique pour que cela vaille le coup.
 242 c             remarque : cela arrive seulement avec du deraffinement.
 243 c
 244       elseif ( etanp1.ge.21 .and. etanp1.le.26 ) then
 245 c
 246         f1hp = filtet(tehnp1)
 247         f1cp = ntesca(f1hp)
 248         f2cp = ntesca(f1hp+1)
 249         prfcap(f1cp) = 1
 250         prfcap(f2cp) = 1
 251         if ( typint.eq.0 ) then
 252           daux1 = unsqu
 253         else
 254           daux1 = unsde
 255         endif
 256         do 22 , nrofon = 1, nbfonc
 257           daux = daux1 * ( vafoen(nrofon,prfcan(f1cn))
 258      >                   + vafoen(nrofon,prfcan(f2cn))
 259      >                   + vafoen(nrofon,prfcan(f3cn))
 260      >                   + vafoen(nrofon,prfcan(f4cn)) )
 261           vafott(nrofon,f1cp) = daux
 262           vafott(nrofon,f2cp) = daux
 263    22   continue
 264 cgn        write(42,7777) f1cp,f2cp
 265 cgn        write(ulsort,7777) f1cn,f2cn,f3cn,f4cn,-1,
 266 cgn     >                     f1cp,f2cp
 267 c
 268 c 2.3. ==> etanp1 = etan : le tetraedre est decoupe en
 269 c                       quatre selon le meme decoupage.
 270 c             c'est ce qui se passe quand un decoupage de conformite
 271 c             est supprime au debut des algorithmes d'adaptation,
 272 c             puis reproduit a la creation du maillage car les faces
 273 c             autour n'ont pas change entre les deux iterations.
 274 c             le fils prend la valeur de la fonction sur l'ancien
 275 c             fils qui etait au meme endroit. comme la procedure de
 276 c             numerotation est la meme (voir cmcdte), le premier fils
 277 c             est toujours le meme, le second egalement. on prendra
 278 c             alors la valeur sur le fils de rang identique a
 279 c             l'iteration n.
 280 c
 281       elseif ( etanp1.eq.etan ) then
 282 c
 283         f1hp = filtet(tehnp1)
 284         f1cp = ntesca(f1hp)
 285         f2cp = ntesca(f1hp+1)
 286         f3cp = ntesca(f1hp+2)
 287         f4cp = ntesca(f1hp+3)
 288         prfcap(f1cp) = 1
 289         prfcap(f2cp) = 1
 290         prfcap(f3cp) = 1
 291         prfcap(f4cp) = 1
 292         do 23 , nrofon = 1, nbfonc
 293           vafott(nrofon,f1cp) = vafoen(nrofon,prfcan(f1cn))
 294           vafott(nrofon,f2cp) = vafoen(nrofon,prfcan(f2cn))
 295           vafott(nrofon,f3cp) = vafoen(nrofon,prfcan(f3cn))
 296           vafott(nrofon,f4cp) = vafoen(nrofon,prfcan(f4cn))
 297    23   continue
 298 cgn        write(43,7777) f1cp,f2cp,f3cp,f4cp
 299 cgn        write(ulsort,7777) f1cn,f2cn,f3cn,f4cn,-1,
 300 cgn     >                     f1cp,f2cp,f3cp,f4cp
 301 c
 302 c 2.4. ==> etanp1 = 41, ..., 47 et different de
 303 c                             etan : le tetraedre est decoupe en quatre
 304 c                             mais par un autre decoupage.
 305 c             c'est ce qui se passe quand un decoupage de conformite
 306 c             est supprime au debut des algorithmes d'adaptation. il y
 307 c             a du deraffinement dans la zone qui induisait le decoupage
 308 c             de conformite et raffinement sur une autre zone.
 309 c             on donne la valeur moyenne de la fonction sur les quatre
 310 c             anciens fils a chaque nouveau fils.
 311 c             remarque : on pourrait certainement faire mieux, avec des
 312 c                        moyennes ponderees en fonction du recouvrement
 313 c                        des anciens et nouveaux fils. c'est trop
 314 c                        complique pour que cela vaille le coup.
 315 c             remarque : cela arrive seulement avec du deraffinement.
 316 c
 317       elseif ( etanp1.ge.41 .and. etanp1.le.47 ) then
 318 c
 319         f1hp = filtet(tehnp1)
 320         f1cp = ntesca(f1hp)
 321         f2cp = ntesca(f1hp+1)
 322         f3cp = ntesca(f1hp+2)
 323         f4cp = ntesca(f1hp+3)
 324         prfcap(f1cp) = 1
 325         prfcap(f2cp) = 1
 326         prfcap(f3cp) = 1
 327         prfcap(f4cp) = 1
 328         do 24 , nrofon = 1, nbfonc
 329           daux = unsqu * ( vafoen(nrofon,prfcan(f1cn))
 330      >                   + vafoen(nrofon,prfcan(f2cn))
 331      >                   + vafoen(nrofon,prfcan(f3cn))
 332      >                   + vafoen(nrofon,prfcan(f4cn)) )
 333           vafott(nrofon,f1cp) = daux
 334           vafott(nrofon,f2cp) = daux
 335           vafott(nrofon,f3cp) = daux
 336           vafott(nrofon,f4cp) = daux
 337    24   continue
 338 cgn        write(44,7777) f1cp,f2cp,f3cp,f4cp
 339 cgn        write(ulsort,7777) f1cn,f2cn,f3cn,f4cn,-1,
 340 cgn     >                     f1cp,f2cp,f3cp,f4cp
 341 c
 342 c 2.5. ==> etanp1 = 85, 86 ou 87 : le tetraedre est
 343 c                      decoupe en huit par une diagonale.
 344 c             c'est ce qui se passe quand un decoupage de conformite
 345 c             est supprime au debut des algorithmes d'adaptation, puis
 346 c             remis differement car l'environnement a change.
 347 c             on donne la valeur moyenne de la fonction sur les quatre
 348 c             anciens fils a chaque nouveau fils.
 349 c             attention : il est possible que les fils sur les bords
 350 c                         soient decoupes par de la conformite. Il faut
 351 c                         alors transmettre la valeur a leurs 2 ou 4
 352 c                         fils.
 353 c             attention : ce n'est pas comme en 2D ; il faut examiner
 354 c                         tous les fils, car par contamination de faces
 355 c                         coupees en 2, les fils centraux peuvent etre
 356 c                         decoupes.
 357 c             remarque : on pourrait certainement faire mieux, avec des
 358 c                        moyennes ponderees en fonction du recouvrement
 359 c                        des anciens et nouveaux fils. c'est trop
 360 c                        complique pour que cela vaille le coup.
 361 c
 362       elseif ( etanp1.ge.85 .and. etanp1.le.87 ) then
 363 c
 364         f1hp = filtet(tehnp1)
 365         if ( typint.eq.0 ) then
 366 c
 367         lglist = 0
 368           do 250 , nrlist = 1 , 8
 369             fihp = f1hp+nrlist-1
 370             iaux = mod(hettet(fihp),100)
 371             if ( iaux.eq.0 ) then
 372               lglist = lglist + 1
 373               list(lglist) = ntesca(fihp)
 374             elseif ( iaux.ge.21 .and. iaux.le.26 ) then
 375               lglist = lglist + 1
 376               list(lglist) = ntesca(filtet(fihp))
 377               lglist = lglist + 1
 378               list(lglist) = ntesca(filtet(fihp)+1)
 379             elseif ( iaux.ge.41 .and. iaux.le.47 ) then
 380               lglist = lglist + 1
 381               list(lglist) = ntesca(filtet(fihp))
 382               lglist = lglist + 1
 383               list(lglist) = ntesca(filtet(fihp)+1)
 384               lglist = lglist + 1
 385               list(lglist) = ntesca(filtet(fihp)+2)
 386               lglist = lglist + 1
 387               list(lglist) = ntesca(filtet(fihp)+3)
 388             else
 389               coderr = 1
 390             endif
 391   250     continue
 392 c
 393           do 260 , nrlist = 1 , lglist
 394             prfcap(list(nrlist)) = 1
 395   260     continue
 396 c
 397           do 270 , nrofon = 1, nbfonc
 398             daux = unsqu * ( vafoen(nrofon,prfcan(f1cn))
 399      >                     + vafoen(nrofon,prfcan(f2cn))
 400      >                     + vafoen(nrofon,prfcan(f3cn))
 401      >                     + vafoen(nrofon,prfcan(f4cn)) )
 402             do 2701 , nrlist = 1 , lglist
 403               vafott(nrofon,list(nrlist)) = daux
 404  2701       continue
 405   270     continue
 406 c
 407         else
 408 c
 409           do 251 , nrlist = 1 , 8
 410 c
 411             fihp = f1hp+nrlist-1
 412             iaux = mod(hettet(fihp),100)
 413             if ( iaux.eq.0 ) then
 414               f1cp = ntesca(fihp)
 415               prfcap(f1cp) = 1
 416               do 2511 , nrofon = 1, nbfonc
 417                 daux = unshu * ( vafoen(nrofon,prfcan(f1cn))
 418      >                         + vafoen(nrofon,prfcan(f2cn))
 419      >                         + vafoen(nrofon,prfcan(f3cn))
 420      >                         + vafoen(nrofon,prfcan(f4cn)) )
 421                 vafott(nrofon,f1cp) = daux
 422  2511         continue
 423             elseif ( iaux.ge.21 .and. iaux.le.26 ) then
 424               f1fcp = ntesca(filtet(fihp))
 425               f2fcp = ntesca(filtet(fihp)+1)
 426               prfcap(f1fcp) = 1
 427               prfcap(f2fcp) = 1
 428               do 2512 , nrofon = 1, nbfonc
 429                 daux = unssz * ( vafoen(nrofon,prfcan(f1cn))
 430      >                         + vafoen(nrofon,prfcan(f2cn))
 431      >                         + vafoen(nrofon,prfcan(f3cn))
 432      >                         + vafoen(nrofon,prfcan(f4cn)) )
 433                 vafott(nrofon,f1fcp) = daux
 434                 vafott(nrofon,f2fcp) = daux
 435  2512         continue
 436             elseif ( iaux.ge.41 .and. iaux.le.47 ) then
 437               f1fcp = ntesca(filtet(fihp))
 438               f2fcp = ntesca(filtet(fihp)+1)
 439               f3fcp = ntesca(filtet(fihp)+2)
 440               f4fcp = ntesca(filtet(fihp)+3)
 441               prfcap(f1fcp) = 1
 442               prfcap(f2fcp) = 1
 443               prfcap(f3fcp) = 1
 444               prfcap(f4fcp) = 1
 445               do 2513 , nrofon = 1, nbfonc
 446                 daux = unstr2 * ( vafoen(nrofon,prfcan(f1cn))
 447      >                          + vafoen(nrofon,prfcan(f2cn))
 448      >                          + vafoen(nrofon,prfcan(f3cn))
 449      >                          + vafoen(nrofon,prfcan(f4cn)) )
 450                 vafott(nrofon,f1fcp) = daux
 451                 vafott(nrofon,f2fcp) = daux
 452                 vafott(nrofon,f3fcp) = daux
 453                 vafott(nrofon,f4fcp) = daux
 454  2513         continue
 455             else
 456               coderr = 1
 457             endif
 458 c
 459   251     continue
 460 c
 461         endif
 462 cgn        write(46,7777) (list(nrlist),nrlist = 1 , lglist)
 463 cgn        write(ulsort,7777) f1cn,f2cn,f3cn,f4cn,-1,
 464 cgn     >                   (list(nrlist),nrlist = 1 , lglist)
 465 c
 466 c 2.6. ==> aucun autre etat sur le tetraedre courant
 467 c                       n'est possible
 468 c
 469       else
 470 c
 471         coderr = 1
 472         write (ulsort,texte(langue,4)) 'n  ', tehn
 473         write (ulsort,texte(langue,5)) 'n  ', etan
 474         write (ulsort,texte(langue,4)) 'n+1', tehnp1
 475         write (ulsort,texte(langue,5)) 'n+1', etanp1
 476 c
 477       endif
 478 c
 479       endif
 480 c
 481 c====
 482 c 3. la fin
 483 c====
 484 c
 485       if ( coderr.ne.0 ) then
 486 c
 487         write (ulsort,texte(langue,1)) 'Sortie', nompro
 488       write (ulsort,texte(langue,2)) coderr
 489       codret = codret + 1
 490 c
 491       endif
 492 c
 493       end