Tpetra Matrix/Vector Services Version of the Day
Tpetra_CrsMatrix_decl.hpp
00001 // @HEADER
00002 // ***********************************************************************
00003 //
00004 //          Tpetra: Templated Linear Algebra Services Package
00005 //                 Copyright (2008) Sandia Corporation
00006 //
00007 // Under the terms of Contract DE-AC04-94AL85000 with Sandia Corporation,
00008 // the U.S. Government retains certain rights in this software.
00009 //
00010 // Redistribution and use in source and binary forms, with or without
00011 // modification, are permitted provided that the following conditions are
00012 // met:
00013 //
00014 // 1. Redistributions of source code must retain the above copyright
00015 // notice, this list of conditions and the following disclaimer.
00016 //
00017 // 2. Redistributions in binary form must reproduce the above copyright
00018 // notice, this list of conditions and the following disclaimer in the
00019 // documentation and/or other materials provided with the distribution.
00020 //
00021 // 3. Neither the name of the Corporation nor the names of the
00022 // contributors may be used to endorse or promote products derived from
00023 // this software without specific prior written permission.
00024 //
00025 // THIS SOFTWARE IS PROVIDED BY SANDIA CORPORATION "AS IS" AND ANY
00026 // EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
00027 // IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
00028 // PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL SANDIA CORPORATION OR THE
00029 // CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
00030 // EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
00031 // PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
00032 // PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY OF
00033 // LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT (INCLUDING
00034 // NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF THIS
00035 // SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
00036 //
00037 // Questions? Contact Michael A. Heroux (maherou@sandia.gov)
00038 //
00039 // ************************************************************************
00040 // @HEADER
00041 
00042 #ifndef TPETRA_CRSMATRIX_DECL_HPP
00043 #define TPETRA_CRSMATRIX_DECL_HPP
00044 
00045 // TODO: row-wise insertion of entries in globalAssemble() may be more efficient
00046 
00047 // TODO: add typeglobs: CrsMatrix<Scalar,typeglob>
00048 // TODO: add template (template) parameter for nonlocal container (this will be part of typeglob)
00049 
00050 #include <Kokkos_DefaultNode.hpp>
00051 #include <Kokkos_DefaultKernels.hpp>
00052 
00053 #include "Tpetra_ConfigDefs.hpp"
00054 #include "Tpetra_RowMatrix.hpp"
00055 #include "Tpetra_DistObject.hpp"
00056 #include "Tpetra_CrsGraph.hpp"
00057 #include "Tpetra_Vector.hpp"
00058 #include "Tpetra_CrsMatrixMultiplyOp_decl.hpp"
00059 
00060 #ifndef DOXYGEN_SHOULD_SKIP_THIS
00061 namespace Tpetra {
00062   // struct for i,j,v triplets
00063   template <class Ordinal, class Scalar>
00064   struct CrsIJV {
00065     CrsIJV();
00066     CrsIJV(Ordinal row, Ordinal col, const Scalar &val);
00067     Ordinal i,j;
00068     Scalar  v;
00069   };
00070 }
00071 
00072 namespace Teuchos {
00073   // SerializationTraits specialization for CrsIJV, using DirectSerialization
00074   template <typename Ordinal, typename Scalar>
00075   class SerializationTraits<int,Tpetra::CrsIJV<Ordinal,Scalar> >
00076   : public DirectSerializationTraits<int,Tpetra::CrsIJV<Ordinal,Scalar> >
00077   {};
00078 }
00079 
00080 namespace std {
00081   template <class Ordinal, class Scalar>
00082   bool operator<(const Tpetra::CrsIJV<Ordinal,Scalar> &ijv1, const Tpetra::CrsIJV<Ordinal,Scalar> &ijv2);
00083 }
00084 #endif
00085 
00086 namespace Tpetra {
00087 
00089 
00180   template <class Scalar,
00181             class LocalOrdinal  = int,
00182             class GlobalOrdinal = LocalOrdinal,
00183             class Node          = Kokkos::DefaultNode::DefaultNodeType,
00184             class LocalMatOps   = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::SparseOps >
00185   class CrsMatrix : public RowMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node>,
00186                     public DistObject<char, LocalOrdinal,GlobalOrdinal,Node> {
00187   public:
00188     typedef Scalar                                scalar_type;
00189     typedef LocalOrdinal                          local_ordinal_type;
00190     typedef GlobalOrdinal                         global_ordinal_type;
00191     typedef Node                                  node_type;
00192     typedef Map<LocalOrdinal,GlobalOrdinal,Node>  map_type;
00193     // backwards compatibility defines both of these
00194     typedef LocalMatOps   mat_vec_type;
00195     typedef LocalMatOps   mat_solve_type;
00196 
00198 
00199 
00217     CrsMatrix (const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& rowMap,
00218                size_t maxNumEntriesPerRow,
00219                ProfileType pftype = DynamicProfile,
00220                const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null);
00221 
00239     CrsMatrix (const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& rowMap,
00240                const ArrayRCP<const size_t>& NumEntriesPerRowToAlloc,
00241                ProfileType pftype = DynamicProfile,
00242                const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null);
00243 
00266     CrsMatrix (const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& rowMap,
00267                const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& colMap,
00268                size_t maxNumEntriesPerRow,
00269                ProfileType pftype = DynamicProfile,
00270                const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null);
00271 
00294     CrsMatrix (const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& rowMap,
00295                const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> >& colMap,
00296                const ArrayRCP<const size_t>& NumEntriesPerRowToAlloc,
00297                ProfileType pftype = DynamicProfile,
00298                const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null);
00299 
00320     explicit CrsMatrix (const Teuchos::RCP<const CrsGraph<LocalOrdinal,GlobalOrdinal,Node,LocalMatOps> >& graph,
00321                         const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null);
00322 
00324     virtual ~CrsMatrix();
00325 
00327 
00328 
00329 
00331 
00340     void insertGlobalValues(GlobalOrdinal globalRow, const ArrayView<const GlobalOrdinal> &cols, const ArrayView<const Scalar> &vals);
00341 
00343 
00354     void insertLocalValues(LocalOrdinal localRow,
00355                            const ArrayView<const LocalOrdinal> &cols,
00356                            const ArrayView<const Scalar> &vals);
00357 
00359 
00364     void replaceGlobalValues(GlobalOrdinal globalRow,
00365                              const ArrayView<const GlobalOrdinal> &cols,
00366                              const ArrayView<const Scalar>        &vals);
00367 
00369 
00371     void replaceLocalValues(LocalOrdinal localRow,
00372                             const ArrayView<const LocalOrdinal> &cols,
00373                             const ArrayView<const Scalar>       &vals);
00374 
00376 
00381     void sumIntoGlobalValues(GlobalOrdinal globalRow,
00382                              const ArrayView<const GlobalOrdinal> &cols,
00383                              const ArrayView<const Scalar>        &vals);
00384 
00385 
00387 
00392     void sumIntoLocalValues(LocalOrdinal globalRow,
00393                             const ArrayView<const LocalOrdinal>  &cols,
00394                             const ArrayView<const Scalar>        &vals);
00395 
00397     void setAllToScalar(const Scalar &alpha);
00398 
00400     void scale(const Scalar &alpha);
00401 
00403 
00404 
00405 
00412     void globalAssemble();
00413 
00422     void resumeFill(const RCP<ParameterList> &params = null);
00423 
00434     void fillComplete(const RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > &domainMap,
00435                       const RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > &rangeMap,
00436                       const RCP<ParameterList> &params = null);
00437 
00450     void fillComplete(const RCP<ParameterList> &params = null);
00451 
00453 
00454 
00455 
00457     const RCP<const Comm<int> > & getComm() const;
00458 
00460     RCP<Node> getNode() const;
00461 
00463     const RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > & getRowMap() const;
00464 
00466     const RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > & getColMap() const;
00467 
00469     RCP<const RowGraph<LocalOrdinal,GlobalOrdinal,Node> > getGraph() const;
00470 
00472     RCP<const CrsGraph<LocalOrdinal,GlobalOrdinal,Node,LocalMatOps> > getCrsGraph() const;
00473 
00493     global_size_t getGlobalNumRows() const;
00494 
00499     global_size_t getGlobalNumCols() const;
00500 
00502     size_t getNodeNumRows() const;
00503 
00505 
00507     size_t getNodeNumCols() const;
00508 
00510     GlobalOrdinal getIndexBase() const;
00511 
00513     global_size_t getGlobalNumEntries() const;
00514 
00516     size_t getNodeNumEntries() const;
00517 
00519 
00520     size_t getNumEntriesInGlobalRow(GlobalOrdinal globalRow) const;
00521 
00523 
00524     size_t getNumEntriesInLocalRow(LocalOrdinal localRow) const;
00525 
00527 
00529     global_size_t getGlobalNumDiags() const;
00530 
00532 
00534     size_t getNodeNumDiags() const;
00535 
00537 
00539     size_t getGlobalMaxNumRowEntries() const;
00540 
00542 
00544     size_t getNodeMaxNumRowEntries() const;
00545 
00547     bool hasColMap() const;
00548 
00550 
00552     bool isLowerTriangular() const;
00553 
00555 
00557     bool isUpperTriangular() const;
00558 
00560     bool isLocallyIndexed() const;
00561 
00563     bool isGloballyIndexed() const;
00564 
00566     bool isFillComplete() const;
00567 
00569     bool isFillActive() const;
00570 
00572 
00578     bool isStorageOptimized() const;
00579 
00581     ProfileType getProfileType() const;
00582 
00584     bool isStaticGraph() const;
00585 
00587 
00593     typename ScalarTraits<Scalar>::magnitudeType getFrobeniusNorm() const;
00594 
00596 
00606     void getGlobalRowCopy(GlobalOrdinal GlobalRow,
00607                           const ArrayView<GlobalOrdinal> &Indices,
00608                           const ArrayView<Scalar> &Values,
00609                           size_t &NumEntries
00610                           ) const;
00611 
00613 
00625     void getLocalRowCopy(LocalOrdinal LocalRow,
00626                          const ArrayView<LocalOrdinal> &Indices,
00627                          const ArrayView<Scalar> &Values,
00628                          size_t &NumEntries
00629                          ) const;
00630 
00632 
00641     void getGlobalRowView(GlobalOrdinal GlobalRow, ArrayView<const GlobalOrdinal> &indices, ArrayView<const Scalar> &values) const;
00642 
00644 
00653     void getLocalRowView(LocalOrdinal LocalRow, ArrayView<const LocalOrdinal> &indices, ArrayView<const Scalar> &values) const;
00654 
00656 
00658     void getLocalDiagCopy(Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node> &diag) const;
00659 
00661     void leftScale(const Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node>& x);
00662 
00664     void rightScale(const Vector<Scalar, LocalOrdinal, GlobalOrdinal, Node>& x);
00665 
00667 
00668 
00669 
00671 
00681     template <class DomainScalar, class RangeScalar>
00682     void
00683     localMultiply (const MultiVector<DomainScalar,LocalOrdinal,GlobalOrdinal,Node>& X,
00684                    MultiVector<RangeScalar,LocalOrdinal,GlobalOrdinal,Node>& Y,
00685                    Teuchos::ETransp trans,
00686                    RangeScalar alpha,
00687                    RangeScalar beta) const;
00688 
00705     template <class DomainScalar, class RangeScalar>
00706     void
00707     localSolve (const MultiVector<RangeScalar,LocalOrdinal,GlobalOrdinal,Node>& Y,
00708                 MultiVector<DomainScalar,LocalOrdinal,GlobalOrdinal,Node>& X,
00709                 Teuchos::ETransp trans) const;
00710 
00712     template <class T>
00713     RCP<CrsMatrix<T,LocalOrdinal,GlobalOrdinal,Node> > convert() const;
00714 
00716 
00717 
00718 
00720 
00723     void apply(const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> & X, MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> &Y,
00724                Teuchos::ETransp mode = Teuchos::NO_TRANS,
00725                Scalar alpha = ScalarTraits<Scalar>::one(),
00726                Scalar beta = ScalarTraits<Scalar>::zero()) const;
00727 
00729     bool hasTransposeApply() const;
00730 
00733     const RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > & getDomainMap() const;
00734 
00737     const RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > & getRangeMap() const;
00738 
00740 
00741 
00742 
00744     std::string description() const;
00745 
00747     void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const;
00748 
00750 
00751 
00752 
00753     bool checkSizes(const DistObject<char, LocalOrdinal,GlobalOrdinal,Node>& source);
00754 
00755     void copyAndPermute(const DistObject<char, LocalOrdinal,GlobalOrdinal,Node>& source,
00756                         size_t numSameIDs,
00757                         const ArrayView<const LocalOrdinal> &permuteToLIDs,
00758                         const ArrayView<const LocalOrdinal> &permuteFromLIDs);
00759 
00760     void packAndPrepare(const DistObject<char, LocalOrdinal,GlobalOrdinal,Node>& source,
00761                         const ArrayView<const LocalOrdinal> &exportLIDs,
00762                         Array<char> &exports,
00763                         const ArrayView<size_t> & numPacketsPerLID,
00764                         size_t& constantNumPackets,
00765                         Distributor &distor);
00766 
00776     void
00777     unpackAndCombine (const Teuchos::ArrayView<const LocalOrdinal> &importLIDs,
00778                       const Teuchos::ArrayView<const char> &imports,
00779                       const Teuchos::ArrayView<size_t> &numPacketsPerLID,
00780                       size_t constantNumPackets,
00781                       Distributor &distor,
00782                       CombineMode combineMode);
00784 
00785   private:
00786     // We forbid copy construction by declaring this method private
00787     // and not implementing it.
00788     CrsMatrix (const CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps> &rhs);
00789 
00790     // We forbid assignment (operator=) by declaring this method
00791     // private and not implementing it.
00792     CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>&
00793     operator= (const CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps> &rhs);
00794 
00802     void
00803     combineGlobalValues (const GlobalOrdinal globalRowIndex,
00804                          const Teuchos::ArrayView<const GlobalOrdinal> columnIndices,
00805                          const Teuchos::ArrayView<const Scalar> values,
00806                          const Tpetra::CombineMode combineMode);
00807 
00822     template<class BinaryFunction>
00823     void
00824     transformGlobalValues (GlobalOrdinal globalRow,
00825                            const Teuchos::ArrayView<const GlobalOrdinal>& indices,
00826                            const Teuchos::ArrayView<const Scalar>        & values,
00827                            BinaryFunction f)
00828     {
00829       typedef LocalOrdinal LO;
00830       typedef GlobalOrdinal GO;
00831       typedef Node NT;
00832       using Teuchos::Array;
00833       using Teuchos::ArrayView;
00834 
00835       TEUCHOS_TEST_FOR_EXCEPTION(values.size() != indices.size(),
00836         std::invalid_argument, "transformGlobalValues: values.size() = "
00837         << values.size() << " != indices.size() = " << indices.size() << ".");
00838 
00839       const LO lrow = this->getRowMap()->getLocalElement(globalRow);
00840 
00841       TEUCHOS_TEST_FOR_EXCEPTION(lrow == LOT::invalid(), std::invalid_argument,
00842         "transformGlobalValues: The given global row index " << globalRow
00843         << " is not owned by the calling process (rank "
00844         << this->getRowMap()->getComm()->getRank() << ").");
00845 
00846       RowInfo rowInfo = staticGraph_->getRowInfo(lrow);
00847       if (indices.size() > 0) {
00848         if (isLocallyIndexed()) {
00849           // Convert global indices to local indices.
00850           const Map<LO, GO, NT> &colMap = *(this->getColMap());
00851           Array<LO> lindices (indices.size());
00852           typename ArrayView<const GO>::iterator gindit = indices.begin();
00853           typename Array<LO>::iterator           lindit = lindices.begin();
00854           while (gindit != indices.end()) {
00855             // No need to filter before asking the column Map to
00856             // convert GID->LID.  If the GID doesn't exist in the
00857             // column Map, the GID will be mapped to invalid(), which
00858             // will not be found in the graph.
00859             *lindit++ = colMap.getLocalElement(*gindit++);
00860           }
00861           typename Graph::SLocalGlobalViews inds_view;
00862           inds_view.linds = lindices();
00863           staticGraph_->template transformValues<LocalIndices>(rowInfo, inds_view, this->getViewNonConst(rowInfo).begin(), values.begin(), f);
00864         }
00865         else if (isGloballyIndexed()) {
00866           typename Graph::SLocalGlobalViews inds_view;
00867           inds_view.ginds = indices;
00868           staticGraph_->template transformValues<GlobalIndices>(rowInfo, inds_view, this->getViewNonConst(rowInfo).begin(), values.begin(), f);
00869         }
00870       }
00871     }
00872 
00873   protected:
00874     // useful typedefs
00875     typedef OrdinalTraits<LocalOrdinal>                     LOT;
00876     typedef OrdinalTraits<GlobalOrdinal>                    GOT;
00877     typedef ScalarTraits<Scalar>                             ST;
00878     typedef typename ST::magnitudeType                Magnitude;
00879     typedef MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>       MV;
00880     typedef Vector<Scalar,LocalOrdinal,GlobalOrdinal,Node>             V;
00881     typedef CrsGraph<LocalOrdinal,GlobalOrdinal,Node,LocalMatOps>  Graph;
00882     typedef typename LocalMatOps::template bind_scalar<Scalar>::other_type                    sparse_ops_type;
00883     typedef typename sparse_ops_type::template graph<LocalOrdinal,Node>::graph_type          local_graph_type;
00884     typedef typename sparse_ops_type::template matrix<Scalar,LocalOrdinal,Node>::matrix_type local_matrix_type;
00885     // Enums
00886     enum GraphAllocationStatus {
00887       GraphAlreadyAllocated,
00888       GraphNotYetAllocated
00889     };
00890 
00907     void allocateValues (ELocalGlobal lg, GraphAllocationStatus gas);
00908 
00914     void sortEntries();
00915 
00921     void mergeRedundantEntries();
00922 
00930     void clearGlobalConstants();
00931 
00940     void computeGlobalConstants();
00941 
00942     // matrix data accessors
00943     ArrayView<const Scalar>    getView(RowInfo rowinfo) const;
00944     ArrayView<      Scalar>    getViewNonConst(RowInfo rowinfo);
00945     // local Kokkos objects
00946     void fillLocalMatrix(const RCP<ParameterList> &params);
00947     void fillLocalGraphAndMatrix(const RCP<ParameterList> &params);
00948     // debugging
00949     void checkInternalState() const;
00950 
00962 
00963     RCP<const Graph> staticGraph_;
00964     RCP<      Graph>     myGraph_;
00966 
00967     RCP<sparse_ops_type>   lclMatOps_;
00968     RCP<local_matrix_type> lclMatrix_;
00969 
00970     // matrix values. before allocation, both are null.
00971     // after allocation, one is null.
00972     // 1D == StaticAllocation, 2D == DynamicAllocation
00973     // The allocation always matches that of graph_, as the graph does the allocation for the matrix.
00974     ArrayRCP<Scalar>                       values1D_;
00975     ArrayRCP<ArrayRCP<Scalar> >            values2D_;
00976     // TODO: these could be allocated at resumeFill() and de-allocated at fillComplete() to make for very fast getView()/getViewNonConst()
00977     // ArrayRCP< typedef ArrayRCP<const Scalar>::iterator > rowPtrs_;
00978     // ArrayRCP< typedef ArrayRCP<      Scalar>::iterator > rowPtrsNC_;
00979 
00981     bool fillComplete_;
00982 
00995     std::map<GlobalOrdinal, Array<std::pair<GlobalOrdinal,Scalar> > > nonlocals_;
00996 
00997     // a wrapper around multiply, for use in apply; it contains a non-owning RCP to *this, therefore, it is not allowed
00998     // to persist past the destruction of *this. therefore, WE MAY NOT SHARE THIS POINTER.
00999     RCP< const CrsMatrixMultiplyOp<Scalar,Scalar,LocalOrdinal,GlobalOrdinal,Node,LocalMatOps> > sameScalarMultiplyOp_;
01000 
01001     // cached frobenius norm: -ST::one() means invalid
01002     mutable Magnitude frobNorm_;
01003 
01005     bool insertGlobalValuesWarnedEfficiency_;
01007     bool insertLocalValuesWarnedEfficiency_;
01008   }; // class CrsMatrix
01009 
01016   template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
01017   Teuchos::RCP<CrsMatrix<Scalar,LocalOrdinal,GlobalOrdinal,Node> >
01018   createCrsMatrix (const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > &map,
01019                    size_t maxNumEntriesPerRow = 0,
01020                    const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null)
01021   {
01022     using Teuchos::rcp;
01023     typedef CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node> matrix_type;
01024 
01025     return rcp (new matrix_type (map, maxNumEntriesPerRow, DynamicProfile, params));
01026   }
01027 
01077   template<class CrsMatrixType>
01078   Teuchos::RCP<CrsMatrixType>
01079   importAndFillCompleteCrsMatrix (const Teuchos::RCP<const CrsMatrixType>& sourceMatrix,
01080                                   const Import<typename CrsMatrixType::local_ordinal_type,
01081                                                typename CrsMatrixType::global_ordinal_type,
01082                                                typename CrsMatrixType::node_type>& importer,
01083                                   const Teuchos::RCP<const Map<typename CrsMatrixType::local_ordinal_type,
01084                                                                typename CrsMatrixType::global_ordinal_type,
01085                                                                typename CrsMatrixType::node_type> >& domainMap = Teuchos::null,
01086                                   const Teuchos::RCP<const Map<typename CrsMatrixType::local_ordinal_type,
01087                                                                typename CrsMatrixType::global_ordinal_type,
01088                                                                typename CrsMatrixType::node_type> >& rangeMap = Teuchos::null,
01089                                   const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null)
01090   {
01091     using Teuchos::as;
01092     using Teuchos::RCP;
01093     using Teuchos::rcp;
01094     typedef Map<typename CrsMatrixType::local_ordinal_type,
01095       typename CrsMatrixType::global_ordinal_type,
01096       typename CrsMatrixType::node_type> map_type;
01097 
01098     // FIXME (mfh 11 Apr 2012) The current implementation of this
01099     // method doesn't actually fuse the Import with fillComplete().
01100     // This will change in the future.
01101     RCP<CrsMatrixType> destMat =
01102       rcp (new CrsMatrixType (importer.getTargetMap (),
01103                               as<size_t> (0),
01104                               DynamicProfile,
01105                               params));
01106     destMat->doImport (*sourceMatrix, importer, INSERT);
01107 
01108     // Use the source matrix's domain Map as the default.
01109     RCP<const map_type> theDomainMap =
01110       domainMap.is_null () ? sourceMatrix->getDomainMap () : domainMap;
01111     // Use the source matrix's range Map as the default.
01112     RCP<const map_type> theRangeMap =
01113       rangeMap.is_null () ? sourceMatrix->getRangeMap () : rangeMap;
01114 
01115     destMat->fillComplete (theDomainMap, theRangeMap);
01116     return destMat;
01117   }
01118 
01152   template<class CrsMatrixType>
01153   Teuchos::RCP<CrsMatrixType>
01154   exportAndFillCompleteCrsMatrix (const Teuchos::RCP<const CrsMatrixType>& sourceMatrix,
01155                                   const Export<typename CrsMatrixType::local_ordinal_type,
01156                                                typename CrsMatrixType::global_ordinal_type,
01157                                                typename CrsMatrixType::node_type>& exporter,
01158                                   const Teuchos::RCP<const Map<typename CrsMatrixType::local_ordinal_type,
01159                                                                typename CrsMatrixType::global_ordinal_type,
01160                                                                typename CrsMatrixType::node_type> >& domainMap = Teuchos::null,
01161                                   const Teuchos::RCP<const Map<typename CrsMatrixType::local_ordinal_type,
01162                                                                typename CrsMatrixType::global_ordinal_type,
01163                                                                typename CrsMatrixType::node_type> >& rangeMap = Teuchos::null,
01164                                   const Teuchos::RCP<Teuchos::ParameterList>& params = Teuchos::null)
01165   {
01166     using Teuchos::as;
01167     using Teuchos::RCP;
01168     using Teuchos::rcp;
01169     typedef Map<typename CrsMatrixType::local_ordinal_type,
01170       typename CrsMatrixType::global_ordinal_type,
01171       typename CrsMatrixType::node_type> map_type;
01172 
01173     // FIXME (mfh 11 Apr 2012) The current implementation of this
01174     // method doesn't actually fuse the Export with fillComplete().
01175     // This will change in the future.
01176     RCP<CrsMatrixType> destMat =
01177       rcp (new CrsMatrixType (exporter.getTargetMap (),
01178                               as<size_t> (0),
01179                               DynamicProfile,
01180                               params));
01181     destMat->doExport (*sourceMatrix, exporter, INSERT);
01182 
01183     // Use the source matrix's domain Map as the default.
01184     RCP<const map_type> theDomainMap =
01185       domainMap.is_null () ? sourceMatrix->getDomainMap () : domainMap;
01186     // Use the source matrix's range Map as the default.
01187     RCP<const map_type> theRangeMap =
01188       rangeMap.is_null () ? sourceMatrix->getRangeMap () : rangeMap;
01189 
01190     destMat->fillComplete (theDomainMap, theRangeMap);
01191     return destMat;
01192   }
01193 } // namespace Tpetra
01194 
01205 #endif
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines