libMesh
Public Member Functions | Public Attributes | List of all members
LaplaceQoI Class Reference

#include <L-qoi.h>

Inheritance diagram for LaplaceQoI:
[legend]

Public Member Functions

 LaplaceQoI ()=default
 
virtual ~LaplaceQoI ()=default
 
virtual void init_qoi_count (System &sys)
 Initialize system qoi. More...
 
virtual void init_context (DiffContext &context)
 Prepares the result of a build_context() call for use. More...
 
virtual void postprocess ()
 
virtual void element_qoi_derivative (DiffContext &context, const QoISet &qois)
 Does any work that needs to be done on elem in a quantity of interest derivative assembly loop, outputting to elem_qoi_derivative. More...
 
virtual void element_qoi (DiffContext &context, const QoISet &qois)
 Does any work that needs to be done on elem in a quantity of interest assembly loop, outputting to elem_qoi. More...
 
virtual std::unique_ptr< DifferentiableQoIclone ()
 Copy of this object. More...
 
virtual void init_qoi (std::vector< Number > &)
 Initialize system qoi. More...
 
void init_qoi (std::vector< Number > &)
 Non-virtual, to try to help deprecated user code catch this change at compile time (if they specified override) More...
 
virtual void clear_qoi ()
 Clear all the data structures associated with the QoI. More...
 
virtual void side_qoi (DiffContext &, const QoISet &)
 Does any work that needs to be done on side of elem in a quantity of interest assembly loop, outputting to elem_qoi. More...
 
virtual void side_qoi_derivative (DiffContext &, const QoISet &)
 Does any work that needs to be done on side of elem in a quantity of interest derivative assembly loop, outputting to elem_qoi_derivative. More...
 
virtual void thread_join (std::vector< Number > &qoi, const std::vector< Number > &other_qoi, const QoISet &qoi_indices)
 Method to combine thread-local qois. More...
 
virtual void parallel_op (const Parallel::Communicator &communicator, std::vector< Number > &sys_qoi, std::vector< Number > &local_qoi, const QoISet &qoi_indices)
 Method to populate system qoi data structure with process-local qoi. More...
 
virtual void finalize_derivative (NumericVector< Number > &derivatives, std::size_t qoi_index)
 Method to finalize qoi derivatives which require more than just a simple sum of element contributions. More...
 

Public Attributes

bool assemble_qoi_sides
 If assemble_qoi_sides is true (it is false by default), the assembly loop for a quantity of interest or its derivatives will loop over domain boundary sides. More...
 
bool assemble_qoi_internal_sides
 If assemble_qoi_internal_sides is true (it is false by default), the assembly loop for a quantity of interest or its derivatives will loop over element sides which do not fall on domain boundaries. More...
 
bool assemble_qoi_elements
 If assemble_qoi_elements is false (it is true by default), the assembly loop for a quantity of interest or its derivatives will skip computing on mesh elements, and will only compute on mesh sides. More...
 

Detailed Description

Definition at line 18 of file L-qoi.h.

Constructor & Destructor Documentation

◆ LaplaceQoI()

LaplaceQoI::LaplaceQoI ( )
default

◆ ~LaplaceQoI()

virtual LaplaceQoI::~LaplaceQoI ( )
virtualdefault

Member Function Documentation

◆ clear_qoi()

virtual void libMesh::DifferentiableQoI::clear_qoi ( )
inlinevirtualinherited

Clear all the data structures associated with the QoI.

Definition at line 91 of file diff_qoi.h.

Referenced by libMesh::DifferentiableSystem::clear().

91 {}

◆ clone()

virtual std::unique_ptr<DifferentiableQoI> LaplaceQoI::clone ( )
inlinevirtual

Copy of this object.

User should override to copy any needed state.

Implements libMesh::DifferentiableQoI.

Definition at line 35 of file L-qoi.h.

36  {
37  return std::make_unique<LaplaceQoI>(*this);
38  }

◆ element_qoi()

void LaplaceQoI::element_qoi ( DiffContext ,
const QoISet  
)
virtual

Does any work that needs to be done on elem in a quantity of interest assembly loop, outputting to elem_qoi.

Only qois included in the supplied QoISet need to be assembled.

Reimplemented from libMesh::DifferentiableQoI.

Definition at line 29 of file L-qoi.C.

References std::abs(), libMesh::FEMContext::get_element_fe(), libMesh::FEMContext::get_element_qrule(), libMesh::DiffContext::get_qois(), libMesh::FEMContext::interior_value(), libMesh::QBase::n_points(), and libMesh::Real.

31 {
32  FEMContext & c = cast_ref<FEMContext &>(context);
33 
34  FEBase * elem_fe = nullptr;
35  c.get_element_fe(0, elem_fe);
36 
37  // Element Jacobian * quadrature weights for interior integration
38  const std::vector<Real> & JxW = elem_fe->get_JxW();
39 
40  const std::vector<Point> & xyz = elem_fe->get_xyz();
41 
42  unsigned int n_qpoints = c.get_element_qrule().n_points();
43 
44  Number dQoI_0 = 0.;
45 
46  // Loop over quadrature points
47  for (unsigned int qp = 0; qp != n_qpoints; qp++)
48  {
49  // Get co-ordinate locations of the current quadrature point
50  const Real xf = xyz[qp](0);
51  const Real yf = xyz[qp](1);
52 
53  // If in the sub-domain omega, add the contribution to the integral R
54  if (std::abs(xf - 0.875) <= 0.125 && std::abs(yf - 0.125) <= 0.125)
55  {
56  // Get the solution value at the quadrature point
57  Number T = c.interior_value(0, qp);
58 
59  // Update the elemental increment dR for each qp
60  dQoI_0 += JxW[qp] * T;
61  }
62  }
63 
64  // Update the computed value of the global functional R, by adding the contribution from this element
65  c.get_qois()[0] = c.get_qois()[0] + dQoI_0;
66 }
const std::vector< Number > & get_qois() const
Const accessor for QoI vector.
Definition: diff_context.h:317
Number interior_value(unsigned int var, unsigned int qp) const
Definition: fem_context.C:404
ADRealEigenVector< T, D, asd > abs(const ADRealEigenVector< T, D, asd > &)
Definition: type_vector.h:57
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
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
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
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

◆ element_qoi_derivative()

void LaplaceQoI::element_qoi_derivative ( DiffContext ,
const QoISet  
)
virtual

Does any work that needs to be done on elem in a quantity of interest derivative assembly loop, outputting to elem_qoi_derivative.

Only qois included in the supplied QoISet need their derivatives assembled.

Reimplemented from libMesh::DifferentiableQoI.

Definition at line 70 of file L-qoi.C.

References std::abs(), libMesh::FEMContext::get_element_fe(), libMesh::FEMContext::get_element_qrule(), libMesh::DiffContext::get_qoi_derivatives(), libMesh::DiffContext::n_dof_indices(), libMesh::QBase::n_points(), and libMesh::Real.

72 {
73  FEMContext & c = cast_ref<FEMContext &>(context);
74 
75  // First we get some references to cell-specific data that
76  // will be used to assemble the linear system.
77  FEBase * elem_fe = nullptr;
78  c.get_element_fe(0, elem_fe);
79 
80  // Element Jacobian * quadrature weights for interior integration
81  const std::vector<Real> & JxW = elem_fe->get_JxW();
82 
83  // The basis functions for the element
84  const std::vector<std::vector<Real>> & phi = elem_fe->get_phi();
85 
86  // The element quadrature points
87  const std::vector<Point > & q_point = elem_fe->get_xyz();
88 
89  // The number of local degrees of freedom in each variable
90  const unsigned int n_T_dofs = c.n_dof_indices(0);
91  unsigned int n_qpoints = c.get_element_qrule().n_points();
92 
93  // Fill the QoI RHS corresponding to this QoI. Since this is the 0th QoI
94  // we fill in the [0][i] subderivatives, i corresponding to the variable index.
95  // Our system has only one variable, so we only have to fill the [0][0] subderivative
97 
98  // Loop over the qps
99  for (unsigned int qp=0; qp != n_qpoints; qp++)
100  {
101  const Real x = q_point[qp](0);
102  const Real y = q_point[qp](1);
103 
104  // If in the sub-domain over which QoI 0 is supported, add contributions
105  // to the adjoint rhs
106  if (std::abs(x - 0.875) <= 0.125 && std::abs(y - 0.125) <= 0.125)
107  for (unsigned int i=0; i != n_T_dofs; i++)
108  Q(i) += JxW[qp]*phi[i][qp];
109  } // end of the quadrature point qp-loop
110 }
unsigned int n_dof_indices() const
Total number of dof indices on the element.
Definition: diff_context.h:395
Defines a dense subvector for use in finite element computations.
ADRealEigenVector< T, D, asd > abs(const ADRealEigenVector< T, D, asd > &)
Definition: type_vector.h:57
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
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
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
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

◆ finalize_derivative()

void libMesh::DifferentiableQoI::finalize_derivative ( NumericVector< Number > &  derivatives,
std::size_t  qoi_index 
)
virtualinherited

Method to finalize qoi derivatives which require more than just a simple sum of element contributions.

Definition at line 53 of file diff_qoi.C.

Referenced by libMesh::FEMSystem::assemble_qoi_derivative().

54 {
55  // by default, do nothing
56 }

◆ init_context()

void LaplaceQoI::init_context ( DiffContext )
virtual

Prepares the result of a build_context() call for use.

FEMSystem-based problems will need to reimplement this in order to call FE::get_*() as their particular QoI requires. Trying to evaluate a QoI without overriding init_context is both inefficient and deprecated.

Reimplemented from libMesh::DifferentiableQoI.

Definition at line 13 of file L-qoi.C.

References libMesh::FEMContext::get_element_fe().

