oomph::OcTreeForest Class Reference

#include <octree.h>

Inheritance diagram for oomph::OcTreeForest:

Public Member Functions
	OcTreeForest (Vector< TreeRoot * > &trees_pt)
	OcTreeForest ()
	Default constructor (empty and broken) More...
	OcTreeForest (const OcTreeForest &dummy)=delete
	Broken copy constructor. More...
void	operator= (const OcTreeForest &)=delete
	Broken assignment operator. More...
virtual	~OcTreeForest ()
void	check_all_neighbours (DocInfo &doc_info)
void	open_hanging_node_files (DocInfo &doc_info, Vector< std::ofstream * > &output_stream)
unsigned	self_test ()
	Self test: Check neighbour finding routine. More...
OcTreeRoot *	octree_pt (const unsigned &i) const
OcTreeRoot *	oc_face_neigh_pt (const unsigned &i, const int &direction)
Vector< TreeRoot * >	oc_edge_neigh_pt (const unsigned &i, const int &direction)
void	construct_up_right_equivalents ()
	Construct the rotation schemes. More...
Public Member Functions inherited from oomph::TreeForest
	TreeForest (Vector< TreeRoot * > &trees_pt)
	TreeForest ()
	Default constructor (empty and broken) More...
	TreeForest (const TreeForest &dummy)=delete
	Broken copy constructor. More...
void	operator= (const TreeForest &)=delete
	Broken assignment operator. More...
virtual	~TreeForest ()
	Kill tree forest: Delete the constituent trees. More...
void	stick_leaves_into_vector (Vector< Tree * > &forest_nodes)
	Traverse forst and stick pointers to leaf "nodes" into Vector. More...
void	stick_all_tree_nodes_into_vector (Vector< Tree * > &all_forest_nodes)
	Traverse forest and stick pointers to all "nodes" into Vector. More...
void	close_hanging_node_files (DocInfo &doc_info, Vector< std::ofstream * > &output_stream)
unsigned	ntree ()
	Number of trees in forest. More...
TreeRoot *	tree_pt (const unsigned &i) const
	Return pointer to i-th tree in forest. More...
void	flush_trees ()
	Flush trees from forest. More...

Private Member Functions
void	find_neighbours ()
	Construct the neighbour scheme. More...

Additional Inherited Members
Protected Attributes inherited from oomph::TreeForest
Vector< TreeRoot * >	Trees_pt
	Vector containing the pointers to the trees. More...

Detailed Description

/////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////////// An OcTreeForest consists of a collection of OcTreeRoots. Each member tree can have neighbours to its L/R/U/D/F/B and DB/UP/... and the orientation of their compasses can differ, allowing for complex, unstructured meshes.

Constructor & Destructor Documentation

OcTreeForest() [1/3]

oomph::OcTreeForest::OcTreeForest

(

Vector< TreeRoot * > &

trees_pt

)

Constructor for OcTree forest: Pass Vector of (pointers to) trees.

/////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////// /////////////////////////////////////////////////////////////////// Constructor for OcTreeForest:

Pass:

• trees_pt[], the Vector of pointers to the constituent trees (OcTreeRoot objects)

                                                        : TreeForest(trees_pt)
  {
#ifdef LEAK_CHECK
    LeakCheckNames::OcTreeForest_build += 1;
#endif
 
    // Don't setup neighbours etc. if forest is empty
    if (trees_pt.size() == 0)
    {
      return;
    }
 
    using namespace OcTreeNames;
 
    // Construct the neighbour and rotation scheme, note that all neighbour
    // pointers must be set before the constructor is called
 
    //  MemoryUsage::doc_memory_usage(
    //   "before find_neighbours in octree forest constr");
 
    // setup the neighbour scheme
    find_neighbours();
 
    //  MemoryUsage::doc_memory_usage(
    //   "after find_neighbours in octree forest constr");
 
    // setup the rotation scheme
    construct_up_right_equivalents();
 
    //  MemoryUsage::doc_memory_usage(
    //   "after construct_up_right_equivalents in octree forest constr");
  }

References construct_up_right_equivalents(), and find_neighbours().

OcTreeForest() [2/3]

oomph::OcTreeForest::OcTreeForest ( )

inline

Default constructor (empty and broken)

    {
      throw OomphLibError("Don't call empty constructor for OcTreeForest!",
                          OOMPH_CURRENT_FUNCTION,
                          OOMPH_EXCEPTION_LOCATION);
    }

References OOMPH_CURRENT_FUNCTION, and OOMPH_EXCEPTION_LOCATION.

OcTreeForest() [3/3]

oomph::OcTreeForest::OcTreeForest ( const OcTreeForest &  dummy )

delete

Broken copy constructor.

~OcTreeForest()

virtual oomph::OcTreeForest::~OcTreeForest ( )

inlinevirtual

Destructor: Delete the constituent octrees (and thus the associated objects!)

{}

Member Function Documentation

check_all_neighbours()

void oomph::OcTreeForest::check_all_neighbours ( DocInfo &  doc_info )

virtual

Document and check all the neighbours of all the nodes in the forest. DocInfo object specifies the output directory and file numbers for the various files. If doc_info.disable_doc() has been called, no output is created.

Document and check all the neighbours of all the nodes in the forest. DocInfo object specifies the output directory and file numbers for the various files. If doc_info.is_doc_enabled()=false no output is created.

