src/expde/grid/hash.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
//

// ------------------------------------------------------------
//
//  hash.h
//
// ------------------------------------------------------------


//  Fuer Hashtable 0: Zellen und Kanten
// -------------------------------------
class Point_hashtable0 {
  friend class Grid_base;
  friend class Grid;
 private:
  Index3D ind;      // Mittelpunkt der Kante oder der Zelle
  int typ;          // typ der Kante bzw. Zelle
  bool isCell_m;      // I.Cell_index() 
  int  depth_m;
  bool send_m;

  // Speicher fuer Zellen
  double *var;

  Point_hashtable0* next;
 public:
  //  I Index des  Mittelpunktes der Zelle oder Kante
  Point_hashtable0(Index3D I) : ind(I) { 
    isCell_m = I.Cell_index(); 
    depth_m  = I.Tiefe();
    typ=noeinfo;
    next = NULL;
    var = NULL;
  };
  Point_hashtable0(Point_hashtable0* point_old);
  Index3D Give_Index() { return ind; }
  double *Give_var() { return var; }
  Celltype Give_cell_typ();
  Point_hashtable0* Next() { return next; }
  bool isCell() { return isCell_m; }
  bool doSend() {  return send_m; }
  void will_be_sent(bool s) { send_m = s; }
  int Give_Tiefe() { return depth_m; }
  ~Point_hashtable0() {
    delete next;
  };
};

inline int hashtable0_function(int x, int y, int z, int length) {
  //  return (x + 179*y + 2137*z) % length;
  return (x + 101*y + 10007*z) % length;
}

//  Fuer Iteration durch Hashtable 0: Grid
// -----------------------------------------
// folgendes define verwendet : point0:
#define iterate_hash0  for(i_iter=0;i_iter<hashtable0_leng;++i_iter) for(point0=hashtable0_start[i_iter];point0!=NULL;point0=point0->next) 

#define Iterate_hash0  for(i_iter=0;i_iter<hashtable0_leng;++i_iter) for(point0=hashtable0_start[i_iter];point0!=NULL;point0=point0->Next()) 


//  Fuer Hashtable 1: Grid
// ------------------------
class Point_hashtable1 {
 friend class Grid_base;
 friend class Grid;
 private:
  Index3D ind;
  Pointtype typ;               // type of the point: interior,  exterior, 
                               //   not responsible points are parallel_p
  Pointtype global_typ;        // type of the point: interior,  exterior, ...
                               //   multigrid is always correct
                               // needed in rers_op.h

  int nummer;                  // Nummer des Punktes (fuer AVS)
  bool label4;                 // Label zum Allocieren von Hashtable4 MG  
  int  depth_m;


  // different level informations
  int level;                   // falls randnaher Punkt: 
                               //       level der lokalen Maschenweite
                               // falls uniformer Punkt:
                               //       level des groebsten 27-Sternes

  int my_finest_parallel_level; // fuer Parallelisierung
                               
  int finest_level;            // feinster moeglicher level fuer multigrid
                               // Punkte und 
                               // andere Punkte; wird auch in Print AVS verw.
  //  bool responsible_mg;         // responsibility for a multigrid point

  // groebster level bezueglich Label fuer multigrid Punkte
  int *label_level;
  types_of_label_send send_label;  // used in: labbo.cc

  Point_hashtable1* next;
 public:
  // I Index des Punktes
  Point_hashtable1(Index3D I);
  Point_hashtable1(Index3D I,Pointtype);
  Point_hashtable1(Point_hashtable1* point_old);
  Index3D Give_Index() { return ind; }
  Pointtype Give_typ() { return typ; }
  int  Give_label(int num);               // label heisst hier ob Dirichlet
  void Put_label_level(int num, int lev); // Rand fuer Level
                                          // groesser gleich lev
  int Give_Tiefe() { return depth_m; }
  ~Point_hashtable1() {
    delete next;
    if(label_level!=NULL) {
      delete label_level;
    }
  };
};

inline int hashtable1_function(int x, int y, int z, int length) {
  //  return (x + 179*y + 2137*z) % length;
  return (x + 101*y + 10007*z) % length;
}


