src/HEXABLOCK/HexQuad.cxx

   1
   2 // C++ : Gestion des Quadrangles
   3
   4 // Copyright (C) 2009-2024  CEA, EDF
   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, or (at your option) any later version.
  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 #include "HexQuad.hxx"
  23
  24 #include "HexDocument.hxx"
  25 #include "HexHexa.hxx"
  26 #include "HexElements.hxx"
  27 #include "HexGlobale.hxx"
  28
  29 #include "HexXmlWriter.hxx"
  30 #include "HexNewShape.hxx"
  31 #include "HexFaceShape.hxx"
  32
  33 BEGIN_NAMESPACE_HEXA
  34
  35 bool db = on_debug();
  36
  37 // ======================================================== Constructeur
  38 Quad::Quad (Vertex* va, Vertex* vb, Vertex* vc, Vertex* vd)
  39     : EltBase (va->dad(), EL_QUAD)
  40 {
  41    q_vertex [E_A] = va;
  42    q_vertex [E_B] = vb;
  43    q_vertex [E_C] = vc;
  44    q_vertex [E_D] = vd;
  45    q_clone        = NULL;
  46    q_orientation  = Q_UNDEFINED;
  47
  48    for (int nro=0 ; nro<QUAD4 ; nro++)
  49        {
  50        q_edge [nro] = new Edge (q_vertex[nro],
  51                                 q_vertex[(nro+1) MODULO QUAD4]);
  52
  53        if (BadElement (q_vertex [nro]) || BadElement (q_edge [nro]))
  54           setError ();
  55        else
  56           for (int nv=nro+1 ; nv<QUAD4 ; nv++)
  57               if (q_vertex[nv] == q_vertex[nro])
  58                  setError ();
  59        }
  60
  61    if (el_root != NULL && el_status==HOK)
  62        el_root->addQuad (this);
  63    majReferences ();
  64 }
  65 // ======================================================== Constructeur bis
  66 Quad::Quad (Edge* ea, Edge* eb, Edge* ec, Edge* ed)
  67     : EltBase (ea->dad(), EL_QUAD)
  68 {
  69    q_edge [E_A] = ea;
  70    q_edge [E_B] = eb;
  71    q_edge [E_C] = ec;
  72    q_edge [E_D] = ed;
  73    q_clone       = NULL;
  74    q_orientation = Q_UNDEFINED;
  75
  76    for (int nro=0 ; nro<QUAD4 ; nro++)
  77        {
  78        if (BadElement (q_edge[nro]))
  79           setError (886);
  80        else
  81           {
  82           for (int nv=nro+1 ; nv<QUAD4 ; nv++)
  83               if (q_edge[nv] == q_edge[nro])
  84                   setError (887);
  85           }
  86        }
  87
  88    if (isValid())
  89       {
  90                 // Cond necessaire : ea disjoint de ec (opposes)
  91       int nc = ea->inter (ec);
  92       if (nc>=0)
  93          {
  94          nc = ea->inter (eb);
  95          if (nc<0)
  96             {
  97             q_edge [E_C] = eb;
  98             q_edge [E_B] = ec;
  99             }
 100          else
 101             {
 102             nc = ea->inter (ed);
 103             if (nc<0)
 104                {
 105                q_edge [E_C] = ed;
 106                q_edge [E_D] = ec;
 107                }
 108             else
 109                setError (880);
 110             }
 111          }
 112       }
 113
 114    for (int nro=0 ; nro<QUAD4 ; nro++)
 115        {
 116        int prec = (nro+1) MODULO QUAD4;
 117        Vertex* node = NULL;
 118        int nc  = q_edge[nro] -> inter (q_edge[prec]);
 119        if (nc>=0)
 120           node = q_edge[nro]->getVertex (nc);
 121        else
 122           setError (888);
 123        q_vertex [prec] = node;
 124        }
 125
 126    if (isBad())
 127       {
 128       printf (" +++++++++++++++++++++++++++++++++++++++++++ \n");
 129       printf (" +++ Quadrangle impossible \n");
 130       printf (" +++++++++++++++++++++++++++++++++++++++++++ \n");
 131       dump ();
 132       printf (" +++++++++++++++++++++++++++++++++++++++++++ \n");
 133       // el_root->dump ();
 134       for (int ned=0; ned<QUAD4; ned++)
 135           {
 136           q_edge[ned]->dumpPlus ();
 137           }
 138       HexDump (q_vertex[0]);
 139       HexDump (q_vertex[1]);
 140       HexDump (q_vertex[2]);
 141       HexDump (q_vertex[3]);
 142
 143       printf (" +++++++++++++++++++++++++++++++++++++++++++ \n");
 144       fatal_error ("Quadrangle impossible");
 145       }
 146
 147    if (el_root != NULL && el_status==HOK)
 148        el_root->addQuad (this);
 149    majReferences ();
 150 }
 151 // ======================================================== Constructeur ter
 152 Quad::Quad (Quad* other)
 153     : EltBase (other->dad(), EL_QUAD)
 154 {
 155    for (int nro=0 ; nro<QUAD4 ; nro++)
 156        {
 157        q_edge   [nro] = NULL;
 158        q_vertex [nro] = NULL;
 159        }
 160    q_orientation = Q_UNDEFINED;
 161    q_clone       = NULL;
 162
 163    if (el_root != NULL && el_status==HOK)
 164        el_root->addQuad (this);
 165 }
 166 // ============================================================  getEdge
 167 Edge* Quad::getEdge (int nro)
 168 {
 169    Edge* elt = NULL;
 170    if (nro >=0 && nro < QUAD4 && el_status == HOK && q_edge [nro]->isValid())
 171       elt = q_edge [nro];
 172
 173    DumpStart  ("getEdge", nro);
 174    DumpReturn (elt);
 175    return elt;
 176 }
 177 // ============================================================  getVertex
 178 Vertex* Quad::getVertex (int nro)
 179 {
 180    Vertex* elt = NULL;
 181    if (nro >=0 && nro < QUAD4 && el_status == HOK && q_vertex [nro]->isValid())
 182       elt = q_vertex [nro];
 183
 184    DumpStart  ("getVertex", nro);
 185    DumpReturn (elt);
 186    return elt;
 187 }
 188 // ========================================================= majReferences
 189 void Quad::majReferences ()
 190 {
 191    for (int nro=0 ; nro<QUAD4 ; nro++)
 192        q_edge [nro] -> addParent (this);
 193 }
 194 // ========================================================= getParent
 195 Hexa* Quad::getParent  (int nro)
 196 {
 197    return static_cast <Hexa*> (getFather (nro));
 198 }
 199 // ======================================================== anaMerge
 200 int Quad::anaMerge (Vertex* v1, Vertex* v2, Vertex* tv1[], Edge* te1[])
 201 {
 202    int orig = NOTHING;
 203    for (int nro=0 ; orig == NOTHING && nro < QUAD4 ; nro++)
 204        if (q_vertex [nro] == v1)
 205            orig = nro;
 206
 207    if (orig==NOTHING)
 208       return HERR;
 209
 210    int nsp1 = (orig+1)       MODULO QUAD4;
 211    int nsm1 = (orig+QUAD4-1) MODULO QUAD4;
 212
 213    if (q_vertex [nsp1] == v2)
 214       {
 215       for (int nro=0 ; nro < QUAD4 ; nro++)
 216           {
 217           tv1 [nro] = q_vertex [(orig+nro) MODULO QUAD4];
 218           te1 [nro] = q_edge   [(orig+nro) MODULO QUAD4];
 219           }
 220       }
 221    else if (q_vertex [nsm1] == v2)
 222       {
 223       for (int nro=0 ; nro < QUAD4 ; nro++)
 224           {
 225           tv1 [nro] = q_vertex [(orig+QUAD4-nro) MODULO QUAD4];
 226           te1 [nro] = q_edge   [(orig+QUAD4-nro) MODULO QUAD4];
 227           }
 228       }
 229    else
 230       return 588;
 231
 232    return HOK;
 233 }
 234 // ======================================================== ordoVertex
 235 int Quad::ordoVertex (Vertex* v1, Vertex* v2, Vertex* tver[])
 236 {
 237    int orig = NOTHING;
 238    for (int nro=0 ; orig == NOTHING && nro < QUAD4 ; nro++)
 239        if (q_vertex [nro] == v1)
 240            orig = nro;
 241
 242    if (orig==NOTHING)
 243       return HERR;
 244
 245    int nsp1 = (orig+1)       MODULO QUAD4;
 246    int nsm1 = (orig+QUAD4-1) MODULO QUAD4;
 247
 248    if (q_vertex [nsp1] == v2)
 249       {
 250       for (int nro=0 ; nro < QUAD4 ; nro++)
 251           tver [nro] = q_vertex [(orig+nro) MODULO QUAD4];
 252       }
 253    else if (q_vertex [nsm1] == v2)
 254       {
 255       for (int nro=0 ; nro < QUAD4 ; nro++)
 256           tver [nro] = q_vertex [(orig+QUAD4-nro) MODULO QUAD4];
 257       }
 258    else
 259       return 588;
 260
 261    return HOK;
 262 }
 263 // ======================================================== ordonner
 264 int Quad::ordonner (Vertex* v1, Vertex* v2, Vertex* tver[], Edge* ted[])
 265 {
 266    tver [0] = tver [1] = tver [2] = tver [3] = NULL;
 267    ted  [0] = ted  [1] = ted  [2] = ted  [3] = NULL;
 268
 269    int ier = ordoVertex (v1, v2, tver);
 270    if (ier != HOK)
 271        return ier;
 272
 273    for (int nro=0 ; nro < QUAD4 ; nro++)
 274        ted [nro] = findEdge (tver[nro], tver [(nro+1) MODULO QUAD4]);
 275
 276    return HOK;
 277 }
 278 // ======================================================== getBrother
 279 Quad* Quad::getBrother (StrOrient* orient)
 280 {
 281 /* *****************************
 282    printf (" getBrother ");
 283    dump ();
 284    printf (" .. Base  : ");
 285    orient->v21->printName();
 286    orient->v22->printName();
 287    printf ("dir=%d, arete=", orient->dir);
 288   ***************************** */
 289
 290    int n21 = indexVertex (orient->v21);
 291    int n22 = indexVertex (orient->v22);
 292
 293    int sens  = n22 - n21;
 294    if (sens >  1) sens -= QUAD4;
 295    if (sens < -1) sens += QUAD4;
 296    if (sens*sens !=1) return NULL;
 297
 298    switch (orient->dir)
 299       {
 300       case OR_LEFT  : n22 = n21 - sens;
 301            break;
 302       case OR_RIGHT : n21 = n22 + sens;
 303            break;
 304       case OR_FRONT : n21 += 2;
 305                       n22 += 2;
 306            break;
 307       default : ;
 308       }
 309
 310    n21 = (n21 + QUAD4) MODULO QUAD4;
 311    n22 = (n22 + QUAD4) MODULO QUAD4;
 312
 313    orient->v21 = q_vertex [n21];
 314    orient->v22 = q_vertex [n22];
 315
 316    Edge* arete  = findEdge (orient->v21, orient->v22);
 317    // arete->printName("\n");
 318
 319    int nbfreres = arete->getNbrParents ();
 320
 321    for (int nq = 0 ; nq < nbfreres ; nq++)
 322        {
 323        Quad* next = arete->getParent (nq);
 324        if (next!=NULL && next != this )
 325           {
 326           int   nbp   = next->getNbrParents();
 327           Hexa* dad   = next->getParent(0);
 328           int   mark  = next->getMark();
 329           int   mark2 = dad ? dad->getMark() : IS_NONE;
 330
 331           if (nbp  == 1  && mark2 != IS_MARRIED && mark == IS_NONE)
 332              return next;
 333           // if (nbp  <= 1  && mark == IS_NONE)
 334              // return next;
 335           }
 336        }
 337    return NULL;
 338 }
 339 // ======================================================== coupler
 340 int Quad::coupler (Quad* other, StrOrient* orient, Elements* table)
 341 {
 342    if (other==NULL)
 343       return HERR;
 344    if (db)
 345       {
 346       std::cout << " Quads::coupler " << el_name << " -> " << other->getName ()
 347                 << std::endl;
 348       }
 349
 350    Hexa* hexa = other->getParent(0);
 351
 352    setMark (IS_MARRIED);
 353    other->setMark (IS_MARRIED);
 354    if (hexa != NULL)
 355        hexa->setMark (IS_MARRIED);
 356
 357    for (int ned = 0 ; ned < QUAD4 ; ned++)
 358        {
 359        Edge* arete  = q_edge[ned];
 360        int nbfreres = arete ->getNbrParents ();
 361        for (int nq = 0 ; nq < nbfreres ; nq++)
 362            {
 363            Quad* next = arete->getParent (nq);
 364            if (next!=NULL && next != this && next->getMark() > 0)
 365               {
 366               StrOrient new_ori (orient);
 367               new_ori.dir = OR_FRONT;
 368               Vertex* va  = arete->getVertex (V_AMONT);
 369               Vertex* vb  = arete->getVertex (V_AVAL);
 370
 371 //    On voit si un point de repere est conserve
 372               if (va == orient->v11)
 373                  {
 374                  new_ori.v12 = vb;
 375                  new_ori.dir += OR_LEFT;
 376                  }
 377               else if (vb == orient->v11)
 378                  {
 379                  new_ori.v12 = va;
 380                  new_ori.dir += OR_LEFT;
 381                  }
 382
 383               if (va == orient->v12)
 384                  {
 385                  new_ori.v11 = vb;
 386                  new_ori.dir += OR_RIGHT;
 387                  }
 388               else if (vb == orient->v12)
 389                  {
 390                  new_ori.v11 = va;
 391                  new_ori.dir += OR_RIGHT;
 392                  }
 393
 394               if (new_ori.dir == OR_FRONT)
 395                  {
 396                  if (definedBy (va, orient->v11))
 397                     {
 398                     new_ori.v11 = va;
 399                     new_ori.v12 = vb;
 400                     }
 401                  else
 402                     {
 403                     new_ori.v11 = vb;
 404                     new_ori.v12 = va;
 405                     }
 406                  }
 407
 408               int nro = next->getMark ();
 409               Quad* beauf = other->getBrother (&new_ori);
 410               int ier = table->coupler (nro, beauf, &new_ori);
 411               if (ier != HOK)
 412                  return ier;
 413               ier = next->coupler (beauf, &new_ori, table);
 414               if (ier != HOK)
 415                  return ier;
 416               }
 417            }
 418        }
 419    return HOK;
 420 }
 421 // ======================================================== getOpposVertex
 422 Vertex* Quad::getOpposVertex (Vertex* start)
 423 {
 424    int  na = indexVertex (start);
 425    if (na==NOTHING)
 426       return NULL;
 427    return  q_vertex [(na+2) MODULO QUAD4];
 428 }
 429 // ======================================================== getOpposEdge
 430 Edge* Quad::getOpposEdge (Edge* start, int& sens)
 431 {
 432    sens = 1;
 433    int  na = indexVertex (start->getVertex (V_AMONT));
 434    int  nb = indexVertex (start->getVertex (V_AVAL));
 435
 436    Vertex* vaprim = q_vertex [(nb+2) MODULO QUAD4];
 437    Vertex* vbprim = q_vertex [(na+2) MODULO QUAD4];
 438
 439    for (int ned = 0 ; ned < QUAD4 ; ned++)
 440        {
 441        if (   q_edge[ned]->getVertex(V_AMONT) == vaprim
 442            && q_edge[ned]->getVertex(V_AVAL ) == vbprim)
 443            {
 444            sens = 1;
 445            return q_edge[ned];
 446            }
 447        else if (   q_edge[ned]->getVertex(V_AMONT) == vbprim
 448                 && q_edge[ned]->getVertex(V_AVAL ) == vaprim)
 449            {
 450            sens = -1;
 451            return q_edge[ned];
 452            }
 453        }
 454    //             TODO : traiter l'erreur
 455    std::cout << " ... Probleme dans Quad::getOpposedEdge :" << std::endl;
 456    HexDisplay (el_name);
 457    PutName (start);
 458    PutName (vaprim);
 459    PutName (vbprim);
 460    HexDisplay (na);
 461    HexDisplay (nb);
 462    dumpPlus ();
 463
 464    for (int ned = 0 ; ned < QUAD4 ; ned++)
 465        q_edge[ned]->dump();
 466
 467    return NULL;
 468 }
 469 // ========================================================= saveXml
 470 void Quad::saveXml (XmlWriter* xml)
 471 {
 472    char buffer[12];
 473    std::string edges;
 474
 475    for (int nro=0 ; nro<QUAD4 ; nro++)
 476        {
 477        if (nro>0) edges += " ";
 478        edges += q_edge[nro]->getName(buffer);
 479        }
 480
 481    xml->openMark     ("Quad");
 482    xml->addAttribute ("id",    getName (buffer));
 483    xml->addAttribute ("edges", edges);
 484    if (el_name!=buffer)
 485        xml->addAttribute ("name", el_name);
 486    xml->closeMark ();
 487
 488    int nbass = tab_assoc.size();
 489    for (int nro=0 ; nro<nbass ; nro++)
 490        if (tab_assoc[nro] != NULL)
 491            tab_assoc[nro]->saveXml (xml);
 492 }
 493 // ======================================================== replaceEdge
 494 void Quad::replaceEdge (Edge* old, Edge* par)
 495 {
 496    for (int nro=0 ; nro<QUAD4 ; nro++)
 497        {
 498        if (q_edge[nro]==old)
 499            {
 500            q_edge[nro] = par;
 501            if (debug())
 502               {
 503               printf (" Dans ");
 504               printName ();
 505               printf (" [%d], ", nro);
 506               old->printName (" est remplace par ");
 507               par->printName ("\n");
 508               }
 509            }
 510        }
 511 }
 512 // ======================================================== replaceVertex
 513 void Quad::replaceVertex (Vertex* old, Vertex* par)
 514 {
 515    for (int nro=0 ; nro<QUAD4 ; nro++)
 516        {
 517        if (q_vertex [nro]==old)
 518           {
 519           q_vertex [nro] = par;
 520            if (debug())
 521               {
 522               printf (" Dans ");
 523               printName ();
 524               printf (" [%d], ", nro);
 525               old->printName (" est remplace par ");
 526               par->printName ("\n");
 527               }
 528           }
 529        }
 530 }
 531 // ======================================================== dump
 532 void Quad::dump ()
 533 {
 534    printName(" = (");
 535    if (NOT isHere ())
 536       {
 537       printf ("*** deleted ***)\n");
 538       return;
 539       }
 540
 541    for (int nro=0 ; nro<QUAD4 ; nro++)
 542         PrintName (q_edge[nro]);
 543    printf (")\n");
 544
 545    printf ("        (");
 546    for (int nro=0 ; nro<QUAD4 ; nro++)
 547         PrintName (q_vertex[nro]);
 548    printf (")");
 549
 550    dumpRef ();
 551    Real3 cg;
 552    getCenter (cg);
 553    printf ("     -> (%g, %g, %g)\n", cg[dir_x], cg[dir_y], cg[dir_z]);
 554 }
 555 // ======================================================== dumpPlus
 556 void Quad::dumpPlus ()
 557 {
 558    dump ();
 559    if (NOT isHere ())
 560       return;
 561
 562    for (int nro=0 ; nro < QUAD4 ; nro++)
 563        {
 564        Vertex* pv = q_vertex[nro];
 565        printf ( "    ");
 566        if (pv!=NULL)
 567           {
 568           pv->printName ("");
 569           printf ( " (%g, %g, %g)\n", pv->getX(),  pv->getY(),  pv->getZ());
 570           }
 571        else
 572           {
 573           printf ( "NULL\n");
 574           }
 575        }
 576 }
 577 // ======================================================== getOpposEdge (2)
 578 Edge* Quad::getOpposEdge (Edge* start)
 579 {
 580    int  na = indexEdge (start);
 581    if (na<0)
 582       return NULL;
 583    return q_edge [(na+2) MODULO QUAD4];
 584 }
 585 // ======================================================== getPerpendicular
 586 Edge* Quad::getPerpendicular (Edge* arete, Vertex* node)
 587 {
 588    int na = indexEdge (arete);
 589    if (na<0)
 590       return NULL;
 591
 592    int nv = arete->index (node);
 593    if (nv<0)
 594       return NULL;
 595    Edge* perp = q_edge [(na+1) MODULO QUAD4];
 596
 597    nv = perp->index (node);
 598    if (nv>=0)
 599       return perp;
 600
 601    perp = q_edge [(na+3) MODULO QUAD4];
 602    nv   = perp->index (node);
 603    if (nv>=0)
 604       return perp;
 605    else
 606       return NULL;
 607 }
 608 static cpchar t_ori[] = {"Q_INSIDE", "Q_DIRECT", "Q_INVERSE", "Q_UNDEF"};
 609 // ======================================================== setOrientation
 610 void Quad::setOrientation (int ori)
 611 {
 612     q_orientation = ori;
 613     if (db && (ori==Q_DIRECT || ori==Q_INVERSE))
 614        printf (" %s = %s\n", el_name.c_str(), t_ori [ q_orientation ]);
 615 }
 616 // ======================================================== setOrientation
 617 int Quad::setOrientation ()
 618 {
 619     q_orientation = Q_INSIDE;
 620     if (getNbrParents() != 1)
 621        return q_orientation;
 622
 623     Real3 cg, orig, pi, pj, vi, vj, vk;
 624
 625     Hexa* hexa = getParent(0);
 626     hexa->getCenter (cg);
 627
 628 /********************************************************************
 629     printName (" = ");
 630
 631     for (int np=0 ; np < QUAD4 ; np++)
 632         {
 633         q_vertex [np        ] -> getPoint (orig);
 634         q_vertex [(np+1) % 4] -> getPoint (pi);
 635         q_vertex [(np+3) % 4] -> getPoint (pj);
 636
 637         calc_vecteur (orig, pi, vi);
 638         calc_vecteur (orig, pj, vj);
 639         calc_vecteur (orig, cg, vk);
 640         double pmixte = prod_mixte (vi, vj, vk);
 641         q_orientation = pmixte > ZEROR ? Q_DIRECT : Q_INVERSE;
 642         if (pmixte>0) printf (">");
 643            else       printf ("<");
 644         }
 645
 646     printf ("\n");
 647     return;
 648   ******************************************************************* */
 649     q_vertex [0] -> getPoint (orig);
 650     q_vertex [1] -> getPoint (pi);
 651     q_vertex [3] -> getPoint (pj);
 652
 653     calc_vecteur (orig, pi, vi);
 654     calc_vecteur (orig, pj, vj);
 655     calc_vecteur (cg, orig, vk);
 656
 657     double pmixte = prod_mixte (vi, vj, vk);
 658     q_orientation = pmixte > ZEROR ? Q_DIRECT : Q_INVERSE;
 659     if (db)
 660        printf (" %s := %s\n", el_name.c_str(), t_ori [ q_orientation ]);
 661     return q_orientation;
 662 }
 663 // ========================================================== clearAssociation
 664 void Quad::clearAssociation ()
 665 {
 666    tab_assoc.clear ();
 667    is_associated = false;
 668 }
 669 // ========================================================== addAssociation
 670 int Quad::addAssociation (NewShape* geom, int subid)
 671 {
 672    if (geom == NULL)
 673       return HERR;
 674
 675    FaceShape* face = geom->findFace (subid);
 676    int ier = addAssociation (face);
 677
 678    return ier;
 679 }
 680 // ========================================================== addAssociation
 681 int Quad::addAssociation (FaceShape* face)
 682 {
 683    if (face == NULL)
 684       return HERR;
 685
 686    face->addAssociation (this);
 687    tab_assoc.push_back (face);
 688    is_associated = true;
 689    return HOK;
 690 }
 691 // ========================================================== getAssociation
 692 FaceShape* Quad::getAssociation (int nro)
 693 {
 694    if (nro < 0 || nro >= (int)tab_assoc.size())
 695       return NULL;
 696
 697    return tab_assoc [nro];
 698 }
 699 // ======================================================== commonEdge
 700 Edge* Quad::commonEdge (Quad* other)
 701 {
 702    for (int ne1=0 ; ne1<QUAD4 ; ne1++)
 703        for (int ne2=0 ; ne2<QUAD4 ; ne2++)
 704            if (q_edge [ne1] == other->q_edge [ne2])
 705               return q_edge [ne1];
 706
 707    return NULL;
 708 }
 709 // ======================================================== Inter
 710 int Quad::inter (Quad* other, int& nother)
 711 {
 712    for (int ne1=0 ; ne1<QUAD4 ; ne1++)
 713        for (int ne2=0 ; ne2<QUAD4 ; ne2++)
 714            if (q_edge [ne1] == other->q_edge [ne2])
 715               {
 716               nother = ne2;
 717               return  ne1;
 718               }
 719
 720    nother = NOTHING;
 721    return NOTHING;
 722 }
 723 // ======================================================== Inter (2)
 724 Edge* Quad::inter (Quad* other)
 725 {
 726    for (int ne1=0 ; ne1<QUAD4 ; ne1++)
 727        for (int ne2=0 ; ne2<QUAD4 ; ne2++)
 728            if (q_edge [ne1] == other->q_edge [ne2])
 729               return  q_edge [ne1];
 730    return NULL;
 731 }
 732 // ============================================================ definedBy (v)
 733 bool Quad::definedBy  (Vertex* v1, Vertex* v2)
 734 {
 735    for (int n1=0 ; n1< QUAD4 ; n1++)
 736        if (v1 == q_vertex[n1] && v2 == q_vertex[(n1+2) MODULO QUAD4])
 737           return true;
 738
 739    return false;
 740 }
 741 // ============================================================ definedBy (e)
 742 bool Quad::definedBy  (Edge* e1, Edge* e2)
 743 {
 744    if (e1==e2 || BadElement (e1) || BadElement (e2))
 745       return false;
 746
 747    bool f1=false, f2=false;
 748    for (int ned=0 ; ned< QUAD4 ; ned++)
 749        if      (e1 == q_edge[ned]) f1 = true;
 750        else if (e2 == q_edge[ned]) f2 = true;
 751     // if (e1 == q_edge[ned] && e2 == q_edge[(ned+2) MODULO QUAD4]) return true;
 752
 753    return f1 && f2;
 754 }
 755 // =============================================================== findEdge
 756 Edge* Quad::findEdge (Vertex* v1, Vertex* v2)
 757 {
 758    for (int nro=0 ; nro< QUAD4 ; nro++)
 759        {
 760        Vertex* va = q_edge[nro]->getVertex(V_AMONT) ;
 761        Vertex* vb = q_edge[nro]->getVertex(V_AVAL) ;
 762        if ((v1==va && v2==vb) || (v1==vb && v2==va))
 763            return q_edge [nro];
 764        }
 765
 766    return NULL;
 767 }
 768 // =============================================================== indexVertex
 769 int Quad::indexVertex  (Vertex* elt)
 770 {
 771    for (int n1=0 ; n1< QUAD4 ; n1++)
 772        if (elt == q_vertex[n1])
 773           return n1;
 774
 775    return NOTHING;
 776 }
 777 // =============================================================== indexEdge
 778 int Quad::indexEdge  (Edge* elt)
 779 {
 780    for (int n1=0 ; n1< QUAD4 ; n1++)
 781        if (elt == q_edge[n1])
 782           return n1;
 783
 784    return NOTHING;
 785 }
 786 // =============================================================== setColor
 787 void Quad::setColor  (double val)
 788 {
 789    for (int n1=0 ; n1< QUAD4 ; n1++)
 790        q_vertex[n1] -> setColor (val);
 791 }
 792 // =============================================================== duplicate
 793 void Quad::duplicate  ()
 794 {
 795    q_orientation  = Q_UNDEFINED;
 796    q_clone = new Quad (GetClone (q_edge [E_A]),
 797                        GetClone (q_edge [E_B]),
 798                        GetClone (q_edge [E_C]),
 799                        GetClone (q_edge [E_D]));
 800    q_clone->tab_assoc  = tab_assoc;
 801 }
 802 // ============================================================ nearestVertex
 803 Vertex* Quad::nearestVertex (Vertex* other)
 804 {
 805    if (BadElement (other))
 806       return NULL;
 807
 808    Vertex* vsol = NULL;
 809    double  dmin = 1e+77;
 810    int nbre = countVertex ();
 811    for (int nro=0 ; nro<nbre ; nro++)
 812        {
 813        Vertex* vert = getVertex (nro);
 814        double  dist = other->dist2 (vert);
 815        if (dist < dmin)
 816           {
 817           dmin = dist;
 818           vsol = vert;
 819           }
 820        }
 821    return vsol;
 822 }
 823 // =============================================================== reorienter
 824 void Quad::reorienter  ()
 825 {
 826    if (q_orientation != Q_INVERSE)
 827          return;
 828
 829    Edge* edge   = q_edge [E_B];
 830    q_edge [E_B] = q_edge [E_D];
 831    q_edge [E_B] = edge;
 832
 833    q_vertex [E_A] = q_edge [E_D]->commonVertex(q_edge [E_A]);
 834    q_vertex [E_B] = q_edge [E_A]->commonVertex(q_edge [E_B]);
 835    q_vertex [E_C] = q_edge [E_B]->commonVertex(q_edge [E_C]);
 836    q_vertex [E_D] = q_edge [E_C]->commonVertex(q_edge [E_D]);
 837
 838    if (db)
 839        printf (" %s est reoriente\n", el_name.c_str());
 840    q_orientation = Q_DIRECT;
 841 }
 842 // =============================================================== getCenter
 843 double* Quad::getCenter (double* center)
 844 {
 845    center[dir_x] = center[dir_y] = center[dir_z] =  0;
 846    if (BadElement (this))
 847       return NULL;
 848
 849    for (int nv=0 ; nv<QUAD4 ; nv++)
 850        {
 851        if (BadElement (q_vertex [nv]))
 852            return NULL;
 853        center [dir_x] += q_vertex[nv]->getX()/4;
 854        center [dir_y] += q_vertex[nv]->getY()/4;
 855        center [dir_z] += q_vertex[nv]->getZ()/4;
 856        }
 857     return center;
 858 }
 859 // =============================================================== dist2
 860 double Quad::dist2 (double* point)
 861 {
 862    Real3 center;
 863    getCenter (center);
 864    double d2 = carre (point[dir_x] - center[dir_x])
 865              + carre (point[dir_y] - center[dir_y])
 866              + carre (point[dir_z] - center[dir_z]) ;
 867    return d2;
 868 }
 869 // =============================================================== opposedHexa
 870 Hexa* Quad::opposedHexa (Hexa* hexa)
 871 {
 872    int nbre = getNbrParents ();
 873    for (int nro=0 ; nro <nbre ; ++nro)
 874        {
 875        Hexa* dad = getParent (nro);
 876        if (dad!= NULL && dad->isValid() && dad != hexa)
 877            return dad;
 878        }
 879    return NULL;
 880 }
 881 END_NAMESPACE_HEXA