Panzer_UniqueGlobalIndexer.hpp

Go to the documentation of this file.

// @HEADER
// ***********************************************************************
//
//           Panzer: A partial differential equation assembly
//       engine for strongly coupled complex multiphysics systems
//                 Copyright (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 Roger P. Pawlowski (rppawlo@sandia.gov) and
// Eric C. Cyr (eccyr@sandia.gov)
// ***********************************************************************
// @HEADER

#ifndef __Panzer_UniqueGlobalIndexer_hpp__
#define __Panzer_UniqueGlobalIndexer_hpp__

#include <vector>
#include <string>

#include <unordered_map> // a hash table for buildLocalIds()

#include "Teuchos_RCP.hpp"
#include "Teuchos_Comm.hpp"

#include "Phalanx_KokkosDeviceTypes.hpp"

namespace panzer {

// Forward declaration.
template <typename LocalOrdinalT> class ConnManagerBase;

class UniqueGlobalIndexerBase {
public:
   virtual ~UniqueGlobalIndexerBase() = 0;

   virtual Teuchos::RCP<Teuchos::Comm<int> > getComm() const = 0;

   virtual int getNumFields() const = 0;

   virtual int getFieldNum(const std::string & str) const = 0;

   virtual void getFieldOrder(std::vector<std::string> & fieldOrder) const = 0;

   virtual const std::string & getFieldString(int num) const = 0;

   virtual void getElementBlockIds(std::vector<std::string> & elementBlockIds) const = 0; 

   virtual bool fieldInBlock(const std::string & field, const std::string & block) const = 0;

   virtual const std::vector<int> & getBlockFieldNumbers(const std::string & blockId) const = 0;

   virtual const std::vector<int> & getGIDFieldOffsets(const std::string & blockId,int fieldNum) const = 0;

   virtual const std::pair<std::vector<int>,std::vector<int> > & 
   getGIDFieldOffsets_closure(const std::string & blockId, int fieldNum,
                              int subcellDim,int subcellId) const = 0;

   virtual int getElementBlockGIDCount(const std::string & blockId) const = 0;

   virtual int getElementBlockGIDCount(const std::size_t & blockIndex) const = 0;
};

template <typename LocalOrdinalT,typename GlobalOrdinalT>
class UniqueGlobalIndexer : public UniqueGlobalIndexerBase {
public:
   virtual ~UniqueGlobalIndexer() = 0;

   virtual Teuchos::RCP<Teuchos::Comm<int> > getComm() const = 0;

   virtual int getNumFields() const = 0;

   virtual int getFieldNum(const std::string & str) const = 0;

   virtual void getFieldOrder(std::vector<std::string> & fieldOrder) const = 0;

   virtual void getElementBlockIds(std::vector<std::string> & elementBlockIds) const = 0; 

   virtual bool fieldInBlock(const std::string & field, const std::string & block) const = 0;

   virtual const std::vector<int> & getBlockFieldNumbers(const std::string & blockId) const = 0;

   virtual const std::vector<int> & getGIDFieldOffsets(const std::string & blockId,int fieldNum) const = 0;

   virtual const std::pair<std::vector<int>,std::vector<int> > & 
   getGIDFieldOffsets_closure(const std::string & blockId, int fieldNum,
                              int subcellDim,int subcellId) const = 0;

   // Methods requiring Local or Global OrdinalT

   virtual void getElementOrientation(LocalOrdinalT localElmtId,std::vector<double> & gidsOrientation) const = 0;

   virtual const std::vector<LocalOrdinalT> & getElementBlock(const std::string & blockId) const = 0;

   virtual void getElementGIDs(LocalOrdinalT localElmtId,std::vector<GlobalOrdinalT> & gids,const std::string & blockIdHint="") const = 0;

   virtual void
   getOwnedIndices(
     std::vector<GlobalOrdinalT>& indices) const = 0;

   virtual void
   getGhostedIndices(
     std::vector<GlobalOrdinalT>& indices) const = 0;

   virtual void
   getOwnedAndGhostedIndices(
     std::vector<GlobalOrdinalT>& indices) const = 0;

   virtual int
   getNumOwned() const = 0;

   virtual int
   getNumGhosted() const = 0;

   virtual int
   getNumOwnedAndGhosted() const = 0;

   virtual void ownedIndices(const std::vector<GlobalOrdinalT> & indices,std::vector<bool> & isOwned) const = 0;

  const Kokkos::View<const LocalOrdinalT*,Kokkos::LayoutRight,PHX::Device> getElementLIDs(LocalOrdinalT localElmtId) const
    { return Kokkos::subview(localIDs_k_, localElmtId, Kokkos::ALL() ); }

  const Kokkos::View<const LocalOrdinalT**,Kokkos::LayoutRight,PHX::Device> getLIDs() const
    {return localIDs_k_;}

   void getElementLIDs(Kokkos::View<const int*,PHX::Device> cellIds,
                       Kokkos::View<LocalOrdinalT**,PHX::Device> lids) const
   { 
     CopyCellLIDsFunctor functor;
     functor.cellIds = cellIds;
     functor.global_lids = localIDs_k_;
     functor.local_lids = lids; // we assume this array is sized correctly!

     Kokkos::parallel_for(cellIds.extent(0),functor);
   }