14 {
15  FEMContext & c = cast_ref<FEMContext &>(context);
16 
17  // Now make sure we have requested all the data
18  // we need to build the linear system.
19  FEBase * elem_fe = nullptr;
20  c.get_element_fe(0, elem_fe);
21  elem_fe->get_JxW();
22  elem_fe->get_phi();
23  elem_fe->get_xyz();
24 }
This class provides all data required for a physics package (e.g.
Definition: fem_context.h:62
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
This class forms the foundation from which generic finite elements may be derived.

◆ init_qoi() [1/2]

virtual void libMesh::DifferentiableQoI::init_qoi ( std::vector< Number > &  )
inlinevirtualinherited

Initialize system qoi.

This version of the function required direct vector access, and is now deprecated.

Definition at line 72 of file diff_qoi.h.

72 {}

◆ init_qoi() [2/2]

void libMesh::DifferentiableQoI::init_qoi ( std::vector< Number > &  )
inlineinherited

Non-virtual, to try to help deprecated user code catch this change at compile time (if they specified override)

Definition at line 78 of file diff_qoi.h.

78 {}

◆ init_qoi_count()

void LaplaceQoI::init_qoi_count ( System )
virtual

Initialize system qoi.

Often this will just call sys.init_qois(some_desired_number_of_qois)

Reimplemented from libMesh::DifferentiableQoI.

Definition at line 7 of file L-qoi.C.

References libMesh::System::init_qois().

8 {
9  // Only 1 qoi to worry about
10  sys.init_qois(1);
11 }

◆ parallel_op()

void libMesh::DifferentiableQoI::parallel_op ( const Parallel::Communicator communicator,
std::vector< Number > &  sys_qoi,
std::vector< Number > &  local_qoi,
const QoISet qoi_indices 
)
virtualinherited

Method to populate system qoi data structure with process-local qoi.

By default, simply sums process qois into system qoi.

Definition at line 41 of file diff_qoi.C.

References communicator.

Referenced by libMesh::FEMSystem::assemble_qoi().

45 {
46  // Sum everything into local_qoi
47  communicator.sum(local_qoi);
48 
49  // Now put into system qoi
50  sys_qoi = local_qoi;
51 }
DIE A HORRIBLE DEATH HERE typedef MPI_Comm communicator

◆ postprocess()

virtual void LaplaceQoI::postprocess ( void  )
inlinevirtual

Definition at line 29 of file L-qoi.h.

29 {}

◆ side_qoi()

virtual void libMesh::DifferentiableQoI::side_qoi ( DiffContext ,
const QoISet  
)
inlinevirtualinherited

Does any work that needs to be done on side of elem in a quantity of interest assembly loop, outputting to elem_qoi.

Only qois included in the supplied QoISet need to be assembled.

Reimplemented in CoupledSystemQoI.

Definition at line 147 of file diff_qoi.h.

149  {}

◆ side_qoi_derivative()

virtual void libMesh::DifferentiableQoI::side_qoi_derivative ( DiffContext ,
const QoISet  
)
inlinevirtualinherited

Does any work that needs to be done on side of elem in a quantity of interest derivative assembly loop, outputting to elem_qoi_derivative.

Only qois included in the supplied QoISet need their derivatives assembled.

Reimplemented in LaplaceSystem, LaplaceSystem, and CoupledSystemQoI.

Definition at line 159 of file diff_qoi.h.

161  {}

◆ thread_join()

void libMesh::DifferentiableQoI::thread_join ( std::vector< Number > &  qoi,
const std::vector< Number > &  other_qoi,
const QoISet qoi_indices 
)
virtualinherited

Method to combine thread-local qois.

By default, simply sums thread qois.

Definition at line 33 of file diff_qoi.C.

References libMesh::index_range().

36 {
37  for (auto i : index_range(qoi))
38  qoi[i] += other_qoi[i];
39 }
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
Definition: int_range.h:111

Member Data Documentation

◆ assemble_qoi_elements

bool libMesh::DifferentiableQoI::assemble_qoi_elements
inherited

If assemble_qoi_elements is false (it is true by default), the assembly loop for a quantity of interest or its derivatives will skip computing on mesh elements, and will only compute on mesh sides.

Definition at line 115 of file diff_qoi.h.

◆ assemble_qoi_internal_sides

bool libMesh::DifferentiableQoI::assemble_qoi_internal_sides
inherited

If assemble_qoi_internal_sides is true (it is false by default), the assembly loop for a quantity of interest or its derivatives will loop over element sides which do not fall on domain boundaries.

Definition at line 107 of file diff_qoi.h.

◆ assemble_qoi_sides

bool libMesh::DifferentiableQoI::assemble_qoi_sides
inherited

If assemble_qoi_sides is true (it is false by default), the assembly loop for a quantity of interest or its derivatives will loop over domain boundary sides.

To add domain interior sides, also set assemble_qoi_internal_sides to true.

Definition at line 99 of file diff_qoi.h.

Referenced by main().


The documentation for this class was generated from the following files: