Sierra Toolkit Version of the Day
QuadFixture.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/QuadFixture.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/BoundaryAnalysis.hpp>
00022 #include <stk_mesh/fem/FEMMetaData.hpp>
00023 #include <stk_mesh/fem/FEMHelpers.hpp>
00024 
00025 namespace stk {
00026 namespace mesh {
00027 namespace fixtures {
00028 
00029 QuadFixture::QuadFixture( stk::ParallelMachine pm ,
00030                           unsigned nx , unsigned ny )
00031   : m_spatial_dimension(2),
00032     m_fem_meta( m_spatial_dimension, fem::entity_rank_names(m_spatial_dimension) ),
00033     m_bulk_data( stk::mesh::fem::FEMMetaData::get_meta_data(m_fem_meta) , pm ),
00034     m_quad_part( fem::declare_part<shards::Quadrilateral<4> >(m_fem_meta, "quad_part" ) ),
00035     m_coord_field( m_fem_meta.declare_field<CoordFieldType>("Coordinates") ),
00036     m_coord_gather_field( m_fem_meta.declare_field<CoordGatherFieldType>("GatherCoordinates") ),
00037     m_nx( nx ),
00038     m_ny( ny )
00039 {
00040   typedef shards::Quadrilateral<4> Quad4 ;
00041   const unsigned nodes_per_elem = Quad4::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 
00051   //put coord-gather-field on all elements:
00052   put_field(
00053       m_coord_gather_field,
00054       m_fem_meta.element_rank(),
00055       m_fem_meta.universal_part(),
00056       nodes_per_elem
00057       );
00058 
00059   // Field relation so coord-gather-field on elements points
00060   // to coord-field of the element's nodes
00061   m_fem_meta.declare_field_relation(
00062       m_coord_gather_field,
00063       fem::element_node_stencil<Quad4, 2>,
00064       m_coord_field
00065       );
00066 }
00067 
00068 void QuadFixture::node_x_y( EntityId entity_id, unsigned &x , unsigned &y ) const
00069 {
00070   entity_id -= 1;
00071 
00072   x = entity_id % (m_nx+1);
00073   entity_id /= (m_nx+1);
00074 
00075   y = entity_id;
00076 }
00077 
00078 void QuadFixture::elem_x_y( EntityId entity_id, unsigned &x , unsigned &y ) const
00079 {
00080   entity_id -= 1;
00081 
00082   x = entity_id % m_nx;
00083   entity_id /= m_nx;
00084 
00085   y = entity_id;
00086 }
00087 
00088 
00089 void QuadFixture::generate_mesh()
00090 {
00091   std::vector<EntityId> element_ids_on_this_processor;
00092 
00093   const unsigned p_size = m_bulk_data.parallel_size();
00094   const unsigned p_rank = m_bulk_data.parallel_rank();
00095   const unsigned num_elems = m_nx * m_ny;
00096 
00097   const EntityId beg_elem = 1 + ( num_elems * p_rank ) / p_size ;
00098   const EntityId end_elem = 1 + ( num_elems * ( p_rank + 1 ) ) / p_size ;
00099 
00100   for ( EntityId i = beg_elem; i != end_elem; ++i) {
00101     element_ids_on_this_processor.push_back(i);
00102   }
00103 
00104   generate_mesh(element_ids_on_this_processor);
00105 }
00106 
00107 void QuadFixture::generate_mesh(std::vector<EntityId> & element_ids_on_this_processor)
00108 {
00109   {
00110     //sort and unique the input elements
00111     std::vector<EntityId>::iterator ib = element_ids_on_this_processor.begin();
00112     std::vector<EntityId>::iterator ie = element_ids_on_this_processor.end();
00113 
00114     std::sort( ib, ie);
00115     ib = std::unique( ib, ie);
00116     element_ids_on_this_processor.erase(ib, ie);
00117   }
00118 
00119   m_bulk_data.modification_begin();
00120 
00121   {
00122     // Declare the elements that belong on this process
00123 
00124     std::vector<EntityId>::const_iterator ib = element_ids_on_this_processor.begin();
00125     const std::vector<EntityId>::const_iterator ie = element_ids_on_this_processor.end();
00126     for (; ib != ie; ++ib) {
00127       EntityId entity_id = *ib;
00128       unsigned ix = 0, iy = 0;
00129       elem_x_y(entity_id, ix, iy);
00130 
00131       stk::mesh::EntityId elem_nodes[4] ;
00132 
00133       elem_nodes[0] = node_id( ix   , iy );
00134       elem_nodes[1] = node_id( ix+1 , iy );
00135       elem_nodes[2] = node_id( ix+1 , iy+1 );
00136       elem_nodes[3] = node_id( ix   , iy+1 );
00137 
00138       stk::mesh::fem::declare_element( m_bulk_data, m_quad_part, elem_id( ix , iy ) , elem_nodes);
00139       for (unsigned i = 0; i<4; ++i) {
00140         stk::mesh::Entity * const node = m_bulk_data.get_entity( fem::FEMMetaData::NODE_RANK , elem_nodes[i] );
00141 
00142         ThrowRequireMsg( node != NULL,
00143           "This process should know about the nodes that make up its element");
00144 
00145         // Compute and assign coordinates to the node
00146         unsigned nx = 0, ny = 0;
00147         node_x_y(elem_nodes[i], nx, ny);
00148 
00149         Scalar * data = stk::mesh::field_data( m_coord_field , *node );
00150 
00151         data[0] = nx ;
00152         data[1] = ny ;
00153       }
00154     }
00155   }
00156 
00157   m_bulk_data.modification_end();
00158 }
00159 
00160 
00161 } // fixtures
00162 } // mesh
00163 } // stk
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends