Adding new layer remapping capaiblity to MeshExtruder

framework/include/meshmodifiers/MeshExtruder.h

@@ -16,6 +16,7 @@
 #define MESHEXTRUDER_H
 
 #include "MeshModifier.h"
+#include "libmesh/mesh_generation.h"
 
 class MeshExtruder;
 
@@ -35,8 +36,33 @@ class MeshExtruder : public MeshModifier
   const unsigned int _num_layers;
   const RealVectorValue _extrusion_vector;
 
+  /**
+   * This class is used during the mesh construction (extrusion) to set element ids as they are created.
+   */
+  class QueryElemSubdomainID : public MeshTools::Generation::QueryElemSubdomainIDBase
+  {
+  public:
+    QueryElemSubdomainID(std::vector<SubdomainID> existing_subdomains,
+                         std::vector<unsigned int> layers,
+                         std::vector<unsigned int> new_ids,
+                         unsigned int num_layers);
+
+    /// The override from the base class for obtaining a new id based on the old (original) element and the specified layer
+    virtual subdomain_id_type get_subdomain_for_layer(const Elem * old_elem, unsigned int layer);
+
+  private:
+    /// Data structure for holding the old -> new id mapping based on the layer number
+    std::map<unsigned int, std::map<SubdomainID, unsigned int> > _layer_data;
+
+    /// The total number of layers in the extrusion
+    unsigned int _num_layers;
+  };
+
 private:
   void changeID(const std::vector<BoundaryName> & names, BoundaryID old_id);
+
+//  MooseSharedPointer<QueryElemSubdomainID> _elem_subdomain_id;
+  bool _map_custom_ids;
 };
 
 #endif /* MESHEXTRUDER_H */

framework/src/meshmodifiers/MeshExtruder.C

@@ -18,6 +18,7 @@
 // libMesh includes
 #include "libmesh/mesh_generation.h"
 #include "libmesh/mesh.h"
+#include "libmesh/elem.h"
 #include "libmesh/boundary_info.h"
 
 template<>
@@ -28,14 +29,27 @@ InputParameters validParams<MeshExtruder>()
   params.addRequiredParam<RealVectorValue>("extrusion_vector", "The direction and length of the extrusion");
   params.addParam<std::vector<BoundaryName> >("bottom_sideset", "The boundary that will be applied to the bottom of the extruded mesh");
   params.addParam<std::vector<BoundaryName> >("top_sideset", "The boundary that will be to the top of the extruded mesh");
+
+  params.addParam<std::vector<SubdomainID> >("existing_subdomains", "The subdomains that will be remapped for specific layers");
+  params.addParam<std::vector<unsigned int> >("layers", "The layers where the \"existing_subdomain\" will be remapped to new ids");
+  params.addParam<std::vector<unsigned int> >("new_ids", "The list of new ids, This list should be either length \"existing_subdomains\" or "
+                                            "\"existing_subdomains\" * layers");
   return params;
 }
 
 MeshExtruder::MeshExtruder(const InputParameters & parameters):
     MeshModifier(parameters),
     _num_layers(getParam<unsigned int>("num_layers")),
-    _extrusion_vector(getParam<RealVectorValue>("extrusion_vector"))
+    _extrusion_vector(getParam<RealVectorValue>("extrusion_vector")),
+    _map_custom_ids(false)
 {
+  if (isParamValid("existing_subdomains") || isParamValid("layers") || isParamValid("new_ids"))
+  {
+    if (isParamValid("existing_subdomains") && isParamValid("layers") && isParamValid("new_ids"))
+      _map_custom_ids = true;
+    else
+      mooseError("If any of the following parameters are defined, all of them must be: \"existing_subdomains\", \"layers\", \"new_ids\"");
+  }
 }
 
 MeshExtruder::~MeshExtruder()
@@ -54,12 +68,18 @@ MeshExtruder::modify()
 
   _mesh_ptr->getMesh().clear();
 
+  QueryElemSubdomainID elem_subdomain_id(getParam<std::vector<SubdomainID> >("existing_subdomains"),
+                                         getParam<std::vector<unsigned int> >("layers"),
+                                         getParam<std::vector<unsigned int> >("new_ids"),
+                                         _num_layers);
+
   // The first argument to build_extrusion() is required to be UnstructuredMesh&, a common
   // base class of both SerialMesh and ParallelMesh, hence the dynamic_cast...
   MeshTools::Generation::build_extrusion(dynamic_cast<libMesh::UnstructuredMesh&>(_mesh_ptr->getMesh()),
                                          source_mesh->getMesh(),
                                          _num_layers,
-                                         _extrusion_vector);
+                                         _extrusion_vector,
+                                         _map_custom_ids ? &elem_subdomain_id : NULL);
 
   // See if the user has requested specific sides for the top and bottom
   const std::set<boundary_id_type> &side_ids = _mesh_ptr->getMesh().get_boundary_info().get_side_boundary_ids();
@@ -95,5 +115,48 @@ MeshExtruder::changeID(const std::vector<BoundaryName> & names, BoundaryID old_i
     _mesh_ptr->getMesh().get_boundary_info().sideset_name(boundary_ids[i]) = names[i];
 }
 
+MeshExtruder::QueryElemSubdomainID::QueryElemSubdomainID(std::vector<SubdomainID> existing_subdomains,
+                                                         std::vector<unsigned int> layers,
+                                                         std::vector<unsigned int> new_ids,
+                                                         unsigned int num_layers) :
+    QueryElemSubdomainIDBase(),
+    _num_layers(num_layers)
+{
+  // Check the length of the vectors
+  if (existing_subdomains.size() != new_ids.size() && existing_subdomains.size() * layers.size() != new_ids.size())
+    mooseError("The length of the \"existing_subdomains\", \"layers\", and \"new_ids\" are not valid");
+
+  // Setup our stride depending on whether the user passed unique sets in new ids or just a single set of new ids
+  const unsigned int zero = 0;
+  unsigned int i = 0;
+  const unsigned int stride = existing_subdomains.size() == new_ids.size() ? zero : existing_subdomains.size();
+
+  // Populate the data structure
+  for (i = 0; i < layers.size(); ++i)
+    for (unsigned int j = 0; j < existing_subdomains.size(); ++j)
+      _layer_data[layers[i]][existing_subdomains[j]] = new_ids[i*stride + j];
+}
 
+subdomain_id_type
+MeshExtruder::QueryElemSubdomainID::get_subdomain_for_layer(const Elem * old_elem, unsigned int layer)
+{
+  mooseAssert(layer < _num_layers, "Access out of bounds: " << layer);
 
+  // First locate the layer if it exists
+  std::map<unsigned int, std::map<SubdomainID, unsigned int> >::const_iterator layer_it = _layer_data.find(layer);
+
+  if (layer_it == _layer_data.end())
+    // If the layer wasn't found, there is no mapping so just return the original subdomain id
+    return old_elem->subdomain_id();
+  else
+  {
+    std::map<SubdomainID, unsigned int>::const_iterator sub_id_it = layer_it->second.find(old_elem->subdomain_id());
+
+    if (sub_id_it == layer_it->second.end())
+      // If the subdomain wasn't found, it won't be remapped, so just return the original subdomain id
+      return old_elem->subdomain_id();
+
+    // Return the remapped id
+    return sub_id_it->second;
+  }
+}

test/tests/mesh_modifiers/mesh_extruder/extrude_remap_layer1.i

@@ -0,0 +1,63 @@
+[Mesh]
+  type = FileMesh
+  file = multiblock.e
+[]
+
+[MeshModifiers]
+  [./extrude]
+    type = MeshExtruder
+    num_layers = 6
+    extrusion_vector = '0 0 2'
+    bottom_sideset = 'new_bottom'
+    top_sideset = 'new_top'
+
+    # Remap layers
+    existing_subdomains = '1 2 5'
+    layers = '1 3 5'
+    new_ids = '10 12 15
+               30 32 35
+               50 52 55'
+  [../]
+[]
+
+[Variables]
+  [./u]
+    order = FIRST
+    family = LAGRANGE
+  [../]
+[]
+
+[Kernels]
+  [./diff]
+    type = Diffusion
+    variable = u
+  [../]
+[]
+
+[BCs]
+  [./bottom]
+    type = DirichletBC
+    variable = u
+    boundary = 'new_bottom'
+    value = 0
+  [../]
+
+  [./top]
+    type = DirichletBC
+    variable = u
+    boundary = 'new_top'
+    value = 1
+  [../]
+[]
+
+[Executioner]
+  type = Steady
+
+  # Preconditioned JFNK (default)
+  solve_type = 'PJFNK'
+[]
+
+[Outputs]
+  output_initial = true
+  exodus = true
+[]

test/tests/mesh_modifiers/mesh_extruder/extrude_remap_layer2.i

@@ -0,0 +1,61 @@
+[Mesh]
+  type = FileMesh
+  file = multiblock.e
+[]
+
+[MeshModifiers]
+  [./extrude]
+    type = MeshExtruder
+    num_layers = 6
+    extrusion_vector = '0 0 2'
+    bottom_sideset = 'new_bottom'
+    top_sideset = 'new_top'
+
+    # Remap layers
+    existing_subdomains = '1 2 5'
+    layers = '1 3 5'
+    new_ids = '10 12 15' # Repeat this remapping for each layer
+  [../]
+[]
+
+[Variables]
+  [./u]
+    order = FIRST
+    family = LAGRANGE
+  [../]
+[]
+
+[Kernels]
+  [./diff]
+    type = Diffusion
+    variable = u
+  [../]
+[]
+
+[BCs]
+  [./bottom]
+    type = DirichletBC
+    variable = u
+    boundary = 'new_bottom'
+    value = 0
+  [../]
+
+  [./top]
+    type = DirichletBC
+    variable = u
+    boundary = 'new_top'
+    value = 1
+  [../]
+[]
+
+[Executioner]
+  type = Steady
+
+  # Preconditioned JFNK (default)
+  solve_type = 'PJFNK'
+[]
+
+[Outputs]
+  output_initial = true
+  exodus = true
+[]

test/tests/mesh_modifiers/mesh_extruder/tests

@@ -22,4 +22,16 @@
     input = 'gen_extrude.i'
     exodiff = 'out_gen.e'
   [../]
+
+  [./extrude_remap_layer1]
+    type = 'Exodiff'
+    input = 'extrude_remap_layer1.i'
+    exodiff = 'extrude_remap_layer1_out.e'
+  [../]
+
+  [./extrude_remap_layer2]
+    type = 'Exodiff'
+    input = 'extrude_remap_layer2.i'
+    exodiff = 'extrude_remap_layer2_out.e'
+  [../]
 []