Reference documentation for deal.II version 9.1.1
\(\newcommand{\dealcoloneq}{\mathrel{\vcenter{:}}=}\)
Public Types | Static Public Attributes | List of all members
Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType > Class Template Reference

#include <deal.II/differentiation/ad/ad_helpers.h>

Inheritance diagram for Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >:
[legend]

Public Types

using scalar_type = typename HelperBase< ADNumberTypeCode, ScalarType >::scalar_type
 
using ad_type = typename HelperBase< ADNumberTypeCode, ScalarType >::ad_type
 
- Public Types inherited from Differentiation::AD::HelperBase< ADNumberTypeCode, ScalarType >
using scalar_type = typename AD::NumberTraits< ScalarType, ADNumberTypeCode >::scalar_type
 
using ad_type = typename AD::NumberTraits< ScalarType, ADNumberTypeCode >::ad_type
 

Public Member Functions

Constructor / destructor
 PointLevelFunctionsBase (const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
 
virtual ~PointLevelFunctionsBase ()=default
 
Operations specific to taped mode: Reusing tapes
void set_independent_variables (const std::vector< scalar_type > &values)
 
template<typename ValueType , typename ExtractorType >
void set_independent_variable (const ValueType &value, const ExtractorType &extractor)
 
- Public Member Functions inherited from Differentiation::AD::HelperBase< ADNumberTypeCode, ScalarType >
 HelperBase (const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
 
virtual ~HelperBase ()=default
 
std::size_t n_independent_variables () const
 
std::size_t n_dependent_variables () const
 
void print (std::ostream &stream) const
 
void print_values (std::ostream &stream) const
 
void print_tape_stats (const typename Types< ad_type >::tape_index tape_index, std::ostream &stream) const
 
bool is_recording () const
 
Types< ad_type >::tape_index active_tape_index () const
 
bool is_registered_tape (const typename Types< ad_type >::tape_index tape_index) const
 
void set_tape_buffer_sizes (const typename Types< ad_type >::tape_buffer_sizes obufsize=64 *1024 *1024, const typename Types< ad_type >::tape_buffer_sizes lbufsize=64 *1024 *1024, const typename Types< ad_type >::tape_buffer_sizes vbufsize=64 *1024 *1024, const typename Types< ad_type >::tape_buffer_sizes tbufsize=64 *1024 *1024)
 
bool start_recording_operations (const typename Types< ad_type >::tape_index tape_index, const bool overwrite_tape=false, const bool keep_independent_values=true)
 
void stop_recording_operations (const bool write_tapes_to_file=false)
 
void activate_recorded_tape (const typename Types< ad_type >::tape_index tape_index)
 
bool recorded_tape_requires_retaping (const typename Types< ad_type >::tape_index tape_index) const
 
bool active_tape_requires_retaping () const
 
void clear_active_tape ()
 

Static Public Attributes

static const unsigned int dimension = dim
 

Independent variables

std::vector< bool > symmetric_independent_variables
 
virtual void reset (const unsigned int n_independent_variables=::numbers::invalid_unsigned_int, const unsigned int n_dependent_variables=::numbers::invalid_unsigned_int, const bool clear_registered_tapes=true) override
 
void register_independent_variables (const std::vector< scalar_type > &values)
 
template<typename ValueType , typename ExtractorType >
void register_independent_variable (const ValueType &value, const ExtractorType &extractor)
 
const std::vector< ad_type > & get_sensitive_variables () const
 
template<typename ExtractorType >
internal::Extractor< dim, ExtractorType >::template tensor_type< ad_typeget_sensitive_variables (const ExtractorType &extractor) const
 
void set_sensitivity_value (const unsigned int index, const bool symmetric_component, const scalar_type &value)
 
bool is_symmetric_independent_variable (const unsigned int index) const
 
unsigned int n_symmetric_independent_variables () const
 

Additional Inherited Members

- Static Public Member Functions inherited from Differentiation::AD::HelperBase< ADNumberTypeCode, ScalarType >
static void configure_tapeless_mode (const unsigned int n_independent_variables, const bool ensure_persistent_setting=true)
 
- Protected Member Functions inherited from Differentiation::AD::HelperBase< ADNumberTypeCode, ScalarType >
void activate_tape (const typename Types< ad_type >::tape_index tape_index, const bool read_mode)
 
void reset_registered_independent_variables ()
 
void set_sensitivity_value (const unsigned int index, const scalar_type &value)
 
void mark_independent_variable (const unsigned int index, ad_type &out) const
 
void finalize_sensitive_independent_variables () const
 
void initialize_non_sensitive_independent_variable (const unsigned int index, ad_type &out) const
 
unsigned int n_registered_independent_variables () const
 
void reset_registered_dependent_variables (const bool flag=false)
 
unsigned int n_registered_dependent_variables () const
 
void register_dependent_variable (const unsigned int index, const ad_type &func)
 
- Protected Attributes inherited from Differentiation::AD::HelperBase< ADNumberTypeCode, ScalarType >
TapedDrivers< ad_type, scalar_typetaped_driver
 
TapelessDrivers< ad_type, scalar_typetapeless_driver
 
std::vector< scalar_typeindependent_variable_values
 
std::vector< ad_typeindependent_variables
 
std::vector< bool > registered_independent_variable_values
 
std::vector< bool > registered_marked_independent_variables
 
std::vector< ad_typedependent_variables
 
std::vector< bool > registered_marked_dependent_variables
 

Detailed Description

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType = double>
class Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >

A base helper class that facilitates the evaluation of point-wise defined functions. This is the point-wise counterpart of the CellLevelBase class, and was conceived for computations at a continuum point, or quadrature point, rather than for finite-element level calculations. That being said, the interface to this and the derived classes are sufficiently generic that the dependent function(s) and their argument(s), that are the independent variables, can be interpreted in any manner that the user may choose.

As it offers a field-based interface, this class would typically be used to compute the derivatives of a constitutive law defined at a quadrature point; however, it may also be used in other contexts, such as to compute the linearization of a set of local nonlinear equations.

Warning
ADOL-C does not support the standard threading models used by deal.II, so this class should not be embedded within a multithreaded function when using ADOL-C number types. It is, however, suitable for use in both serial and MPI routines.
Author
Jean-Paul Pelteret, 2016, 2017, 2018

Definition at line 2647 of file ad_helpers.h.

Member Typedef Documentation

◆ scalar_type

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType = double>
using Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::scalar_type = typename HelperBase<ADNumberTypeCode, ScalarType>::scalar_type

Type definition for the floating point number type that is used in, and results from, all computations.

Definition at line 2662 of file ad_helpers.h.

◆ ad_type

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType = double>
using Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::ad_type = typename HelperBase<ADNumberTypeCode, ScalarType>::ad_type

Type definition for the auto-differentiation number type that is used in all computations.

Definition at line 2669 of file ad_helpers.h.

Constructor & Destructor Documentation

◆ PointLevelFunctionsBase()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType >
Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::PointLevelFunctionsBase ( const unsigned int  n_independent_variables,
const unsigned int  n_dependent_variables 
)

The constructor for the class.

Parameters
[in]n_independent_variablesThe number of independent variables that will be used in the definition of the functions that it is desired to compute the sensitivities of. In the computation of \(\mathbf{f}(\mathbf{X})\), this will be the number of inputs \(\mathbf{X}\), i.e., the dimension of the domain space.
[in]n_dependent_variablesThe number of scalar functions to be defined that will have a sensitivity to the given independent variables. In the computation of \(\mathbf{f}(\mathbf{X})\), this will be the number of outputs \(\mathbf{f}\), i.e., the dimension of the image space.

Definition at line 1164 of file ad_helpers.cc.

◆ ~PointLevelFunctionsBase()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType = double>
virtual Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::~PointLevelFunctionsBase ( )
virtualdefault

Destructor

Member Function Documentation

◆ reset()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType >
void Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::reset ( const unsigned int  n_independent_variables = ::numbers::invalid_unsigned_int,
const unsigned int  n_dependent_variables = ::numbers::invalid_unsigned_int,
const bool  clear_registered_tapes = true 
)
overridevirtual

Reset the state of the helper class.

When an instance of an HelperBase is stored as a class member object with the intention to reuse its instance, it may be necessary to reset the state of the object before use. This is because, internally, there is error checking performed to ensure that the correct auto-differentiable data is being tracked and used only when appropriate. This function clears all member data and, therefore, allows the state of all internal flags to be safely reset to their initial state.

In the rare case that the number of independent or dependent variables has changed, this can also reconfigured by passing in the appropriate arguments to the function.

Parameters
[in]n_independent_variablesThe number of independent variables that will be used in the definition of the functions that it is desired to compute the sensitivities of. In the computation of \(\mathbf{f}(\mathbf{X})\), this will be the number of inputs \(\mathbf{X}\), i.e., the dimension of the domain space.
[in]n_dependent_variablesThe number of scalar functions to be defined that will have a sensitivity to the given independent variables. In the computation of \(\mathbf{f}(\mathbf{X})\), this will be the number of outputs \(\mathbf{f}\), i.e., the dimension of the image space.
[in]clear_registered_tapesA flag that indicates the that list of registered_tapes must be cleared. If set to true then the data structure that tracks which tapes have been recorded is cleared as well. It is then expected that any preexisting tapes be re-recorded.
Note
This also resets the active tape number to an invalid number, and deactivates the recording mode for taped variables.

Reimplemented from Differentiation::AD::HelperBase< ADNumberTypeCode, ScalarType >.

Definition at line 1177 of file ad_helpers.cc.

◆ register_independent_variables()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType >
void Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::register_independent_variables ( const std::vector< scalar_type > &  values)

Register the complete set of independent variables \(\mathbf{X}\).

Parameters
[in]valuesA field that defines the values of all independent variables. When considering taped AD numbers with branching functions, to avoid potential issues with branch switching it may be a good idea to choose these values close or equal to those that will be later evaluated and differentiated around.
Note
The input value type must correspond to this class's scalar_type. Depending on the selected ADNumberTypeCode, this may or may not correspond with the ScalarType prescribed as a template argument.
For taped AD numbers, this operation is only valid in recording mode.

Definition at line 1229 of file ad_helpers.cc.

◆ register_independent_variable()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType = double>
template<typename ValueType , typename ExtractorType >
void Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::register_independent_variable ( const ValueType &  value,
const ExtractorType &  extractor 
)

Register the subset of independent variables \(\mathbf{A} \subset \mathbf{X}\).

Parameters
[in]valueA field that defines a number of independent variables. When considering taped AD numbers with branching functions, to avoid potential issues with branch switching it may be a good idea to choose these values close or equal to those that will be later evaluated and differentiated around.
[in]extractorAn extractor associated with the input field variables. This effectively defines which components of the global set of independent variables this field is associated with.
Note
The input value type must correspond to this class's scalar_type. Depending on the selected ADNumberTypeCode, this may or may not correspond with the ScalarType prescribed as a template argument.
The input extractor must correspond to the input ValueType. So, for example, if a value is a rank-1 tensor (i.e. of type Tensor<1,dim,scalar_type>), then the extractor must be an FEValuesExtractors::Vector or FEValuesExtractors::Tensor<1>.
This function may be repeatedly used until a call to finalize_sensitive_independent_variables() or get_sensitive_variables() is made.
For taped AD numbers, this operation is only valid in recording mode.

◆ get_sensitive_variables()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType >
const std::vector< typename PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::ad_type > & Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::get_sensitive_variables ( ) const

Return the complete set of independent variables as represented by auto-differentiable numbers. These are the independent variables \(\mathbf{X}\) at which the dependent values are evaluated and differentiated.

This function indicates to the AD library that implements the auto-differentiable number type that operations performed on these numbers are to be tracked so they are considered "sensitive" variables. This is, therefore, the set of variables with which one would then perform computations, and based on which one can then extract both the value of the function and its derivatives with the member functions below. The values of the components of the returned object are initialized to the values set with register_independent_variable().

Returns
An array of auto-differentiable type numbers.
Note
For taped AD numbers, this operation is only valid in recording mode.

Definition at line 1255 of file ad_helpers.cc.

◆ set_independent_variables()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType >
void Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::set_independent_variables ( const std::vector< scalar_type > &  values)

Set the values for the independent variables \(\mathbf{X}\).

Parameters
[in]valuesA vector that defines the values of all independent variables.
Note
The input value type must correspond to this class's scalar_type. Depending on the selected ADNumberTypeCode, this may or may not correspond with the ScalarType prescribed as a template argument.
If the keep_independent_values flag has been set when HelperBase::start_recording_operations() is called then the tape is immediately usable after creation, and the values of the independent variables set by register_independent_variables() are those at which the function is to be evaluated. In this case, a separate call to this function is not strictly necessary.

Definition at line 1305 of file ad_helpers.cc.

◆ set_independent_variable()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType = double>
template<typename ValueType , typename ExtractorType >
void Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::set_independent_variable ( const ValueType &  value,
const ExtractorType &  extractor 
)

Set the values for a subset of independent variables \(\mathbf{A} \subset \mathbf{X}\).

Parameters
[in]valueA field that defines the values of a number of independent variables.
[in]extractorAn extractor associated with the input field variables. This effectively defines which components of the global set of independent variables this field is associated with.
Note
The input value type must correspond to this class's scalar_type. Depending on the selected ADNumberTypeCode, this may or may not correspond with the ScalarType prescribed as a template argument.
The input extractor must correspond to the input ValueType. So, for example, if a value is a rank-1 tensor (i.e. of type Tensor<1,dim,scalar_type>), then the extractor must be an FEValuesExtractors::Vector or FEValuesExtractors::Tensor<1>.
If the keep_independent_values flag has been set when HelperBase::start_recording_operations() is called then the tape is immediately usable after creation, and the values of the independent variables set by register_independent_variable() are those at which the function is to be evaluated. In this case, a separate call to this function is not strictly necessary.

◆ set_sensitivity_value()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType >
void Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::set_sensitivity_value ( const unsigned int  index,
const bool  symmetric_component,
const scalar_type value 
)
protected

Set the actual value of the independent variable \(X_{i}\).

Parameters
[in]indexThe index in the vector of independent variables.
[in]symmetric_componentMark whether this index relates to a component of a field that has a symmetric counterpart (e.g. if index represents an off-diagonal entry in a symmetric tensor).
[in]valueThe value to set the index'd independent variable to.

Definition at line 1283 of file ad_helpers.cc.

◆ is_symmetric_independent_variable()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType >
bool Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::is_symmetric_independent_variable ( const unsigned int  index) const
protected

Return whether the index'th independent variables is one for which we must take into account symmetry when extracting their gradient or Hessian values.

Definition at line 1201 of file ad_helpers.cc.

◆ n_symmetric_independent_variables()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType >
unsigned int Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::n_symmetric_independent_variables ( ) const
protected

Return the number of independent variables that have been marked as being components of a symmetric field.

Definition at line 1215 of file ad_helpers.cc.

Member Data Documentation

◆ dimension

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType = double>
const unsigned int Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::dimension = dim
static

Type definition for the dimension of the associated input and output tensor types.

Definition at line 2655 of file ad_helpers.h.

◆ symmetric_independent_variables

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType = double>
std::vector<bool> Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >::symmetric_independent_variables
private

The independent variables for which we must take into account symmetry when extracting their gradient or Hessian values.

Definition at line 2927 of file ad_helpers.h.


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