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

// ------------------------------------------------------------
//
// boundy.h
//
// ------------------------------------------------------------

#ifndef BOUNDY_H_
#define BOUNDY_H_
        
class Bo_description;

// 1. Tetraeder
// 2. Beschreibung einer Randzelle;
// 3. Speicher fuer Randzellfreiheitsgrad
// 4. ( lokale Randzellen-Matrizen ) existiert nicht mehr
// 5. Beschreibung einer Randzelle;
// 6. inline Implementierung einiger Memberfunktionenen


// 1. Tetraeder


class Tetraeder_storage {
 public:
  // Constructor for interior tetrahedrals
  inline Tetraeder_storage(Bo_description* desc,
                    int P0, int P1, int P2, int P3,
                    Tetraeder_storage* Next);
  // Constructor for tetrahedrals with boundary face
  Tetraeder_storage(Tetraeder_storage* Next) {
    next = Next;
    var = NULL;
  }
  
  // Alles zur Listenverwaltung; Next
  void Put_Next(Tetraeder_storage* Next) {
    next = Next;
  }
  Tetraeder_storage* Next() { return next; }

  // Alles zur Determinante
  inline void Put_det(double deter) { det = deter; }
  inline double Det() { return det; }  

  // Nummer der Eckpunkte
  int N0() { return num0; }
  int N1() { return num1; }
  int N2() { return num2; }
  int N3() { return num3; }

  // Speicher fuer lokale Steifigkeitsmatrizen
  double Xm0, Ym0, Zm0;
  double Xm1, Ym1, Zm1;
  double Xm2, Ym2, Zm2;
  double Xm3, Ym3, Zm3;

  // Fuer Tetraeder, die eine Randflaeche haben
  virtual bool Boundary_tet() { return false; }
  virtual double Give_det_surface() { return 0.0; }
  virtual void    Put_det_surface(double value) { };

  // test interior edge and face angles
  inline bool Check_angles(Bo_description* desc);

  // for cell variable
  double *Give_variable() { return var; };
  void    Set_var(double *v) { var = v; };

  virtual ~Tetraeder_storage() {
    delete next;
    if(var!=NULL) delete var;
  }

 protected:
  int num0;  int num1;  int num2;  int num3;
  double det;
  Tetraeder_storage* next;

  // for cell variable
  double *var;
};

// Tetraeder, die eine Randflaeche haben
class Boundary_tetraeder_storage : public Tetraeder_storage {
 public:
  Boundary_tetraeder_storage(Bo_description* desc, 
                             int P0, int P1, int P2, int P3,
                             Tetraeder_storage* Next) : 
    Tetraeder_storage(Next) {
    num0 = P0; num1 = P1; num2 = P2; num3 = P3;
  };

  bool Boundary_tet() { return true; }
  double Give_det_surface()             { return det_surface; }
  void    Put_det_surface(double value) { det_surface = value; }

  ~Boundary_tetraeder_storage() {
    delete next;
    if(var!=NULL) delete var;
  }
 private:
  double det_surface;
};


// 2. Beschreibung einer Randzelle;
// dies wird jedoch nicht abgespeichert,
// sondern nur temporaer erzeugt

class BoCeData;

class Bo_description {
 public:
  Bo_description(); // Konstruktor, setzt alles false bzw. 0.0
  void Put_zero();  // setzt alles false bzw. 0.0, wie Konstruktor

  // Abfragen der Informationen ueber die Zelle
  bool corner_point(dir_sons e);          // Ist Eckpunkt e im Gebiet?
  bool edge_point(dir_sons e,dir_3D d);   // Ist Kantenpunkt mit Ecke e und
                                          // Richtung d im Gebiet?
                                          // (Hier ist gemeint EW, NS, TD) 
  bool edge_point(Edge_Corner_point poi);
  double h(dir_sons e,dir_3D d);          // Falls edge_point(e,d)==true :
                                          // Abstand von Ecke e in Richtung d
                                          // zum Randpunkt
  bool Edge(Edges_cell ed);               // Ist ed innere Kante?
  double Meshsize();                      // Gib Maschenweite der Zelle
  int  edge_info();                       // Gib die ganze Information 
                                          // ueber Kanten

  // Koordinaten der Punkte:  bezueglich [0,h^3]
  D3vector coord(dir_sons e);             // Koordinate eines Eckpunktes
  D3vector coord(dir_sons e,dir_3D);      // Koordinate eines Kantenpunktes
  D3vector coord(Edge_Corner_point poi);

