src/tool/AP_Conversion/pcs2h5.F

   1       subroutine pcs2h5 ( nbfop2, profho, vap2ho,
   2      >                    somare, np2are,
   3      >                    hetqua,
   4      >                    quahex,
   5      >                    lehexa, listso, listno,
   6      >                    areint,
   7      >                    ulsort, langue, codret )
   8 c ______________________________________________________________________
   9 c
  10 c                             H O M A R D
  11 c
  12 c Outil de Maillage Adaptatif par Raffinement et Deraffinement d'EDF R&D
  13 c
  14 c Version originale enregistree le 18 juin 1996 sous le numero 96036
  15 c aupres des huissiers de justice Simart et Lavoir a Clamart
  16 c Version 11.2 enregistree le 13 fevrier 2015 sous le numero 2015/014
  17 c aupres des huissiers de justice
  18 c Lavoir, Silinski & Cherqui-Abrahmi a Clamart
  19 c
  20 c    HOMARD est une marque deposee d'Electricite de France
  21 c
  22 c Copyright EDF 1996
  23 c Copyright EDF 1998
  24 c Copyright EDF 2002
  25 c Copyright EDF 2020
  26 c ______________________________________________________________________
  27 c
  28 c    aPres adaptation - Conversion de Solution -
  29 c     -                 -             -
  30 c    interpolation p2 sur les noeuds - decoupage Hexaedres - 5
  31 c                   -                            -           -
  32 c    D'un milieu de face a un sommet (par face)
  33 c ______________________________________________________________________
  34 c .        .     .        .                                            .
  35 c .  nom   . e/s . taille .           description                      .
  36 c .____________________________________________________________________.
  37 c . nbfop2 . e   .    1   . nombre de fonctions P2                     .
  38 c . profho . es  .   *    . pour chaque entite en numerotation homard :.
  39 c .        .     .        . 0 : l'entite est absente du profil         .
  40 c .        .     .        . 1 : l'entite est presente dans le profil   .
  41 c . vap2ho . es  . nbfop2*. variables p2 numerotation homard           .
  42 c .        .     . nbnoto .                                            .
  43 c . somare . e   .2*nbarto. numeros des extremites d'arete             .
  44 c . np2are . e   . nbarto . numero des noeuds p2 milieux d'aretes      .
  45 c . hetqua . e   . nbquto . historique de l'etat des quadrangles       .
  46 c . quahex . e   .nbhecf*6. numeros des 6 quadrangles des hexaedres    .
  47 c . lehexa . e   .    1   . hexaedre a traiter                         .
  48 c . listso . e   .    8   . liste des sommets de l'hexaedre            .
  49 c . listno . e   .   12   . liste des noeuds de l'hexaedre             .
  50 c . areint . e   .   *    . numeros globaux des aretes internes        .
  51 c . ulsort . e   .   1    . numero d'unite logique de la liste standard.
  52 c . langue . e   .    1   . langue des messages                        .
  53 c .        .     .        . 1 : francais, 2 : anglais                  .
  54 c . codret . es  .    1   . code de retour des modules                 .
  55 c .        .     .        . 0 : pas de probleme                        .
  56 c .        .     .        . 1 : probleme                               .
  57 c ______________________________________________________________________
  58 c
  59 c====
  60 c 0. declarations et dimensionnement
  61 c====
  62 c
  63 c 0.1. ==> generalites
  64 c
  65       implicit none
  66       save
  67 c
  68 #include "fractf.h"
  69 #include "fractg.h"
  70 #include "fracth.h"
  71 c
  72       character*6 nompro
  73       parameter ( nompro = 'PCS2H5' )
  74 c
  75 #include "nblang.h"
  76 c
  77 c 0.2. ==> communs
  78 c
  79 #include "envex1.h"
  80 c
  81 #include "nombar.h"
  82 #include "nombqu.h"
  83 #include "nombhe.h"
  84 c
  85 c 0.3. ==> arguments
  86 c
  87       integer nbfop2
  88       integer profho(*)
  89       integer somare(2,nbarto), np2are(nbarto)
  90       integer hetqua(nbquto)
  91       integer quahex(nbhecf,6)
  92       integer lehexa
  93       integer listso(8), listno(12)
  94       integer areint(4)
  95 c
  96       double precision vap2ho(nbfop2,*)
  97 c
  98       integer ulsort, langue, codret
  99 c
 100 c 0.4. ==> variables locales
 101 c
 102       integer iaux, jaux
 103       integer nuface, etahex
 104       integer nusomm, laface, larete
 105       integer listns(20)
 106       integer sm, nuv
 107       integer iaux1, iaux2, iaux3, iaux4
 108 c
 109       integer nbmess
 110       parameter ( nbmess = 100 )
 111       character*80 texte(nblang,nbmess)
 112 c
 113 #include "impr03.h"
 114 c ______________________________________________________________________
 115 c
 116 #include "impr01.h"
 117 c
 118 cgn        write (ulsort,texte(langue,1)) 'Entree', nompro
 119 #ifdef _DEBUG_HOMARD_
 120       write (ulsort,90002) 'listso', listso
 121       write (ulsort,90002) 'listno  1-8', (listno(iaux),iaux=1,8)
 122       write (ulsort,90002) 'listno 9-12', (listno(iaux),iaux=9,12)
 123 #endif
 124 c
 125 c====
 126 c 0. Reperage de la face coupees
 127 c====
 128 c
 129       do 11 , nuface = 1 , 6
 130 c
 131         laface = quahex(lehexa,nuface)
 132         if ( mod(hetqua(laface),100).eq.4 ) then
 133           etahex = 40 + nuface
 134         endif
 135 c
 136    11 continue
 137 c
 138       if ( etahex.lt.41 .or. etahex.gt.46 ) then
 139         write(ulsort,*) 'Pb. Dans pcs2h5, etahex =',etahex
 140         codret = 11
 141       endif
 142 c
 143 c    On passe en revue tous les sommets
 144 c    Ils sont parcourus dans l'ordre des aretes a1, a2, a3 et a4 de
 145 c    la pyramide de base.
 146 c
 147       do 10 , nusomm = 1 , 4
 148 cgn      write(6,*) 'Dans pcs2h5, nusomm =',nusomm
 149 c
 150 c====
 151 c 1. Les 2 sommets les plus proches et les 2 les plus eloignes
 152 c    Remarques :
 153 c    . L'un des sommets les plus proches est celui ou aboutit l'arete
 154 c      du noeud milieu a interpoler (iaux1)
 155 c    . L'autre sommet le plus proche est la 2nde extremite de l'arete
 156 c      de l'hexaedre qui relie ce sommet a la face coupee (iaux2)
 157 c    . Les sommets eloignes sont deduits par la regle :
 158 c          somme des numeros locaux de deux sommets opposes = 9
 159 c====
 160 c
 161       if ( etahex.eq.41 ) then
 162         if ( nusomm.eq.1 ) then
 163           iaux1 = 6
 164           iaux2 = 1
 165         elseif ( nusomm.eq.2 ) then
 166           iaux1 = 5
 167           iaux2 = 2
 168         elseif ( nusomm.eq.3 ) then
 169           iaux1 = 8
 170           iaux2 = 3
 171         else
 172           iaux1 = 7
 173           iaux2 = 4
 174         endif
 175       elseif ( etahex.eq.42 ) then
 176         if ( nusomm.eq.1 ) then
 177           iaux1 = 3
 178           iaux2 = 2
 179         elseif ( nusomm.eq.2 ) then
 180           iaux1 = 4
 181           iaux2 = 1
 182         elseif ( nusomm.eq.3 ) then
 183           iaux1 = 7
 184           iaux2 = 6
 185         else
 186           iaux1 = 8
 187           iaux2 = 5
 188         endif
 189       elseif ( etahex.eq.43 ) then
 190         if ( nusomm.eq.1 ) then
 191           iaux1 = 2
 192           iaux2 = 1
 193         elseif ( nusomm.eq.2 ) then
 194           iaux1 = 3
 195           iaux2 = 4
 196         elseif ( nusomm.eq.3 ) then
 197           iaux1 = 8
 198           iaux2 = 7
 199         else
 200           iaux1 = 5
 201           iaux2 = 6
 202         endif
 203       elseif ( etahex.eq.44 ) then
 204         if ( nusomm.eq.1 ) then
 205           iaux1 = 1
 206           iaux2 = 2
 207         elseif ( nusomm.eq.2 ) then
 208           iaux1 = 6
 209           iaux2 = 5
 210         elseif ( nusomm.eq.3 ) then
 211           iaux1 = 7
 212           iaux2 = 8
 213         else
 214           iaux1 = 4
 215           iaux2 = 3
 216         endif
 217       elseif ( etahex.eq.45 ) then
 218         if ( nusomm.eq.1 ) then
 219           iaux1 = 1
 220           iaux2 = 4
 221         elseif ( nusomm.eq.2 ) then
 222           iaux1 = 2
 223           iaux2 = 3
 224         elseif ( nusomm.eq.3 ) then
 225           iaux1 = 5
 226           iaux2 = 8
 227         else
 228           iaux1 = 6
 229           iaux2 = 7
 230         endif
 231       else
 232         if ( nusomm.eq.1 ) then
 233           iaux1 = 2
 234           iaux2 = 5
 235         elseif ( nusomm.eq.2 ) then
 236           iaux1 = 1
 237           iaux2 = 6
 238         elseif ( nusomm.eq.3 ) then
 239           iaux1 = 4
 240           iaux2 = 7
 241         else
 242           iaux1 = 3
 243           iaux2 = 8
 244         endif
 245       endif
 246 c
 247       listns(1) = listso(iaux1)
 248       listns(2) = listso(iaux2)
 249       listns(7) = listso(9-iaux1)
 250       listns(8) = listso(9-iaux2)
 251 c
 252 c====
 253 c 2. Les sommets intermediaires
 254 c    Il suffit d'en noter 2 sur la meme arete de l'hexaedre, les deux
 255 c    autres sont deduits par la regle :
 256 c          somme des numeros locaux de deux sommets opposes = 9
 257 c====
 258 c
 259       if ( etahex.eq.41 ) then
 260         if ( nusomm.eq.1 .or. nusomm.eq.3 ) then
 261           iaux1 = 2
 262           iaux2 = 5
 263         else
 264           iaux1 = 1
 265           iaux2 = 6
 266         endif
 267       elseif ( etahex.eq.42 ) then
 268         if ( nusomm.eq.1 .or. nusomm.eq.3 ) then
 269           iaux1 = 1
 270           iaux2 = 4
 271         else
 272           iaux1 = 2
 273           iaux2 = 3
 274         endif
 275       elseif ( etahex.eq.43 .or. etahex.eq.44 ) then
 276         if ( nusomm.eq.1 .or. nusomm.eq.3 ) then
 277           iaux1 = 3
 278           iaux2 = 4
 279         else
 280           iaux1 = 1
 281           iaux2 = 2
 282         endif
 283       elseif ( etahex.eq.45 ) then
 284         if ( nusomm.eq.1 .or. nusomm.eq.3 ) then
 285           iaux1 = 2
 286           iaux2 = 3
 287         else
 288           iaux1 = 1
 289           iaux2 = 4
 290         endif
 291       else
 292         if ( nusomm.eq.1 .or. nusomm.eq.3 ) then
 293           iaux1 = 1
 294           iaux2 = 6
 295         else
 296           iaux1 = 2
 297           iaux2 = 5
 298         endif
 299       endif
 300 c
 301       listns(3) = listso(iaux1)
 302       listns(4) = listso(iaux2)
 303       listns(5) = listso(9-iaux1)
 304       listns(6) = listso(9-iaux2)
 305 c
 306 c====
 307 c 3. Le noeud le plus proche, le plus eloigne et les 2 coplanaires
 308 c    dans un plan parallelle a la face coupee
 309 c    Remarques :
 310 c    . Le noeud le plus proche est celui au milieu de l'arete
 311 c      de l'hexaedre qui relie ce sommet a la face coupee (iaux1)
 312 c    . Le noeud eloigne est deduit par la regle :
 313 c          somme des numeros locaux de deux noeuds opposes = 13
 314 c    . Le dernier noeud est deduit par la regle :
 315 c          somme des numeros locaux de deux noeuds opposes = 13
 316 c====
 317 c
 318       if ( etahex.eq.41 ) then
 319         if ( nusomm.eq.1 ) then
 320           iaux1 = 5
 321           iaux2 = 6
 322         elseif ( nusomm.eq.2 ) then
 323           iaux1 = 6
 324           iaux2 = 5
 325         elseif ( nusomm.eq.3 ) then
 326           iaux1 = 8
 327           iaux2 = 6
 328         else
 329           iaux1 = 7
 330           iaux2 = 5
 331         endif
 332       elseif ( etahex.eq.42 ) then
 333         if ( nusomm.eq.1 ) then
 334           iaux1 = 3
 335           iaux2 = 2
 336         elseif ( nusomm.eq.2 ) then
 337           iaux1 = 2
 338           iaux2 = 3
 339         elseif ( nusomm.eq.3 ) then
 340           iaux1 = 10
 341           iaux2 = 2
 342         else
 343           iaux1 = 11
 344           iaux2 = 3
 345         endif
 346       elseif ( etahex.eq.43 ) then
 347         if ( nusomm.eq.1 ) then
 348           iaux1 = 1
 349           iaux2 = 4
 350         elseif ( nusomm.eq.2 ) then
 351           iaux1 = 4
 352           iaux2 = 1
 353         elseif ( nusomm.eq.3 ) then
 354           iaux1 = 12
 355           iaux2 = 4
 356         else
 357           iaux1 = 9
 358           iaux2 = 1
 359         endif
 360       elseif ( etahex.eq.44 ) then
 361         if ( nusomm.eq.1 ) then
 362           iaux1 = 1
 363           iaux2 = 4
 364         elseif ( nusomm.eq.2 ) then
 365           iaux1 = 9
 366           iaux2 = 1
 367         elseif ( nusomm.eq.3 ) then
 368           iaux1 = 12
 369           iaux2 = 4
 370         else
 371           iaux1 = 4
 372           iaux2 = 1
 373         endif
 374       elseif ( etahex.eq.45 ) then
 375         if ( nusomm.eq.1 ) then
 376           iaux1 = 2
 377           iaux2 = 3
 378         elseif ( nusomm.eq.2 ) then
 379           iaux1 = 3
 380           iaux2 = 2
 381         elseif ( nusomm.eq.3 ) then
 382           iaux1 = 11
 383           iaux2 = 3
 384         else
 385           iaux1 = 10
 386           iaux2 = 2
 387         endif
 388       else
 389         if ( nusomm.eq.1 ) then
 390           iaux1 = 6
 391           iaux2 = 5
 392         elseif ( nusomm.eq.2 ) then
 393           iaux1 = 5
 394           iaux2 = 6
 395         elseif ( nusomm.eq.3 ) then
 396           iaux1 = 7
 397           iaux2 = 5
 398         else
 399           iaux1 = 8
 400           iaux2 = 6
 401         endif
 402       endif
 403 c
 404       listns( 9) = listno(iaux1)
 405       listns(20) = listno(13-iaux1)
 406       listns(18) = listno(iaux2)
 407       listns(19) = listno(13-iaux2)
 408 c
 409 c====
 410 c 4. Les 4 noeuds intermediaires les plus proches
 411 c    et les 4 noeuds intermediaires les plus eloignes
 412 c====
 413 c
 414       if ( etahex.eq.41 ) then
 415         if ( nusomm.eq.1 ) then
 416           iaux1 = 1
 417           iaux2 = 2
 418           iaux3 = 9
 419           iaux4 = 10
 420         elseif ( nusomm.eq.2 ) then
 421           iaux1 = 1
 422           iaux2 = 3
 423           iaux3 = 9
 424           iaux4 = 11
 425         elseif ( nusomm.eq.3 ) then
 426           iaux1 = 3
 427           iaux2 = 4
 428           iaux3 = 11
 429           iaux4 = 12
 430         else
 431           iaux1 = 2
 432           iaux2 = 4
 433           iaux3 = 10
 434           iaux4 = 12
 435         endif
 436       elseif ( etahex.eq.42 ) then
 437         if ( nusomm.eq.1 ) then
 438           iaux1 = 1
 439           iaux2 = 4
 440           iaux3 = 6
 441           iaux4 = 8
 442         elseif ( nusomm.eq.2 ) then
 443           iaux1 = 1
 444           iaux2 = 4
 445           iaux3 = 5
 446           iaux4 = 7
 447         elseif ( nusomm.eq.3 ) then
 448           iaux1 = 5
 449           iaux2 = 7
 450           iaux3 = 9
 451           iaux4 = 12
 452         else
 453           iaux1 = 6
 454           iaux2 = 8
 455           iaux3 = 9
 456           iaux4 = 12
 457         endif
 458       elseif ( etahex.eq.43 ) then
 459         if ( nusomm.eq.1 ) then
 460           iaux1 = 2
 461           iaux2 = 3
 462           iaux3 = 5
 463           iaux4 = 6
 464         elseif ( nusomm.eq.2 ) then
 465           iaux1 = 2
 466           iaux2 = 3
 467           iaux3 = 7
 468           iaux4 = 8
 469         elseif ( nusomm.eq.3 ) then
 470           iaux1 = 7
 471           iaux2 = 8
 472           iaux3 = 10
 473           iaux4 = 11
 474         else
 475           iaux1 = 5
 476           iaux2 = 6
 477           iaux3 = 10
 478           iaux4 = 11
 479         endif
 480       elseif ( etahex.eq.44 ) then
 481         if ( nusomm.eq.1 ) then
 482           iaux1 = 2
 483           iaux2 = 3
 484           iaux3 = 5
 485           iaux4 = 6
 486         elseif ( nusomm.eq.2 ) then
 487           iaux1 = 5
 488           iaux2 = 6
 489           iaux3 = 10
 490           iaux4 = 11
 491         elseif ( nusomm.eq.3 ) then
 492           iaux1 = 7
 493           iaux2 = 8
 494           iaux3 = 10
 495           iaux4 = 11
 496         else
 497           iaux1 = 2
 498           iaux2 = 3
 499           iaux3 = 7
 500           iaux4 = 8
 501         endif
 502       elseif ( etahex.eq.45 ) then
 503         if ( nusomm.eq.1 ) then
 504           iaux1 = 1
 505           iaux2 = 4
 506           iaux3 = 5
 507           iaux4 = 7
 508         elseif ( nusomm.eq.2 ) then
 509           iaux1 = 1
 510           iaux2 = 6
 511           iaux3 = 4
 512           iaux4 = 8
 513         elseif ( nusomm.eq.3 ) then
 514           iaux1 = 6
 515           iaux2 = 8
 516           iaux3 = 9
 517           iaux4 = 12
 518         else
 519           iaux1 = 5
 520           iaux2 = 7
 521           iaux3 = 9
 522           iaux4 = 12
 523         endif
 524       else
 525         if ( nusomm.eq.1 ) then
 526           iaux1 = 1
 527           iaux2 = 3
 528           iaux3 = 9
 529           iaux4 = 11
 530         elseif ( nusomm.eq.2 ) then
 531           iaux1 = 1
 532           iaux2 = 2
 533           iaux3 = 9
 534           iaux4 = 10
 535         elseif ( nusomm.eq.3 ) then
 536           iaux1 = 2
 537           iaux2 = 4
 538           iaux3 = 10
 539           iaux4 = 12
 540         else
 541           iaux1 = 3
 542           iaux2 = 4
 543           iaux3 = 11
 544           iaux4 = 12
 545         endif
 546       endif
 547 c
 548       listns(10) = listno(iaux1)
 549       listns(11) = listno(iaux2)
 550       listns(12) = listno(iaux3)
 551       listns(13) = listno(iaux4)
 552 c
 553       listns(14) = listno(13-iaux1)
 554       listns(15) = listno(13-iaux2)
 555       listns(16) = listno(13-iaux3)
 556       listns(17) = listno(13-iaux4)
 557 c
 558 c====
 559 c 5. L'arete concernee : celle des aretes internes qui demarrent
 560 c    ou finissent sur le sommet en cours
 561 c====
 562 c
 563 cgn      write (ulsort,90002) 'listns(1)', listns(1)
 564       do 62 , iaux = 1 , 4
 565         larete = areint(iaux)
 566         if ( ( somare(1,larete).eq.listns(1) ) .or.
 567      >       ( somare(2,larete).eq.listns(1) ) ) then
 568           sm = np2are(larete)
 569           goto 620
 570         endif
 571   62  continue
 572       write(ulsort,*) nompro//' - aucune arete interne ne correspond ?'
 573       codret = 62
 574 c
 575   620 continue
 576 c
 577 c====
 578 c 7. Interpolation
 579 c====
 580 c
 581       if ( codret.eq.0 ) then
 582 c
 583 cgn      write (ulsort,90002) 'sm', sm
 584       profho(sm) = 1
 585 c
 586 cgn        write (ulsort,90002) 'listns  1- 8',(listns(jaux),jaux=1,8)
 587 cgn        write (ulsort,90002) 'listns  9-16',(listns(jaux),jaux=9,16)
 588 cgn        write (ulsort,90002) 'listns 17-20',(listns(jaux),jaux=17,20)
 589 c
 590       do 71, nuv = 1 , nbfop2
 591 cgn          do 711 , jaux =1 ,20
 592 cgn        write (ulsort,90014) listns(jaux), vap2ho(nuv,listns(jaux))
 593 cgn  711 continue
 594 c
 595           vap2ho(nuv,sm) = - nfstr2 * ( vap2ho(nuv,listns(1))
 596      >                                + vap2ho(nuv,listns(2)) )
 597      >                     - trssz  * ( vap2ho(nuv,listns(3))
 598      >                                + vap2ho(nuv,listns(4))
 599      >                                + vap2ho(nuv,listns(5))
 600      >                                + vap2ho(nuv,listns(6)) )
 601      >                     - trstr2 * ( vap2ho(nuv,listns(7))
 602      >                                + vap2ho(nuv,listns(8)) )
 603      >                     + nessz  *   vap2ho(nuv,listns(9))
 604      >                     + nfstr2 * ( vap2ho(nuv,listns(10))
 605      >                                + vap2ho(nuv,listns(11))
 606      >                                + vap2ho(nuv,listns(12))
 607      >                                + vap2ho(nuv,listns(13)) )
 608      >                     + trstr2 * ( vap2ho(nuv,listns(14))
 609      >                                + vap2ho(nuv,listns(15))
 610      >                                + vap2ho(nuv,listns(16))
 611      >                                + vap2ho(nuv,listns(17)) )
 612      >                     + trssz  * ( vap2ho(nuv,listns(18))
 613      >                                + vap2ho(nuv,listns(19)) )
 614      >                     + unssz  *   vap2ho(nuv,listns(20))
 615 c
 616 cgn        write (ulsort,*) 'vap2ho(nuv,',sm,') =',vap2ho(nuv,sm)
 617    71 continue
 618 c
 619       endif
 620 c
 621    10 continue
 622 c
 623 c====
 624 c 8. La fin
 625 c====
 626 c
 627       if ( codret.ne.0 ) then
 628 c
 629 #include "envex2.h"
 630 c
 631       write (ulsort,texte(langue,1)) 'Sortie', nompro
 632       write (ulsort,texte(langue,2)) codret
 633 c
 634       endif
 635 c
 636 #ifdef _DEBUG_HOMARD_
 637       write (ulsort,texte(langue,1)) 'Sortie', nompro
 638       call dmflsh (iaux)
 639 #endif
 640 c
 641       end