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

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

Inheritance diagram for Differentiation::AD::VectorFunction< 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::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >
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
 VectorFunction (const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
 
virtual ~VectorFunction ()=default
 
- Public Member Functions inherited from Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >
 PointLevelFunctionsBase (const unsigned int n_independent_variables, const unsigned int n_dependent_variables)
 
virtual ~PointLevelFunctionsBase ()=default
 
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)
 
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
 
- 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 ()
 

Dependent variables

void register_dependent_variables (const std::vector< ad_type > &funcs)
 
template<typename ValueType , typename ExtractorType >
void register_dependent_variable (const ValueType &funcs, const ExtractorType &extractor)
 
void compute_values (Vector< scalar_type > &values) const
 
void compute_jacobian (FullMatrix< scalar_type > &jacobian) const
 
template<typename ExtractorType_Row >
static internal::VectorFieldValue< dim, scalar_type, ExtractorType_Row >::type extract_value_component (const Vector< scalar_type > &values, const ExtractorType_Row &extractor_row)
 
template<typename ExtractorType_Row , typename ExtractorType_Col >
static internal::VectorFieldJacobian< dim, scalar_type, ExtractorType_Row, ExtractorType_Col >::type extract_jacobian_component (const FullMatrix< scalar_type > &jacobian, const ExtractorType_Row &extractor_row, const ExtractorType_Col &extractor_col)
 
static Tensor< 0, dim, scalar_typeextract_jacobian_component (const FullMatrix< scalar_type > &jacobian, const FEValuesExtractors::Scalar &extractor_row, const FEValuesExtractors::Scalar &extractor_col)
 
static SymmetricTensor< 4, dim, scalar_typeextract_jacobian_component (const FullMatrix< scalar_type > &jacobian, const FEValuesExtractors::SymmetricTensor< 2 > &extractor_row, const FEValuesExtractors::SymmetricTensor< 2 > &extractor_col)
 

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)
 
- Static Public Attributes inherited from Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >
static const unsigned int dimension = dim
 
- Protected Member Functions inherited from Differentiation::AD::PointLevelFunctionsBase< dim, ADNumberTypeCode, ScalarType >
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
 
- 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::VectorFunction< dim, ADNumberTypeCode, ScalarType >

A helper class that facilitates the evaluation of a vector of functions, typically one that represents a collection of coupled, multi-dimensional fields. This class would typically be used to compute the linearization of a set of kinetic field variables defined at the quadrature point level.

An example of its usage in the case of linearizing the kinetic variables derived from a multi-field constitutive law might be as follows:

