44 #ifndef STOKHOS_SPARSE_3_TENSOR_UTILITIES_HPP 45 #define STOKHOS_SPARSE_3_TENSOR_UTILITIES_HPP 48 #include "Epetra_Comm.h" 49 #include "Epetra_CrsGraph.h" 50 #include "Epetra_BlockMap.h" 51 #include "Epetra_LocalMap.h" 52 #include "Epetra_CrsMatrix.h" 53 #include "EpetraExt_RowMatrixOut.h" 65 template <
typename ordinal_type,
typename value_type>
66 Teuchos::RCP<Epetra_CrsGraph>
70 const Epetra_Comm& comm)
78 Epetra_LocalMap map(num_rows, 0, comm);
81 Teuchos::RCP<Epetra_CrsGraph> graph =
82 Teuchos::rcp(
new Epetra_CrsGraph(Copy, map, 0));
86 for (
typename Cijk_type::k_iterator k_it=Cijk.
k_begin();
87 k_it!=Cijk.
k_end(); ++k_it) {
88 for (
typename Cijk_type::kj_iterator j_it = Cijk.
j_begin(k_it);
89 j_it != Cijk.
j_end(k_it); ++j_it) {
91 for (
typename Cijk_type::kji_iterator i_it = Cijk.
i_begin(j_it);
92 i_it != Cijk.
i_end(j_it); ++i_it) {
94 graph->InsertGlobalIndices(i, 1, &
j);
100 graph->FillComplete();
113 template <
typename ordinal_type,
typename value_type>
114 Teuchos::RCP<Epetra_CrsGraph>
117 const Epetra_BlockMap& map)
122 Teuchos::RCP<Epetra_CrsGraph> graph =
123 Teuchos::rcp(
new Epetra_CrsGraph(Copy, map, 0));
127 for (
typename Cijk_type::k_iterator k_it=Cijk.
k_begin();
128 k_it!=Cijk.
k_end(); ++k_it) {
129 for (
typename Cijk_type::kj_iterator j_it = Cijk.
j_begin(k_it);
130 j_it != Cijk.
j_end(k_it); ++j_it) {
132 for (
typename Cijk_type::kji_iterator i_it = Cijk.
i_begin(j_it);
133 i_it != Cijk.
i_end(j_it); ++i_it) {
135 graph->InsertGlobalIndices(i, 1, &
j);
141 graph->FillComplete();
146 template <
typename ordinal_type,
typename value_type>
151 const Epetra_Comm& comm,
152 const std::string& file)
154 Teuchos::RCP<Epetra_CrsGraph> graph =
156 Epetra_CrsMatrix mat(Copy, *graph);
159 EpetraExt::RowMatrixToMatrixMarketFile(file.c_str(), mat);
162 template <
typename ordinal_type,
typename value_type>
166 const Epetra_BlockMap& map,
167 const std::string& file)
169 Teuchos::RCP<Epetra_CrsGraph> graph =
171 Epetra_CrsMatrix mat(Copy, *graph);
174 EpetraExt::RowMatrixToMatrixMarketFile(file.c_str(), mat);
179 #endif // SPARSE_3_TENSOR_UTILITIES_HPP
kji_iterator i_end(const kj_iterator &j) const
Iterator pointing to last i entry for given j and k.
Data structure storing a sparse 3-tensor C(i,j,k) in a a compressed format.
void sparse3Tensor2MatrixMarket(const Stokhos::OrthogPolyBasis< ordinal_type, value_type > &basis, const Stokhos::Sparse3Tensor< ordinal_type, value_type > &Cijk, const Epetra_Comm &comm, const std::string &file)
Abstract base class for multivariate orthogonal polynomials.
Top-level namespace for Stokhos classes and functions.
Stokhos::Sparse3Tensor< int, double > Cijk_type
Teuchos::RCP< Epetra_CrsGraph > sparse3Tensor2CrsGraph(const Stokhos::OrthogPolyBasis< ordinal_type, value_type > &basis, const Stokhos::Sparse3Tensor< ordinal_type, value_type > &Cijk, const Epetra_Comm &comm)
Build an Epetra_CrsGraph from a sparse 3 tensor.
k_iterator k_end() const
Iterator pointing to last k entry.
kj_iterator j_begin(const k_iterator &k) const
Iterator pointing to first j entry for given k.
kj_iterator j_end(const k_iterator &k) const
Iterator pointing to last j entry for given k.
kji_iterator i_begin(const kj_iterator &j) const
Iterator pointing to first i entry for given j and k.
virtual ordinal_type size() const =0
Return total size of basis.
k_iterator k_begin() const
Iterator pointing to first k entry.