Sierra Toolkit Version of the Day
HexFixture.cpp
00001 /*------------------------------------------------------------------------*/
00002 /*                 Copyright 2010, 2011 Sandia Corporation.                     */
00003 /*  Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive   */
00004 /*  license for use of this work by or on behalf of the U.S. Government.  */
00005 /*  Export of this program may require a license from the                 */
00006 /*  United States Government.                                             */
00007 /*------------------------------------------------------------------------*/
00008 
00009 #include <algorithm>
00010 
00011 #include <stk_util/environment/ReportHandler.hpp>
00012 
00013 #include <stk_mesh/fixtures/HexFixture.hpp>
00014 
00015 #include <stk_mesh/base/FieldData.hpp>
00016 #include <stk_mesh/base/Types.hpp>
00017 #include <stk_mesh/base/Entity.hpp>
00018 #include <stk_mesh/base/BulkModification.hpp>
00019 
00020 #include <stk_mesh/fem/Stencils.hpp>
00021 #include <stk_mesh/fem/FEMHelpers.hpp>
00022 #include <stk_mesh/fem/BoundaryAnalysis.hpp>
00023 
00024 namespace stk {
00025 namespace mesh {
00026 namespace fixtures {
00027 
00028 HexFixture::HexFixture(stk::ParallelMachine pm, unsigned nx, unsigned ny, unsigned nz)
00029   : m_spatial_dimension(3),
00030     m_nx(nx),
00031     m_ny(ny),
00032     m_nz(nz),
00033     m_fem_meta( m_spatial_dimension, fem::entity_rank_names(m_spatial_dimension) ),
00034     m_bulk_data( stk::mesh::fem::FEMMetaData::get_meta_data(m_fem_meta) , pm ),
00035     m_hex_part( fem::declare_part<shards::Hexahedron<8> >(m_fem_meta, "hex_part") ),
00036     m_node_part( m_fem_meta.declare_part("node_part", m_fem_meta.node_rank() ) ),
00037     m_coord_field( m_fem_meta.declare_field<CoordFieldType>("Coordinates") ),
00038     m_coord_gather_field(m_fem_meta.declare_field<CoordGatherFieldType>("GatherCoordinates"))
00039 {
00040   typedef shards::Hexahedron<8> Hex8 ;
00041   const unsigned nodes_per_elem = Hex8::node_count;
00042 
00043   //put coord-field on all nodes:
00044   put_field(
00045     m_coord_field,
00046     fem::FEMMetaData::NODE_RANK,
00047     m_fem_meta.universal_part(),
00048     m_spatial_dimension);
00049 
00050   //put coord-gather-field on all elements:
00051   put_field(
00052     m_coord_gather_field,
00053     m_fem_meta.element_rank(),
00054     m_fem_meta.universal_part(),
00055     nodes_per_elem);
00056 
00057   // Field relation so coord-gather-field on elements points
00058   // to coord-field of the element's nodes
00059   m_fem_meta.declare_field_relation( m_coord_gather_field,
00060                                      fem::element_node_stencil<Hex8, 3>,
00061                                      m_coord_field);
00062 }
00063 
00064 void HexFixture::generate_mesh()
00065 {
00066   std::vector<EntityId> element_ids_on_this_processor;
00067 
00068   const unsigned p_size = m_bulk_data.parallel_size();
00069   const unsigned p_rank = m_bulk_data.parallel_rank();
00070   const unsigned num_elems = m_nx * m_ny * m_nz ;
00071 
00072   const EntityId beg_elem = 1 + ( num_elems * p_rank ) / p_size ;
00073   const EntityId end_elem = 1 + ( num_elems * ( p_rank + 1 ) ) / p_size ;
00074 
00075   for ( EntityId i = beg_elem; i != end_elem; ++i) {
00076     element_ids_on_this_processor.push_back(i);
00077   }
00078 
00079   generate_mesh(element_ids_on_this_processor);
00080 }
00081 
00082 void HexFixture::node_x_y_z( EntityId entity_id, unsigned &x , unsigned &y , unsigned &z ) const
00083 {
00084   entity_id -= 1;
00085 
00086   x = entity_id % (m_nx+1);
00087   entity_id /= (m_nx+1);
00088 
00089   y = entity_id % (m_ny+1);
00090   entity_id /= (m_ny+1);
00091 
00092   z = entity_id;
00093 }
00094 
00095 void HexFixture::elem_x_y_z( EntityId entity_id, unsigned &x , unsigned &y , unsigned &z ) const
00096 {
00097   entity_id -= 1;
00098 
00099   x = entity_id % m_nx;
00100   entity_id /= m_nx;
00101 
00102   y = entity_id % m_ny;
00103   entity_id /= m_ny;
00104 
00105   z = entity_id;
00106 }
00107 
00108 void HexFixture::generate_mesh(std::vector<EntityId> & element_ids_on_this_processor)
00109 {
00110   {
00111     //sort and unique the input elements
00112     std::vector<EntityId>::iterator ib = element_ids_on_this_processor.begin();
00113     std::vector<EntityId>::iterator ie = element_ids_on_this_processor.end();
00114 
00115     std::sort( ib, ie);
00116     ib = std::unique( ib, ie);
00117     element_ids_on_this_processor.erase(ib, ie);
00118   }
00119 
00120   m_bulk_data.modification_begin();
00121 
00122   {
00123     // Declare the elements that belong on this process
00124 
00125     PartVector add_parts;
00126     add_parts.push_back(&m_node_part);
00127 
00128     std::vector<EntityId>::iterator ib = element_ids_on_this_processor.begin();
00129     const std::vector<EntityId>::iterator ie = element_ids_on_this_processor.end();
00130     for (; ib != ie; ++ib) {
00131       EntityId entity_id = *ib;
00132       unsigned ix = 0, iy = 0, iz = 0;
00133       elem_x_y_z(entity_id, ix, iy, iz);
00134 
00135       stk::mesh::EntityId elem_node[8] ;
00136 
00137       elem_node[0] = node_id( ix   , iy   , iz   );
00138       elem_node[1] = node_id( ix+1 , iy   , iz   );
00139       elem_node[2] = node_id( ix+1 , iy   , iz+1 );
00140       elem_node[3] = node_id( ix   , iy   , iz+1 );
00141       elem_node[4] = node_id( ix   , iy+1 , iz   );
00142       elem_node[5] = node_id( ix+1 , iy+1 , iz   );
00143       elem_node[6] = node_id( ix+1 , iy+1 , iz+1 );
00144       elem_node[7] = node_id( ix   , iy+1 , iz+1 );
00145 
00146       stk::mesh::fem::declare_element( m_bulk_data, m_hex_part, elem_id( ix , iy , iz ) , elem_node);
00147 
00148       for (unsigned i = 0; i<8; ++i) {
00149         stk::mesh::Entity * const node = m_bulk_data.get_entity( fem::FEMMetaData::NODE_RANK , elem_node[i] );
00150         m_bulk_data.change_entity_parts(*node, add_parts);
00151 
00152         ThrowRequireMsg( node != NULL,
00153           "This process should know about the nodes that make up its element");
00154 
00155         // Compute and assign coordinates to the node
00156         unsigned nx = 0, ny = 0, nz = 0;
00157         node_x_y_z(elem_node[i], nx, ny, nz);
00158 
00159         Scalar * data = stk::mesh::field_data( m_coord_field , *node );
00160 
00161         data[0] = nx ;
00162         data[1] = ny ;
00163         data[2] = -(Scalar)nz ;
00164       }
00165     }
00166   }
00167 
00168   m_bulk_data.modification_end();
00169 }
00170 
00171 } // fixtures
00172 } // mesh
00173 } // stk
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines