refineable_brick_mesh.h

Go to the documentation of this file.

//LIC// ====================================================================
//LIC// This file forms part of oomph-lib, the object-oriented, 
//LIC// multi-physics finite-element library, available 
//LIC// at http://www.oomph-lib.org.
//LIC// 
//LIC//    Version 1.0; svn revision $LastChangedRevision$
//LIC//
//LIC// $LastChangedDate$
//LIC// 
//LIC// Copyright (C) 2006-2016 Matthias Heil and Andrew Hazel
//LIC// 
//LIC// This library is free software; you can redistribute it and/or
//LIC// modify it under the terms of the GNU Lesser General Public
//LIC// License as published by the Free Software Foundation; either
//LIC// version 2.1 of the License, or (at your option) any later version.
//LIC// 
//LIC// This library is distributed in the hope that it will be useful,
//LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
//LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
//LIC// Lesser General Public License for more details.
//LIC// 
//LIC// You should have received a copy of the GNU Lesser General Public
//LIC// License along with this library; if not, write to the Free Software
//LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
//LIC// 02110-1301  USA.
//LIC// 
//LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
//LIC// 
//LIC//====================================================================
#ifndef QMESH3OOMPH_D_HEADER
#define QMESH3OOMPH_D_HEADER


// Config header generated by autoconfig
#ifdef HAVE_CONFIG_H
  #include <oomph-lib-config.h>
#endif

#include<limits.h>


// ooomph-lib includes
#include "brick_mesh.h"
#include "refineable_mesh.h"
#include "refineable_brick_element.h"
//Include to fill in additional_setup_shared_node_scheme() function
#include "refineable_mesh.template.cc"

namespace oomph
{

//=======================================================================
/// Intermediate mesh class that implements the mesh adaptation functions
/// specified in the TreeBasedRefineableMesh class for meshes that contain the
/// refineable variant of QElement s [The class ELEMENT provided
/// as the template parameter must be of type 
/// RefineableQElement<3>].
/// 
/// Mesh adaptation/refinement is implemented by OcTree 
/// procedures and any concrete implementation of this class needs to
/// provide a OcTreeForest representation of the initial (coarse) mesh.
//=======================================================================
template <class ELEMENT>
class RefineableBrickMesh : public TreeBasedRefineableMesh<ELEMENT>, 
 public virtual BrickMeshBase
{

public:


 /// Constructor: Setup static octree data
 RefineableBrickMesh()
  {
   // OcTree static data needs to be setup before octree-based mesh 
   // refinement works
   OcTree::setup_static_data();
  }

 /// Broken copy constructor
 RefineableBrickMesh(const RefineableBrickMesh& dummy) 
  { 
   BrokenCopy::broken_copy("RefineableBrickMesh");
  } 
 
 /// Broken assignment operator
 void operator=(const RefineableBrickMesh&) 
  {
   BrokenCopy::broken_assign("RefineableBrickMesh");
  }

 /// Destructor:
 virtual ~RefineableBrickMesh() {}

 /// \short Set up the tree forest associated with the Mesh. 
 /// Forwards call to setup_octree_forest()
 virtual void setup_tree_forest()
  {
   setup_octree_forest();
  }

 /// Do what it says...
 void setup_octree_forest()
  {
   if (this->Forest_pt!=0)
    {
     // Get all the tree nodes
     Vector<Tree*> all_tree_nodes_pt;
     this->Forest_pt->stick_all_tree_nodes_into_vector(all_tree_nodes_pt);

     // Get min and max refinement level from the tree
     unsigned local_min_ref=0;
     unsigned local_max_ref=0;
     this->get_refinement_levels(local_min_ref,local_max_ref);

     unsigned min_ref=local_min_ref;
#ifdef OOMPH_HAS_MPI
     if (Comm_pt!=0)
      {
       // Reconcile between processors: If (e.g. following distribution/pruning)
       // the mesh has no elements on this processor) then ignore its
       // contribution to the poll of max/min refinement levels
       int int_local_min_ref=local_min_ref;       
       if (this->nelement()==0)
        {
         int_local_min_ref=INT_MAX;
        }
       int int_min_ref=0;
       MPI_Allreduce(&int_local_min_ref,&int_min_ref,1,
                     MPI_INT,MPI_MIN,
                     Comm_pt->mpi_comm());       
       min_ref=int_min_ref;
      }
#endif

     // If we have no elements there's nothing more to be done --
     // we only came in here to participate in the communication
     if (this->nelement()==0)
      {
       // Flush the Forest's current trees
       this->Forest_pt->flush_trees();
       
       // Delete the old Forest
       delete this->Forest_pt;

       // Empty dummy vector to build empty forest
       Vector<TreeRoot*> trees_pt;

       // Make a new (empty) Forest
       this->Forest_pt = new OcTreeForest(trees_pt);

       return;
      }


     // Vector to store trees for new Forest
     Vector<TreeRoot*> trees_pt;

     // Loop over tree nodes (e.g. elements)
     unsigned n_tree_nodes=all_tree_nodes_pt.size();
     for (unsigned e=0;e<n_tree_nodes;e++)
      {
       Tree* tree_pt=all_tree_nodes_pt[e];

       // If the object_pt has been flushed then we don't want to keep
       // this tree
       if (tree_pt->object_pt()!=0)
        {
         // Get the refinement level of the current tree node
         RefineableElement* el_pt=dynamic_cast<RefineableElement*>
          (tree_pt->object_pt());
         unsigned level=el_pt->refinement_level();

         // If we are below the minimum refinement level, remove tree
         if (level<min_ref)
          {
           // Flush sons for this tree
           tree_pt->flush_sons();

           // Delete the tree (no recursion)
           delete tree_pt;

           // Delete the element
           delete el_pt;
          }
         else if (level==min_ref)
          {
           // Get the sons (if there are any) and store them
           unsigned n_sons=tree_pt->nsons();
           Vector<Tree*> backed_up_sons(n_sons);
           for (unsigned i_son=0;i_son<n_sons;i_son++)
            {
             backed_up_sons[i_son]=tree_pt->son_pt(i_son);
            }

           // Make the element into a new treeroot
           OcTreeRoot* tree_root_pt=new OcTreeRoot(el_pt);


           // Pass sons
           tree_root_pt->set_son_pt(backed_up_sons);

           // Loop over sons and make the new treeroot their father
           for (unsigned i_son=0;i_son<n_sons;i_son++)
            {
             Tree* son_pt=backed_up_sons[i_son];

             // Tell the son about its new father (which is also the root)
             son_pt->set_father_pt(tree_root_pt);
             son_pt->root_pt()=tree_root_pt;

             // ...and then tell all the descendants too
             Vector<Tree*> all_sons_pt;
             son_pt->stick_all_tree_nodes_into_vector(all_sons_pt);
             unsigned n=all_sons_pt.size();
             for (unsigned i=0;i<n;i++)
              {
               all_sons_pt[i]->root_pt()=tree_root_pt;
              }
            }

           // Add tree root to the trees_pt vector
           trees_pt.push_back(tree_root_pt);

           // Now kill the original (non-root) tree: First
           // flush sons for this tree
           tree_pt->flush_sons();

           // ...then delete the tree (no recursion)
           delete tree_pt;
          }
        }
       else // tree_pt->object_pt() is null, so delete tree
        {
         // Flush sons for this tree
         tree_pt->flush_sons();

         // Delete the tree (no recursion)
         delete tree_pt;
        }
      }

     // Flush the Forest's current trees
     this->Forest_pt->flush_trees();

     // Delete the old Forest
     delete this->Forest_pt;

     // Make a new Forest with the trees_pt roots created earlier
     this->Forest_pt = new OcTreeForest(trees_pt);
    }
   else // Create a new Forest from scratch in the "usual" uniform way
    {
     //Turn elements into individual octrees and plant in forest
     Vector<TreeRoot*> trees_pt;
     unsigned nel=nelement();
     for (unsigned iel=0;iel<nel;iel++)
      {
       // Get pointer to full element type 
       ELEMENT* el_pt=dynamic_cast<ELEMENT*>(element_pt(iel));
     
       // Build associated octree(root) -- pass pointer to corresponding
       // finite element and add the pointer to vector of octree (roots):
       OcTreeRoot* octree_root_pt=new OcTreeRoot(el_pt);
       trees_pt.push_back(octree_root_pt);
      }
     // Plant OcTreeRoots in OcTreeForest
     this->Forest_pt = new OcTreeForest(trees_pt);
    }
  }

  protected:

};

}

#endif