src/expde/formulas/subdivII.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
//

// ------------------------------------------------------------
//
// subdivII.h
//
// ------------------------------------------------------------

#ifndef SUBDIVII_H_
#define SUBDIVII_H_

void Act_on_subdivision_partII(Bo_description* Description, 
                               Calculator_FE* Act,
                               bool* corners, bool& stop) {
  int i;

  Edge_Corner_point poi_0,   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,  poi_11,  poi_12,  poi_13;
  Edge_Corner_point poi_A,   poi_B,   poi_C,   poi_D;
  Edge_Corner_point poi_A2,  poi_B2,  poi_C2,  poi_D2;
  Edge_Corner_point poi_boA, poi_boB, poi_boC, poi_boD;

  dir_3D    dirA, dirB, dirC;
  dir_3D sur, suropp;
  bool later_stop, LA, LB, LC, LD;

  int num_max;
  double h_6, h_7, h_8, h_9, h_10, H, h_max;
  f_type x_face_type;
  f_type y_face_type;

  f_type A_face_type, B_face_type, C_face_type;


  later_stop = false;

  // Bemerkung:  richtig? (nicht so wichtig!)
    // Ab hier gibt es:
    // a) ein oder zwei Dreibeine
    // b) ein oder zwei Quader
    // c) ein oder zwei Typ 2.3
    // d) genau einen weitere Randtyp

  // stop = true;
  // stop = false;
  //  Boundary cells : 1.
  // ---------------------
  //  Ganze Flaeche im Gebiet
  for(i=0;i<6 && !stop;++i) {
    if(Description->Edge(Edge_surface((dir_3D)i,Edir2D)) & 
       Description->Edge(Edge_surface((dir_3D)i,Wdir2D)) &
       Description->Edge(Edge_surface((dir_3D)i,Ndir2D)) &
       Description->Edge(Edge_surface((dir_3D)i,Sdir2D))) {
      // Flaeche ist im Gebiet
      stop=true;

      poi_A.edge_point = corner_poi_typ;
      poi_B.edge_point = corner_poi_typ;
      poi_C.edge_point = corner_poi_typ;
      poi_D.edge_point = corner_poi_typ;

      poi_A2.edge_point = corner_poi_typ;
      poi_B2.edge_point = corner_poi_typ;
      poi_C2.edge_point = corner_poi_typ;
      poi_D2.edge_point = corner_poi_typ;

      poi_boA.edge_point = edge_poi_typ;
      poi_boB.edge_point = edge_poi_typ;
      poi_boC.edge_point = edge_poi_typ;
      poi_boD.edge_point = edge_poi_typ;

      sur    = (dir_3D)i;
      suropp = opposite3D((dir_3D)i);
      if(i==Wdir || i==Ndir || i==Ddir) {       
        poi_A.corner = corner_sur(sur,SWdir);
        poi_B.corner = corner_sur(sur,SEdir);
        poi_C.corner = corner_sur(sur,NWdir);
        poi_D.corner = corner_sur(sur,NEdir);
        
        poi_A2.corner = corner_sur(suropp,SWdir);
        poi_B2.corner = corner_sur(suropp,SEdir);
        poi_C2.corner = corner_sur(suropp,NWdir);
        poi_D2.corner = corner_sur(suropp,NEdir);
        
        poi_boA.corner = poi_A.corner;   poi_boA.d = opposite3D((dir_3D)i);
        LA = Description->edge_point(poi_boA);
        poi_boB.corner = poi_B.corner;   poi_boB.d = opposite3D((dir_3D)i);
        LB = Description->edge_point(poi_boB);
        poi_boC.corner = poi_C.corner;   poi_boC.d = opposite3D((dir_3D)i);
        LC = Description->edge_point(poi_boC);
        poi_boD.corner = poi_D.corner;   poi_boD.d = opposite3D((dir_3D)i);
        LD = Description->edge_point(poi_boD);
      }
      else {
        poi_A.corner = corner_sur(sur,SWdir);
        poi_B.corner = corner_sur(sur,NWdir);
        poi_C.corner = corner_sur(sur,SEdir);
        poi_D.corner = corner_sur(sur,NEdir);
        
        poi_A2.corner = corner_sur(suropp,SWdir);
        poi_B2.corner = corner_sur(suropp,NWdir);
        poi_C2.corner = corner_sur(suropp,SEdir);
        poi_D2.corner = corner_sur(suropp,NEdir);
        
        poi_boA.corner = poi_A.corner;   poi_boA.d = opposite3D((dir_3D)i);
        LA = Description->edge_point(poi_boA);
        poi_boB.corner = poi_B.corner;   poi_boB.d = opposite3D((dir_3D)i);
        LB = Description->edge_point(poi_boB);
        poi_boC.corner = poi_C.corner;   poi_boC.d = opposite3D((dir_3D)i);
        LC = Description->edge_point(poi_boC);
        poi_boD.corner = poi_D.corner;   poi_boD.d = opposite3D((dir_3D)i);
        LD = Description->edge_point(poi_boD);
      }

      //   Quader von oben:
      //   A C     
      //   B D

      //  Boundary cell : 1.1.
      // ----------------------
      // Ein Quader
      //      if(LA && LB && LC && LD) {
      if(LA && LB && LC && LD) {
        stop=false; 
        later_stop = true;

          // C D     3 2     7 6
          // A B     0 1     4 5 
        poi_3 = poi_boC;     poi_2 = poi_boD;
        poi_0 = poi_boA;     poi_1 = poi_boB;

        poi_7 = poi_C;       poi_6 = poi_D;
        poi_4 = poi_A;       poi_5 = poi_B;

        dirA = direction_from_to(poi_A.corner,poi_B.corner);
        dirB = direction_from_to(poi_A.corner,poi_C.corner);
        dirC = poi_0.d;

        dirC = opposite3D(dirC);
        A_face_type = Calc_f_type(dirA,dirC);
        B_face_type = Calc_f_type(dirB,dirC);
        C_face_type = Calc_f_type(dirA,dirB);
        if(C_face_type == an) {
          Act_on_Prism(Description,Act,A_face_type,B_face_type,
                       poi_0,poi_1,poi_3,poi_4,poi_5,poi_7); 
          Act_on_Prism(Description,Act,opposite(A_face_type),
                       opposite(B_face_type),
                       poi_2,poi_3,poi_1,poi_6,poi_7,poi_5); 
        }
        else {
          Act_on_Prism(Description,Act,A_face_type,B_face_type,
                       poi_1,poi_2,poi_0,poi_5,poi_6,poi_4); 
          Act_on_Prism(Description,Act,A_face_type,B_face_type,
                       poi_3,poi_0,poi_2,poi_7,poi_4,poi_6); 
        }
        
      }
      //  Boundary cell : 1.2.
      // ----------------------
      else if((LA && LB && LC && !LD) || (LA && !LB && LC && LD) ||
              (LA && LB && !LC && LD) || (!LA && LB && LC && LD)) {
        poi_1.edge_point = poi_2.edge_point = corner_poi_typ;
        poi_3.edge_point = poi_4.edge_point = corner_poi_typ;
        poi_5.edge_point = corner_poi_typ;
        poi_6.edge_point = edge_poi_typ;
        poi_7.edge_point = poi_8.edge_point = edge_poi_typ;
        poi_9.edge_point = poi_10.edge_point = edge_poi_typ;
        if(LA && LB && LC && !LD) {
          // B A
          // D C
          poi_2 = poi_B;          poi_3 = poi_A;
          poi_1 = poi_D;          poi_4 = poi_C;

          poi_8 = poi_boB;        poi_9  = poi_boA;
          poi_5 = poi_D2;         poi_10 = poi_boC;

          poi_7.corner = poi_5.corner;   
          poi_7.d = direction_from_to(poi_5.corner,poi_C2.corner);
          poi_6.corner = poi_5.corner;   
          poi_6.d = direction_from_to(poi_5.corner,poi_B2.corner);
        }
        if(LA && !LB && LC && LD) {
          // A C
          // B D
          poi_2 = poi_A;          poi_3 = poi_C;
          poi_1 = poi_B;          poi_4 = poi_D;

          poi_8 = poi_boA;        poi_9  = poi_boC;
          poi_5 = poi_B2;         poi_10 = poi_boD;

          poi_7.corner = poi_5.corner;   
          poi_7.d = direction_from_to(poi_5.corner,poi_D2.corner);
          poi_6.corner = poi_5.corner;   
          poi_6.d = direction_from_to(poi_5.corner,poi_A2.corner);
        }            
        if(LA && LB && !LC && LD) {
          // D B
          // C A
          poi_2 = poi_D;          poi_3 = poi_B;
          poi_1 = poi_C;          poi_4 = poi_A;

          poi_8 = poi_boD;        poi_9  = poi_boB;
          poi_5 = poi_C2;         poi_10 = poi_boA;

          poi_7.corner = poi_5.corner;   
          poi_7.d = direction_from_to(poi_5.corner,poi_A2.corner);
          poi_6.corner = poi_5.corner;   
          poi_6.d = direction_from_to(poi_5.corner,poi_D2.corner);
        }
        if(!LA && LB && LC && LD) {
          // C D
          // A B
          poi_2 = poi_C;          poi_3 = poi_D;
          poi_1 = poi_A;          poi_4 = poi_B;

          poi_8 = poi_boC;        poi_9  = poi_boD;
          poi_5 = poi_A2;         poi_10 = poi_boB;

          poi_7.corner = poi_5.corner;   
          poi_7.d = direction_from_to(poi_5.corner,poi_B2.corner);
          poi_6.corner = poi_5.corner;   
          poi_6.d = direction_from_to(poi_5.corner,poi_C2.corner);
        }
        
        
        Act->Tetraeder(poi_1,poi_5,poi_6,poi_7);

        dirA = direction_from_to(poi_4.corner,poi_1.corner);
        dirB = direction_from_to(poi_4.corner,poi_3.corner);
        if(Calc_f_type(dirA,dirB)==an) { // Bei Punkt 4 an
          Act->Tetraeder(poi_1,poi_6,poi_9,poi_7);
          Act->Tetraeder(poi_1,poi_7,poi_9,poi_10);
          Act->Tetraeder(poi_1,poi_6,poi_8,poi_9);
          // Anstelle von Pyr 1, 10, 9, 3, 4
          dirA = poi_10.d;
          dirB = direction_from_to(poi_4.corner,poi_3.corner);
          if(Calc_f_type(dirA,dirB)==an) {
            Act->Tetraeder(poi_1,poi_10,poi_3,poi_4);
            Act->Tetraeder(poi_1,poi_10,poi_9,poi_3);
          }
          else {
            Act->Tetraeder(poi_1,poi_10,poi_9,poi_4);
            Act->Tetraeder(poi_1,poi_9,poi_3,poi_4);
          }
          // Anstelle von Pyr 1, 9, 8, 2, 3
          dirA = poi_9.d;
          dirB = direction_from_to(poi_3.corner,poi_2.corner);
          if(Calc_f_type(dirA,dirB)==an) {
            Act->Tetraeder(poi_1,poi_9,poi_2,poi_3);
            Act->Tetraeder(poi_1,poi_8,poi_2,poi_9);
          }
          else {
            Act->Tetraeder(poi_1,poi_8,poi_3,poi_9);
            Act->Tetraeder(poi_1,poi_8,poi_2,poi_3);
          }
        }
        else { // Bei Punkt 4 spitz
          // Anstelle von:  Pyr 1, 6, 8, 10, 7
          if(parallel_version)
            Tell_about_round_of_problem(Description->h(poi_7.corner, poi_7.d),
                                        Description->h(poi_6.corner, poi_6.d));
          if(Description->h(poi_7.corner, poi_7.d) <
             Description->h(poi_6.corner, poi_6.d)) {
            Act->Tetraeder(poi_1,poi_7,poi_6,poi_10);
            Act->Tetraeder(poi_1,poi_6,poi_8,poi_10);
          }
          else {
            Act->Tetraeder(poi_1,poi_7,poi_6,poi_8);
            Act->Tetraeder(poi_1,poi_7,poi_8,poi_10);
          }
          dirA = opposite3D(poi_9.d);
          dirB = direction_from_to(poi_3.corner,poi_4.corner);
          dirC = direction_from_to(poi_3.corner,poi_2.corner);

          x_face_type = Calc_f_type(dirA,dirC);
          y_face_type = Calc_f_type(dirA,dirB);
          
          A_face_type = Act_on_Prism(Description,Act,x_face_type,y_face_type,
                                     poi_9,poi_8,poi_10,poi_3,poi_2,poi_4); 

          if(A_face_type==an) {
            Act->Tetraeder(poi_1,poi_8,poi_2,poi_10);
            Act->Tetraeder(poi_1,poi_10,poi_2,poi_4);
          }
          else {
            Act->Tetraeder(poi_1,poi_8,poi_4,poi_10);
            Act->Tetraeder(poi_1,poi_8,poi_2,poi_4);
          }
        }
        
      }
      //  Boundary cell : 1.3.
      // ----------------------
      else if((LA && LB && !LC && !LD) || (LA && !LB && LC && !LD) ||
              (!LA && !LB && LC && LD) || (!LA && LB && !LC && LD)) {
        poi_1.edge_point = poi_2.edge_point = corner_poi_typ;
        poi_3.edge_point = poi_4.edge_point = corner_poi_typ;
        poi_5.edge_point = poi_6.edge_point = corner_poi_typ;
        poi_7.edge_point = poi_8.edge_point = edge_poi_typ;
        poi_9.edge_point = poi_10.edge_point = edge_poi_typ;
        if(!LA && !LB && LC && LD) {
          // A C
          // B D
          poi_2 = poi_A;          poi_3 = poi_C;
          poi_1 = poi_B;          poi_4 = poi_D;

          poi_6 = poi_A2;         poi_10 = poi_boC;
          poi_5 = poi_B2;         poi_9  = poi_boD;

          poi_7.corner = poi_6.corner;
          poi_7.d = direction_from_to(poi_6.corner,poi_C2.corner);
          poi_8.corner = poi_5.corner;
          poi_8.d = direction_from_to(poi_5.corner,poi_D2.corner);
        }
        if(LA && LB && !LC && !LD) {
          // D B
          // C A
          poi_2 = poi_D;          poi_3 = poi_B;
          poi_1 = poi_C;          poi_4 = poi_A;

          poi_6 = poi_D2;         poi_10 = poi_boB;
          poi_5 = poi_C2;         poi_9  = poi_boA;

          poi_7.corner = poi_6.corner;
          poi_7.d = direction_from_to(poi_6.corner,poi_B2.corner);
          poi_8.corner = poi_5.corner;
          poi_8.d = direction_from_to(poi_5.corner,poi_A2.corner);
        }
        if(LA && !LB && LC && !LD) {
          // B A
          // D C
          poi_2 = poi_B;          poi_3 = poi_A;
          poi_1 = poi_D;          poi_4 = poi_C;

          poi_6 = poi_B2;         poi_10 = poi_boA;
          poi_5 = poi_D2;         poi_9  = poi_boC;

          poi_7.corner = poi_6.corner;
          poi_7.d = direction_from_to(poi_6.corner,poi_A2.corner);
          poi_8.corner = poi_5.corner;
          poi_8.d = direction_from_to(poi_5.corner,poi_C2.corner);
        }
        if(!LA && LB && !LC && LD) {
          // C D
          // A B
          poi_2 = poi_C;          poi_3 = poi_D;
          poi_1 = poi_A;          poi_4 = poi_B;

          poi_6 = poi_C2;         poi_10 = poi_boD;
          poi_5 = poi_A2;         poi_9  = poi_boB;

          poi_7.corner = poi_6.corner;
          poi_7.d = direction_from_to(poi_6.corner,poi_D2.corner);
          poi_8.corner = poi_5.corner;
          poi_8.d = direction_from_to(poi_5.corner,poi_B2.corner);
        }

        h_9  = Description->h(poi_9.corner,  poi_9.d);
        h_7  = Description->h(poi_7.corner,  poi_7.d);
        h_8  = Description->h(poi_8.corner,  poi_8.d);
        h_10 = Description->h(poi_10.corner, poi_10.d);

        H = Description->Meshsize();

        num_max = 9;
        h_max   = h_9;
        if(parallel_version) Tell_about_round_of_problem(h_7 , h_max);
        if(h_7 > h_max) {
          num_max = 7;
          h_max   = h_7;
        }
        if(parallel_version) Tell_about_round_of_problem(h_8 , h_max);
        if(h_8 > h_max) {
          num_max = 8;
          h_max   = h_8;
        }
        if(parallel_version) Tell_about_round_of_problem(h_10 , h_max);
        if(h_10 > h_max) {
          num_max = 10;
          h_max   = h_10;
        }

        if(parallel_version) Tell_about_round_of_problem(h_max,0.5*H);
        if(h_max > 0.5*H) {
          // version with additional point
          poi_11.edge_point = cell_poi_typ;
          dirA = direction_from_to(poi_1.corner,poi_2.corner);

          if(dirA==Wdir || dirA==Edir)
            Act->Add_Bo_freedom(D3vector( 0.0,
                                          -0.1+0.2*y1DCoord(poi_1.corner),
                                          -0.1+0.2*z1DCoord(poi_1.corner)));
          if(dirA==Ndir || dirA==Sdir)
            Act->Add_Bo_freedom(D3vector(-0.1+0.2*x1DCoord(poi_1.corner),
                                         0.0,
                                         -0.1+0.2*z1DCoord(poi_1.corner)));
          if(dirA==Tdir || dirA==Ddir)
            Act->Add_Bo_freedom(D3vector(-0.1+0.2*x1DCoord(poi_1.corner),
                                         -0.1+0.2*y1DCoord(poi_1.corner),
                                         0.0));
          
          // East
          // Anstelle von Pyr 11, 10, 9, 3, 4
          dirA = poi_9.d;
          dirB = direction_from_to(poi_4.corner,poi_3.corner);
          if(Calc_f_type(dirA,dirB)==an) {
            Act->Tetraeder(poi_11,poi_9,poi_3,poi_4);
            Act->Tetraeder(poi_11,poi_9,poi_10,poi_3);
          }
          else {
            Act->Tetraeder(poi_11,poi_9,poi_10,poi_4);
            Act->Tetraeder(poi_11,poi_10,poi_3,poi_4);
          }
          // West
          // Anstelle von Pyr 11, 5, 1, 2, 6
          dirA = poi_9.d;
          dirB = direction_from_to(poi_1.corner,poi_2.corner);
          if(Calc_f_type(dirA,dirB)==an) {
            Act->Tetraeder(poi_11,poi_5,poi_1,poi_2);
            Act->Tetraeder(poi_11,poi_6,poi_5,poi_2);
          }
          else {
            Act->Tetraeder(poi_11,poi_6,poi_1,poi_2);
            Act->Tetraeder(poi_11,poi_6,poi_5,poi_1);
          }
          // Top
          // Anstelle von Pyr 11, 8, 5, 6, 7
          dirA = poi_8.d;
          dirB = direction_from_to(poi_1.corner,poi_2.corner);
          if(Calc_f_type(dirA,dirB)==an) {
            Act->Tetraeder(poi_11,poi_8,poi_5,poi_6);
            Act->Tetraeder(poi_11,poi_8,poi_6,poi_7);
          }
          else {
            Act->Tetraeder(poi_11,poi_8,poi_5,poi_7);
            Act->Tetraeder(poi_11,poi_5,poi_6,poi_7);
          }
          // Bottom
          // Anstelle von Pyr 11, 2, 1, 4, 3
          dirA = poi_8.d;
          dirB = direction_from_to(poi_1.corner,poi_2.corner);
          if(Calc_f_type(dirA,dirB)==an) {
            Act->Tetraeder(poi_11,poi_2,poi_1,poi_4);
            Act->Tetraeder(poi_11,poi_2,poi_4,poi_3);
          }
          else {
            Act->Tetraeder(poi_11,poi_3,poi_1,poi_4);
            Act->Tetraeder(poi_11,poi_3,poi_2,poi_1);
          }
          // South
          Act->Tetraeder(poi_11,poi_8,poi_1,poi_5);
          Act->Tetraeder(poi_11,poi_8,poi_9,poi_1);
          Act->Tetraeder(poi_11,poi_9,poi_4,poi_1);
          
          // North
          Act->Tetraeder(poi_11,poi_7,poi_6,poi_2);
          Act->Tetraeder(poi_11,poi_7,poi_2,poi_10);
          Act->Tetraeder(poi_11,poi_10,poi_2,poi_3);
          // Oberflaeche
          // Anstelle von Pyr 11, 8, 7, 10, 9
          if(num_max==9 || num_max==7) {
            Act->Tetraeder(poi_11,poi_9,poi_8,poi_7);
            Act->Tetraeder(poi_11,poi_9,poi_7,poi_10);
          }
          else {
            Act->Tetraeder(poi_11,poi_9,poi_8,poi_10);
            Act->Tetraeder(poi_11,poi_8,poi_7,poi_10);
          }
        }
        else {
          // version with prisms
          // Anstelle: Act->Prism(poi_5,poi_1,poi_8,poi_6,poi_2,poi_7);
          dirA = direction_from_to(poi_5.corner,poi_1.corner);
          dirB = direction_from_to(poi_5.corner,poi_6.corner);
          x_face_type=Calc_f_type(dirA,dirB);
          dirA = poi_8.d;
          y_face_type=Calc_f_type(dirA,dirB);
          A_face_type = Act_on_Prism_long(Description,Act,
                                          x_face_type,y_face_type,
                                          poi_5,poi_1,poi_8,poi_6,poi_2,poi_7);
          
          // Anstelle: Act->Prism(poi_4,poi_9,poi_1,poi_3,poi_10,poi_2);
          dirA = poi_9.d;
          dirB = direction_from_to(poi_4.corner,poi_3.corner);
          x_face_type=Calc_f_type(dirA,dirB);
          dirA = direction_from_to(poi_4.corner,poi_1.corner);
          y_face_type=Calc_f_type(dirA,dirB);
          B_face_type = Act_on_Prism_short(Description,Act,
                                           x_face_type,y_face_type,
                                           poi_4,poi_9,poi_1,
                                           poi_3,poi_10,poi_2);
          // Anstelle: Act->Prism(poi_8,poi_1,poi_9,poi_7,poi_2,poi_10); 
          Act_on_Prism(Description,Act,opposite_face_type(B_face_type),
                       A_face_type,
                       poi_1,poi_9,poi_8,poi_2,poi_10,poi_7);
        }
      }
     //  Boundary cell : 1.4.
     // ----------------------
      else if((LA && !LB && !LC && LD) || (!LA && LB && LC && !LD)) {
        poi_1.edge_point = poi_2.edge_point = corner_poi_typ;
        poi_3.edge_point = poi_4.edge_point = corner_poi_typ;
        poi_5.edge_point = poi_8.edge_point = corner_poi_typ;
        poi_6.edge_point = poi_7.edge_point = edge_poi_typ;
        poi_9.edge_point = poi_10.edge_point = edge_poi_typ;
        poi_11.edge_point = poi_12.edge_point = edge_poi_typ;
        if(!LA && LB && LC && !LD) {
          // A C
          // B D
          poi_2 = poi_A;          poi_3 = poi_C;
          poi_1 = poi_B;          poi_4 = poi_D;

          poi_8  = poi_A2;        poi_11 = poi_boC;  
          poi_12 = poi_boB;       poi_5  = poi_D2;
        }
        if(LA && !LB && !LC && LD) {
          // B A
          // D C
          poi_2 = poi_B;          poi_3 = poi_A;
          poi_1 = poi_D;          poi_4 = poi_C;

          poi_8  = poi_B2;        poi_11 = poi_boA;  
          poi_12 = poi_boD;       poi_5  = poi_C2;
        }
        poi_9.corner  = poi_8.corner;
        poi_10.corner = poi_8.corner;
        poi_9.d  = direction_from_to(poi_2.corner,poi_3.corner);
        poi_10.d = direction_from_to(poi_2.corner,poi_1.corner);
        poi_6.corner = poi_5.corner;
        poi_7.corner = poi_5.corner;
        poi_6.d = direction_from_to(poi_4.corner,poi_3.corner);
        poi_7.d = direction_from_to(poi_4.corner,poi_1.corner);
        
        poi_13.edge_point = cell_poi_typ;
        Act->Add_Bo_freedom(D3vector(0.0,0.0,0.0));

        // Ecke 5 und 8
        Act->Tetraeder(poi_5,poi_6,poi_7,poi_4);
        Act->Tetraeder(poi_8,poi_10,poi_9,poi_2);

        // Verbindung zu poi_13
        Act->Tetraeder(poi_13,poi_10,poi_2,poi_9);
        Act->Tetraeder(poi_13,poi_7,poi_6,poi_4);

        // Dreieck in der Mitte einer c)-Zelle zu poi_13
        Act->Tetraeder(poi_13,poi_7,poi_4,poi_12);
        Act->Tetraeder(poi_13,poi_6,poi_11,poi_4);
        Act->Tetraeder(poi_13,poi_11,poi_9,poi_2);
        Act->Tetraeder(poi_13,poi_10,poi_12,poi_2);

        // Restliche Dreiecke  einer c)-Zelle zu poi_13
        Act->Tetraeder(poi_13,poi_12,poi_4,poi_1);
        Act->Tetraeder(poi_13,poi_12,poi_1,poi_2);
        Act->Tetraeder(poi_13,poi_11,poi_2,poi_3);
        Act->Tetraeder(poi_13,poi_11,poi_3,poi_4);
                
        // Anstelle von Pyr 13, 2, 1, 4, 3
        dirA = poi_6.d;
        dirB = poi_7.d;
        if(Calc_f_type(dirA,dirB)==an) {
          Act->Tetraeder(poi_13,poi_3,poi_1,poi_4);
          Act->Tetraeder(poi_13,poi_3,poi_2,poi_1);
        }
        else {
          Act->Tetraeder(poi_13,poi_2,poi_4,poi_3);
          Act->Tetraeder(poi_13,poi_2,poi_1,poi_4);
        }       

        // Oberflaeche:
        Act->Tetraeder(poi_13,poi_7,poi_12,poi_10);
        Act->Tetraeder(poi_13,poi_6,poi_9,poi_11);
        
        h_6 = Description->h(poi_6.corner, poi_6.d);
        h_7 = Description->h(poi_7.corner, poi_7.d);
        h_9 = Description->h(poi_9.corner, poi_9.d);
        h_10 = Description->h(poi_10.corner, poi_10.d);

        double H; 
        double angle_6;
        double angle_7;
        double angle_9;
        double angle_10;
        
        H = Description->Meshsize();

        angle_6  = calc_special_angle(h_7,h_6,H-h_9,H-h_6);
        angle_7  = calc_special_angle(H-h_10,H-h_7,h_6,h_7);
        angle_10 = calc_special_angle(h_9,h_10,H-h_7,H-h_10);
        angle_9  = calc_special_angle(H-h_6,H-h_9,h_10,h_9);

        // Anstelle von Pyr 13, 7, 10, 9, 6
        //      if(num_min==6 || num_min==10) {
        if(parallel_version) 
          Tell_about_round_of_problem(MAX(angle_7,angle_9),
                                      MAX(angle_10,angle_6));
        if(MAX(angle_7,angle_9) > MAX(angle_10,angle_6)) {
          Act->Tetraeder(poi_13,poi_7,poi_10,poi_9);
          Act->Tetraeder(poi_13,poi_7,poi_9,poi_6);
        }
        else {
          Act->Tetraeder(poi_13,poi_7,poi_10,poi_6);
          Act->Tetraeder(poi_13,poi_6,poi_10,poi_9);
        }
      }
      //  Boundary cell : 1.5.
      // ----------------------
      else if((!LA && !LB && !LC && LD) || (!LA && !LB && LC && !LD) ||
              (!LA && LB && !LC && !LD) || (LA && !LB && !LC && !LD)) {
        poi_1.edge_point = poi_2.edge_point = corner_poi_typ;
        poi_3.edge_point = poi_4.edge_point = corner_poi_typ;
        poi_5.edge_point = poi_6.edge_point = corner_poi_typ;
        poi_7.edge_point = corner_poi_typ;
        poi_8.edge_point = edge_poi_typ;
        poi_9.edge_point = poi_10.edge_point = edge_poi_typ;

        if(!LA && !LB && !LC && LD) {
          // A C
          // B D
          poi_2 = poi_A;          poi_3 = poi_C;
          poi_1 = poi_B;          poi_4 = poi_D;

          poi_6 = poi_A2;         poi_7 = poi_C2;
          poi_5 = poi_B2;

          poi_8.corner = poi_7.corner;
          poi_8.d = direction_from_to(poi_7.corner,poi_D2.corner);
          poi_9.corner = poi_5.corner;
          poi_9.d = direction_from_to(poi_5.corner,poi_D2.corner);
          poi_10.corner = poi_4.corner;
          poi_10.d = direction_from_to(poi_4.corner,poi_D2.corner);
        }
        if(LA && !LB && !LC && !LD) {
          // D B
          // C A
          poi_2 = poi_D;          poi_3 = poi_B;
          poi_1 = poi_C;          poi_4 = poi_A;

          poi_6 = poi_D2;        poi_7 = poi_B2;
          poi_5 = poi_C2;

          poi_8.corner = poi_7.corner;
          poi_8.d = direction_from_to(poi_7.corner,poi_A2.corner);
          poi_9.corner = poi_5.corner;
          poi_9.d = direction_from_to(poi_5.corner,poi_A2.corner);
          poi_10.corner = poi_4.corner;
          poi_10.d = direction_from_to(poi_4.corner,poi_A2.corner);
        }
        if(!LA && !LB && LC && !LD) {
          // B A
          // D C
          poi_2 = poi_B;          poi_3 = poi_A;
          poi_1 = poi_D;          poi_4 = poi_C;

          poi_6 = poi_B2;         poi_7 = poi_A2;         
          poi_5 = poi_D2;         

          poi_8.corner = poi_7.corner;
          poi_8.d = direction_from_to(poi_7.corner,poi_C2.corner);
          poi_9.corner = poi_5.corner;
          poi_9.d = direction_from_to(poi_5.corner,poi_C2.corner);
          poi_10.corner = poi_4.corner;
          poi_10.d = direction_from_to(poi_4.corner,poi_C2.corner);
        }
        if(!LA && LB && !LC && !LD) {
          // C D
          // A B
          poi_2 = poi_C;          poi_3 = poi_D;
          poi_1 = poi_A;          poi_4 = poi_B;

          poi_6 = poi_C2;         poi_7 = poi_D2;
          poi_5 = poi_A2;

          poi_8.corner = poi_7.corner;
          poi_8.d = direction_from_to(poi_7.corner,poi_B2.corner);
          poi_9.corner = poi_5.corner;
          poi_9.d = direction_from_to(poi_5.corner,poi_B2.corner);
          poi_10.corner = poi_4.corner;
          poi_10.d = direction_from_to(poi_4.corner,poi_B2.corner);
        }

        poi_11.edge_point = cell_poi_typ;
        /*
        Act->Add_Bo_freedom(D3vector(-0.25+0.5*x1DCoord(poi_2.corner),
                                     -0.25+0.5*y1DCoord(poi_2.corner),
                                     -0.25+0.5*z1DCoord(poi_2.corner)));
        */
        Act->Add_Bo_freedom(D3vector(0.0,0.0,0.0));

        // Anstelle: Act->Pyramid(poi_11,poi_5,poi_1,poi_2,poi_6);
        dirA = direction_from_to(poi_5.corner,poi_6.corner);
        dirB = direction_from_to(poi_5.corner,poi_1.corner);
        if(Calc_f_type(dirA,dirB)==an) {
          Act->Tetraeder(poi_11,poi_5,poi_1,poi_6);
          Act->Tetraeder(poi_11,poi_6,poi_1,poi_2);
        }
        else {
          Act->Tetraeder(poi_11,poi_6,poi_5,poi_2);
          Act->Tetraeder(poi_11,poi_5,poi_1,poi_2);
        }
        // Anstelle: Act->Pyramid(poi_11,poi_4,poi_3,poi_2,poi_1);
        dirA = direction_from_to(poi_4.corner,poi_1.corner);
        dirB = direction_from_to(poi_4.corner,poi_3.corner);
        if(Calc_f_type(dirA,dirB)==an) {
          Act->Tetraeder(poi_11,poi_3,poi_1,poi_4);
          Act->Tetraeder(poi_11,poi_2,poi_1,poi_3);
        }
        else {
          Act->Tetraeder(poi_11,poi_2,poi_1,poi_4);
          Act->Tetraeder(poi_11,poi_2,poi_4,poi_3);
        }
        // Anstelle: Act->Pyramid(poi_11,poi_6,poi_2,poi_3,poi_7);
        dirA = direction_from_to(poi_7.corner,poi_6.corner);
        dirB = direction_from_to(poi_7.corner,poi_3.corner);
        if(Calc_f_type(dirA,dirB)==an) {
          Act->Tetraeder(poi_11,poi_6,poi_3,poi_7);
          Act->Tetraeder(poi_11,poi_6,poi_2,poi_3);
        }
        else {
          Act->Tetraeder(poi_11,poi_6,poi_2,poi_7);
          Act->Tetraeder(poi_11,poi_7,poi_2,poi_3);
        }
        
        Act->Tetraeder(poi_11,poi_1,poi_5,poi_9);
        Act->Tetraeder(poi_11,poi_1,poi_10,poi_4);
        Act->Tetraeder(poi_11,poi_1,poi_9,poi_10);
        
        Act->Tetraeder(poi_11,poi_6,poi_9,poi_5);
        Act->Tetraeder(poi_11,poi_6,poi_8,poi_9);
        Act->Tetraeder(poi_11,poi_6,poi_7,poi_8);

        Act->Tetraeder(poi_11,poi_3,poi_8,poi_7);
        Act->Tetraeder(poi_11,poi_3,poi_10,poi_8);
        Act->Tetraeder(poi_11,poi_3,poi_4,poi_10);

        Act->Tetraeder(poi_11,poi_9,poi_8,poi_10);

        }
    }
}

if(later_stop) stop = true; 
}


#endif
   

---