Opm::AutoDiffBlock< Scalar > Class Template Reference

#include <AutoDiffBlock.hpp>

Inheritance diagram for Opm::AutoDiffBlock< Scalar >:
Inheritance graph

Public Types

typedef Eigen::Array< Scalar,
Eigen::Dynamic, 1 > 
V
 Underlying type for values. More...
 
typedef AutoDiffMatrix M
 Underlying type for jacobians. More...
 

Public Member Functions

AutoDiffBlockoperator+= (const AutoDiffBlock &rhs)
 Elementwise operator +=. More...
 
AutoDiffBlockoperator-= (const AutoDiffBlock &rhs)
 Elementwise operator -=. More...
 
AutoDiffBlock operator+ (const AutoDiffBlock &rhs) const
 Elementwise operator +. More...
 
AutoDiffBlock operator- (const AutoDiffBlock &rhs) const
 Elementwise operator -. More...
 
AutoDiffBlock operator* (const AutoDiffBlock &rhs) const
 Elementwise operator *. More...
 
AutoDiffBlock operator/ (const AutoDiffBlock &rhs) const
 Elementwise operator /. More...
 
template<class Ostream >
Ostream & print (Ostream &os) const
 I/O. More...
 
void swap (AutoDiffBlock &other)
 Efficient swap function. More...
 
int size () const
 Number of elements. More...
 
int numBlocks () const
 Number of Jacobian blocks. More...
 
std::vector< int > blockPattern () const
 Sizes (number of columns) of Jacobian blocks. More...
 
const Vvalue () const
 Function value. More...
 
const std::vector< M > & derivative () const
 Function derivatives. More...
 

Static Public Member Functions

static AutoDiffBlock null ()
 Construct an empty AutoDiffBlock. More...
 
static AutoDiffBlock constant (V &&val)
 
static AutoDiffBlock constant (const V &val)
 
static AutoDiffBlock constant (const V &val, const std::vector< int > &blocksizes)
 
static AutoDiffBlock variable (const int index, V &&val, const std::vector< int > &blocksizes)
 
static AutoDiffBlock variable (const int index, const V &val, const std::vector< int > &blocksizes)
 
static AutoDiffBlock function (V &&val, std::vector< M > &&jac)
 
static AutoDiffBlock function (const V &val, const std::vector< M > &jac)
 
static std::vector< AutoDiffBlockvariables (const std::vector< V > &initial_values)
 

Detailed Description

template<typename Scalar>
class Opm::AutoDiffBlock< Scalar >

A class for forward-mode automatic differentiation with vector values and sparse jacobian matrices.

The class contains a (column) vector of values and multiple sparse matrices representing its partial derivatives. Each such matrix has a number of rows equal to the number of rows in the value vector, and a number of columns equal to the number of discrete variables we want to compute the derivatives with respect to. The reason to have multiple such jacobians instead of just one is to allow simpler grouping of variables, making it easier to implement various preconditioning schemes. Only basic arithmetic operators are implemented for this class, reflecting our needs so far.

The class is built on the Eigen library, using an Eigen array type to contain the values and Eigen sparse matrices for the jacobians. The overloaded operators are intended to behave in a similar way to Eigen arrays, meaning for example that the * operator is elementwise multiplication. The only exception is multiplication with a sparse matrix from the left, which is treated as an Eigen matrix operation.

There are no public constructors, instead we use the Named Constructor pattern. In general, one needs to know which variables one wants to compute the derivatives with respect to before constructing an AutoDiffBlock. Some of the constructors require you to pass a block pattern. This should be a vector containing the number of columns you want for each jacobian matrix.

For example: you want the derivatives with respect to three different variables p, r and s. Assuming that there are 10 elements in p, and 20 in each of r and s, the block pattern is { 10, 20, 20 }. When creating the variables p, r and s in your program you have two options:

  • Use the variable() constructor three times, passing the index (0 for p, 1 for r and 2 for s), initial value of each variable and the block pattern.
  • Use the variables() constructor passing only the initial values of each variable. The block pattern will be inferred from the size of the initial value vectors. This is usually the simplest option if you have multiple variables. Note that this constructor returns a vector of all three variables, so you need to use index access (operator[]) to get the individual variables (that is p, r and s).

After this, the r variable for example will have a size() of 20 and three jacobian matrices. The first is a 20 by 10 zero matrix, the second is a 20 by 20 identity matrix, and the third is a 20 by 20 zero matrix.

Member Typedef Documentation

template<typename Scalar>
typedef AutoDiffMatrix Opm::AutoDiffBlock< Scalar >::M

Underlying type for jacobians.

template<typename Scalar>
typedef Eigen::Array<Scalar, Eigen::Dynamic, 1> Opm::AutoDiffBlock< Scalar >::V

Underlying type for values.

Member Function Documentation

template<typename Scalar>
std::vector<int> Opm::AutoDiffBlock< Scalar >::blockPattern ( ) const
inline
template<typename Scalar>
static AutoDiffBlock Opm::AutoDiffBlock< Scalar >::constant ( const V val)
inlinestatic

Create an AutoDiffBlock representing a constant.

Parameters
[in]valvalues
template<typename Scalar>
static AutoDiffBlock Opm::AutoDiffBlock< Scalar >::constant ( const V val,
const std::vector< int > &  blocksizes 
)
inlinestatic

Create an AutoDiffBlock representing a constant. This variant requires specifying the block sizes used for the Jacobians even though the Jacobian matrices themselves will be zero.

Parameters
[in]valvalues
[in]blocksizesblock pattern
template<typename Scalar>
static AutoDiffBlock Opm::AutoDiffBlock< Scalar >::function ( V &&  val,
std::vector< M > &&  jac 
)
inlinestatic

Create an AutoDiffBlock by directly specifying values and jacobians. This version of function() moves its arguments and is therefore quite efficient, but leaves the argument variables empty (but valid).

Parameters
[in]valvalues
[in]jacvector of jacobians

Referenced by Opm::operator*().

template<typename Scalar>
static AutoDiffBlock Opm::AutoDiffBlock< Scalar >::function ( const V val,
const std::vector< M > &  jac 
)
inlinestatic

Create an AutoDiffBlock by directly specifying values and jacobians. This version of function() copies its arguments and is therefore less efficient than the other (moving) overload.

Parameters
[in]valvalues
[in]jacvector of jacobians
template<typename Scalar>
static AutoDiffBlock Opm::AutoDiffBlock< Scalar >::null ( )
inlinestatic

Construct an empty AutoDiffBlock.

template<typename Scalar>
AutoDiffBlock Opm::AutoDiffBlock< Scalar >::operator* ( const AutoDiffBlock< Scalar > &  rhs) const
inline

Elementwise operator *.

template<typename Scalar>
AutoDiffBlock Opm::AutoDiffBlock< Scalar >::operator+ ( const AutoDiffBlock< Scalar > &  rhs) const
inline

Elementwise operator +.

template<typename Scalar>
AutoDiffBlock& Opm::AutoDiffBlock< Scalar >::operator+= ( const AutoDiffBlock< Scalar > &  rhs)
inline

Elementwise operator +=.

template<typename Scalar>
AutoDiffBlock Opm::AutoDiffBlock< Scalar >::operator- ( const AutoDiffBlock< Scalar > &  rhs) const
inline

Elementwise operator -.

template<typename Scalar>
AutoDiffBlock& Opm::AutoDiffBlock< Scalar >::operator-= ( const AutoDiffBlock< Scalar > &  rhs)
inline

Elementwise operator -=.

template<typename Scalar>
AutoDiffBlock Opm::AutoDiffBlock< Scalar >::operator/ ( const AutoDiffBlock< Scalar > &  rhs) const
inline

Elementwise operator /.

template<typename Scalar>
template<class Ostream >
Ostream& Opm::AutoDiffBlock< Scalar >::print ( Ostream &  os) const
inline

I/O.

template<typename Scalar>
void Opm::AutoDiffBlock< Scalar >::swap ( AutoDiffBlock< Scalar > &  other)
inline

Efficient swap function.

template<typename Scalar>
const V& Opm::AutoDiffBlock< Scalar >::value ( ) const
inline

Function value.

Referenced by Opm::BlackoilSolventModel< Grid >::addWellContributionToMassBalanceEq(), Opm::BlackoilModelBase< Grid, Implementation >::applyThresholdPressures(), Opm::BlackoilSolventModel< Grid >::assemble(), Opm::collapseJacs(), Opm::BlackoilSolventModel< Grid >::computeMassFlux(), Opm::BlackoilModelBase< Grid, Implementation >::computeMassFlux(), Opm::BlackoilSolventModel< Grid >::computeWellConnectionPressures(), Opm::BlackoilModelBase< Grid, Implementation >::computeWellConnectionPressures(), Opm::BlackoilModelBase< Grid, Implementation >::computeWellFlux(), Opm::BlackoilModelBase< Grid, Implementation >::convergenceReduction(), Opm::BlackoilModelBase< Grid, Implementation >::fluidRsSat(), Opm::BlackoilModelBase< Grid, Implementation >::fluidRvSat(), Opm::BlackoilModelBase< Grid, Implementation >::getConvergence(), Opm::BlackoilModelBase< Grid, Implementation >::getWellConvergence(), Opm::detail::infinityNorm(), Opm::detail::infinityNormWell(), Opm::BlackoilModelBase< Grid, Implementation >::makeConstantState(), Opm::operator*(), Opm::AutoDiffBlock< double >::operator+=(), Opm::AutoDiffBlock< double >::operator-=(), Opm::BlackoilModelBase< Grid, Implementation >::poroMult(), Opm::BlackoilModelBase< Grid, Implementation >::solveWellEq(), Opm::subset(), Opm::BlackoilModelBase< Grid, Implementation >::transMult(), Opm::BlackoilSolventModel< Grid >::updateEquationsScaling(), Opm::BlackoilModelBase< Grid, Implementation >::updateEquationsScaling(), Opm::AutoDiffBlock< double >::variable(), and Opm::detail::zeroIfNan().

template<typename Scalar>
static AutoDiffBlock Opm::AutoDiffBlock< Scalar >::variable ( const int  index,
V &&  val,
const std::vector< int > &  blocksizes 
)
inlinestatic

Create an AutoDiffBlock representing a single variable block.

Parameters
[in]indexindex of the variable you are constructing
[in]valvalues
[in]blocksizesblock pattern The resulting object will have size() equal to block_pattern[index]. Its jacobians will all be zero, except for derivative()[index], which will be an identity matrix.

Referenced by Opm::AutoDiffBlock< double >::variable(), and Opm::AutoDiffBlock< double >::variables().

template<typename Scalar>
static AutoDiffBlock Opm::AutoDiffBlock< Scalar >::variable ( const int  index,
const V val,
const std::vector< int > &  blocksizes 
)
inlinestatic

Create an AutoDiffBlock representing a single variable block.

Parameters
[in]indexindex of the variable you are constructing
[in]valvalues
[in]blocksizesblock pattern The resulting object will have size() equal to block_pattern[index]. Its jacobians will all be zero, except for derivative()[index], which will be an identity matrix.
template<typename Scalar>
static std::vector<AutoDiffBlock> Opm::AutoDiffBlock< Scalar >::variables ( const std::vector< V > &  initial_values)
inlinestatic

Construct a set of primary variables, each initialized to a given vector.


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