libMesh
element_qoi_derivative.C
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2024 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 
20 // General libMesh includes
21 #include "libmesh/libmesh_common.h"
22 #include "libmesh/elem.h"
23 #include "libmesh/fe_base.h"
24 #include "libmesh/fem_context.h"
25 #include "libmesh/point.h"
26 #include "libmesh/quadrature.h"
27 
28 // Local includes
29 #include "heatsystem.h"
30 
31 // Bring in everything from the libMesh namespace
32 using namespace libMesh;
33 
34 // We only have one QoI, so we don't bother checking the qois argument
35 // to see if it was requested from us
37  const QoISet & /* qois */)
38 {
39  FEMContext & c = cast_ref<FEMContext &>(context);
40 
41  // First we get some references to cell-specific data that
42  // will be used to assemble the linear system.
43  FEBase * elem_fe = nullptr;
44  c.get_element_fe(0, elem_fe);
45 
46  // Element Jacobian * quadrature weights for interior integration
47  const std::vector<Real> & JxW = elem_fe->get_JxW();
48 
49  // The basis functions for the element
50  const std::vector<std::vector<Real>> & phi = elem_fe->get_phi();
51 
52  // The number of local degrees of freedom in each variable
53  const unsigned int n_T_dofs = c.n_dof_indices(0);
54  unsigned int n_qpoints = c.get_element_qrule().n_points();
55 
56  // Fill the QoI RHS corresponding to this QoI. Since this is the 0th QoI
57  // we fill in the [0][i] subderivatives, i corresponding to the variable index.
58  // Our system has only one variable, so we only have to fill the [0][0] subderivative
60 
61  // A reference to the system context is built with
62  const System & sys = c.get_system();
63 
64  // Get a pointer to the adjoint solution vector
65  NumericVector<Number> & adjoint_solution = const_cast<System &>(sys).get_adjoint_solution(0);
66 
67  // Get the previous adjoint solution values at all the qps
68 
69  std::vector<Number> old_adjoint (n_qpoints, 0);
70 
71  c.interior_values<Number>(0, adjoint_solution, old_adjoint);
72 
73  // Our QoI depends solely on the final time, so there are no QoI contributions.
74  // However, there is a contribution from the adjoint solution timestep, for the
75  // time part of the residual of the adjoint problem
76  // Loop over the qps
77  for (unsigned int qp=0; qp != n_qpoints; qp++)
78  for (unsigned int i=0; i != n_T_dofs; i++)
79  Q(i) += -JxW[qp] * old_adjoint[qp] * phi[i][qp] * (1./(dynamic_cast<const HeatSystem &>(sys).deltat));
80 }
This class provides all data required for a physics package (e.g.
Definition: diff_context.h:55
Data structure for specifying which Quantities of Interest should be calculated in an adjoint or a pa...
Definition: qoi_set.h:45
unsigned int n_dof_indices() const
Total number of dof indices on the element.
Definition: diff_context.h:395
The libMesh namespace provides an interface to certain functionality in the library.
void interior_values(unsigned int var, const NumericVector< Number > &_system_vector, std::vector< OutputType > &interior_values_vector) const
Fills a vector of values of the _system_vector at the all the quadrature points in the current elemen...
Definition: fem_context.C:424
Defines a dense subvector for use in finite element computations.
const System & get_system() const
Accessor for associated system.
Definition: diff_context.h:106
Manages consistently variables, degrees of freedom, and coefficient vectors.
Definition: system.h:96
This class provides all data required for a physics package (e.g.
Definition: fem_context.h:62
unsigned int n_points() const
Definition: quadrature.h:123
void get_element_fe(unsigned int var, FEGenericBase< OutputShape > *&fe) const
Accessor for interior finite element object for variable var for the largest dimension in the mesh...
Definition: fem_context.h:277
virtual void element_qoi_derivative(DiffContext &context, const QoISet &)
Does any work that needs to be done on elem in a quantity of interest derivative assembly loop...
const std::vector< DenseVector< Number > > & get_qoi_derivatives() const
Const accessor for QoI derivatives.
Definition: diff_context.h:329
This class forms the foundation from which generic finite elements may be derived.
const QBase & get_element_qrule() const
Accessor for element interior quadrature rule for the dimension of the current _elem.
Definition: fem_context.h:802