libMesh
Public Types | Public Member Functions | Public Attributes | Protected Attributes | Private Attributes | List of all members
libMesh::DiffContext Class Reference

This class provides all data required for a physics package (e.g. More...

#include <diff_context.h>

Inheritance diagram for libMesh::DiffContext:
[legend]

Public Types

typedef std::map< const NumericVector< Number > *, std::pair< DenseVector< Number >, std::vector< DenseSubVector< Number > > > >::iterator localized_vectors_iterator
 Typedef for the localized_vectors iterator. More...
 

Public Member Functions

 DiffContext (const System &, bool allocate_local_matrices=true)
 Constructor. More...
 
virtual ~DiffContext ()
 Destructor. More...
 
virtual void elem_reinit (Real)
 Gives derived classes the opportunity to reinitialize data (FE objects in FEMSystem, for example) needed for an interior integration at a new point within a timestep. More...
 
virtual void elem_side_reinit (Real)
 Gives derived classes the opportunity to reinitialize data needed for a side integration at a new point within a timestep. More...
 
virtual void elem_edge_reinit (Real)
 Gives derived classes the opportunity to reinitialize data needed for an edge integration at a new point within a timestep. More...
 
virtual void nonlocal_reinit (Real)
 Gives derived classes the opportunity to reinitialize data needed for nonlocal calculations at a new point within a timestep. More...
 
unsigned int n_vars () const
 Number of variables in solution. More...
 
const Systemget_system () const
 Accessor for associated system. More...
 
const DenseVector< Number > & get_elem_solution () const
 Accessor for element solution. More...
 
DenseVector< Number > & get_elem_solution ()
 Non-const accessor for element solution. More...
 
const DenseSubVector< Number > & get_elem_solution (unsigned int var) const
 Accessor for element solution of a particular variable corresponding to the variable index argument. More...
 
DenseSubVector< Number > & get_elem_solution (unsigned int var)
 Accessor for element solution of a particular variable corresponding to the variable index argument. More...
 
const DenseVector< Number > & get_elem_solution_rate () const
 Accessor for element solution rate of change w.r.t. More...
 
DenseVector< Number > & get_elem_solution_rate ()
 Non-const accessor for element solution rate of change w.r.t. More...
 
const DenseSubVector< Number > & get_elem_solution_rate (unsigned int var) const
 Accessor for element solution rate for a particular variable corresponding to the variable index argument. More...
 
DenseSubVector< Number > & get_elem_solution_rate (unsigned int var)
 Accessor for element solution rate for a particular variable corresponding to the variable index argument. More...
 
const DenseVector< Number > & get_elem_solution_accel () const
 Accessor for element solution accel of change w.r.t. More...
 
DenseVector< Number > & get_elem_solution_accel ()
 Non-const accessor for element solution accel of change w.r.t. More...
 
const DenseSubVector< Number > & get_elem_solution_accel (unsigned int var) const
 Accessor for element solution accel for a particular variable corresponding to the variable index argument. More...
 
DenseSubVector< Number > & get_elem_solution_accel (unsigned int var)
 Accessor for element solution accel for a particular variable corresponding to the variable index argument. More...
 
const DenseVector< Number > & get_elem_fixed_solution () const
 Accessor for element fixed solution. More...
 
DenseVector< Number > & get_elem_fixed_solution ()
 Non-const accessor for element fixed solution. More...
 
const DenseSubVector< Number > & get_elem_fixed_solution (unsigned int var) const
 Accessor for element fixed solution of a particular variable corresponding to the variable index argument. More...
 
DenseSubVector< Number > & get_elem_fixed_solution (unsigned int var)
 Accessor for element fixed solution of a particular variable corresponding to the variable index argument. More...
 
const DenseVector< Number > & get_elem_residual () const
 Const accessor for element residual. More...
 
DenseVector< Number > & get_elem_residual ()
 Non-const accessor for element residual. More...
 
const DenseSubVector< Number > & get_elem_residual (unsigned int var) const
 Const accessor for element residual of a particular variable corresponding to the variable index argument. More...
 
DenseSubVector< Number > & get_elem_residual (unsigned int var)
 Non-const accessor for element residual of a particular variable corresponding to the variable index argument. More...
 
const DenseMatrix< Number > & get_elem_jacobian () const
 Const accessor for element Jacobian. More...
 
DenseMatrix< Number > & get_elem_jacobian ()
 Non-const accessor for element Jacobian. More...
 
const DenseSubMatrix< Number > & get_elem_jacobian (unsigned int var1, unsigned int var2) const
 Const accessor for element Jacobian of particular variables corresponding to the variable index arguments. More...
 
DenseSubMatrix< Number > & get_elem_jacobian (unsigned int var1, unsigned int var2)
 Non-const accessor for element Jacobian of particular variables corresponding to the variable index arguments. More...
 
const std::vector< Number > & get_qois () const
 Const accessor for QoI vector. More...
 
std::vector< Number > & get_qois ()
 Non-const accessor for QoI vector. More...
 
const std::vector< DenseVector< Number > > & get_qoi_derivatives () const
 Const accessor for QoI derivatives. More...
 
std::vector< DenseVector< Number > > & get_qoi_derivatives ()
 Non-const accessor for QoI derivatives. More...
 
const DenseSubVector< Number > & get_qoi_derivatives (std::size_t qoi, unsigned int var) const
 Const accessor for QoI derivative of a particular qoi and variable corresponding to the index arguments. More...
 
DenseSubVector< Number > & get_qoi_derivatives (std::size_t qoi, unsigned int var)
 Non-const accessor for QoI derivative of a particular qoi and variable corresponding to the index arguments. More...
 
const std::vector< dof_id_type > & get_dof_indices () const
 Accessor for element dof indices. More...
 
std::vector< dof_id_type > & get_dof_indices ()
 Non-const accessor for element dof indices. More...
 
const std::vector< dof_id_type > & get_dof_indices (unsigned int var) const
 Accessor for element dof indices of a particular variable corresponding to the index argument. More...
 
std::vector< dof_id_type > & get_dof_indices (unsigned int var)
 Accessor for element dof indices of a particular variable corresponding to the index argument. More...
 
unsigned int n_dof_indices () const
 Total number of dof indices on the element. More...
 
unsigned int n_dof_indices (unsigned int var) const
 Total number of dof indices of the particular variable corresponding to the index argument. More...
 
Real get_system_time () const
 Accessor for the time variable stored in the system class. More...
 
Real get_time () const
 Accessor for the time for which the current nonlinear_solution is defined. More...
 
void set_time (Real time_in)
 Set the time for which the current nonlinear_solution is defined. More...
 
Real get_elem_solution_derivative () const
 The derivative of the current elem_solution w.r.t. More...
 
Real get_elem_solution_rate_derivative () const
 The derivative of the current elem_solution_rate w.r.t. More...
 
Real get_elem_solution_accel_derivative () const
 The derivative of the current elem_solution_accel w.r.t. More...
 
Real get_fixed_solution_derivative () const
 The derivative of the current fixed_elem_solution w.r.t. More...
 
bool is_adjoint () const
 Accessor for querying whether we need to do a primal or adjoint solve. More...
 
bool & is_adjoint ()
 Accessor for setting whether we need to do a primal or adjoint solve. More...
 
void set_deltat_pointer (Real *dt)
 Points the _deltat member of this class at a timestep value stored in the creating System, for example DiffSystem::deltat. More...
 
Real get_deltat_value ()
 
void add_localized_vector (NumericVector< Number > &localized_vector, const System &sys)
 Adds a vector to the map of localized vectors. More...
 
DenseVector< Number > & get_localized_vector (const NumericVector< Number > &localized_vector)
 Return a reference to DenseVector localization of localized_vector contained in the _localized_vectors map. More...
 
const DenseVector< Number > & get_localized_vector (const NumericVector< Number > &localized_vector) const
 const accessible version of get_localized_vector function More...
 
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 the _localized_vectors map. More...
 
const DenseSubVector< Number > & get_localized_subvector (const NumericVector< Number > &localized_vector, unsigned int var) const
 const accessible version of get_localized_subvector function More...
 

Public Attributes

Real time
 For time-dependent problems, this is the time t for which the current nonlinear_solution is defined. More...
 
const Real system_time
 This is the time stored in the System class at the time this context was created, i.e. More...
 
Real elem_solution_derivative
 The derivative of elem_solution with respect to the current nonlinear solution. More...
 
Real elem_solution_rate_derivative
 The derivative of elem_solution_rate with respect to the current nonlinear solution, for use by systems with non default mass_residual terms. More...
 
Real elem_solution_accel_derivative
 The derivative of elem_solution_accel with respect to the current nonlinear solution, for use by systems with non default mass_residual terms. More...
 
Real fixed_solution_derivative
 The derivative of elem_fixed_solution with respect to the nonlinear solution, for use by systems constructing jacobians with elem_fixed_solution based methods. More...
 

Protected Attributes

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 localized versions of that vector and per variable views. More...
 
const bool _have_local_matrices
 Whether we have local matrices allocated/initialized. More...
 
DenseVector< Number_elem_solution
 Element by element components of nonlinear_solution as adjusted by a time_solver. More...
 
std::vector< DenseSubVector< Number > > _elem_subsolutions
 
DenseVector< Number_elem_solution_rate
 Element by element components of du/dt as adjusted by a time_solver. More...
 
std::vector< DenseSubVector< Number > > _elem_subsolution_rates
 
DenseVector< Number_elem_solution_accel
 Element by element components of du/dt as adjusted by a time_solver. More...
 
std::vector< DenseSubVector< Number > > _elem_subsolution_accels
 
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. More...
 
std::vector< DenseSubVector< Number > > _elem_fixed_subsolutions
 
DenseVector< Number_elem_residual
 Element residual vector. More...
 
DenseMatrix< Number_elem_jacobian
 Element jacobian: derivatives of elem_residual with respect to elem_solution. More...
 
std::vector< Number_elem_qoi
 Element quantity of interest contributions. More...
 
std::vector< DenseVector< Number > > _elem_qoi_derivative
 Element quantity of interest derivative contributions. More...
 
std::vector< std::vector< DenseSubVector< Number > > > _elem_qoi_subderivatives
 
std::vector< DenseSubVector< Number > > _elem_subresiduals
 Element residual subvectors and (if _have_local_matrices) Jacobian submatrices. More...
 
std::vector< std::vector< DenseSubMatrix< Number > > > _elem_subjacobians
 
std::vector< dof_id_type_dof_indices
 Global Degree of freedom index lists. More...
 
std::vector< std::vector< dof_id_type > > _dof_indices_var
 

Private Attributes

Real_deltat
 Defaults to nullptr, can optionally be used to point to a timestep value in the System-derived class responsible for creating this DiffContext. More...
 
const System_system
 A reference to the system this context is constructed with. More...
 
bool _is_adjoint
 Is this context to be used for a primal or adjoint solve? More...
 

Detailed Description

This class provides all data required for a physics package (e.g.

a DifferentiableSystem subclass) to perform local element residual and jacobian integrations.

This class is part of the new DifferentiableSystem framework, which is still experimental. Users of this framework should beware of bugs and future API changes.

Author
Roy H. Stogner
Date
2009

Definition at line 55 of file diff_context.h.

Member Typedef Documentation

◆ localized_vectors_iterator

typedef std::map<const NumericVector<Number> *, std::pair<DenseVector<Number>, std::vector<DenseSubVector<Number> > > >::iterator libMesh::DiffContext::localized_vectors_iterator

Typedef for the localized_vectors iterator.

Definition at line 540 of file diff_context.h.

Constructor & Destructor Documentation

◆ DiffContext()

libMesh::DiffContext::DiffContext ( const System sys,
bool  allocate_local_matrices = true 
)
explicit

Constructor.

Optionally initializes required data structures.

Definition at line 32 of file diff_context.C.

References _elem_fixed_solution, _elem_fixed_subsolutions, _elem_jacobian, _elem_qoi, _elem_qoi_derivative, _elem_qoi_subderivatives, _elem_residual, _elem_solution, _elem_solution_accel, _elem_solution_rate, _elem_subjacobians, _elem_subresiduals, _elem_subsolution_accels, _elem_subsolution_rates, _elem_subsolutions, libMesh::DifferentiablePhysics::get_second_order_vars(), libMesh::DifferentiableSystem::get_time_solver(), libMesh::TimeSolver::is_steady(), libMesh::System::n_qois(), libMesh::System::n_vars(), libMesh::UnsteadySolver::time_order(), and libMesh::System::use_fixed_solution.

33  :
34  time(sys.time),
35  system_time(sys.time),
40  _have_local_matrices(allocate_local_matrices),
41  _dof_indices_var(sys.n_vars()),
42  _deltat(nullptr),
43  _system(sys),
44  _is_adjoint(false)
45 {
46  // Finally initialize solution/residual/jacobian data structures
47  unsigned int nv = sys.n_vars();
48 
49  _elem_subsolutions.reserve(nv);
50  _elem_subresiduals.reserve(nv);
51  _elem_subsolution_rates.reserve(nv);
52  _elem_subsolution_accels.reserve(nv);
53 
54  if (sys.use_fixed_solution)
55  _elem_fixed_subsolutions.reserve(nv);
56 
57  if (allocate_local_matrices)
58  _elem_subjacobians.resize(nv);
59 
60  // If the user resizes sys.qoi, it will invalidate us
61  std::size_t n_qoi = sys.n_qois();
62  _elem_qoi.resize(n_qoi);
63  _elem_qoi_derivative.resize(n_qoi);
64  _elem_qoi_subderivatives.resize(n_qoi);
65  for (std::size_t q=0; q != n_qoi; ++q)
66  _elem_qoi_subderivatives[q].reserve(nv);
67 
68  for (unsigned int i=0; i != nv; ++i)
69  {
70  _elem_subsolutions.emplace_back(_elem_solution);
71  _elem_subresiduals.emplace_back(_elem_residual);
72  for (std::size_t q=0; q != n_qoi; ++q)
74  if (allocate_local_matrices)
75  _elem_subjacobians[i].reserve(nv);
76 
77  // Only make space for these if we're using DiffSystem
78  // This is assuming *only* DiffSystem is using elem_solution_rate/accel
79  const DifferentiableSystem * diff_system = dynamic_cast<const DifferentiableSystem *>(&sys);
80  if (diff_system)
81  {
82  // Now, we only need these if the solver is unsteady
83  if (!diff_system->get_time_solver().is_steady())
84  {
86 
87  // We only need accel space if the TimeSolver is second order
88  const UnsteadySolver & time_solver = cast_ref<const UnsteadySolver &>(diff_system->get_time_solver());
89 
90  if (time_solver.time_order() >= 2 || !diff_system->get_second_order_vars().empty())
92  }
93  }
94 
95  if (sys.use_fixed_solution)
97 
98  if (allocate_local_matrices)
99  for (unsigned int j=0; j != nv; ++j)
100  _elem_subjacobians[i].emplace_back(_elem_jacobian);
101  }
102 }
std::vector< std::vector< dof_id_type > > _dof_indices_var
Definition: diff_context.h:640
DenseVector< Number > _elem_solution
Element by element components of nonlinear_solution as adjusted by a time_solver. ...
Definition: diff_context.h:582
std::vector< Number > _elem_qoi
Element quantity of interest contributions.
Definition: diff_context.h:621
Real fixed_solution_derivative
The derivative of elem_fixed_solution with respect to the nonlinear solution, for use by systems cons...
Definition: diff_context.h:517
std::vector< DenseSubVector< Number > > _elem_fixed_subsolutions
Definition: diff_context.h:605
std::vector< DenseSubVector< Number > > _elem_subsolution_accels
Definition: diff_context.h:597
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:604
std::vector< DenseVector< Number > > _elem_qoi_derivative
Element quantity of interest derivative contributions.
Definition: diff_context.h:626
DenseMatrix< Number > _elem_jacobian
Element jacobian: derivatives of elem_residual with respect to elem_solution.
Definition: diff_context.h:616
Real * _deltat
Defaults to nullptr, can optionally be used to point to a timestep value in the System-derived class ...
Definition: diff_context.h:655
std::vector< std::vector< DenseSubVector< Number > > > _elem_qoi_subderivatives
Definition: diff_context.h:627
const bool _have_local_matrices
Whether we have local matrices allocated/initialized.
Definition: diff_context.h:576
std::vector< DenseSubVector< Number > > _elem_subsolution_rates
Definition: diff_context.h:590
Real elem_solution_rate_derivative
The derivative of elem_solution_rate with respect to the current nonlinear solution, for use by systems with non default mass_residual terms.
Definition: diff_context.h:503
Real elem_solution_derivative
The derivative of elem_solution with respect to the current nonlinear solution.
Definition: diff_context.h:496
const System & _system
A reference to the system this context is constructed with.
Definition: diff_context.h:660
const Real system_time
This is the time stored in the System class at the time this context was created, i...
Definition: diff_context.h:490
std::vector< DenseSubVector< Number > > _elem_subsolutions
Definition: diff_context.h:583
std::vector< std::vector< DenseSubMatrix< Number > > > _elem_subjacobians
Definition: diff_context.h:634
DenseVector< Number > _elem_residual
Element residual vector.
Definition: diff_context.h:610
bool _is_adjoint
Is this context to be used for a primal or adjoint solve?
Definition: diff_context.h:665
DenseVector< Number > _elem_solution_accel
Element by element components of du/dt as adjusted by a time_solver.
Definition: diff_context.h:596
Real elem_solution_accel_derivative
The derivative of elem_solution_accel with respect to the current nonlinear solution, for use by systems with non default mass_residual terms.
Definition: diff_context.h:510
std::vector< DenseSubVector< Number > > _elem_subresiduals
Element residual subvectors and (if _have_local_matrices) Jacobian submatrices.
Definition: diff_context.h:633
Real time
For time-dependent problems, this is the time t for which the current nonlinear_solution is defined...
Definition: diff_context.h:481
DenseVector< Number > _elem_solution_rate
Element by element components of du/dt as adjusted by a time_solver.
Definition: diff_context.h:589

◆ ~DiffContext()

libMesh::DiffContext::~DiffContext ( )
virtualdefault

Destructor.

Member Function Documentation

◆ add_localized_vector()

void libMesh::DiffContext::add_localized_vector ( NumericVector< Number > &  localized_vector,
const System sys 
)

Adds a vector to the map of localized vectors.

We can later evaluate interior_values, interior_gradients and side_values for these fields these vectors represent.

Definition at line 125 of file diff_context.C.

References _localized_vectors, and libMesh::System::n_vars().

126 {
127  // Make an empty pair keyed with a reference to this _localized_vector
128  _localized_vectors[&localized_vector] = std::make_pair(DenseVector<Number>(), std::vector<DenseSubVector<Number>>());
129 
130  unsigned int nv = sys.n_vars();
131 
132  _localized_vectors[&localized_vector].second.reserve(nv);
133 
134  // Fill the DenseSubVector with nv copies of DenseVector
135  for (unsigned int i=0; i != nv; ++i)
136  _localized_vectors[&localized_vector].second.emplace_back(_localized_vectors[&localized_vector].first);
137 }
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:571

◆ elem_edge_reinit()

virtual void libMesh::DiffContext::elem_edge_reinit ( Real  )
inlinevirtual

Gives derived classes the opportunity to reinitialize data needed for an edge integration at a new point within a timestep.

Reimplemented in libMesh::FEMContext.

Definition at line 89 of file diff_context.h.

89 {}

◆ elem_reinit()

virtual void libMesh::DiffContext::elem_reinit ( Real  )
inlinevirtual

Gives derived classes the opportunity to reinitialize data (FE objects in FEMSystem, for example) needed for an interior integration at a new point within a timestep.

Reimplemented in libMesh::FEMContext.

Definition at line 77 of file diff_context.h.

Referenced by libMesh::EulerSolver::element_residual(), libMesh::Euler2Solver::element_residual(), and libMesh::NewmarkSolver::element_residual().

77 {}

◆ elem_side_reinit()

virtual void libMesh::DiffContext::elem_side_reinit ( Real  )
inlinevirtual

Gives derived classes the opportunity to reinitialize data needed for a side integration at a new point within a timestep.

Reimplemented in libMesh::FEMContext.

Definition at line 83 of file diff_context.h.

Referenced by libMesh::EulerSolver::side_residual(), libMesh::Euler2Solver::side_residual(), and libMesh::NewmarkSolver::side_residual().

83 {}

◆ get_deltat_value()

Real libMesh::DiffContext::get_deltat_value ( )
Returns
The value currently pointed to by this class's _deltat member

Definition at line 117 of file diff_context.C.

References _deltat, and libMesh::libmesh_assert().

Referenced by libMesh::FEMContext::_update_time_from_system().

118 {
120 
121  return *_deltat;
122 }
Real * _deltat
Defaults to nullptr, can optionally be used to point to a timestep value in the System-derived class ...
Definition: diff_context.h:655
libmesh_assert(ctx)

◆ get_dof_indices() [1/4]

const std::vector<dof_id_type>& libMesh::DiffContext::get_dof_indices ( ) const
inline

Accessor for element dof indices.

Definition at line 363 of file diff_context.h.

References _dof_indices.

Referenced by libMesh::EulerSolver::_general_residual(), libMesh::Euler2Solver::_general_residual(), libMesh::NewmarkSolver::_general_residual(), libMesh::RBConstruction::add_scaled_matrix_and_vector(), assembly_with_dg_fem_context(), AssemblyA0::boundary_assembly(), AssemblyA1::boundary_assembly(), AssemblyF0::boundary_assembly(), AssemblyF1::boundary_assembly(), AssemblyA2::boundary_assembly(), AssemblyF2::boundary_assembly(), A2::boundary_assembly(), A3::boundary_assembly(), F0::boundary_assembly(), Output0::boundary_assembly(), libMesh::FirstOrderUnsteadySolver::compute_second_order_eqns(), SecondOrderScalarSystemSecondOrderTimeSolverBase::damping_residual(), SecondOrderScalarSystemFirstOrderTimeSolverBase::damping_residual(), FirstOrderScalarSystemBase::element_time_derivative(), SecondOrderScalarSystemFirstOrderTimeSolverBase::element_time_derivative(), libMesh::FEMContext::fixed_point_gradient(), libMesh::FEMContext::fixed_point_hessian(), libMesh::FEMContext::fixed_point_value(), A0::interior_assembly(), B::interior_assembly(), M0::interior_assembly(), A1::interior_assembly(), AssemblyA0::interior_assembly(), AcousticsInnerProduct::interior_assembly(), AssemblyA1::interior_assembly(), A2::interior_assembly(), AssemblyA2::interior_assembly(), F0::interior_assembly(), OutputAssembly::interior_assembly(), EIM_IP_assembly::interior_assembly(), EIM_F::interior_assembly(), AssemblyEIM::interior_assembly(), InnerProductAssembly::interior_assembly(), AssemblyF0::interior_assembly(), AssemblyF1::interior_assembly(), Ex6InnerProduct::interior_assembly(), libMesh::FEMContext::interior_curl(), libMesh::FEMContext::interior_div(), libMesh::FEMContext::interior_gradients(), libMesh::FEMContext::interior_hessians(), libMesh::FEMContext::interior_values(), libMesh::FEMPhysics::mass_residual(), FirstOrderScalarSystemBase::mass_residual(), SecondOrderScalarSystemSecondOrderTimeSolverBase::mass_residual(), SecondOrderScalarSystemFirstOrderTimeSolverBase::mass_residual(), libMesh::FEMSystem::mesh_position_get(), libMesh::DGFEMContext::neighbor_side_fe_reinit(), libMesh::DifferentiablePhysics::nonlocal_mass_residual(), libMesh::FEMSystem::numerical_jacobian(), libMesh::FEMContext::point_curl(), libMesh::FEMContext::point_gradient(), libMesh::FEMContext::point_hessian(), libMesh::FEMContext::point_value(), libMesh::FEMContext::pre_fe_reinit(), OverlappingCouplingGhostingTest::run_sparsity_pattern_test(), libMesh::FEMContext::side_gradient(), libMesh::FEMContext::side_gradients(), libMesh::FEMContext::side_hessians(), libMesh::FEMContext::side_values(), libMesh::FEMContext::some_gradient(), libMesh::FEMContext::some_hessian(), and libMesh::FEMContext::some_value().

364  { return _dof_indices; }
std::vector< dof_id_type > _dof_indices
Global Degree of freedom index lists.
Definition: diff_context.h:639

◆ get_dof_indices() [2/4]

std::vector<dof_id_type>& libMesh::DiffContext::get_dof_indices ( )
inline

Non-const accessor for element dof indices.

Definition at line 369 of file diff_context.h.

References _dof_indices.

370  { return _dof_indices; }
std::vector< dof_id_type > _dof_indices
Global Degree of freedom index lists.
Definition: diff_context.h:639

◆ get_dof_indices() [3/4]

const std::vector<dof_id_type>& libMesh::DiffContext::get_dof_indices ( unsigned int  var) const
inline

Accessor for element dof indices of a particular variable corresponding to the index argument.

Definition at line 376 of file diff_context.h.

References _dof_indices_var.

377  {
378  libmesh_assert_greater(_dof_indices_var.size(), var);
379  return _dof_indices_var[var];
380  }
std::vector< std::vector< dof_id_type > > _dof_indices_var
Definition: diff_context.h:640

◆ get_dof_indices() [4/4]

std::vector<dof_id_type>& libMesh::DiffContext::get_dof_indices ( unsigned int  var)
inline

Accessor for element dof indices of a particular variable corresponding to the index argument.

Definition at line 386 of file diff_context.h.

References _dof_indices_var.

387  {
388  libmesh_assert_greater(_dof_indices_var.size(), var);
389  return _dof_indices_var[var];
390  }
std::vector< std::vector< dof_id_type > > _dof_indices_var
Definition: diff_context.h:640

◆ get_elem_fixed_solution() [1/4]

const DenseVector<Number>& libMesh::DiffContext::get_elem_fixed_solution ( ) const
inline

◆ get_elem_fixed_solution() [2/4]

DenseVector<Number>& libMesh::DiffContext::get_elem_fixed_solution ( )
inline

Non-const accessor for element fixed solution.

Definition at line 216 of file diff_context.h.

References _elem_fixed_solution.

217  { return _elem_fixed_solution; }
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:604

◆ get_elem_fixed_solution() [3/4]

const DenseSubVector<Number>& libMesh::DiffContext::get_elem_fixed_solution ( unsigned int  var) const
inline

Accessor for element fixed solution of a particular variable corresponding to the variable index argument.

Definition at line 223 of file diff_context.h.

References _elem_fixed_subsolutions.

224  {
225  libmesh_assert_greater(_elem_fixed_subsolutions.size(), var);
226  return _elem_fixed_subsolutions[var];
227  }
std::vector< DenseSubVector< Number > > _elem_fixed_subsolutions
Definition: diff_context.h:605

◆ get_elem_fixed_solution() [4/4]

DenseSubVector<Number>& libMesh::DiffContext::get_elem_fixed_solution ( unsigned int  var)
inline

Accessor for element fixed solution of a particular variable corresponding to the variable index argument.

Definition at line 233 of file diff_context.h.

References _elem_fixed_subsolutions.

234  {
235  libmesh_assert_greater(_elem_fixed_subsolutions.size(), var);
236  return _elem_fixed_subsolutions[var];
237  }
std::vector< DenseSubVector< Number > > _elem_fixed_subsolutions
Definition: diff_context.h:605

◆ get_elem_jacobian() [1/4]

const DenseMatrix<Number>& libMesh::DiffContext::get_elem_jacobian ( ) const
inline

Const accessor for element Jacobian.

Definition at line 274 of file diff_context.h.

References _elem_jacobian, _have_local_matrices, and libMesh::libmesh_assert().

Referenced by libMesh::EulerSolver::_general_residual(), libMesh::Euler2Solver::_general_residual(), libMesh::NewmarkSolver::_general_residual(), libMesh::RBConstruction::add_scaled_matrix_and_vector(), libMesh::FEMSystem::assembly(), assembly_with_dg_fem_context(), AssemblyA0::boundary_assembly(), AssemblyA1::boundary_assembly(), AssemblyA2::boundary_assembly(), A2::boundary_assembly(), A3::boundary_assembly(), libMesh::FirstOrderUnsteadySolver::compute_second_order_eqns(), SecondOrderScalarSystemSecondOrderTimeSolverBase::damping_residual(), SecondOrderScalarSystemFirstOrderTimeSolverBase::damping_residual(), NavierSystem::element_constraint(), CoupledSystem::element_constraint(), libMesh::EigenTimeSolver::element_residual(), NavierSystem::element_time_derivative(), SolidSystem::element_time_derivative(), LaplaceSystem::element_time_derivative(), PoissonSystem::element_time_derivative(), CurlCurlSystem::element_time_derivative(), ElasticitySystem::element_time_derivative(), L2System::element_time_derivative(), CoupledSystem::element_time_derivative(), HeatSystem::element_time_derivative(), B::interior_assembly(), A0::interior_assembly(), M0::interior_assembly(), A1::interior_assembly(), AssemblyA0::interior_assembly(), AcousticsInnerProduct::interior_assembly(), AssemblyA1::interior_assembly(), A2::interior_assembly(), AssemblyA2::interior_assembly(), EIM_IP_assembly::interior_assembly(), AssemblyEIM::interior_assembly(), InnerProductAssembly::interior_assembly(), Ex6InnerProduct::interior_assembly(), NavierSystem::mass_residual(), ElasticitySystem::mass_residual(), libMesh::FEMPhysics::mass_residual(), FirstOrderScalarSystemBase::mass_residual(), SecondOrderScalarSystemSecondOrderTimeSolverBase::mass_residual(), SecondOrderScalarSystemFirstOrderTimeSolverBase::mass_residual(), libMesh::DifferentiablePhysics::nonlocal_mass_residual(), libMesh::EigenTimeSolver::nonlocal_residual(), libMesh::FEMSystem::numerical_jacobian(), libMesh::FEMContext::pre_fe_reinit(), OverlappingCouplingGhostingTest::run_sparsity_pattern_test(), NavierSystem::side_constraint(), LaplaceSystem::side_constraint(), libMesh::EigenTimeSolver::side_residual(), SolidSystem::side_time_derivative(), and CurlCurlSystem::side_time_derivative().

275  {
277  return _elem_jacobian;
278  }
DenseMatrix< Number > _elem_jacobian
Element jacobian: derivatives of elem_residual with respect to elem_solution.
Definition: diff_context.h:616
const bool _have_local_matrices
Whether we have local matrices allocated/initialized.
Definition: diff_context.h:576
libmesh_assert(ctx)

◆ get_elem_jacobian() [2/4]

DenseMatrix<Number>& libMesh::DiffContext::get_elem_jacobian ( )
inline

Non-const accessor for element Jacobian.

Definition at line 283 of file diff_context.h.

References _elem_jacobian, _have_local_matrices, and libMesh::libmesh_assert().

284  {
286  return _elem_jacobian;
287  }
DenseMatrix< Number > _elem_jacobian
Element jacobian: derivatives of elem_residual with respect to elem_solution.
Definition: diff_context.h:616
const bool _have_local_matrices
Whether we have local matrices allocated/initialized.
Definition: diff_context.h:576
libmesh_assert(ctx)

◆ get_elem_jacobian() [3/4]

const DenseSubMatrix<Number>& libMesh::DiffContext::get_elem_jacobian ( unsigned int  var1,
unsigned int  var2 
) const
inline

Const accessor for element Jacobian of particular variables corresponding to the variable index arguments.

Definition at line 293 of file diff_context.h.

References _elem_subjacobians, _have_local_matrices, and libMesh::libmesh_assert().

294  {
296  libmesh_assert_greater(_elem_subjacobians.size(), var1);
297  libmesh_assert_greater(_elem_subjacobians[var1].size(), var2);
298  return _elem_subjacobians[var1][var2];
299  }
const bool _have_local_matrices
Whether we have local matrices allocated/initialized.
Definition: diff_context.h:576
libmesh_assert(ctx)
std::vector< std::vector< DenseSubMatrix< Number > > > _elem_subjacobians
Definition: diff_context.h:634

◆ get_elem_jacobian() [4/4]

DenseSubMatrix<Number>& libMesh::DiffContext::get_elem_jacobian ( unsigned int  var1,
unsigned int  var2 
)
inline

Non-const accessor for element Jacobian of particular variables corresponding to the variable index arguments.

Only available if _have_local_matrices

Definition at line 306 of file diff_context.h.

References _elem_subjacobians, _have_local_matrices, and libMesh::libmesh_assert().

307  {
309  libmesh_assert_greater(_elem_subjacobians.size(), var1);
310  libmesh_assert_greater(_elem_subjacobians[var1].size(), var2);
311  return _elem_subjacobians[var1][var2];
312  }
const bool _have_local_matrices
Whether we have local matrices allocated/initialized.
Definition: diff_context.h:576
libmesh_assert(ctx)
std::vector< std::vector< DenseSubMatrix< Number > > > _elem_subjacobians
Definition: diff_context.h:634

◆ get_elem_residual() [1/4]

const DenseVector<Number>& libMesh::DiffContext::get_elem_residual ( ) const
inline

Const accessor for element residual.

Definition at line 242 of file diff_context.h.

References _elem_residual.

Referenced by libMesh::Euler2Solver::_general_residual(), libMesh::RBConstruction::add_scaled_matrix_and_vector(), libMesh::FEMSystem::assembly(), assembly_with_dg_fem_context(), AssemblyF0::boundary_assembly(), AssemblyF1::boundary_assembly(), AssemblyF2::boundary_assembly(), F0::boundary_assembly(), Output0::boundary_assembly(), libMesh::FirstOrderUnsteadySolver::compute_second_order_eqns(), SecondOrderScalarSystemSecondOrderTimeSolverBase::damping_residual(), SecondOrderScalarSystemFirstOrderTimeSolverBase::damping_residual(), NavierSystem::element_constraint(), CoupledSystem::element_constraint(), NavierSystem::element_time_derivative(), SolidSystem::element_time_derivative(), LaplaceSystem::element_time_derivative(), PoissonSystem::element_time_derivative(), CurlCurlSystem::element_time_derivative(), ElasticitySystem::element_time_derivative(), L2System::element_time_derivative(), CoupledSystem::element_time_derivative(), HeatSystem::element_time_derivative(), FirstOrderScalarSystemBase::element_time_derivative(), SecondOrderScalarSystemFirstOrderTimeSolverBase::element_time_derivative(), F0::interior_assembly(), OutputAssembly::interior_assembly(), EIM_F::interior_assembly(), AssemblyF0::interior_assembly(), AssemblyF1::interior_assembly(), NavierSystem::mass_residual(), ElasticitySystem::mass_residual(), libMesh::FEMPhysics::mass_residual(), FirstOrderScalarSystemBase::mass_residual(), SecondOrderScalarSystemSecondOrderTimeSolverBase::mass_residual(), SecondOrderScalarSystemFirstOrderTimeSolverBase::mass_residual(), libMesh::DifferentiablePhysics::nonlocal_mass_residual(), libMesh::FEMSystem::numerical_jacobian(), libMesh::FEMContext::pre_fe_reinit(), NavierSystem::side_constraint(), LaplaceSystem::side_constraint(), SolidSystem::side_time_derivative(), CurlCurlSystem::side_time_derivative(), and ElasticitySystem::side_time_derivative().

243  { return _elem_residual; }
DenseVector< Number > _elem_residual
Element residual vector.
Definition: diff_context.h:610

◆ get_elem_residual() [2/4]

DenseVector<Number>& libMesh::DiffContext::get_elem_residual ( )
inline

Non-const accessor for element residual.

Definition at line 248 of file diff_context.h.

References _elem_residual.

249  { return _elem_residual; }
DenseVector< Number > _elem_residual
Element residual vector.
Definition: diff_context.h:610

◆ get_elem_residual() [3/4]

const DenseSubVector<Number>& libMesh::DiffContext::get_elem_residual ( unsigned int  var) const
inline

Const accessor for element residual of a particular variable corresponding to the variable index argument.

Definition at line 255 of file diff_context.h.

References _elem_subresiduals.

256  {
257  libmesh_assert_greater(_elem_subresiduals.size(), var);
258  return _elem_subresiduals[var];
259  }
std::vector< DenseSubVector< Number > > _elem_subresiduals
Element residual subvectors and (if _have_local_matrices) Jacobian submatrices.
Definition: diff_context.h:633

◆ get_elem_residual() [4/4]

DenseSubVector<Number>& libMesh::DiffContext::get_elem_residual ( unsigned int  var)
inline

Non-const accessor for element residual of a particular variable corresponding to the variable index argument.

Definition at line 265 of file diff_context.h.

References _elem_subresiduals.

266  {
267  libmesh_assert_greater(_elem_subresiduals.size(), var);
268  return _elem_subresiduals[var];
269  }
std::vector< DenseSubVector< Number > > _elem_subresiduals
Element residual subvectors and (if _have_local_matrices) Jacobian submatrices.
Definition: diff_context.h:633

◆ get_elem_solution() [1/4]

const DenseVector<Number>& libMesh::DiffContext::get_elem_solution ( ) const
inline

◆ get_elem_solution() [2/4]

DenseVector<Number>& libMesh::DiffContext::get_elem_solution ( )
inline

Non-const accessor for element solution.

Definition at line 118 of file diff_context.h.

References _elem_solution.

119  { return _elem_solution; }
DenseVector< Number > _elem_solution
Element by element components of nonlinear_solution as adjusted by a time_solver. ...
Definition: diff_context.h:582

◆ get_elem_solution() [3/4]

const DenseSubVector<Number>& libMesh::DiffContext::get_elem_solution ( unsigned int  var) const
inline

Accessor for element solution of a particular variable corresponding to the variable index argument.

Definition at line 125 of file diff_context.h.

References _elem_subsolutions.

126  {
127  libmesh_assert_greater(_elem_subsolutions.size(), var);
128  return _elem_subsolutions[var];
129  }
std::vector< DenseSubVector< Number > > _elem_subsolutions
Definition: diff_context.h:583

◆ get_elem_solution() [4/4]

DenseSubVector<Number>& libMesh::DiffContext::get_elem_solution ( unsigned int  var)
inline

Accessor for element solution of a particular variable corresponding to the variable index argument.

Definition at line 135 of file diff_context.h.

References _elem_subsolutions.

136  {
137  libmesh_assert_greater(_elem_subsolutions.size(), var);
138  return _elem_subsolutions[var];
139  }
std::vector< DenseSubVector< Number > > _elem_subsolutions
Definition: diff_context.h:583

◆ get_elem_solution_accel() [1/4]

const DenseVector<Number>& libMesh::DiffContext::get_elem_solution_accel ( ) const
inline

Accessor for element solution accel of change w.r.t.

time.

Definition at line 177 of file diff_context.h.

References _elem_solution_accel.

Referenced by libMesh::NewmarkSolver::_general_residual(), libMesh::FEMContext::interior_accel(), libMesh::FEMContext::pre_fe_reinit(), libMesh::FirstOrderUnsteadySolver::prepare_accel(), and libMesh::FEMContext::side_accel().

178  { return _elem_solution_accel; }
DenseVector< Number > _elem_solution_accel
Element by element components of du/dt as adjusted by a time_solver.
Definition: diff_context.h:596

◆ get_elem_solution_accel() [2/4]

DenseVector<Number>& libMesh::DiffContext::get_elem_solution_accel ( )
inline

Non-const accessor for element solution accel of change w.r.t.

time.

Definition at line 184 of file diff_context.h.

References _elem_solution_accel.

185  { return _elem_solution_accel; }
DenseVector< Number > _elem_solution_accel
Element by element components of du/dt as adjusted by a time_solver.
Definition: diff_context.h:596

◆ get_elem_solution_accel() [3/4]

const DenseSubVector<Number>& libMesh::DiffContext::get_elem_solution_accel ( unsigned int  var) const
inline

Accessor for element solution accel for a particular variable corresponding to the variable index argument.

Definition at line 191 of file diff_context.h.

References _elem_subsolution_accels.

192  {
193  libmesh_assert_greater(_elem_subsolution_accels.size(), var);
194  return _elem_subsolution_accels[var];
195  }
std::vector< DenseSubVector< Number > > _elem_subsolution_accels
Definition: diff_context.h:597

◆ get_elem_solution_accel() [4/4]

DenseSubVector<Number>& libMesh::DiffContext::get_elem_solution_accel ( unsigned int  var)
inline

Accessor for element solution accel for a particular variable corresponding to the variable index argument.

Definition at line 201 of file diff_context.h.

References _elem_subsolution_accels.

202  {
203  libmesh_assert_greater(_elem_subsolution_accels.size(), var);
204  return _elem_subsolution_accels[var];
205  }
std::vector< DenseSubVector< Number > > _elem_subsolution_accels
Definition: diff_context.h:597

◆ get_elem_solution_accel_derivative()

Real libMesh::DiffContext::get_elem_solution_accel_derivative ( ) const
inline

The derivative of the current elem_solution_accel w.r.t.

the unknown solution. Corresponding Jacobian contributions should be multiplied by this amount, or may be skipped if get_elem_solution_accel_derivative() is 0.

Definition at line 450 of file diff_context.h.

References elem_solution_accel_derivative.

Referenced by ElasticitySystem::mass_residual(), SecondOrderScalarSystemSecondOrderTimeSolverBase::mass_residual(), and SecondOrderScalarSystemFirstOrderTimeSolverBase::mass_residual().

Real elem_solution_accel_derivative
The derivative of elem_solution_accel with respect to the current nonlinear solution, for use by systems with non default mass_residual terms.
Definition: diff_context.h:510

◆ get_elem_solution_derivative()

Real libMesh::DiffContext::get_elem_solution_derivative ( ) const
inline

The derivative of the current elem_solution w.r.t.

the unknown solution. Corresponding Jacobian contributions should be multiplied by this amount, or may be skipped if get_elem_solution_derivative() is 0.

Definition at line 432 of file diff_context.h.

References elem_solution_derivative.

Referenced by libMesh::FirstOrderUnsteadySolver::compute_second_order_eqns(), NavierSystem::element_constraint(), ElasticitySystem::element_time_derivative(), L2System::element_time_derivative(), HeatSystem::element_time_derivative(), NavierSystem::mass_residual(), and NavierSystem::side_constraint().

433  { return elem_solution_derivative; }
Real elem_solution_derivative
The derivative of elem_solution with respect to the current nonlinear solution.
Definition: diff_context.h:496

◆ get_elem_solution_rate() [1/4]

const DenseVector<Number>& libMesh::DiffContext::get_elem_solution_rate ( ) const
inline

◆ get_elem_solution_rate() [2/4]

DenseVector<Number>& libMesh::DiffContext::get_elem_solution_rate ( )
inline

Non-const accessor for element solution rate of change w.r.t.

time.

Definition at line 151 of file diff_context.h.

References _elem_solution_rate.

152  { return _elem_solution_rate; }
DenseVector< Number > _elem_solution_rate
Element by element components of du/dt as adjusted by a time_solver.
Definition: diff_context.h:589

◆ get_elem_solution_rate() [3/4]

const DenseSubVector<Number>& libMesh::DiffContext::get_elem_solution_rate ( unsigned int  var) const
inline

Accessor for element solution rate for a particular variable corresponding to the variable index argument.

Definition at line 158 of file diff_context.h.

References _elem_subsolution_rates.

159  {
160  libmesh_assert_greater(_elem_subsolution_rates.size(), var);
161  return _elem_subsolution_rates[var];
162  }
std::vector< DenseSubVector< Number > > _elem_subsolution_rates
Definition: diff_context.h:590

◆ get_elem_solution_rate() [4/4]

DenseSubVector<Number>& libMesh::DiffContext::get_elem_solution_rate ( unsigned int  var)
inline

Accessor for element solution rate for a particular variable corresponding to the variable index argument.

Definition at line 168 of file diff_context.h.

References _elem_subsolution_rates.

169  {
170  libmesh_assert_greater(_elem_subsolution_rates.size(), var);
171  return _elem_subsolution_rates[var];
172  }
std::vector< DenseSubVector< Number > > _elem_subsolution_rates
Definition: diff_context.h:590

◆ get_elem_solution_rate_derivative()

Real libMesh::DiffContext::get_elem_solution_rate_derivative ( ) const
inline

The derivative of the current elem_solution_rate w.r.t.

the unknown solution. Corresponding Jacobian contributions should be multiplied by this amount, or may be skipped if get_elem_solution_rate_derivative() is 0.

Definition at line 441 of file diff_context.h.

References elem_solution_rate_derivative.

Referenced by libMesh::FirstOrderUnsteadySolver::compute_second_order_eqns(), SecondOrderScalarSystemSecondOrderTimeSolverBase::damping_residual(), SecondOrderScalarSystemFirstOrderTimeSolverBase::damping_residual(), FirstOrderScalarSystemBase::mass_residual(), and libMesh::FirstOrderUnsteadySolver::prepare_accel().

Real elem_solution_rate_derivative
The derivative of elem_solution_rate with respect to the current nonlinear solution, for use by systems with non default mass_residual terms.
Definition: diff_context.h:503

◆ get_fixed_solution_derivative()

Real libMesh::DiffContext::get_fixed_solution_derivative ( ) const
inline

The derivative of the current fixed_elem_solution w.r.t.

the unknown solution. Corresponding Jacobian contributions should be multiplied by this amount, or may be skipped if get_fixed_elem_solution_derivative() is 0.

Definition at line 459 of file diff_context.h.

References fixed_solution_derivative.

460  { return fixed_solution_derivative; }
Real fixed_solution_derivative
The derivative of elem_fixed_solution with respect to the nonlinear solution, for use by systems cons...
Definition: diff_context.h:517

◆ get_localized_subvector() [1/2]

DenseSubVector< Number > & libMesh::DiffContext::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 the _localized_vectors map.

Definition at line 154 of file diff_context.C.

References _localized_vectors.

Referenced by libMesh::FEMContext::interior_gradients(), libMesh::FEMContext::interior_hessians(), libMesh::FEMContext::interior_values(), libMesh::FEMContext::side_gradients(), libMesh::FEMContext::side_hessians(), and libMesh::FEMContext::side_values().

155 {
156  return _localized_vectors[&localized_vector].second[var];
157 }
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:571

◆ get_localized_subvector() [2/2]

const DenseSubVector< Number > & libMesh::DiffContext::get_localized_subvector ( const NumericVector< Number > &  localized_vector,
unsigned int  var 
) const

const accessible version of get_localized_subvector function

Definition at line 160 of file diff_context.C.

References _localized_vectors, and libMesh::libmesh_assert().

161 {
162  auto localized_vectors_it = _localized_vectors.find(&localized_vector);
163  libmesh_assert(localized_vectors_it != _localized_vectors.end());
164  return localized_vectors_it->second.second[var];
165 }
libmesh_assert(ctx)
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:571

◆ get_localized_vector() [1/2]

DenseVector< Number > & libMesh::DiffContext::get_localized_vector ( const NumericVector< Number > &  localized_vector)

Return a reference to DenseVector localization of localized_vector contained in the _localized_vectors map.

Definition at line 140 of file diff_context.C.

References _localized_vectors.

141 {
142  return _localized_vectors[&localized_vector].first;
143 }
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:571

◆ get_localized_vector() [2/2]

const DenseVector< Number > & libMesh::DiffContext::get_localized_vector ( const NumericVector< Number > &  localized_vector) const

const accessible version of get_localized_vector function

Definition at line 146 of file diff_context.C.

References _localized_vectors, and libMesh::libmesh_assert().

147 {
148  auto localized_vectors_it = _localized_vectors.find(&localized_vector);
149  libmesh_assert(localized_vectors_it != _localized_vectors.end());
150  return localized_vectors_it->second.first;
151 }
libmesh_assert(ctx)
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:571

◆ get_qoi_derivatives() [1/4]

const std::vector<DenseVector<Number> >& libMesh::DiffContext::get_qoi_derivatives ( ) const
inline

Const accessor for QoI derivatives.

Definition at line 329 of file diff_context.h.

References _elem_qoi_derivative.

Referenced by LaplaceQoI::element_qoi_derivative(), LaplaceSystem::element_qoi_derivative(), HeatSystem::element_qoi_derivative(), libMesh::FEMContext::pre_fe_reinit(), CoupledSystemQoI::side_qoi_derivative(), and LaplaceSystem::side_qoi_derivative().

330  { return _elem_qoi_derivative; }
std::vector< DenseVector< Number > > _elem_qoi_derivative
Element quantity of interest derivative contributions.
Definition: diff_context.h:626

◆ get_qoi_derivatives() [2/4]

std::vector<DenseVector<Number> >& libMesh::DiffContext::get_qoi_derivatives ( )
inline

Non-const accessor for QoI derivatives.

Definition at line 335 of file diff_context.h.

References _elem_qoi_derivative.

336  { return _elem_qoi_derivative; }
std::vector< DenseVector< Number > > _elem_qoi_derivative
Element quantity of interest derivative contributions.
Definition: diff_context.h:626

◆ get_qoi_derivatives() [3/4]

const DenseSubVector<Number>& libMesh::DiffContext::get_qoi_derivatives ( std::size_t  qoi,
unsigned int  var 
) const
inline

Const accessor for QoI derivative of a particular qoi and variable corresponding to the index arguments.

Definition at line 342 of file diff_context.h.

References _elem_qoi_subderivatives.

343  {
344  libmesh_assert_greater(_elem_qoi_subderivatives.size(), qoi);
345  libmesh_assert_greater(_elem_qoi_subderivatives[qoi].size(), var);
346  return _elem_qoi_subderivatives[qoi][var];
347  }
std::vector< std::vector< DenseSubVector< Number > > > _elem_qoi_subderivatives
Definition: diff_context.h:627

◆ get_qoi_derivatives() [4/4]

DenseSubVector<Number>& libMesh::DiffContext::get_qoi_derivatives ( std::size_t  qoi,
unsigned int  var 
)
inline

Non-const accessor for QoI derivative of a particular qoi and variable corresponding to the index arguments.

Definition at line 353 of file diff_context.h.

References _elem_qoi_subderivatives.

354  {
355  libmesh_assert_greater(_elem_qoi_subderivatives.size(), qoi);
356  libmesh_assert_greater(_elem_qoi_subderivatives[qoi].size(), var);
357  return _elem_qoi_subderivatives[qoi][var];
358  }
std::vector< std::vector< DenseSubVector< Number > > > _elem_qoi_subderivatives
Definition: diff_context.h:627

◆ get_qois() [1/2]

const std::vector<Number>& libMesh::DiffContext::get_qois ( ) const
inline

Const accessor for QoI vector.

Definition at line 317 of file diff_context.h.

References _elem_qoi.

Referenced by LaplaceQoI::element_qoi(), HeatSystem::element_qoi(), and CoupledSystemQoI::side_qoi().

318  { return _elem_qoi; }
std::vector< Number > _elem_qoi
Element quantity of interest contributions.
Definition: diff_context.h:621

◆ get_qois() [2/2]

std::vector<Number>& libMesh::DiffContext::get_qois ( )
inline

Non-const accessor for QoI vector.

Definition at line 323 of file diff_context.h.

References _elem_qoi.

324  { return _elem_qoi; }
std::vector< Number > _elem_qoi
Element quantity of interest contributions.
Definition: diff_context.h:621

◆ get_system()

const System& libMesh::DiffContext::get_system ( ) const
inline

◆ get_system_time()

Real libMesh::DiffContext::get_system_time ( ) const
inline

Accessor for the time variable stored in the system class.

Definition at line 411 of file diff_context.h.

References system_time.

Referenced by libMesh::FEMContext::_update_time_from_system().

412  { return system_time; }
const Real system_time
This is the time stored in the System class at the time this context was created, i...
Definition: diff_context.h:490

◆ get_time()

Real libMesh::DiffContext::get_time ( ) const
inline

Accessor for the time for which the current nonlinear_solution is defined.

Definition at line 417 of file diff_context.h.

References time.

418  { return time; }
Real time
For time-dependent problems, this is the time t for which the current nonlinear_solution is defined...
Definition: diff_context.h:481

◆ is_adjoint() [1/2]

bool libMesh::DiffContext::is_adjoint ( ) const
inline

Accessor for querying whether we need to do a primal or adjoint solve.

Definition at line 466 of file diff_context.h.

References _is_adjoint.

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

467  { return _is_adjoint; }
bool _is_adjoint
Is this context to be used for a primal or adjoint solve?
Definition: diff_context.h:665

◆ is_adjoint() [2/2]

bool& libMesh::DiffContext::is_adjoint ( )
inline

Accessor for setting whether we need to do a primal or adjoint solve.

Definition at line 473 of file diff_context.h.

References _is_adjoint.

474  { return _is_adjoint; }
bool _is_adjoint
Is this context to be used for a primal or adjoint solve?
Definition: diff_context.h:665

◆ n_dof_indices() [1/2]

unsigned int libMesh::DiffContext::n_dof_indices ( ) const
inline

◆ n_dof_indices() [2/2]

unsigned int libMesh::DiffContext::n_dof_indices ( unsigned int  var) const
inline

Total number of dof indices of the particular variable corresponding to the index argument.

Definition at line 402 of file diff_context.h.

References _dof_indices_var.

403  {
404  libmesh_assert_greater(_dof_indices_var.size(), var);
405  return cast_int<unsigned int>(_dof_indices_var[var].size());
406  }
std::vector< std::vector< dof_id_type > > _dof_indices_var
Definition: diff_context.h:640

◆ n_vars()

unsigned int libMesh::DiffContext::n_vars ( ) const
inline

◆ nonlocal_reinit()

virtual void libMesh::DiffContext::nonlocal_reinit ( Real  )
inlinevirtual

Gives derived classes the opportunity to reinitialize data needed for nonlocal calculations at a new point within a timestep.

Reimplemented in libMesh::FEMContext.

Definition at line 95 of file diff_context.h.

Referenced by libMesh::EulerSolver::nonlocal_residual(), libMesh::Euler2Solver::nonlocal_residual(), and libMesh::NewmarkSolver::nonlocal_residual().

95 {}

◆ set_deltat_pointer()

void libMesh::DiffContext::set_deltat_pointer ( Real dt)

Points the _deltat member of this class at a timestep value stored in the creating System, for example DiffSystem::deltat.

Definition at line 109 of file diff_context.C.

References _deltat.

Referenced by libMesh::FEMSystem::build_context(), libMesh::DifferentiableSystem::build_context(), and libMesh::FEMSystem::init_context().

110 {
111  // We may actually want to be able to set this pointer to nullptr, so
112  // don't report an error for that.
113  _deltat = dt;
114 }
Real * _deltat
Defaults to nullptr, can optionally be used to point to a timestep value in the System-derived class ...
Definition: diff_context.h:655

◆ set_time()

void libMesh::DiffContext::set_time ( Real  time_in)
inline

Set the time for which the current nonlinear_solution is defined.

Definition at line 423 of file diff_context.h.

References time.

Referenced by libMesh::FEMContext::_update_time_from_system().

424  { time = time_in; }
Real time
For time-dependent problems, this is the time t for which the current nonlinear_solution is defined...
Definition: diff_context.h:481

Member Data Documentation

◆ _deltat

Real* libMesh::DiffContext::_deltat
private

Defaults to nullptr, can optionally be used to point to a timestep value in the System-derived class responsible for creating this DiffContext.

In DiffSystem's build_context() function, is assigned to point to the deltat member of that class.

Accessible via public get_deltat()/set_deltat() methods of this class.

Always test for nullptr before using!

Definition at line 655 of file diff_context.h.

Referenced by get_deltat_value(), and set_deltat_pointer().

◆ _dof_indices

std::vector<dof_id_type> libMesh::DiffContext::_dof_indices
protected

Global Degree of freedom index lists.

Definition at line 639 of file diff_context.h.

Referenced by get_dof_indices(), and n_dof_indices().

◆ _dof_indices_var

std::vector<std::vector<dof_id_type> > libMesh::DiffContext::_dof_indices_var
protected

Definition at line 640 of file diff_context.h.

Referenced by get_dof_indices(), n_dof_indices(), and n_vars().

◆ _elem_fixed_solution

DenseVector<Number> libMesh::DiffContext::_elem_fixed_solution
protected

Element by element components of nonlinear_solution at a fixed point in a timestep, for optional use by e.g.

stabilized methods

Definition at line 604 of file diff_context.h.

Referenced by DiffContext(), and get_elem_fixed_solution().

◆ _elem_fixed_subsolutions

std::vector<DenseSubVector<Number> > libMesh::DiffContext::_elem_fixed_subsolutions
protected

◆ _elem_jacobian

DenseMatrix<Number> libMesh::DiffContext::_elem_jacobian
protected

Element jacobian: derivatives of elem_residual with respect to elem_solution.

Only initialized if _have_local_matrices

Definition at line 616 of file diff_context.h.

Referenced by DiffContext(), and get_elem_jacobian().

◆ _elem_qoi

std::vector<Number> libMesh::DiffContext::_elem_qoi
protected

Element quantity of interest contributions.

Definition at line 621 of file diff_context.h.

Referenced by DiffContext(), and get_qois().

◆ _elem_qoi_derivative

std::vector<DenseVector<Number> > libMesh::DiffContext::_elem_qoi_derivative
protected

Element quantity of interest derivative contributions.

Definition at line 626 of file diff_context.h.

Referenced by DiffContext(), and get_qoi_derivatives().

◆ _elem_qoi_subderivatives

std::vector<std::vector<DenseSubVector<Number> > > libMesh::DiffContext::_elem_qoi_subderivatives
protected

◆ _elem_residual

DenseVector<Number> libMesh::DiffContext::_elem_residual
protected

Element residual vector.

Definition at line 610 of file diff_context.h.

Referenced by DiffContext(), and get_elem_residual().

◆ _elem_solution

DenseVector<Number> libMesh::DiffContext::_elem_solution
protected

Element by element components of nonlinear_solution as adjusted by a time_solver.

Definition at line 582 of file diff_context.h.

Referenced by DiffContext(), and get_elem_solution().

◆ _elem_solution_accel

DenseVector<Number> libMesh::DiffContext::_elem_solution_accel
protected

Element by element components of du/dt as adjusted by a time_solver.

Definition at line 596 of file diff_context.h.

Referenced by DiffContext(), and get_elem_solution_accel().

◆ _elem_solution_rate

DenseVector<Number> libMesh::DiffContext::_elem_solution_rate
protected

Element by element components of du/dt as adjusted by a time_solver.

Definition at line 589 of file diff_context.h.

Referenced by DiffContext(), and get_elem_solution_rate().

◆ _elem_subjacobians

std::vector<std::vector<DenseSubMatrix<Number> > > libMesh::DiffContext::_elem_subjacobians
protected

Definition at line 634 of file diff_context.h.

Referenced by DiffContext(), and get_elem_jacobian().

◆ _elem_subresiduals

std::vector<DenseSubVector<Number> > libMesh::DiffContext::_elem_subresiduals
protected

Element residual subvectors and (if _have_local_matrices) Jacobian submatrices.

Definition at line 633 of file diff_context.h.

Referenced by DiffContext(), and get_elem_residual().

◆ _elem_subsolution_accels

std::vector<DenseSubVector<Number> > libMesh::DiffContext::_elem_subsolution_accels
protected

Definition at line 597 of file diff_context.h.

Referenced by DiffContext(), and get_elem_solution_accel().

◆ _elem_subsolution_rates

std::vector<DenseSubVector<Number> > libMesh::DiffContext::_elem_subsolution_rates
protected

Definition at line 590 of file diff_context.h.

Referenced by DiffContext(), and get_elem_solution_rate().

◆ _elem_subsolutions

std::vector<DenseSubVector<Number> > libMesh::DiffContext::_elem_subsolutions
protected

◆ _have_local_matrices

const bool libMesh::DiffContext::_have_local_matrices
protected

Whether we have local matrices allocated/initialized.

Definition at line 576 of file diff_context.h.

Referenced by get_elem_jacobian(), and libMesh::FEMContext::pre_fe_reinit().

◆ _is_adjoint

bool libMesh::DiffContext::_is_adjoint
private

Is this context to be used for a primal or adjoint solve?

Definition at line 665 of file diff_context.h.

Referenced by is_adjoint().

◆ _localized_vectors

std::map<const NumericVector<Number> *, std::pair<DenseVector<Number>, std::vector<DenseSubVector<Number> > > > libMesh::DiffContext::_localized_vectors
protected

Contains pointers to vectors the user has asked to be localized, keyed with pairs of element localized versions of that vector and per variable views.

Definition at line 571 of file diff_context.h.

Referenced by add_localized_vector(), get_localized_subvector(), get_localized_vector(), and libMesh::FEMContext::pre_fe_reinit().

◆ _system

const System& libMesh::DiffContext::_system
private

A reference to the system this context is constructed with.

Definition at line 660 of file diff_context.h.

Referenced by get_system().

◆ elem_solution_accel_derivative

Real libMesh::DiffContext::elem_solution_accel_derivative

The derivative of elem_solution_accel with respect to the current nonlinear solution, for use by systems with non default mass_residual terms.

Definition at line 510 of file diff_context.h.

Referenced by libMesh::NewmarkSolver::_general_residual(), get_elem_solution_accel_derivative(), and libMesh::FirstOrderUnsteadySolver::prepare_accel().

◆ elem_solution_derivative

Real libMesh::DiffContext::elem_solution_derivative

◆ elem_solution_rate_derivative

Real libMesh::DiffContext::elem_solution_rate_derivative

◆ fixed_solution_derivative

Real libMesh::DiffContext::fixed_solution_derivative

The derivative of elem_fixed_solution with respect to the nonlinear solution, for use by systems constructing jacobians with elem_fixed_solution based methods.

Definition at line 517 of file diff_context.h.

Referenced by libMesh::EulerSolver::_general_residual(), libMesh::Euler2Solver::_general_residual(), libMesh::SteadySolver::_general_residual(), and get_fixed_solution_derivative().

◆ system_time

const Real libMesh::DiffContext::system_time

This is the time stored in the System class at the time this context was created, i.e.

the time at the beginning of the current timestep. This value gets set in the constructor and unlike DiffContext::time, is not tweaked mid-timestep by transient solvers: it remains equal to the value it was assigned at construction.

Definition at line 490 of file diff_context.h.

Referenced by get_system_time().

◆ time

Real libMesh::DiffContext::time

For time-dependent problems, this is the time t for which the current nonlinear_solution is defined.

FIXME - this needs to be tweaked mid-timestep by all transient solvers!

Definition at line 481 of file diff_context.h.

Referenced by get_time(), and set_time().


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