  // Setzen der Informationen ueber die Zelle
  void Put_corner_point(dir_sons e);                 // Setze einen Eckpunkt
  void Put_edge_point(dir_sons e,dir_3D d,double h); // Setze einen Kantenpunkt
  void Put_Meshsize(double H);                       // Setze Maschenweite 
  void Put_edge(Edges_cell ed);                      // Setze innere Kante

  // Funktionen zum Testen:  muss spaeter wohl weg:
  void Test_Initialisierung();

  // Nummern der Eck und Kantenpunkte, die man spaeter setzen kann
  void Set_num_corner(dir_sons e, int num);
  void Set_num_edge_poi(dir_sons e, dir_3D d, int num);
  void Set_num_cellpoi(int num);
  int Num_corner(dir_sons e);
  int Num_edge_poi(dir_sons e,dir_3D d);
  int Num(Edge_Corner_point poi);
  
  // For calc of max interior angle
  void Put_bocedata(BoCeData* boce) { bocedat=boce; };
  BoCeData* Give_bocedata() { return bocedat; };

  // Only for some testing:
  double * Give_dis() { return dis; }

 private:
  bool corner[8];    // Info ueber Eckpunkte                
  bool edge_poi[24]; // Info ueber  Kantenpunkte, Implementierung: e*3+(d>>1) 
  double dis[24];    // Abstand bei Kantenpunkte, Implementierung: e*3+(d>>1) 
  int  edge;         // Info ueber Kanten
  double H;        // Maschenweite

  // Hilfsspeicher, in diesen stehen spaeter die Nummern der Punkte
  int num_corner[8];
  int num_edge_poi[24];
  int num_cellpoi; 

  BoCeData* bocedat; 
};


// 3. Speicher fuer Randzellfreiheitsgrad
// Die lokale Koordinate local_coord
// bezieht sich auf den Mittelpunkt der Randzelle
// Randzelle [-H/2,H/2]^3 Zelle
class Bocellpoint {
 public:
  Bocellpoint(D3vector l_coord) { 
    local_coord=l_coord;
    label_level=NULL;
  };
  D3vector Give_local_coord() { return local_coord; };
  void Set_number_avs(int avs_num) { number_avs = avs_num; };
  int  Number_avs() { return number_avs; };

  int *label_level; // for boundary labels
 private:
  D3vector local_coord;
  int number_avs;
};




// 5. Beschreibung einer Randzelle;
// und Informationen ueber FE-Tetraeder
// und Information ueber Randzellfreiheitsgrad
// dies wird abgespeicher


class BoCeData {
 public:
  BoCeData() { bocellpoint=NULL; 
               Tets_boundary=NULL;
               Tets=NULL;
               vars = NULL;                //neu
               info = NULL;                //neu
               hhh  = NULL;
  };


  // Additionals:
  int  Give_total_number_points() const;
  bool Is_point_clipped_point(int num) const;

  void Set_h(double *h, Bo_description& desc);


  // Funktionen zu FE-Tetraeder:
  void Calc_FE_elements(Bo_description& desc);   // Berechne FE-elemente
  int Give_number_points() const; // Gib Anzahl N der Punkte (Eck- und Kanten-)
                                  // Nummerierung: 0,...,N-1
                                  // N ist Nummer des  Randzellfreiheitsgrad
                                  // falls vorhanden
  void Add_Tetraeder(Bo_description* desc, int P0, int P1, int P2, int P3);
  void Add_Tetraeder_bo(Bo_description* desc, int P0, int P1, int P2, int P3);

  // Koordinaten der Punkte bezueglich [0,H]^3 Zelle
  D3vector coord(int num, double* h);       
  double *hhh;            // neu fuer Andreas

  // Funktionen die Informationen ueber Punktes mit Nummer num angeben:
  inline dir_sons corner(int num) const;   // Ecke des Punktes mit Nummer num
  inline Cell_type_points edge_point(int num) const; 
                        // Ist Punkt mit Nummer num
                        // Kantenpunkt, oder Eckpunkt oder Zellfreiheitsgrad?
  inline dir_3D edge_dir(int num); // Falls Punkt mit Nummer num Kantenpunkt:
                                   // In welcher Richtung liegt die Kante?
  bool Edge(Edges_cell ed);        // Ist ed innere Kante?

  inline double Meshsize();        // Gib Maschenweite der Zelle

  // Funktionen die die Implementierung der lokalen Bilinearform vereinfachen:
  bool boundary_EW(int num);       // Gibt es in EW von Ecke num aus Randpunkt?
  bool boundary_NS(int num);       // Gibt es in NS von Ecke num aus Randpunkt?
  bool boundary_TD(int num);       // Gibt es in TD von Ecke num aus Randpunkt?

  // Anfang der Liste der Tetraeder:
  Tetraeder_storage* Give_tets() const { return Tets; } 
  Boundary_tetraeder_storage* Give_boundary_tets() const
    { return Tets_boundary; } 

  // Randzellfreiheitsgrad
  bool Exists_bocellpoint() const { return bocellpoint!=NULL; }
  int* Give_label_level() const { return bocellpoint->label_level; }
  void Set_label_level(int* ll) { bocellpoint->label_level = ll; }
  void Set_Dirichlet(int num)   { bocellpoint->label_level[num-1]=true; }
  bool Is_Dirichlet(int num) const { return bocellpoint->label_level[num-1]; }

  void Set_number_avs(int avs_num);
  int Number_avs();
  void Add_Bo_freedom(D3vector local_coord); 
  D3vector Local_coord_bocellpoint(); // bezueglich [-H/2,H/2]^3 Zelle
  inline D3vector Coord_bocellpoint(); // bezueglich [0,H]^3 Zelle
  inline D3vector Coord_bocellpoint_normal() const; // bezueglich [0,1]^3 Zelle


  // Zum setzen der Zeiger auf die Variablen
  void Set_variable_pointer(int num, double* poi);
  // Zeiger auf die Variablen der Punkte
  double** vars; // zerst Nummer des Punktes, dann Nummer der Variable

  // Zur Fehlersuche: Informationen ausdrucken:
  void Print();
  void Analyse_angles(Bo_description& desc);   // Analysiere Winkel

  ~BoCeData() {
    if(bocellpoint!=NULL)
      delete bocellpoint;
    if(Tets!=NULL)
      delete Tets;
    if(vars!=NULL)                //neu
      delete vars;                //neu
    if(info!=NULL)                //neu
      delete info;                //neu
  };

private:
  // Setzt Informationen ueber Punkte und H entsprechend desc
  void Initialize(Bo_description& desc);
  // Diese Funktion verwendet folgende 
  // Funktionen zur Eingabe:
  void Put_corner_point(int num, dir_sons c);       // Setze Eckpunkt
  void Put_edge_point(int num,dir_sons c,dir_3D d); // Setze Kantenpunkt
  void Allocate(int n_points);              // Stelle Speicher fuer n_points
                                            // bereit
  void Put_edge_info(int ed);
  void Put_Meshsize(double h);                      // Setze Maschenweite 

  int number_points;     // Anzahl der Ecken
  int*    info;           // Zeiger auf den Vektor der Informationen 
                          // ueber dir Ecken
  int  edge;              // Info ueber Kanten
  double H;               // Maschenweite



  // Anfang der Liste der Tetraeder:
  Tetraeder_storage* Tets; 
  Boundary_tetraeder_storage* Tets_boundary; 

  // Zeiger auf eventuell vorhandenen Randzellfreiheitsgrad
  Bocellpoint* bocellpoint;

  // Erklaerung zu info (Zur Anwendung unwichtig!!)
  // -----------------------------------------------------------------
  // Bit:    1-3       | 4            | 5 - 7                        |
  // Typ     dir_sons  | bool         | dir_3D                       |
  // Inhalt: Ecke      | Kantenpunkt? | Richtung der Kante           |
  // -----------------------------------------------------------------


  // Einfache Printfunktionen:
  void Print(dir_3D d);
  void Print(dir_sons d);
  void Print(Edges_cell ed);
};


// 6. inline Implementierung einiger Memberfunktionenen


inline Bo_description::Bo_description() {
  static int i;
  for(i=0;i<8;++i)  corner[i] = false;
  for(i=0;i<24;++i) edge_poi[i]   = false;
  for(i=0;i<24;++i) dis[i]    = 0.0;
  edge = 0;
}

inline void Bo_description::Put_zero() {
  static int i;
  for(i=0;i<8;++i)  corner[i] = false;
  for(i=0;i<24;++i) edge_poi[i]   = false;
  for(i=0;i<24;++i) dis[i]    = 0.0;
  edge = 0;
}

inline bool Bo_description::corner_point(dir_sons e) {
  return corner[e];
};  

inline bool Bo_description::edge_point(dir_sons e,dir_3D d) {
  return edge_poi[3*e + (d>>1)];
}

inline void Bo_description::Set_num_corner(dir_sons e, int num) {
  num_corner[e]=num;
}
inline void Bo_description::Set_num_edge_poi(dir_sons e, dir_3D d, int num) {
  num_edge_poi[3*e + (d>>1)] = num;
};
inline void Bo_description::Set_num_cellpoi(int num) {
  num_cellpoi = num;
};

inline int Bo_description::Num_corner(dir_sons e) {
  /*
 if(e<0 || e>10) { // yet another uninitialized obj
   e=(dir_sons)0;
   printf("%4d ../boundary.h : Ne=%d\n",
          __LINE__,e);
 }
 fflush(stdout);                  
  */
  return num_corner[e];
}
inline int Bo_description::Num_edge_poi(dir_sons e, dir_3D d) {
  return num_edge_poi[3*e + (d>>1)];
};

inline int Bo_description::Num(Edge_Corner_point poi) {
  if(poi.edge_point == edge_poi_typ)   return Num_edge_poi(poi.corner,poi.d);
  if(poi.edge_point == corner_poi_typ) return Num_corner(poi.corner);
  return num_cellpoi; 
}


inline bool Bo_description::edge_point(Edge_Corner_point poi) {
  return edge_point(poi.corner,poi.d);
};

inline double Bo_description::h(dir_sons e,dir_3D d) {
  return dis[3*e + (d>>1)];
}

inline double Bo_description::Meshsize() {
  return H;
}

inline void Bo_description::Put_Meshsize(double h) {
  H=h;
}

inline void Bo_description::Put_corner_point(dir_sons e) {
  corner[e]=true;
}

inline void Bo_description::Put_edge_point(dir_sons e,dir_3D d,double h) {
  edge_poi[3*e + (d>>1)] = true;
  dis[3*e + (d>>1)] = h;
  // New:
  Put_corner_point(e);
}

inline bool Bo_description::Edge(Edges_cell ed) {
  return (bool)((edge>>ed)&1);
}

inline int Bo_description::edge_info() { return edge; };

inline void Bo_description::Put_edge(Edges_cell ed) {
  edge = edge|(1<<ed);
  // New:
  Put_corner_point(Transform(ed,Ld));
  Put_corner_point(Transform(ed,Rd));
}

inline D3vector Bo_description::coord(dir_sons e) {
  return D3vector((e&1)*H,((e>>1)&1)*H,(e>>2)*H);
};

inline D3vector Bo_description::coord(dir_sons e,dir_3D d) {
  if(d < Sdir)
     return D3vector((e&1)*H +  h(e,d) *(2*d-1), 
                     ((e>>1)&1)*H,(e>>2)*H);
  if(d > Ndir)
    return D3vector((e&1)*H,((e>>1)&1)*H,
                    (e>>2)*H + h(e,d) * (2*d-9));
  return D3vector((e&1)*H, 
                  ((e>>1)&1)*H + h(e,d) * (2*d-5),(e>>2)*H);

};

inline D3vector Bo_description::coord(Edge_Corner_point poi) {
  if(poi.edge_point == edge_poi_typ) return coord(poi.corner,poi.d);
  if(poi.edge_point == corner_poi_typ) return coord(poi.corner);
  return Give_bocedata()->Coord_bocellpoint();
}


inline D3vector BoCeData::Local_coord_bocellpoint() {
  if(developer_version) {
    if(bocellpoint == NULL)
      cout << " \n Fehler in BoCeData::Local_coord_bocellpoint!" << endl;
  }
  return bocellpoint->Give_local_coord();
}

inline D3vector BoCeData::Coord_bocellpoint() {
  if(developer_version) {
    if(bocellpoint == NULL)
      cout << " \n Fehler in BoCeData::Local_coord_bocellpoint!" << endl;
  }
  return bocellpoint->Give_local_coord() + D3vector(0.5,0.5,0.5)*H;
}

inline D3vector BoCeData::Coord_bocellpoint_normal() const {
  return bocellpoint->Give_local_coord() / H + D3vector(0.5,0.5,0.5);
}

inline int BoCeData::Number_avs() {
  if(developer_version) {
    if(bocellpoint == NULL)
      cout << " \n Fehler in BoCeData::Number_avs!" << endl;
  }
  return bocellpoint->Number_avs();
}

inline void BoCeData::Set_number_avs(int avs_num) {
  if(developer_version) {
    if(bocellpoint == NULL)
      cout << " \n Fehler in BoCeData::Set_number_avs!" << endl;
  }
  bocellpoint->Set_number_avs(avs_num);
}

