Merge branch 'develop'

Author: MaxSagebaum

Makefile

@@ -146,6 +146,9 @@ doc:
 	@mkdir -p $(BUILD_DIR)/documentation
 	CODI_VERSION=$(CODI_VERSION) doxygen
 
+single_header:
+	quom --include_directory include include/codi.hpp $(BUILD_DIR)/codi_single.hpp
+
 .PHONY: format
 format:
 	find include tests/general/include tests/general/src tests/events/include tests/events/src -type f -exec $(CLANG_FORMAT) -i {} \;

documentation/Changelog.md

@@ -1,6 +1,22 @@
 Changelog {#Changelog}
 ===========================
 
+### v3.?.? - ????-??-??
+ - Features:
+   * It is now possible to define custom tape evaluators for all CoDiPack tapes. The evaluators have access to the full
+     statement data and functionality for low level functions.
+   * Low level function and external functions can now iterate over their input and output identifiers.
+   * New tool for optimizing the cache access of reuse index tapes. See \ref Example_29_Tape_cache_optimization.
+
+ - Internal:
+   * Restructure of per value tape data handling. Each tape and each index manager in CoDiPack can now define data that
+     is stored in each value. This is a breaking interface change but it will not affect the default CoDiPack tapes.
+     This change is mostly used for debugging.
+
+ - Bugfix:
+   * Explicitly set the language of CoDiPack to C++ for CMake.
+   * Add missing setter functions to complex numbers.
+
 ### v3.0.0 - 2025-07-08
  - General:
    * Raised default cpp version of CoDiPack to 17. If you require a lower cpp version please use CoDiPack 2.*.

documentation/CoDiLayout.xml

@@ -46,6 +46,7 @@
       <tab type="user" title="Example 26 - Jacobian tape readers" url="@ref Example_26_Jacobian_Tape_Readers"/>
       <tab type="user" title="Example 27 - Primal tape readers" url="@ref Example_27_Primal_Tape_Readers"/>
       <tab type="user" title="Example 28 - Complex numbers" url="@ref Example_28_Complex_numbers"/>
+      <tab type="user" title="Example 29 - Tape cache optimization" url="@ref Example_29_Tape_cache_optimization"/>
     </tab>
     <tab type="user" title="Papers" url="@ref Papers"/>
     <tab type="user" title="Taping strategies" url="@ref TapingStrategy"/>

documentation/developer/simpleTape.cpp

@@ -46,6 +46,7 @@ struct SimpleTape : public codi::ReverseTapeInterface<double, double, int> {
     using Real = double;
     using Gradient = double;
     using Identifier = int;
+    using ActiveTypeTapeData = int;
 
 //! [Data stream - Type definition]
     using OperatorData = codi::ChunkedData<codi::Chunk1<OperatorCode>>;
@@ -176,17 +177,25 @@ struct SimpleTape : public codi::ReverseTapeInterface<double, double, int> {
 
     static bool constexpr AllowJacobianOptimization = false; // If certain operations can be hidden from the tape.
 
-//! [Identifiers - Initialization]
+//! [Identifiers - Initialization and handling]
     template<typename Real>
-    void initIdentifier(Real& value, Identifier& identifier) {
-      identifier = 0; // Initialize with zero we perform an online activity analysis.
+    void initTapeData(Real& value, ActiveTypeTapeData& data) {
+      data = 0; // Initialize with zero we perform an online activity analysis.
     }
 
     template<typename Real>
-    void destroyIdentifier(Real& value, Identifier& identifier) {
+    void destroyTapeData(Real& value, ActiveTypeTapeData& data) {
       // Do nothing: Identifiers are not reused.
     }
-//! [Identifiers - Initialization]
+
+    Identifier const& getIdentifier(ActiveTypeTapeData const& data) {
+      return data;
+    }
+
+    Identifier& getIdentifier(ActiveTypeTapeData& data) {
+      return data;
+    }
+//! [Identifiers - Initialization and handling]
 
 //! [Storing - Entry]
     template<typename Lhs, typename Rhs>

documentation/developer/simpleTape.md

@@ -97,11 +97,13 @@ left hand side identifier from the data stream in a linear index management sche
 \snippet developer/simpleTape.cpp Identifiers - Registration
 
 The identifiers are stored in the AD type provided by CoDiPack. The initialization of the identifier in
-the AD value is done by the function `initIdentifier` required by the codi::InternalStatementRecordingTapeInterface. We
+the AD value is done by the function `initTapeData` required by the codi::InternalStatementRecordingTapeInterface. We
 implement an online activity analysis in this tape. Therefore, all identifiers in the AD values can be initialized with
 zero. The zero identifier is used in our implementation to track _passive_ values. These are values that do not depend
-on the input values. How this is done is explained in the next section.
-\snippet developer/simpleTape.cpp Identifiers - Initialization
+on the input values. How this is done is explained in the next section. In addition the
+codi::IdentifierInformationTapeInterface requires the function 'getIdentifier' for a const and non-const argument. Since
+we do not have any tape specific data beside the identifier, these function are the identity.
+\snippet developer/simpleTape.cpp Identifiers - Initialization and handling
 
 #### Storing of expressions/operators
 

documentation/examples/Example_23_OpenMP_Parallel_Codes.cpp

@@ -1,7 +1,6 @@
 #if CODI_EnableOpDiLib
 //! [Example 23 - OpenMP Parallel Codes]
 #include <codi.hpp>
-#include <opdi.hpp>
 #include <iostream>
 
 #include <opdi/backend/macro/macroBackend.hpp>
@@ -19,6 +18,7 @@ int main(int nargs, char** args) {
   opdi::logic = new opdi::OmpLogic;
   opdi::logic->init();
   opdi::tool = new CoDiOpDiLibTool<Real>;
+  opdi::tool->init();
 
   // usual AD workflow in the serial parts of the code
 
@@ -55,14 +55,16 @@ int main(int nargs, char** args) {
   y.setGradient(1.0);
 
   opdi::logic->prepareEvaluate();  // prepare OpDiLib for evaluation
-
   tape.evaluate();
+  opdi::logic->postEvaluate();  // OpDiLib-specific postprocessing
 
   std::cout << "f(" << x << ") = " << y << std::endl;
   std::cout << "df/dx(" << x << ") = " << x.getGradient() << std::endl;
 
   // finalize OpDiLib
 
+  opdi::tool->finalize();
+  opdi::logic->finalize();
   opdi::backend->finalize();
   delete opdi::backend;
   delete opdi::logic;

documentation/examples/Example_29_Tape_Cache_Optimization.cpp

@@ -0,0 +1,92 @@
+//! [Example 29 - Tape cache optimization]
+#include <codi.hpp>
+#include <iostream>
+
+//! [Function]
+template<typename Real>
+void func(const Real* x, size_t l, Real* y) {
+  y[0] = 0.0;
+  y[1] = 1.0;
+  for(size_t i = 0; i < l; ++i) {
+    y[0] += x[i];
+    y[1] *= x[i];
+  }
+}
+//! [Function]
+
+int main(int nargs, char** args) {
+
+  using Real = codi::RealReverseIndex;
+  using Identifier = typename Real::Identifier;
+  using Tape = typename Real::Tape;
+
+  Real x[5];
+  Real y[2];
+  x[0] = 1.0;
+  x[1] = 2.0;
+  x[2] = 3.0;
+  x[3] = 4.0;
+  x[4] = 5.0;
+
+  // Step 1: Record the tape.
+  Tape& tape = Real::getTape();
+  tape.setActive();
+
+  for(size_t i = 0; i < 5; ++i) {
+    tape.registerInput(x[i]);
+  }
+
+  func(x, 5, y);
+
+  tape.registerOutput(y[0]);
+  tape.registerOutput(y[1]);
+
+  tape.setPassive();
+
+  // Step 2: Gather the input and output identifiers.
+  Identifier xIds[5];
+  Identifier yIds[2];
+  for(int i = 0; i < 5; i += 1) {
+    xIds[i] = x[i].getIdentifier();
+  }
+  for(int i = 0; i < 2; i += 1) {
+    yIds[i] = y[i].getIdentifier();
+  }
+
+  // Step 3: Define the input and output iterators.
+  auto iterX = [&xIds](auto&& func) {
+    for(size_t i = 0; i < 5; ++i) {
+      func(xIds[i]);
+    }
+  };
+  auto iterY = [&yIds](auto&& func) {
+    for(size_t i = 0; i < 2; ++i) {
+      func(yIds[i]);
+    }
+  };
+
+  // Step 4: Apply the optimization.
+  codi::IdentifierCacheOptimizerHotCold<Tape> co{tape};
+  co.eval(iterX, iterY);
+
+  // Step 5: Do a tape evaluation with the translated ids.
+  codi::Jacobian<double> jacobian(2,5);
+  for(size_t curY = 0; curY < 2; curY += 1) {
+    tape.gradient(yIds[curY]) = 1.0;
+    tape.evaluate();
+
+    for(size_t curX = 0; curX < 5; curX += 1) {
+      jacobian(curY,curX) = tape.gradient(xIds[curX]);
+      tape.gradient(xIds[curX]) = 0.0;
+    }
+  }
+
+  std::cout << "Reverse Jacobian:" << std::endl;
+  std::cout << "f(1 .. 5) = (" << y[0] << ", " << y[1] << ")" << std::endl;
+  std::cout << "df/dx (1 .. 5) = \n" << jacobian << std::endl;
+
+  tape.reset();
+
+  return 0;
+}
+//! [Example 29 - Tape cache optimization]

documentation/examples/Example_29_Tape_Cache_Optimization.md

@@ -0,0 +1,22 @@
+Example 29 - Tape cache optimization {#Example_29_Tape_cache_optimization}
+=======
+
+**Goal:** Demonstrate the use of std::complex with CoDiPack.
+
+**Prerequisite:** \ref Tutorial_02_Reverse_mode_AD, \ref Example_02_Custom_adjoint_vector_evaluation
+
+**Function:**
+\snippet examples/Example_28_Complex_numbers.cpp Function implementations
+
+**Full code:**
+\snippet examples/Example_29_Tape_Cache_Optimization.cpp Example 29 - Tape cache optimization
+
+**Additional information:**
+The cache optimizer performs a lifetime analysis of the identifiers on the tape. The identifiers are redistributed
+for more cache performance during the reverse evaluation.
+
+Since identifiers are redistributed, identifiers from values like `x` should not be used. Instead, the identifier of
+`x` should be stored and then the stored value should be given to the optimizer.
+
+The cache optimization is only meaningfully if the tape is evaluated at least 10 times, e.g., as in a reverse
+accumulation process.

documentation/user/Tutorials.md

@@ -48,6 +48,7 @@ or pointers to other features.
 | \subpage Example_26_Jacobian_Tape_Readers "" | Reading Jacobian tapes from disk. |
 | \subpage Example_27_Primal_Tape_Readers "" | Rading primal value tapes from disk. |
 | \subpage Example_28_Complex_numbers "" | How to use complex numbers in CoDiPack. |
+| \subpage Example_29_Tape_cache_optimization "" | Applying a cache optimimization for faster reverse evaluations to the tape.|
 
 The graph shows how the tutorials and examples are connected. Usually it is better to understand first the prerequisites
 of a tutorial/example before reading the actual example.

documentation/user/TutorialsGraph.dot

@@ -76,6 +76,8 @@ digraph Tutorials {
 
   E28 [label="E28 - Complex numbers"];
 
+  E29 [label="E29 - Tape cache optimization"];
+
   // Edges (sorted)
   E02:e -> E08:w;
   E02:e -> E09:w;
@@ -102,6 +104,7 @@ digraph Tutorials {
   T02:e -> E23:w;
   T02:e -> E25:w;
   T02:e -> E28:w;
+  T02:e -> E29:w;
   T02:e -> T03:w;
   T02:e -> T04:w;
   T02:e -> T05:w;