src/tool/AP_Conversion/pcsqug.F

   1       subroutine pcsqug ( nbfonc, ngauss, nbnorf, typgeo, deraff,
   2      >                    prfcan, prfcap,
   3      >                    hetqua, ancqua,
   4      >                    filqua,
   5      >                    nbanqu, anfiqu, anhequ,
   6      >                    nqueca, nqusca,
   7      >                    ntreca, ntrsca,
   8      >                    vafoen, vafott,
   9      >                    conorf, copgrf, wipg,
  10      >                    ulsort, langue, codret )
  11 c ______________________________________________________________________
  12 c
  13 c                             H O M A R D
  14 c
  15 c Outil de Maillage Adaptatif par Raffinement et Deraffinement d'EDF R&D
  16 c
  17 c Version originale enregistree le 18 juin 1996 sous le numero 96036
  18 c aupres des huissiers de justice Simart et Lavoir a Clamart
  19 c Version 11.2 enregistree le 13 fevrier 2015 sous le numero 2015/014
  20 c aupres des huissiers de justice
  21 c Lavoir, Silinski & Cherqui-Abrahmi a Clamart
  22 c
  23 c    HOMARD est une marque deposee d'Electricite de France
  24 c
  25 c Copyright EDF 1996
  26 c Copyright EDF 1998
  27 c Copyright EDF 2002
  28 c Copyright EDF 2020
  29 c ______________________________________________________________________
  30 c
  31 c    aPres adaptation - Conversion de Solution -
  32 c     -                 -             -
  33 c                       QUadrangles a plusieurs points de Gauss
  34 c                       --                                -
  35 c ______________________________________________________________________
  36 c .        .     .        .                                            .
  37 c .  nom   . e/s . taille .           description                      .
  38 c .____________________________________________________________________.
  39 c . nbfonc . e   .    1   . nombre de fonctions aux points de Gauss    .
  40 c . ngauss . e   .   1    . nbre de points de Gauss des fonctions pg   .
  41 c . nbnorf . e   .   1    . nbre de noeuds de l'element de reference   .
  42 c . typgeo . e   .   1    . type geometrique au sens MED               .
  43 c . deraff . e   .    1   . vrai, s'il y a eu du deraffinement en      .
  44 c .        .     .        . passant de l'iteration n a n+1 ; faux sinon.
  45 c . prfcan . e   .   *    . En numero du calcul a l'iteration n :      .
  46 c .        .     .        . 0 : l'entite est absente du profil         .
  47 c .        .     .        . i : l'entite est au rang i dans le profil  .
  48 c . prfcap . es  .   *    . En numero du calcul a l'iteration n+1 :    .
  49 c .        .     .        . 0 : l'entite est absente du profil         .
  50 c .        .     .        . 1 : l'entite est presente dans le profil   .
  51 c . hetqua . e   . nbquto . historique de l'etat des quadrangles       .
  52 c . filqua . e   . nbquto . premier fils des quadrangles               .
  53 c . nbanqu . e   .   1    . nombre de quadrangles decoupes par         .
  54 c .        .     .        . conformite sur le maillage avant adaptation.
  55 c . anfiqu . e   . nbanqu . tableau filqua du maillage de l'iteration n
  56 c . anhequ . e   . nbanqu . tableau hetqua du maillage de l'iteration n
  57 c . nqueca . e   .   *    . nro des quadrangles dans le calcul en ent. .
  58 c . nqusca . e   . rsquto . numero des quadrangles du calcul           .
  59 c . ntreca . e   .   *    . nro des triangles dans le calcul en entree .
  60 c . ntrsca . e   . rstrto . numero des triangles du calcul             .
  61 c . vafoen .   e . nbfonc*. variables en entree de l'adaptation        .
  62 c .        .     .    *   .                                            .
  63 c . vafott .   a . nbfonc*. tableau temporaire de la solution          .
  64 c .        .     .    *   .                                            .
  65 c . conorf . e   .  sdim* . coordonnees des noeuds de l'element de     .
  66 c .        .     . nbnorf . reference                                  .
  67 c . copgrf . e   .  sdim* . coordonnees des points de Gauss            .
  68 c .        .     . ngauss . de l'element de reference                  .
  69 c . wipg   .  a  . nbnorf*. fonctions de forme exprimees aux points de .
  70 c .        .     . ngauss . Gauss                                      .
  71 c . ulsort . e   .   1    . numero d'unite logique de la liste standard.
  72 c . langue . e   .    1   . langue des messages                        .
  73 c .        .     .        . 1 : francais, 2 : anglais                  .
  74 c . codret . es  .    1   . code de retour des modules                 .
  75 c .        .     .        . 0 : pas de probleme                        .
  76 c .        .     .        . 1 : probleme                               .
  77 c ______________________________________________________________________
  78 c
  79 c====
  80 c 0. declarations et dimensionnement
  81 c====
  82 c
  83 c 0.1. ==> generalites
  84 c
  85       implicit none
  86       save
  87 c
  88       character*6 nompro
  89       parameter ( nompro = 'PCSQUG' )
  90 c
  91 #include "nblang.h"
  92 #include "fractc.h"
  93 c
  94 c 0.2. ==> communs
  95 c
  96 #include "envex1.h"
  97 c
  98 #include "envca1.h"
  99 #include "nombqu.h"
 100 #include "nombsr.h"
 101 #include "nomber.h"
 102 c
 103 c 0.3. ==> arguments
 104 c
 105       integer nbfonc
 106       integer ngauss, nbnorf, typgeo
 107       integer prfcan(*), prfcap(*)
 108       integer hetqua(nbquto), ancqua(*)
 109       integer filqua(nbquto)
 110       integer nbanqu, anfiqu(nbanqu), anhequ(nbanqu)
 111       integer nqueca(requto), nqusca(rsquto)
 112       integer ntreca(retrto), ntrsca(rstrto)
 113 c
 114       double precision vafoen(nbfonc,ngauss,*)
 115       double precision vafott(nbfonc,ngauss,*)
 116       double precision conorf(sdim,nbnorf), copgrf(sdim,ngauss)
 117       double precision wipg(nbnorf,ngauss)
 118 c
 119       logical deraff
 120 c
 121       integer ulsort, langue, codret
 122 c
 123 c 0.4. ==> variables locales
 124 c
 125       integer iaux
 126 c
 127 c     qucn   = QUadrangle courant en numerotation Calcul a l'it. N
 128 c     qucnp1 = QUadrangle courant en numerotation Calcul a l'it. N+1
 129 c     quhn   = QUadrangle courant en numerotation Homard a l'it. N
 130 c     quhnp1 = QUadrangle courant en numerotation Homard a l'it. N+1
 131 c
 132       integer qucn, qucnp1, quhn, quhnp1
 133 c
 134 c     f1hp = Fils 1er du quadrangle en numerotation Homard a l'it. N+1
 135 c     f1cp = Fils 1er du quadrangle en numerotation Calcul a l'it. N+1
 136 c     f2cp = Fils 2eme du quadrangle en numerotation Calcul a l'it. N+1
 137 c     f3cp = Fils 3eme du quadrangle en numerotation Calcul a l'it. N+1
 138 c     f4cp = Fils 4eme du quadrangle en numerotation Calcul a l'it. N+1
 139 c
 140       integer f1hp
 141       integer f1cp, f2cp, f3cp, f4cp
 142 c
 143 c     f1hn = Fils 1er du quadrangle en numerotation Homard a l'it. N
 144 c     f1cn = Fils 1er du quadrangle en numerotation Calcul a l'it. N
 145 c     f2cn = Fils 2eme du quadrangle en numerotation Calcul a l'it. N
 146 c     f3cn = Fils 3eme du quadrangle en numerotation Calcul a l'it. N
 147 c     f4cn = Fils 4eme du quadrangle en numerotation Calcul a l'it. N
 148 c
 149       integer f1hn
 150       integer f1cn, f2cn, f3cn, f4cn
 151 c
 152 c     etan   = ETAt du quadrangle a l'iteration N
 153 c     etanp1 = ETAt du quadrangle a l'iteration N+1
 154 c
 155       integer etan, etanp1
 156 c
 157       integer nrofon, nugaus
 158 c
 159       double precision daux
 160 c
 161       integer nbmess
 162       parameter ( nbmess = 10 )
 163       character*80 texte(nblang,nbmess)
 164 c
 165 c 0.5. ==> initialisations
 166 c ______________________________________________________________________
 167 c
 168 c====
 169 c 1. initialisations
 170 c====
 171 c
 172 #include "pcimp0.h"
 173 #include "impr01.h"
 174 c
 175 #ifdef _DEBUG_HOMARD_
 176       write (ulsort,texte(langue,1)) 'Entree', nompro
 177       call dmflsh (iaux)
 178       write(ulsort,*) 'nbfonc, ngauss, nbquto = ',nbfonc, ngauss, nbquto
 179 #endif
 180 c
 181       texte(1,4) =
 182      >'(/,''Quad. en cours : nro a l''''iteration '',a3,'' :'',i10)'
 183       texte(1,5) =
 184      > '(  ''                etat a l''''iteration '',a3,''   : '',i4)'
 185       texte(1,6) =
 186      >'(/,''Quadrangle decoupe en triangles : on ne sait pas faire.'')'
 187 c
 188       texte(2,4) =
 189      >'(/,''Current quadrangle : # at iteration '',a3,''     : '',i10)'
 190       texte(2,5) =
 191      > '(  ''              status at iteration '',a3,'' : '',i4)'
 192       texte(2,6) =
 193      >'(/,''Quadrangle cut into triangles : not available.'')'
 194 c
 195       codret = 0
 196 c
 197 c====
 198 c 2. calcul des valeurs des fonctions de forme aux points de Gauss
 199 c    Inutile pour le moment
 200 c====
 201 c
 202 c====
 203 c 3. on boucle sur tous les quadrangles du maillage HOMARD n+1
 204 c    on trie en fonction de l'etat du quadrangle dans le maillage n
 205 c    on numerote les paragraphes en fonction de la documentation, a
 206 c    savoir : le paragraphe doc.n.p traite de la mise a jour de solution
 207 c    pour un quadrangle dont l'etat est passe de n a p.
 208 c    les autres paragraphes sont numerotes classiquement
 209 c====
 210 c
 211 #ifdef _DEBUG_HOMARD_
 212       write (ulsort,*) '3. Boucle ; codret = ', codret
 213 #endif
 214 c
 215       if ( codret.eq.0 ) then
 216 c
 217       if ( nbfonc.ne.0 ) then
 218 c
 219       do 30 , quhnp1 = 1 , nbquto
 220 c
 221 c 3.1. ==> caracteristiques du quadrangle :
 222 c 3.1.1. ==> son numero homard dans le maillage precedent
 223 c
 224         if ( deraff ) then
 225           quhn = ancqua(quhnp1)
 226         else
 227           quhn = quhnp1
 228         endif
 229 c
 230 c 3.1.2. ==> l'historique de son etat
 231 c          On rappelle que l'etat vaut :
 232 c          etan = 0 : le quadrangle etait actif
 233 c          etan = 4 : le quadrangle etait coupe en 4 ; il y a eu
 234 c                     deraffinement.
 235 c          etan = 55 : le quadrangle n'existait pas ; il a ete produit
 236 c                      par un decoupage.
 237 c          etan = 31, 32, 33, 34 : le quadrangle etait coupe en 3
 238 c                                  quadrangles ; il y a eu deraffinement.
 239 c
 240         etanp1 = mod(hetqua(quhnp1),100)
 241         etan   = (hetqua(quhnp1)-etanp1) / 100
 242 c
 243 cgn        write (ulsort,1792) 'Quadrangle', quhn, etan, quhnp1, etanp1
 244 c
 245 c=======================================================================
 246 c doc.0.p. ==> etan = 0 : le quadrangle etait actif
 247 c=======================================================================
 248 c
 249         if ( etan.eq.0 ) then
 250 c
 251 c         on repere son ancien numero dans le calcul
 252 c
 253           qucn = nqueca(quhn)
 254 c
 255           if ( prfcan(qucn).gt.0 ) then
 256 c
 257 cgn          print 1789,(vafoen(nrofon,nugaus,qucn), nugaus = 1 , ngauss)
 258 cgn 1789 format('  Valeurs anciennes  :  ',5g12.5)
 259 c
 260 c doc.0.0. ===> etanp1 = 0 : le quadrangle etait actif et l'est encore ;
 261 c               il est inchange
 262 c             c'est le cas le plus simple : on prend la valeur de la
 263 c             fonction associee a l'ancien numero du quadrangle.
 264 c
 265 c           .................         .................
 266 c           .               .         .               .
 267 c           .               .         .               .
 268 c           .               .         .               .
 269 c           .               .  ===>   .               .
 270 c           .               .         .               .
 271 c           .               .         .               .
 272 c           .               .         .               .
 273 c           .................         .................
 274 c
 275           if ( etanp1.eq.0 ) then
 276 c
 277             qucnp1 = nqusca(quhnp1)
 278             prfcap(qucnp1) = 1
 279 c
 280             do 321 , nrofon = 1 , nbfonc
 281 cgn      write(ulsort,7778)
 282 cgn     > (vafoen(nrofon,nugaus,prfcan(qucn)),nugaus=1,ngauss)
 283               do 3211 , nugaus = 1 , ngauss
 284                 vafott(nrofon,nugaus,qucnp1) =
 285      >                                vafoen(nrofon,nugaus,prfcan(qucn))
 286  3211         continue
 287   321       continue
 288 cgn        write(31,7777) qucnp1
 289 cgn        write(ulsort,7777) qucn,-1,qucnp1
 290 cgn7777   format(I3)
 291 cgn7778   format(8g14.7)
 292 c
 293 c doc.0.1/2/3 ==> etanp1 = 31, 32, 33 ou 34 : le quadrangle etait actif
 294 c                 et est decoupe en 3 triangles.
 295 c                 les trois fils prennent la valeur de la fonction sur
 296 c                 le pere
 297 c           .................         .................
 298 c           .               .         .      . .      .
 299 c           .               .         .     .   .     .
 300 c           .               .         .    .     .    .
 301 c           .               .  ===>   .   .       .   .
 302 c           .               .         .  .         .  .
 303 c           .               .         . .           . .
 304 c           .               .         ..             ..
 305 c           .................         .................
 306 c
 307           elseif ( etanp1.ge.31 .and. etanp1.le.34 ) then
 308 c
 309             f1hp = -filqua(quhnp1)
 310             f1cp = nqusca(f1hp)
 311             f2cp = nqusca(f1hp+1)
 312             f3cp = nqusca(f1hp+2)
 313             prfcap(f1cp) = 1
 314             prfcap(f2cp) = 1
 315             prfcap(f3cp) = 1
 316             write (ulsort,texte(langue,6))
 317             codret = codret + 1
 318 c
 319 c doc.0.4/6/7/8. ==> etanp1 = 4 : le quadrangle etait actif et
 320 c                    est decoupe en 4.
 321 c                    les quaque fils prennent la valeur de la fonction
 322 c                    sur le pere
 323 c           .................         .................
 324 c           .               .         .       .       .
 325 c           .               .         .       .       .
 326 c           .               .         .       .       .
 327 c           .               .  ===>   .................
 328 c           .               .         .       .       .
 329 c           .               .         .       .       .
 330 c           .               .         .       .       .
 331 c           .................         .................
 332 c
 333           elseif ( etanp1.eq.4 ) then
 334 c
 335             f1hp = filqua(quhnp1)
 336             f1cp = nqusca(f1hp)
 337             f2cp = nqusca(f1hp+1)
 338             f3cp = nqusca(f1hp+2)
 339             f4cp = nqusca(f1hp+3)
 340             prfcap(f1cp) = 1
 341             prfcap(f2cp) = 1
 342             prfcap(f3cp) = 1
 343             prfcap(f4cp) = 1
 344             do 323 , nrofon = 1 , nbfonc
 345               do 3231 , nugaus = 1 , ngauss
 346                 daux = vafoen(nrofon,nugaus,prfcan(qucn))
 347                 vafott(nrofon,nugaus,f1cp) = daux
 348                 vafott(nrofon,nugaus,f2cp) = daux
 349                 vafott(nrofon,nugaus,f3cp) = daux
 350                 vafott(nrofon,nugaus,f4cp) = daux
 351  3231         continue
 352   323       continue
 353 c
 354 c doc.0.erreur. ==> aucun autre etat sur le quadrangle courant n'est
 355 c                   possible
 356 c
 357           else
 358 c
 359             codret = codret + 1
 360             write (ulsort,texte(langue,4)) 'n  ', quhn
 361             write (ulsort,texte(langue,5)) 'n  ', etan
 362             write (ulsort,texte(langue,4)) 'n+1', quhnp1
 363             write (ulsort,texte(langue,5)) 'n+1', etanp1
 364 c
 365           endif
 366 c
 367           endif
 368 c
 369 c=======================================================================
 370 c doc.1/2/3.p. ==> etan = 31, 32, 33 ou 34 : le quadrangle etait coupe
 371 c                  en 3 triangles
 372 c=======================================================================
 373 c
 374         elseif ( etan.ge.31 .and. etan.le.34 ) then
 375 c
 376 c         on repere les numeros dans le calcul pour ses trois fils a
 377 c         l'iteration n
 378 c
 379           f1hn = -anfiqu(quhn)
 380           f1cn = ntreca(f1hn)
 381           f2cn = ntreca(f1hn+1)
 382           f3cn = ntreca(f1hn+2)
 383 c
 384           if ( prfcan(f1cn).gt.0 .and. prfcan(f2cn).gt.0 .and.
 385      >         prfcan(f3cn).gt.0 ) then
 386 c
 387           write (ulsort,texte(langue,6))
 388           codret = codret + 1
 389 c
 390 c doc.1/2/3.0. ===> etanp1 = 0 : le quadrangle est actif. il est
 391 c                   reactive.
 392 c                   on lui attribue la valeur moyenne sur les trois
 393 c                   anciens fils.
 394 c        remarque : cela arrive seulement avec du deraffinement.
 395 c           .................         .................
 396 c           .      . .      .         .               .
 397 c           .     .   .     .         .               .
 398 c           .    .     .    .         .               .
 399 c           .   .       .   .  ===>   .               .
 400 c           .  .         .  .         .               .
 401 c           . .           . .         .               .
 402 c           ..             ..         .               .
 403 c           .................         .................
 404 c
 405           if ( etanp1.eq.0 ) then
 406 c
 407             qucnp1 = nqusca(quhnp1)
 408             prfcap(qucnp1) = 1
 409 c
 410 c doc.1/2/3.1/2/3. ===> etanp1 = etan : le quadrangle est decoupe en
 411 c                       trois selon le meme decoupage.
 412 c             c'est ce qui se passe quand un decoupage de conformite
 413 c             est supprime au debut des algorithmes d'adaptation,
 414 c             puis reproduit a la creation du maillage car les
 415 c             quadrangles autour n'ont pas change entre les deux
 416 c             iterations.
 417 c             le fils prend la valeur de la fonction sur l'ancien
 418 c             fils qui etait au meme endroit. comme la procedure de
 419 c             numerotation est la meme (voir cmcdqu), le premier fils
 420 c             est toujours le meme, le 2eme et le 3eme egalement.
 421 c             on prendra alors la valeur sur le fils de rang identique
 422 c             a l'iteration n.
 423 c           .................         .................
 424 c           .      . .      .         .      . .      .
 425 c           .     .   .     .         .     .   .     .
 426 c           .    .     .    .         .    .     .    .
 427 c           .   .       .   .  ===>   .   .       .   .
 428 c           .  .         .  .         .  .         .  .
 429 c           . .           . .         . .           . .
 430 c           ..             ..         ..             ..
 431 c           .................         .................
 432 c
 433           elseif ( etanp1.eq.etan ) then
 434 c
 435             f1hp = -filqua(quhnp1)
 436             f1cp = ntrsca(f1hp)
 437             f2cp = ntrsca(f1hp+1)
 438             f3cp = ntrsca(f1hp+2)
 439             prfcap(f1cp) = 1
 440             prfcap(f2cp) = 1
 441             prfcap(f3cp) = 1
 442 c
 443 c doc.1/2/3.perm(1/2/3). ===> etanp1 = 31, 32, 33 ou 34 et different de
 444 c                             etan : le quadrangle est encore decoupe
 445 c                             en trois, mais par un autre decoupage.
 446 c             c'est ce qui se passe quand un decoupage de conformite
 447 c             est supprime au debut des algorithmes d'adaptation. il y
 448 c             a du deraffinement dans la zone qui induisait le decoupage
 449 c             de conformite et raffinement sur une autre zone.
 450 c             On donne la valeur moyenne de la fonction sur les trois
 451 c             anciens fils a chaque nouveau fils.
 452 c             remarque : on pourrait certainement faire mieux, avec des
 453 c                        moyennes ponderees en fonction du recouvrement
 454 c                        des anciens et nouveaux fils. c'est trop
 455 c                        complique pour que cela vaille le coup.
 456 c             remarque : cela arrive seulement avec du deraffinement.
 457 c
 458 c           .................         .................
 459 c           .      . .      .         ..             ..
 460 c           .     .   .     .         . .           . .
 461 c           .    .     .    .         .  .         .  .
 462 c           .   .       .   .  ===>   .   .       .   .
 463 c           .  .         .  .         .    .     .    .
 464 c           . .           . .         .     .   .     .
 465 c           ..             ..         .      . .      .
 466 c           .................         .................
 467 c
 468           elseif ( etanp1.ge.31 .and. etanp1.le.34 ) then
 469 c
 470             f1hp = -filqua(quhnp1)
 471             f1cp = ntrsca(f1hp)
 472             f2cp = ntrsca(f1hp+1)
 473             f3cp = ntrsca(f1hp+2)
 474             prfcap(f1cp) = 1
 475             prfcap(f2cp) = 1
 476             prfcap(f3cp) = 1
 477 c
 478 c doc.1/2/3.4/6/7/8. ===> etanp1 = 4 : le quadrangle est
 479 c                         decoupe en quatre.
 480 c             c'est ce qui se passe quand un decoupage de conformite
 481 c             est supprime au debut des algorithmes d'adaptation. il y
 482 c             a ensuite raffinement du quadrangle. qui plus est, par
 483 c             suite de la regle des ecarts de niveau, on peut avoir
 484 c             induit un decoupage de conformite sur l'un des fils.
 485 c
 486 c             On donne la valeur moyenne de la fonction sur les trois
 487 c             anciens fils a chaque nouveau fils.
 488 c             remarque : on pourrait certainement faire mieux, avec des
 489 c                        moyennes ponderees en fonction du recouvrement
 490 c                        des anciens et nouveaux fils. c'est trop
 491 c                        complique pour que cela vaille le coup.
 492 c
 493 c           .................         .................
 494 c           .      . .      .         .       .       .
 495 c           .     .   .     .         .       .       .
 496 c           .    .     .    .         .       .       .
 497 c           .   .       .   .  ===>   .................
 498 c           .  .         .  .         .       .       .
 499 c           . .           . .         .       .       .
 500 c           ..             ..         .       .       .
 501 c           .................         .................
 502 c
 503 c
 504           elseif ( etanp1.eq.4 ) then
 505 c
 506             f1hp = filqua(quhnp1)
 507             f1cp = nqusca(f1hp)
 508             f2cp = nqusca(f1hp+1)
 509             f3cp = nqusca(f1hp+2)
 510             f4cp = nqusca(f1hp+3)
 511             prfcap(f1cp) = 1
 512             prfcap(f2cp) = 1
 513             prfcap(f3cp) = 1
 514             prfcap(f4cp) = 1
 515 c
 516 c doc.1/2/3.erreur. ==> aucun autre etat sur le quadrangle courant
 517 c                       n'est possible
 518 c
 519           else
 520 c
 521             codret = codret + 1
 522             write (ulsort,texte(langue,4)) 'n  ', quhn
 523             write (ulsort,texte(langue,5)) 'n  ', etan
 524             write (ulsort,texte(langue,4)) 'n+1', quhnp1
 525             write (ulsort,texte(langue,5)) 'n+1', etanp1
 526 c
 527           endif
 528 c
 529           endif
 530 c
 531 c=======================================================================
 532 c doc.4. ==> le quadrangle etait coupe en 4 :
 533 c=======================================================================
 534 c
 535         elseif ( etan.eq.4 ) then
 536 c
 537 c         on repere les numeros dans le calcul pour ses quatre fils
 538 c         a l'iteration n
 539 c
 540           f1hn = anfiqu(quhn)
 541           f1cn = nqueca(f1hn)
 542           f2cn = nqueca(f1hn+1)
 543           f3cn = nqueca(f1hn+2)
 544           f4cn = nqueca(f1hn+3)
 545 c
 546           if ( prfcan(f1cn).gt.0 .and. prfcan(f2cn).gt.0 .and.
 547      >         prfcan(f3cn).gt.0 .and. prfcan(f4cn).gt.0 ) then
 548 c
 549 c doc.4.0. ===> etanp1 = 0 : le quadrangle est actif ; il est reactive.
 550 c             on lui attribue la valeur moyenne sur les quatre anciens
 551 c             fils.
 552 c             remarque : cela arrive seulement avec du deraffinement.
 553 c           .................         .................
 554 c           .       .       .         .               .
 555 c           .       .       .         .               .
 556 c           .       .       .         .               .
 557 c           .................  ===>   .               .
 558 c           .       .       .         .               .
 559 c           .       .       .         .               .
 560 c           .       .       .         .               .
 561 c           .................         .................
 562 c
 563           if ( etanp1.eq.0 ) then
 564 c
 565             qucnp1 = nqusca(quhnp1)
 566             prfcap(qucnp1) = 1
 567 c
 568             do 341 , nrofon = 1 , nbfonc
 569               do 3411 , nugaus = 1 , ngauss
 570                 vafott(nrofon,nugaus,qucnp1) =
 571      >                   unsqu * ( vafoen(nrofon,nugaus,prfcan(f1cn))
 572      >                           + vafoen(nrofon,nugaus,prfcan(f2cn))
 573      >                           + vafoen(nrofon,nugaus,prfcan(f3cn))
 574      >                           + vafoen(nrofon,nugaus,prfcan(f4cn)) )
 575  3411         continue
 576   341       continue
 577 c
 578 c doc.4.1/2/3. ===> etanp1 = 31, 32, 33 ou 34 : le quadrangle est
 579 c             decoupe en trois.
 580 c             on attribue la valeur moyenne sur les quatre anciens
 581 c             fils a chacune des trois nouveaux fils.
 582 c             remarque : on pourrait certainement faire mieux, avec des
 583 c                        moyennes ponderees en fonction du recouvrement
 584 c                        des anciens et nouveaux fils. c'est trop
 585 c                        complique pour que cela vaille le coup.
 586 c             remarque : cela arrive seulement avec du deraffinement.
 587 c           .................         .................
 588 c           .       .       .         .      . .      .
 589 c           .       .       .         .     .   .     .
 590 c           .       .       .         .    .     .    .
 591 c           .................  ===>   .   .       .   .
 592 c           .       .       .         .  .         .  .
 593 c           .       .       .         . .           . .
 594 c           .       .       .         ..             ..
 595 c           .................         .................
 596 c
 597           elseif ( etanp1.ge.31 .and. etanp1.le.34 ) then
 598 c
 599             f1hp = -filqua(quhnp1)
 600             f1cp = ntrsca(f1hp)
 601             f2cp = ntrsca(f1hp+1)
 602             f3cp = ntrsca(f1hp+2)
 603             prfcap(f1cp) = 1
 604             prfcap(f2cp) = 1
 605             prfcap(f3cp) = 1
 606             write (ulsort,texte(langue,6))
 607             codret = codret + 1
 608 c
 609           endif
 610 c
 611           endif
 612 c
 613 c=======================================================================
 614 c
 615         endif
 616 c
 617    30 continue
 618 c
 619       endif
 620 c
 621       endif
 622 c
 623 c====
 624 c 4. la fin
 625 c====
 626 c
 627 cgn      do 922 , nrofon = 1 , nbfonc
 628 cgn        print *,'fonction numero ', nrofon
 629 cgn        iaux = 0
 630 cgn        do 9222 , quhnp1 = 1 , nbtrto
 631 cgn           if ( mod(hettri(quhnp1),100).eq.0 ) then
 632 cgn           iaux = iaux+1
 633 cgn          print 1788,quhnp1,
 634 cgn     >     (vafott(nrofon,nugaus,iaux), nugaus = 1 , ngauss)
 635 cgn           endif
 636 cgn 9222   continue
 637 cgn  922 continue
 638       if ( codret.ne.0 ) then
 639 c
 640 #include "envex2.h"
 641 c
 642       write (ulsort,texte(langue,1)) 'Sortie', nompro
 643       write (ulsort,texte(langue,2)) codret
 644 c
 645       endif
 646 c
 647 #ifdef _DEBUG_HOMARD_
 648       write (ulsort,texte(langue,1)) 'Sortie', nompro
 649       call dmflsh (iaux)
 650 #endif
 651 c
 652       end