inline void BoCeData::Add_Bo_freedom(D3vector local_coord) {
  if(developer_version) {
    if(bocellpoint != NULL)
      cout << " \n Randzellfreiheitsgrad existiert schon! " << endl;
  }
  bocellpoint = new Bocellpoint(local_coord);
}


inline void BoCeData::Add_Tetraeder(Bo_description* desc,
                                    int P0, int P1, int P2, int P3) {
  Tets = new Tetraeder_storage(desc,P0,P1,P2,P3,Tets);
}

inline void BoCeData::Add_Tetraeder_bo(Bo_description* desc,
                                       int P0, int P1, int P2, int P3) {
  Tets_boundary =
    new Boundary_tetraeder_storage(desc,P0,P1,P2,P3,Tets_boundary);
}


inline D3vector BoCeData::coord(int num, double *h) {
  dir_sons e;
  dir_3D d;
  if(developer_version) {
    if(num<0 || num > number_points) {
      cout << " \n Fehler in BoCeData::coord " << endl;
      return D3vector(0,0,0);
    }
  }
  if(num == number_points) {
    return bocellpoint->Give_local_coord() + D3vector(0.5,0.5,0.5)*H;
  }
  e = corner(num);
  if(edge_point(num)) {
    d = edge_dir(num);
    if(d < Sdir)
      return D3vector((e&1)*H +  h[num] *(2*d-1), 
                      ((e>>1)&1)*H,(e>>2)*H);
    if(d > Ndir)
      return D3vector((e&1)*H,((e>>1)&1)*H,
                      (e>>2)*H + h[num] * (2*d-9));
    return D3vector((e&1)*H, 
                    ((e>>1)&1)*H + h[num] * (2*d-5),(e>>2)*H);
  }
  else {
    return D3vector((e&1)*H,((e>>1)&1)*H,(e>>2)*H);
  }
}

                           

inline bool BoCeData::Edge(Edges_cell ed) {
  return (bool)((edge>>ed)&1);
}

inline void BoCeData::Put_edge_info(int ed) { edge=ed; };

inline double BoCeData::Meshsize() {
  return H;
}

inline void BoCeData::Put_Meshsize(double h) {
  H=h;
}



inline int BoCeData::Give_number_points() const {
  return number_points;
};


inline dir_sons BoCeData::corner(int num) const {
  if(developer_version) {
    if(num<0 || num >= number_points) {
      cout << " \n Fehler in BoCeData::corner " << endl;
      return WSDd;
    }
  }
  return (dir_sons)(info[num]&7);
}; 

inline Cell_type_points BoCeData::edge_point(int num) const {
  if(developer_version) {
    if(num<0 || num > number_points) {
      cout << " \n Fehler in BoCeData::edge_point " << endl;
      return corner_poi_typ;
    }
  }
  if( num == number_points) return cell_poi_typ;
  else return (Cell_type_points)((info[num]>>3)&1);
};


inline dir_3D BoCeData::edge_dir(int num) {
  if(developer_version) {
    if(num<0 || num >= number_points) {
      cout << " \n Fehler in BoCeData::edge_dir " << endl;
      return Wdir;
    }
  }
  return (dir_3D)((info[num]>>4)&7);
}


inline bool BoCeData::boundary_EW(int num) {
  if(developer_version) {
    if(num<0 || num >= number_points) {
      cout << " \n Fehler in BoCeData::boundary_EW " << endl;
      return false;
    }
  }
  return (bool)((info[num]>>7)&1);
};  


inline bool BoCeData::boundary_NS(int num) {
  if(developer_version) {
    if(num<0 || num >= number_points) {
      cout << " \n Fehler in BoCeData::boundary_NS " << endl;
      return false;
    }
  }
  return (bool)((info[num]>>8)&1);
};  


inline bool BoCeData::boundary_TD(int num) {
  if(developer_version) {
    if(num<0 || num >= number_points) {
      cout << " \n Fehler in BoCeData::boundary_TD " << endl;
      return false;
    }
  }
  return (bool)((info[num]>>9)&1);
};  


inline void BoCeData::Put_corner_point(int num, dir_sons c) {
  if(developer_version) {
    if(num<0 || num >= number_points) {
      cout << " \n Fehler in BoCeData::Put_corner " << endl;
      return;
    }
  }
  info[num] = c; 
}