//  Fuer Iteration durch Hashtable 1: Grid
// -----------------------------------------
// folgendes define verwendet : point1:
#define iterate_hash1  for(i_iter=0;i_iter<hashtable1_leng;++i_iter) for(point1=hashtable1_start[i_iter];point1!=NULL;point1=point1->next) 



//  Fuer Hashtable 2: Randzellen, Boundary Cell
//  Beachte obiges: #define Point_hashtable2 BoCell 
// -------------------------------------------------
class BoCell : public BoCeData {
 friend class Grid_base;
 friend class Grid;
 private:
  Index3D ind;

  // Speicher fuer Freiheitsgrad in Zelle; existiert nicht immer
  double *var;

  Point_hashtable2* next;
 public:
  BoCell(Index3D I, Grid_base* gitter);
  Index3D Give_Index() { return ind; }
  double *Give_var() { return var; }
  Point_hashtable2* Next() { return next; }
  ~BoCell() {
    delete next;
  };
};

inline int hashtable2_function(int x, int y, int z, int length) {
  //  return (x + 179*y + 2137*z) % length;
  return (x + 101*y + 10007*z) % length;
}


//  Fuer Iteration durch Hashtable 2: Randzellen
// -----------------------------------------------
// folgendes define verwendet : bocell:
#define iterate_hash2  for(i_iter=0;i_iter<hashtable2_leng;++i_iter) for(bocell=hashtable2_start[i_iter];bocell!=NULL;bocell=bocell->next) 

#define Iterate_hash2  for(i_iter=0;i_iter<hashtable2_leng;++i_iter) for(bocell=hashtable2_start[i_iter];bocell!=NULL;bocell=bocell->Next()) 




//  Fuer Hashtable 3: Randpunkte 
// --------------------------------
class Bo2Point {
 friend class Grid_base;
 friend class Grid;
 private:
  Index3D ind;
  dir_3D  direction;
  int level;
  double h;                // Abstand zum Rand gerastert
  double h_origional;      // Exacter Abstand zum Rand
  double *var;
  int nummer;              // for avs

  // Label fuer Randpunkte
  TypLabel label;

  Point_hashtable3* next;
 public:
  inline Bo2Point(Index3D I, dir_3D d, int l, Grid_base* gitter,
                  int number_var);
  void Set_number_avs(int numavs) { nummer = numavs; };
  // Koordinate des Punktes gerastert
  inline D3vector transform_coord(D3vector A, double H);
  // Koordinate des Orginalpunktes verschoben
  inline D3vector transform_coord_org_shift(D3vector A, double h, 
                                            double H, double shift);
  int Nummer() { return nummer; }

  void Put_label(bool lab, int num);    // label heisst hier
  bool Give_label(int num);             // Punkt im Gebiet oder nicht
  IndexBo Give_boindex() const { return IndexBo(ind,direction); }
  ~Bo2Point() {
    delete next;
  };
};

inline int hashtable3_function(int x, int y, int z, dir_3D  d, int length) {
  //  return (d + 10*x + 1079*y + 20137*z) % length;
  return (d + 11*x + 101*y + 10007*z) % length;
}


//  Fuer Iteration durch Hashtable 3: Randpunkte 2.Teil
// ----------------------------------------------------
// folgendes define verwendet : bo2point:
#define iterate_hash3  for(i_iter=0;i_iter<hashtable3_leng;++i_iter) for(bo2point=hashtable3_start[i_iter];bo2point!=NULL;bo2point=bo2point->next) 



//  Fuer Hashtable 4: Datenspeicher
// ----------------------------------
class Data {
 friend class Grid_base;
 private:
  Index3D ind;
  int level;
  double *var;

  Point_hashtable4* next;
 public:
  double *Get_var() { return var; };
  Data(Index3D I, int l, int number_var);
  ~Data() {
    delete next;
  };
};

// hier Achtung, dass keine negative Zahl entsteht:
// Parameter_MG_Variable is Parameter aus parameter.h
inline int hashtable4_function(int x, int y, int z, int ebene, int length) {
  // return ((100 + Parameter_MG_Variable) + 11*x + 1079*y + 20137*z) % length;
  return ((100 + Parameter_MG_Variable) + 11*x + 101*y + 10007*z) % length;
}


