CoDiPack  2.2.0
A Code Differentiation Package
SciComp TU Kaiserslautern
Loading...
Searching...
No Matches
jacobianLinearTape.hpp
1/*
2 * CoDiPack, a Code Differentiation Package
3 *
4 * Copyright (C) 2015-2024 Chair for Scientific Computing (SciComp), University of Kaiserslautern-Landau
5 * Homepage: http://scicomp.rptu.de
6 * Contact: Prof. Nicolas R. Gauger (codi@scicomp.uni-kl.de)
7 *
8 * Lead developers: Max Sagebaum, Johannes Blühdorn (SciComp, University of Kaiserslautern-Landau)
9 *
10 * This file is part of CoDiPack (http://scicomp.rptu.de/software/codi).
11 *
12 * CoDiPack is free software: you can redistribute it and/or
13 * modify it under the terms of the GNU General Public License
14 * as published by the Free Software Foundation, either version 3 of the
15 * License, or (at your option) any later version.
16 *
17 * CoDiPack is distributed in the hope that it will be useful,
18 * but WITHOUT ANY WARRANTY; without even the implied warranty
19 * of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.
20 *
21 * See the GNU General Public License for more details.
22 * You should have received a copy of the GNU
23 * General Public License along with CoDiPack.
24 * If not, see <http://www.gnu.org/licenses/>.
25 *
26 * For other licensing options please contact us.
27 *
28 * Authors:
29 * - SciComp, University of Kaiserslautern-Landau:
30 * - Max Sagebaum
31 * - Johannes Blühdorn
32 * - Former members:
33 * - Tim Albring
34 */
35#pragma once
36
37#include <algorithm>
38#include <type_traits>
39
40#include "../config.h"
41#include "../expressions/lhsExpressionInterface.hpp"
42#include "../expressions/logic/compileTimeTraversalLogic.hpp"
43#include "../expressions/logic/traversalLogic.hpp"
44#include "../misc/macros.hpp"
45#include "../traits/expressionTraits.hpp"
46#include "data/chunk.hpp"
47#include "indices/linearIndexManager.hpp"
48#include "interfaces/reverseTapeInterface.hpp"
49#include "jacobianBaseTape.hpp"
50
52namespace codi {
53
61 template<typename T_TapeTypes>
62 struct JacobianLinearTape : public JacobianBaseTape<T_TapeTypes, JacobianLinearTape<T_TapeTypes>> {
63 public:
64
65 using TapeTypes = CODI_DD(T_TapeTypes,
68
70 friend Base;
71
72 using Real = typename TapeTypes::Real;
73 using Gradient = typename TapeTypes::Gradient;
74 using IndexManager = typename TapeTypes::IndexManager;
75 using Identifier = typename TapeTypes::Identifier;
76 using Position = typename Base::Position;
77
78 CODI_STATIC_ASSERT(IndexManager::IsLinear, "This class requires an index manager with a linear scheme.");
79
82
84
86 void clearAdjoints(Position const& start, Position const& end,
88 if (AdjointsManagement::Automatic == adjointsManagement) {
89 this->adjoints.beginUse();
90 }
91
92 using IndexPosition = CODI_DD(typename IndexManager::Position, int);
93 IndexPosition startIndex = this->llfByteData.template extractPosition<IndexPosition>(start);
94 IndexPosition endIndex = this->llfByteData.template extractPosition<IndexPosition>(end);
95
96 startIndex = std::min(startIndex, (IndexPosition)this->adjoints.size() - 1);
97 endIndex = std::min(endIndex, (IndexPosition)this->adjoints.size() - 1);
98
99 for (IndexPosition curPos = endIndex + 1; curPos <= startIndex; curPos += 1) {
100 this->adjoints[curPos] = Gradient();
101 }
102
103 if (AdjointsManagement::Automatic == adjointsManagement) {
104 this->adjoints.endUse();
105 }
106 }
107
108 protected:
109
112 CODI_INLINE void pushStmtData(Identifier const& index, Config::ArgumentSize const& numberOfArguments) {
113 CODI_UNUSED(index);
114
115 this->statementData.pushData(numberOfArguments);
116 }
117
119 template<typename Adjoint>
121 /* data from call */
122 JacobianLinearTape& tape, Adjoint* adjointVector,
123 /* data from low level function byte data vector */
124 size_t& curLLFByteDataPos, size_t const& endLLFByteDataPos, char* dataPtr,
125 /* data from low level function info data vector */
126 size_t& curLLFInfoDataPos, size_t const& endLLFInfoDataPos, Config::LowLevelFunctionToken* const tokenPtr,
127 Config::LowLevelFunctionDataSize* const dataSizePtr,
128 /* data from jacobian vector */
129 size_t& curJacobianPos, size_t const& endJacobianPos, Real const* const rhsJacobians,
130 Identifier const* const rhsIdentifiers,
131 /* data from statement vector */
132 size_t& curStmtPos, size_t const& endStmtPos, Config::ArgumentSize const* const numberOfJacobians,
133 /* data from index handler */
134 size_t const& startAdjointPos, size_t const& endAdjointPos) {
135 CODI_UNUSED(endJacobianPos, endStmtPos, endLLFByteDataPos, endLLFInfoDataPos);
136
137 typename Base::template VectorAccess<Adjoint> vectorAccess(adjointVector);
138
139 size_t curAdjointPos = startAdjointPos;
140
141 while (curAdjointPos < endAdjointPos) CODI_Likely {
142 curAdjointPos += 1;
143
144 Config::ArgumentSize const argsSize = numberOfJacobians[curStmtPos];
145
147 Base::template callLowLevelFunction<LowLevelFunctionEntryCallKind::Forward>(
148 tape, true, curLLFByteDataPos, dataPtr, curLLFInfoDataPos, tokenPtr, dataSizePtr, &vectorAccess);
149 } else if (Config::StatementInputTag == argsSize) CODI_Unlikely {
150 // Do nothing.
151 } else CODI_Likely {
152 Adjoint lhsAdjoint = Adjoint();
153 Base::incrementTangents(adjointVector, lhsAdjoint, argsSize, curJacobianPos, rhsJacobians, rhsIdentifiers);
154 adjointVector[curAdjointPos] = lhsAdjoint;
155
157 tape, curAdjointPos, GradientTraits::dim<Adjoint>(), GradientTraits::toArray(lhsAdjoint).data());
158 }
159
160 curStmtPos += 1;
161 }
162 }
163
165 template<typename Adjoint>
167 /* data from call */
168 JacobianLinearTape& tape, Adjoint* adjointVector,
169 /* data from low level function byte data vector */
170 size_t& curLLFByteDataPos, size_t const& endLLFByteDataPos, char* dataPtr,
171 /* data from low level function info data vector */
172 size_t& curLLFInfoDataPos, size_t const& endLLFInfoDataPos, Config::LowLevelFunctionToken* const tokenPtr,
173 Config::LowLevelFunctionDataSize* const dataSizePtr,
174 /* data from jacobianData */
175 size_t& curJacobianPos, size_t const& endJacobianPos, Real const* const rhsJacobians,
176 Identifier const* const rhsIdentifiers,
177 /* data from statementData */
178 size_t& curStmtPos, size_t const& endStmtPos, Config::ArgumentSize const* const numberOfJacobians,
179 /* data from index handler */
180 size_t const& startAdjointPos, size_t const& endAdjointPos) {
181 CODI_UNUSED(endJacobianPos, endStmtPos, endLLFByteDataPos, endLLFInfoDataPos);
182
183 typename Base::template VectorAccess<Adjoint> vectorAccess(adjointVector);
184
185 size_t curAdjointPos = startAdjointPos;
186
187 while (curAdjointPos > endAdjointPos) CODI_Likely {
188 curStmtPos -= 1;
189 Config::ArgumentSize const argsSize = numberOfJacobians[curStmtPos];
190
192 Base::template callLowLevelFunction<LowLevelFunctionEntryCallKind::Reverse>(
193 tape, false, curLLFByteDataPos, dataPtr, curLLFInfoDataPos, tokenPtr, dataSizePtr, &vectorAccess);
194 } else if (Config::StatementInputTag == argsSize) CODI_Unlikely {
195 // Do nothing.
196 } else CODI_Likely {
197 // No input value, perform regular statement evaluation.
198
199 Adjoint const lhsAdjoint = adjointVector[curAdjointPos]; // We do not use the zero index, decrement of
200 // curAdjointPos at the end of the loop.
201
203 tape, (Identifier)curAdjointPos, GradientTraits::dim<Adjoint>(), GradientTraits::toArray(lhsAdjoint).data());
204
206 adjointVector[curAdjointPos] = Adjoint();
207 }
208
209 Base::incrementAdjoints(adjointVector, lhsAdjoint, argsSize, curJacobianPos, rhsJacobians, rhsIdentifiers);
210 }
211
212 curAdjointPos -= 1;
213 }
214 }
215 };
216}
#define CODI_Unlikely
Declare unlikely evaluation of an execution path.
Definition config.h:399
#define CODI_INLINE
See codi::Config::ForcedInlines.
Definition config.h:457
#define CODI_Likely
Declare likely evaluation of an execution path.
Definition config.h:397
#define CODI_DD(Type, Default)
Abbreviation for CODI_DECLARE_DEFAULT.
Definition macros.hpp:94
#define CODI_STATIC_ASSERT(cond, message)
Static assert in CoDiPack.
Definition macros.hpp:123
#define CODI_T(...)
Abbreviation for CODI_TEMPLATE.
Definition macros.hpp:111
uint16_t LowLevelFunctionDataSize
Size store type for a low level function.
Definition config.h:98
uint16_t LowLevelFunctionToken
Token type for low level functions in the tapes.
Definition config.h:108
size_t constexpr StatementLowLevelFunctionTag
Statement tag for low level functions.
Definition config.h:126
bool constexpr ReversalZeroesAdjoints
With a linear index management, control if adjoints are set to zero during reversal.
Definition config.h:289
uint8_t ArgumentSize
Type for the number of arguments in statements.
Definition config.h:117
size_t constexpr StatementInputTag
Tag for statements that are inputs. Used in linear index management context.
Definition config.h:123
std::array< AtomicTraits::RemoveAtomic< typename TraitsImplementation< Gradient >::Real >, TraitsImplementation< Gradient >::dim > toArray(Gradient const &gradient)
Converts the (possibly multi-component) gradient to an array of Reals.
Definition gradientTraits.hpp:116
CoDiPack - Code Differentiation Package.
Definition codi.hpp:90
void CODI_UNUSED(Args const &...)
Disable unused warnings for an arbitrary number of arguments.
Definition macros.hpp:46
AdjointsManagement
Policies for management of the tape's interal adjoints.
Definition tapeParameters.hpp:98
@ Automatic
Manage internal adjoints automatically, including locking, bounds checking, and resizing.
Implementation of VectorAccessInterface for adjoint vectors.
Definition adjointVectorAccess.hpp:59
Data is stored chunk-wise in this DataInterface implementation. If a chunk runs out of space,...
Definition chunkedData.hpp:64
LowLevelFunctionByteData llfByteData
Byte data for low level functions.
Definition commonTapeImplementation.hpp:155
static void notifyStatementEvaluateListeners(Tape &tape, Identifier const &lhsIdentifier, size_t sizeLhsAdjoint, Real const *lhsAdjoint)
Invoke callbacks for StatementEvaluate events.
Definition eventSystem.hpp:712
Indices enable the mapping of primal values to their adjoint counterparts.
Definition indexManagerInterface.hpp:78
Base class for all standard Jacobian tape implementations.
Definition jacobianBaseTape.hpp:125
typename Base::Position Position
See TapeTypesInterface.
Definition jacobianBaseTape.hpp:149
StatementData statementData
Data stream for statement specific data.
Definition jacobianBaseTape.hpp:170
Adjoints adjoints
Evaluation vector for AD.
Definition jacobianBaseTape.hpp:173
void clearAdjoints(AdjointsManagement adjointsManagement=AdjointsManagement::Automatic)
Clear all adjoint values, that is, set them to zero.
Definition jacobianBaseTape.hpp:468
static void incrementAdjoints(Adjoint *adjointVector, Adjoint const &lhsAdjoint, Config::ArgumentSize const &numberOfArguments, size_t &curJacobianPos, Real const *const rhsJacobians, Identifier const *const rhsIdentifiers)
Performs the AD reverse equation for a statement.
Definition jacobianBaseTape.hpp:520
static void incrementTangents(Adjoint const *const adjointVector, Adjoint &lhsAdjoint, Config::ArgumentSize const &numberOfArguments, size_t &curJacobianPos, Real const *const rhsJacobians, Identifier const *const rhsIdentifiers)
Performs the AD forward equation for a statement.
Definition jacobianBaseTape.hpp:542
Final implementation for a Jacobian tape with a linear index management.
Definition jacobianLinearTape.hpp:62
typename TapeTypes::Real Real
See TapeTypesInterface.
Definition jacobianLinearTape.hpp:72
void clearAdjoints(Position const &start, Position const &end, AdjointsManagement adjointsManagement=AdjointsManagement::Automatic)
Clear all adjoints that would be set in a tape evaluation from start to end. It has to hold start >= ...
Definition jacobianLinearTape.hpp:86
typename TapeTypes::IndexManager IndexManager
See TapeTypesInterface.
Definition jacobianLinearTape.hpp:74
static void internalEvaluateForward_EvalStatements(JacobianLinearTape &tape, Adjoint *adjointVector, size_t &curLLFByteDataPos, size_t const &endLLFByteDataPos, char *dataPtr, size_t &curLLFInfoDataPos, size_t const &endLLFInfoDataPos, Config::LowLevelFunctionToken *const tokenPtr, Config::LowLevelFunctionDataSize *const dataSizePtr, size_t &curJacobianPos, size_t const &endJacobianPos, Real const *const rhsJacobians, Identifier const *const rhsIdentifiers, size_t &curStmtPos, size_t const &endStmtPos, Config::ArgumentSize const *const numberOfJacobians, size_t const &startAdjointPos, size_t const &endAdjointPos)
Perform a forward evaluation of the tape. Arguments are from the recursive eval methods of the DataIn...
Definition jacobianLinearTape.hpp:120
JacobianLinearTape()
Constructor.
Definition jacobianLinearTape.hpp:81
typename TapeTypes::Gradient Gradient
See TapeTypesInterface.
Definition jacobianLinearTape.hpp:73
T_TapeTypes TapeTypes
See JacobianLinearTape.
Definition jacobianLinearTape.hpp:65
typename TapeTypes::Identifier Identifier
See TapeTypesInterface.
Definition jacobianLinearTape.hpp:75
friend Base
Allow the base class to call protected and private methods.
Definition jacobianLinearTape.hpp:70
static void internalEvaluateReverse_EvalStatements(JacobianLinearTape &tape, Adjoint *adjointVector, size_t &curLLFByteDataPos, size_t const &endLLFByteDataPos, char *dataPtr, size_t &curLLFInfoDataPos, size_t const &endLLFInfoDataPos, Config::LowLevelFunctionToken *const tokenPtr, Config::LowLevelFunctionDataSize *const dataSizePtr, size_t &curJacobianPos, size_t const &endJacobianPos, Real const *const rhsJacobians, Identifier const *const rhsIdentifiers, size_t &curStmtPos, size_t const &endStmtPos, Config::ArgumentSize const *const numberOfJacobians, size_t const &startAdjointPos, size_t const &endAdjointPos)
Perform a reverse evaluation of the tape. Arguments are from the recursive eval methods of the DataIn...
Definition jacobianLinearTape.hpp:166
void pushStmtData(Identifier const &index, Config::ArgumentSize const &numberOfArguments)
Add statement specific data to the data streams.
Definition jacobianLinearTape.hpp:112
typename Base::Position Position
See TapeTypesInterface.
Definition jacobianLinearTape.hpp:76
Type definitions for the Jacobian tapes.
Definition jacobianBaseTape.hpp:75