Use connected lists for tape data structure

Instead of a single list use mutliple connected lists to store elements.
This allows to dynamically increase the size of tape without the need of
relocating elements.

benchmark/AlgorithmicComplexity.cpp

@@ -19,7 +19,7 @@ static void BM_NumericGausP(benchmark::State& state) {
   double p[] = {1, 2, 3, 4, 5};
   double dx[5] = {0, 0, 0, 0, 0};
   double dp[5] = {0, 0, 0, 0, 0};
-  clad::tape<clad::array_ref<double>> results = {};
+  clad::old_tape<clad::array_ref<double>> results = {};
   int dim = 5;
   results.emplace_back(dx, dim);
   results.emplace_back(dp, dim);

demos/CustomTypeNumDiff.cpp

@@ -134,7 +134,7 @@ int main() {
   // This is how we return the derivative with respect to all arguments.
   // The order of being placed in this tape should be the same as the order of
   // the arguments being passed to the function.
-  clad::tape<clad::array_ref<
+  clad::old_tape<clad::array_ref<
       double /*This should be the return value of the function you want to differentiate.*/>>
       grad = {};
   // Place the l-value reference of the variables in the tape.

include/clad/Differentiator/Differentiator.h

@@ -14,6 +14,7 @@
 #include "DynamicGraph.h"
 #include "FunctionTraits.h"
 #include "Matrix.h"
+#include "NewTape.h"
 #include "NumericalDiff.h"
 #include "Tape.h"
 
@@ -43,15 +44,16 @@ inline CUDA_HOST_DEVICE unsigned int GetLength(const char* code) {
 }
 
   /// Tape type used for storing values in reverse-mode AD inside loops.
-  template <typename T>
-  using tape = tape_impl<T>;
+template <typename T> using tape = new_tape_impl<T>;
 
-  /// Add value to the end of the tape, return the same value.
-  template <typename T, typename... ArgsT>
-  CUDA_HOST_DEVICE T push(tape<T>& to, ArgsT... val) {
-    to.emplace_back(std::forward<ArgsT>(val)...);
-    return to.back();
-  }
+template <typename T> using old_tape = tape_impl<T>;
+
+/// Add value to the end of the tape, return the same value.
+template <typename T, typename... ArgsT>
+CUDA_HOST_DEVICE T push(tape<T>& to, ArgsT... val) {
+  to.emplace_back(std::forward<ArgsT>(val)...);
+  return to.back();
+}
 
   /// Add value to the end of the tape, return the same value.
   /// A specialization for clad::array_ref types to use in reverse mode.

include/clad/Differentiator/NewTape.h

@@ -0,0 +1,87 @@
+#ifndef CLAD_NEWTAPE_H
+#define CLAD_NEWTAPE_H
+
+#include <cassert>
+#include <cstdio>
+#include <memory>
+#include <type_traits>
+#include <utility>
+
+namespace clad {
+
+template <typename T> class Block {
+public:
+  T* _data;
+  Block<T>* next;
+  Block<T>* prev;
+
+  Block(std::size_t& _capacity)
+      : next(nullptr), prev(nullptr),
+        _data(static_cast<T*>(
+            ::operator new(_capacity * sizeof(T), std::nothrow))) {}
+
+  ~Block() = default;
+  Block(const Block& other) = delete;
+  Block& operator=(const Block& other) = delete;
+
+  Block(Block&& other) = delete;
+  Block& operator=(const Block&& other) = delete;
+};
+
+template <typename T> class new_tape_impl {
+  Block<T>* m_cur_block = nullptr;
+  std::size_t _capacity = 32;
+  std::size_t m_size = 0;
+  using pointer = T*;
+  using reference = T&;
+
+  using iterator = pointer;
+
+public:
+  template <typename... ArgsT>
+
+  void emplace_back(ArgsT&&... args) {
+    if (!m_cur_block || m_size >= _capacity) {
+      Block<T>* prev_block = m_cur_block;
+      m_cur_block = new Block<T>(_capacity);
+      if (prev_block) {
+        prev_block->next = m_cur_block;
+        m_cur_block->prev = prev_block;
+      }
+      m_size = 0;
+    }
+    m_size += 1;
+    ::new (const_cast<void*>(static_cast<const volatile void*>(block_end())))
+        T(std::forward<ArgsT>(args)...);
+  }
+
+  [[nodiscard]] std::size_t size() const { return m_size; }
+
+  iterator block_begin() {
+    return reinterpret_cast<iterator>(m_cur_block->_data);
+  }
+
+  iterator block_end() {
+    return reinterpret_cast<iterator>(m_cur_block->_data + m_size - 1);
+  }
+
+  reference back() {
+    assert(m_size || m_cur_block->prev);
+    return block_begin()[m_size - 1];
+  }
+
+  void pop_back() {
+    assert(m_size || m_cur_block->prev);
+    block_end()->~T();
+    m_size -= 1;
+    if (m_size == 0) {
+      Block<T>* temp = m_cur_block;
+      m_cur_block = m_cur_block->prev;
+      temp->~Block<T>();
+      m_size = _capacity;
+    }
+  }
+};
+} // namespace clad
+
+#endif // CLAD_NEWTAPE_H

test/NumericalDiff/PureCentralDiffCalls.C

@@ -75,7 +75,7 @@ int main() { // expected-no-diagnostics
   printf("Result is = %f\n", func1_res); // CHECK-EXEC: Result is = 2.000000
 
   // Gradients, derivative wrt all args
-  clad::tape<clad::array_ref<double>> grad = {};
+  clad::old_tape<clad::array_ref<double>> grad = {};
   grad.emplace_back(dx, 3);
   grad.emplace_back(&dy);
   grad.emplace_back(&dz);