Xpetra_StridedMap.hpp

Go to the documentation of this file.

// @HEADER
//
// ***********************************************************************
//
//             Xpetra: A linear algebra interface package
//                  Copyright 2012 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
//                    Jonathan Hu       (jhu@sandia.gov)
//                    Andrey Prokopenko (aprokop@sandia.gov)
//                    Ray Tuminaro      (rstumin@sandia.gov)
//
// ***********************************************************************
//
// @HEADER
 
// WARNING: This code is experimental. Backwards compatibility should not be expected.
 
#ifndef XPETRA_STRIDEDMAP_HPP
#define XPETRA_STRIDEDMAP_HPP
 
/* this file is automatically generated - do not edit (see script/interfaces.py) */
 
#include <Kokkos_DefaultNode.hpp>
 
#include <Teuchos_Describable.hpp>
#include <Teuchos_OrdinalTraits.hpp>
 
#include "Xpetra_ConfigDefs.hpp"
#include "Xpetra_Exceptions.hpp"
 
#include "Xpetra_Map.hpp"
#include "Xpetra_MapFactory.hpp"
 
namespace Xpetra {
 
  template <class LocalOrdinal = Map<>::local_ordinal_type,
            class GlobalOrdinal = typename Map<LocalOrdinal>::global_ordinal_type,
            class Node = typename Map<LocalOrdinal,GlobalOrdinal>::node_type>
  class StridedMap : public virtual Map<LocalOrdinal, GlobalOrdinal, Node> {
  public:
    typedef LocalOrdinal local_ordinal_type;
    typedef GlobalOrdinal global_ordinal_type;
    typedef Node node_type;
 
  private:
 
  typedef Xpetra::MapFactory<LocalOrdinal, GlobalOrdinal, Node> MapFactory_t;
#undef XPETRA_STRIDEDMAP_SHORT
#include "Xpetra_UseShortNamesOrdinal.hpp"
 
  public:
 
 
 
