libMesh
unsteady_solver.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 #include "libmesh/diff_solver.h"
20 #include "libmesh/diff_system.h"
21 #include "libmesh/dof_map.h"
22 #include "libmesh/numeric_vector.h"
23 #include "libmesh/unsteady_solver.h"
24 
25 namespace libMesh
26 {
27 
28 
29 
31  : TimeSolver(s),
32  old_local_nonlinear_solution (NumericVector<Number>::build(s.comm())),
33  first_solve (true),
34  first_adjoint_step (true)
35 {
36 }
37 
38 
39 
41 {
42 }
43 
44 
45 
47 {
49 
50  _system.add_vector("_old_nonlinear_solution");
51 }
52 
53 
54 
56 {
58 
59 #ifdef LIBMESH_ENABLE_GHOSTED
62  GHOSTED);
63 #else
65 #endif
66 }
67 
68 
69 
71 {
73 
74 #ifdef LIBMESH_ENABLE_GHOSTED
77  GHOSTED);
78 #else
80 #endif
81 
82  // localize the old solution
83  NumericVector<Number> & old_nonlinear_soln =
84  _system.get_vector("_old_nonlinear_solution");
85 
86  old_nonlinear_soln.localize
89 }
90 
91 
92 
94 {
95  if (first_solve)
96  {
98  first_solve = false;
99  }
100 
101  unsigned int solve_result = _diff_solver->solve();
102 
103  // If we requested the UnsteadySolver to attempt reducing dt after a
104  // failed DiffSolver solve, check the results of the solve now.
106  {
107  bool backtracking_failed =
109 
110  bool max_iterations =
112 
113  if (backtracking_failed || max_iterations)
114  {
115  // Cut timestep in half
116  for (unsigned int nr=0; nr<reduce_deltat_on_diffsolver_failure; ++nr)
117  {
118  _system.deltat *= 0.5;
119  libMesh::out << "Newton backtracking failed. Trying with smaller timestep, dt="
120  << _system.deltat << std::endl;
121 
122  solve_result = _diff_solver->solve();
123 
124  // Check solve results with reduced timestep
125  bool backtracking_still_failed =
127 
128  bool backtracking_max_iterations =
130 
131  if (!backtracking_still_failed && !backtracking_max_iterations)
132  {
133  if (!quiet)
134  libMesh::out << "Reduced dt solve succeeded." << std::endl;
135  return;
136  }
137  }
138 
139  // If we made it here, we still couldn't converge the solve after
140  // reducing deltat
141  libMesh::out << "DiffSolver::solve() did not succeed after "
142  << reduce_deltat_on_diffsolver_failure
143  << " attempts." << std::endl;
144  libmesh_convergence_failure();
145 
146  } // end if (backtracking_failed || max_iterations)
147  } // end if (reduce_deltat_on_diffsolver_failure)
148 }
149 
150 
151 
153 {
154  if (!first_solve)
155  {
156  // Store the solution, does nothing by default
157  // User has to attach appropriate solution_history object for this to
158  // actually store anything anywhere
159  solution_history->store();
160 
162  }
163 
164  NumericVector<Number> & old_nonlinear_soln =
165  _system.get_vector("_old_nonlinear_solution");
166  NumericVector<Number> & nonlinear_solution =
167  *(_system.solution);
168 
169  old_nonlinear_soln = nonlinear_solution;
170 
171  old_nonlinear_soln.localize
174 }
175 
176 
177 
179 {
180  // On the first call of this function, we dont save the adjoint solution or
181  // decrement the time, we just call the retrieve function below
182  if (!first_adjoint_step)
183  {
184  // Call the store function to store the last adjoint before decrementing the time
185  solution_history->store();
186  // Decrement the system time
188  }
189  else
190  {
191  first_adjoint_step = false;
192  }
193 
194  // Retrieve the primal solution vectors at this time using the
195  // solution_history object
196  solution_history->retrieve();
197 
198  // Dont forget to localize the old_nonlinear_solution !
199  _system.get_vector("_old_nonlinear_solution").localize
202 }
203 
205 {
206  // Retrieve all the stored vectors at the current time
207  solution_history->retrieve();
208 
209  // Dont forget to localize the old_nonlinear_solution !
210  _system.get_vector("_old_nonlinear_solution").localize
213 }
214 
215 
217  const
218 {
219  libmesh_assert_less (global_dof_number, _system.get_dof_map().n_dofs());
220  libmesh_assert_less (global_dof_number, old_local_nonlinear_solution->size());
221 
222  return (*old_local_nonlinear_solution)(global_dof_number);
223 }
224 
225 
226 
227 Real UnsteadySolver::du(const SystemNorm & norm) const
228 {
229 
230  UniquePtr<NumericVector<Number>> solution_copy =
231  _system.solution->clone();
232 
233  solution_copy->add(-1., _system.get_vector("_old_nonlinear_solution"));
234 
235  solution_copy->close();
236 
237  return _system.calculate_norm(*solution_copy, norm);
238 }
239 
240 } // namespace libMesh
bool quiet
Print extra debugging information if quiet == false.
Definition: time_solver.h:191
Real time
For time-dependent problems, this is the time t at the beginning of the current timestep.
Definition: system.h:1545
UnsteadySolver(sys_type &s)
Constructor.
This is a generic class that defines a solver to handle time integration of DifferentiableSystems.
Definition: time_solver.h:58
virtual void init_data()
The data initialization function.
Definition: time_solver.C:77
virtual void reinit()
The reinitialization function.
Definition: time_solver.C:48
Numeric vector.
Definition: dof_map.h:66
This class defines a norm/seminorm to be applied to a NumericVector which contains coefficients in a ...
Definition: system_norm.h:43
The libMesh namespace provides an interface to certain functionality in the library.
virtual void advance_timestep() libmesh_override
This method advances the solution to the next timestep, after a solve() has been performed.
virtual void retrieve_timestep() libmesh_override
This method retrieves all the stored solutions at the current system.time.
dof_id_type n_dofs() const
Definition: dof_map.h:510
virtual void init() libmesh_override
The initialization function.
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
This class provides a specific system class.
Definition: diff_system.h:53
The DiffSolver reached the maximum allowed number of nonlinear iterations before satisfying any conve...
Definition: diff_solver.h:270
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:257
Number old_nonlinear_solution(const dof_id_type global_dof_number) const
NumericVector< Number > & add_vector(const std::string &vec_name, const bool projections=true, const ParallelType type=PARALLEL)
Adds the additional vector vec_name to this system.
Definition: system.C:681
const DofMap & get_dof_map() const
Definition: system.h:2030
virtual ~UnsteadySolver()
Destructor.
virtual void init()
The initialization function.
Definition: time_solver.C:64
Real deltat
For time-dependent problems, this is the amount delta t to advance the solution in time...
Definition: diff_system.h:248
UniquePtr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition: system.h:1523
virtual void adjoint_advance_timestep() libmesh_override
This method advances the adjoint solution to the previous timestep, after an adjoint_solve() has been...
unsigned int reduce_deltat_on_diffsolver_failure
This value (which defaults to zero) is the number of times the TimeSolver is allowed to halve deltat ...
Definition: time_solver.h:220
const NumericVector< Number > & get_vector(const std::string &vec_name) const
Definition: system.C:794
virtual void reinit() libmesh_override
The reinitialization function.
Real calculate_norm(const NumericVector< Number > &v, unsigned int var, FEMNormType norm_type, std::set< unsigned int > *skip_dimensions=libmesh_nullptr) const
Definition: system.C:1383
dof_id_type n_local_dofs() const
Definition: system.C:185
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual Real du(const SystemNorm &norm) const libmesh_override
Computes the size of ||u^{n+1} - u^{n}|| in some norm.
OStreamProxy out
dof_id_type n_dofs() const
Definition: system.C:148
The DiffSolver failed to find a descent direction by backtracking (See newton_solver.C)
Definition: diff_solver.h:276
UniquePtr< DiffSolver > _diff_solver
An implicit linear or nonlinear solver to use at each timestep.
Definition: time_solver.h:247
UniquePtr< SolutionHistory > solution_history
A UniquePtr to a SolutionHistory object.
Definition: time_solver.h:264
const std::vector< dof_id_type > & get_send_list() const
Definition: dof_map.h:394
virtual void localize(std::vector< T > &v_local) const =0
Creates a copy of the global vector in the local vector v_local.
virtual void init_data() libmesh_override
The data initialization function.
virtual void solve() libmesh_override
This method solves for the solution at the next timestep.
bool first_adjoint_step
A bool that will be true the first time adjoint_advance_timestep() is called, (when the primal soluti...
uint8_t dof_id_type
Definition: id_types.h:64
bool first_solve
A bool that will be true the first time solve() is called, and false thereafter.
UniquePtr< NumericVector< Number > > old_local_nonlinear_solution
Serial vector of _system.get_vector("_old_nonlinear_solution")