// Define some extractors that will help us set independent variables
// and later get the computed values related to the dependent
// variables. Each of these extractors is related to the gradient of a
// component of the solution field (in this case, displacement and
// magnetic scalar potential). Here "C" is the right Cauchy-Green
// tensor and "H" is the magnetic field.
const unsigned int n_independent_variables =
// Declare how many dependent variables we expect to compute.
// In this case, we will be computing a stress field (a symmetric
// rank-2 tensor) and the magnetic induction (a vector field).
// At the same time we define some additional extractors associated
// with these kinetic fields. In general, these need not be of the same
// layout as the independent variables.
const unsigned int n_dependent_variables =
// Define the helper that we will use in the AD computations for our
// scalar energy function. Note that we expect it to return values of
// type double.
VectorFunction<dim,double> ad_helper (n_independent_variables,
using ADNumberType = typename ADHelper::ad_type;
// Compute the fields that provide the independent values.
// When the tape is being replayed, these should be set to something
// meaningful.
const Tensor<1,dim> H = ...;
// If using a taped AD number, then at this point we would initiate
// taping of the expression for the material stored energy function
// for this particular set of material parameters:
// Select a tape number to record to
const typename Types<ADNumberType>::tape_index tape_index = ...;
// Indicate that we are about to start tracing the operations for
// function evaluation on the tape. If this tape has already been
// used (i.e. the operations are already recorded) then we
// (optionally) load the tape and reuse this data.
const bool is_recording
= ad_helper.start_recording_operations(tape_index);
// The steps that follow in the recording phase are required for
// tapeless methods as well.
if (is_recording == true)
{
// This is the "recording" phase of the operations.
// First, we set the values for all fields.
// These could happily be set to anything, unless the function will
// be evaluated along a branch not otherwise traversed during later
// use. For this reason, in this example instead of using some dummy
// values, we'll actually map out the function at the same point
// around which we'll later linearize it.
ad_helper.register_independent_variable(H, H_dofs);
ad_helper.register_independent_variable(C, C_dofs);
// NOTE: We have to extract the sensitivities in the order we wish to
// introduce them. So this means we have to do it by logical order
// of the extractors that we've created.
ad_helper.get_sensitive_variables(C_dofs);
ad_helper.get_sensitive_variables(H_dofs);
// Here we define the stress and magnetic induction in terms
// of the independent values C_AD and H_AD.
const Tensor<1, dim, ad_type> B_AD = ...;
// Register the definition of the kinetic fields. The second
// argument to the function provides a non-overlapping ordering
// of the
ad_helper.register_dependent_variable(S_AD, S_dofs);
ad_helper.register_dependent_variable(B_AD, B_dofs);
// Indicate that we have completed tracing the operations onto
// the tape.
ad_helper.stop_recording_operations(false); // write_tapes_to_file
}
else
{
// This is the "tape reuse" phase of the operations.
// Here we will leverage the already traced operations that reside
// on a tape, and simply re-evaluate the tape at a different point
// to get the function values and their derivatives.
// Load the existing tape to be reused
ad_helper.activate_recorded_tape(tape_index);
// Set the new values of the independent variables where the
// recorded dependent functions are to be evaluated (and
// differentiated around).
ad_helper.set_independent_variable(C, C_dofs);
ad_helper.set_independent_variable(H, H_dofs);
}
// Play the tape and store the output function value, its gradient and
// linearization. These are expensive to compute, so we'll do this once
// and extract the desired values from these intermediate outputs.
Vector<double> values (ad_helper.n_dependent_variables());
FullMatrix<double> jacobian (ad_helper.n_dependent_variables(),
ad_helper.n_independent_variables());
ad_helper.compute_values(values);
ad_helper.compute_jacobian(jacobian);
// Extract the desired components of the value vector and its Jacobian
// matrix. In this example, we use them to compute the Piola-Kirchhoff
// stress tensor S and its associated tangent, defined
// as HH = 2*dS/dC...
ad_helper.extract_value_component(values,S_dofs);
2.0*ad_helper.extract_jacobian_component(jacobian,S_dofs,C_dofs);
// ... the magnetic induction B and its associated tangent defined
// as BB = dB/dH...
const Tensor<1,dim> B =
ad_helper.extract_value_component(values,H_dofs);
symmetrize(ad_helper.extract_jacobian_component(jacobian,B_dofs,H_dofs));
// ... and finally the magnetoelastic coupling tangent, defined
// as PP = -dS/dH. Here the order of the extractor arguments is
// especially important, as it dictates the field that is being
// differentiated, and which the directional derivatives are being
// computed.
-ad_helper.extract_jacobian_component(jacobian,S_dofs,H_dofs)
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 3487 of file ad_helpers.h.

Member Typedef Documentation

◆ scalar_type

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType = double>
using Differentiation::AD::VectorFunction< 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 3496 of file ad_helpers.h.

◆ ad_type

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType = double>
using Differentiation::AD::VectorFunction< 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 3503 of file ad_helpers.h.

Constructor & Destructor Documentation

◆ VectorFunction()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType >
Differentiation::AD::VectorFunction< dim, ADNumberTypeCode, ScalarType >::VectorFunction ( 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 1655 of file ad_helpers.cc.

◆ ~VectorFunction()

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

Destructor.

Member Function Documentation

◆ register_dependent_variables()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType >
void Differentiation::AD::VectorFunction< dim, ADNumberTypeCode, ScalarType >::register_dependent_variables ( const std::vector< ad_type > &  funcs)

Register the definition of the vector field \(\boldsymbol{\Psi}(\mathbf{X})\).

Parameters
[in]funcsA vector of recorded functions that defines the dependent variables.
Note
For this class that expects only vector field of dependent variables, this function must only be called once per tape.
For taped AD numbers, this operation is only valid in recording mode.

Definition at line 1670 of file ad_helpers.cc.

◆ register_dependent_variable()

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

Register the definition of the vector field \(\hat{\mathbf{g}}(\mathbf{X}) \subset \boldsymbol{\Psi}(\mathbf{X})\) that may represent a subset of the dependent variables.

Parameters
[in]funcsThe recorded functions that define a set of dependent variables.
[in]extractorAn extractor associated with the input field variables. This effectively defines which components of the global set of dependent variables this field is associated with.
Note
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>.
For taped AD numbers, this operation is only valid in recording mode.

◆ compute_values()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType >
void Differentiation::AD::VectorFunction< dim, ADNumberTypeCode, ScalarType >::compute_values ( Vector< scalar_type > &  values) const

Compute the value of the vector field \(\boldsymbol{\Psi}(\mathbf{X})\).

Parameters
[out]valuesA Vector object with the value for each component of the vector field evaluated at the point defined by the independent variable values. The output values vector has a length corresponding to n_dependent_variables.

Definition at line 1686 of file ad_helpers.cc.

◆ compute_jacobian()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType >
void Differentiation::AD::VectorFunction< dim, ADNumberTypeCode, ScalarType >::compute_jacobian ( FullMatrix< scalar_type > &  jacobian) const

Compute the Jacobian (first derivative) of the vector field with respect to all independent variables, i.e.

\[ \mathbf{J}(\boldsymbol{\Psi}) = \frac{\partial\boldsymbol{\Psi}(\mathbf{X})}{\partial\mathbf{X}} \]

Parameters
[out]jacobianA FullMatrix with the gradient of each component of the vector field evaluated at the point defined by the independent variable values. The output jacobian matrix has dimensions corresponding to n_dependent_variables \(\times\)n_independent_variables.

Definition at line 1742 of file ad_helpers.cc.

◆ extract_value_component()

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType = double>
template<typename ExtractorType_Row >
static internal:: VectorFieldValue<dim, scalar_type, ExtractorType_Row>::type Differentiation::AD::VectorFunction< dim, ADNumberTypeCode, ScalarType >::extract_value_component ( const Vector< scalar_type > &  values,
const ExtractorType_Row &  extractor_row 
)
static

Extract the set of functions' values for a subset of dependent variables \(\mathbf{g} \subset \boldsymbol{\Psi}(\mathbf{X})\).

Parameters
[in]valuesA Vector object with the value for each component of the vector field evaluated at the point defined by the independent variable values.
[in]extractor_rowAn extractor associated with the input field variables. This effectively defines which components of the global set of dependent variables this field is associated with.

◆ extract_jacobian_component() [1/3]

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType = double>
template<typename ExtractorType_Row , typename ExtractorType_Col >
static internal::VectorFieldJacobian<dim, scalar_type, ExtractorType_Row, ExtractorType_Col>::type Differentiation::AD::VectorFunction< dim, ADNumberTypeCode, ScalarType >::extract_jacobian_component ( const FullMatrix< scalar_type > &  jacobian,
const ExtractorType_Row &  extractor_row,
const ExtractorType_Col &  extractor_col 
)
static

Extract the Jacobian of the subset of dependent functions \(\mathbf{g} \subset \boldsymbol{\Psi}(\mathbf{X})\) for a subset of independent variables \(\mathbf{A} \subset \mathbf{X}\), i.e.

\[ \mathbf{J}(\mathbf{g}) = \frac{\partial\mathbf{g}(\mathbf{X})}{\partial\mathbf{A}} \]

The first index of the Jacobian matrix \(\mathbf{J}(\mathbf{g})\) relates to the dependent variables, while the second index relates to the independent variables.

Parameters
[in]jacobianThe Jacobian of the vector function with respect to all independent variables, i.e., that returned by compute_jacobian().
[in]extractor_rowAn extractor associated with the input field variables for which the first index of the Jacobian is extracted. This effectively defines the correspondence between components of the global set of dependent variables and the field (representing a subset of dependent functions) associated with the extractor.
[in]extractor_colAn extractor associated with the input field variables for which the second index of the Jacobian is extracted. This effectively defines the correspondence between components of the global set of independent variables and the field (representing a subset of independent variables) associated with the extractor.
Returns
A Tensor or SymmetricTensor with its rank and symmetries determined by the extractor_row and extractor_col . This corresponds to subsetting a whole set of rows and columns of the Jacobian matrix, scaling those entries to take account of tensor symmetries, and then reshaping the (sub-)matrix so obtained into a tensor, the final result. For example, if extractor_row is a FEValuesExtractors::Vector and extractor_col is a FEValuesExtractors::Tensor, then the returned object is a Tensor of rank 3, with its first index associated with the field corresponding to the row extractor and the second and third indices associated with the field corresponding to the column extractor. Similarly, if extractor_row is a FEValuesExtractors::SymmetricTensor and extractor_col is a FEValuesExtractors::SymmetricTensor, then the returned object is a SymmetricTensor of rank 4, with its first two indices associated with the field corresponding to the row extractor and the last two indices associated with the field corresponding to the column extractor.

◆ extract_jacobian_component() [2/3]

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType >
Tensor< 0, dim, typename VectorFunction< dim, ADNumberTypeCode, ScalarType >::scalar_type > Differentiation::AD::VectorFunction< dim, ADNumberTypeCode, ScalarType >::extract_jacobian_component ( const FullMatrix< scalar_type > &  jacobian,
const FEValuesExtractors::Scalar extractor_row,
const FEValuesExtractors::Scalar extractor_col 
)
static

Extract the Jacobian of the subset of dependent functions \(\mathbf{g} \subset \boldsymbol{\Psi}(\mathbf{X})\) for a subset of independent variables \(\mathbf{A} \subset \mathbf{X}\), i.e.

\[ \mathbf{J}(\mathbf{g}) = \frac{\partial\mathbf{g}(\mathbf{X})}{\partial\mathbf{A}} \]

This function is a specialization of the above for rank-0 tensors (scalars). This corresponds to extracting a single entry of the Jacobian matrix because both extractors imply selection of just a single row or column of the matrix.

Definition at line 1821 of file ad_helpers.cc.

◆ extract_jacobian_component() [3/3]

template<int dim, enum AD::NumberTypes ADNumberTypeCode, typename ScalarType >
SymmetricTensor< 4, dim, typename VectorFunction< dim, ADNumberTypeCode, ScalarType >::scalar_type > Differentiation::AD::VectorFunction< dim, ADNumberTypeCode, ScalarType >::extract_jacobian_component ( const FullMatrix< scalar_type > &  jacobian,
const FEValuesExtractors::SymmetricTensor< 2 > &  extractor_row,
const FEValuesExtractors::SymmetricTensor< 2 > &  extractor_col 
)
static

Extract the Jacobian of the subset of dependent functions \(\mathbf{g} \subset \boldsymbol{\Psi}(\mathbf{X})\) for a subset of independent variables \(\mathbf{A} \subset \mathbf{X}\), i.e.

\[ \mathbf{J}(\mathbf{g}) = \frac{\partial\mathbf{g}(\mathbf{X})}{\partial\mathbf{A}} \]

This function is a specialization of the above for rank-4 symmetric tensors.

Definition at line 1859 of file ad_helpers.cc.


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