    StridedMap (UnderlyingLib xlib,
                global_size_t numGlobalElements,
                GlobalOrdinal indexBase,
                std::vector<size_t>& stridingInfo,
                const Teuchos::RCP< const Teuchos::Comm< int > >& comm,
                LocalOrdinal stridedBlockId = -1,  // FIXME (mfh 03 Sep 2014) This breaks for unsigned LocalOrdinal
                GlobalOrdinal offset = 0,
                LocalGlobal lg = GloballyDistributed,
                const Teuchos::RCP< Node >& node = Teuchos::rcp(new Node))
      : stridingInfo_ (stridingInfo),
        stridedBlockId_ (stridedBlockId),
        offset_ (offset),
        indexBase_ (indexBase)
    {
      size_t blkSize = getFixedBlockSize ();
      TEUCHOS_TEST_FOR_EXCEPTION(
        stridingInfo.size() == 0, Exceptions::RuntimeError,
        "StridedMap::StridedMap: stridingInfo not valid: stridingInfo.size() = 0?");
      TEUCHOS_TEST_FOR_EXCEPTION(
        numGlobalElements == Teuchos::OrdinalTraits<global_size_t>::invalid (),
        std::invalid_argument,
        "StridedMap::StridedMap: numGlobalElements is invalid");
      TEUCHOS_TEST_FOR_EXCEPTION(
        numGlobalElements % blkSize != 0, Exceptions::RuntimeError,
        "StridedMap::StridedMap: stridingInfo not valid: getFixedBlockSize "
        "is not an integer multiple of numGlobalElements.");
      if (stridedBlockId != -1)
        TEUCHOS_TEST_FOR_EXCEPTION(
          stridingInfo.size() < static_cast<size_t> (stridedBlockId),
          Exceptions::RuntimeError, "StridedTpetraMap::StridedTpetraMap: "
          "stridedBlockId > stridingInfo.size()");
 
      // Try to create a shortcut
      if (blkSize != 1 || offset_ != 0) {
        // check input data and reorganize map
        global_size_t numGlobalNodes = numGlobalElements / blkSize;
 
        // build an equally distributed node map
        RCP<Map> nodeMap = MapFactory_t::Build(xlib, numGlobalNodes, indexBase, comm, lg, node);
        global_size_t numLocalNodes = nodeMap->getNodeNumElements();
 
        // translate local node ids to local dofs
        size_t nStridedOffset = 0;
        size_t nDofsPerNode = blkSize; // dofs per node for local striding block
        if (stridedBlockId > -1) {
          for (int j = 0; j < stridedBlockId; j++)
            nStridedOffset += stridingInfo_[j];
 
          nDofsPerNode = stridingInfo_[stridedBlockId];
          numGlobalElements = numGlobalNodes * Teuchos::as<global_size_t>(nDofsPerNode);
        }
        size_t numLocalElements = numLocalNodes * Teuchos::as<size_t>(nDofsPerNode);
 
        std::vector<GlobalOrdinal> dofgids(numLocalElements);
        for (LocalOrdinal i = 0; i < Teuchos::as<LocalOrdinal>(numLocalNodes); i++) {
          GlobalOrdinal nodeGID = nodeMap->getGlobalElement(i);
 
          for (size_t j = 0; j < nDofsPerNode; j++)
            dofgids[i*nDofsPerNode + j] = indexBase_ + offset_ + (nodeGID - indexBase_)*Teuchos::as<GlobalOrdinal>(blkSize) + Teuchos::as<GlobalOrdinal>(nStridedOffset + j);
        }
 
        map_ = MapFactory_t::Build(xlib, numGlobalElements, dofgids, indexBase, comm, node);
 
        if (stridedBlockId == -1) {
          TEUCHOS_TEST_FOR_EXCEPTION(getNodeNumElements() != Teuchos::as<size_t>(nodeMap->getNodeNumElements()*nDofsPerNode), Exceptions::RuntimeError,
                                     "StridedTpetraMap::StridedTpetraMap: wrong distribution of dofs among processors.");
          TEUCHOS_TEST_FOR_EXCEPTION(getGlobalNumElements() != Teuchos::as<size_t>(nodeMap->getGlobalNumElements()*nDofsPerNode), Exceptions::RuntimeError,
                                     "StridedTpetraMap::StridedTpetraMap: wrong distribution of dofs among processors.");
 
        } else {
          size_t nDofsInStridedBlock = stridingInfo[stridedBlockId];
          TEUCHOS_TEST_FOR_EXCEPTION(getNodeNumElements() != Teuchos::as<size_t>(nodeMap->getNodeNumElements()*nDofsInStridedBlock), Exceptions::RuntimeError,
                                     "StridedTpetraMap::StridedTpetraMap: wrong distribution of dofs among processors.");
          TEUCHOS_TEST_FOR_EXCEPTION(getGlobalNumElements() != Teuchos::as<size_t>(nodeMap->getGlobalNumElements()*nDofsInStridedBlock), Exceptions::RuntimeError,
                                     "StridedTpetraMap::StridedTpetraMap: wrong distribution of dofs among processors.");
        }
      } else {
        map_ = MapFactory_t::Build(xlib, numGlobalElements, indexBase, comm, lg, node);
      }
 
      TEUCHOS_TEST_FOR_EXCEPTION(CheckConsistency() == false, Exceptions::RuntimeError, "StridedTpetraMap::StridedTpetraMap: CheckConsistency() == false");
    }
 
 
    StridedMap(UnderlyingLib xlib, global_size_t numGlobalElements, size_t numLocalElements, GlobalOrdinal indexBase, std::vector<size_t>& stridingInfo,
               const Teuchos::RCP< const Teuchos::Comm< int > > &comm, LocalOrdinal stridedBlockId = -1, GlobalOrdinal offset = 0,
               const Teuchos::RCP< Node > &node = Teuchos::rcp(new Node))
    : stridingInfo_(stridingInfo), stridedBlockId_(stridedBlockId), offset_(offset), indexBase_(indexBase)
    {
      size_t blkSize = getFixedBlockSize();
      TEUCHOS_TEST_FOR_EXCEPTION(stridingInfo.size() == 0, Exceptions::RuntimeError,
                                 "StridedMap::StridedMap: stridingInfo not valid: stridingInfo.size() = 0?");
      if (numGlobalElements != Teuchos::OrdinalTraits<global_size_t>::invalid()) {
        TEUCHOS_TEST_FOR_EXCEPTION(numGlobalElements % blkSize != 0, Exceptions::RuntimeError,
                                   "StridedMap::StridedMap: stridingInfo not valid: getFixedBlockSize is not an integer multiple of numGlobalElements.");
#ifdef HAVE_XPETRA_DEBUG
        // We have to do this check ourselves, as we don't necessarily construct the full Tpetra map
        global_size_t sumLocalElements;
        Teuchos::reduceAll(*comm, Teuchos::REDUCE_SUM, Teuchos::as<global_size_t>(numLocalElements), Teuchos::outArg(sumLocalElements));
        TEUCHOS_TEST_FOR_EXCEPTION(sumLocalElements != numGlobalElements, std::invalid_argument,
                                   "StridedMap::StridedMap: sum of numbers of local elements is different from the provided number of global elements.");
#endif
      }
      TEUCHOS_TEST_FOR_EXCEPTION(numLocalElements % blkSize != 0, Exceptions::RuntimeError,
                                 "StridedMap::StridedMap: stridingInfo not valid: getFixedBlockSize is not an integer multiple of numLocalElements.");
      if (stridedBlockId != -1)
        TEUCHOS_TEST_FOR_EXCEPTION(stridingInfo.size() < Teuchos::as<size_t>(stridedBlockId), Exceptions::RuntimeError,
                                   "StridedTpetraMap::StridedTpetraMap: stridedBlockId > stridingInfo.size()");
 
      // Try to create a shortcut
      if (blkSize != 1 || offset_ != 0) {
        // check input data and reorganize map
        global_size_t numGlobalNodes = Teuchos::OrdinalTraits<global_size_t>::invalid();
        if (numGlobalElements != Teuchos::OrdinalTraits<global_size_t>::invalid())
          numGlobalNodes = numGlobalElements / blkSize;
        global_size_t numLocalNodes = numLocalElements / blkSize;
 
        // build an equally distributed node map
        RCP<Map> nodeMap = MapFactory_t::Build(xlib, numGlobalNodes, numLocalNodes, indexBase, comm, node);
 
        // translate local node ids to local dofs
        size_t nStridedOffset = 0;
        size_t nDofsPerNode = blkSize; // dofs per node for local striding block
        if (stridedBlockId > -1) {
          for (int j = 0; j < stridedBlockId; j++)
            nStridedOffset += stridingInfo_[j];
 
          nDofsPerNode = stridingInfo_[stridedBlockId];
          numGlobalElements = nodeMap->getGlobalNumElements() * Teuchos::as<global_size_t>(nDofsPerNode);
        }
        numLocalElements = numLocalNodes * Teuchos::as<size_t>(nDofsPerNode);
 
        std::vector<GlobalOrdinal> dofgids(numLocalElements);
        for (LocalOrdinal i = 0; i < Teuchos::as<LocalOrdinal>(numLocalNodes); i++) {
          GlobalOrdinal nodeGID = nodeMap->getGlobalElement(i);
 
          for (size_t j = 0; j < nDofsPerNode; j++)
            dofgids[i*nDofsPerNode + j] = indexBase_ + offset_ + (nodeGID - indexBase_)*Teuchos::as<GlobalOrdinal>(blkSize) + Teuchos::as<GlobalOrdinal>(nStridedOffset + j);
        }
 
        map_ = MapFactory_t::Build(xlib, numGlobalElements, dofgids, indexBase, comm, node);
 
        if (stridedBlockId == -1) {
          TEUCHOS_TEST_FOR_EXCEPTION(getNodeNumElements() != Teuchos::as<size_t>(nodeMap->getNodeNumElements()*nDofsPerNode), Exceptions::RuntimeError,
                                     "StridedTpetraMap::StridedTpetraMap: wrong distribution of dofs among processors.");
          TEUCHOS_TEST_FOR_EXCEPTION(getGlobalNumElements() != Teuchos::as<size_t>(nodeMap->getGlobalNumElements()*nDofsPerNode), Exceptions::RuntimeError,
                                     "StridedTpetraMap::StridedTpetraMap: wrong distribution of dofs among processors.");
 
        } else {
          int nDofsInStridedBlock = stridingInfo[stridedBlockId];
          TEUCHOS_TEST_FOR_EXCEPTION(getNodeNumElements() != Teuchos::as<size_t>(nodeMap->getNodeNumElements()*nDofsInStridedBlock), Exceptions::RuntimeError,
                                     "StridedTpetraMap::StridedTpetraMap: wrong distribution of dofs among processors.");
          TEUCHOS_TEST_FOR_EXCEPTION(getGlobalNumElements() != Teuchos::as<size_t>(nodeMap->getGlobalNumElements()*nDofsInStridedBlock), Exceptions::RuntimeError,
                                     "StridedTpetraMap::StridedTpetraMap: wrong distribution of dofs among processors.");
        }
 
      } else {
        map_ = MapFactory_t::Build(xlib, numGlobalElements, numLocalElements, indexBase, comm, node);
      }
 
      TEUCHOS_TEST_FOR_EXCEPTION(CheckConsistency() == false, Exceptions::RuntimeError, "StridedTpetraMap::StridedTpetraMap: CheckConsistency() == false");
    }
 
