src/HEXABLOCK/HexElements_piq.cxx

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