Implements oomph::TreeForest.

  {
    // Create vector of all elements in the tree
    Vector<Tree*> all_tree_nodes_pt;
    this->stick_all_tree_nodes_into_vector(all_tree_nodes_pt);
 
    // Face neighbours
    //----------------
    {
      // Create storage for the files
      std::ofstream neigh_file;
      std::ofstream neigh_txt_file;
      // If we are documenting, then do so
      if (doc_info.is_doc_enabled())
      {
        std::ostringstream fullname;
        fullname << doc_info.directory() << "/neighbours" << doc_info.number()
                 << ".dat";
        oomph_info << "opened " << fullname.str() << " to doc neighbours"
                   << std::endl;
        neigh_file.open(fullname.str().c_str());
        fullname.str("");
        fullname << doc_info.directory() << "/neighbours" << doc_info.number()
                 << ".txt";
        oomph_info << "opened " << fullname.str() << " to doc neighbours"
                   << std::endl;
        neigh_txt_file.open(fullname.str().c_str());
      }
 
      // Call the static member of OcTree function
      double max_error = 0.0;
      OcTree::doc_face_neighbours(
        all_tree_nodes_pt, neigh_file, neigh_txt_file, max_error);
      if (max_error > Tree::max_neighbour_finding_tolerance())
      {
        std::ostringstream error_stream;
        error_stream << "\nMax. error in octree neighbour finding: "
                     << max_error << " is too big" << std::endl;
        error_stream
          << "i.e. bigger than Tree::max_neighbour_finding_tolerance()="
          << Tree::max_neighbour_finding_tolerance() << std::endl;
 
        // Close the files if they were opened
        if (doc_info.is_doc_enabled())
        {
          neigh_file.close();
          neigh_txt_file.close();
        }
 
        throw OomphLibError(
          error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
      }
      else
      {
        oomph_info << "\nMax. error in octree neighbour finding: " << max_error
                   << " is OK" << std::endl;
        oomph_info
          << "i.e. less than OcTree::max_neighbour_finding_tolerance()="
          << OcTree::max_neighbour_finding_tolerance() << std::endl;
      }
 
      // Close the documentation files if they were opened
      if (doc_info.is_doc_enabled())
      {
        neigh_file.close();
        neigh_txt_file.close();
      }
    }
 
    // Edge neighbours
    //----------------
    {
      // Create storage for the files
      std::ofstream neigh_file;
      std::ofstream no_true_edge_file;
      std::ofstream neigh_txt_file;
      // If we are documenting, then do so
      if (doc_info.is_doc_enabled())
      {
        std::ostringstream fullname;
        fullname << doc_info.directory() << "/edge_neighbours"
                 << doc_info.number() << ".dat";
        neigh_file.open(fullname.str().c_str());
        fullname.str("");
        fullname << doc_info.directory() << "/no_true_edge" << doc_info.number()
                 << ".dat";
        no_true_edge_file.open(fullname.str().c_str());
        fullname.str("");
        fullname << doc_info.directory() << "/edge_neighbours"
                 << doc_info.number() << ".txt";
        neigh_txt_file.open(fullname.str().c_str());
      }
 
      // Call the static member of OcTree function
      double max_error = 0.0;
      // Get the maximum error between edge neighbours
      OcTree::doc_true_edge_neighbours(all_tree_nodes_pt,
                                       neigh_file,
                                       no_true_edge_file,
                                       neigh_txt_file,
                                       max_error);
      if (max_error > Tree::max_neighbour_finding_tolerance())
      {
        std::ostringstream error_stream;
        error_stream << "Max. error in octree edge neighbour finding: "
                     << max_error << " is too big" << std::endl;
        error_stream
          << "i.e. bigger than Tree::max_neighbour_finding_tolerance()="
          << Tree::max_neighbour_finding_tolerance() << std::endl;
 
        // Close the files if they were opened
        if (doc_info.is_doc_enabled())
        {
          neigh_file.close();
          no_true_edge_file.close();
          neigh_txt_file.close();
        }
 
        throw OomphLibError(
          error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
      }
      else
      {
        oomph_info << "Max. error in octree edge neighbour finding: "
                   << max_error << " is OK" << std::endl;
        oomph_info
          << "i.e. less than OcTree::max_neighbour_finding_tolerance()="
          << OcTree::max_neighbour_finding_tolerance() << std::endl;
      }
 
      // Close the documentation files if they were opened
      if (doc_info.is_doc_enabled())
      {
        neigh_file.close();
        no_true_edge_file.close();
        neigh_txt_file.close();
      }
    }
  }

References oomph::DocInfo::directory(), oomph::OcTree::doc_face_neighbours(), oomph::OcTree::doc_true_edge_neighbours(), oomph::DocInfo::is_doc_enabled(), MeshRefinement::max_error, oomph::Tree::max_neighbour_finding_tolerance(), oomph::DocInfo::number(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION, oomph::oomph_info, and oomph::TreeForest::stick_all_tree_nodes_into_vector().

construct_up_right_equivalents()

void oomph::OcTreeForest::construct_up_right_equivalents

(

)

Construct the rotation schemes.

Construct the rotation scheme for the octree forest. Note that all pointers to neighbours must have been allocated for this to work.

  {
    using namespace OcTreeNames;
 
    // Number of trees in forest
    unsigned numtrees = ntree();
 
    // nnode1d from first element (if it exists!)
    int nnode1d = 0;
    if (numtrees > 0)
    {
      nnode1d = Trees_pt[0]->object_pt()->nnode_1d();
    }
 
    OcTreeRoot* neigh_pt = 0;
 
    // Loop over all the trees
    //------------------------
    for (unsigned i = 0; i < numtrees; i++)
    {
      // Find the pointer to the Up neighbour
      //------------------------------------
      neigh_pt = oc_face_neigh_pt(i, U);
 
      // If there is a neighbour?
      if (neigh_pt != 0)
      {
        // Find the direction of the present tree, as viewed from the neighbour
        int direction = neigh_pt->direction_of_neighbour(octree_pt(i));
 
        // If up neighbour has a pointer to this tree
        if (direction != Tree::OMEGA)
        {
          // The direction to the element in the neighbour
          // must be equivalent to the down direction in this element
          // Hence, the up equivalent is the reflection of that direction
          octree_pt(i)->set_up_equivalent(neigh_pt,
                                          OcTree::Reflect_face[direction]);
 
          // The right equivalent is the direction to the equivalent of
          // the right face in the present element, which will
          // be connected to the UR edge of the present element,
          // but will not be the face adjacent to the U boundary
          // i.e. the "direction" face.
          // We find the local node numbers, in the neighbour, of
          // the UR edge of the present element
          int nod1 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt(nnode1d * nnode1d * nnode1d -
                                               1));
          int nod2 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt(nnode1d * nnode1d - 1));
 
          // Now get the other face connected to that edge in the
          // neighbour. It is the right equivalent
          octree_pt(i)->set_right_equivalent(
            neigh_pt,
            OcTree::get_the_other_face(nod1, nod2, nnode1d, direction));
        }
        else
        // If U neighbour does not have pointer to this tree, die
        {
          std::ostringstream error_stream;
          error_stream << "Tree " << i
                       << "'s Up neighbour has no neighbour pointer to Tree "
                       << i << std::endl;
 
          throw OomphLibError(error_stream.str(),
                              OOMPH_CURRENT_FUNCTION,
                              OOMPH_EXCEPTION_LOCATION);
        }
      }
 
 
      // Find the pointer to the Right neighbour
      //---------------------------------------
      neigh_pt = oc_face_neigh_pt(i, R);
 
      // If there is a neighbour?
      if (neigh_pt != 0)
      {
        // Find the direction of the present tree, as viewed from the neighbour
        int direction = neigh_pt->direction_of_neighbour(octree_pt(i));
 
        // If the neighbour has a pointer to this tree
        if (direction != Tree::OMEGA)
        {
          // The direction to the element in the neighbour
          // must be equivalent to the left direction in this element
          // Hence, the right equivalent is the reflection of that direction
          octree_pt(i)->set_right_equivalent(neigh_pt,
                                             OcTree::Reflect_face[direction]);
 
          // The up equivalent will be the direction to the equivalent of
          // the up face in the neighbour, which will be connected to the
          // UR edge of the present element, but will not be the face
          // adjacent to the R boundary, i.e. the "direction" face
          // We find the local node numbers, in the neighbour, of the
          // UR edge of the present element
          int nod1 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt(nnode1d * nnode1d * nnode1d -
                                               1));
          int nod2 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt(nnode1d * nnode1d - 1));
 
          // Now get the other face connected to that edge in the
          // neighbour. It is the up equivalent
          octree_pt(i)->set_up_equivalent(
            neigh_pt,
            OcTree::get_the_other_face(nod1, nod2, nnode1d, direction));
        }
        else
        {
          // If R neighbour does not have pointer to this tree, die
          std::ostringstream error_stream;
          error_stream << "Tree " << i
                       << "'s Right neighbour has no neighbour pointer to Tree "
                       << i << std::endl;
 
          throw OomphLibError(error_stream.str(),
                              OOMPH_CURRENT_FUNCTION,
                              OOMPH_EXCEPTION_LOCATION);
        }
      }
 
 
      // Find the pointer to the Down neighbour
      //--------------------------------------
      neigh_pt = oc_face_neigh_pt(i, D);
 
      // If there is a neighbour?
      if (neigh_pt != 0)
      {
        // Find the direction of the present tree, as viewed from the neighbour
        int direction = neigh_pt->direction_of_neighbour(octree_pt(i));
 
        // If the neighbour has a pointer to this element
        if (direction != Tree::OMEGA)
        {
          // The direction to the element in the neighbour must be
          // equivalent to the up direction
          octree_pt(i)->set_up_equivalent(neigh_pt, direction);
 
          // The right equivalent is the direction to the equivalent of
          // the right face in the present element, which will be
          // connected to the DR edge of the present element, but will
          // not be the face adjacent to the D boundary,
          // i.e. the "direction" face.
          // We find the local node numbers, in the neighbour, of
          // the RD edge of the present element
          int nod1 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt(nnode1d - 1));
          int nod2 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt((nnode1d * nnode1d + 1) *
                                               (nnode1d - 1)));
 
          // Now get the other face connected to that edge in the neighbour.
          // It is the right equivalent
          octree_pt(i)->set_right_equivalent(
            neigh_pt,
            OcTree::get_the_other_face(nod1, nod2, nnode1d, direction));
        }
        else
        {
          // If D neighbour does not have pointer to this tree, die
          std::ostringstream error_stream;
          error_stream << "Tree " << i
                       << "'s Down neighbour has no neighbour pointer to Tree "
                       << i << std::endl;
 
          throw OomphLibError(error_stream.str(),
                              OOMPH_CURRENT_FUNCTION,
                              OOMPH_EXCEPTION_LOCATION);
        }
      }
 
 
      // Find the pointer to the Left neighbour
      //--------------------------------------
      neigh_pt = oc_face_neigh_pt(i, L);
 
      // If there is a neighbour
      if (neigh_pt != 0)
      {
        // Find the direction of the present tree, as viewed from the neighbour
        int direction = neigh_pt->direction_of_neighbour(octree_pt(i));
 
        // If the neighbour has a pointer to the present element
        if (direction != Tree::OMEGA)
        {
          // The direction to the element in the neighbour is
          // must be equivalent to the right direction
          octree_pt(i)->set_right_equivalent(neigh_pt, direction);
 
          // The up equivalent is the direction to the equivalent of the
          // up face in the present element, which will be connected to
          // the UL edge of the present element, but will not
          // be the face adjacent to the L boundary, i.e. the "direction"
          // face.
          // We find the local node numbers, in the neighbour, of the UL
          // edge
          int nod1 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt((nnode1d - 1) * nnode1d));
          int nod2 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt((nnode1d * nnode1d - 1) *
                                               nnode1d));
 
          // Now get the other face connected to that edge in the
          // neighbour.It is the up equivalent
          octree_pt(i)->set_up_equivalent(
            neigh_pt,
            OcTree::get_the_other_face(nod1, nod2, nnode1d, direction));
        }
        else
        {
          // If L neighbour does not have pointer to this tree, die
          std::ostringstream error_stream;
          error_stream << "Tree " << i
                       << "'s Left neighbour has no neighbour pointer to Tree "
                       << i << std::endl;
 
          throw OomphLibError(error_stream.str(),
                              OOMPH_CURRENT_FUNCTION,
                              OOMPH_EXCEPTION_LOCATION);
        }
      }
 
 
      // Find the pointer to the Front neighbour
      //---------------------------------------
      neigh_pt = oc_face_neigh_pt(i, F);
 
      // If there is a neighbour?
      if (neigh_pt != 0)
      {
        // Find the direction of the present tree, as viewed from the neighbour
        int direction = neigh_pt->direction_of_neighbour(octree_pt(i));
 
        // If the neighbour has a pointer to the present element
        if (direction != Tree::OMEGA)
        {
          // The direction to the up face will be given by the
          // other face connected to the UF edge in the neighbour
          // i.e. the face that is not given by direction
          // We obtain the local node numbers, in the neighbour,
          // of the UF edge in the present element
          int nod1 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt(nnode1d * nnode1d * nnode1d -
                                               1));
          int nod2 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt((nnode1d * nnode1d - 1) *
                                               nnode1d));
 
          // We now get the other face connected to that edge
          // It is the up equivalent.
          octree_pt(i)->set_up_equivalent(
            neigh_pt,
            OcTree::get_the_other_face(nod1, nod2, nnode1d, direction));
 
          // The direction to the right face will be given by the
          // other face connected to the RF edge in the neighbour
          // i.e. the face that is not given by the direction
          // We get the local node numbers, in the neighbour,
          // of the RF edge in the present element
          nod1 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt(nnode1d * nnode1d * nnode1d -
                                               1));
          nod2 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt((nnode1d * nnode1d + 1) *
                                               (nnode1d - 1)));
 
          // We now get the other face connected to that edge
          // It is the right equivalent
          octree_pt(i)->set_right_equivalent(
            neigh_pt,
            OcTree::get_the_other_face(nod1, nod2, nnode1d, direction));
        }
        else
        {
          // If F neighbour does not have pointer to this tree, die
          std::ostringstream error_stream;
          error_stream << "Tree " << i
                       << "'s Front neighbour has no neighbour pointer to Tree "
                       << i << std::endl;
 
          throw OomphLibError(error_stream.str(),
                              OOMPH_CURRENT_FUNCTION,
                              OOMPH_EXCEPTION_LOCATION);
        }
      }
 
 
      // Find the pointer to the Back neighbour
      //--------------------------------------
      neigh_pt = oc_face_neigh_pt(i, B);
 
      // If there is a neighbour?
      if (neigh_pt != 0)
      {
        // Find the direction of the present tree, as viewed from the neighbour
        int direction = neigh_pt->direction_of_neighbour(octree_pt(i));
 
        // If the neighbour has a pointer to the present element
        if (direction != Tree::OMEGA)
        {
          // The direction to the up face will be given by the
          // other face connected to the UB edge of the present element
          // i.e. not the "direction" face
          // We obtain the local node numbers, in the neighbour,
          // of the UB edge in the present element
          int nod1 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt((nnode1d - 1) * nnode1d));
          int nod2 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt(nnode1d * nnode1d - 1));
 
          // We now get the other face connected to that edge
          // It is the up equivalent
          octree_pt(i)->set_up_equivalent(
            neigh_pt,
            OcTree::get_the_other_face(nod1, nod2, nnode1d, direction));
 
          // The direction to the right face will be given by
          // the other face connected to the RB edge of the present
          // element, i.e. not the direction face.
          // We obtain local node numbers, in the neighbour,
          // of the RB edge in the present element
          nod1 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt(nnode1d * nnode1d - 1));
          nod2 = neigh_pt->object_pt()->get_node_number(
            octree_pt(i)->object_pt()->node_pt(nnode1d - 1));
 
          // We now get the other face connected to that edge
          // It is the right equivalent
          octree_pt(i)->set_right_equivalent(
            neigh_pt,
            OcTree::get_the_other_face(nod1, nod2, nnode1d, direction));
        }
        else
        {
          // If B neighbour does not have pointer to this tree, die
          std::ostringstream error_stream;
          error_stream << "Tree " << i
                       << "'s Back neighbour has no neighbour pointer to Tree "
                       << i << std::endl;
 
          throw OomphLibError(error_stream.str(),
                              OOMPH_CURRENT_FUNCTION,
                              OOMPH_EXCEPTION_LOCATION);
        }
      }
 
 
      // EDGE NEIGHBOURS
      //----------------
 
      // Loop over all edges:
      for (int edge = LB; edge <= UF; edge++)
      {
        // Get vector to pointers of edge neighbours
        Vector<TreeRoot*> edge_neigh_pt = oc_edge_neigh_pt(i, edge);
 
        // Loop over all edge neighbours
        unsigned n_neigh = edge_neigh_pt.size();
        for (unsigned e = 0; e < n_neigh; e++)
        {
          // If there is a neighbour
          if (edge_neigh_pt[e] != 0)
          {
            // Here are the two vertices at the two ends of the edge
            // in the present element
            int vertex = OcTree::Vertex_at_end_of_edge[edge][0];
            int other_vertex = OcTree::Vertex_at_end_of_edge[edge][1];
 
            // Get node numbers of the two vertices on edge
            unsigned nod1 = OcTree::vertex_to_node_number(vertex, nnode1d);
            unsigned nod2 =
              OcTree::vertex_to_node_number(other_vertex, nnode1d);
 
 
            // Here are the local node numbers (in the neighbouring element)
            // of the start/end nodes of present element's edge:
            unsigned neighb_nod1 =
              edge_neigh_pt[e]->object_pt()->get_node_number(
                octree_pt(i)->object_pt()->node_pt(nod1));
 
            unsigned neighb_nod2 =
              edge_neigh_pt[e]->object_pt()->get_node_number(
                octree_pt(i)->object_pt()->node_pt(nod2));
 
            // Convert to vertices
            int neighb_vertex =
              OcTree::node_number_to_vertex(neighb_nod1, nnode1d);
            int neighb_other_vertex =
              OcTree::node_number_to_vertex(neighb_nod2, nnode1d);
 
            // Up equivalent is stored first in the pair that's returned
            // by the lookup table
            octree_pt(i)->set_up_equivalent(
              edge_neigh_pt[e],
              OcTree::Up_and_right_equivalent_for_pairs_of_vertices
                [std::make_pair(
                   std::make_pair(neighb_vertex, vertex),
                   std::make_pair(neighb_other_vertex, other_vertex))]
                  .first);
 
            // Right equivalent is stored second in the pair that's returned
            // by the lookup table
            octree_pt(i)->set_right_equivalent(
              edge_neigh_pt[e],
              OcTree::Up_and_right_equivalent_for_pairs_of_vertices
                [std::make_pair(
                   std::make_pair(neighb_vertex, vertex),
                   std::make_pair(neighb_other_vertex, other_vertex))]
                  .second);
          }
        }
      }
    } // end of loop over all trees.
  }