    StridedMap(UnderlyingLib xlib, global_size_t numGlobalElements, const Teuchos::ArrayView< const GlobalOrdinal > &elementList, GlobalOrdinal indexBase,
               std::vector<size_t>& stridingInfo, const Teuchos::RCP< const Teuchos::Comm< int > > &comm, LocalOrdinal stridedBlockId = -1,
               const Teuchos::RCP< Node > &node = Teuchos::rcp(new Node))
    : stridingInfo_(stridingInfo), stridedBlockId_(stridedBlockId), indexBase_(indexBase)
    {
      size_t blkSize = getFixedBlockSize();
 
      TEUCHOS_TEST_FOR_EXCEPTION(stridingInfo.size() == 0, Exceptions::RuntimeError,
                                 "StridedMap::StridedMap: stridingInfo not valid: stridingInfo.size() = 0?");
      if (stridedBlockId != -1)
        TEUCHOS_TEST_FOR_EXCEPTION(stridingInfo.size() < Teuchos::as<size_t>(stridedBlockId), Exceptions::RuntimeError,
                                   "StridedTpetraMap::StridedTpetraMap: stridedBlockId > stridingInfo.size()");
      if (numGlobalElements != Teuchos::OrdinalTraits<global_size_t>::invalid()) {
        TEUCHOS_TEST_FOR_EXCEPTION(numGlobalElements % blkSize != 0, Exceptions::RuntimeError,
                                   "StridedMap::StridedMap: stridingInfo not valid: getFixedBlockSize is not an integer multiple of numGlobalElements.");
#ifdef HAVE_XPETRA_DEBUG
        // We have to do this check ourselves, as we don't necessarily construct the full Tpetra map
        global_size_t sumLocalElements, numLocalElements = elementList.size();
        Teuchos::reduceAll(*comm, Teuchos::REDUCE_SUM, numLocalElements, Teuchos::outArg(sumLocalElements));
        TEUCHOS_TEST_FOR_EXCEPTION(sumLocalElements != numGlobalElements, std::invalid_argument,
                                   "StridedMap::StridedMap: sum of numbers of local elements is different from the provided number of global elements.");
#endif
      }
 
      if (stridedBlockId == -1) {
        // numGlobalElements can be -1! FIXME
        // TEUCHOS_TEST_FOR_EXCEPTION(numGlobalElements  % blkSize != 0, Exceptions::RuntimeError,
                                   // "StridedMap::StridedMap: stridingInfo not valid: getFixedBlockSize is not an integer multiple of numGlobalElements.");
        TEUCHOS_TEST_FOR_EXCEPTION(elementList.size() % blkSize != 0, Exceptions::RuntimeError,
                                   "StridedMap::StridedMap: stridingInfo not valid: getFixedBlockSize is not an integer multiple of elementList.size().");
 
      } else {
        // numGlobalElements can be -1! FIXME
        // TEUCHOS_TEST_FOR_EXCEPTION(numGlobalElements  % stridingInfo[stridedBlockId] != 0, Exceptions::RuntimeError,
                                   // "StridedMap::StridedMap: stridingInfo not valid: stridingBlockInfo[stridedBlockId] is not an integer multiple of numGlobalElements.");
        TEUCHOS_TEST_FOR_EXCEPTION(elementList.size() % stridingInfo[stridedBlockId] != 0, Exceptions::RuntimeError,
                                   "StridedMap::StridedMap: stridingInfo not valid: stridingBlockInfo[stridedBlockId] is not an integer multiple of elementList.size().");
      }
 
      map_ = MapFactory_t::Build(xlib, numGlobalElements, elementList, indexBase, comm, node);
 
      // calculate offset_
 
      // find minimum GID over all procs
      GlobalOrdinal minGidOnCurProc = Teuchos::OrdinalTraits<GlobalOrdinal>::max();
      for (Teuchos_Ordinal k = 0; k < elementList.size(); k++) // TODO fix occurence of Teuchos_Ordinal
        if (elementList[k] < minGidOnCurProc)
          minGidOnCurProc = elementList[k];
 
      Teuchos::reduceAll(*comm, Teuchos::REDUCE_MIN, minGidOnCurProc, Teuchos::outArg(offset_));
 
      // calculate striding index
      size_t nStridedOffset = 0;
      for (int j = 0; j < stridedBlockId; j++)
        nStridedOffset += stridingInfo[j];
      const GlobalOrdinal goStridedOffset = Teuchos::as<GlobalOrdinal>(nStridedOffset);
 
      // adapt offset_
      offset_ -= goStridedOffset + indexBase_;
 
      TEUCHOS_TEST_FOR_EXCEPTION(CheckConsistency() == false, Exceptions::RuntimeError, "StridedTpetraMap::StridedTpetraMap: CheckConsistency() == false");
    }
 
