Tpetra Matrix/Vector Services Version of the Day
TpetraExt_MMHelpers_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_MMHELPERS_DEF_HPP
00030 #define TPETRA_MMHELPERS_DEF_HPP
00031 
00032 #include "Teuchos_VerboseObject.hpp"
00033 #include "Tpetra_ConfigDefs.hpp"
00034 #include "Tpetra_CrsMatrix.hpp"
00035 #include "TpetraExt_MMHelpers_decl.hpp"
00036 
00041 namespace Tpetra {
00042 
00043 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class SpMatOps>
00044 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps>::CrsMatrixStruct()
00045 {
00046 }
00047 
00048 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class SpMatOps>
00049 CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps>::~CrsMatrixStruct()
00050 {
00051   deleteContents();
00052 }
00053 
00054 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class SpMatOps>
00055 void CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps>::deleteContents()
00056 {
00057   numRows = 0;
00058   numEntriesPerRow.clear();
00059   indices.clear();
00060   values.clear();
00061   remote.clear();
00062   numRemote = 0;
00063   importMatrix.reset();
00064 }
00065 
00066 template <class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class SpMatOps>
00067 int dumpCrsMatrixStruct(const CrsMatrixStruct<Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps>& M)
00068 {
00069   std::cout << "proc " << M.rowMap->Comm().MyPID()<<std::endl;
00070   std::cout << "numRows: " << M.numRows<<std::endl;
00071   for(LocalOrdinal i=0; i<M.numRows; ++i) {
00072     for(LocalOrdinal j=0; j<M.numEntriesPerRow[i]; ++j) {
00073       if (M.remote[i]) {
00074         std::cout << "  *"<<M.rowMap->GID(i)<<"   "
00075              <<M.importColMap->GID(M.indices[i][j])<<"   "<<M.values[i][j]<<std::endl;
00076       }
00077       else {
00078         std::cout << "   "<<M.rowMap->GID(i)<<"   "
00079              <<M.colMap->GID(M.indices[i][j])<<"   "<<M.values[i][j]<<std::endl;
00080       }
00081     }
00082   }
00083   return(0);
00084 }
00085 
00086 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class SpMatOps>
00087 CrsWrapper_CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps>::CrsWrapper_CrsMatrix(CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps>& crsmatrix)
00088  : crsmat_(crsmatrix)
00089 {
00090 }
00091 
00092 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class SpMatOps>
00093 CrsWrapper_CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps>::~CrsWrapper_CrsMatrix()
00094 {
00095 }
00096 
00097 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class SpMatOps>
00098 Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> >
00099 CrsWrapper_CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps>::getRowMap() const
00100 {
00101   return crsmat_.getRowMap();
00102 }
00103 
00104 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class SpMatOps>
00105 bool CrsWrapper_CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps>::isFillComplete()
00106 {
00107   return crsmat_.isFillComplete();
00108 }
00109 
00110 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class SpMatOps>
00111 void
00112 CrsWrapper_CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps>::insertGlobalValues(GlobalOrdinal globalRow, const Teuchos::ArrayView<const GlobalOrdinal> &indices, const Teuchos::ArrayView<const Scalar> &values)
00113 {
00114   crsmat_.insertGlobalValues(globalRow, indices, values);
00115 }
00116 
00117 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class SpMatOps>
00118 void 
00119 CrsWrapper_CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps>::sumIntoGlobalValues(GlobalOrdinal globalRow, const Teuchos::ArrayView<const GlobalOrdinal> &indices, const Teuchos::ArrayView<const Scalar> &values)
00120 {
00121   crsmat_.sumIntoGlobalValues(globalRow, indices, values);
00122 }
00123 
00124 
00125 //------------------------------------
00126 
00127 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
00128 CrsWrapper_GraphBuilder<Scalar, LocalOrdinal, GlobalOrdinal, Node>::CrsWrapper_GraphBuilder(const Teuchos::RCP<const Map<LocalOrdinal, GlobalOrdinal, Node> >& map)
00129  : graph_(),
00130    rowmap_(map),
00131    max_row_length_(0)
00132 {
00133   Teuchos::ArrayView<const GlobalOrdinal> rows= map->getNodeElementList();
00134 
00135   for(int i=0; i<rows.size(); ++i) {
00136     graph_[rows[i]] = new std::set<GlobalOrdinal>;
00137   }
00138 }
00139 
00140 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
00141 CrsWrapper_GraphBuilder<Scalar, LocalOrdinal, GlobalOrdinal, Node>::~CrsWrapper_GraphBuilder()
00142 {
00143   typename std::map<GlobalOrdinal,std::set<GlobalOrdinal>*>::iterator
00144     iter = graph_.begin(), iter_end = graph_.end();
00145   for(; iter!=iter_end; ++iter) {
00146     delete iter->second;
00147   }
00148 
00149   graph_.clear();
00150 }
00151 
00152 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
00153 bool CrsWrapper_GraphBuilder<Scalar, LocalOrdinal, GlobalOrdinal, Node>::isFillComplete()
00154 {
00155   return false;
00156 }
00157 
00158 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
00159 void
00160 CrsWrapper_GraphBuilder<Scalar, LocalOrdinal, GlobalOrdinal, Node>::insertGlobalValues(GlobalOrdinal globalRow, const Teuchos::ArrayView<const GlobalOrdinal> &indices, const Teuchos::ArrayView<const Scalar> &values)
00161 {
00162   /*std::cout << "inserting" << std::endl;
00163   std::cout << "  indices: " << indices << std::endl;
00164   std::cout << "  values: " << values << std::endl;*/
00165   typename std::map<GlobalOrdinal,std::set<GlobalOrdinal>*>::iterator
00166     iter = graph_.find(globalRow);
00167 
00168   TEST_FOR_EXCEPTION(iter == graph_.end(), std::runtime_error,
00169   Teuchos::typeName(*this) << "::insertGlobalValues could not find row " << globalRow << " in the graph. Super bummer man. Hope you figure it out.\n");
00170 
00171   std::set<GlobalOrdinal>& cols = *(iter->second);
00172 
00173   for(int i=0; i<indices.size(); ++i) {
00174     cols.insert(indices[i]);
00175   }
00176 
00177   global_size_t row_length = cols.size();
00178   if (row_length > max_row_length_) max_row_length_ = row_length;
00179 
00180   
00181 }
00182 
00183 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
00184 void
00185 CrsWrapper_GraphBuilder<Scalar, LocalOrdinal, GlobalOrdinal, Node>::sumIntoGlobalValues(GlobalOrdinal globalRow, const Teuchos::ArrayView<const GlobalOrdinal> &indices, const Teuchos::ArrayView<const Scalar> &values)
00186 {
00187   insertGlobalValues(globalRow, indices, values);
00188 }
00189 
00190 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node>
00191 std::map<GlobalOrdinal,std::set<GlobalOrdinal>*>&
00192 CrsWrapper_GraphBuilder<Scalar, LocalOrdinal, GlobalOrdinal, Node>::get_graph()
00193 {
00194   return graph_;
00195 }
00196 
00197 template<class Scalar, class LocalOrdinal, class GlobalOrdinal, class Node, class SpMatOps>
00198 void insert_matrix_locations(CrsWrapper_GraphBuilder<Scalar, LocalOrdinal, GlobalOrdinal, Node> & graphbuilder,
00199   CrsMatrix<Scalar, LocalOrdinal, GlobalOrdinal, Node, SpMatOps>& C)
00200 {
00201   global_size_t max_row_length = graphbuilder.get_max_row_length();
00202   if (max_row_length < 1) return;
00203 
00204   Array<GlobalOrdinal> indices(max_row_length);
00205   Array<Scalar> zeros(max_row_length, 0.0);
00206 
00207   typedef std::map<GlobalOrdinal,std::set<GlobalOrdinal>*> Graph;
00208   typedef typename Graph::iterator    GraphIter;
00209   typedef std::set<GlobalOrdinal>     Set;
00210   typedef typename Set::iterator      SetIter;
00211   Graph& graph = graphbuilder.get_graph();
00212 
00213   const GraphIter iter_end = graph.end();
00214   for(GraphIter iter=graph.begin(); iter!=iter_end; ++iter) {
00215     const GlobalOrdinal row = iter->first;
00216     const std::set<GlobalOrdinal> &cols = *(iter->second);
00217     // "copy" entries out of set into contiguous array storage
00218     const size_t num_entries = std::copy(cols.begin(), cols.end(), indices.begin()) - indices.begin();
00219     // insert zeros into the result matrix at the appropriate locations
00220     C.insertGlobalValues(row, indices(0,num_entries), zeros(0,num_entries));
00221   }
00222 }
00223 
00224 } //namespace Tpetra
00225 
00226 //
00227 // Explicit instantiation macro
00228 //
00229 // Must be expanded from within the Tpetra namespace!
00230 //
00231 
00232 #define TPETRA_CRSMATRIXSTRUCT_INSTANT(SCALAR,LO,GO,NODE) \
00233   \
00234   template class CrsMatrixStruct< SCALAR , LO , GO , NODE >;
00235 
00236 #define TPETRA_CRSWRAPPER_INSTANT(SCALAR,LO,GO,NODE) \
00237   \
00238   template class CrsWrapper< SCALAR , LO , GO , NODE >;
00239 
00240 #define TPETRA_CRSWRAPPER_CRSMATRIX_INSTANT(SCALAR,LO,GO,NODE) \
00241   \
00242   template class CrsWrapper_CrsMatrix< SCALAR , LO , GO , NODE >;
00243 
00244 #define TPETRA_CRSWRAPPER_GRAPHBUILDER_INSTANT(SCALAR,LO,GO,NODE) \
00245   \
00246   template class CrsWrapper_GraphBuilder< SCALAR , LO , GO , NODE >;
00247 
00248 #endif // TPETRA_MMHELPERS_DEF_HPP
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Defines