Shards_CellTopology.hpp

/*
//@HEADER
// ************************************************************************
//
//                Shards : Shared Discretization Tools
//                 Copyright 2008, 2011 Sandia Corporation
//
// Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
// the U.S. Government retains certain rights in this software.
//
// Redistribution and use in source and binary forms, with or without
// modification, are permitted provided that the following conditions are
// met:
//
// 1. Redistributions of source code must retain the above copyright
// notice, this list of conditions and the following disclaimer.
//
// 2. Redistributions in binary form must reproduce the above copyright
// notice, this list of conditions and the following disclaimer in the
// documentation and/or other materials provided with the distribution.
//
// 3. Neither the name of the Corporation nor the names of the
// contributors may be used to endorse or promote products derived from
// this software without specific prior written permission.
//
// THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
// EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
// IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
// PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
// CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
// EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
// PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
// PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
// LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
// NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
// SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
//
// Questions? Contact Carter Edwards (hcedwar@sandia.gov),
//                    Pavel Bochev (pbboche@sandia.gov), or
//                    Denis Ridzal (dridzal@sandia.gov).
//
// ************************************************************************
//@HEADER
*/

#ifndef Shards_CellTopology_hpp
#define Shards_CellTopology_hpp

#ifdef HAVE_SHARDS_DEBUG
#define SHARDS_REQUIRE( S )  S
#else
#define SHARDS_REQUIRE( S ) /* empty */
#endif

#include <string>
#include <vector>
#include <Shards_CellTopologyData.h>
#include <Shards_BasicTopologies.hpp>

namespace shards {

/*------------------------------------------------------------------------*/

class CellTopology ;

//std::ostream & operator << ( std::ostream & , const CellTopologyData & );


std::ostream & operator << ( std::ostream & , const CellTopology & );


enum ECellType {
  ALL_CELLS = 0,
  STANDARD_CELL,
  NONSTANDARD_CELL
};

inline std::string ECellTypeToString(ECellType cellType) {
  std::string retString;
  switch(cellType){
    case ALL_CELLS:         retString = "All";           break;
    case STANDARD_CELL:     retString = "Standard";      break;
    case NONSTANDARD_CELL:  retString = "Nonstandard";   break;
    default:                retString = "Invalid Cell";
  }
  return retString;
}


enum ETopologyType {
  ALL_TOPOLOGIES,
  BASE_TOPOLOGY,
  EXTENDED_TOPOLOGY
};

inline std::string ETopologyTypeToString(ETopologyType topologyType) {
  std::string retString;
  switch(topologyType){
    case ALL_TOPOLOGIES:      retString = "All";            break;
    case BASE_TOPOLOGY:       retString = "Base";           break;
    case EXTENDED_TOPOLOGY:   retString = "Extended";       break;
    default:                  retString = "Invalid Topology";
  }
  return retString;
}


void getTopologies(std::vector<shards::CellTopology>& topologies,
                   const unsigned       cellDim   = 4,
                   const ECellType      cellType      = ALL_CELLS,
                   const ETopologyType  topologyType  = ALL_TOPOLOGIES);



int isPredefinedCell(const CellTopology &  cell);



/*------------------------------------------------------------------------*/
class CellTopology {
private:
  
  void requireCell() const ;
  
  
  void requireDimension( const unsigned subcellDim ) const ;
  
  
  void requireSubcell( const unsigned subcellDim ,
                       const unsigned subcellOrd ) const ;
  
  
  void requireNodeMap( const unsigned subcellDim ,
                       const unsigned subcellOrd ,
                       const unsigned nodeOrd ) const ;

  void requireNodePermutation( const unsigned permutationOrd ,
                               const unsigned nodeOrd ) const ;
  
  const CellTopologyData    * m_cell ;

public:

  /*------------------------------------------------------------------*/
  unsigned getDimension() const
    {
      SHARDS_REQUIRE( requireCell() );
      return m_cell->dimension ;
    }
  
        
  unsigned getKey() const
    {
      SHARDS_REQUIRE( requireCell() );
      return m_cell->key ;
    }

        
  unsigned getBaseKey() const
    {
      SHARDS_REQUIRE( requireCell() );
      return m_cell->base->key ;
    }
        
        
        