References D, oomph::OcTreeRoot::direction_of_neighbour(), e(), oomph::OcTreeNames::F, oomph::FiniteElement::get_node_number(), oomph::OcTree::get_the_other_face(), i, L, oomph::OcTreeNames::LB, oomph::OcTree::node_number_to_vertex(), oomph::TreeForest::ntree(), oomph::Tree::object_pt(), oc_edge_neigh_pt(), oc_face_neigh_pt(), octree_pt(), oomph::Tree::OMEGA, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION, R, oomph::OcTree::Reflect_face, oomph::OcTreeRoot::set_right_equivalent(), oomph::OcTreeRoot::set_up_equivalent(), oomph::TreeForest::Trees_pt, RachelsAdvectionDiffusion::U, oomph::OcTreeNames::UF, oomph::OcTree::Up_and_right_equivalent_for_pairs_of_vertices, oomph::OcTree::Vertex_at_end_of_edge, and oomph::OcTree::vertex_to_node_number().

Referenced by OcTreeForest().

find_neighbours()

void oomph::OcTreeForest::find_neighbours ( )

private

Construct the neighbour scheme.

setup the neighbour scheme : tells to each element in the forest which (if any) element is its {R,L,U,D,B,F} face neighbour and which is its {LB,RB,...,UF} edge neighbour.

  {
    using namespace OcTreeNames;
    unsigned numtrees = ntree();
    int n = 0; // to store nnode1d
    if (numtrees > 0)
    {
      n = Trees_pt[0]->object_pt()->nnode_1d();
    }
    else
    {
      throw OomphLibError(
        "Trying to setup the neighbour scheme for an empty forest",
        OOMPH_CURRENT_FUNCTION,
        OOMPH_EXCEPTION_LOCATION);
    }
 
 
    // Number of vertex nodes: 8
    unsigned n_vertex_node = 8;
 
    // Find potentially connected trees by identifying
    // those whose associated elements share a common vertex node
    std::map<Node*, std::set<unsigned>> tree_assoc_with_vertex_node;
 
    // Loop over all trees
    for (unsigned i = 0; i < numtrees; i++)
    {
      // Loop over the vertex nodes of the associated element
      for (unsigned j = 0; j < n_vertex_node; j++)
      {
        Node* nod_pt = dynamic_cast<BrickElementBase*>(Trees_pt[i]->object_pt())
                         ->vertex_node_pt(j);
        tree_assoc_with_vertex_node[nod_pt].insert(i);
      }
    }
 
 
    // For each tree we store a set of potentially neighbouring trees
    // i.e. trees that share at least one node
    Vector<std::set<unsigned>> potentially_neighb_tree(numtrees);
 
    // Loop over vertex nodes
    for (std::map<Node*, std::set<unsigned>>::iterator it =
           tree_assoc_with_vertex_node.begin();
         it != tree_assoc_with_vertex_node.end();
         it++)
    {
      // Loop over connected elements twice
      for (std::set<unsigned>::iterator it_el1 = it->second.begin();
           it_el1 != it->second.end();
           it_el1++)
      {
        unsigned i = (*it_el1);
        for (std::set<unsigned>::iterator it_el2 = it->second.begin();
             it_el2 != it->second.end();
             it_el2++)
        {
          unsigned j = (*it_el2);
          // These two elements are potentially connected
          if (i != j)
          {
            potentially_neighb_tree[i].insert(j);
          }
        }
      }
    }
 
 
    // Loop over all trees
    for (unsigned i = 0; i < numtrees; i++)
    {
      // Cast to OcTreeRoot
      OcTreeRoot* octree_root_i_pt = dynamic_cast<OcTreeRoot*>(Trees_pt[i]);
 
      // Loop over their potential neighbours
      for (std::set<unsigned>::iterator it = potentially_neighb_tree[i].begin();
           it != potentially_neighb_tree[i].end();
           it++)
      {
        unsigned j = (*it);
 
        // So far, we haven't identified the candidate element as a face
        // neighbour of element/tree i
        bool is_face_neighbour = false;
 
        // is it the Up neighbour ?
        bool is_Up_neighbour =
          ((Trees_pt[j]->object_pt()->get_node_number(
              Trees_pt[i]->object_pt()->node_pt(n * (n * n - 1))) != -1) &&
           (Trees_pt[j]->object_pt()->get_node_number(
              Trees_pt[i]->object_pt()->node_pt(n * n - 1)) != -1));
 
        // is it the Down neighbour ?
        bool is_Down_neighbour =
          ((Trees_pt[j]->object_pt()->get_node_number(
              Trees_pt[i]->object_pt()->node_pt(n * n * (n - 1))) != -1) &&
           (Trees_pt[j]->object_pt()->get_node_number(
              Trees_pt[i]->object_pt()->node_pt(n - 1)) != -1));
 
        // is it the Right neighbour ?
        bool is_Right_neighbour =
          ((Trees_pt[j]->object_pt()->get_node_number(
              Trees_pt[i]->object_pt()->node_pt(n * n * n - 1)) != -1) &&
           (Trees_pt[j]->object_pt()->get_node_number(
              Trees_pt[i]->object_pt()->node_pt(n - 1)) != -1));
 
        // is it the Left neighbour ?
        bool is_Left_neighbour =
          ((Trees_pt[j]->object_pt()->get_node_number(
              Trees_pt[i]->object_pt()->node_pt(n * (n * n - 1))) != -1) &&
           (Trees_pt[j]->object_pt()->get_node_number(
              Trees_pt[i]->object_pt()->node_pt(0)) != -1));
 
        // is it the Back neighbour ?
        bool is_Back_neighbour =
          ((Trees_pt[j]->object_pt()->get_node_number(
              Trees_pt[i]->object_pt()->node_pt(n * n - 1)) != -1) &&
           (Trees_pt[j]->object_pt()->get_node_number(
              Trees_pt[i]->object_pt()->node_pt(0)) != -1));
 
        // is it the Front neighbour ?
        bool is_Front_neighbour =
          ((Trees_pt[j]->object_pt()->get_node_number(
              Trees_pt[i]->object_pt()->node_pt(n * n * n - 1)) != -1) &&
           (Trees_pt[j]->object_pt()->get_node_number(
              Trees_pt[i]->object_pt()->node_pt(n * n * (n - 1))) != -1));
 
 
        if (is_Down_neighbour)
        {
          is_face_neighbour = true;
          Trees_pt[i]->neighbour_pt(D) = Trees_pt[j];
        }
        if (is_Up_neighbour)
        {
          is_face_neighbour = true;
          Trees_pt[i]->neighbour_pt(U) = Trees_pt[j];
        }
        if (is_Right_neighbour)
        {
          is_face_neighbour = true;
          Trees_pt[i]->neighbour_pt(R) = Trees_pt[j];
        }
        if (is_Left_neighbour)
        {
          is_face_neighbour = true;
          Trees_pt[i]->neighbour_pt(L) = Trees_pt[j];
        }
        if (is_Back_neighbour)
        {
          is_face_neighbour = true;
          Trees_pt[i]->neighbour_pt(B) = Trees_pt[j];
        }
        if (is_Front_neighbour)
        {
          is_face_neighbour = true;
          Trees_pt[i]->neighbour_pt(F) = Trees_pt[j];
        }
 
 
        // If it's not a face neighbour, it may still be an edge
        // neighbour. We check this by checking if the
        // vertex nodes coincide. Note: This test would also
        // evaluate to true for face neighbours but we've already
        // determined that the element is not a face neighbour!
        if (!is_face_neighbour)
        {
          // is it the left back neighbour ?
          bool is_left_back_neighbour =
            ((Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt(0)) != -1) &&
             (Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt(n * (n - 1))) != -1));
 
          // is it the right back neighbour ?
          bool is_right_back_neighbour =
            ((Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt(n - 1)) != -1) &&
             (Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt(n * n - 1)) != -1));
 
 
          // is it the down back neighbour ?
          bool is_down_back_neighbour =
            ((Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt(n - 1)) != -1) &&
             (Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt(0)) != -1));
 
          // is it the up back neighbour ?
          bool is_up_back_neighbour =
            ((Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt(n * (n - 1))) != -1) &&
             (Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt(n * n - 1)) != -1));
 
 
          // is it the left down neighbour ?
          bool is_left_down_neighbour =
            ((Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt(n * n * (n - 1))) != -1) &&
             (Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt(0)) != -1));
 
 
          // is it the right down neighbour ?
          bool is_right_down_neighbour =
            ((Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt((n * n + 1) * (n - 1))) !=
              -1) &&
             (Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt(n - 1)) != -1));
 
 
          // is it the left up neighbour ?
          bool is_left_up_neighbour =
            ((Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt((n * n * n - n))) != -1) &&
             (Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt(n * (n - 1))) != -1));
 
 
          // is it the right up neighbour ?
          bool is_right_up_neighbour =
            ((Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt((n * n * n - 1))) != -1) &&
             (Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt(n * n - 1)) != -1));
 
 
          // is it the left front neighbour ?
          bool is_left_front_neighbour =
            ((Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt((n * n * n - n))) != -1) &&
             (Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt(n * n * (n - 1))) != -1));
 
 
          // is it the right front neighbour ?
          bool is_right_front_neighbour =
            ((Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt((n * n * n - 1))) != -1) &&
             (Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt((n * n + 1) * (n - 1))) !=
              -1));
 
 
          // is it the down front neighbour ?
          bool is_down_front_neighbour =
            ((Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt(n * n * (n - 1))) != -1) &&
             (Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt((n * n + 1) * (n - 1))) !=
              -1));
 
 
          // is it the up front neighbour ?
          bool is_up_front_neighbour =
            ((Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt((n * n * n - n))) != -1) &&
             (Trees_pt[j]->object_pt()->get_node_number(
                Trees_pt[i]->object_pt()->node_pt((n * n * n - 1))) != -1));
 
 
          // Add to storage scheme for edge neighbours (only!)
 
          if (is_left_back_neighbour)
          {
            // Trees_pt[i]->neighbour_pt(LB)=Trees_pt[j];
            octree_root_i_pt->add_edge_neighbour_pt(Trees_pt[j], LB);
          }
          if (is_right_back_neighbour)
          {
            // Trees_pt[i]->neighbour_pt(RB)=Trees_pt[j];
            octree_root_i_pt->add_edge_neighbour_pt(Trees_pt[j], RB);
          }
          if (is_down_back_neighbour)
          {
            // Trees_pt[i]->neighbour_pt(DB)=Trees_pt[j];
            octree_root_i_pt->add_edge_neighbour_pt(Trees_pt[j], DB);
          }
          if (is_up_back_neighbour)
          {
            // Trees_pt[i]->neighbour_pt(UB)=Trees_pt[j];
            octree_root_i_pt->add_edge_neighbour_pt(Trees_pt[j], UB);
          }
 
 
          if (is_left_down_neighbour)
          {
            // Trees_pt[i]->neighbour_pt(LD)=Trees_pt[j];
            octree_root_i_pt->add_edge_neighbour_pt(Trees_pt[j], LD);
          }
          if (is_right_down_neighbour)
          {
            // Trees_pt[i]->neighbour_pt(RD)=Trees_pt[j];
            octree_root_i_pt->add_edge_neighbour_pt(Trees_pt[j], RD);
          }
          if (is_left_up_neighbour)
          {
            // Trees_pt[i]->neighbour_pt(LU)=Trees_pt[j];
            octree_root_i_pt->add_edge_neighbour_pt(Trees_pt[j], LU);
          }
          if (is_right_up_neighbour)
          {
            // Trees_pt[i]->neighbour_pt(RU)=Trees_pt[j];
            octree_root_i_pt->add_edge_neighbour_pt(Trees_pt[j], RU);
          }
 
 
          if (is_left_front_neighbour)
          {
            // Trees_pt[i]->neighbour_pt(LF)=Trees_pt[j];
            octree_root_i_pt->add_edge_neighbour_pt(Trees_pt[j], LF);
          }
          if (is_right_front_neighbour)
          {
            // Trees_pt[i]->neighbour_pt(RF)=Trees_pt[j];
            octree_root_i_pt->add_edge_neighbour_pt(Trees_pt[j], RF);
          }
          if (is_down_front_neighbour)
          {
            // Trees_pt[i]->neighbour_pt(DF)=Trees_pt[j];
            octree_root_i_pt->add_edge_neighbour_pt(Trees_pt[j], DF);
          }
          if (is_up_front_neighbour)
          {
            // Trees_pt[i]->neighbour_pt(UF)=Trees_pt[j];
            octree_root_i_pt->add_edge_neighbour_pt(Trees_pt[j], UF);
          }
        }
      }
    }
  }

