libMesh
diff_context.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_context.h"
20 #include "libmesh/diff_system.h"
21 #include "libmesh/diff_system.h"
22 #include "libmesh/unsteady_solver.h"
23 
24 namespace libMesh
25 {
26 
27 
28 
30  time(sys.time),
31  system_time(sys.time),
32  elem_solution_derivative(1.),
33  elem_solution_rate_derivative(1.),
34  elem_solution_accel_derivative(1.),
35  fixed_solution_derivative(0.),
36  _dof_indices_var(sys.n_vars()),
37  _deltat(libmesh_nullptr),
38  _system(sys),
39  _is_adjoint(false)
40 {
41  // Finally initialize solution/residual/jacobian data structures
42  unsigned int nv = sys.n_vars();
43 
44  _elem_subsolutions.reserve(nv);
45  _elem_subresiduals.reserve(nv);
46  _elem_subjacobians.resize(nv);
47  _elem_subsolution_rates.reserve(nv);
48  _elem_subsolution_accels.reserve(nv);
49  if (sys.use_fixed_solution)
50  _elem_fixed_subsolutions.reserve(nv);
51 
52  // If the user resizes sys.qoi, it will invalidate us
53 
54  std::size_t n_qoi = sys.qoi.size();
55  _elem_qoi.resize(n_qoi);
56  _elem_qoi_derivative.resize(n_qoi);
57  _elem_qoi_subderivatives.resize(n_qoi);
58  for (std::size_t q=0; q != n_qoi; ++q)
59  _elem_qoi_subderivatives[q].reserve(nv);
60 
61  for (unsigned int i=0; i != nv; ++i)
62  {
65  for (std::size_t q=0; q != n_qoi; ++q)
67  _elem_subjacobians[i].reserve(nv);
68 
69  // Only make space for these if we're using DiffSystem
70  // This is assuming *only* DiffSystem is using elem_solution_rate/accel
71  const DifferentiableSystem * diff_system = dynamic_cast<const DifferentiableSystem *>(&sys);
72  if (diff_system)
73  {
74  // Now, we only need these if the solver is unsteady
75  if (!diff_system->get_time_solver().is_steady())
76  {
78 
79  // We only need accel space if the TimeSolver is second order
80  const UnsteadySolver & time_solver = cast_ref<const UnsteadySolver &>(diff_system->get_time_solver());
81 
82  if (time_solver.time_order() >= 2 || !diff_system->get_second_order_vars().empty())
84  }
85  }
86 
87  if (sys.use_fixed_solution)
88  _elem_fixed_subsolutions.push_back
90 
91  for (unsigned int j=0; j != nv; ++j)
92  {
93  _elem_subjacobians[i].push_back
95  }
96  }
97 }
98 
99 
100 
102 {
103  for (std::size_t i=0; i != _elem_subsolutions.size(); ++i)
104  {
105  delete _elem_subsolutions[i];
106  delete _elem_subresiduals[i];
107  for (std::size_t q=0; q != _elem_qoi_subderivatives.size(); ++q)
108  delete _elem_qoi_subderivatives[q][i];
109  if (!_elem_subsolution_rates.empty())
110  delete _elem_subsolution_rates[i];
111  if (!_elem_subsolution_accels.empty())
112  delete _elem_subsolution_accels[i];
113  if (!_elem_fixed_subsolutions.empty())
114  delete _elem_fixed_subsolutions[i];
115 
116  for (std::size_t j=0; j != _elem_subjacobians[i].size(); ++j)
117  delete _elem_subjacobians[i][j];
118  }
119 
120  // We also need to delete all the DenseSubVectors from the localized_vectors map
121  // localized_vectors iterators
122  std::map<const NumericVector<Number> *, std::pair<DenseVector<Number>, std::vector<DenseSubVector<Number> *>>>::iterator localized_vectors_it = _localized_vectors.begin();
123  std::map<const NumericVector<Number> *, std::pair<DenseVector<Number>, std::vector<DenseSubVector<Number> *>>>::iterator localized_vectors_end = _localized_vectors.end();
124 
125  // Loop over every localized_vector
126  for (; localized_vectors_it != localized_vectors_end; ++localized_vectors_it)
127  {
128  // Grab the DenseSubVector to be deleted
129  std::vector<DenseSubVector<Number> *> & localized_vector_dsv = localized_vectors_it->second.second;
130 
131  // Loop over that vector and delete each entry
132  for (std::size_t i=0; i != localized_vector_dsv.size(); ++i)
133  delete localized_vector_dsv[i];
134  }
135 }
136 
137 
139 {
140  // We may actually want to be able to set this pointer to NULL, so
141  // don't report an error for that.
142  _deltat = dt;
143 }
144 
145 
147 {
149 
150  return *_deltat;
151 }
152 
153 
155 {
156  // Make an empty pair keyed with a reference to this _localized_vector
157  _localized_vectors[&localized_vector] = std::make_pair(DenseVector<Number>(), std::vector<DenseSubVector<Number> *>());
158 
159  unsigned int nv = sys.n_vars();
160 
161  _localized_vectors[&localized_vector].second.reserve(nv);
162 
163  // Fill the DenseSubVector with nv copies of DenseVector
164  for (unsigned int i=0; i != nv; ++i)
165  _localized_vectors[&localized_vector].second.push_back(new DenseSubVector<Number>(_localized_vectors[&localized_vector].first));
166 }
167 
168 
170 {
171  return _localized_vectors[&localized_vector].first;
172 }
173 
174 
176 {
177  std::map<const NumericVector<Number> *, std::pair<DenseVector<Number>, std::vector<DenseSubVector<Number> *>>>::const_iterator
178  localized_vectors_it = _localized_vectors.find(&localized_vector);
179  libmesh_assert(localized_vectors_it != _localized_vectors.end());
180  return localized_vectors_it->second.first;
181 }
182 
183 
185 {
186  return *_localized_vectors[&localized_vector].second[var];
187 }
188 
189 
190 const DenseSubVector<Number> & DiffContext::get_localized_subvector (const NumericVector<Number> & localized_vector, unsigned int var) const
191 {
192  std::map<const NumericVector<Number> *, std::pair<DenseVector<Number>, std::vector<DenseSubVector<Number> *>>>::const_iterator
193  localized_vectors_it = _localized_vectors.find(&localized_vector);
194  libmesh_assert(localized_vectors_it != _localized_vectors.end());
195  return *localized_vectors_it->second.second[var];
196 }
197 
198 } // namespace libMesh
std::vector< DenseSubVector< Number > * > _elem_subsolution_rates
Definition: diff_context.h:572
DenseVector< Number > _elem_solution
Element by element components of nonlinear_solution as adjusted by a time_solver. ...
Definition: diff_context.h:564
const std::set< unsigned int > & get_second_order_vars() const
Definition: diff_physics.h:529
std::vector< Number > _elem_qoi
Element quantity of interest contributions.
Definition: diff_context.h:603
std::map< const NumericVector< Number > *, std::pair< DenseVector< Number >, std::vector< DenseSubVector< Number > * > > > _localized_vectors
Contains pointers to vectors the user has asked to be localized, keyed with pairs of element localize...
Definition: diff_context.h:558
std::vector< std::vector< DenseSubMatrix< Number > * > > _elem_subjacobians
Definition: diff_context.h:615
DenseVector< Number > _elem_fixed_solution
Element by element components of nonlinear_solution at a fixed point in a timestep, for optional use by e.g.
Definition: diff_context.h:586
ImplicitSystem & sys
std::vector< DenseSubVector< Number > * > _elem_subsolution_accels
Definition: diff_context.h:579
std::vector< std::vector< DenseSubVector< Number > * > > _elem_qoi_subderivatives
Definition: diff_context.h:609
std::vector< DenseVector< Number > > _elem_qoi_derivative
Element quantity of interest derivative contributions.
Definition: diff_context.h:608
DenseMatrix< Number > _elem_jacobian
Element jacobian: derivatives of elem_residual with respect to elem_solution.
Definition: diff_context.h:598
Real * _deltat
Default NULL, can optionally be used to point to a timestep value in the System-derived class respons...
Definition: diff_context.h:636
const class libmesh_nullptr_t libmesh_nullptr
const unsigned int n_vars
Definition: tecplot_io.C:68
The libMesh namespace provides an interface to certain functionality in the library.
std::vector< DenseSubVector< Number > * > _elem_subresiduals
Element residual subvectors and Jacobian submatrices.
Definition: diff_context.h:614
DiffContext(const System &)
Constructor.
Definition: diff_context.C:29
libmesh_assert(j)
Defines a dense subvector for use in finite element computations.
virtual ~DiffContext()
Destructor.
Definition: diff_context.C:101
This class provides a specific system class.
Definition: diff_system.h:53
void add_localized_vector(NumericVector< Number > &localized_vector, const System &sys)
Adds a vector to the map of localized vectors.
Definition: diff_context.C:154
virtual bool is_steady() const =0
Is this effectively a steady-state solver?
virtual unsigned int time_order() const =0
std::vector< DenseSubVector< Number > * > _elem_fixed_subsolutions
Definition: diff_context.h:587
std::vector< Number > qoi
Values of the quantities of interest.
Definition: system.h:1553
bool use_fixed_solution
A boolean to be set to true by systems using elem_fixed_solution, for optional use by e...
Definition: system.h:1493
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
Defines a dense submatrix for use in Finite Element-type computations.
DenseSubVector< Number > & get_localized_subvector(const NumericVector< Number > &localized_vector, unsigned int var)
Return a reference to DenseSubVector localization of localized_vector at variable var contained in th...
Definition: diff_context.C:184
This is a generic class that defines a solver to handle time integration of DifferentiableSystems.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
DenseVector< Number > _elem_residual
Element residual vector.
Definition: diff_context.h:592
DenseVector< Number > & get_localized_vector(const NumericVector< Number > &localized_vector)
Return a reference to DenseVector localization of localized_vector contained in the _localized_vector...
Definition: diff_context.C:169
void set_deltat_pointer(Real *dt)
Points the _deltat member of this class at a timestep value stored in the creating System...
Definition: diff_context.C:138
unsigned int n_vars() const
Definition: system.h:2086
std::vector< DenseSubVector< Number > * > _elem_subsolutions
Definition: diff_context.h:565
DenseVector< Number > _elem_solution_accel
Element by element components of du/dt as adjusted by a time_solver.
Definition: diff_context.h:578
TimeSolver & get_time_solver()
Definition: diff_system.h:402
DenseVector< Number > _elem_solution_rate
Element by element components of du/dt as adjusted by a time_solver.
Definition: diff_context.h:571