src/HEXABLOCK/HexQuad.cxx

   1
   2 // C++ : Gestion des Quadrangles
   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 #include "HexQuad.hxx"
  23
  24 #include "HexDocument.hxx"
  25 #include "HexHexa.hxx"
  26 #include "HexElements.hxx"
  27
  28 #include "HexXmlWriter.hxx"
  29 #include "HexShape.hxx"
  30
  31 BEGIN_NAMESPACE_HEXA
  32
  33 // ======================================================== Constructeur
  34 Quad::Quad (Vertex* va, Vertex* vb, Vertex* vc, Vertex* vd)
  35     : EltBase (va->dad(), EL_QUAD)
  36 {
  37    q_vertex [E_A] = va;
  38    q_vertex [E_B] = vb;
  39    q_vertex [E_C] = vc;
  40    q_vertex [E_D] = vd;
  41    q_clone        = NULL;
  42    q_orientation  = Q_UNDEFINED;
  43
  44    for (int nro=0 ; nro<QUAD4 ; nro++)
  45        {
  46        q_edge [nro] = new Edge (q_vertex[nro],
  47                                 q_vertex[(nro+1) MODULO QUAD4]);
  48
  49        if (BadElement (q_vertex [nro]) || BadElement (q_edge [nro]))
  50           setError ();
  51        else
  52           for (int nv=nro+1 ; nv<QUAD4 ; nv++)
  53               if (q_vertex[nv] == q_vertex[nro])
  54                  setError ();
  55        }
  56
  57    majReferences ();
  58 }
  59 // ======================================================== Constructeur bis
  60 Quad::Quad (Edge* ea, Edge* eb, Edge* ec, Edge* ed)
  61     : EltBase (ea->dad(), EL_QUAD)
  62 {
  63    q_edge [E_A] = ea;
  64    q_edge [E_B] = eb;
  65    q_edge [E_C] = ec;
  66    q_edge [E_D] = ed;
  67    q_clone       = NULL;
  68    q_orientation = Q_UNDEFINED;
  69
  70    for (int nro=0 ; nro<QUAD4 ; nro++)
  71        {
  72        int prec = (nro+1) MODULO QUAD4;
  73        Vertex* node = NULL;
  74
  75        if (BadElement (q_edge[nro]))
  76           setError ();
  77        else
  78           {
  79           for (int nv=nro+1 ; nv<QUAD4 ; nv++)
  80               if (q_edge[nv] == q_edge[nro])
  81                   setError ();
  82           int nc  = q_edge[nro] -> inter (q_edge[prec]);
  83           if (nc>=0)
  84              node = q_edge[nro]->getVertex (nc);
  85           else
  86              setError (888);
  87           }
  88        q_vertex [prec] = node;
  89        }
  90
  91    if (isBad())
  92       {
  93       printf (" +++++++++++++++++++++++++++++++++++++++++++ \n");
  94       printf (" +++ Quadrangle impossible \n");
  95       printf (" +++++++++++++++++++++++++++++++++++++++++++ \n");
  96       dump ();
  97       printf (" +++++++++++++++++++++++++++++++++++++++++++ \n");
  98       // el_root->dump ();
  99       for (int ned=0; ned<QUAD4; ned++)
 100           {
 101           q_edge[ned]->dumpPlus ();
 102           }
 103       HexDump (q_vertex[0]);
 104       HexDump (q_vertex[1]);
 105       HexDump (q_vertex[2]);
 106       HexDump (q_vertex[3]);
 107
 108       printf (" +++++++++++++++++++++++++++++++++++++++++++ \n");
 109       fatal_error ("Quadrangle impossible");
 110       }
 111
 112    majReferences ();
 113 }
 114 // ======================================================== Constructeur bis
 115 Quad::Quad (Quad* other)
 116     : EltBase (other->dad(), EL_QUAD)
 117 {
 118    for (int nro=0 ; nro<QUAD4 ; nro++)
 119        {
 120        q_edge   [nro] = NULL;
 121        q_vertex [nro] = NULL;
 122        }
 123    q_orientation = Q_UNDEFINED;
 124    q_clone       = NULL;
 125 }
 126 // ========================================================= majReferences
 127 void Quad::majReferences ()
 128 {
 129    for (int nro=0 ; nro<QUAD4 ; nro++)
 130        q_edge [nro] -> addParent (this);
 131 }
 132 // ========================================================= getParent
 133 Hexa* Quad::getParent  (int nro)
 134 {
 135    return static_cast <Hexa*> (getFather (nro));
 136 }
 137 // ======================================================== anaMerge
 138 int Quad::anaMerge (Vertex* v1, Vertex* v2, Vertex* tv1[], Edge* te1[])
 139 {
 140    int orig = NOTHING;
 141    for (int nro=0 ; orig == NOTHING && nro < QUAD4 ; nro++)
 142        if (q_vertex [nro] == v1)
 143            orig = nro;
 144
 145    if (orig==NOTHING)
 146       return HERR;
 147
 148    int nsp1 = (orig+1)       MODULO QUAD4;
 149    int nsm1 = (orig+QUAD4-1) MODULO QUAD4;
 150
 151    if (q_vertex [nsp1] == v2)
 152       {
 153       for (int nro=0 ; nro < QUAD4 ; nro++)
 154           {
 155           tv1 [nro] = q_vertex [(orig+nro) MODULO QUAD4];
 156           te1 [nro] = q_edge   [(orig+nro) MODULO QUAD4];
 157           }
 158       }
 159    else if (q_vertex [nsm1] == v2)
 160       {
 161       for (int nro=0 ; nro < QUAD4 ; nro++)
 162           {
 163           tv1 [nro] = q_vertex [(orig+QUAD4-nro) MODULO QUAD4];
 164           te1 [nro] = q_edge   [(orig+QUAD4-nro) MODULO QUAD4];
 165           }
 166       }
 167    else
 168       return 588;
 169
 170    return HOK;
 171 }
 172 // ======================================================== ordoVertex
 173 int Quad::ordoVertex (Vertex* v1, Vertex* v2, Vertex* tv1[])
 174 {
 175    int orig = NOTHING;
 176    for (int nro=0 ; orig == NOTHING && nro < QUAD4 ; nro++)
 177        if (q_vertex [nro] == v1)
 178            orig = nro;
 179
 180    if (orig==NOTHING)
 181       return HERR;
 182
 183    int nsp1 = (orig+1)       MODULO QUAD4;
 184    int nsm1 = (orig+QUAD4-1) MODULO QUAD4;
 185
 186    if (q_vertex [nsp1] == v2)
 187       {
 188       for (int nro=0 ; nro < QUAD4 ; nro++)
 189           tv1 [nro] = q_vertex [(orig+nro) MODULO QUAD4];
 190       }
 191    else if (q_vertex [nsm1] == v2)
 192       {
 193       for (int nro=0 ; nro < QUAD4 ; nro++)
 194           tv1 [nro] = q_vertex [(orig+QUAD4-nro) MODULO QUAD4];
 195       }
 196    else
 197       return 588;
 198
 199    return HOK;
 200 }
 201 // ======================================================== getBrother
 202 Quad* Quad::getBrother (StrOrient* orient)
 203 {
 204 /* *****************************
 205    printf (" getBrother ");
 206    dump ();
 207    printf (" .. Base  : ");
 208    orient->v21->printName();
 209    orient->v22->printName();
 210    printf ("dir=%d, arete=", orient->dir);
 211   ***************************** */
 212
 213    int n21 = indexVertex (orient->v21);
 214    int n22 = indexVertex (orient->v22);
 215
 216    int sens  = n22 - n21;
 217    if (sens >  1) sens -= QUAD4;
 218    if (sens < -1) sens += QUAD4;
 219    if (sens*sens !=1) return NULL;
 220
 221    switch (orient->dir)
 222       {
 223       case OR_LEFT  : n22 = n21 - sens;
 224            break;
 225       case OR_RIGHT : n21 = n22 + sens;
 226            break;
 227       case OR_FRONT : n21 += 2;
 228                       n22 += 2;
 229            break;
 230       default : ;
 231       }
 232
 233    n21 = (n21 + QUAD4) MODULO QUAD4;
 234    n22 = (n22 + QUAD4) MODULO QUAD4;
 235
 236    orient->v21 = q_vertex [n21];
 237    orient->v22 = q_vertex [n22];
 238
 239    Edge* arete  = findEdge (orient->v21, orient->v22);
 240    // arete->printName("\n");
 241
 242    int nbfreres = arete->getNbrParents ();
 243
 244    for (int nq = 0 ; nq < nbfreres ; nq++)
 245        {
 246        Quad* next = arete->getParent (nq);
 247        if (next!=NULL && next != this )
 248           {
 249           int   nbp   = next->getNbrParents();
 250           Hexa* dad   = next->getParent(0);
 251           int   mark  = next->getMark();
 252           int   mark2 = dad ? dad->getMark() : IS_NONE;
 253
 254           if (nbp  <= 1  && mark2 != IS_MARRIED && mark == IS_NONE)
 255              return next;
 256           // if (nbp  <= 1  && mark == IS_NONE)
 257              // return next;
 258           }
 259        }
 260    return NULL;
 261 }
 262 // ======================================================== coupler
 263 int Quad::coupler (Quad* other, StrOrient* orient, Elements* table)
 264 {
 265    if (other==NULL)
 266       return HERR;
 267
 268    Hexa* hexa = other->getParent(0);
 269
 270    setMark (IS_MARRIED);
 271    other->setMark (IS_MARRIED);
 272    if (hexa != NULL)
 273        hexa->setMark (IS_MARRIED);
 274
 275    for (int ned = 0 ; ned < QUAD4 ; ned++)
 276        {
 277        Edge* arete  = q_edge[ned];
 278        int nbfreres = arete ->getNbrParents ();
 279        for (int nq = 0 ; nq < nbfreres ; nq++)
 280            {
 281            Quad* next = arete->getParent (nq);
 282            if (next!=NULL && next != this && next->getMark() > 0)
 283               {
 284               StrOrient new_ori (orient);
 285               new_ori.dir = OR_FRONT;
 286               Vertex* va  = arete->getVertex (V_AMONT);
 287               Vertex* vb  = arete->getVertex (V_AVAL);
 288
 289 //    On voit si un point de repere est conserve
 290               if (va == orient->v11)
 291                  {
 292                  new_ori.v12 = vb;
 293                  new_ori.dir += OR_LEFT;
 294                  }
 295               else if (vb == orient->v11)
 296                  {
 297                  new_ori.v12 = va;
 298                  new_ori.dir += OR_LEFT;
 299                  }
 300
 301               if (va == orient->v12)
 302                  {
 303                  new_ori.v11 = vb;
 304                  new_ori.dir += OR_RIGHT;
 305                  }
 306               else if (vb == orient->v12)
 307                  {
 308                  new_ori.v11 = va;
 309                  new_ori.dir += OR_RIGHT;
 310                  }
 311
 312               if (new_ori.dir == OR_FRONT)
 313                  {
 314                  if (definedBy (va, orient->v11))
 315                     {
 316                     new_ori.v11 = va;
 317                     new_ori.v12 = vb;
 318                     }
 319                  else
 320                     {
 321                     new_ori.v11 = vb;
 322                     new_ori.v12 = va;
 323                     }
 324                  }
 325
 326               int nro = next->getMark ();
 327               Quad* beauf = other->getBrother (&new_ori);
 328               int ier = table->coupler (nro, beauf, &new_ori);
 329               if (ier != HOK)
 330                  return ier;
 331               ier = next->coupler (beauf, &new_ori, table);
 332               if (ier != HOK)
 333                  return ier;
 334               }
 335            }
 336        }
 337    return HOK;
 338 }
 339 // ======================================================== getOpposVertex
 340 Vertex* Quad::getOpposVertex (Vertex* start)
 341 {
 342    int  na = indexVertex (start);
 343    if (na==NOTHING)
 344       return NULL;
 345    return  q_vertex [(na+2) MODULO QUAD4];
 346 }
 347 // ======================================================== getOpposEdge
 348 Edge* Quad::getOpposEdge (Edge* start, int& sens)
 349 {
 350    sens = 1;
 351    int  na = indexVertex (start->getVertex (V_AMONT));
 352    int  nb = indexVertex (start->getVertex (V_AVAL));
 353
 354    Vertex* vaprim = q_vertex [(nb+2) MODULO QUAD4];
 355    Vertex* vbprim = q_vertex [(na+2) MODULO QUAD4];
 356
 357    for (int ned = 0 ; ned < QUAD4 ; ned++)
 358        {
 359        if (   q_edge[ned]->getVertex(V_AMONT) == vaprim
 360            && q_edge[ned]->getVertex(V_AVAL ) == vbprim)
 361            {
 362            sens = 1;
 363            return q_edge[ned];
 364            }
 365        else if (   q_edge[ned]->getVertex(V_AMONT) == vbprim
 366                 && q_edge[ned]->getVertex(V_AVAL ) == vaprim)
 367            {
 368            sens = -1;
 369            return q_edge[ned];
 370            }
 371        }
 372    //             TODO : traiter l'erreur
 373    cout << " ... Probleme dans Quad::getOpposedEdge :" << endl;
 374    PutName (start);
 375    PutName (vaprim);
 376    PutName (vbprim);
 377    HexDisplay (na);
 378    HexDisplay (nb);
 379    dumpPlus ();
 380
 381    for (int ned = 0 ; ned < QUAD4 ; ned++)
 382        q_edge[ned]->dump();
 383
 384    return NULL;
 385 }
 386 // ========================================================= saveXml
 387 void Quad::saveXml (XmlWriter* xml)
 388 {
 389    char buffer[12];
 390    string edges;
 391
 392    for (int nro=0 ; nro<QUAD4 ; nro++)
 393        {
 394        if (nro>0) edges += " ";
 395        edges += q_edge[nro]->getName(buffer);
 396        }
 397
 398    xml->openMark     ("Quad");
 399    xml->addAttribute ("id",    getName (buffer));
 400    xml->addAttribute ("edges", edges);
 401    if (el_name!=buffer)
 402        xml->addAttribute ("name", el_name);
 403    xml->closeMark ();
 404
 405    int nbass = tab_assoc.size();
 406    for (int nro=0 ; nro<nbass ; nro++)
 407        if (tab_assoc[nro] != NULL)
 408            tab_assoc[nro]->saveXml (xml);
 409 }
 410 // ======================================================== replaceEdge
 411 void Quad::replaceEdge (Edge* old, Edge* par)
 412 {
 413    for (int nro=0 ; nro<QUAD4 ; nro++)
 414        {
 415        if (q_edge[nro]==old)
 416            {
 417            q_edge[nro] = par;
 418            if (debug())
 419               {
 420               printf (" Dans ");
 421               printName ();
 422               printf (" [%d], ", nro);
 423               old->printName (" est remplace par ");
 424               par->printName ("\n");
 425               }
 426            }
 427        }
 428 }
 429 // ======================================================== replaceVertex
 430 void Quad::replaceVertex (Vertex* old, Vertex* par)
 431 {
 432    for (int nro=0 ; nro<QUAD4 ; nro++)
 433        {
 434        if (q_vertex [nro]==old)
 435           {
 436           q_vertex [nro] = par;
 437            if (debug())
 438               {
 439               printf (" Dans ");
 440               printName ();
 441               printf (" [%d], ", nro);
 442               old->printName (" est remplace par ");
 443               par->printName ("\n");
 444               }
 445           }
 446        }
 447 }
 448 // ======================================================== dump
 449 void Quad::dump ()
 450 {
 451    printName(" = (");
 452    if (NOT isHere ())
 453       {
 454       printf ("*** deleted ***)\n");
 455       return;
 456       }
 457
 458    for (int nro=0 ; nro<QUAD4 ; nro++)
 459         PrintName (q_edge[nro]);
 460    printf (")\n");
 461
 462    printf ("        (");
 463    for (int nro=0 ; nro<QUAD4 ; nro++)
 464         PrintName (q_vertex[nro]);
 465    printf (")");
 466
 467    dumpRef ();
 468 }
 469 // ======================================================== dumpPlus
 470 void Quad::dumpPlus ()
 471 {
 472    dump ();
 473    if (NOT isHere ())
 474       return;
 475
 476    for (int nro=0 ; nro < QUAD4 ; nro++)
 477        {
 478        Vertex* pv = q_vertex[nro];
 479        printf ( "    ");
 480        if (pv!=NULL)
 481           {
 482           pv->printName ("");
 483           printf ( " (%g, %g, %g)\n", pv->getX(),  pv->getY(),  pv->getZ());
 484           }
 485        else
 486           {
 487           printf ( "NULL\n");
 488           }
 489        }
 490 }
 491 // ======================================================== getOpposEdge (2)
 492 Edge* Quad::getOpposEdge (Edge* start)
 493 {
 494    int  na = indexEdge (start);
 495    if (na<0)
 496       return NULL;
 497    return q_edge [(na+2) MODULO QUAD4];
 498 }
 499 // ======================================================== getPerpendicular
 500 Edge* Quad::getPerpendicular (Edge* arete, Vertex* node)
 501 {
 502    int na = indexEdge (arete);
 503    if (na<0)
 504       return NULL;
 505
 506    int nv = arete->index (node);
 507    if (nv<0)
 508       return NULL;
 509    Edge* perp = q_edge [(na+1) MODULO QUAD4];
 510
 511    nv = perp->index (node);
 512    if (nv>=0)
 513       return perp;
 514
 515    perp = q_edge [(na+3) MODULO QUAD4];
 516    nv   = perp->index (node);
 517    if (nv>=0)
 518       return perp;
 519    else
 520       return NULL;
 521 }
 522 static cpchar t_ori[] = {"Q_INSIDE", "Q_DIRECT", "Q_INVERSE", "Q_UNDEF"};
 523 // ======================================================== setOrientation
 524 void Quad::setOrientation (int ori)
 525 {
 526     q_orientation = ori;
 527     if (ori==Q_DIRECT || ori==Q_INVERSE)
 528        printf (" %s = %s\n", el_name.c_str(), t_ori [ q_orientation ]);
 529 }
 530 // ======================================================== setOrientation
 531 int Quad::setOrientation ()
 532 {
 533     q_orientation = Q_INSIDE;
 534     if (getNbrParents() != 1)
 535        return q_orientation;
 536
 537     Real3 cg, orig, pi, pj, vi, vj, vk;
 538
 539     Hexa* hexa = getParent(0);
 540     hexa->getCenter (cg);
 541
 542 /********************************************************************
 543     printName (" = ");
 544
 545     for (int np=0 ; np < QUAD4 ; np++)
 546         {
 547         q_vertex [np        ] -> getPoint (orig);
 548         q_vertex [(np+1) % 4] -> getPoint (pi);
 549         q_vertex [(np+3) % 4] -> getPoint (pj);
 550
 551         calc_vecteur (orig, pi, vi);
 552         calc_vecteur (orig, pj, vj);
 553         calc_vecteur (orig, cg, vk);
 554         double pmixte = prod_mixte (vi, vj, vk);
 555         q_orientation = pmixte > ZEROR ? Q_DIRECT : Q_INVERSE;
 556         if (pmixte>0) printf (">");
 557            else       printf ("<");
 558         }
 559
 560     printf ("\n");
 561     return;
 562   ******************************************************************* */
 563     q_vertex [0] -> getPoint (orig);
 564     q_vertex [1] -> getPoint (pi);
 565     q_vertex [3] -> getPoint (pj);
 566
 567     calc_vecteur (orig, pi, vi);
 568     calc_vecteur (orig, pj, vj);
 569     calc_vecteur (cg, orig, vk);
 570
 571     double pmixte = prod_mixte (vi, vj, vk);
 572     q_orientation = pmixte > ZEROR ? Q_DIRECT : Q_INVERSE;
 573     printf (" %s = %s\n", el_name.c_str(), t_ori [ q_orientation ]);
 574     return q_orientation;
 575 }
 576 // ========================================================== setAssociation
 577 void Quad::setAssociation (Shape* forme)
 578 {
 579    clearAssociation ();
 580    addAssociation (forme);
 581 }
 582 END_NAMESPACE_HEXA
 583