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(
00045     const Teuchos::RCP<const Tpetra::BlockMap<LocalOrdinal,GlobalOrdinal,Node> >& blockMap,
00046     const Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> >& ptColMap)
00047 {
00048   global_size_t numGlobalBlocks = Teuchos::OrdinalTraits<global_size_t>::invalid();
00049   Teuchos::ArrayView<const GlobalOrdinal> blockIDs = ptColMap->getNodeElementList();
00050   Teuchos::ArrayRCP<const LocalOrdinal> firstPoints = blockMap->getNodeFirstPointInBlocks();
00051   Teuchos::Array<GlobalOrdinal> points(firstPoints.size()-1);
00052   Teuchos::Array<LocalOrdinal> blockSizes(firstPoints.size()-1);
00053 
00054   typedef typename Teuchos::ArrayView<const LocalOrdinal>::size_type Tsize_t;
00055   for(Tsize_t i=0; i<blockSizes.size(); ++i) {
00056     points[i] = blockMap->getFirstGlobalPointInLocalBlock(i);
00057     blockSizes[i] = firstPoints[i+1]-firstPoints[i];
00058   }
00059 
00060   //We will create a block-map where each block corresponds to a point in
00061   //the input-map 'ptColMap', and where each block's size is obtained from
00062   //the input-block-map 'blockMap'.
00063 
00064   //if we are running on a single processor, then it's easy because
00065   //we know that blockMap is distributed the same as ptColMap:
00066   int numProcessors = ptColMap->getComm()->getSize();
00067   if (numProcessors == 1) {
00068     return Teuchos::rcp(new Tpetra::BlockMap<LocalOrdinal,GlobalOrdinal,Node>(
00069         numGlobalBlocks, blockIDs, points(), blockSizes(),
00070         ptColMap->getIndexBase(), ptColMap->getComm(), ptColMap->getNode()
00071     ));
00072   }
00073 
00074   //If we get to here, we're running on multiple processors, and blockMap
00075   //is probably not distributed the same as ptColMap, so we have to do
00076   //some communication to get the block-sizes from blockMap corresponding
00077   //to the blockIDs we got from ptColMap.
00078   //We also have to do communication to get the global first-points-in-block
00079   //for the blocks in the new block-column-map.
00080   //
00081   //I think the simplest way to do this is to create vectors where the values
00082   //are block-sizes (or first-points-in-block), and import one to the other.
00083   typedef Tpetra::Vector<GlobalOrdinal,LocalOrdinal,GlobalOrdinal,Node> GOVec;
00084   typedef Tpetra::Vector<LocalOrdinal,LocalOrdinal,GlobalOrdinal,Node> LOVec;
00085 
00086   Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > blkPtMap =
00087       convertBlockMapToPointMap(*blockMap);
00088 
00089   LOVec source_sizes(blkPtMap, blockSizes());
00090   GOVec source_points(blkPtMap, points());
00091   LOVec target_sizes(ptColMap);
00092   GOVec target_points(ptColMap);
00093 
00094   Tpetra::Import<LocalOrdinal,GlobalOrdinal,Node> importer(blkPtMap, ptColMap);
00095   target_sizes.doImport(source_sizes, importer, REPLACE);
00096   target_points.doImport(source_points, importer, REPLACE);
00097 
00098   //now we can create our block-column-map:
00099   return Teuchos::rcp(new Tpetra::BlockMap<LocalOrdinal,GlobalOrdinal,Node>(
00100     numGlobalBlocks, blockIDs, target_points.get1dView()(), target_sizes.get1dView()(),
00101     ptColMap->getIndexBase(), ptColMap->getComm(), ptColMap->getNode()
00102   ));
00103 }
00104 
00105 template<class LocalOrdinal,class GlobalOrdinal,class Node>
00106 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::BlockCrsGraph(
00107    const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> >& blkRowMap,
00108    size_t maxNumEntriesPerRow, ProfileType pftype
00109     )
00110  : ptGraph_(),
00111    blkRowMap_(blkRowMap),
00112    blkColMap_(),
00113    blkDomainMap_(),
00114    blkRangeMap_()
00115 {
00116   ptGraph_ = Teuchos::rcp(new CrsGraph<LocalOrdinal,GlobalOrdinal,Node>(
00117           convertBlockMapToPointMap(*blkRowMap), maxNumEntriesPerRow, pftype));
00118 }
00119 
00120 //-------------------------------------------------------------------
00121 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00122 void
00123 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::insertGlobalIndices(GlobalOrdinal row, const Teuchos::ArrayView<const GlobalOrdinal> &indices)
00124 {
00125   ptGraph_->insertGlobalIndices(row, indices);
00126 }
00127 
00128 //-------------------------------------------------------------------
00129 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00130 Teuchos::ArrayRCP<const size_t>
00131 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getNodeRowOffsets() const
00132 {
00133   return ptGraph_->getNodeRowBegs();
00134 }
00135 
00136 //-------------------------------------------------------------------
00137 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00138 Teuchos::ArrayRCP<const LocalOrdinal>
00139 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getNodePackedIndices() const
00140 {
00141   return ptGraph_->getNodePackedIndices();
00142 }
00143 
00144 //-------------------------------------------------------------------
00145 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00146 void
00147 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::globalAssemble()
00148 {
00149   ptGraph_->globalAssemble();
00150 }
00151 
00152 //-------------------------------------------------------------------
00153 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00154 void
00155 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)
00156 {
00157   blkDomainMap_ = blkDomainMap;
00158   blkRangeMap_  = blkRangeMap;
00159 
00160   ptGraph_->fillComplete(convertBlockMapToPointMap(*blkDomainMap),
00161                          convertBlockMapToPointMap(*blkRangeMap), os);
00162 
00163   //Now we need to take the point-column-map from ptGraph_ and create a
00164   //corresponding block-column-map.
00165   //Our block-column-map will have the same distribution as
00166   //blkGraph_->getColMap, and we'll get block-size info from the
00167   //blkDomainMap_. This will require some communication in cases
00168   //where blkDomainMap is distributed differently.
00169   blkColMap_ = makeBlockColumnMap(blkDomainMap_, ptGraph_->getColMap());
00170 }
00171 
00172 //-------------------------------------------------------------------
00173 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00174 void
00175 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::fillComplete(OptimizeOption os)
00176 {
00177   fillComplete(getBlockRowMap(), getBlockRowMap(), os);
00178 }
00179 
00180 //-------------------------------------------------------------------
00181 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00182 void
00183 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::optimizeStorage()
00184 {
00185   if (ptGraph_->isFillComplete()) {
00186     ptGraph_->resumeFill();
00187   }
00188   ptGraph_->fillComplete(DoOptimizeStorage);
00189 }
00190 
00191 //-------------------------------------------------------------------
00192 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00193 bool
00194 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::isFillComplete() const
00195 {
00196   return ptGraph_->isFillComplete();
00197 }
00198 
00199 //-------------------------------------------------------------------
00200 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00201 bool
00202 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::isUpperTriangular() const
00203 {
00204   return ptGraph_->isUpperTriangular();
00205 }
00206 
00207 //-------------------------------------------------------------------
00208 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00209 bool
00210 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::isLowerTriangular() const
00211 {
00212   return ptGraph_->isLowerTriangular();
00213 }
00214 
00215 //-------------------------------------------------------------------
00216 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00217 bool
00218 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::isLocallyIndexed() const
00219 {
00220   return ptGraph_->isLocallyIndexed();
00221 }
00222 
00223 //-------------------------------------------------------------------
00224 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00225 size_t
00226 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getNodeNumBlockEntries() const
00227 {
00228   return ptGraph_->getNodeNumEntries();
00229 }
00230 
00231 //-------------------------------------------------------------------
00232 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00233 size_t
00234 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getGlobalBlockRowLength(GlobalOrdinal row) const
00235 {
00236   return ptGraph_->getNumEntriesInGlobalRow(row);
00237 }
00238 
00239 //-------------------------------------------------------------------
00240 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00241 void
00242 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getGlobalBlockRowView(GlobalOrdinal row,
00243                              Teuchos::ArrayView<const GlobalOrdinal>& blockCols) const
00244 {
00245   ptGraph_->getGlobalRowView(row, blockCols);
00246 }
00247 
00248 //-------------------------------------------------------------------
00249 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00250 void
00251 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getLocalBlockRowView(LocalOrdinal row,
00252                              Teuchos::ArrayView<const LocalOrdinal>& blockCols) const
00253 {
00254   ptGraph_->getLocalRowView(row, blockCols);
00255 }
00256 
00257 //-------------------------------------------------------------------
00258 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00259 size_t
00260 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getNodeNumBlockRows() const
00261 {
00262   return ptGraph_->getNodeNumRows();
00263 }
00264 
00265 //-------------------------------------------------------------------
00266 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00267 size_t
00268 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getGlobalNumBlockRows() const
00269 {
00270   return ptGraph_->getGlobalNumRows();
00271 }
00272 
00273 //-------------------------------------------------------------------
00274 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00275 size_t
00276 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getNodeNumBlockDiags() const
00277 {
00278   return ptGraph_->getNodeNumDiags();
00279 }
00280 
00281 //-------------------------------------------------------------------
00282 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00283 const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > &
00284 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getBlockRowMap() const
00285 {
00286   return blkRowMap_;
00287 }
00288 
00289 //-------------------------------------------------------------------
00290 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00291 const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > &
00292 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getBlockColMap() const
00293 {
00294   return blkColMap_;
00295 }
00296 
00297 //-------------------------------------------------------------------
00298 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00299 const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > &
00300 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getBlockDomainMap() const
00301 {
00302   return blkDomainMap_;
00303 }
00304 
00305 //-------------------------------------------------------------------
00306 template<class LocalOrdinal, class GlobalOrdinal, class Node>
00307 const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > &
00308 BlockCrsGraph<LocalOrdinal,GlobalOrdinal,Node>::getBlockRangeMap() const
00309 {
00310   return blkRangeMap_;
00311 }
00312 
00313 //
00314 // Explicit instantiation macro
00315 //
00316 // Must be expanded from within the Tpetra namespace!
00317 //
00318 
00319 #define TPETRA_BLOCKCRSGRAPH_INSTANT(LO,GO,NODE) \
00320   \
00321   template class BlockCrsGraph< LO , GO , NODE >;
00322 
00323 
00324 }//namespace Tpetra
00325 
00326 #endif
00327 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines