Adjoint API

Overview

This page describes the API for adjoint-based sensitivities in Kratos. It is the reference for developers creating new adjoint elements and objective functions.

Adjoint Variables

An adjoint variable is defined for each primal variable by prefixing the variable name with ADJOINT. For example, the adjoint variable of VELOCITY is ADJOINT_VELOCITY.

Adjoint Elements

Adjoint elements are created in the same way as regular elements. The element functions and their meaning for the adjoint problem are given below.

Solving the adjoint problem

The following functions are used in the solution of the adjoint problem:

GetDofList(DofsVectorType& rElementalDofList, ProcessInfo& rCurrentProcessInfo)

The list of adjoint dofs.

CalculateLeftHandSide(MatrixType& rLeftHandSideMatrix, ProcessInfo& rCurrentProcessInfo)

The adjoint matrix obtained from the gradient of the residual with respect to the primal variable (e.g., DISPLACEMENT). This is zero if the primal variable is a rate (e.g., VELOCITY).

CalculateFirstDerivativesLHS(MatrixType& rLeftHandSideMatrix, ProcessInfo& rCurrentProcessInfo)

The adjoint matrix obtained from the gradient of the residual with respect to first derivatives of the primal variable (e.g., VELOCITY).

CalculateSecondDerivativesLHS(MatrixType& rLeftHandSideMatrix, ProcessInfo& rCurrentProcessInfo)

The adjoint matrix obtained from the gradient of the residual with respect to second derivatives of the primal variable (e.g., ACCELERATION).

GetValuesVector(Vector& values, int Step = 0)

Get the vector of adjoint values for the primal variable (e.g., ADJOINT_DISPLACEMENT).

GetFirstDerivativesVector(Vector& values, int Step = 0)

Get the vector of adjoint values for first derivatives of primal variable (e.g., ADJOINT_VELOCITY).

GetSecondDerivativesVector(Vector& values, int Step = 0)

Get the vector of adjoint values for second derivatives of primal variable (e.g., ADJOINT_ACCELERATION).

Calculating sensitivities

The following functions are used after the solution of the adjoint problem to calculate sensitivities:

Calculate(const Variable<double>& rVariable, Matrix& Output, const ProcessInfo& rCurrentProcessInfo)
Calculate(const Variable<array_1d<double,3>>& rVariable, Matrix& Output, const ProcessInfo& rCurrentProcessInfo)

The transpose of the gradient of the element residual with respect to the design variable.

The sensitivities depend on the design variable which may be a double (e.g., THICKNESS) or array_1d<double,3> (e.g., NODAL_COORDINATES). The row dimension of the output is determined from the dimension of design variable (scalar or vector) and its type (nodal or elemental).

Adjoint Schemes

Calculating Sensitivities

The sensitivities are calculated in the update function.

Matrix pRpXt;
Vector adjoint, scalar_sensitivity, vector_sensitivity;
vector<Variable <double>*> scalar_design_variables;
vector<Variable < array_1d<double,3>>*> vector_design_variables;
...
for (auto it : model_part.Elements())
{
    // get the adjoint vector from the element (e.g., it->GetAdjointValuesVector(adjoint))
    ...
    for (auto p_current_variable : scalar_design_variables)
    {
        it->CalculateAdjointSensitivityContributions(pRpXt, *p_current_variable, pinfo);
        scalar_sensitivity = prod(A, adjoint);
        for (unsigned int i=0; i < it->GetGeometry().size(); ++i)
            it->GetGeometry()[i].FastGetSolutionStepValue(*p_current_variable) += scalar_sensitivity[i];
    }

    for (auto p_current_variable : vector_design_variables)
    {
        it->CalculateAdjointSensitivityContributions(pRpXt, *p_current_variable, pinfo);
        vector_sensitivity = prod(pRpXt, adjoint);
        unsigned k = 0;
        for (unsigned int i=0; i < it->GetGeometry().size(); ++i)
        {
            array_1d<double>& var = it->GetGeometry()[i].FastGetSolutionStepValue(*p_current_variable);
            for (unsigned d=0; d<3; ++d)
                var[d] += vector_sensitivity[k++];
        }
    }
}