libMesh
adaptive_time_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 #include "libmesh/adaptive_time_solver.h"
19 #include "libmesh/diff_system.h"
20 #include "libmesh/numeric_vector.h"
21 
22 namespace libMesh
23 {
24 
25 
26 
29  core_time_solver(),
30  target_tolerance(1.e-3),
31  upper_tolerance(0.0),
32  max_deltat(0.),
33  min_deltat(0.),
34  max_growth(0.),
35  global_tolerance(true)
36 {
37  // the child class must populate core_time_solver
38  // with whatever actual time solver is to be used
39 
40  // As an UnsteadySolver, we have an old_local_nonlinear_solution, but we're
41  // going to drop it and use our core time solver's instead.
43 }
44 
45 
46 
48 {
49  // As an UnsteadySolver, we have an old_local_nonlinear_solution, but it
50  // is being managed by our core_time_solver. Make sure we don't delete
51  // it out from under them, in case the user wants to keep using the core
52  // solver after they're done with us.
54 }
55 
56 
57 
59 {
61 
62  // We override this because our core_time_solver is the one that
63  // needs to handle new vectors, diff_solver->init(), etc
64  core_time_solver->init();
65 
66  // As an UnsteadySolver, we have an old_local_nonlinear_solution, but it
67  // isn't pointing to the right place - fix it
68  //
69  // This leaves us with two UniquePtrs holding the same pointer - dangerous
70  // for future use. Replace with shared_ptr?
72  UniquePtr<NumericVector<Number>>(core_time_solver->old_local_nonlinear_solution.get());
73 }
74 
75 
76 
78 {
80 
81  // We override this because our core_time_solver is the one that
82  // needs to handle new vectors, diff_solver->reinit(), etc
83  core_time_solver->reinit();
84 }
85 
86 
88 {
89  NumericVector<Number> & old_nonlinear_soln =
90  _system.get_vector("_old_nonlinear_solution");
91  NumericVector<Number> & nonlinear_solution =
92  *(_system.solution);
93  // _system.get_vector("_nonlinear_solution");
94 
95  old_nonlinear_soln = nonlinear_solution;
96 
97  if (!first_solve)
99 }
100 
101 
102 
104 {
106 
107  return core_time_solver->error_order();
108 }
109 
110 
111 
112 bool AdaptiveTimeSolver::element_residual (bool request_jacobian,
113  DiffContext & context)
114 {
116 
117  return core_time_solver->element_residual(request_jacobian, context);
118 }
119 
120 
121 
122 bool AdaptiveTimeSolver::side_residual (bool request_jacobian,
123  DiffContext & context)
124 {
126 
127  return core_time_solver->side_residual(request_jacobian, context);
128 }
129 
130 
131 
132 bool AdaptiveTimeSolver::nonlocal_residual (bool request_jacobian,
133  DiffContext & context)
134 {
136 
137  return core_time_solver->nonlocal_residual(request_jacobian, context);
138 }
139 
140 
141 
143 {
144  return core_time_solver->diff_solver();
145 }
146 
147 
148 
150 {
151  return core_time_solver->linear_solver();
152 }
153 
154 
155 
158 {
159  return s.calculate_norm(v, component_norm);
160 }
161 
162 } // namespace libMesh
virtual Real calculate_norm(System &, NumericVector< Number > &)
A helper function to calculate error norms.
Real time
For time-dependent problems, this is the time t at the beginning of the current timestep.
Definition: system.h:1545
virtual bool side_residual(bool get_jacobian, DiffContext &) libmesh_override
This method is passed on to the core_time_solver.
This class provides all data required for a physics package (e.g.
Definition: diff_context.h:54
SystemNorm component_norm
Error calculations are done in this norm, DISCRETE_L2 by default.
Real last_deltat
We need to store the value of the last deltat used, so that advance_timestep() will increment the sys...
The libMesh namespace provides an interface to certain functionality in the library.
virtual UniquePtr< LinearSolver< Number > > & linear_solver() libmesh_override
An implicit linear solver to use for adjoint and sensitivity problems.
virtual void advance_timestep() libmesh_override
This method advances the solution to the next timestep, after a solve() has been performed.
libmesh_assert(j)
AdaptiveTimeSolver(sys_type &s)
Constructor.
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual bool element_residual(bool get_jacobian, DiffContext &) libmesh_override
This method is passed on to the core_time_solver.
This class provides a specific system class.
Definition: diff_system.h:53
sys_type & _system
A reference to the system we are solving.
Definition: time_solver.h:257
Generic class from which first order UnsteadySolvers should subclass.
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
virtual void init() libmesh_override
The initialization function.
UniquePtr< NumericVector< Number > > solution
Data structure to hold solution values.
Definition: system.h:1523
const NumericVector< Number > & get_vector(const std::string &vec_name) const
Definition: system.C:794
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
virtual void reinit() libmesh_override
The reinitialization function.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual ~AdaptiveTimeSolver()
Destructor.
UniquePtr< UnsteadySolver > core_time_solver
This object is used to take timesteps.
virtual UniquePtr< DiffSolver > & diff_solver() libmesh_override
An implicit linear or nonlinear solver to use at each timestep.
virtual Real error_order() const libmesh_override
This method is passed on to the core_time_solver.
virtual bool nonlocal_residual(bool get_jacobian, DiffContext &) libmesh_override
This method is passed on to the core_time_solver.
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")