   virtual int getElementBlockGIDCount(const std::string & blockId) const = 0;

   virtual int getElementBlockGIDCount(const std::size_t & blockIndex) const = 0;

   virtual Teuchos::RCP<const ConnManagerBase<LocalOrdinalT> > getConnManagerBase() const = 0;

   class CopyCellLIDsFunctor {
   public:
     typedef typename PHX::Device execution_space;

     Kokkos::View<const int*,PHX::Device> cellIds;
     Kokkos::View<const LocalOrdinalT**,Kokkos::LayoutRight,PHX::Device> global_lids;
     Kokkos::View<LocalOrdinalT**,PHX::Device> local_lids;

     KOKKOS_INLINE_FUNCTION
     void operator()(const int cell) const
     {
       for(int i=0;i<static_cast<int>(local_lids.extent(1));i++) 
         local_lids(cell,i) = global_lids(cellIds(cell),i);
     }
     
   };

protected:

   void buildLocalIds()
   { 
     // this method is implmented as two steps to ensure
     // that setLocalIds works, it would be better to simply
     // call:
     //   buildLocalIdsFromOwnedElements(localIDs_); 

     std::vector<std::vector<LocalOrdinalT> > localIDs; 
     buildLocalIdsFromOwnedElements(localIDs); 
     setLocalIds(localIDs);
   }

   void buildLocalIdsFromOwnedElements(std::vector<std::vector<LocalOrdinalT> > & localIDs) const ; 

   void setLocalIds(const std::vector<std::vector<LocalOrdinalT> > & localIDs)
   {  
     // determine the maximium second dimension of the local IDs
     std::size_t max = 0;
     for(std::size_t i=0;i<localIDs.size();i++)
       max = localIDs[i].size() > max ? localIDs[i].size() : max;

     // allocate for the kokkos size
     Kokkos::View<LocalOrdinalT**,Kokkos::LayoutRight,PHX::Device> localIDs_k 
       = Kokkos::View<LocalOrdinalT**,Kokkos::LayoutRight,PHX::Device>("ugi:localIDs_",localIDs.size(),max);
     for(std::size_t i=0;i<localIDs.size();i++) {
       for(std::size_t j=0;j<localIDs[i].size();j++)
         localIDs_k(i,j) = localIDs[i][j];
     }

     // store in Kokkos type
     localIDs_k_ = localIDs_k;
   }

   void shareLocalIDs(const UniqueGlobalIndexer<LocalOrdinalT,GlobalOrdinalT> & src)
   {
     localIDs_k_ = src.localIDs_k_;
   }

private:
  Kokkos::View<const LocalOrdinalT**,Kokkos::LayoutRight,PHX::Device> localIDs_k_;
};

// prevents a warning because a destructor does not exist
inline UniqueGlobalIndexerBase::~UniqueGlobalIndexerBase() {}

// prevents a warning because a destructor does not exist
template <typename LocalOrdinalT,typename GlobalOrdinalT>
inline UniqueGlobalIndexer<LocalOrdinalT,GlobalOrdinalT>::~UniqueGlobalIndexer() {}

template <typename LocalOrdinalT,typename GlobalOrdinalT>
inline void UniqueGlobalIndexer<LocalOrdinalT,GlobalOrdinalT>::
buildLocalIdsFromOwnedElements(std::vector<std::vector<LocalOrdinalT> > & localIDs) const
{
  std::vector<GlobalOrdinalT> ownedAndGhosted;
  this->getOwnedAndGhostedIndices(ownedAndGhosted);
   
  // build global to local hash map (temporary and used only once)
  std::unordered_map<GlobalOrdinalT,LocalOrdinalT> hashMap;
  for(std::size_t i=0;i<ownedAndGhosted.size();i++)
    hashMap[ownedAndGhosted[i]] = i;

  std::vector<std::string> elementBlocks;
  this->getElementBlockIds(elementBlocks);
 
  // compute total number of elements
  std::size_t numElmts = 0;
  for(std::size_t eb=0;eb<elementBlocks.size();eb++)
    numElmts += this->getElementBlock(elementBlocks[eb]).size();
  localIDs.resize(numElmts); // allocate local ids

  // perform computation of local ids
  for(std::size_t eb=0;eb<elementBlocks.size();eb++) {
    std::vector<GlobalOrdinalT> gids;
    const std::vector<LocalOrdinalT> & elmts = this->getElementBlock(elementBlocks[eb]);

    for(std::size_t e=0;e<elmts.size();e++) {
      this->getElementGIDs(elmts[e],gids,elementBlocks[eb]);
      std::vector<LocalOrdinalT> & lids = localIDs[elmts[e]];
      lids.resize(gids.size());
 
      for(std::size_t g=0;g<gids.size();g++)
        lids[g] = hashMap[gids[g]];
    }
  } 
}

}

#endif