    StridedMap(const RCP<const Map>& map, std::vector<size_t>& stridingInfo, GlobalOrdinal indexBase, LocalOrdinal stridedBlockId = -1, GlobalOrdinal offset = 0)
    : stridingInfo_(stridingInfo), stridedBlockId_(stridedBlockId), offset_(offset), indexBase_(map->getIndexBase())
    {
      // TAW: 11/24/15
      //      A strided map never can be built from a strided map. getMap always returns the underlying
      //      Xpetra::Map object which contains the data (either in a Xpetra::EpetraMapT or Xpetra::TpetraMap
      //      object)
      if(Teuchos::rcp_dynamic_cast<const StridedMap>(map) == Teuchos::null)
        map_ = map; // if map is not a strided map, just store it (standard case)
      else
        map_ = map->getMap(); // if map is also a strided map, store the underlying plain Epetra/Tpetra Xpetra map object
    }
 
 
    virtual ~StridedMap() { }
 
 
 
 
    std::vector<size_t> getStridingData() const             { return stridingInfo_; }
 
    void setStridingData(std::vector<size_t> stridingInfo)  { stridingInfo_ = stridingInfo; }
 
    size_t getFixedBlockSize() const {
      size_t blkSize = 0;
      for (std::vector<size_t>::const_iterator it = stridingInfo_.begin(); it != stridingInfo_.end(); ++it)
        blkSize += *it;
      return blkSize;
    }
 
