Tpetra Matrix/Vector Services Version of the Day
Tpetra_BlockCrsGraph_def.hpp
00001 // @HEADER
00002 // ***********************************************************************
00003 // 
00004 //          Tpetra: Templated Linear Algebra Services Package
00005 //                 Copyright (2008) Sandia Corporation
00006 // 
00007 // Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive
00008 // license for use of this work by or on behalf of the U.S. Government.
00009 // 
00010 // This library is free software; you can redistribute it and/or modify
00011 // it under the terms of the GNU Lesser General Public License as
00012 // published by the Free Software Foundation; either version 2.1 of the
00013 // License, or (at your option) any later version.
00014 //  
00015 // This library is distributed in the hope that it will be useful, but
00016 // WITHOUT ANY WARRANTY; without even the implied warranty of
00017 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
00018 // Lesser General Public License for more details.
00019 //  
00020 // You should have received a copy of the GNU Lesser General Public
00021 // License along with this library; if not, write to the Free Software
00022 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
00023 // USA
00024 // Questions? Contact Michael A. Heroux (maherou@sandia.gov) 
00025 // 
00026 // ***********************************************************************
00027 // @HEADER
00028 
00029 #ifndef TPETRA_BLOCKCRSGRAPH_DEF_HPP
00030 #define TPETRA_BLOCKCRSGRAPH_DEF_HPP
00031 
00032 #include "Tpetra_ConfigDefs.hpp"
00033 #include "Tpetra_Vector.hpp"
00034 
00035 #ifdef DOXYGEN_USE_ONLY
00036 #include "Tpetra_BlockCrsGraph_decl.hpp"
00037 #endif
00038 
00039 namespace Tpetra {
00040 
00041 //-----------------------------------------------------------------
00042 template<class LocalOrdinal,class GlobalOrdinal,class Node>
00043 Teuchos::RCP<const Tpetra::BlockMap<LocalOrdinal,GlobalOrdinal,Node> >
00044 makeBlockColumnMap(const Teuchos::RCP<const Tpetra::BlockMap<LocalOrdinal,GlobalOrdinal,Node> >& blockMap, const Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> >& ptColMap)
00045 {
00046   global_size_t numGlobalBlocks = Teuchos::OrdinalTraits<global_size_t>::invalid();
00047   Teuchos::ArrayView<const GlobalOrdinal> blockIDs = ptColMap->getNodeElementList();
00048   Teuchos::ArrayRCP<const LocalOrdinal> firstPoints = blockMap->getNodeFirstPointInBlocks();
00049   Teuchos::Array<GlobalOrdinal> points(firstPoints.size()-1);
00050   Teuchos::Array<LocalOrdinal> blockSizes(firstPoints.size()-1);
00051 
00052   typedef typename Teuchos::ArrayView<const LocalOrdinal>::size_type Tsize_t;
00053   for(Tsize_t i=0; i<blockSizes.size(); ++i) {
00054     points[i] = blockMap->getFirstGlobalPointInLocalBlock(i);
00055     blockSizes[i] = firstPoints[i+1]-firstPoints[i];
00056   }
00057 
00058   //We will create a block-map where each block corresponds to a point in
00059   //the input-map 'ptColMap', and where each block's size is obtained from
00060   //the input-block-map 'blockMap'.
00061 
00062   //if we are running on a single processor, then it's easy because
00063   //we know that blockMap is distributed the same as ptColMap:
00064   int numProcessors = ptColMap->getComm()->getSize();
00065   if (numProcessors == 1) {
00066     return Teuchos::rcp(new Tpetra::BlockMap<LocalOrdinal,GlobalOrdinal,Node>(
00067         numGlobalBlocks, blockIDs, points(), blockSizes(),
00068         ptColMap->getIndexBase(), ptColMap->getComm(), ptColMap->getNode()
00069     ));
00070   }
00071 
00072   //If we get to here, we're running on multiple processors, and blockMap
00073   //is probably not distributed the same as ptColMap, so we have to do
00074   //some communication to get the block-sizes from blockMap corresponding
00075   //to the blockIDs we got from ptColMap.
00076   //We also have to do communication to get the global first-points-in-block
00077   //for the blocks in the new block-column-map.
00078   //
00079   //I think the simplest way to do this is to create vectors where the values
00080   //are block-sizes (or first-points-in-block), and import one to the other.
00081   typedef Tpetra::Vector<GlobalOrdinal,LocalOrdinal,GlobalOrdinal,Node> GOVec;
00082   typedef Tpetra::Vector<LocalOrdinal,LocalOrdinal,GlobalOrdinal,Node> LOVec;
00083 
00084   Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > blkPtMap =
00085       convertBlockMapToPointMap(blockMap);
00086 
00087   LOVec source_sizes(blkPtMap, blockSizes());
00088   GOVec source_points(blkPtMap, points());
00089   LOVec target_sizes(ptColMap);
00090   GOVec target_points(ptColMap);
00091 
00092   Tpetra::Import<LocalOrdinal,GlobalOrdinal,Node> importer(blkPtMap, ptColMap);
00093   target_sizes.doImport(source_sizes, importer, REPLACE);
00094   target_points.doImport(source_points, importer, REPLACE);
00095 
00096   //now we can create our block-column-map:
00097   return Teuchos::rcp(new Tpetra::BlockMap<LocalOrdinal,GlobalOrdinal,Node>(
00098     numGlobalBlocks, blockIDs, target_points.get1dView()(), target_sizes.get1dView()(),
00099     ptColMap->getIndexBase(), ptColMap->getComm(), ptColMap->getNode()
00100   ));
00101 }
00102 
00103 template<class LocalOrdinal,class GlobalOrdinal,class Node>
00104 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::BlockCrsGraph(
00105    const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> >& blkRowMap,
00106    size_t maxNumEntriesPerRow, ProfileType pftype
00107     )
00108  : ptGraph_(),
00109    blkRowMap_(blkRowMap),
00110    blkColMap_(),
00111    blkDomainMap_(),
00112    blkRangeMap_()
00113 {
00114   ptGraph_ = Teuchos::rcp(new CrsGraph<LocalOrdinal,GlobalOrdinal,Node>(
00115           convertBlockMapToPointMap(blkRowMap), maxNumEntriesPerRow, pftype));
00116 }
00117 
00118 //-------------------------------------------------------------------
00119 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00120 void
00121 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::insertGlobalIndices(GlobalOrdinal row, const Teuchos::ArrayView<const GlobalOrdinal> &indices)
00122 {
00123   ptGraph_->insertGlobalIndices(row, indices);
00124 }
00125 
00126 //-------------------------------------------------------------------
00127 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00128 Teuchos::ArrayRCP<const size_t>
00129 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getNodeRowOffsets() const
00130 {
00131   return ptGraph_->getNodeRowBegs();
00132 }
00133 
00134 //-------------------------------------------------------------------
00135 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00136 Teuchos::ArrayRCP<const LocalOrdinal>
00137 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getNodePackedIndices() const
00138 {
00139   return ptGraph_->getNodePackedIndices();
00140 }
00141 
00142 //-------------------------------------------------------------------
00143 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00144 void
00145 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::globalAssemble()
00146 {
00147   ptGraph_->globalAssemble();
00148 }
00149 
00150 //-------------------------------------------------------------------
00151 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00152 void
00153 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::fillComplete(const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > &blkDomainMap, const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > &blkRangeMap, OptimizeOption os)
00154 {
00155   blkDomainMap_ = blkDomainMap;
00156   blkRangeMap_  = blkRangeMap;
00157 
00158   ptGraph_->fillComplete(convertBlockMapToPointMap(blkDomainMap),
00159                          convertBlockMapToPointMap(blkRangeMap), os);
00160 
00161   //Now we need to take the point-column-map from ptGraph_ and create a
00162   //corresponding block-column-map.
00163   //Our block-column-map will have the same distribution as
00164   //blkGraph_->getColMap, and we'll get block-size info from the
00165   //blkDomainMap_. This will require some communication in cases
00166   //where blkDomainMap is distributed differently.
00167   blkColMap_ = makeBlockColumnMap(blkDomainMap_, ptGraph_->getColMap());
00168 }
00169 
00170 //-------------------------------------------------------------------
00171 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00172 void
00173 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::fillComplete(OptimizeOption os)
00174 {
00175   fillComplete(getBlockRowMap(), getBlockRowMap(), os);
00176 }
00177 
00178 //-------------------------------------------------------------------
00179 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00180 void
00181 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::optimizeStorage()
00182 {
00183   if (ptGraph_->isFillComplete()) {
00184     ptGraph_->resumeFill();
00185   }
00186   ptGraph_->fillComplete(DoOptimizeStorage);
00187 }
00188 
00189 //-------------------------------------------------------------------
00190 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00191 bool
00192 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::isFillComplete() const
00193 {
00194   return ptGraph_->isFillComplete();
00195 }
00196 
00197 //-------------------------------------------------------------------
00198 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00199 bool
00200 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::isUpperTriangular() const
00201 {
00202   return ptGraph_->isUpperTriangular();
00203 }
00204 
00205 //-------------------------------------------------------------------
00206 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00207 bool
00208 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::isLowerTriangular() const
00209 {
00210   return ptGraph_->isLowerTriangular();
00211 }
00212 
00213 //-------------------------------------------------------------------
00214 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00215 bool
00216 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::isLocallyIndexed() const
00217 {
00218   return ptGraph_->isLocallyIndexed();
00219 }
00220 
00221 //-------------------------------------------------------------------
00222 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00223 size_t
00224 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getNodeNumEntries() const
00225 {
00226   return ptGraph_->getNodeNumEntries();
00227 }
00228 
00229 //-------------------------------------------------------------------
00230 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00231 size_t
00232 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getNodeNumRows() const
00233 {
00234   return ptGraph_->getNodeNumRows();
00235 }
00236 
00237 //-------------------------------------------------------------------
00238 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00239 size_t
00240 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getNodeNumDiags() const
00241 {
00242   return ptGraph_->getNodeNumDiags();
00243 }
00244 
00245 //-------------------------------------------------------------------
00246 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00247 const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > &
00248 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getBlockRowMap() const
00249 {
00250   return blkRowMap_;
00251 }
00252 
00253 //-------------------------------------------------------------------
00254 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00255 const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > &
00256 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getBlockColMap() const
00257 {
00258   return blkColMap_;
00259 }
00260 
00261 //-------------------------------------------------------------------
00262 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00263 const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > &
00264 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getBlockDomainMap() const
00265 {
00266   return blkDomainMap_;
00267 }
00268 
00269 //-------------------------------------------------------------------
00270 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00271 const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > &
00272 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getBlockRangeMap() const
00273 {
00274   return blkRangeMap_;
00275 }
00276 
00277 //
00278 // Explicit instantiation macro
00279 //
00280 // Must be expanded from within the Tpetra namespace!
00281 //
00282 
00283 #define TPETRA_BLOCKCRSGRAPH_INSTANT(LO,GO,NODE) \
00284   \
00285   template class BlockCrsGraph< LO , GO , NODE >;
00286 
00287 
00288 }//namespace Tpetra
00289 
00290 #endif
00291 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines