src/HEXABLOCK/HexElements_piq.cxx

   1
   2 // C++ : Table d'hexaedres (Evol Versions 3)
   3
   4 // Copyright (C) 2009-2013  CEA/DEN, EDF R&D
   5 //
   6 // This library is free software; you can redistribute it and/or
   7 // modify it under the terms of the GNU Lesser General Public
   8 // License as published by the Free Software Foundation; either
   9 // version 2.1 of the License.
  10 //
  11 // This library is distributed in the hope that it will be useful,
  12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
  13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14 // Lesser General Public License for more details.
  15 //
  16 // You should have received a copy of the GNU Lesser General Public
  17 // License along with this library; if not, write to the Free Software
  18 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  19 //
  20 // See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  21 //
  22
  23 #include "HexElements.hxx"
  24
  25 #include "HexVector.hxx"
  26 #include "HexVertex.hxx"
  27 #include "HexEdge.hxx"
  28 #include "HexDiagnostics.hxx"
  29 #include "HexDocument.hxx"
  30 #include "HexHexa.hxx"
  31
  32 #include <cmath>
  33 #include <map>
  34
  35 BEGIN_NAMESPACE_HEXA
  36
  37 static bool db = false;
  38
  39 // --------------------------------------------------------
  40 struct PatQuad
  41 {
  42    int   q_edge [QUAD4];
  43    Quad* refer;
  44 };
  45 // --------------------------------------------------------
  46 struct PatVertex
  47 {
  48    double  v_x, v_y;
  49    Vertex* refer;
  50 };
  51 // --------------------------------------------------------
  52 struct PatEdge
  53 {
  54    int   v_amont, v_aval;
  55    Edge* refer;
  56    int   nbr_refer;
  57 };
  58
  59 // --------------------------------------------------------
  60 class Pattern
  61 {
  62 friend class Elements;
  63 public  :
  64     int initialize (Vertex* v1, Vertex* v2, Vertex* v3);
  65     int addQuad    (Quad* quad);
  66     int verify     (int &nbed, int &nbver);
  67
  68 private :
  69     int addEdge   (Edge*   edge);
  70     int addVertex (Vertex* vertex);
  71
  72 private :
  73    enum EnumProj { ProjXY, ProjYZ, ProjZX };
  74
  75    vector <PatVertex> pat_vertex;
  76    vector <PatEdge  > pat_edge;
  77    vector <PatQuad  > pat_quad;
  78
  79    int      nbr_vertex, nbr_edges, nbr_quads;
  80    double   determinant;
  81    EnumProj projection;
  82    Real3    base_i, base_j, origine;
  83
  84    int pos_vertex4;
  85    int pos_edge3, pos_edge4;
  86 };
  87
  88 // ====================================================== initialize
  89 int Pattern::initialize (Vertex* v1, Vertex* v2, Vertex* v3)
  90 {
  91    nbr_vertex  = nbr_edges = nbr_quads = 0;
  92    projection  = ProjXY;
  93    determinant = 1;
  94    base_i[0]   = base_i[1] =  base_i[2] = 0;
  95    base_j[0]   = base_j[1] =  base_j[2] = 0;
  96
  97    if (v1==NULL || v2==NULL || v3==NULL)
  98       return HERR;
  99
 100    addVertex (v1);
 101    addVertex (v2);
 102    addVertex (v3);
 103    if (nbr_vertex!=3)
 104       return HERR;
 105
 106    Document* doc = v1->dad ();
 107    Edge* edc1  = doc->findEdge (v1, v2);
 108    Edge* edc2  = doc->findEdge (v2, v3);
 109    if (edc1==NULL || edc2==NULL)
 110       return HERR;
 111
 112    addEdge (edc1);
 113    addEdge (edc2);
 114    if (nbr_edges!=2)
 115       return HERR;
 116
 117    pat_edge[0].nbr_refer = pat_edge[1].nbr_refer = 0;
 118
 119    pat_vertex[0].v_x = pat_vertex[1].v_x = 0;
 120    pat_vertex[2].v_x = 1;
 121    pat_vertex[0].v_y = 1;
 122    pat_vertex[1].v_y = pat_vertex[2].v_y = 0;
 123
 124    Real3 pb, pc;
 125    v1->getPoint (pc);
 126    v2->getPoint (origine);
 127    v3->getPoint (pb);
 128
 129    calc_vecteur (origine, pb, base_i);
 130    calc_vecteur (origine, pc, base_j);
 131
 132    if (db)
 133       {
 134       PutCoord (origine);
 135       PutCoord (base_i);
 136       PutCoord (base_j);
 137       PutData (determinant);
 138       }
 139
 140 /* ******************************
 141  * AB = pu.vI + pv.vJ
 142  *
 143  * vx = pu.vI[0] + pv.vJ[0]   (1) *+vJ[1] ) };
 144  *
 145  * (1 & 2 ) : pu =  (vx.vJ[1] - vy.vJ[0]) / detxy
 146  *            pv = -(vx.vI[1] - vy.vI[0]) / detxy
 147  *x
 148  * (2 & 3 ) : pu =  (vy.vJ[2] - vz.vJ[1]) / detyz
 149  *            pv = -(vy.vI[2] - vz.vI[1]) / detyz
 150  *
 151  * (3 & 1 ) : pu =  (vz.vJ[0] - vx.vJ[2]) / detzx
 152  *            pv = -(vz.vI[0] - vx.vI[2]) / detzx
 153  *
 154  *  Les 3 systemes d'equations sont valides.
 155  *  On va choisir celui dont la valeur absolue du determinant est maximale
 156    ****************************** */
 157
 158    double detxy = base_i[dir_x]*base_j[dir_y] - base_i[dir_y]*base_j[dir_x];
 159    double detyz = base_i[dir_y]*base_j[dir_z] - base_i[dir_z]*base_j[dir_y];
 160    double detzx = base_i[dir_z]*base_j[dir_x] - base_i[dir_x]*base_j[dir_z];
 161
 162    determinant = detxy;
 163    projection  = ProjXY;
 164
 165    if (fabs (detyz) > fabs (determinant))
 166       {
 167       determinant = detyz;
 168       projection  = ProjYZ;
 169       }
 170    if (fabs (detzx) > fabs (determinant))
 171       {
 172       determinant = detzx;
 173       projection  = ProjZX;
 174       }
 175
 176    return HOK;
 177 }
 178 // ====================================================== verify
 179 int Pattern::verify (int &nbed, int &nbver)
 180 {
 181    nbed  = nbr_edges;
 182    nbver = nbr_vertex;
 183    pos_edge3 = pos_edge4 = pos_vertex4 = NOTHING;
 184
 185    if (pat_edge[0].nbr_refer!=1 || pat_edge[1].nbr_refer!=1)
 186       return HERR;
 187
 188    for (int nro=2 ; nro<nbr_edges; nro++)
 189        {
 190        if (pat_edge[nro].nbr_refer==1)
 191           {
 192           int pv1 = pat_edge[nro].v_amont;
 193           int pv2 = pat_edge[nro].v_aval;
 194
 195           if (pv1==2)
 196              {
 197              pos_edge3   = nro;
 198              if (pos_vertex4 == NOTHING)
 199                  pos_vertex4 =  pv2;
 200              else if (pos_vertex4 != pv2)
 201                  {
 202                  return HERR;
 203                  }
 204              }
 205           else if (pv2==2)
 206              {
 207              pos_edge3   = nro;
 208              if (pos_vertex4 == NOTHING)
 209                  pos_vertex4 =  pv1;
 210              else if (pos_vertex4 != pv1)
 211                  {
 212                  return HERR;
 213                  }
 214              }
 215           else if (pv1==0)
 216              {
 217              pos_edge4   = nro;
 218              if (pos_vertex4 == NOTHING)
 219                  pos_vertex4 =  pv2;
 220              else if (pos_vertex4 != pv2)
 221                  {
 222                  return HERR;
 223                  }
 224              }
 225           else if (pv2==0)
 226              {
 227              pos_edge4   = nro;
 228              if (pos_vertex4 == NOTHING)
 229                  pos_vertex4 =  pv1;
 230              else if (pos_vertex4 != pv1)
 231                  {
 232                  return HERR;
 233                  }
 234              }
 235           }
 236        }
 237
 238    if (pos_edge3==NOTHING || pos_edge4==NOTHING)
 239       return HERR;
 240
 241    return HOK;
 242 }
 243 // ====================================================== addQuad
 244 int Pattern::addQuad (Quad* elt)
 245 {
 246    if (elt==NULL)
 247       return HERR;
 248
 249    PatQuad quad;
 250    quad.refer = elt;
 251    for (int nro=0; nro<QUAD4 ; nro++)
 252        {
 253        Edge* edge = elt->getEdge (nro);
 254        quad.q_edge [nro] = addEdge (edge);
 255        }
 256
 257    pat_quad.push_back (quad);
 258    nbr_quads++;
 259    return HOK;
 260 }
 261 // ====================================================== addEdge
 262 int Pattern::addEdge (Edge* elt)
 263 {
 264    for (int nro=0; nro<nbr_edges ; nro++)
 265        {
 266        if (elt==pat_edge [nro].refer)
 267           {
 268           pat_edge [nro].nbr_refer++;
 269           return nro;
 270           }
 271        }
 272
 273    PatEdge edge;
 274    edge.nbr_refer = 1;
 275    edge.refer    = elt;
 276    edge.v_amont  = addVertex (elt->getVertex (V_AMONT));
 277    edge.v_aval   = addVertex (elt->getVertex (V_AVAL));
 278
 279    pat_edge.push_back (edge);
 280    nbr_edges++;
 281    return nbr_edges-1;
 282 }
 283 // ====================================================== addVertex
 284 int Pattern::addVertex (Vertex* elt)
 285 {
 286    for (int nro=0; nro<nbr_vertex ; nro++)
 287        {
 288        if (elt==pat_vertex [nro].refer)
 289           return nro;
 290        }
 291
 292    PatVertex vertex;
 293    vertex.refer = elt;
 294    double vx = elt->getX() - origine [dir_x];
 295    double vy = elt->getY() - origine [dir_y];
 296    double vz = elt->getZ() - origine [dir_z];
 297    switch (projection)
 298       {
 299       case ProjXY :  default :
 300            vertex.v_x =  (vx*base_j[dir_y] - vy*base_j[dir_x]) / determinant;
 301            vertex.v_y = -(vx*base_i[dir_y] - vy*base_i[dir_x]) / determinant;
 302            break;
 303
 304       case ProjYZ :
 305            vertex.v_y =  (vy*base_j[dir_z] - vz*base_j[dir_y]) / determinant;
 306            vertex.v_x = -(vy*base_i[dir_z] - vz*base_i[dir_y]) / determinant;
 307            break;
 308
 309       case ProjZX :
 310            vertex.v_x =  (vz*base_j[dir_x] - vx*base_j[dir_z]) / determinant;
 311            vertex.v_y = -(vz*base_i[dir_x] - vx*base_i[dir_z]) / determinant;
 312            break;
 313       }
 314
 315    if (db)
 316       printf (" Vertex nro %d : (%g, %g, %g)\t -> (%g, %g)\n",
 317                         nbr_vertex, vx, vy, vz, vertex.v_x, vertex.v_y);
 318    pat_vertex.push_back (vertex);
 319    nbr_vertex++;
 320    return nbr_vertex-1;
 321 }
 322 // -------------------------------------------------------------------
 323 // -------------------------------------------------------------------
 324 // -------------------------------------------------------------------
 325 // -------------------------------------------------------------------
 326 // ====================================================== replaceHexas
 327 int Elements::replaceHexas (Quads& liste, Vertex* p1, Vertex* c1,
 328                             Vertex* p2, Vertex* c2, Vertex* p3, Vertex* c3)
 329 {
 330     Edge* edp1  = el_root->findEdge (p1, p2);
 331     Edge* edp2  = el_root->findEdge (p2, p3);
 332     Edge* edc1  = el_root->findEdge (c1, c2);
 333     Edge* edc2  = el_root->findEdge (c2, c3);
 334
 335     Quad* quadc = el_root->findQuad (c1, c3);
 336
 337     if (edp1==NULL || edp2==NULL || edc1==NULL || edc2==NULL || quadc==NULL)
 338        {
 339        printf ("... Error in HexaBlock function \n");
 340        printf ("... doc.replace (lquads, p1,c1, p2,c2, p3,c3)\n");
 341        if (edp1==NULL)
 342            printf ("Vertices p1 and p2 don't define an edge\n");
 343        else if (edp2==NULL)
 344            printf ("Vertices p2 and p3 don't define an edge\n");
 345        else if (edc1==NULL)
 346            printf ("Vertices c1 and c2 don't define an edge\n");
 347        else if (edc2==NULL)
 348            printf ("Vertices c2 and c3 don't define an edge\n");
 349        else if (quadc==NULL)
 350            printf ("Vertices c1 and c3 don't define a quad\n");
 351        return HERR;
 352        }
 353
 354     int np = quadc->getNbrParents ();
 355     Hexa* hexac = quadc->getParent (0);
 356
 357     if (np!=1 || hexac==NULL)
 358        {
 359        printf ("... Error in HexaBlock function \n");
 360        printf ("... doc.replace (lquads, p1,c1, p2,c2, p3,c3)\n");
 361        printf ("Quad (c1,c2,c3) is not an external quad\n");
 362        return HERR;
 363        }
 364                               // Analyse du pattern
 365     int nbquads = liste.size();
 366     Pattern   pattern;
 367     int ier = pattern.initialize (p1, p2, p3);
 368     if (ier!=HOK)
 369        {
 370        printf ("... Error in HexaBlock function \n");
 371        printf ("... doc.replace (lquads, p1,c1, p2,c2, p3,c3)\n");
 372        printf ("Vertices (p1,p2,p3) don't define a virtual quad\n");
 373        return HERR;
 374        }
 375
 376     for (int nq=0 ; nq<nbquads ; nq++)
 377         {
 378         ier = pattern.addQuad (liste[nq]);
 379         if (ier!=HOK)
 380            {
 381            printf ("... Error in HexaBlock function \n");
 382            printf ("... doc.replace (lquads, p1,c1, p2,c2, p3,c3)\n");
 383            printf ("Quad nr %d of the list is NULL\n", nq+1);
 384            return HERR;
 385            }
 386         }
 387     int nbv  = 0;
 388     int nbed = 0;
 389     int nbh  = 0;
 390     ier = pattern.verify (nbed, nbv);
 391     if (ier!=HOK)
 392        {
 393        printf ("... Error in HexaBlock function \n");
 394        printf ("... doc.replace (lquads, p1,c1, p2,c2, p3,c3)\n");
 395        printf ("Vertices (p1,p2,p3) don't define a virtual quad\n");
 396        return HERR;
 397        }
 398                               // Analyse de la cible
 399     Quads pil_quad;
 400     Hexas pil_hexa;
 401
 402     pil_hexa.push_back (hexac);
 403     pil_quad.push_back (quadc);
 404     bool more = true;
 405                               // Constitution de la pile des hexaedres
 406     while (more)
 407           {
 408           nbh ++;
 409           Quad* oppos = hexac->getOpposedQuad (quadc);
 410           pil_quad.push_back (oppos);
 411           more = oppos->getNbrParents() == 2;
 412
 413           if (more)
 414              {
 415              if (oppos->getParent (0)==hexac)
 416                  hexac = oppos->getParent(1);
 417              else
 418                  hexac = oppos->getParent(0);
 419              pil_hexa.push_back (hexac);
 420              }
 421           quadc = oppos;
 422           }
 423
 424     resize (GR_REPLACE, nbquads, nbh, nbv, nbed);
 425
 426                 // 1) Constitution des 4 coins de la cible
 427                 // 2) Quadriller la face du haut
 428                 // 3) Trouver les 4 coins de la face du bas
 429                 // 4) Quadriller la face du bas
 430                 // 6) Decouper l'hexaedre defini par ces 2 faces
 431                 // 7) S'il y a un hexaedre dessous goto 3)
 432
 433     Vertex* tvert[QUAD4] = { c1, c2, c3, pil_quad[0]->getOpposVertex(c2) };
 434     replaceQuad (1, &pattern, pil_quad[0], tvert);
 435     for (int nh=1 ; nh<gr_hauteur ; nh++)
 436         {
 437         for (int nro=0 ; nro<QUAD4 ; nro++)
 438             {
 439             Edge*  edv  = pil_hexa[nh-1]->getPerpendicularEdge (pil_quad[nh-1],
 440                                                                 tvert[nro]);
 441             tvert [nro] = edv ->opposedVertex(tvert[nro]);
 442             }
 443         replaceQuad (nh+1, &pattern, pil_quad[nh], tvert);
 444         replaceHexa (nh,   &pattern, pil_hexa[nh-1]);
 445         }
 446
 447     for (int nh=0 ; nh<=nbh ; nh++)
 448         pil_quad[nh]->remove ();
 449
 450     extrudeQuad (&pattern);
 451     replaceHexa (0, &pattern, NULL);
 452     return HOK;
 453 }
 454 // ====================================================== repVertex
 455 void Elements::repVertex (int nh, int nro, Vertex* elt)
 456 {
 457    int addr = nh * pat_nbvertex + nro;
 458    if (tab_vertex[addr] == NULL)
 459        tab_vertex[addr]  = elt;
 460 }
 461 // ====================================================== repVertex
 462 Vertex* Elements::repVertex (int nh, int nro, double px, double py, double pz)
 463 {
 464    int addr = nh * pat_nbvertex + nro;
 465
 466    if (tab_vertex[addr] == NULL)
 467        tab_vertex[addr] = el_root->addVertex (px, py, pz);
 468
 469    return tab_vertex[addr];
 470 }
 471 // ====================================================== repVertex
 472 Vertex* Elements::repVertex (int nh, int nro)
 473 {
 474    int addr = nh * pat_nbvertex + nro;
 475    return tab_vertex[addr];
 476 }
 477 // ====================================================== repEdgeH
 478 void Elements::repEdgeH (int nh, int nro, Edge* elt)
 479 {
 480    int addr = nh * (pat_nbedges + pat_nbvertex) + nro;
 481    if (tab_edge[addr] == NULL)
 482        tab_edge[addr]  = elt;
 483 }
 484 // ====================================================== repEdgeH
 485 Edge* Elements::repEdgeH (int nh, int nro, Vertex* v1, Vertex* v2)
 486 {
 487    int addr = nh * (pat_nbedges + pat_nbvertex) + nro;
 488    if (tab_edge[addr] == NULL)
 489        tab_edge[addr]  = el_root->addEdge (v1, v2);
 490    return tab_edge[addr];
 491 }
 492 // ====================================================== repEdgeH
 493 Edge* Elements::repEdgeH (int nh, int nro)
 494 {
 495    int addr = nh * (pat_nbedges + pat_nbvertex) + nro;
 496    return tab_edge[addr];
 497 }
 498 // ====================================================== repEdgeV
 499 void Elements::repEdgeV (int nh, int nro, Edge* elt)
 500 {
 501    int addr = nh * (pat_nbedges + pat_nbvertex) + pat_nbedges + nro;
 502    if (tab_edge[addr] == NULL)
 503        tab_edge[addr]  = elt;
 504 }
 505 // ====================================================== repEdgeV
 506 Edge* Elements::repEdgeV (int nh, int nro, Vertex* v1, Vertex* v2)
 507 {
 508    int addr = nh * (pat_nbedges + pat_nbvertex) + pat_nbedges + nro;
 509    if (tab_edge[addr] == NULL)
 510        tab_edge[addr]  = el_root->addEdge (v1, v2);
 511    return tab_edge[addr];
 512 }
 513 // ====================================================== repEdgeV
 514 Edge* Elements::repEdgeV (int nh, int nro)
 515 {
 516    int addr = nh * (pat_nbedges + pat_nbvertex) + pat_nbedges + nro;
 517    return tab_edge[addr];
 518 }
 519 // ====================================================== repQuadH
 520 void Elements::repQuadH (int nh, int nro, Quad* elt)
 521 {
 522    int addr = nh * (nbr_orig + pat_nbedges) + nro;
 523    if (tab_quad[addr] == NULL)
 524        tab_quad[addr]  = elt;
 525 }
 526 // ====================================================== repQuadH
 527 Quad* Elements::repQuadH (int nh, int nro)
 528 {
 529    int addr = nh * (nbr_orig + pat_nbedges) + nro;
 530    return tab_quad [addr];
 531 }
 532 // ====================================================== repQuadH
 533 Quad* Elements::repQuadH (int nh, int nro, Edge* e1, Edge* e2, Edge* e3,
 534                                            Edge* e4)
 535 {
 536    int addr = nh * (nbr_orig + pat_nbedges) + nro;
 537    if (tab_quad[addr] == NULL)
 538        tab_quad[addr]  = el_root->addQuad (e1, e2, e3, e4);
 539    return tab_quad [addr];
 540 }
 541 // ====================================================== repQuadV
 542 void Elements::repQuadV (int nh, int nro, Quad* elt)
 543 {
 544    int addr = nh * (nbr_orig + pat_nbedges) + nbr_orig + nro;
 545    if (tab_quad[addr] == NULL)
 546        tab_quad[addr]  = elt;
 547 }
 548 // ====================================================== repQuadV
 549 Quad* Elements::repQuadV (int nh, int nro)
 550 {
 551    int addr = nh * (nbr_orig + pat_nbedges) + nbr_orig + nro;
 552    return tab_quad [addr];
 553 }
 554 // ====================================================== repQuadV
 555 Quad* Elements::repQuadV (int nh, int nro, Edge* e1, Edge* e2, Edge* e3,
 556                                            Edge* e4)
 557 {
 558    int addr = nh * (nbr_orig + pat_nbedges) + nbr_orig + nro;
 559    if (tab_quad[addr] == NULL)
 560        tab_quad[addr]  = el_root->addQuad (e1, e2, e3, e4);
 561    return tab_quad [addr];
 562 }
 563 // ====================================================== repHexa
 564 Hexa* Elements::repHexa (int nh, int nro, Quad* qa, Quad* qb, Quad* qc,
 565                                           Quad* qd, Quad* qe, Quad* qf)
 566 {
 567    int addr = nh * nbr_orig  + nro;
 568    if (tab_hexa[addr] == NULL)
 569        tab_hexa[addr]  = el_root->addHexa (qa, qb, qc, qd, qe, qf);
 570    return tab_hexa [addr];
 571 }
 572 // ====================================================== replaceQuad
 573 // ==== Creation des quads horizontaux
 574 int Elements::replaceQuad (int nh, Pattern* pat, Quad* quad, Vertex* tvert[])
 575 {
 576    Vertex* pvert[QUAD4] = { tvert[0],  tvert[1],  tvert[2],  tvert[3] };
 577
 578    int vnro [4] = { 0, 1, 2, pat->pos_vertex4};
 579    int ednro[4] = { 0, 1, pat->pos_edge3, pat->pos_edge4 };
 580
 581                               // Enregistrement des vertex & edges existant
 582    for (int nro=0 ; nro<QUAD4 ; nro++)
 583        {
 584        Vertex* vh    = pvert [nro];
 585        Vertex* vh1   = pvert [(nro+1) MODULO QUAD4];
 586        Edge*   edh   = quad->findEdge (vh, vh1);
 587
 588        repVertex (nh,   vnro  [nro], vh);
 589        repEdgeH  (nh,   ednro [nro], edh);
 590        }
 591
 592    Real3 orig, ib, jb, pb, pc;
 593                               // Creation des vertex
 594    tvert[0]->getPoint (pc);
 595    tvert[1]->getPoint (orig);
 596    tvert[2]->getPoint (pb);
 597    calc_vecteur (orig, pb, ib);
 598    calc_vecteur (orig, pc, jb);
 599
 600    for (int nro=0 ; nro<pat->nbr_vertex ; nro++)
 601        {
 602        double lambda = pat->pat_vertex [nro].v_x;
 603        double mu     = pat->pat_vertex [nro].v_y;
 604        double px = orig[dir_x] + lambda*ib[dir_x] + mu*jb[dir_x];
 605        double py = orig[dir_y] + lambda*ib[dir_y] + mu*jb[dir_y];
 606        double pz = orig[dir_z] + lambda*ib[dir_z] + mu*jb[dir_z];
 607        repVertex (nh, nro, px, py, pz);
 608        }
 609                               // Creation des edges horizontaux
 610    for (int nro=0 ; nro<pat->nbr_edges ; nro++)
 611        {
 612        int nv1 = pat->pat_edge [nro].v_amont;
 613        int nv2 = pat->pat_edge [nro].v_aval;
 614        Vertex* v1 = repVertex (nh, nv1);
 615        Vertex* v2 = repVertex (nh, nv2);
 616        repEdgeH (nh, nro, v1, v2);
 617        }
 618                               // Creation des quads horizontaux
 619    for (int nro=0 ; nro<pat->nbr_quads ; nro++)
 620        {
 621        Edge* eda = repEdgeH (nh, pat->pat_quad [nro].q_edge [0]);
 622        Edge* edb = repEdgeH (nh, pat->pat_quad [nro].q_edge [1]);
 623        Edge* edc = repEdgeH (nh, pat->pat_quad [nro].q_edge [2]);
 624        Edge* edd = repEdgeH (nh, pat->pat_quad [nro].q_edge [3]);
 625
 626        repQuadH (nh,  nro, eda, edb, edc, edd);
 627        }
 628    return HOK;
 629 }
 630 // ====================================================== extrudeQuad
 631 // ==== Creation des quads horizontaux
 632 int Elements::extrudeQuad (Pattern* pat)
 633 {
 634                               // Creation des vertex de niveau 0
 635    for (int nro=0 ; nro<pat->nbr_vertex ; nro++)
 636        {
 637        Vertex* v1 = repVertex (1, nro);
 638        Vertex* v2 = repVertex (2, nro);
 639        double px = 2*v1->getX() - v2->getX();
 640        double py = 2*v1->getY() - v2->getY();
 641        double pz = 2*v1->getZ() - v2->getZ();
 642        repVertex (0, nro, px, py, pz);
 643        }
 644                               // Creation des edges horizontaux
 645    for (int nro=0 ; nro<pat->nbr_edges ; nro++)
 646        {
 647        int nv1 = pat->pat_edge [nro].v_amont;
 648        int nv2 = pat->pat_edge [nro].v_aval;
 649        Vertex* v1 = repVertex (0, nv1);
 650        Vertex* v2 = repVertex (0, nv2);
 651        repEdgeH (0, nro, v1, v2);
 652        }
 653                               // Creation des quads horizontaux
 654    for (int nro=0 ; nro<pat->nbr_quads ; nro++)
 655        {
 656        Edge* eda = repEdgeH (0, pat->pat_quad [nro].q_edge [0]);
 657        Edge* edb = repEdgeH (0, pat->pat_quad [nro].q_edge [1]);
 658        Edge* edc = repEdgeH (0, pat->pat_quad [nro].q_edge [2]);
 659        Edge* edd = repEdgeH (0, pat->pat_quad [nro].q_edge [3]);
 660
 661        repQuadH (0,  nro, eda, edb, edc, edd);
 662        }
 663    return HOK;
 664 }
 665 // ====================================================== replaceHexa
 666 int Elements::replaceHexa (int nh, Pattern* pat, Hexa* hexa)
 667 {
 668    int vnro [4] = { 0, 1, 2, pat->pos_vertex4};
 669    int ednro[4] = { 0, 1, pat->pos_edge3, pat->pos_edge4 };
 670
 671                               // Enregistrement des edges & quads existants
 672    if (hexa!=NULL)
 673        {
 674        for (int nro=0 ; nro<QUAD4 ; nro++)
 675            {
 676            Vertex* vh  = repVertex (nh,   vnro [nro]);
 677            Vertex* vh1 = repVertex (nh+1, vnro [nro]);
 678            Edge*   edv = hexa->findEdge (vh, vh1);
 679            repEdgeV (nh, vnro [nro], edv);
 680
 681            Edge*   edh   = repEdgeH (nh,   ednro [nro]);
 682            Edge*   edh1  = repEdgeH (nh+1, ednro [nro]);
 683            Quad*   quadv = hexa->findQuad (edh, edh1);
 684            repQuadV (nh, ednro [nro], quadv);
 685            }
 686        }
 687                               // Creation des edges verticaux
 688    for (int nro=0 ; nro<pat->nbr_vertex ; nro++)
 689        {
 690        Vertex* v1 = repVertex (nh,   nro);
 691        Vertex* v2 = repVertex (nh+1, nro);
 692        repEdgeV (nh, nro, v1, v2);
 693        }
 694                               // Creation des quads verticaux
 695    for (int nro=0 ; nro<pat->nbr_edges ; nro++)
 696        {
 697        int nv1 = pat->pat_edge [nro].v_amont;
 698        int nv2 = pat->pat_edge [nro].v_aval;
 699
 700        Edge* eda = repEdgeH (nh,   nro);
 701        Edge* edb = repEdgeV (nh,   nv1);
 702        Edge* edc = repEdgeH (nh+1, nro);
 703        Edge* edd = repEdgeV (nh,   nv2);
 704        repQuadV (nh,  nro, eda, edb, edc, edd);
 705        }
 706                               // Creation des hexaedres
 707    for (int nro=0 ; nro<pat->nbr_quads ; nro++)
 708        {
 709        Quad* qa = repQuadH (nh,   nro);
 710        Quad* qb = repQuadH (nh+1, nro);
 711        Quad* qc = repQuadV (nh, pat->pat_quad [nro].q_edge [0]);
 712        Quad* qd = repQuadV (nh, pat->pat_quad [nro].q_edge [2]);
 713        Quad* qe = repQuadV (nh, pat->pat_quad [nro].q_edge [1]);
 714        Quad* qf = repQuadV (nh, pat->pat_quad [nro].q_edge [3]);
 715        repHexa (nh, nro, qa, qb, qc, qd, qe, qf);
 716        }
 717
 718    return HOK;
 719 }
 720 END_NAMESPACE_HEXA