inline void BoCeData::Put_edge_point(int num,dir_sons c, dir_3D d) {
  if(developer_version) {
    if(num<0 || num >= number_points) {
      cout << " \n Fehler in BoCeData::Put_edge_point " << endl;
      return;
    }
  }
  info[num] = c; 
  info[num] =  info[num] | (1<<3); 
  info[num] = (info[num] & (~(7<<4))) | (d<<4); 
}

inline void BoCeData::Allocate(int n_points) {
  static int i;
  number_points = n_points;
  info   = new int[number_points];
  for(i=0;i<number_points;i++) info[i]=0;
  Tets   = NULL; 
}

inline void BoCeData::Set_variable_pointer(int num, double* poi) {
  vars[num]=poi;
};

inline Tetraeder_storage::Tetraeder_storage(Bo_description* desc,
                                     int P0, int P1, int P2, int P3,
                                     Tetraeder_storage* Next) {
    next = Next; num0 = P0; num1 = P1; num2 = P2; num3 = P3; 
    var= NULL;

    if(developer_version) {
      Check_angles(desc);
    }
};

inline bool Tetraeder_storage::Check_angles(Bo_description* desc) {
  bool oh_key;

  oh_key = true;
  // 1. Calc Edge_Corner_point of Tet
  Edge_Corner_point ec0;
  Edge_Corner_point ec1;
  Edge_Corner_point ec2;
  Edge_Corner_point ec3;
  
  ec0.edge_point=desc->Give_bocedata()->edge_point(num0);
  if(ec0.edge_point!=cell_poi_typ)
    ec0.corner    =desc->Give_bocedata()->corner(num0);
  if(ec0.edge_point==edge_poi_typ)
    ec0.d = desc->Give_bocedata()->edge_dir(num0);
  
  
  ec1.edge_point=desc->Give_bocedata()->edge_point(num1);
  if(ec1.edge_point!=cell_poi_typ)
    ec1.corner    =desc->Give_bocedata()->corner(num1);
  if(ec1.edge_point==edge_poi_typ)
    ec1.d = desc->Give_bocedata()->edge_dir(num1);
  
  
  ec2.edge_point=desc->Give_bocedata()->edge_point(num2);
  if(ec2.edge_point!=cell_poi_typ)
    ec2.corner    =desc->Give_bocedata()->corner(num2);
  if(ec2.edge_point==edge_poi_typ)
    ec2.d = desc->Give_bocedata()->edge_dir(num2);
  
  ec3.edge_point=desc->Give_bocedata()->edge_point(num3);
  if(ec3.edge_point!=cell_poi_typ)
    ec3.corner    =desc->Give_bocedata()->corner(num3);
  if(ec3.edge_point==edge_poi_typ)
    ec3.d = desc->Give_bocedata()->edge_dir(num3);
  
  // Berechnung der Winkel
  double max_face_ang;
  double max_edge_ang;
  max_face_ang = calc_maximal_face_angle(desc->coord(ec0),
                                         desc->coord(ec1),
                                         desc->coord(ec2),
                                         desc->coord(ec3));
  if(max_face_ang > 163) {
    cout << "\n Error: Too large interior angle!" <<
      " Tell Christoph Pflaum this error! ";
    cout << "\n Maximal interior face angle: " << max_face_ang << endl;
    desc->coord(ec0).Print();  cout << endl;
    desc->coord(ec1).Print();  cout << endl;
    desc->coord(ec2).Print();  cout << endl;
    desc->coord(ec3).Print();  cout << endl;
    oh_key = false;
  }
  max_edge_ang = calc_maximal_edge_angle(desc->coord(ec0),
                                         desc->coord(ec1),
                                         desc->coord(ec2),
                                         desc->coord(ec3));
  if(max_edge_ang > 144) {
    cout << "\n Error: Too large interior angle!" <<
      " Tell Christoph Pflaum this error! ";
    cout << "\n Maximal interior edge angle: " << max_edge_ang << endl;
    desc->coord(ec0).Print();  cout << endl;
    desc->coord(ec1).Print();  cout << endl;
    desc->coord(ec2).Print();  cout << endl;
    desc->coord(ec3).Print();  cout << endl;
    oh_key = false;
  }
  return oh_key;
};



// For Ralf
inline int  BoCeData::Give_total_number_points() const {
  if(Exists_bocellpoint()) return number_points+1;
  else                     return number_points;
}

inline bool BoCeData::Is_point_clipped_point(int num) const {
  if(num<number_points) return true;
  return false;
}

#endif

---