Tpetra_VbrMatrix_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_VBRMATRIX_DEF_HPP
00030 #define TPETRA_VBRMATRIX_DEF_HPP
00031 
00032 #include <algorithm>
00033 #include <sstream>
00034 
00035 #include <Kokkos_NodeHelpers.hpp>
00036 #include <Teuchos_Array.hpp>
00037 #include <Teuchos_CommHelpers.hpp>
00038 #include <Teuchos_ScalarTraits.hpp>
00039 #include <Teuchos_OrdinalTraits.hpp>
00040 #include <Teuchos_TypeNameTraits.hpp>
00041 #include <Tpetra_Vector.hpp>
00042 
00043 #ifdef DOXYGEN_USE_ONLY
00044 #include "Tpetra_VbrMatrix_decl.hpp"
00045 #endif
00046 
00047 namespace Tpetra {
00048 
00049 //-----------------------------------------------------------------
00050 template<class LocalOrdinal,class GlobalOrdinal,class Node>
00051 Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> >
00052 convertBlockMapToPointMap(const Teuchos::RCP<const Tpetra::BlockMap<LocalOrdinal,GlobalOrdinal,Node> >& blockMap)
00053 {
00054   global_size_t numGlobalElems = Teuchos::OrdinalTraits<global_size_t>::invalid();
00055   GlobalOrdinal indexBase = blockMap->getPointMap()->getIndexBase();
00056   const Teuchos::RCP<const Teuchos::Comm<int> >& comm = blockMap->getPointMap()->getComm();
00057   const Teuchos::RCP<Node>& node = blockMap->getPointMap()->getNode();
00058 
00059   //Create a point-entry map where each point
00060   //corresponds to a block in the block-map:
00061   return Teuchos::rcp(new Map<LocalOrdinal,GlobalOrdinal,Node>(numGlobalElems, blockMap->getNodeBlockIDs(), indexBase, comm, node));
00062 }
00063 
00064 //-----------------------------------------------------------------
00065 template<class LocalOrdinal,class GlobalOrdinal,class Node>
00066 Teuchos::RCP<const Tpetra::BlockMap<LocalOrdinal,GlobalOrdinal,Node> >
00067 makeBlockColumnMap(const Teuchos::RCP<const Tpetra::BlockMap<LocalOrdinal,GlobalOrdinal,Node> >& blockMap, const Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> >& ptColMap)
00068 {
00069   global_size_t numGlobalBlocks = Teuchos::OrdinalTraits<global_size_t>::invalid();
00070   Teuchos::ArrayView<const GlobalOrdinal> blockIDs = ptColMap->getNodeElementList();
00071   Teuchos::ArrayRCP<const LocalOrdinal> firstPoints = blockMap->getNodeFirstPointInBlocks();
00072   Teuchos::Array<GlobalOrdinal> points(firstPoints.size()-1);
00073   Teuchos::Array<LocalOrdinal> blockSizes(firstPoints.size()-1);
00074 
00075   typedef typename Teuchos::ArrayView<const LocalOrdinal>::size_type Tsize_t;
00076   for(Tsize_t i=0; i<blockSizes.size(); ++i) {
00077     points[i] = blockMap->getFirstGlobalPointInLocalBlock(i);
00078     blockSizes[i] = firstPoints[i+1]-firstPoints[i];
00079   }
00080 
00081   //We will create a block-map where each block corresponds to a point in
00082   //the input-map 'ptColMap', and where each block's size is obtained from
00083   //the input-block-map 'blockMap'.
00084 
00085   //if we are running on a single processor, then it's easy because
00086   //we know that blockMap is distributed the same as ptColMap:
00087   int numProcessors = ptColMap->getComm()->getSize();
00088   if (numProcessors == 1) {
00089     return Teuchos::rcp(new Tpetra::BlockMap<LocalOrdinal,GlobalOrdinal,Node>(
00090         numGlobalBlocks, blockIDs, points(), blockSizes(),
00091         ptColMap->getIndexBase(), ptColMap->getComm(), ptColMap->getNode()
00092     ));
00093   }
00094 
00095   //If we get to here, we're running on multiple processors, and blockMap
00096   //is probably not distributed the same as ptColMap, so we have to do
00097   //some communication to get the block-sizes from blockMap corresponding
00098   //to the blockIDs we got from ptColMap.
00099   //We also have to do communication to get the global first-points-in-block
00100   //for the blocks in the new block-column-map.
00101   //
00102   //I think the simplest way to do this is to create vectors where the values
00103   //are block-sizes (or first-points-in-block), and import one to the other.
00104   typedef Tpetra::Vector<GlobalOrdinal,LocalOrdinal,GlobalOrdinal,Node> GOVec;
00105   typedef Tpetra::Vector<LocalOrdinal,LocalOrdinal,GlobalOrdinal,Node> LOVec;
00106 
00107   Teuchos::RCP<const Tpetra::Map<LocalOrdinal,GlobalOrdinal,Node> > blkPtMap =
00108       convertBlockMapToPointMap(blockMap);
00109 
00110   LOVec source_sizes(blkPtMap, blockSizes());
00111   GOVec source_points(blkPtMap, points());
00112   LOVec target_sizes(ptColMap);
00113   GOVec target_points(ptColMap);
00114 
00115   Tpetra::Import<LocalOrdinal,GlobalOrdinal,Node> importer(blkPtMap, ptColMap);
00116   target_sizes.doImport(source_sizes, importer, REPLACE);
00117   target_points.doImport(source_points, importer, REPLACE);
00118 
00119   //now we can create our block-column-map:
00120   return Teuchos::rcp(new Tpetra::BlockMap<LocalOrdinal,GlobalOrdinal,Node>(
00121     numGlobalBlocks, blockIDs, target_points.get1dView()(), target_sizes.get1dView()(),
00122     ptColMap->getIndexBase(), ptColMap->getComm(), ptColMap->getNode()
00123   ));
00124 }
00125 
00126 template<class Scalar,class Node>
00127 void fill_device_ArrayRCP(Teuchos::RCP<Node>& node, Teuchos::ArrayRCP<Scalar>& ptr, Scalar value)
00128 {
00129   Kokkos::ReadyBufferHelper<Node> rbh(node);
00130   Kokkos::InitOp<Scalar> wdp;
00131   wdp.alpha = value;
00132   rbh.begin();
00133   wdp.x = rbh.addNonConstBuffer(ptr);
00134   rbh.end();
00135   node->template parallel_for<Kokkos::InitOp<Scalar> >(0, ptr.size(), wdp);
00136 }
00137 
00138 //-------------------------------------------------------------------
00139 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00140 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::VbrMatrix(const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > &blkRowMap, size_t maxNumEntriesPerRow, ProfileType pftype)
00141  : blkRowMap_(blkRowMap),
00142    blkColMap_(Teuchos::null),
00143    blkDomainMap_(Teuchos::null),
00144    blkRangeMap_(Teuchos::null),
00145    blkGraph_(),
00146    lclMatrix_(blkRowMap->getNodeNumBlocks(), blkRowMap->getPointMap()->getNode()),
00147    pbuf_values1D_(),
00148    pbuf_bptr_(),
00149    pbuf_bindx_(),
00150    pbuf_indx_(),
00151    lclMatVec_(blkRowMap->getPointMap()->getNode()),
00152    importer_(),
00153    exporter_(),
00154    importedVec_(),
00155    exportedVec_(),
00156    col_ind_2D_global_(Teuchos::rcp(new Teuchos::Array<std::map<GlobalOrdinal,Teuchos::ArrayRCP<Scalar> > >)),
00157    values2D_(Teuchos::rcp(new Teuchos::Array<Teuchos::Array<Teuchos::ArrayRCP<Scalar> > >)),
00158    is_fill_completed_(false),
00159    is_storage_optimized_(false)
00160 {
00161   //The graph of this VBR matrix will be a CrsGraph where each entry in the
00162   //graph corresponds to a block-entry in the matrix.
00163   //That is, you can think of a VBR matrix as a Crs matrix of dense
00164   //submatrices...
00165 
00166   //First create a point-entry map where each point corresponds to
00167   //a block in our block-row-map:
00168   Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >
00169       blkPointMap = convertBlockMapToPointMap(blkRowMap);
00170 
00171   //now create the graph:
00172   blkGraph_ = Teuchos::rcp(new CrsGraph<LocalOrdinal,GlobalOrdinal,Node>(blkPointMap, maxNumEntriesPerRow, pftype));
00173 }
00174 
00175 //-------------------------------------------------------------------
00176 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00177 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::~VbrMatrix()
00178 {
00179 }
00180 
00181 //-------------------------------------------------------------------
00182 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00183 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > &
00184 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::getDomainMap() const
00185 {
00186   return blkDomainMap_->getPointMap();
00187 }
00188 
00189 //-------------------------------------------------------------------
00190 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00191 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > &
00192 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::getRangeMap() const
00193 {
00194   return blkRangeMap_->getPointMap();
00195 }
00196 
00197 //-------------------------------------------------------------------
00198 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00199 bool
00200 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::isFillComplete() const
00201 {
00202   return is_fill_completed_;
00203 }
00204 
00205 //-------------------------------------------------------------------
00206 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00207 template <class DomainScalar, class RangeScalar>
00208 void
00209 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::multiply(const MultiVector<DomainScalar,LocalOrdinal,GlobalOrdinal,Node> & X, MultiVector<RangeScalar,LocalOrdinal,GlobalOrdinal,Node> &Y, Teuchos::ETransp trans, RangeScalar alpha, RangeScalar beta) const
00210 {
00211   TEST_FOR_EXCEPTION(!isFillComplete(), std::runtime_error,
00212     "Tpetra::VbrMatrix::multiply ERROR, multiply may only be called after fillComplete has been called.");
00213 
00214   const Kokkos::MultiVector<Scalar,Node> *lclX = &X.getLocalMV();
00215   Kokkos::MultiVector<Scalar,Node>        *lclY = &Y.getLocalMVNonConst();
00216 
00217   lclMatVec_.template multiply<Scalar,Scalar>(trans,alpha,*lclX,beta,*lclY);
00218 }
00219 
00220 //-------------------------------------------------------------------
00221 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00222 void
00223 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::updateImport(const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& X) const
00224 {
00225   if (importer_ != Teuchos::null) {
00226     if (importedVec_ == Teuchos::null || importedVec_->getNumVectors() != X.getNumVectors()) {
00227       importedVec_ = Teuchos::rcp(new MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>(blkColMap_->getPointMap(), X.getNumVectors()));
00228     }
00229 
00230     importedVec_->doImport(X, *importer_, REPLACE);
00231   }
00232 }
00233 
00234 //-------------------------------------------------------------------
00235 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00236 void
00237 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::updateExport(const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Y) const
00238 {
00239   if (exporter_ != Teuchos::null) {
00240     if (exportedVec_ == Teuchos::null || exportedVec_->getNumVectors() != Y.getNumVectors()) {
00241       exportedVec_ = Teuchos::rcp(new MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>(blkColMap_->getPointMap(), Y.getNumVectors()));
00242     }
00243 
00244     exportedVec_->doExport(Y, *exporter_, REPLACE);
00245   }
00246 }
00247 
00248 //-------------------------------------------------------------------
00249 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00250 void
00251 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::apply(
00252          const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> &X,
00253                MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> &Y,
00254                Teuchos::ETransp trans,
00255                Scalar alpha,
00256                Scalar beta) const
00257 {
00258   if (trans == Teuchos::NO_TRANS) {
00259     updateImport(X);
00260     const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Xref = importedVec_ == Teuchos::null ? X : *importedVec_;
00261     this->template multiply<Scalar,Scalar>(Xref, Y, trans, alpha, beta);
00262   }
00263   else if (trans == Teuchos::TRANS || trans == Teuchos::CONJ_TRANS) {
00264     updateImport(Y);
00265     MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Yref = importedVec_ == Teuchos::null ? Y : *importedVec_;
00266     this->template multiply<Scalar,Scalar>(X, Yref, trans, alpha, beta);
00267     if (importedVec_ != Teuchos::null) {
00268       Y.doExport(*importedVec_, *importer_, ADD);
00269     }
00270   }
00271 }
00272 
00273 //-------------------------------------------------------------------
00274 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00275 bool
00276 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::hasTransposeApply() const
00277 {
00278   return false;
00279 }
00280 
00281 //-------------------------------------------------------------------
00282 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00283 const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > &
00284 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::getBlockRowMap() const
00285 {
00286   return blkRowMap_;
00287 }
00288 
00289 //-------------------------------------------------------------------
00290 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00291 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > &
00292 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::getPointRowMap() const
00293 {
00294   return blkRowMap_->getPointMap();
00295 }
00296 
00297 //-------------------------------------------------------------------
00298 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00299 const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > &
00300 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::getBlockColMap() const
00301 {
00302   return blkColMap_;
00303 }
00304 
00305 //-------------------------------------------------------------------
00306 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00307 const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > &
00308 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::getPointColMap() const
00309 {
00310   return blkColMap_->getPointMap();
00311 }
00312 
00313 //-------------------------------------------------------------------
00314 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00315 void
00316 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::getGlobalBlockEntryViewNonConst(
00317     GlobalOrdinal globalBlockRow,
00318     GlobalOrdinal globalBlockCol,
00319     LocalOrdinal& numPtRows,
00320     LocalOrdinal& numPtCols,
00321     Teuchos::ArrayRCP<Scalar>& blockEntry)
00322 {
00323   //Return a non-const, read-write view of a block-entry (as an ArrayRCP),
00324   //creating/allocating the block-entry if it doesn't already exist, (but only
00325   //if fillComplete hasn't been called yet, and if the arguments numPtRows
00326   //and numPtCols have been set on input).
00327 
00328   if (is_storage_optimized_) {
00329     TEST_FOR_EXCEPTION(!isFillComplete(), std::runtime_error,
00330       "Tpetra::VbrMatrix::getGlobalBlockEntryViewNonConst internal ERROR, storage is optimized but isFillComplete() is false.");
00331 
00332     LocalOrdinal localBlockRow = blkRowMap_->getLocalBlockID(globalBlockRow);
00333     LocalOrdinal localBlockCol = blkColMap_->getLocalBlockID(globalBlockCol);
00334 
00335     getLocalBlockEntryViewNonConst(localBlockRow, localBlockCol,
00336                                    numPtRows, numPtCols, blockEntry);
00337     return;
00338   }
00339 
00340   //If we get to here, fillComplete hasn't been called yet, and the matrix data
00341   //is stored in un-packed '2D' form.
00342 
00343   if (values2D_->size() == 0) {
00344     LocalOrdinal numBlockRows = blkGraph_->getNodeNumRows();
00345     values2D_->resize(numBlockRows);
00346     col_ind_2D_global_->resize(numBlockRows);
00347   }
00348 
00349   //this acts as a range-check for globalBlockRow:
00350   LocalOrdinal localBlockRow = blkRowMap_->getLocalBlockID(globalBlockRow);
00351   TEST_FOR_EXCEPTION( localBlockRow == Teuchos::OrdinalTraits<LocalOrdinal>::invalid(),
00352      std::runtime_error,
00353      "Tpetra::VbrMatrix::getGlobalBlockEntryViewNonConst, globalBlockRow not on the local processor.");
00354 
00355   MapGlobalArrayRCP& blkrow = (*col_ind_2D_global_)[localBlockRow];
00356   typename MapGlobalArrayRCP::iterator col_iter = blkrow.find(globalBlockCol);
00357 
00358   if (col_iter != blkrow.end()) {
00359     blockEntry = col_iter->second;
00360   }
00361   else {
00362     //blockEntry doesn't already exist, so we will create it.
00363 
00364     //make sure block-size is specified:
00365     TEST_FOR_EXCEPTION(numPtRows==0 || numPtCols==0, std::runtime_error,
00366       "Tpetra::VbrMatrix::getGlobalBlockEntryViewNonConst ERROR: creating block-entry, but numPtRows and/or numPtCols is 0.");
00367 
00368     Teuchos::RCP<Node> node = getNode();
00369     size_t blockSize = numPtRows*numPtCols;
00370     blockEntry = Teuchos::arcp(new Scalar[blockSize], 0, blockSize);
00371     std::fill(blockEntry.begin(), blockEntry.end(), 0);
00372     (*values2D_)[localBlockRow].push_back(blockEntry);
00373     blkrow.insert(std::make_pair(globalBlockCol, blockEntry));
00374     blkGraph_->insertGlobalIndices(globalBlockRow, Teuchos::ArrayView<GlobalOrdinal>(&globalBlockCol, 1));
00375   }
00376 }
00377 
00378 //-------------------------------------------------------------------
00379 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00380 void
00381 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::getGlobalBlockEntryView(
00382       GlobalOrdinal globalBlockRow,
00383       GlobalOrdinal globalBlockCol,
00384       LocalOrdinal& numPtRows,
00385       LocalOrdinal& numPtCols,
00386       Teuchos::ArrayRCP<const Scalar>& blockEntry) const
00387 {
00388   //This method returns a const-view of a block-entry (as an ArrayRCP).
00389   //Throws an exception if the block-entry doesn't already exist.
00390 
00391   if (is_storage_optimized_) {
00392     TEST_FOR_EXCEPTION(!isFillComplete(), std::runtime_error,
00393       "Tpetra::VbrMatrix::getGlobalBlockEntryView internal ERROR, storage is optimized but isFillComplete() is false.");
00394 
00395     LocalOrdinal localBlockRow = blkRowMap_->getLocalBlockID(globalBlockRow);
00396     LocalOrdinal localBlockCol = blkColMap_->getLocalBlockID(globalBlockCol);
00397     getLocalBlockEntryView(localBlockRow, localBlockCol, numPtRows, numPtCols, blockEntry);
00398     return;
00399   }
00400 
00401   TEST_FOR_EXCEPTION(values2D_->size() == 0, std::runtime_error,
00402     "Tpetra::VbrMatrix::getGlobalBlockEntryView ERROR, matrix storage not yet allocated, can't return a const view.");
00403 
00404   //this acts as a range-check for globalBlockRow:
00405   LocalOrdinal localBlockRow = blkRowMap_->getLocalBlockID(globalBlockRow);
00406   TEST_FOR_EXCEPTION( localBlockRow == Teuchos::OrdinalTraits<LocalOrdinal>::invalid(),
00407      std::runtime_error,
00408      "Tpetra::VbrMatrix::getGlobalBlockEntryView, globalBlockRow not on the local processor.");
00409 
00410   MapGlobalArrayRCP& blkrow = (*col_ind_2D_global_)[localBlockRow];
00411   typename MapGlobalArrayRCP::iterator col_iter = blkrow.find(globalBlockCol);
00412 
00413   if (col_iter == blkrow.end()) {
00414     throw std::runtime_error("Tpetra::VbrMatrix::getGlobalBlockEntryView ERROR, specified block-entry doesn't exist.");
00415   }
00416 
00417   numPtRows = blkRowMap_->getLocalBlockSize(localBlockRow);
00418 
00419   TEST_FOR_EXCEPTION(numPtRows == 0, std::runtime_error,
00420     "Tpetra::VbrMatrix::getGlobalBlockEntryView ERROR, numPtRows == 0.");
00421 
00422   blockEntry = col_iter->second;
00423   numPtCols = blockEntry.size() / numPtRows;
00424 }
00425 
00426 //-------------------------------------------------------------------
00427 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00428 void
00429 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::getLocalBlockEntryViewNonConst(
00430     LocalOrdinal localBlockRow,
00431     LocalOrdinal localBlockCol,
00432     LocalOrdinal& numPtRows,                     
00433     LocalOrdinal& numPtCols,
00434     Teuchos::ArrayRCP<Scalar>& blockEntry)
00435 {
00436   typedef Teuchos::ArrayRCP<const LocalOrdinal> Host_View;
00437   typedef typename Host_View::iterator ITER;
00438   //This method returns a non-constant view of a block-entry (as an ArrayRCP).
00439 
00440   TEST_FOR_EXCEPTION(isFillComplete() == false, std::runtime_error,
00441    "Tpetra::VbrMatrix::getLocalBlockEntryViewNonConst ERROR, this method can only be called after fillComplete() has been called.");
00442 
00443   TEST_FOR_EXCEPTION(is_storage_optimized_ == false, std::runtime_error,
00444    "Tpetra::VbrMatrix::getLocalBlockEntryViewNonConst ERROR, this method can only be called if storage is optimized.");
00445 
00446   Teuchos::RCP<Node> node = getNode();
00447 
00448   Host_View bptr = node->template viewBuffer<LocalOrdinal>(2,pbuf_bptr_+localBlockRow);
00449   LocalOrdinal bindx_offset = bptr[0];
00450   LocalOrdinal length = bptr[1] - bindx_offset;
00451 
00452   TEST_FOR_EXCEPTION( length < 1, std::runtime_error,
00453     "Tpetra::VbrMatrix::getLocalBlockEntryViewNonConst ERROR, specified localBlockCol not found in localBlockRow.");
00454 
00455   Host_View bindx = node->template viewBuffer<LocalOrdinal>(length, pbuf_bindx_+bindx_offset);
00456   ITER bindx_beg = bindx.begin(), bindx_end = bindx.end();
00457   ITER it = std::lower_bound(bindx_beg, bindx_end, localBlockCol);
00458 
00459   TEST_FOR_EXCEPTION(it == bindx_end || *it != localBlockCol, std::runtime_error,
00460     "Tpetra::VbrMatrix::getLocalBlockEntryViewNonConst ERROR, specified localBlockCol not found.");
00461 
00462   numPtRows = blkRowMap_->getLocalBlockSize(localBlockRow);
00463   numPtCols = blkColMap_->getLocalBlockSize(localBlockCol);
00464 
00465   const LocalOrdinal blkSize = numPtRows*numPtCols;
00466   Host_View indx = node->template viewBuffer<LocalOrdinal>(1,pbuf_indx_+bptr[0]+(it-bindx_beg));
00467   const LocalOrdinal offset = indx[0];
00468   blockEntry = node->template viewBufferNonConst<Scalar>(Kokkos::ReadWrite, blkSize, pbuf_values1D_ + offset);
00469 }
00470 
00471 //-------------------------------------------------------------------
00472 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00473 void
00474 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::getLocalBlockEntryView(
00475       LocalOrdinal localBlockRow,
00476       LocalOrdinal localBlockCol,
00477       LocalOrdinal& numPtRows,
00478       LocalOrdinal& numPtCols,
00479       Teuchos::ArrayRCP<const Scalar>& blockEntry) const
00480 {
00481   typedef Teuchos::ArrayRCP<const LocalOrdinal> Host_View;
00482   typedef typename Host_View::iterator ITER;
00483   //This method returns a constant view of a block-entry (as an ArrayRCP).
00484 
00485   TEST_FOR_EXCEPTION(isFillComplete() == false, std::runtime_error,
00486    "Tpetra::VbrMatrix::getLocalBlockEntryView ERROR, this method can only be called after fillComplete() has been called.");
00487 
00488   TEST_FOR_EXCEPTION(is_storage_optimized_ == false, std::runtime_error,
00489    "Tpetra::VbrMatrix::getLocalBlockEntryView ERROR, this method can only be called if storage is optimized.");
00490 
00491   Teuchos::RCP<Node> node = getNode();
00492 
00493   Host_View bptr = node->template viewBuffer<LocalOrdinal>(2, pbuf_bptr_+localBlockRow);
00494   LocalOrdinal bindx_offset = bptr[0];
00495   LocalOrdinal length = bptr[1] - bindx_offset;
00496 
00497   Host_View bindx = node->template viewBuffer<LocalOrdinal>(length, pbuf_bindx_+bindx_offset);
00498   ITER bindx_beg = bindx.begin(), bindx_end = bindx.end();
00499   ITER it = std::lower_bound(bindx_beg, bindx_end, localBlockCol);
00500 
00501   TEST_FOR_EXCEPTION(it == bindx_end || *it != localBlockCol, std::runtime_error,
00502     "Tpetra::VbrMatrix::getLocalBlockEntryView ERROR, specified localBlockCol not found.");
00503 
00504   numPtRows = blkRowMap_->getLocalBlockSize(localBlockRow);
00505   numPtCols = blkRowMap_->getLocalBlockSize(localBlockCol);
00506 
00507   const LocalOrdinal blkSize = numPtRows*numPtCols;
00508   Host_View indx = node->template viewBuffer<LocalOrdinal>(1, pbuf_indx_+bptr[0]+(it-bindx_beg));
00509   const LocalOrdinal offset = indx[0];
00510   blockEntry = node->template viewBuffer<Scalar>(blkSize, pbuf_values1D_ + offset);
00511 }
00512 
00513 //-------------------------------------------------------------------
00514 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00515 Teuchos::RCP<Node>
00516 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::getNode() const
00517 {
00518   return blkRowMap_->getPointMap()->getNode();
00519 }
00520 
00521 //-------------------------------------------------------------------
00522 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00523 void
00524 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::putScalar(Scalar s)
00525 {
00526   Teuchos::RCP<Node> node = getNode();
00527 
00528   if (isFillComplete()) {
00529     fill_device_ArrayRCP(node, pbuf_values1D_, s);
00530   }
00531   else {
00532     typedef typename Teuchos::Array<Teuchos::Array<Teuchos::ArrayRCP<Scalar> > >::size_type Tsize_t;
00533     Teuchos::Array<Teuchos::Array<Teuchos::ArrayRCP<Scalar> > >& pbuf_vals2D = *values2D_;
00534     LocalOrdinal numBlockRows = blkGraph_->getNodeNumRows();
00535     for(LocalOrdinal r=0; r<numBlockRows; ++r) {
00536       for(Tsize_t c=0; c<pbuf_vals2D[r].size(); ++c) {
00537         std::fill(pbuf_vals2D[r][c].begin(), pbuf_vals2D[r][c].end(), s);
00538       }
00539     }
00540   }
00541 }
00542 
00543 //-------------------------------------------------------------------
00544 template<class ArrayType1,class ArrayType2,class Ordinal>
00545 void set_array_values(ArrayType1& dest, ArrayType2& src, Ordinal rows, Ordinal cols, Ordinal stride, Tpetra::CombineMode mode)
00546 {
00547   size_t offset = 0;
00548   size_t src_offset = 0;
00549 
00550   if (mode == ADD) {
00551     for(Ordinal col=0; col<cols; ++col) {
00552       for(Ordinal row=0; row<rows; ++row) {
00553         dest[offset++] += src[src_offset + row];
00554       }
00555       src_offset += stride;
00556     }
00557   }
00558   else {
00559     for(Ordinal col=0; col<cols; ++col) {
00560       for(Ordinal row=0; row<rows; ++row) {
00561         dest[offset++] = src[src_offset + row];
00562       }
00563       src_offset += stride;
00564     }
00565   }
00566 }
00567 
00568 //-------------------------------------------------------------------
00569 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00570 void
00571 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::setGlobalBlockEntry(GlobalOrdinal globalBlockRow, GlobalOrdinal globalBlockCol, const Teuchos::SerialDenseMatrix<GlobalOrdinal,Scalar>& blockEntry)
00572 {
00573   //first get an ArrayRCP for the internal storage for this block-entry:
00574   Teuchos::ArrayRCP<Scalar> internalBlockEntry;
00575   LocalOrdinal blkRowSize = blockEntry.numRows();
00576   LocalOrdinal blkColSize = blockEntry.numCols();
00577   getGlobalBlockEntryViewNonConst(globalBlockRow,globalBlockCol, blkRowSize, blkColSize, internalBlockEntry);
00578 
00579   //now copy the incoming block-entry into internal storage:
00580   const Scalar* inputvalues = blockEntry.values();
00581   set_array_values(internalBlockEntry, inputvalues, blkRowSize, blkColSize, blockEntry.stride(), REPLACE);
00582 }
00583 
00584 //-------------------------------------------------------------------
00585 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00586 void
00587 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::sumIntoGlobalBlockEntry(GlobalOrdinal globalBlockRow, GlobalOrdinal globalBlockCol, const Teuchos::SerialDenseMatrix<GlobalOrdinal,Scalar>& blockEntry)
00588 {
00589   //first get an ArrayRCP for the internal storage for this block-entry:
00590   Teuchos::ArrayRCP<Scalar> internalBlockEntry;
00591   LocalOrdinal blkRowSize = blockEntry.numRows();
00592   LocalOrdinal blkColSize = blockEntry.numCols();
00593   getGlobalBlockEntryViewNonConst(globalBlockRow,globalBlockCol, blkRowSize, blkColSize, internalBlockEntry);
00594 
00595   //now sum (add) the incoming block-entry into internal storage:
00596   const Scalar* inputvalues = blockEntry.values();
00597   set_array_values(internalBlockEntry, inputvalues, blkRowSize, blkColSize, blockEntry.stride(), ADD);
00598 }
00599 
00600 //-------------------------------------------------------------------
00601 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00602 void
00603 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::setGlobalBlockEntry(GlobalOrdinal globalBlockRow, GlobalOrdinal globalBlockCol, LocalOrdinal blkRowSize, LocalOrdinal blkColSize, LocalOrdinal LDA, const Teuchos::ArrayView<const Scalar>& blockEntry)
00604 {
00605   //first get an ArrayRCP for the internal storage for this block-entry:
00606   Teuchos::ArrayRCP<Scalar> internalBlockEntry;
00607   getGlobalBlockEntryViewNonConst(globalBlockRow,globalBlockCol, blkRowSize, blkColSize, internalBlockEntry);
00608 
00609   LocalOrdinal blk_size = blockEntry.size();
00610   TEST_FOR_EXCEPTION(blkColSize*LDA > blk_size, std::runtime_error,
00611     "Tpetra::VbrMatrix::setGlobalBlockEntry ERROR, inconsistent block-entry sizes.");
00612 
00613   //copy the incoming block-entry into internal storage:
00614   set_array_values(internalBlockEntry, blockEntry, blkRowSize, blkColSize, LDA, REPLACE);
00615 }
00616 
00617 //-------------------------------------------------------------------
00618 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00619 void
00620 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::sumIntoGlobalBlockEntry(GlobalOrdinal globalBlockRow, GlobalOrdinal globalBlockCol, LocalOrdinal blkRowSize, LocalOrdinal blkColSize, LocalOrdinal LDA, const Teuchos::ArrayView<const Scalar>& blockEntry)
00621 {
00622   //first get an ArrayRCP for the internal storage for this block-entry:
00623   Teuchos::ArrayRCP<Scalar> internalBlockEntry;
00624   getGlobalBlockEntryViewNonConst(globalBlockRow,globalBlockCol, blkRowSize, blkColSize, internalBlockEntry);
00625 
00626   LocalOrdinal blk_size = blockEntry.size();
00627   TEST_FOR_EXCEPTION(blkColSize*LDA > blk_size, std::runtime_error,
00628     "Tpetra::VbrMatrix::setGlobalBlockEntry ERROR, inconsistent block-entry sizes.");
00629 
00630   //copy the incoming block-entry into internal storage:
00631   set_array_values(internalBlockEntry, blockEntry, blkRowSize, blkColSize, LDA, ADD);
00632 }
00633 
00634 //-------------------------------------------------------------------
00635 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00636 void
00637 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::optimizeStorage()
00638 {
00639   if (is_storage_optimized_ == true) return;
00640 
00641   size_t num_block_rows = blkGraph_->getNodeNumRows();
00642   size_t num_block_nonzeros = blkGraph_->getNodeNumEntries();
00643 
00644   //need to count the number of point-entries:
00645   size_t num_point_entries = 0;
00646   typedef typename Teuchos::Array<Teuchos::Array<Teuchos::ArrayRCP<Scalar> > >::size_type Tsize_t;
00647   for(Tsize_t r=0; r<values2D_->size(); ++r) {
00648     Tsize_t rlen = (*values2D_)[r].size();
00649     for(Tsize_t c=0; c<rlen; ++c) {
00650       num_point_entries += (*values2D_)[r][c].size();
00651     }
00652   }
00653 
00654   const Teuchos::RCP<Node>& node = blkRowMap_->getPointMap()->getNode();
00655 
00656   //Don't change these allocation sizes unless you know what you're doing!
00657   //The '+1's are in the right place. (Trust me, I'm a trained professional.)
00658   pbuf_bptr_ = node->template allocBuffer<LocalOrdinal>(num_block_rows+1);
00659   pbuf_bindx_= node->template allocBuffer<LocalOrdinal>(num_block_nonzeros);
00660   pbuf_indx_ = node->template allocBuffer<LocalOrdinal>(num_block_nonzeros+1);
00661   pbuf_values1D_ = node->template allocBuffer<Scalar>(num_point_entries);
00662 
00663   Teuchos::ArrayRCP<LocalOrdinal> v_bptr = node->template viewBufferNonConst<LocalOrdinal>(Kokkos::WriteOnly, num_block_rows+1, pbuf_bptr_);
00664   Teuchos::ArrayRCP<LocalOrdinal> v_bindx = node->template viewBufferNonConst<LocalOrdinal>(Kokkos::WriteOnly, num_block_nonzeros, pbuf_bindx_);
00665   Teuchos::ArrayRCP<LocalOrdinal> v_indx = node->template viewBufferNonConst<LocalOrdinal>(Kokkos::WriteOnly, num_block_nonzeros+1, pbuf_indx_);
00666   Teuchos::ArrayRCP<Scalar> v_values1D = node->template viewBufferNonConst<Scalar>(Kokkos::WriteOnly, num_point_entries, pbuf_values1D_);
00667 
00668   size_t roffset = 0;
00669   size_t ioffset = 0;
00670   size_t offset = 0;
00671   Teuchos::Array<Teuchos::Array<Teuchos::ArrayRCP<Scalar> > >& pbuf_vals2D = *values2D_;
00672   for(Tsize_t r=0; r<pbuf_vals2D.size(); ++r) {
00673     LocalOrdinal rsize = blkRowMap_->getLocalBlockSize(r);
00674 
00675     v_bptr[r] = roffset;
00676     Tsize_t rlen = pbuf_vals2D[r].size();
00677     roffset += rlen;
00678 
00679     Teuchos::ArrayRCP<const LocalOrdinal> blk_row_inds = blkGraph_->getLocalRowView(r);
00680     MapGlobalArrayRCP& blk_row = (*col_ind_2D_global_)[r];
00681 
00682     for(Tsize_t c=0; c<rlen; ++c) {
00683       v_bindx[ioffset] = blk_row_inds[c];
00684       v_indx[ioffset++] = offset;
00685 
00686       LocalOrdinal csize = blkColMap_->getLocalBlockSize(blk_row_inds[c]);
00687       Tsize_t blkSize = rsize*csize;
00688 
00689       GlobalOrdinal global_col = blkGraph_->getColMap()->getGlobalElement(blk_row_inds[c]);
00690       typename MapGlobalArrayRCP::iterator iter = blk_row.find(global_col);
00691       TEST_FOR_EXCEPTION(iter == blk_row.end(), std::runtime_error,
00692         "Tpetra::VbrMatrix::optimizeStorage ERROR, global_col not found in row.");
00693 
00694       Teuchos::ArrayRCP<Scalar> vals = iter->second;
00695       for(Tsize_t n=0; n<blkSize; ++n) {
00696         v_values1D[offset+n] = vals[n];
00697       }
00698       offset += blkSize;
00699     }
00700   }
00701   v_bptr[num_block_rows] = roffset;
00702   v_indx[ioffset] = offset;
00703 
00704   //Final step: release memory for the "2D" storage:
00705   col_ind_2D_global_ = Teuchos::null;
00706   values2D_ = Teuchos::null;
00707 
00708   is_storage_optimized_ = true;
00709 }
00710 
00711 //-------------------------------------------------------------------
00712 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00713 void VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::fillLocalMatrix()
00714 {
00715   //We insist that optimzeStorage has already been called.
00716   //We don't care whether this function (fillLocalMatrix()) is being
00717   //called for the first time or not.
00718   TEST_FOR_EXCEPTION(is_storage_optimized_ != true, std::runtime_error,
00719     "Tpetra::VbrMatrix::fillLocalMatrix ERROR, optimizeStorage is required to have already been called.");
00720 
00721   lclMatrix_.setPackedValues(pbuf_values1D_,
00722                              blkRowMap_->getNodeFirstPointInBlocks_Device(),
00723                              blkColMap_->getNodeFirstPointInBlocks_Device(),
00724                              pbuf_bptr_,
00725                              pbuf_bindx_,
00726                              pbuf_indx_);
00727 }
00728 
00729 //-------------------------------------------------------------------
00730 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00731 void VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::fillLocalMatVec()
00732 {
00733   //We insist that optimzeStorage has already been called.
00734   //We don't care whether this function (fillLocalMatVec()) is being
00735   //called for the first time or not.
00736   TEST_FOR_EXCEPTION(is_storage_optimized_ != true, std::runtime_error,
00737     "Tpetra::VbrMatrix::fillLocalMatrix ERROR, optimizeStorage is required to have already been called.");
00738 
00739   lclMatVec_.initializeValues(lclMatrix_);
00740 }
00741 
00742 //-------------------------------------------------------------------
00743 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00744 void
00745 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::fillComplete(const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> >& blockDomainMap, const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> >& blockRangeMap)
00746 {
00747   if (isFillComplete()) return;
00748 
00749   blkDomainMap_ = blockDomainMap;
00750   blkRangeMap_  = blockRangeMap;
00751 
00752   blkGraph_->fillComplete(convertBlockMapToPointMap(blockDomainMap),
00753                           convertBlockMapToPointMap(blockRangeMap));
00754 
00755   //Now we need to take the point-column-map from blkGraph_ and create a
00756   //corresponding block-column-map.
00757   //Our block-column-map will have the same distribution as
00758   //blkGraph_->getColMap, and we'll get block-size info from the
00759   //blkDomainMap_. This will require some communication in cases
00760   //where blkDomainMap is distributed differently.
00761   blkColMap_ = makeBlockColumnMap(blkDomainMap_, blkGraph_->getColMap());
00762 
00763   if (!blkDomainMap_->getPointMap()->isSameAs(*blkColMap_->getPointMap())) {
00764     importer_ = Teuchos::rcp(new Tpetra::Import<LocalOrdinal,GlobalOrdinal,Node>(blkDomainMap_->getPointMap(), blkColMap_->getPointMap()));
00765   }
00766   if (!blkRangeMap_->getPointMap()->isSameAs(*blkRowMap_->getPointMap())) {
00767     exporter_ = Teuchos::rcp(new Tpetra::Export<LocalOrdinal,GlobalOrdinal,Node>(blkRowMap_->getPointMap(), blkRangeMap_->getPointMap()));
00768   }
00769 
00770   is_fill_completed_ = true;
00771 
00772   optimizeStorage();
00773 
00774   fillLocalMatrix();
00775   fillLocalMatVec();
00776 }
00777 
00778 //-------------------------------------------------------------------
00779 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00780 void VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::fillComplete()
00781 {
00782   //In this case, our block-row-map will also be our domain-map and
00783   //range-map.
00784 
00785   fillComplete(getBlockRowMap(), getBlockRowMap());
00786 }
00787 
00788 //-------------------------------------------------------------------
00789 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00790 std::string
00791 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::description() const
00792 {
00793   std::ostringstream oss;
00794   oss << Teuchos::Describable::description();
00795   if (isFillComplete()) {
00796     oss << "{status = fill complete, global num block rows = " << blkRowMap_->getGlobalNumBlocks() << "}";
00797   }
00798   else {
00799     oss << "{status = fill not complete, global num block rows = " << blkRowMap_->getGlobalNumBlocks() << "}";
00800   }
00801   return oss.str();
00802 }
00803 
00804 //-------------------------------------------------------------------
00805 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class LocalMatOps>
00806 void
00807 VbrMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>::describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel) const
00808 {
00809   Teuchos::EVerbosityLevel vl = verbLevel;
00810   if (vl == Teuchos::VERB_DEFAULT) vl = Teuchos::VERB_LOW;
00811 
00812   if (vl == Teuchos::VERB_NONE) return;
00813 
00814   Teuchos::RCP<const Teuchos::Comm<int> > comm = blkRowMap_->getPointMap()->getComm();
00815   const int myProc = comm->getRank();
00816 
00817   if (myProc == 0) out << "VbrMatrix::describe is under construction..." << std::endl;
00818 }
00819 
00820 }//namespace Tpetra
00821 
00822 //
00823 // Explicit instantiation macro
00824 //
00825 // Must be expanded from within the Tpetra namespace!
00826 //
00827 
00828 #define TPETRA_VBRMATRIX_INSTANT(SCALAR,LO,GO,NODE) \
00829   \
00830   template class VbrMatrix< SCALAR , LO , GO , NODE >;
00831 
00832 #endif //TPETRA_VBRMATRIX_DEF_HPP
00833 
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines
Generated on Wed Apr 13 10:21:41 2011 for Tpetra Matrix/Vector Services by  doxygen 1.6.3