  const char* getName() const
    {
        SHARDS_REQUIRE( requireCell() );
        return m_cell->name ;
    }

        
  const char* getBaseName() const
    {
      SHARDS_REQUIRE( requireCell() );
      return m_cell->base->name ;
    }
        
        
  unsigned getNodeCount() const
    {
      SHARDS_REQUIRE( requireCell() );
      return m_cell->node_count ;
    }
  
        
  unsigned getVertexCount() const
    {
      SHARDS_REQUIRE( requireCell() );
      return m_cell->vertex_count ;
    }
  
        
  unsigned getEdgeCount() const
    {
      SHARDS_REQUIRE( requireCell() );
      return m_cell->edge_count ;
    }
  
  unsigned getFaceCount() const
    {
      SHARDS_REQUIRE( requireCell() );
      return m_cell->dimension == 3 ? m_cell->side_count : 0 ;
    }
  
        
  unsigned getSideCount() const
    {
      SHARDS_REQUIRE( requireCell() );
      return m_cell->side_count ;
    }
  
        
  bool isValid() const
    { return m_cell != 0 ; }

        
  const CellTopologyData * getCellTopologyData() const
    { return m_cell ; }

        
  const CellTopologyData * getBaseCellTopologyData() const
    {
      SHARDS_REQUIRE( requireCell() );
      return m_cell->base ;
    }

        
  const CellTopologyData * getCellTopologyData( const unsigned subcell_dim ,
                                                const unsigned subcell_ord ) const
    {
      SHARDS_REQUIRE( requireCell() );
      SHARDS_REQUIRE( requireDimension(subcell_dim) );
      SHARDS_REQUIRE( requireSubcell(subcell_dim,subcell_ord) );
      return m_cell->subcell[subcell_dim][subcell_ord].topology ;
    }

        
  const CellTopologyData * getBaseCellTopologyData( const unsigned subcell_dim ,
                                                    const unsigned subcell_ord ) const
    {
      return getCellTopologyData(subcell_dim,subcell_ord)->base ;
    }

        
  unsigned getKey( const unsigned subcell_dim ,
                   const unsigned subcell_ord ) const
    {
      return getCellTopologyData(subcell_dim,subcell_ord)->key ;
    }


        
  const char * getName(const unsigned subcell_dim,   
                       const unsigned subcell_ord) const
    {
      return getCellTopologyData(subcell_dim,subcell_ord) -> name;
    }
        
        
  unsigned getNodeCount( const unsigned subcell_dim ,
                         const unsigned subcell_ord ) const
    {
      return getCellTopologyData(subcell_dim,subcell_ord)->node_count ;
    }

        
  unsigned getVertexCount( const unsigned subcell_dim ,
                           const unsigned subcell_ord ) const
    {
      return getCellTopologyData(subcell_dim,subcell_ord)->vertex_count ;
    }

        
  unsigned getEdgeCount( const unsigned subcell_dim ,
                         const unsigned subcell_ord ) const
    {
      return getCellTopologyData(subcell_dim,subcell_ord)->edge_count ;
    }
  
        
  unsigned getSideCount( const unsigned subcell_dim ,
                         const unsigned subcell_ord ) const
    {
      return getCellTopologyData(subcell_dim,subcell_ord)->side_count ;
    }

        
  unsigned getSubcellCount( const unsigned subcell_dim ) const
    {
      SHARDS_REQUIRE( requireCell() );
      SHARDS_REQUIRE( requireDimension(subcell_dim) );
      return m_cell->subcell_count[subcell_dim] ;
    }
  
        
  bool getSubcellHomogeneity( const unsigned subcell_dim ) const
    {
      SHARDS_REQUIRE( requireCell() );
      SHARDS_REQUIRE( requireDimension(subcell_dim) );
      return 0 != m_cell->subcell_homogeneity[subcell_dim] ;
    }
  
        
  unsigned getNodeMap( const unsigned  subcell_dim ,
                       const unsigned  subcell_ord ,
                       const unsigned  subcell_node_ord ) const
    {
      SHARDS_REQUIRE( requireCell() );
      SHARDS_REQUIRE( requireDimension(subcell_dim) );
      SHARDS_REQUIRE( requireSubcell(subcell_dim,subcell_ord) );
      SHARDS_REQUIRE( requireNodeMap(subcell_dim,subcell_ord,subcell_node_ord));
      return m_cell->subcell[subcell_dim][subcell_ord].node[subcell_node_ord];
    }


  unsigned getNodePermutationCount() const
    {
      SHARDS_REQUIRE(requireCell());
      return m_cell->permutation_count ;
    }

  unsigned getNodePermutation( const unsigned permutation_ord ,
                               const unsigned node_ord ) const
    {
      SHARDS_REQUIRE(requireCell());
      SHARDS_REQUIRE(requireNodePermutation(permutation_ord,node_ord));
      return m_cell->permutation[permutation_ord].node[node_ord];
    }
  
  unsigned getNodePermutationPolarity( const unsigned permutation_ord ) const
    {
      SHARDS_REQUIRE(requireCell());
      SHARDS_REQUIRE(requireNodePermutation(permutation_ord,0));
      return m_cell->permutation[permutation_ord].polarity;
    }
  
  unsigned getNodePermutationInverse( const unsigned permutation_ord ,
                                      const unsigned node_ord ) const
    {
      SHARDS_REQUIRE(requireCell());
      SHARDS_REQUIRE(requireNodePermutation(permutation_ord,node_ord));
      return m_cell->permutation_inverse[permutation_ord].node[node_ord];
    }
  
  /*------------------------------------------------------------------*/
  CellTopology( const CellTopologyData * cell )
    : m_cell( cell )
  {}
  
        
  CellTopology( const std::string & name,
                const unsigned      nodeCount);
  
        
  CellTopology( const std::string                             & name,
                const unsigned                                  vertex_count,
                const unsigned                                  node_count,
                const std::vector< const CellTopologyData * > & edges ,
                const std::vector< unsigned >                 & edge_node_map ,
                const CellTopologyData                        * base = NULL );

  
  CellTopology( const std::string                             & name,
                const unsigned                                  vertex_count,
                const unsigned                                  node_count,
                const std::vector< const CellTopologyData * > & edges ,
                const std::vector< unsigned >                 & edge_node_map ,
                const std::vector< const CellTopologyData * > & faces ,
                const std::vector< unsigned >                 & face_node_map ,
                const CellTopologyData                        * base = NULL );
  
  
  CellTopology& operator = (const CellTopology& right);

  CellTopology( const CellTopology& right );

  CellTopology();

  ~CellTopology();
  
}; // class CellTopology

/*------------------------------------------------------------------------*/
/* \brief  Find the permutation from the expected nodes to the actual nodes,
 *
 *  Find permutation 'p' such that:
 *    actual_node[j] == expected_node[ top.permutation[p].node[j] ]
 *  for all vertices.
 */
template< typename id_type >
int findPermutation( const CellTopologyData & top ,
                     const id_type * const expected_node ,
                     const id_type * const actual_node )
{
  const int nv = top.vertex_count ;
  const int np = top.permutation_count ;
  int p = 0 ;
  for ( ; p < np ; ++p ) {
    const unsigned * const perm_node = top.permutation[p].node ;
    int j = 0 ;
    for ( ; j < nv && actual_node[j] == expected_node[ perm_node[j] ] ; ++j );
    if ( nv == j ) break ;
  }
  if ( np == p ) p = -1 ;
  return p ;
}

template< typename id_type >
int findPermutation( const CellTopology & top ,
                     const id_type * const expected_node ,
                     const id_type * const actual_node )
{
  return findPermutation( * top.getCellTopologyData() , expected_node , actual_node );
}

/*------------------------------------------------------------------------*/
void badCellTopologyKey( const unsigned dimension ,
                         const unsigned face_count ,
                         const unsigned edge_count ,
                         const unsigned vertex_count ,
                         const unsigned node_count );
  

inline
unsigned cellTopologyKey( const unsigned dimension ,
                          const unsigned face_count ,
                          const unsigned edge_count ,
                          const unsigned vertex_count ,
                          const unsigned node_count )
{
  const bool bad = ( dimension    >> 3 ) ||
                   ( face_count   >> 6 ) ||
                   ( edge_count   >> 6 ) ||
                   ( vertex_count >> 6 ) ||
                   ( node_count   >> 10 );

  if ( bad ) {
    badCellTopologyKey( dimension , 
                        face_count ,
                        edge_count , 
                        vertex_count , 
                        node_count );
  }

  const unsigned key = ( dimension    << 28  ) |
                       ( face_count   << 22  ) |
                       ( edge_count   << 16  ) |
                       ( vertex_count << 10  ) |
                       ( node_count          ) ;

  return key ;
}

inline
bool operator==(const CellTopology &left, const CellTopology &right)
{
  return left.getCellTopologyData() == right.getCellTopologyData();
}

inline
bool operator<(const CellTopology &left, const CellTopology &right)
{
  return left.getCellTopologyData() < right.getCellTopologyData();
// FIXME: Should is be this?  
//  return left.getKey() < right.getKey(); 
}

inline
bool operator!=(const CellTopology &left, const CellTopology &right) {
  return !(left == right);
}


} // namespace shards

#undef SHARDS_REQUIRE

#endif // Shards_CellTopology_hpp