References oomph::OcTreeRoot::add_edge_neighbour_pt(), D, oomph::OcTreeNames::DB, oomph::OcTreeNames::DF, oomph::OcTreeNames::F, i, j, L, oomph::OcTreeNames::LB, oomph::OcTreeNames::LD, oomph::OcTreeNames::LF, oomph::OcTreeNames::LU, n, oomph::TreeForest::ntree(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION, R, oomph::OcTreeNames::RB, oomph::OcTreeNames::RD, oomph::OcTreeNames::RF, oomph::OcTreeNames::RU, oomph::TreeForest::Trees_pt, RachelsAdvectionDiffusion::U, oomph::OcTreeNames::UB, and oomph::OcTreeNames::UF.

Referenced by OcTreeForest().

oc_edge_neigh_pt()

Vector<TreeRoot*> oomph::OcTreeForest::oc_edge_neigh_pt ( const unsigned &  i, const int &  direction  )

inline

Given the number i of the root octree in this forest, return the vector of pointers to the true edge neighbours in the specified (edge) direction.

    {
      // Note: paranoia check is done in edge_neighbour_pt
      return dynamic_cast<OcTreeRoot*>(Trees_pt[i])
        ->edge_neighbour_pt(direction);
    }

References oomph::OcTreeRoot::edge_neighbour_pt(), and i.

Referenced by construct_up_right_equivalents().

oc_face_neigh_pt()

OcTreeRoot* oomph::OcTreeForest::oc_face_neigh_pt ( const unsigned &  i, const int &  direction  )

inline

Given the number i of the root octree in this forest, return pointer to its face neighbour in the specified direction. NULL if neighbour doesn't exist. (This does the dynamic cast from a TreeRoot to a OcTreeRoot internally).

    {
#ifdef PARANOID
      using namespace OcTreeNames;
      if ((direction != U) && (direction != D) && (direction != F) &&
          (direction != B) && (direction != L) && (direction != R))
      {
        std::ostringstream error_stream;
        error_stream << "Wrong edge_direction: "
                     << OcTree::Direct_string[direction] << std::endl;
        throw OomphLibError(
          error_stream.str(), OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION);
      }
#endif
      return dynamic_cast<OcTreeRoot*>(Trees_pt[i]->neighbour_pt(direction));
    }

References oomph::OcTreeNames::D, oomph::OcTree::Direct_string, oomph::OcTreeNames::F, i, L, OOMPH_CURRENT_FUNCTION, OOMPH_EXCEPTION_LOCATION, R, and oomph::OcTreeNames::U.

Referenced by construct_up_right_equivalents().

octree_pt()

OcTreeRoot* oomph::OcTreeForest::octree_pt ( const unsigned &  i ) const

inline

Return pointer to i-th OcTree in forest (Performs a dynamic cast from the TreeRoot to a OcTreeRoot).

    {
      return dynamic_cast<OcTreeRoot*>(Trees_pt[i]);
    }

References i.

Referenced by construct_up_right_equivalents().

open_hanging_node_files()

void oomph::OcTreeForest::open_hanging_node_files ( DocInfo &  doc_info, Vector< std::ofstream * > &  output_stream  )

virtual

Open output files that will store any hanging nodes in the forest and return a vector of the streams.

Open output files that will stored any hanging nodes that are

created in the mesh refinement process.

Implements oomph::TreeForest.

  {
    // In 3D, there will be six output files
    for (unsigned i = 0; i < 6; i++)
    {
      output_stream.push_back(new std::ofstream);
    }
 
    // If we are documenting the output, open the files
    if (doc_info.is_doc_enabled())
    {
      std::ostringstream fullname;
      fullname << doc_info.directory() << "/hang_nodes_u" << doc_info.number()
               << ".dat";
      output_stream[0]->open(fullname.str().c_str());
      fullname.str("");
      fullname << doc_info.directory() << "/hang_nodes_d" << doc_info.number()
               << ".dat";
      output_stream[1]->open(fullname.str().c_str());
      fullname.str("");
      fullname << doc_info.directory() << "/hang_nodes_l" << doc_info.number()
               << ".dat";
      output_stream[2]->open(fullname.str().c_str());
      fullname.str("");
      fullname << doc_info.directory() << "/hang_nodes_r" << doc_info.number()
               << ".dat";
      output_stream[3]->open(fullname.str().c_str());
      fullname.str("");
      fullname << doc_info.directory() << "/hang_nodes_b" << doc_info.number()
               << ".dat";
      output_stream[4]->open(fullname.str().c_str());
      fullname.str("");
      fullname << doc_info.directory() << "/hang_nodes_f" << doc_info.number()
               << ".dat";
      output_stream[5]->open(fullname.str().c_str());
    }
  }

References oomph::DocInfo::directory(), i, oomph::DocInfo::is_doc_enabled(), and oomph::DocInfo::number().

operator=()

void oomph::OcTreeForest::operator= ( const OcTreeForest &  )

delete

Broken assignment operator.

self_test()

unsigned oomph::OcTreeForest::self_test

(

)

Self test: Check neighbour finding routine.

Self-test: Check all neighbours. Return success (0) if the max. distance between corresponding points in the neighbours is less than the tolerance specified in the static value Tree::Max_neighbour_finding_tolerance.

  {
    // Stick pointers to all nodes into Vector and number elements
    // in the process
    Vector<Tree*> all_nodes_pt;
    stick_all_tree_nodes_into_vector(all_nodes_pt);
    long int count = 0;
    unsigned long num_nodes = all_nodes_pt.size();
    for (unsigned long i = 0; i < num_nodes; i++)
    {
      all_nodes_pt[i]->object_pt()->set_number(++count);
    }
 
    // Check neighbours vertices and their opposite points
    std::ofstream neighbours_file;
    std::ofstream no_true_edge_file;
    std::ofstream neighbours_txt_file;
 
    double max_error_face = 0.0;
    OcTree::doc_face_neighbours(
      all_nodes_pt, neighbours_file, neighbours_txt_file, max_error_face);
 
    double max_error_edge = 0.0;
    OcTree::doc_true_edge_neighbours(all_nodes_pt,
                                     neighbours_file,
                                     no_true_edge_file,
                                     neighbours_txt_file,
                                     max_error_edge);
 
    bool failed = false;
    if (max_error_face > OcTree::max_neighbour_finding_tolerance())
    {
      oomph_info << "\n\n Failed self_test() for OcTreeForest because of "
                    "faces: Max. error "
                 << max_error_face << std::endl
                 << std::endl;
      failed = true;
    }
 
    if (max_error_edge > OcTree::max_neighbour_finding_tolerance())
    {
      oomph_info << "\n\n Failed self_test() for OcTreeForest because of "
                    "edges: Max. error "
                 << max_error_edge << std::endl
                 << std::endl;
      failed = true;
    }
 
    double max_error = max_error_face;
    if (max_error_edge > max_error) max_error = max_error_edge;
 
    if (failed)
    {
      return 1;
    }
    else
    {
      oomph_info << "\nPassed self_test() for OcTreeForest: Max. error "
                 << max_error << std::endl;
      return 0;
    }
  }

References oomph::OcTree::doc_face_neighbours(), oomph::OcTree::doc_true_edge_neighbours(), i, MeshRefinement::max_error, oomph::Tree::max_neighbour_finding_tolerance(), oomph::oomph_info, and oomph::TreeForest::stick_all_tree_nodes_into_vector().

---

The documentation for this class was generated from the following files:

• octree.h
• octree.cc