src/tool/AP_Conversion/pcequ1.F

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