Hermes2D  3.0
linear_solver.cpp
Go to the documentation of this file.
1 // This file is part of HermesCommon
2 //
3 // Copyright (c) 2009 hp-FEM group at the University of Nevada, Reno (UNR).
4 // Email: hpfem-group@unr.edu, home page: http://www.hpfem.org/.
5 //
6 // Hermes2D is free software; you can redistribute it and/or modify
7 // it under the terms of the GNU General Public License as published
8 // by the Free Software Foundation; either version 2 of the License,
9 // or (at your option) any later version.
10 //
11 // Hermes2D is distributed in the hope that it will be useful,
12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14 // GNU General Public License for more details.
15 //
16 // You should have received a copy of the GNU General Public License
17 // along with Hermes2D; if not, write to the Free Software
18 // Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
22 #include "solver/linear_solver.h"
23 #include "solvers/matrix_solver.h"
24 
25 using namespace Hermes::Algebra;
26 
27 namespace Hermes
28 {
29  namespace Hermes2D
30  {
31  template<typename Scalar>
32  LinearSolver<Scalar>::LinearSolver(bool force_use_direct_solver) : Solver<Scalar>(false), Hermes::Solvers::MatrixSolver<Scalar>(force_use_direct_solver)
33  {
34  this->dp = new DiscreteProblem<Scalar>(true);
35  this->own_dp = true;
36  }
37 
38  template<typename Scalar>
39  LinearSolver<Scalar>::LinearSolver(DiscreteProblem<Scalar>* dp, bool force_use_direct_solver) : Solver<Scalar>(dp), Hermes::Solvers::MatrixSolver<Scalar>(force_use_direct_solver)
40  {
41  }
42 
43  template<typename Scalar>
44  LinearSolver<Scalar>::LinearSolver(WeakFormSharedPtr<Scalar> wf, SpaceSharedPtr<Scalar> space, bool force_use_direct_solver) : Solver<Scalar>(false), Hermes::Solvers::MatrixSolver<Scalar>(force_use_direct_solver)
45  {
46  this->dp = new DiscreteProblem<Scalar>(wf, space, true);
47  this->own_dp = true;
48  }
49 
50  template<typename Scalar>
51  LinearSolver<Scalar>::LinearSolver(WeakFormSharedPtr<Scalar> wf, std::vector<SpaceSharedPtr<Scalar> > spaces, bool force_use_direct_solver) : Solver<Scalar>(false), Hermes::Solvers::MatrixSolver<Scalar>(force_use_direct_solver)
52  {
53  this->dp = new DiscreteProblem<Scalar>(wf, spaces, true);
54  this->own_dp = true;
55  }
56 
57  template<typename Scalar>
58  LinearSolver<Scalar>::~LinearSolver()
59  {
60  }
61 
62  template<typename Scalar>
64  {
65  return this->sln_vector;
66  }
67 
68  template<typename Scalar>
70  {
71  return Solver<Scalar>::isOkay();
72  }
73 
74  template<typename Scalar>
76  {
78  this->jacobian_reusable = false;
79  }
80 
81  template<typename Scalar>
83  {
85  this->jacobian_reusable = false;
86  }
87 
88  template<typename Scalar>
90  {
91  Hermes::Solvers::MatrixSolver<Scalar>::set_verbose_output(to_set);
92  this->dp->set_verbose_output(to_set);
93  }
94 
95  template<typename Scalar>
96  void LinearSolver<Scalar>::solve(Scalar* coeff_vec)
97  {
98  this->check();
99 
100  this->on_initialization();
101 
102  this->tick();
103 
104  // Extremely important.
105  Space<Scalar>::assign_dofs(this->dp->get_spaces());
106 
107  // Assemble the residual always and the Matrix when necessary (nonconstant jacobian, not reusable, ...).
108  if (this->jacobian_reusable && this->constant_jacobian)
109  {
110  this->info("\tLinearSolver: assembling... [reusing matrix, assembling rhs].");
111  this->dp->assemble(coeff_vec, this->get_residual());
112  this->linear_matrix_solver->set_reuse_scheme(Hermes::Solvers::HERMES_REUSE_MATRIX_STRUCTURE_COMPLETELY);
113  }
114  else
115  {
116  if (this->jacobian_reusable)
117  this->info("\tLinearSolver: assembling... [re-assembling with a reusable matrix structure].");
118  else
119  this->info("\tLinearSolver: assembling... [assembling the matrix and rhs anew].");
120  this->dp->assemble(coeff_vec, this->get_jacobian(), this->get_residual());
121  this->linear_matrix_solver->set_reuse_scheme(Hermes::Solvers::HERMES_CREATE_STRUCTURE_FROM_SCRATCH);
122  }
123 
124  this->process_matrix_output(this->get_jacobian(), 1);
125  this->process_vector_output(this->get_residual(), 1);
126 
127  this->tick();
128  this->info("\tLinearSolver: assembling done in %s. Solving...", this->last_str().c_str());
129  this->tick();
130 
131  // Solve, if the solver is iterative, give him the initial guess.
132  this->linear_matrix_solver->solve(coeff_vec);
133 
134  this->sln_vector = this->linear_matrix_solver->get_sln_vector();
135 
136  this->on_finish();
137 
138  this->tick();
139  this->info("\tLinearSolver: solving done in %s.", this->last_str().c_str());
140  }
141 
142  template class HERMES_API LinearSolver < double > ;
143  template class HERMES_API LinearSolver < std::complex<double> > ;
144  }
145 }
Definition: adapt.h:24
Used to pass the instances of Space around.
Definition: space.h:34
Used to pass the instances of WeakForm around.
Definition: weakform.h:55
Represents a finite element space over a domain.
Definition: api2d.h:34