Sierra Toolkit  Version of the Day
LinearSystem.cpp
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 
10 #include <stk_linsys/LinearSystem.hpp>
11 #include <stk_linsys/ImplDetails.hpp>
12 #include <stk_mesh/base/GetBuckets.hpp>
13 
14 #include <stk_linsys/LinsysFunctions.hpp>
15 
16 #include <fei_Trilinos_Helpers.hpp>
17 
18 namespace stk_classic {
19 namespace linsys {
20 
21 LinearSystem::LinearSystem(MPI_Comm comm, fei::SharedPtr<fei::Factory> factory)
22  : m_fei_factory(factory),
23  m_dof_mapper(comm),
24  m_fei_mgraph(new fei::MatrixGraph_Impl2(m_dof_mapper.get_fei_VectorSpace(), fei::SharedPtr<fei::VectorSpace>())),
25  m_fei_linearsystem(),
26  m_param_set()
27 {
28 }
29 
31 {
32 }
33 
34 void
35 LinearSystem::set_parameters(Teuchos::ParameterList& paramlist)
36 {
37  Trilinos_Helpers::copy_parameterlist(paramlist, m_param_set);
38  m_dof_mapper.get_fei_VectorSpace()->setParameters(m_param_set);
39  if (m_fei_factory.get() != NULL) {
40  m_fei_factory->parameters(m_param_set);
41  }
42  if (m_fei_mgraph.get() != NULL) {
43  m_fei_mgraph->setParameters(m_param_set);
44  }
45 }
46 
47 void
49 {
50  m_fei_mgraph->initComplete();
51  m_dof_mapper.finalize();
52 }
53 
54 void
56 {
57  m_fei_linearsystem = m_fei_factory->createLinearSystem(m_fei_mgraph);
58  m_fei_linearsystem->parameters(m_param_set);
59 
60  fei::SharedPtr<fei::Matrix> mtx = m_fei_factory->createMatrix(m_fei_mgraph);
61  mtx->parameters(m_param_set);
62  m_fei_linearsystem->setMatrix(mtx);
63  fei::SharedPtr<fei::Vector> rhs = m_fei_factory->createVector(m_fei_mgraph);
64  m_fei_linearsystem->setRHS(rhs);
65  fei::SharedPtr<fei::Vector> soln = m_fei_factory->createVector(m_fei_mgraph,true);
66  m_fei_linearsystem->setSolutionVector(soln);
67 }
68 
69 void
71 {
72  m_fei_linearsystem->loadComplete();
73 }
74 
75 const DofMapper&
77 {
78  return m_dof_mapper;
79 }
80 
81 DofMapper&
83 {
84  return m_dof_mapper;
85 }
86 
87 void
88 LinearSystem::reset_to_zero()
89 {
90  fei::SharedPtr<fei::Matrix> mtx = m_fei_linearsystem->getMatrix();
91  if (mtx.get() != NULL) {
92  mtx->putScalar(0);
93  }
94 
95  fei::SharedPtr<fei::Vector> rhs = m_fei_linearsystem->getRHS();
96  if (rhs.get() != NULL) {
97  rhs->putScalar(0);
98  }
99 }
100 
101 const fei::SharedPtr<fei::MatrixGraph>
103 {
104  return m_fei_mgraph;
105 }
106 
107 fei::SharedPtr<fei::MatrixGraph>
109 {
110  return m_fei_mgraph;
111 }
112 
113 const fei::SharedPtr<fei::LinearSystem>
115 {
116  return m_fei_linearsystem;
117 }
118 
119 fei::SharedPtr<fei::LinearSystem>
121 {
122  return m_fei_linearsystem;
123 }
124 
125 void
126 LinearSystem::write_files(const std::string& base_name) const
127 {
128  fei::SharedPtr<fei::Matrix> A = m_fei_linearsystem->getMatrix();
129  if (A.get() != NULL) {
130  std::ostringstream ossA;
131  ossA << "A_" << base_name << ".mtx";
132  std::string Aname = ossA.str();
133  A->writeToFile(Aname.c_str());
134  }
135  fei::SharedPtr<fei::Vector> b = m_fei_linearsystem->getRHS();
136  if (b.get() != NULL) {
137  std::ostringstream ossb;
138  ossb << "b_" << base_name << ".vec";
139  std::string bname = ossb.str();
140  b->writeToFile(bname.c_str());
141  }
142 }
143 
144 int
145 LinearSystem::solve(int &status, const Teuchos::ParameterList & params )
146 {
147  return fei_solve(status, *m_fei_linearsystem, params);
148 }
149 
150 }//namespace linsys
151 }//namespace stk_classic
152 
int fei_solve(int &status, fei::LinearSystem &fei_ls, const Teuchos::ParameterList &params)
const fei::SharedPtr< fei::MatrixGraph > get_fei_MatrixGraph() const
const fei::SharedPtr< fei::VectorSpace > get_fei_VectorSpace() const
Definition: DofMapper.hpp:121
const fei::SharedPtr< fei::LinearSystem > get_fei_LinearSystem() const
const DofMapper & get_DofMapper() const
int solve(int &status, const Teuchos::ParameterList &params)
Sierra Toolkit.
LinearSystem(MPI_Comm comm, fei::SharedPtr< fei::Factory > factory)