//  Fuer Iteration durch Hashtable 4: Datenspeicher
// -------------------------------------------------
// folgendes define verwendet : varpoint:
#define iterate_hash4  for(i_iter=0;i_iter<hashtable4_leng;++i_iter) for(varpoint=hashtable4_start[i_iter];varpoint!=NULL;varpoint=varpoint->next) 


// Implementierung einiger Memberfunktionen

// allgemein
inline double *newdouble(int number) {
  static int i;
  static double* sp;
  sp = new double[number];
  for(i=0;i<number;++i) sp[i]=0.0;
  return sp;
}

// hash 0

inline Point_hashtable0* Grid_base::hashtable0_point(Index3D I) const {
  Point_hashtable0* point;
  point = hashtable0_start
    [hashtable0_function(I.ind_x.index,I.ind_y.index,I.ind_z.index,
                         hashtable0_leng)];
  while(point!=NULL) {
    if(point->ind==I) return point;
    point = point->next;
  }
  return NULL;
}

inline Point_hashtable0::Point_hashtable0(Point_hashtable0* point_old) {
  ind = point_old->ind;
  typ = point_old->typ;
  var = point_old->var;

  isCell_m = point_old->isCell_m;
  depth_m  = point_old->depth_m;

  next=NULL;
}

inline Celltype Point_hashtable0::Give_cell_typ() {
  if(developer_version) {
    if(!ind.Cell_index()) {
      cout << "\n Error 1 in Point_hashtable0::Give_cell_typ!" << endl;
      cout << " Point: "; ind.coordinate().Print(); cout << endl;
    }
  }
  return (Celltype)typ;
}

// hash 1

inline Point_hashtable1::Point_hashtable1(Index3D I) {
  ind = I;
  typ = exterior;
  next=NULL;
  //  responsible_mg = false;

  depth_m  = I.Tiefe();

  label_level=NULL;

  finest_level=0;
  my_finest_parallel_level=0;
  label4 = false;
}

inline Point_hashtable1::Point_hashtable1(Index3D I, Pointtype type) {
  ind = I;
  typ = type;
  next=NULL;
  //  responsible_mg = false;

  depth_m  = I.Tiefe();

  label_level=NULL;

  finest_level=0;
  my_finest_parallel_level=0;
  label4 = false;
}

inline void Point_hashtable1::Put_label_level(int num, int lev) {
  label_level[num-1] = lev; 
}

inline int Point_hashtable1::Give_label(int num) {
  return label_level[num-1]; 
}


inline Point_hashtable1::Point_hashtable1(Point_hashtable1* point_old) {
  ind = point_old->ind;
  typ = point_old->typ;

  label_level=NULL;

  nummer = point_old->nummer;
  label4 = point_old->label4;
  finest_level=0;
  my_finest_parallel_level=0;
  next=NULL;

  depth_m  = point_old->depth_m;
}


inline Point_hashtable1* Grid_base::hashtable1_point(Index3D I) const {
  Point_hashtable1* point;
  point = hashtable1_start
    [hashtable1_function(I.ind_x.index,I.ind_y.index,I.ind_z.index,
                         hashtable1_leng)];
  while(point!=NULL) {
    if(point->ind==I) return point;
    point = point->next;
  }
  return NULL;
}

// hash 2

inline Point_hashtable2::Point_hashtable2(Index3D I, Grid_base* gitter) {
  ind = I;  
  var = NULL;
  next=NULL;
}

inline Point_hashtable2* Grid_base::hashtable2_point(Index3D I) const {
  static Point_hashtable2* point;
  point = hashtable2_start
    [hashtable2_function(I.ind_x.index,I.ind_y.index,I.ind_z.index,
                         hashtable2_leng)];
  while(point!=NULL) {
    if(point->ind==I) return point;
    point = point->next;
  }
  return NULL;
}


// hash 3