    LocalOrdinal getStridedBlockId() const                  { return stridedBlockId_; }
 
    bool isStrided() const                                  { return stridingInfo_.size() > 1 ? true : false; }
 
    bool isBlocked() const                                  { return getFixedBlockSize() > 1 ? true : false; }
 
    GlobalOrdinal getOffset() const                         { return offset_; }
 
    void setOffset(GlobalOrdinal offset)                    { offset_ = offset; }
 
    // returns number of strided block id which gid belongs to.
    size_t GID2StridingBlockId(GlobalOrdinal gid) const {
      GlobalOrdinal tgid = gid - offset_ - indexBase_;
      tgid = tgid % getFixedBlockSize();
 
      size_t nStridedOffset = 0;
      size_t stridedBlockId = 0;
      for (size_t j = 0; j < stridingInfo_.size(); j++) {
        nStridedOffset += stridingInfo_[j];
        if (Teuchos::as<size_t>(tgid) < nStridedOffset) {
          stridedBlockId = j;
          break;
        }
      }
      return stridedBlockId;
    }
 
 
 
    RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> > getMap() const { return map_; }
 
#ifdef HAVE_XPETRA_KOKKOS_REFACTOR
#ifdef HAVE_XPETRA_TPETRA
    using local_map_type = typename Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node>::local_map_type;
    local_map_type getLocalMap () const {
      return map_->getLocalMap();
    }
#else
#ifdef __GNUC__
#warning "Xpetra Kokkos interface for CrsMatrix is enabled (HAVE_XPETRA_KOKKOS_REFACTOR) but Tpetra is disabled. The Kokkos interface needs Tpetra to be enabled, too."
#endif
#endif
#endif
 
 
    /* // function currently not needed but maybe useful
    std::vector<GlobalOrdinal> NodeId2GlobalDofIds(GlobalOrdinal nodeId) const {
      TEUCHOS_TEST_FOR_EXCEPTION(stridingInfo_.size() == 0, Exceptions::RuntimeError, "StridedMap::NodeId2GlobalDofIds: stridingInfo not valid: stridingInfo.size() = 0?");
      std::vector<GlobalOrdinal> dofs;
      if(stridedBlockId_ > -1) {
          TEUCHOS_TEST_FOR_EXCEPTION(stridingInfo_[stridedBlockId_] == 0, Exceptions::RuntimeError, "StridedMap::NodeId2GlobalDofIds: stridingInfo not valid: stridingInfo[stridedBlockId] = 0?");
 
          // determine nStridedOffset
          size_t nStridedOffset = 0;
          for(int j=0; j<stridedBlockId_; j++) {
            nStridedOffset += stridingInfo_[j];
          }
 
          for(size_t i = 0; i<stridingInfo_[stridedBlockId_]; i++) {
            GlobalOrdinal gid =
                nodeId * Teuchos::as<GlobalOrdinal>(getFixedBlockSize()) +
                offset_ +
                Teuchos::as<GlobalOrdinal>(nStridedOffset) +
                Teuchos::as<GlobalOrdinal>(i);
            dofs.push_back(gid);
          }
      } else {
        for(size_t i = 0; i<getFixedBlockSize(); i++) {
          GlobalOrdinal gid =
              nodeId * Teuchos::as<GlobalOrdinal>(getFixedBlockSize()) +
              offset_ +
              Teuchos::as<GlobalOrdinal>(i);
          dofs.push_back(gid);
        }
      }
      return dofs;
    }*/
 
  private:
    virtual bool CheckConsistency() {
#ifndef HAVE_XPETRA_DEBUG
      return true;
#else
      if (getStridedBlockId() == -1) {
        // Strided map contains the full map
        if (getNodeNumElements()   % getFixedBlockSize() != 0 ||    // number of local  elements is not a multiple of block size
            getGlobalNumElements() % getFixedBlockSize() != 0)      // number of global    -//-
          return false;
 
      } else {
        // Strided map contains only the partial map
        Teuchos::ArrayView<const GlobalOrdinal> dofGids = getNodeElementList();
        // std::sort(dofGids.begin(), dofGids.end());
 
        if (dofGids.size() == 0)  // special treatment for empty processors
          return true;
 
        if (dofGids.size() % stridingInfo_[stridedBlockId_] != 0)
          return false;
 
 
        // Calculate nStridedOffset
        size_t nStridedOffset = 0;
        for (int j = 0; j < stridedBlockId_; j++)
          nStridedOffset += stridingInfo_[j];
 
        const GlobalOrdinal goStridedOffset = Teuchos::as<GlobalOrdinal>(nStridedOffset);
        const GlobalOrdinal goZeroOffset    = (dofGids[0] - nStridedOffset - offset_ - indexBase_) / Teuchos::as<GlobalOrdinal>(getFixedBlockSize());
 
        GlobalOrdinal cnt = 0;
        for (size_t i = 0; i < Teuchos::as<size_t>(dofGids.size())/stridingInfo_[stridedBlockId_]; i += stridingInfo_[stridedBlockId_]) {
          const GlobalOrdinal first_gid = dofGids[i];
 
          // We expect this to be the same for all DOFs of the same node
          cnt = (first_gid - goStridedOffset - offset_ - indexBase_) / Teuchos::as<GlobalOrdinal>(getFixedBlockSize()) - goZeroOffset;
 
          // Loop over all DOFs that belong to current node
          for (size_t j = 0; j < stridingInfo_[stridedBlockId_]; j++) {
            const GlobalOrdinal gid = dofGids[i+j];
            const GlobalOrdinal r   = (gid - Teuchos::as<GlobalOrdinal>(j) - goStridedOffset - offset_ - indexBase_) /
                                      Teuchos::as<GlobalOrdinal>(getFixedBlockSize()) - goZeroOffset - cnt;
            // TAW 1/18/2016: We cannot use Teuchos::OrdinalTraits<GlobalOrdinal>::zero() ) here,
            //                If, e.g., GO=long long is disabled, OrdinalTraits<long long> is not available.
            //                But we instantiate stubs on GO=long long which might contain StridedMaps.
            //                These lead to compilation errors, then.
            if (r != 0 ) {
              std::cout << "goZeroOffset   : " <<  goZeroOffset << std::endl
                        << "dofGids[0]     : " <<  dofGids[0] << std::endl
                        << "stridedOffset  : " <<  nStridedOffset << std::endl
                        << "offset_        : " <<  offset_ << std::endl
                        << "goStridedOffset: " <<  goStridedOffset << std::endl
                        << "getFixedBlkSize: " <<  getFixedBlockSize() << std::endl
                        << "gid: " << gid << " GID: " << r << std::endl;
 
              return false;
            }
          }
        }
      }
 
      return true;
#endif
    }
 
