Sierra Toolkit  Version of the Day
AggregateLinearSystem.hpp
1 /*------------------------------------------------------------------------*/
2 /* Copyright 2010 Sandia Corporation. */
3 /* Under terms of Contract DE-AC04-94AL85000, there is a non-exclusive */
4 /* license for use of this work by or on behalf of the U.S. Government. */
5 /* Export of this program may require a license from the */
6 /* United States Government. */
7 /*------------------------------------------------------------------------*/
8 
9 #ifndef stk_linsys_AggregateLinearSystem_hpp
10 #define stk_linsys_AggregateLinearSystem_hpp
11 
12 #include <stk_linsys/FeiBaseIncludes.hpp>
13 #include <stk_linsys/DofMapper.hpp>
14 #include <stk_linsys/LinearSystemInterface.hpp>
15 #include <stk_linsys/LinearSystem.hpp>
16 
17 #include <Teuchos_ParameterList.hpp>
18 
19 namespace stk_classic {
20 namespace linsys {
21 
32 class AggregateLinearSystem : public LinearSystemInterface {
33  public:
35  AggregateLinearSystem(MPI_Comm comm, fei::SharedPtr<fei::Factory> factory, size_t num_matrices=1, size_t num_rhsvecs=1);
36 
38  virtual ~AggregateLinearSystem();
39 
40  void set_parameters(Teuchos::ParameterList& paramlist);
41 
43  void set_num_matrices_rhsvecs(size_t num_matrices, size_t num_rhsvecs);
44 
51 
57 
60  fei::SharedPtr<fei::Matrix> get_matrix(size_t index);
61 
64  fei::SharedPtr<fei::Vector> get_rhsvec(size_t index);
65 
70  void aggregate_system(const std::vector<double>& mat_scalars,
71  const std::vector<double>& rhs_scalars);
72 
80  void finalize_assembly();
81 
83  const DofMapper& get_DofMapper() const;
84 
87 
88  void reset_to_zero();
89 
91  const fei::SharedPtr<fei::MatrixGraph> get_fei_MatrixGraph() const;
92 
94  fei::SharedPtr<fei::MatrixGraph> get_fei_MatrixGraph();
95 
97  const fei::SharedPtr<fei::LinearSystem> get_fei_LinearSystem() const;
98 
100  fei::SharedPtr<fei::LinearSystem> get_fei_LinearSystem();
101 
102  void write_files(const std::string& base_name) const;
103 
122  int solve(int & status, const Teuchos::ParameterList & params);
123 
124  private:
125 
126  fei::SharedPtr<fei::Factory> m_fei_factory;
127  stk_classic::linsys::LinearSystem m_linear_system;
128 
129  std::vector<fei::SharedPtr<fei::Matrix> > m_matrices;
130  std::vector<fei::SharedPtr<fei::Vector> > m_rhsvecs;
131 };//class AggregateLinearSystem
132 
133 }//namespace linsys
134 }//namespace stk_classic
135 
136 #endif
137 
fei::SharedPtr< fei::Vector > get_rhsvec(size_t index)
int solve(int &status, const Teuchos::ParameterList &params)
void aggregate_system(const std::vector< double > &mat_scalars, const std::vector< double > &rhs_scalars)
fei::SharedPtr< fei::Matrix > get_matrix(size_t index)
const fei::SharedPtr< fei::LinearSystem > get_fei_LinearSystem() const
AggregateLinearSystem(MPI_Comm comm, fei::SharedPtr< fei::Factory > factory, size_t num_matrices=1, size_t num_rhsvecs=1)
const fei::SharedPtr< fei::MatrixGraph > get_fei_MatrixGraph() const
Sierra Toolkit.
void set_num_matrices_rhsvecs(size_t num_matrices, size_t num_rhsvecs)