Tpetra_VbrMatrix_decl.hpp

Go to the documentation of this file.
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_DECL_HPP
00030 #define TPETRA_VBRMATRIX_DECL_HPP
00031 
00032 #include <Kokkos_DefaultNode.hpp>
00033 #include <Kokkos_DefaultKernels.hpp>
00034 #include <Kokkos_VbrMatrix.hpp>
00035 
00036 #include <Teuchos_ScalarTraits.hpp>
00037 #include <Teuchos_OrdinalTraits.hpp>
00038 #include <Teuchos_SerializationTraits.hpp>
00039 
00040 #include "Tpetra_ConfigDefs.hpp"
00041 #include "Tpetra_Operator.hpp"
00042 #include "Tpetra_BlockMap.hpp"
00043 #include "Tpetra_CrsGraph.hpp"
00044 
00049 namespace Tpetra {
00050 
00052 
00080 template <class Scalar, 
00081           class LocalOrdinal  = int, 
00082           class GlobalOrdinal = LocalOrdinal, 
00083           class Node          = Kokkos::DefaultNode::DefaultNodeType, 
00084           class LocalMatOps   = typename Kokkos::DefaultKernels<Scalar,LocalOrdinal,Node>::BlockSparseOps >
00085 class VbrMatrix : public Tpetra::Operator<Scalar,LocalOrdinal,GlobalOrdinal,Node> {
00086  public:
00087   typedef Scalar        scalar_type;
00088   typedef LocalOrdinal  local_ordinal_type;
00089   typedef GlobalOrdinal global_ordinal_type;
00090   typedef Node          node_type;
00091   typedef LocalMatOps   mat_vec_type;
00092 
00094 
00095 
00097 
00102   VbrMatrix(const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > &blkRowMap, size_t maxNumEntriesPerRow, ProfileType pftype = DynamicProfile);
00103 
00105 
00122 //  VbrMatrix(const Teuchos::RCP<const CrsGraph<LocalOrdinal,GlobalOrdinal,Node> >& blkGraph, const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> >& blkDomainMap, const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> >& blkRangeMap);
00123 
00125   virtual ~VbrMatrix();
00126 
00128 
00130 
00132 
00137   template <class DomainScalar, class RangeScalar>
00138       void multiply(const MultiVector<DomainScalar,LocalOrdinal,GlobalOrdinal,Node> & X, MultiVector<RangeScalar,LocalOrdinal,GlobalOrdinal,Node> &Y, Teuchos::ETransp trans, RangeScalar alpha, RangeScalar beta) const;
00139 
00141 
00143 
00144 
00146 
00148     const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > & getDomainMap() const;
00149 
00151 
00153     const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > & getRangeMap() const;
00154 
00156 
00161     void apply(const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> &X,
00162                MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> &Y,
00163                Teuchos::ETransp trans = Teuchos::NO_TRANS,
00164                Scalar alpha = Teuchos::ScalarTraits<Scalar>::one(),
00165                Scalar beta = Teuchos::ScalarTraits<Scalar>::zero()) const;
00166 
00168     bool hasTransposeApply() const;
00169 
00171 
00173 
00174 
00176   const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > & getBlockRowMap() const;
00177 
00179   const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > & getBlockColMap() const;
00180 
00182   const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > & getPointRowMap() const;
00183 
00185   const Teuchos::RCP<const Map<LocalOrdinal,GlobalOrdinal,Node> > & getPointColMap() const;
00186 
00188   bool isFillComplete() const;
00190 
00192 
00193 
00195 
00198   void putScalar(Scalar s);
00199 
00201 
00210   void setGlobalBlockEntry(GlobalOrdinal globalBlockRow, GlobalOrdinal globalBlockCol, const Teuchos::SerialDenseMatrix<GlobalOrdinal,Scalar>& blockEntry);
00211 
00213 
00222   void sumIntoGlobalBlockEntry(GlobalOrdinal globalBlockRow, GlobalOrdinal globalBlockCol, const Teuchos::SerialDenseMatrix<GlobalOrdinal,Scalar>& blockEntry);
00223 
00225 
00234   void setGlobalBlockEntry(GlobalOrdinal globalBlockRow, GlobalOrdinal globalBlockCol, LocalOrdinal blkRowSize, LocalOrdinal blkColSize, LocalOrdinal LDA, const Teuchos::ArrayView<const Scalar>& blockEntry);
00235 
00237 
00246   void sumIntoGlobalBlockEntry(GlobalOrdinal globalBlockRow, GlobalOrdinal globalBlockCol, LocalOrdinal blkRowSize, LocalOrdinal blkColSize, LocalOrdinal LDA, const Teuchos::ArrayView<const Scalar>& blockEntry);
00247 
00249 
00251 
00252 
00253   void fillComplete(const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> >& blockDomainMap, const Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> >& blockRangeMap);
00254 
00255   void fillComplete();
00257 
00259 
00260 
00262 
00270   void getGlobalBlockEntryView(GlobalOrdinal globalBlockRow,
00271                                GlobalOrdinal globalBlockCol,
00272                                LocalOrdinal& numPtRows,
00273                                LocalOrdinal& numPtCols,
00274                                Teuchos::ArrayRCP<const Scalar>& blockEntry) const;
00275 
00277 
00287   void getGlobalBlockEntryViewNonConst(GlobalOrdinal globalBlockRow,
00288                                        GlobalOrdinal globalBlockCol,
00289                                        LocalOrdinal& numPtRows,
00290                                        LocalOrdinal& numPtCols,
00291                                        Teuchos::ArrayRCP<Scalar>& blockEntry);
00292 
00294 
00304   void getLocalBlockEntryView(LocalOrdinal localBlockRow,
00305                               LocalOrdinal localBlockCol,
00306                               LocalOrdinal& numPtRows, 
00307                               LocalOrdinal& numPtCols,
00308                               Teuchos::ArrayRCP<const Scalar>& blockEntry) const;
00309 
00311 
00327   void getLocalBlockEntryViewNonConst(LocalOrdinal localBlockRow,
00328                                       LocalOrdinal localBlockCol,
00329                                       LocalOrdinal& numPtRows,
00330                                       LocalOrdinal& numPtCols,
00331                                       Teuchos::ArrayRCP<Scalar>& blockEntry);
00332 
00334 
00336 
00337   std::string description() const;
00338 
00341   void describe(Teuchos::FancyOStream &out, const Teuchos::EVerbosityLevel verbLevel=Teuchos::Describable::verbLevel_default) const;
00343 
00344  private:
00345   //private methods:
00346 
00347   Teuchos::RCP<Node> getNode() const;
00348 
00349   void updateImport(const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& X) const;
00350   void updateExport(const MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node>& Y) const;
00351 
00352   void optimizeStorage();
00353   void fillLocalMatrix();
00354   void fillLocalMatVec();
00355 
00356   //private data members:
00357 
00358   Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > blkRowMap_;
00359   Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > blkColMap_;
00360   Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > blkDomainMap_;
00361   Teuchos::RCP<const BlockMap<LocalOrdinal,GlobalOrdinal,Node> > blkRangeMap_;
00362 
00363   Teuchos::RCP<CrsGraph<LocalOrdinal,GlobalOrdinal,Node> > blkGraph_;
00364   Kokkos::VbrMatrix<Scalar,LocalOrdinal,Node> lclMatrix_;
00365 
00366   //It takes 6 arrays to adequately represent a variable-block-row
00367   //matrix in packed (contiguous storage) form. For a description of these
00368   //arrays, see the text at the bottom of this file.
00369   //(Note that 2 of those arrays, rptr and cptr, are represented by arrays in the
00370   //blkRowMap_ and blkColMap_ objects.)
00371   //
00372   //These arrays are handled as if they may point to memory that resides on
00373   //a separate device (e.g., a GPU). In other words, when the contents of these
00374   //arrays are manipulated, we use views or buffers obtained from the Node object.
00375   Teuchos::ArrayRCP<Scalar> pbuf_values1D_;
00376   Teuchos::ArrayRCP<LocalOrdinal> pbuf_bptr_;
00377   Teuchos::ArrayRCP<LocalOrdinal> pbuf_bindx_;
00378   Teuchos::ArrayRCP<LocalOrdinal> pbuf_indx_;
00379 
00380   LocalMatOps lclMatVec_;
00381   Teuchos::RCP<Tpetra::Import<LocalOrdinal,GlobalOrdinal,Node> > importer_;
00382   Teuchos::RCP<Tpetra::Export<LocalOrdinal,GlobalOrdinal,Node> > exporter_;
00383   mutable Teuchos::RCP<Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > importedVec_;
00384   mutable Teuchos::RCP<Tpetra::MultiVector<Scalar,LocalOrdinal,GlobalOrdinal,Node> > exportedVec_;
00385 
00386   typedef typename std::map<GlobalOrdinal,Teuchos::ArrayRCP<Scalar> > MapGlobalArrayRCP;
00387   typedef typename std::map<LocalOrdinal,Teuchos::ArrayRCP<Scalar> > MapLocalArrayRCP;
00388 
00389   //We use 2 arrays (well, array-of-maps, array-of-array-of-arrays...) to
00390   //represent the variable-block-row matrix in un-packed '2D' form.
00391   //
00392   //Note that these arrays are assumed to be resident in CPU (host) memory.
00393   //It doesn't make sense to copy this kind of data back and forth to a separate
00394   //compute device (e.g., a GPU), since we don't support doing matrix-vector
00395   //products until after fillComplete is called, at which time contiguous
00396   //arrays are allocated on the device and matrix data is copied into them.
00397   Teuchos::RCP<Teuchos::Array<MapGlobalArrayRCP> > col_ind_2D_global_;
00398   Teuchos::RCP<Teuchos::Array<Teuchos::Array<Teuchos::ArrayRCP<Scalar> > > > values2D_;
00399 
00400   bool is_fill_completed_;
00401   bool is_storage_optimized_;
00402 };//class VbrMatrix
00403 
00404 }//namespace Tpetra
00405 
00406 //----------------------------------------------------------------------------
00407 // Description of arrays representing the VBR format:
00408 //
00409 // (For more description, see this URL (valid as of 5/26/2010):
00410 // http://docs.sun.com/source/817-0086-10/prog-sparse-support.html)
00411 // ...and of course more can be found using google...
00412 // The old Aztec manual was a great resource for this but I can't
00413 // find a copy of that these days...
00414 //
00415 // rptr: length num_block_rows + 1
00416 //       rptr[i]: the pt-row corresponding to the i-th block-row
00417 //
00418 // cptr: length num_distinct_block_cols + 1
00419 //       cptr[j]: the pt-col corresponding to the j-th block-col
00420 //
00421 // bptr: length num_block_rows + 1
00422 //       bptr[i]: location in bindx of the first nonzero block-entry
00423 //                of the i-th block-row
00424 //
00425 // bindx: length num-nonzero-block-entries
00426 //        bindx[j]: block-col-index of j-th block-entry
00427 //
00428 // indx: length num-nonzero-block-entries + 1
00429 //       indx[j] location in vals of the beginning of the j-th
00430 //       block-entry
00431 //
00432 // vals: length num-nonzero-scalar-entries
00433 //
00434 //
00435 // Some example look-ups:
00436 //
00437 // int nbr = num_block_rows;
00438 // int total_num_block_nonzeros = bptr[nbr];
00439 // int total_num_scalar_nonzeros = indx[num_block_nonzeros];
00440 // 
00441 // //get arrays for i-th block-row:
00442 // int* bindx_i = &bindx[bptr[i]];
00443 // double* vals_i = &val[indx[bptr[i]]];
00444 // int num_block_nonzeros_in_row_i = bptr[i+1]-bptr[i];
00445 // 
00446 //----------------------------------------------------------------------------
00447 
00448 #endif //TPETRA_VBRMATRIX_DECL_HPP
00449 
 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