elem.h

一个用来实现偏微分方程中网格的计算库

  /**
   * Make the classes that need to access our build
   * member friends.  These classes do not really fit
   * the profile of what a "friend" should be, but
   * if we are going to protect the constructor and
   * the build method, there's no way around it.
   *
   * Do we *really* need to protect the build member?
   * It would seem that we are just getting around it
   * by using friends!
   */
  friend class MeshRefinement;    // (Elem::nullify_neighbors)

 private:
  /**
   * This function is used internally for node key generation.
   * It handles casting of pointers on various architectures.
   */
  unsigned int _cast_node_address_to_unsigned_int(const unsigned int n);

  // Prime numbers used for computing node keys.  These are the same
  // for every instance of the Elem class.
  static const unsigned int _bp1;
  static const unsigned int _bp2;
};




// ------------------------------------------------------------
// Elem class member functions
inline
Elem::Elem(const unsigned int nn,
	   const unsigned int ns,
	   Elem* p) :
  _parent(p)
#ifdef ENABLE_AMR
  , _p_level(0)
#endif
{
  this->subdomain_id() = 0;
  this->processor_id() = DofObject::invalid_processor_id;
  
  // Initialize the nodes data structure
  _nodes = NULL;
  
  if (nn != 0)
    {
      _nodes = new Node*[nn]; 
      
      for (unsigned int n=0; n<nn; n++)
	_nodes[n] = NULL;
    }
  
  // Initialize the neighbors data structure
  _neighbors = NULL;
  
  if (ns != 0)
    {
      _neighbors = new Elem*[ns]; 
      
      for (unsigned int n=0; n<ns; n++)
	_neighbors[n] = NULL;
    }

  // Optionally initialize data from the parent
  if (this->parent() != NULL)
    {
      this->subdomain_id() = this->parent()->subdomain_id();
      this->processor_id() = this->parent()->processor_id();
    }  

#ifdef ENABLE_AMR
  
  _children = NULL;

  this->set_refinement_flag(Elem::DO_NOTHING);

  this->set_p_refinement_flag(Elem::DO_NOTHING);

  if (this->parent())
    this->set_p_level(parent()->p_level());
  else
    this->set_p_level(0);

#endif
}



inline
Elem::~Elem() 
{
  // Delete my node storage
  if (_nodes != NULL)
    delete [] _nodes;
  _nodes = NULL;

  // Delete my neighbor storage
  if (_neighbors != NULL)
    delete [] _neighbors;
  _neighbors = NULL;

#ifdef ENABLE_AMR

  // Delete my children's storage
  if (_children != NULL)
    delete [] _children;
  _children = NULL;
  
#endif
}



inline
const Point & Elem::point (const unsigned int i) const
{
  libmesh_assert (i < this->n_nodes());
  libmesh_assert (_nodes[i] != NULL);
  libmesh_assert (_nodes[i]->id() != Node::invalid_id);

  return *_nodes[i];
}



inline
Point & Elem::point (const unsigned int i)
{
  libmesh_assert (i < this->n_nodes());

  return *_nodes[i];
}



inline
unsigned int Elem::node (const unsigned int i) const
{
  libmesh_assert (i < this->n_nodes());
  libmesh_assert (_nodes[i] != NULL);
  libmesh_assert (_nodes[i]->id() != Node::invalid_id);

  return _nodes[i]->id();
}



inline
Node* Elem::get_node (const unsigned int i) const
{
  libmesh_assert (i < this->n_nodes());
  libmesh_assert (_nodes[i] != NULL);

  return _nodes[i];
}



inline
Node* & Elem::set_node (const unsigned int i)
{
  libmesh_assert (i < this->n_nodes());

  return _nodes[i];
}



inline
unsigned char Elem::subdomain_id () const
{
  return _sbd_id;
}



inline
unsigned char & Elem::subdomain_id ()
{
  return _sbd_id;
}



inline
Elem* Elem::neighbor (const unsigned int i) const
{
  libmesh_assert (i < this->n_neighbors());

  return _neighbors[i];
}



inline
void Elem::set_neighbor (const unsigned int i, Elem* n)
{
  libmesh_assert (i < this->n_neighbors());
  
  _neighbors[i] = n;
}



inline
bool Elem::has_neighbor (const Elem* elem) const
{
  for (unsigned int n=0; n<this->n_neighbors(); n++)
    if (this->neighbor(n) == elem)
      return true;

  return false;
}



inline
Elem* Elem::child_neighbor (Elem* elem) const
{
  for (unsigned int n=0; n<elem->n_neighbors(); n++)
    if (elem->neighbor(n) &&
	elem->neighbor(n)->parent() == this)
      return elem->neighbor(n);

  return NULL;
}



inline
const Elem* Elem::child_neighbor (const Elem* elem) const
{
  for (unsigned int n=0; n<elem->n_neighbors(); n++)
    if (elem->neighbor(n) &&
	elem->neighbor(n)->parent() == this)
      return elem->neighbor(n);

  return NULL;
}



inline
bool Elem::on_boundary () const
{
  // By convention, the element is on the boundary
  // if it has a NULL neighbor.
  return this->has_neighbor(NULL);
}



inline
unsigned int Elem::which_neighbor_am_i (const Elem* e) const
{
  libmesh_assert (e != NULL);

  const Elem* eparent = e;

  while (eparent->level() > this->level())
    {
      eparent = eparent->parent();
      libmesh_assert(eparent);
    }
  
  for (unsigned int s=0; s<this->n_neighbors(); s++)
    if (this->neighbor(s) == eparent)
      return s;
    
  return libMesh::invalid_uint;
}



inline
bool Elem::active() const
{
#ifdef ENABLE_AMR
  if ((this->refinement_flag() == INACTIVE) ||
      (this->refinement_flag() == COARSEN_INACTIVE))
    return false;
  else
    return true;
#else
  return true;
#endif
}





inline
bool Elem::subactive() const
{
#ifdef ENABLE_AMR
  if (this->active())
    return false;
  if (!this->has_children())
    return true;
  return this->child(0)->subactive();
#else
  return false;
#endif
}



inline
bool Elem::has_children() const
{
#ifdef ENABLE_AMR
  if (_children == NULL)
    return false;
  else
    return true;
#else
  return false;
#endif
}


inline
bool Elem::has_ancestor_children() const
{
#ifdef ENABLE_AMR
  if (_children == NULL)
    return false;
  else
    for (unsigned int c=0; c != this->n_children(); c++)
      if (this->child(c)->has_children())
	return true;
#endif
  return false;
}



#ifdef ENABLE_AMR
inline
bool Elem::is_ancestor_of(const Elem *descendant) const
{
  const Elem *e = descendant;
  while (e)
    {
      if (this == e)
        return true;
      e = e->parent();
    }
#else
inline
bool Elem::is_ancestor_of(const Elem *) const
{
#endif
  return false;
}



inline
const Elem* Elem::parent () const
{
  return _parent;
}



inline
Elem* Elem::parent ()
{
  return _parent;
}



inline
void Elem::set_parent (Elem *p)
{
  _parent = p;
}



inline
const Elem* Elem::top_parent () const
{
  const Elem* tp = this;

  // Keep getting the element's parent
  // until that parent is at level-0
  while (tp->parent() != NULL)
    tp = tp->parent();
  
  libmesh_assert (tp != NULL);
  libmesh_assert (tp->level() == 0);

  return tp;  
}



inline
unsigned int Elem::level() const
{
#ifdef ENABLE_AMR

  // if I don't have a parent I was
  // created directly from file
  // or by the user, so I am a
  // level-0 element
  if (this->parent() == NULL)
    return 0;

  // otherwise we are at a level one
  // higher than our parent
  return (this->parent()->level() + 1);

#else

  // Without AMR all elements are
  // at level 0.
  return 0;
  
#endif
}



inline
unsigned int Elem::p_level() const
{
#ifdef ENABLE_AMR
  return _p_level;
#else
  return 0;
#endif
}



#ifdef ENABLE_AMR

inline
Elem* Elem::child (const unsigned int i) const
{
  libmesh_assert (_children    != NULL);
  libmesh_assert (_children[i] != NULL);
  
  return _children[i];
}



inline
unsigned int Elem::which_child_am_i (const Elem* e) const
{
  libmesh_assert (e != NULL);
  libmesh_assert (this->has_children());

  for (unsigned int c=0; c<this->n_children(); c++)
    if (this->child(c) == e)
      return c;

  std::cerr << "ERROR:  which_child_am_i() was called with a non-child!" 
	    << std::endl;

  libmesh_error();

  return libMesh::invalid_uint;
}



inline
Elem::RefinementState Elem::refinement_flag () const
{
  return static_cast<RefinementState>(_rflag);
}



inline
void Elem::set_refinement_flag(RefinementState rflag)
{
#ifdef DEBUG
  if (rflag != static_cast<RefinementState>(static_cast<unsigned char>(rflag)))
    {
      std::cerr << "ERROR: unsigned char too small to hold Elem::_rflag!"
		<< std::endl
		<< "Recompile with Elem:_*flag set to something bigger!"
		<< std::endl;
      libmesh_error();
    }
#endif

  _rflag = rflag;
}