  private:
    RCP<const Xpetra::Map<LocalOrdinal, GlobalOrdinal, Node> >      map_;
 
    std::vector<size_t> stridingInfo_;      
    LocalOrdinal        stridedBlockId_;    
                                            //     stridedBlock == -1: the full map (with all strided block dofs)
                                            //     stridedBlock  > -1: only dofs of strided block with index "stridedBlockId" are stored in this map
    GlobalOrdinal       offset_;                    
    GlobalOrdinal       indexBase_;         
 
  public:
 
 
 
    global_size_t getGlobalNumElements() const { return map_->getGlobalNumElements(); }
 
    size_t getNodeNumElements() const { return map_->getNodeNumElements(); }
 
    GlobalOrdinal getIndexBase() const { return map_->getIndexBase(); }
 
    LocalOrdinal getMinLocalIndex() const { return map_->getMinLocalIndex(); }
 
    LocalOrdinal getMaxLocalIndex() const { return map_->getMaxLocalIndex(); }
 
    GlobalOrdinal getMinGlobalIndex() const { return map_->getMinGlobalIndex(); }
 
    GlobalOrdinal getMaxGlobalIndex() const { return map_->getMaxGlobalIndex(); }
 
    GlobalOrdinal getMinAllGlobalIndex() const { return map_->getMinAllGlobalIndex(); }
 
    GlobalOrdinal getMaxAllGlobalIndex() const { return map_->getMaxAllGlobalIndex(); }
 
    LocalOrdinal getLocalElement(GlobalOrdinal globalIndex) const { return map_->getLocalElement(globalIndex); }
 
    GlobalOrdinal getGlobalElement(LocalOrdinal localIndex) const { return map_->getGlobalElement(localIndex); }
 
    LookupStatus getRemoteIndexList(const Teuchos::ArrayView<const GlobalOrdinal> &GIDList, const Teuchos::ArrayView<int> &nodeIDList, const Teuchos::ArrayView<LocalOrdinal> &LIDList) const {
      return map_->getRemoteIndexList(GIDList, nodeIDList, LIDList);
    }
 
    LookupStatus getRemoteIndexList(const Teuchos::ArrayView<const GlobalOrdinal> &GIDList, const Teuchos::ArrayView< int > &nodeIDList) const {
      return map_->getRemoteIndexList(GIDList, nodeIDList);
    }
 
    Teuchos::ArrayView< const GlobalOrdinal > getNodeElementList() const { return map_->getNodeElementList(); }
 
    bool isNodeLocalElement(LocalOrdinal localIndex) const { return map_->isNodeLocalElement(localIndex); }
 
    bool isNodeGlobalElement(GlobalOrdinal globalIndex) const { return map_->isNodeGlobalElement(globalIndex); }
 
    bool isContiguous() const { return map_->isContiguous(); }
 
    bool isDistributed() const { return map_->isDistributed(); }
 
 
    bool isCompatible(const Map& map) const { return map_->isCompatible(map); }
 
    bool isSameAs(const Map& map) const { return map_->isSameAs(map); }
 
    Teuchos::RCP< const Teuchos::Comm< int > > getComm() const { return map_->getComm(); }
 
    Teuchos::RCP<Node>  getNode() const { return map_->getNode(); }
 
    RCP<const Map> removeEmptyProcesses  () const { return map_->removeEmptyProcesses(); }
    RCP<const Map> replaceCommWithSubset (const Teuchos::RCP<const Teuchos::Comm<int> >& newComm) const { return map_->replaceCommWithSubset(newComm); }
 
    std::string description() const { return map_->description(); }
 
    void describe(Teuchos::FancyOStream& out, const Teuchos::EVerbosityLevel verbLevel = Teuchos::Describable::verbLevel_default) const { map_->describe(out, verbLevel); }
 
    UnderlyingLib lib() const { return map_->lib(); }
 
  }; // StridedMap class
 
} // Xpetra namespace
 
#define XPETRA_STRIDEDMAP_SHORT
#endif // XPETRA_STRIDEDMAP_HPP