TopologyHelpers.hpp

/*------------------------------------------------------------------------*/
/*                 Copyright 2010 Sandia Corporation.                     */
/*  Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive   */
/*  license for use of this work by or on behalf of the U.S. Government.  */
/*  Export of this program may require a license from the                 */
/*  United States Government.                                             */
/*------------------------------------------------------------------------*/

#ifndef stk_mesh_TopologyHelpers_hpp
#define stk_mesh_TopologyHelpers_hpp

#include <sstream>
#include <stdexcept>
#include <Shards_CellTopologyTraits.hpp>
#include <stk_mesh/base/Types.hpp>
#include <stk_mesh/fem/FEMTypes.hpp>
#include <stk_mesh/fem/EntityRanks.hpp>
#include <stk_mesh/base/BulkData.hpp>
#include <stk_mesh/fem/TopologicalMetaData.hpp>

#ifndef SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS
#include <stk_mesh/fem/TopologyHelpersDeprecated.hpp>
#endif // SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS

namespace stk {
namespace mesh {





//----------------------------------------------------------------------
template< class Traits >
void get_parts_with_topology(stk::mesh::BulkData& mesh,
                             stk::mesh::PartVector& parts)
{
  parts.clear();

  const stk::mesh::PartVector& all_parts = mesh.mesh_meta_data().get_parts();

  stk::mesh::PartVector::const_iterator
    iter = all_parts.begin(),
    iter_end = all_parts.end();

  const CellTopologyData* topology = shards::getCellTopologyData<Traits>();

  for(; iter!=iter_end; ++iter) {
    stk::mesh::Part* part =  *iter;
#ifndef SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS
    if (get_cell_topology(*part) == topology) {
      parts.push_back(part);
    }
#else // SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS
    if (TopologicalMetaData::get_cell_topology(*part) == topology) {
      parts.push_back(part);
    }
#endif // SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS
  }
}

//----------------------------------------------------------------------
template< typename IdType >
inline
Entity & declare_element( BulkData & mesh ,
                          Part & part ,
                          const IdType elem_id ,
                          const IdType node_id[] )
{
#ifndef SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS
  const CellTopologyData * const top = get_cell_topology( part );
#else // SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS
  const CellTopologyData * const top = TopologicalMetaData::get_cell_topology( part );
#endif // SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS

  if ( top == NULL ) {
    std::ostringstream msg ;
    msg << "stk::mesh::declare_element( mesh , " ;
    msg << part.name();
    msg << " , " ;
    msg << elem_id ;
    msg << " , node_id[] ) ERROR, Part does not have a local topology" ;
    throw std::runtime_error( msg.str() );
  }

  PartVector empty ;
  PartVector add( 1 ); add[0] = & part ;

#ifndef SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS
  const EntityRank entity_rank = element_rank_deprecated(part.mesh_meta_data());
#else // SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS
  const EntityRank entity_rank = top->dimension;
#endif // SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS

  Entity & elem = mesh.declare_entity( entity_rank, elem_id, add );

  for ( unsigned i = 0 ; i < top->node_count ; ++i ) {
    //declare node if it doesn't already exist
    Entity * node = mesh.get_entity( NodeRank , node_id[i]);
    if ( NULL == node) {
      node = & mesh.declare_entity( NodeRank , node_id[i], empty );
    }

    mesh.declare_relation( elem , *node , i );
  }
  return elem ;
}

template< typename IdType >
Entity & declare_element( BulkData & mesh ,
                          Part & part ,
                          const IdType elem_id ,
                          Entity * node[] )
{
#ifndef SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS
  const CellTopologyData * const top = get_cell_topology( part );
#else // SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS
  const CellTopologyData * const top = TopologicalMetaData::get_cell_topology( part );
#endif // SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS

  if ( top == NULL ) {
    std::ostringstream msg ;
    msg << "stk::mesh::declare_element( mesh , " ;
    msg << part.name();
    msg << " , " ;
    msg << elem_id ;
    msg << " , node[] ) ERROR, Part does not have a local topology" ;
    throw std::runtime_error( msg.str() );
  }

  PartVector add( 1 ); add[0] = & part ;

#ifndef SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS
  const EntityRank entity_rank = element_rank_deprecated(part.mesh_meta_data());
#else // SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS
  const EntityRank entity_rank = top->dimension;
#endif // SKIP_DEPRECATED_STK_MESH_TOPOLOGY_HELPERS

  Entity & elem = mesh.declare_entity( entity_rank, elem_id, add );

  for ( unsigned i = 0 ; i < top->node_count ; ++i ) {
    mesh.declare_relation( elem , *node[i] , i );
  }
  return elem ;
}

//----------------------------------------------------------------------
Entity & declare_element_side( BulkData & mesh ,
                               const stk::mesh::EntityId global_side_id ,
                               Entity & elem ,
                               const unsigned local_side_id ,
                               Part * part = NULL );

Entity & declare_element_side( Entity & elem ,
                               Entity & side ,
                               const unsigned local_side_id ,
                               Part * part = NULL );

bool element_side_polarity( const Entity & elem ,
                            const Entity & side , int local_side_id = -1 );


int element_local_side_id( const Entity & elem ,
                           const CellTopologyData * side_topology,
                           const std::vector<Entity*>& side_nodes );

//----------------------------------------------------------------------

}//namespace mesh
}//namespace stk

#endif