src/expde/formulas/subdiIII.h

Go to the documentation of this file.

//    expde: expression templates for partial differential equations.
//    Copyright (C) 2001  Christoph Pflaum
//    This program is free software; you can redistribute it and/or modify
//    it under the terms of the GNU General Public License as published by
//    the Free Software Foundation; either version 2 of the License, or
//    (at your option) any later version.
//
//    This program is distributed in the hope that it will be useful,
//    but WITHOUT ANY WARRANTY; without even the implied warranty of
//    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
//    GNU General Public License for more details.
//
//                 SEE  Notice1.doc made by 
//                 LAWRENCE LIVERMORE NATIONAL LABORATORY
//

// ------------------------------------------------------------
//
// subdiIII.h
//
// ------------------------------------------------------------

#ifndef SUBDIIII_H_
#define SUBDIIII_H_

void Act_on_subdivision_partIII(Bo_description* Description, 
                                Calculator_FE* Act,
                                bool* corners, bool& stop) {

  int i;

  Edge_Corner_point poi_1,   poi_2,   poi_3,   poi_4;
  Edge_Corner_point poi_5,   poi_6,   poi_7,   poi_8;
  Edge_Corner_point poi_9,   poi_10;
  Edge_Corner_point poi_A,   poi_B,   poi_C,   poi_D;
  Edge_Corner_point poi_A2,  poi_B2,  poi_C2;
  Edge_Corner_point poi_oppA,   poi_oppB,   poi_oppC;

  dir_sons    oppA, oppB, oppC;
  bool later_stop, LA, LB, LC;
  // int ori;

  // leere Initialisierung
  // ori = 0;






  later_stop = false;




  // stop=true;
  //  Boundary cell : 2.2.
  // ----------------------
  // Dreibein
  // Bemerkung: Test ob Kanten im Gebiet nicht notwendig wegen corners
  for(i=0;i<8 && !stop;++i) {
    // Punkt am Eck des Dreibeins
    poi_D.edge_point = corner_poi_typ;
    poi_D.corner = (dir_sons)i;

    // Drei Ecken des Dreibeins
    poi_oppA.edge_point = poi_oppB.edge_point = 
      poi_oppC.edge_point = corner_poi_typ;

    oppA = Next_corner((dir_sons)i,opposite3D(xCoord((dir_sons)i)));
    oppB = Next_corner((dir_sons)i,opposite3D(yCoord((dir_sons)i)));
    oppC = Next_corner((dir_sons)i,opposite3D(zCoord((dir_sons)i)));

    poi_oppA.corner = oppA;
    poi_oppB.corner = oppB;
    poi_oppC.corner = oppC;

    // Kantenpunkte des Dreibeins
    poi_A.edge_point  = poi_B.edge_point  = poi_C.edge_point  = edge_poi_typ;
    poi_A2.edge_point = poi_B2.edge_point = poi_C2.edge_point = edge_poi_typ;

    poi_A.corner  = oppA; poi_A.d  = opposite3D(yCoord((dir_sons)i));
    poi_A2.corner = oppA; poi_A2.d = opposite3D(zCoord((dir_sons)i));
    poi_B.corner  = oppB; poi_B.d  = opposite3D(xCoord((dir_sons)i));
    poi_B2.corner = oppB; poi_B2.d = opposite3D(zCoord((dir_sons)i));
    poi_C.corner  = oppC; poi_C.d  = opposite3D(xCoord((dir_sons)i));
    poi_C2.corner = oppC; poi_C2.d = opposite3D(yCoord((dir_sons)i));

    if(Description->edge_point(poi_A) && Description->edge_point(poi_A2) &&
       Description->edge_point(poi_B) && Description->edge_point(poi_B2) &&
       Description->edge_point(poi_C) && Description->edge_point(poi_C2)
       && corners[i] && corners[oppA] && corners[oppB] && corners[oppC]) {
      corners[i] = corners[oppA] = corners[oppB] = corners[oppC] = false;
      //  stop = true;
      later_stop = true;
      if(((((i>>2)&1) + ((i>>1)&1) + ((i>>0)&1))&1)==0) { // Orientierung 
        //  1. Version
        poi_10 =  poi_D;

        poi_7  =  poi_oppA;
        poi_6  =  poi_A;
        poi_8  =  poi_A2;

        poi_4  =  poi_oppB;
        poi_3  =  poi_B2;
        poi_5  =  poi_B;

        poi_1  =  poi_oppC;
        poi_2  =  poi_C2;
        poi_9  =  poi_C;
      }
      else {
        //  2. Version
        poi_10 =  poi_D;

        poi_7  =  poi_oppA;
        poi_6  =  poi_A2;
        poi_8  =  poi_A;

        poi_4  =  poi_oppC;
        poi_3  =  poi_C2;
        poi_5  =  poi_C;

        poi_1  =  poi_oppB;
        poi_2  =  poi_B2;
        poi_9  =  poi_B;
      }

      Act->Tetraeder(poi_10,poi_2,poi_9,poi_1);
      Act->Tetraeder(poi_10,poi_6,poi_7,poi_8);
      Act->Tetraeder(poi_10,poi_4,poi_5,poi_3);

     LA = Description->h(poi_9.corner, poi_9.d) <
          Description->h(poi_2.corner, poi_2.d);
     LB = Description->h(poi_3.corner, poi_3.d) <
          Description->h(poi_5.corner, poi_5.d);
     LC = Description->h(poi_6.corner, poi_6.d) <
          Description->h(poi_8.corner, poi_8.d);


     if(LA && LB && LC) {
       Act->Tetraeder(poi_10,poi_8,poi_9,poi_2);
       Act->Tetraeder(poi_10,poi_2,poi_3,poi_5);
       Act->Tetraeder(poi_10,poi_5,poi_6,poi_8);
       Act->Tetraeder(poi_10,poi_8,poi_2,poi_5);
     }
     else if(!LA && !LB && !LC) {
       Act->Tetraeder(poi_10,poi_3,poi_9,poi_2);
       Act->Tetraeder(poi_10,poi_5,poi_6,poi_3);
       Act->Tetraeder(poi_10,poi_6,poi_8,poi_9);
       Act->Tetraeder(poi_10,poi_9,poi_3,poi_6);
     }
     else {
       if(LA && !LC) {
         poi_7 = poi_2;  poi_1 = poi_3;
         poi_2 = poi_5;  poi_3 = poi_6;
         poi_5 = poi_8;  poi_6 = poi_9;
         poi_8 = poi_7;  poi_9 = poi_1;
       }
       if(LB && !LA) {
         poi_4 = poi_9;  poi_1 = poi_8;
         poi_9 = poi_6;  poi_8 = poi_5;
         poi_6 = poi_3;  poi_5 = poi_2;
         poi_3 = poi_4;  poi_2 = poi_1;
       }
       Act->Tetraeder(poi_10,poi_5,poi_8,poi_3);
       Act->Tetraeder(poi_10,poi_6,poi_8,poi_5);
       if(Description->h(poi_9.corner, poi_9.d) <
          Description->h(poi_2.corner, poi_2.d)) {
         Act->Tetraeder(poi_10,poi_8,poi_9,poi_2);
         Act->Tetraeder(poi_10,poi_3,poi_8,poi_2);
       }
       else {
         Act->Tetraeder(poi_10,poi_3,poi_9,poi_2);
         Act->Tetraeder(poi_10,poi_3,poi_8,poi_9);
       }
     }
    }
}

if(later_stop) stop = true; 
}

#endif
   

---