inline Point_hashtable3::Point_hashtable3(Index3D I,  dir_3D d, int l,
                                   Grid_base* gitter, int number_var) {
  double h_min, h_max;
  ind = I;
  direction = d;
  level = l;
  label=0;

  // Minimaler Abstand am Rand 
  //^^^^^parallel^^^^^^^^^ deswegen braucht man hier -1
  // h_min = gitter->H_mesh() / Zweipotenz(l*2-1) * factor_boundary_dis; //old

  h_min =  gitter->h_min_for_boundary_points;  //new

  //  h_min = gitter->finest_mesh_size() / 4.0; // (for testing)
  //  h_min = 0.0;  // Test
  // Maximaler Abstand am Rand
  //  h_max = gitter->H_mesh() / Zweipotenz(l) - h_min;//old
  h_max = gitter->finest_mesh_size() - h_min; //new
  h = gitter->distanceD(I.coordinate(),d);
  // Test 
  //  h = gitter->finest_mesh_size() * 0.25;

  h_origional = h;


  if(h<h_min)      h = h_min;
  else if(h>h_max) h = h_max;

  //  h_max = gitter->H_mesh() / Zweipotenz(l); //old
  h_max = gitter->finest_mesh_size();       //new
  h_min = 0.0;
  if(h_origional<h_min)      h_origional = h_min;
  else if(h_origional>h_max) h_origional = h_max;
  if(developer_version==true) {
    if(h>h_max*1.2) 
      cout << "\n Fehler bei Randabstand! (zu gross)! " << endl;  
    if(h<0.0)
      cout << "\n Fehler bei Randabstand! (negativ)! " << endl;  
  }  

  next=NULL;
  if(number_var>0)    var   = newdouble(number_var);
  else var   = NULL;
}



inline D3vector Point_hashtable3::transform_coord(D3vector A,
                                                  double H) {
  A.x = A.x + ind.coordinate_x() * H;
  A.y = A.y + ind.coordinate_y() * H;
  A.z = A.z + ind.coordinate_z() * H;

  if(direction==Wdir)
    return D3vector(A.x-h,A.y,A.z);
  if(direction==Edir)
    return D3vector(A.x+h,A.y,A.z);
  if(direction==Ndir)
    return D3vector(A.x,A.y+h,A.z);
  if(direction==Sdir)
    return D3vector(A.x,A.y-h,A.z);
  if(direction==Tdir)
    return D3vector(A.x,A.y,A.z+h);
  //  if(direction==Ddir)
  return D3vector(A.x,A.y,A.z-h);
};


// Orginal Randpunkt wird leucht verschoben
inline D3vector Point_hashtable3::transform_coord_org_shift(D3vector A,
                                                            double h_s,
                                                            double H, 
                                                            double shift) {
  double h_shift;
  // Rand_punkt wird hier leicht verschoben
  h_shift = h_origional + h_s * shift;

  A.x = A.x + ind.coordinate_x() * H;
  A.y = A.y + ind.coordinate_y() * H;
  A.z = A.z + ind.coordinate_z() * H;

  if(direction==Wdir)
    return D3vector(A.x-h_shift,A.y,A.z);
  if(direction==Edir)
    return D3vector(A.x+h_shift,A.y,A.z);
  if(direction==Ndir)
    return D3vector(A.x,A.y+h_shift,A.z);
  if(direction==Sdir)
    return D3vector(A.x,A.y-h_shift,A.z);
  if(direction==Tdir)
    return D3vector(A.x,A.y,A.z+h_shift);
  //  if(direction==Ddir)
  return D3vector(A.x,A.y,A.z-h_shift);
};



inline void Point_hashtable3::Put_label(bool lab, int num) {
  int i, z;
  if(developer_version) {
    if(num>Max_label_number || num<=0) {
      cout << " Label number wrong in Point_hashtable3::Put_label! " << endl;
      cout << " num is: " << num << endl;
      num=0;
    } 
  }
  z=0;
  for(i=0;i<num-1;++i) {
    z = (z<<1) + 1; 
  }

  label = ((((label>>num)<<1)+lab)<<(num-1)) + (label&z);
};

inline bool Point_hashtable3::Give_label(int num) {
  if(developer_version) {
    if(num>Max_label_number || num<=0) {
      cout << " Label number wrong in Point_hashtable3::Give_label! " << endl;
      num=0;
    } 
  }
  return (bool)((label>>(num-1))&1);
};


inline Point_hashtable3* Grid_base::hashtable3_point(IndexBo IBo) const {
  static Point_hashtable3* point;
  point = hashtable3_start
    [hashtable3_function(IBo.I.ind_x.index,IBo.I.ind_y.index,IBo.
                         I.ind_z.index, IBo.d, hashtable3_leng)];
  while(point!=NULL) {
    if(point->ind==IBo.I && point->direction==IBo.d) return point;
    point = point->next;
  }
  return NULL;
}

inline Point_hashtable3* Grid_base::hashtable3_point(Index3D Ind, dir_3D dir)
     const {
  static Point_hashtable3* point;
  point = hashtable3_start
    [hashtable3_function(Ind.ind_x.index,Ind.ind_y.index,
                         Ind.ind_z.index, dir, hashtable3_leng)];
  while(point!=NULL) {
    if(point->ind==Ind && point->direction==dir) return point;
    point = point->next;
  }
  return NULL;
}


inline double* Grid_base::Give_variable(Index3D I, int ebene) const {
  static Point_hashtable4* point;
  point = hashtable4_start
    [hashtable4_function(I.ind_x.index,I.ind_y.index,I.ind_z.index, ebene,
                         hashtable4_leng)];
  if(developer_version) {
    if(point==NULL) {
      cout << "\n Error 1 in Give_variable (H4) "
           << my_rank 
           << " x: " << I.I_x().get()
           << " y: " << I.I_y().get()
           << " z: " << I.I_z().get()
           << endl;
      I.Print();
      I.coordinate().Print(); cout << " ebene: " << ebene << endl;
      return NULL;
    }
  }
  while(true) {
    if(point->ind==I && point->level==Parameter_MG_Variable) return point->var;
    point = point->next;
    if(developer_version) {
      if(point==NULL) {
        cout << "\n Error 2 in Give_variable (H4) " 
             << my_rank 
             << " x: " << I.I_x().get()
             << " y: " << I.I_y().get()
             << " z: " << I.I_z().get()
             << endl;
        I.coordinate().Print(); cout << " ebene: " << ebene << endl;
        return NULL;
      }
    }
  }
}

inline double* Grid_base::Give_variable_slow(const Index3D I, int ebene) const {
  static Point_hashtable4* point;
  point = hashtable4_start
    [hashtable4_function(I.ind_x.index,I.ind_y.index,I.ind_z.index, ebene, 
                         hashtable4_leng)];
  while(point!=NULL) {
    if(point->ind==I && point->level==Parameter_MG_Variable) return point->var;
    point = point->next;
  }
  return NULL;
}


// Elementares fuer hashtable4
// ----------------------------

inline Point_hashtable4::Point_hashtable4(Index3D I, int ebene, 
                                          int number_var) {
  ind = I;
  // Parameter_MG_Variable is Parameter aus parameter.h
  level = Parameter_MG_Variable;
  var = newdouble(number_var);
  next=NULL;
}


inline Point_hashtable4 *Grid_base::hashtable4_point(Index3D I, int ebene) {
  static Point_hashtable4* point;
  point = hashtable4_start
    [hashtable4_function(I.ind_x.index,I.ind_y.index,
                         I.ind_z.index, ebene, hashtable4_leng)];
  while(point!=NULL) {
    // Parameter_MG_Variable is Parameter aus parameter.h
    if(point->ind==I && point->level==Parameter_MG_Variable) return point;
    point = point->next;
  }
  return NULL;
}

inline Point_hashtable4* Grid_base::hashtable4_point(Index3D I) {
  static Point_hashtable4* point;
  point = hashtable4_start
    [hashtable4_function(I.ind_x.index,I.ind_y.index,
                         I.ind_z.index, 0, hashtable4_leng)];
  while(point!=NULL) {
    if(point->ind==I && point->level==0) return point;
    point = point->next;
  }
  return NULL;
}

---