Merge branch '19-flexible-stack' into 'master'

Resolve "Flexible Stack"

Closes #19

See merge request !73

CMakeModules/CorsikaUtilities.cmake

@@ -86,6 +86,7 @@ endmacro(CORSIKA_ADD_FILES_ABSOLUTE)
 #
 
 function (CORSIKA_ADD_TEST name)
+  target_include_directories (${name} PRIVATE ${CMAKE_CURRENT_SOURCE_DIR})
   file (MAKE_DIRECTORY ${PROJECT_BINARY_DIR}/test_outputs/)
   add_test (NAME ${name} COMMAND ${name} -o ${PROJECT_BINARY_DIR}/test_outputs/junit-${name}.xml -r junit)
 endfunction (CORSIKA_ADD_TEST)

COAST/ParticleConversion.h

@@ -191,19 +191,19 @@ namespace corsika::coast {
       {CoastCode::Tritium, corsika::particles::Code::Tritium},
       {CoastCode::He3, corsika::particles::Code::Helium3},
       {CoastCode::Helium, corsika::particles::Code::Helium}, // 402
-      {CoastCode::Lithium, corsika::particles::Code::Lithium},
-      {CoastCode::Beryllium, corsika::particles::Code::Beryllium},
-      {CoastCode::Boron, corsika::particles::Code::Boron},
-      {CoastCode::Carbon, corsika::particles::Code::Carbon},
-      {CoastCode::Carbon13, corsika::particles::Code::Carbon13},
-      {CoastCode::Nitrogen, corsika::particles::Code::Nitrogen},
-      {CoastCode::Fluor, corsika::particles::Code::Fluor},
-      {CoastCode::Oxygen, corsika::particles::Code::Oxygen},
-      {CoastCode::Neon21, corsika::particles::Code::Neon21},
-      {CoastCode::Neon, corsika::particles::Code::Neon},
-      {CoastCode::Argon, corsika::particles::Code::Argon},
-      {CoastCode::Xenon, corsika::particles::Code::Xenon},
-      {CoastCode::Radon, corsika::particles::Code::Radon},
+      //{CoastCode::Lithium, corsika::particles::Code::Lithium},
+      //{CoastCode::Beryllium, corsika::particles::Code::Beryllium},
+      //{CoastCode::Boron, corsika::particles::Code::Boron},
+      //{CoastCode::Carbon, corsika::particles::Code::Carbon},
+      //{CoastCode::Carbon13, corsika::particles::Code::Carbon13},
+      //{CoastCode::Nitrogen, corsika::particles::Code::Nitrogen},
+      //{CoastCode::Fluor, corsika::particles::Code::Fluor},
+      //{CoastCode::Oxygen, corsika::particles::Code::Oxygen},
+      //{CoastCode::Neon21, corsika::particles::Code::Neon21},
+      //{CoastCode::Neon, corsika::particles::Code::Neon},
+      //{CoastCode::Argon, corsika::particles::Code::Argon},
+      //{CoastCode::Xenon, corsika::particles::Code::Xenon},
+      //{CoastCode::Radon, corsika::particles::Code::Radon},
       {CoastCode::Iron, corsika::particles::Code::Iron}, // 5628
 
       //{CoastCode::, corsika::particles::Code::},

Framework/Cascade/Cascade.h

@@ -29,6 +29,8 @@
 namespace corsika::cascade {
 
   /**
+   * \class Cascade
+   *
    * The Cascade class is constructed from template arguments making
    * it very versatile. Via the template arguments physics models are
    * plugged into the cascade simulation.
@@ -53,6 +55,7 @@ namespace corsika::cascade {
   class Cascade {
     using Particle = typename Stack::ParticleType;
 
+    // we only want fully configured objects
     Cascade() = delete;
 
   public:

Framework/Particles/ParticleProperties.cc

@@ -30,9 +30,10 @@ namespace corsika::particles {
   
   Code ConvertFromPDG(PDGCode p) {
       static_assert(detail::conversionArray.size() % 2 == 1);
-      auto constexpr maxPDG{(detail::conversionArray.size() - 1) >> 1};
+      // this will fail, for the strange case where the maxPDG is negative...
+      unsigned int constexpr maxPDG{(detail::conversionArray.size() - 1) >> 1};
       auto k = static_cast<PDGCodeType>(p);
-      if (abs(k) <= maxPDG) {
+      if ((unsigned int)abs(k) <= maxPDG) {
           return detail::conversionArray[k + maxPDG];
       } else {
           return detail::conversionMap.at(p);

Framework/StackInterface/CMakeLists.txt

 set (
   CORSIKAstackinterface_HEADERS
+  ParticleBase.h
+  StackIteratorInterface.h
   Stack.h
   SecondaryView.h
-  StackIteratorInterface.h
-  ParticleBase.h
+  CombinedStack.h
   )
 
 set (
@@ -36,3 +37,11 @@ install (
 add_executable (testStackInterface testStackInterface.cc)
 target_link_libraries (testStackInterface CORSIKAstackinterface CORSIKAthirdparty) # for catch2
 CORSIKA_ADD_TEST(testStackInterface)
+
+add_executable (testSecondaryView testSecondaryView.cc)
+target_link_libraries (testSecondaryView CORSIKAstackinterface CORSIKAthirdparty) # for catch2
+CORSIKA_ADD_TEST(testSecondaryView)
+
+add_executable (testCombinedStack testCombinedStack.cc)
+target_link_libraries (testCombinedStack CORSIKAstackinterface CORSIKAthirdparty) # for catch2
+CORSIKA_ADD_TEST(testCombinedStack)

Framework/StackInterface/CombinedStack.h

+
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
+#ifndef _include_stack_combinedstack_h_
+#define _include_stack_combinedstack_h_
+
+#include <corsika/particles/ParticleProperties.h>
+#include <corsika/stack/Stack.h>
+#include <corsika/units/PhysicalUnits.h>
+
+namespace corsika::stack {
+
+  /**
+   *
+   *
+   */
+  template <template <typename> typename ParticleInterface,
+            template <typename> typename ParticleInterfaceAdd, typename StackIterator>
+  class CombinedParticleInterface
+      : public ParticleInterfaceAdd<ParticleInterface<StackIterator>> {
+
+    using C =
+        CombinedParticleInterface<ParticleInterface, ParticleInterfaceAdd, StackIterator>;
+    using T = ParticleInterfaceAdd<ParticleInterface<StackIterator>>;
+    using I = ParticleInterface<StackIterator>;
+
+  protected:
+    using T::GetIndex;
+    using T::GetStackData;
+
+  public:
+    template <typename... Args1>
+    void SetParticleData(const std::tuple<Args1...> vA) {
+      I::SetParticleData(vA);
+      T::SetParticleData();
+    }
+    template <typename... Args1, typename... Args2>
+    void SetParticleData(const std::tuple<Args1...> vA, const std::tuple<Args2...> vB) {
+      I::SetParticleData(vA);
+      T::SetParticleData(vB);
+    }
+
+    template <typename... Args1>
+    void SetParticleData(C& p, const std::tuple<Args1...> vA) {
+      // static_assert(MT<I>::has_not, "error");
+      I::SetParticleData(static_cast<I&>(p), vA);
+      T::SetParticleData(static_cast<T&>(p));
+    }
+    template <typename... Args1, typename... Args2>
+    void SetParticleData(C& p, const std::tuple<Args1...> vA,
+                         const std::tuple<Args2...> vB) {
+      I::SetParticleData(static_cast<I&>(p), vA);
+      T::SetParticleData(static_cast<T&>(p), vB);
+    }
+  };
+
+  /**
+   * Memory implementation of the most simple (stupid) particle stack object.
+   */
+  template <typename Stack1Impl, typename Stack2Impl>
+  class CombinedStackImpl : public Stack1Impl, public Stack2Impl {
+
+  public:
+    void Init() {
+      Stack1Impl::Init();
+      Stack2Impl::Init();
+    }
+
+    void Clear() {
+      Stack1Impl::Clear();
+      Stack2Impl::Clear();
+    }
+
+    unsigned int GetSize() const { return Stack1Impl::GetSize(); }
+    unsigned int GetCapacity() const { return Stack1Impl::GetCapacity(); }
+
+    /**
+     *   Function to copy particle at location i1 in stack to i2
+     */
+    void Copy(const unsigned int i1, const unsigned int i2) {
+      if (i1 >= GetSize() || i2 >= GetSize()) {
+        std::ostringstream err;
+        err << "CombinedStack: trying to access data beyond size of stack!";
+        throw std::runtime_error(err.str());
+      }
+      Stack1Impl::Copy(i1, i2);
+      Stack2Impl::Copy(i1, i2);
+    }
+
+    /**
+     *   Function to copy particle at location i2 in stack to i1
+     */
+    void Swap(const unsigned int i1, const unsigned int i2) {
+      if (i1 >= GetSize() || i2 >= GetSize()) {
+        std::ostringstream err;
+        err << "CombinedStack: trying to access data beyond size of stack!";
+        throw std::runtime_error(err.str());
+      }
+      Stack1Impl::Swap(i1, i2);
+      Stack2Impl::Swap(i1, i2);
+    }
+
+    void IncrementSize() {
+      Stack1Impl::IncrementSize();
+      Stack2Impl::IncrementSize();
+    }
+
+    void DecrementSize() {
+      Stack1Impl::DecrementSize();
+      Stack2Impl::DecrementSize();
+    }
+
+  private:
+    /// the actual memory to store particle data
+
+  }; // end class CombinedStackImpl
+
+  template <typename Stack1Impl, typename Stack2Impl, template <typename> typename _PI>
+  using CombinedStack = Stack<CombinedStackImpl<Stack1Impl, Stack2Impl>, _PI>;
+
+} // namespace corsika::stack
+
+#endif

Framework/StackInterface/ParticleBase.h

@@ -87,9 +87,9 @@ namespace corsika::stack {
       return static_cast<const StackIterator&>(*this);
     }
 
-  protected:
+    //  protected:
     /**
-        @name Access to underlying stack data, these are service
+        @name Access to underlying stack dfata, these are service
         function for user classes. User code can only rely on GetIndex
         and GetStackData to retrieve data
         @{
@@ -106,6 +106,25 @@ namespace corsika::stack {
     ///@}
   };
 
+  template <typename T>
+  class ParticleBaseAdd {
+
+  public:
+    ParticleBaseAdd() = default;
+
+    using T::GetIndex;
+    using T::GetStackData;
+
+  public:
+    template <typename... Args1, typename... Args2>
+    void SetParticleData(Args1... args1, Args2... args2) {
+      T::SetParticleData(args1...);
+    }
+
+    template <typename... Args1, typename... Args2>
+    void SetParticleData(T& p, Args1... args1, Args2... args2) {}
+  };
+
 } // namespace corsika::stack
 
 #endif

Framework/StackInterface/SecondaryView.h

@@ -3,7 +3,6 @@
 
 #include <corsika/stack/Stack.h>
 
-#include <algorithm>
 #include <vector>
 
 namespace corsika::stack {
@@ -176,6 +175,11 @@ namespace corsika::stack {
     bool IsEmpty() { return GetSize() == 0; }
 
   protected:
+    /**
+     * We only want to 'see' secondaries indexed in fIndices. In this
+     * function the conversion form iterator-index to stack-index is
+     * performed.
+     */
     unsigned int GetIndexFromIterator(const unsigned int vI) const {
       if (vI == 0) return fProjectileIndex;
       return fIndices[vI - 1];

Framework/StackInterface/Stack.h

@@ -49,7 +49,7 @@ namespace corsika::stack {
      <b>Important:</b> ParticleInterface must inherit from ParticleBase !
    */
 
-  template <typename>
+  template <typename>      //, bool>
   class ParticleInterface; // forward decl
 
   /**
@@ -83,7 +83,11 @@ namespace corsika::stack {
         delete; ///< since Stack can be very big, we don't want to copy it
 
   public:
-    //    template<typename = std::enable_if_t<std::is_reference<StackDataType>{}>>
+    /**
+     * if StackDataType is a reference member we *have* to initialize
+     * it in the constructor, this is typically needed for SecondaryView
+     */
+    template <typename = std::enable_if_t<std::is_reference<StackDataType>{}>>
     Stack(StackDataType vD)
         : fData(vD) {}
 
@@ -96,15 +100,13 @@ namespace corsika::stack {
               typename = std::enable_if_t<!std::is_reference<StackDataType>{}>>
     Stack(Args... args)
         : fData(args...) {}
-    // , typename std::enable_if<!std::is_reference<StackDataType>::value,
-    // std::nullptr_t>::type = nullptr)
 
   public:
     typedef StackDataType
         StackImpl; ///< this is the type of the user-provided data structure
 
-    template <typename SI>
-    using PIType = ParticleInterface<SI>;
+    template <typename SI>                //, bool IsBase>
+    using PIType = ParticleInterface<SI>; //, IsBase>;
 
     /**
      * Via the StackIteratorInterface and ConstStackIteratorInterface
@@ -125,18 +127,21 @@ namespace corsika::stack {
      */
     typedef typename StackIterator::ParticleInterfaceType ParticleType;
 
+    // friends are needed since they need access to protected members
     friend class StackIteratorInterface<
         typename std::remove_reference<StackDataType>::type, ParticleInterface,
         StackType>;
-
     friend class ConstStackIteratorInterface<
         typename std::remove_reference<StackDataType>::type, ParticleInterface,
         StackType>;
 
   public:
+    /**
+     * @name Most generic proxy methods for StackDataType fData
+     * @{
+     */
     unsigned int GetCapacity() const { return fData.GetCapacity(); }
     unsigned int GetSize() const { return fData.GetSize(); }
-
     template <typename... Args>
     auto Init(Args... args) {
       return fData.Init(args...);
@@ -145,6 +150,7 @@ namespace corsika::stack {
     auto Clear(Args... args) {
       return fData.Clear(args...);
     }
+    ///@}
 
   public:
     /**
@@ -231,12 +237,23 @@ namespace corsika::stack {
     // typename std::enable_if<HasGetIndexFromIterator<T>::value, unsigned int>::type
     // typename std::enable_if<std::is_base_of<decltype(*this)>,
     // SecondaryView<StackDataType, ParticleInterface>>::value, unsigned int>::type
+    /**
+     * Function to perform eventual transformation from
+     * StackIterator::GetIndex() to index in data stored in
+     * StackDataType fData. By default (and in almost all cases) this
+     * should just be identiy. See class SecondaryView for an alternative implementation.
+     */
     unsigned int GetIndexFromIterator(const unsigned int vI) const { return vI; }
 
+    /**
+     * @name Return reference to StackDataType object fData for data access
+     * @{
+     */
     typename std::remove_reference<StackDataType>::type& GetStackData() { return fData; }
     const typename std::remove_reference<StackDataType>::type& GetStackData() const {
       return fData;
     }
+    ///@}
   };
 
 } // namespace corsika::stack

Framework/StackInterface/StackIteratorInterface.h

@@ -42,31 +42,38 @@ namespace corsika::stack {
      the StackIteratorInterface. In addition to Stack the iterator only knows
      the index fIndex in the Stack data.
 
-     The template argument Particles acts as a policy to provide
-     readout function of Particle data from the stack. The Particle
+     The template argument `ParticleInterface` acts as a policy to provide
+     readout function of Particle data from the stack. The ParticleInterface
      class must know how to retrieve information from the Stack data
      for a particle entry at any index fIndex.
+
+     The ParticleInterface class must be written and provided by the
+     user, it contains methods like <code> auto GetData() const {
+     return GetStackData().GetData(GetIndex()); }</code>, where
+     StackIteratorInterface::GetStackData() return a reference to the
+     object storing the particle data of type StackDataType. And
+     StackIteratorInterface::GetIndex() provides the iterator index to
+     be readout. The StackDataType is another user-provided class to
+     store data and must implement functions compatible with
+     ParticleInterface, in this example StackDataType::GetData(const unsigned int
+     vIndex).
+
+     For an example see stack_example.cc, or the
+     corsika::processes::sibyll::SibStack class
   */
 
   template <typename StackDataType, template <typename> typename ParticleInterface,
-            typename StackType = Stack<StackDataType, ParticleInterface>>
+            typename StackType =
+                Stack<StackDataType, ParticleInterface>> //, bool IsBase=true >
   class StackIteratorInterface
-      : public ParticleInterface<
-            StackIteratorInterface<StackDataType, ParticleInterface, StackType>> {
+      : public ParticleInterface<StackIteratorInterface<StackDataType, ParticleInterface,
+                                                        StackType>> { //,IsBase> {
 
   public:
-    typedef StackDataType SD;
-    // typedef Stack<StackDataType, ParticleInterface> StackType;
-    /*typedef
-        typename std::conditional<std::is_const<StackDataType>::value,
-                                  const Stack<const StackDataType, ParticleInterface>&,
-                                  Stack<StackDataType, ParticleInterface>&>::type
-       StackType;*/
-
-    typedef ParticleInterface<
-        StackIteratorInterface<StackDataType, ParticleInterface, StackType>>
-        ParticleInterfaceType;
+    using ParticleInterfaceType = ParticleInterface<
+        StackIteratorInterface<StackDataType, ParticleInterface, StackType>>; //,IsBase>;
 
+    // friends are needed for access to protected methods
     friend class Stack<StackDataType,
                        ParticleInterface>; // for access to GetIndex for Stack
     friend class Stack<StackDataType&, ParticleInterface>; // for access to GetIndex
@@ -123,7 +130,8 @@ namespace corsika::stack {
         : fIndex(index)
         , fData(&data) {
       ParticleInterfaceType& p = **this;
-      p.SetParticleData(*parent, args...);
+      ParticleInterfaceType& pa = *parent;
+      p.SetParticleData(pa, args...);
     }
 
   public:
@@ -192,12 +200,12 @@ namespace corsika::stack {
   template <typename StackDataType, template <typename> typename ParticleInterface,
             typename StackType = Stack<StackDataType, ParticleInterface>>
   class ConstStackIteratorInterface
-      : public ParticleInterface<
-            ConstStackIteratorInterface<StackDataType, ParticleInterface, StackType>> {
+      : public ParticleInterface<ConstStackIteratorInterface<
+            StackDataType, ParticleInterface, StackType>> { //,IsBase> {
 
   public:
-    typedef ParticleInterface<
-        ConstStackIteratorInterface<StackDataType, ParticleInterface, StackType>>
+    typedef ParticleInterface<ConstStackIteratorInterface<
+        StackDataType, ParticleInterface, StackType>> //,IsBase>
         ParticleInterfaceType;
 
     friend class Stack<StackDataType, ParticleInterface>;   // for access to GetIndex

Framework/StackInterface/comp

+TestParticleInterface2<TestParticleInterface<corsika::stack::StackIteratorInterface<corsika::stack::CombinedStackImpl<corsika::stack::CombinedStackImpl<TestStackData, TestStackData2>, TestStackData3>, CombinedTestInterfaceType2, corsika::stack::Stack<corsika::stack::CombinedStackImpl<corsika::stack::CombinedStackImpl<TestStackData, TestStackData2>, TestStackData3>, CombinedTestInterfaceType2> > > >::content
+

Framework/StackInterface/testCombinedStack.cc

+
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
+#include <corsika/stack/CombinedStack.h>
+#include <corsika/stack/Stack.h>
+
+#include <testTestStack.h> // for testing: simple stack. This is a
+// test-build, and inluce file is obtained from CMAKE_CURRENT_SOURCE_DIR
+
+#include <boost/type_index.hpp>
+#include <type_traits>
+using boost::typeindex::type_id_with_cvr;
+
+#include <iomanip>
+#include <iostream>
+#include <vector>
+
+#define CATCH_CONFIG_MAIN // This tells Catch to provide a main() - only do this in one
+                          // cpp file
+#include <catch2/catch.hpp>
+
+using namespace corsika::stack;
+using namespace std;
+
+// definition of stack-data object
+class TestStackData2 {
+
+public:
+  // these functions are needed for the Stack interface
+  void Init() {}
+  void Clear() { fData2.clear(); }
+  unsigned int GetSize() const { return fData2.size(); }
+  unsigned int GetCapacity() const { return fData2.size(); }
+  void Copy(const int i1, const int i2) { fData2[i2] = fData2[i1]; }
+  void Swap(const int i1, const int i2) {
+    double tmp0 = fData2[i1];
+    fData2[i1] = fData2[i2];
+    fData2[i2] = tmp0;
+  }
+
+  // custom data access function
+  void SetData2(const int i, const double v) { fData2[i] = v; }
+  double GetData2(const int i) const { return fData2[i]; }
+
+  // these functions are also needed by the Stack interface
+  void IncrementSize() { fData2.push_back(0.); }
+  void DecrementSize() {
+    if (fData2.size() > 0) { fData2.pop_back(); }
+  }
+
+  // custom private data section
+private:
+  std::vector<double> fData2;
+};
+
+// defintion of a stack-readout object, the iteractor dereference
+// operator will deliver access to these function
+template <typename T>
+class TestParticleInterface2
+    : public T { // corsika::stack::ParticleBaseAdd<T> {//public T {
+
+public:
+  using T::GetIndex;
+  using T::GetStackData;
+  using T::SetParticleData;
+
+  // default version for particle-creation from input data
+  void SetParticleData(const std::tuple<double> v = {0.}) { SetData2(std::get<0>(v)); }
+  void SetParticleData(TestParticleInterface2<T>& parent,
+                       const std::tuple<double> v = {0.}) {
+    SetData2(parent.GetData2() + std::get<0>(v));
+  }
+  void SetData2(const double v) { GetStackData().SetData2(GetIndex(), v); }
+  double GetData2() const { return GetStackData().GetData2(GetIndex()); }
+};
+
+// combined stack
+template <typename StackIter>
+using CombinedTestInterfaceType =
+    corsika::stack::CombinedParticleInterface<TestParticleInterface,
+                                              TestParticleInterface2, StackIter>;
+
+using StackTest = CombinedStack<TestStackData, TestStackData2, CombinedTestInterfaceType>;
+typedef StackTest::ParticleType Particle;
+
+TEST_CASE("Combined Stack", "[stack]") {
+
+  // helper function for sum over stack data
+  auto sum = [](const StackTest& stack) {
+    double v = 0;
+    for (const auto& p : stack) v += p.GetData();
+    return v;
+  };
+  auto sum2 = [](const StackTest& stack) {
+    double v = 0;
+    for (const auto& p : stack) v += p.GetData2();
+    return v;
+  };
+
+  SECTION("StackInterface") {
+
+    // construct a valid Stack object
+    StackTest s;
+    s.Init();
+    s.Clear();
+    s.AddParticle(std::tuple{0.});
+    s.Copy(s.cbegin(), s.begin());
+    s.Swap(s.begin(), s.begin());
+    REQUIRE(s.GetSize() == 1);
+  }
+
+  SECTION("construct") {
+
+    // construct a valid, empty Stack object
+    StackTest s;
+  }
+
+  SECTION("write and read") {
+
+    StackTest s;
+    s.AddParticle(std::tuple{9.9});
+    REQUIRE(sum2(s) == 0.);
+    REQUIRE(sum(s) == 9.9);
+  }
+
+  SECTION("delete from stack") {
+
+    StackTest s;
+    REQUIRE(s.GetSize() == 0);
+    StackTest::StackIterator p =
+        s.AddParticle(std::tuple{0.}); // valid way to access particle data
+    p.SetData(8.9);
+    p.SetData2(3.);
+    REQUIRE(sum2(s) == 3.);
+    REQUIRE(sum(s) == 8.9);
+    REQUIRE(s.GetSize() == 1);
+    s.Delete(p);
+    REQUIRE(s.GetSize() == 0);
+  }
+
+  SECTION("delete particle") {
+
+    StackTest s;
+    REQUIRE(s.GetSize() == 0);
+    auto p = s.AddParticle(
+        std::tuple{9.9}); // also valid way to access particle data, identical to above
+    REQUIRE(s.GetSize() == 1);
+    p.Delete();
+    REQUIRE(s.GetSize() == 0);
+  }
+
+  SECTION("create secondaries") {
+    StackTest s;
+    REQUIRE(s.GetSize() == 0);
+    auto iter = s.AddParticle(std::tuple{9.9});
+    Particle& p = *iter; // also this is valid to access particle data
+    p.SetData2(2);
+    REQUIRE(s.GetSize() == 1);
+    p.AddSecondary(std::tuple{4.4});
+    REQUIRE(s.GetSize() == 2);
+    // p.AddSecondary(3.3, 2.2, 1.);
+    // REQUIRE(s.GetSize() == 3);
+    double v = 0;
+    for (const auto& i : s) {
+      v += i.GetData();
+      REQUIRE(i.GetData2() == 2);
+    }
+    REQUIRE(v == 9.9 + 4.4);
+  }
+
+  SECTION("get next particle") {
+    StackTest s;
+    REQUIRE(s.GetSize() == 0);
+    auto p1 = s.AddParticle(std::tuple{9.9});
+    auto p2 = s.AddParticle(std::tuple{8.8});
+    p1.SetData2(20.2);
+    p2.SetData2(20.3);
+    auto particle = s.GetNextParticle(); // first particle
+    REQUIRE(particle.GetData() == 8.8);
+    REQUIRE(particle.GetData2() == 20.3);
+
+    particle.Delete();
+    auto particle2 = s.GetNextParticle(); // first particle
+    REQUIRE(particle2.GetData() == 9.9);
+    REQUIRE(particle2.GetData2() == 20.2);
+    particle2.Delete();
+
+    REQUIRE(s.GetSize() == 0);
+  }
+}
+
+////////////////////////////////////////////////////////////
+
+// definition of stack-data object
+class TestStackData3 {
+
+public:
+  // these functions are needed for the Stack interface
+  void Init() {}
+  void Clear() { fData3.clear(); }
+  unsigned int GetSize() const { return fData3.size(); }
+  unsigned int GetCapacity() const { return fData3.size(); }
+  void Copy(const int i1, const int i2) { fData3[i2] = fData3[i1]; }
+  void Swap(const int i1, const int i2) {
+    double tmp0 = fData3[i1];
+    fData3[i1] = fData3[i2];
+    fData3[i2] = tmp0;
+  }
+
+  // custom data access function
+  void SetData3(const int i, const double v) { fData3[i] = v; }
+  double GetData3(const int i) const { return fData3[i]; }
+
+  // these functions are also needed by the Stack interface
+  void IncrementSize() { fData3.push_back(0.); }
+  void DecrementSize() {
+    if (fData3.size() > 0) { fData3.pop_back(); }
+  }
+
+  // custom private data section
+private:
+  std::vector<double> fData3;
+};
+
+// defintion of a stack-readout object, the iteractor dereference
+// operator will deliver access to these function
+template <typename T>
+class TestParticleInterface3 : public T {
+
+public:
+  using T::GetIndex;
+  using T::GetStackData;
+  using T::SetParticleData;
+
+  // default version for particle-creation from input data
+  void SetParticleData(const std::tuple<double> v = {0.}) { SetData3(std::get<0>(v)); }
+  void SetParticleData(TestParticleInterface3<T>& parent,
+                       const std::tuple<double> v = {0.}) {
+    SetData3(parent.GetData3() + std::get<0>(v));
+  }
+  void SetData3(const double v) { GetStackData().SetData3(GetIndex(), v); }
+  double GetData3() const { return GetStackData().GetData3(GetIndex()); }
+};
+
+// double combined stack
+// combined stack
+template <typename StackIter>
+using CombinedTestInterfaceType2 =
+    corsika::stack::CombinedParticleInterface<CombinedTestInterfaceType,
+                                              TestParticleInterface3, StackIter>;
+
+using StackTest2 = CombinedStack<typename StackTest::StackImpl, TestStackData3,
+                                 CombinedTestInterfaceType2>;
+typedef StackTest2::ParticleType Particle2;
+
+TEST_CASE("Combined Stack - multi", "[stack]") {
+
+  SECTION("create secondaries") {
+
+    StackTest2 s;
+    REQUIRE(s.GetSize() == 0);
+    auto p1 = s.AddParticle(std::tuple{9.9});
+    auto p2 = s.AddParticle(std::tuple{8.8}, std::tuple{0.1});
+    p2.SetData2(0.1); // not clear why this is needed, need to check
+		      // SetParticleData workflow for more complicated
+		      // settings
+    // auto p3 = s.AddParticle( std::tuple {8.8}, std::tuple{1.}, std::tuple{0.1} );
+    p1.SetData3(20.2);
+    p2.SetData3(10.3);
+
+    REQUIRE(p1.GetData() == 9.9);
+    REQUIRE(p1.GetData2() == 0.);
+    p1.SetData2(10.2);
+    REQUIRE(p1.GetData2() == 10.2);
+    REQUIRE(p1.GetData3() == 20.2);
+
+    REQUIRE(p2.GetData() == 8.8);
+    REQUIRE(p2.GetData2() == 0.1);
+    REQUIRE(p2.GetData3() == 10.3);
+
+    auto particle = s.GetNextParticle(); // first particle
+    REQUIRE(particle.GetData() == 8.8);
+    REQUIRE(particle.GetData2() == 0.1);
+    REQUIRE(particle.GetData3() == 10.3);
+
+    REQUIRE(s.GetSize() == 2);
+    auto sec = particle.AddSecondary(std::tuple{4.4});
+    REQUIRE(s.GetSize() == 3);
+    REQUIRE(sec.GetData() == 4.4);
+    REQUIRE(sec.GetData2() == 0.1);
+    REQUIRE(sec.GetData3() == 10.3);
+
+    sec.Delete();
+    s.DeleteLast();
+    s.GetNextParticle().Delete();
+    REQUIRE(s.GetSize() == 0);
+  }
+}

Framework/StackInterface/testSecondaryView.cc

+
+/*
+ * (c) Copyright 2018 CORSIKA Project, corsika-project@lists.kit.edu
+ *
+ * See file AUTHORS for a list of contributors.
+ *
+ * This software is distributed under the terms of the GNU General Public
+ * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
+ * the license.
+ */
+
+#include <corsika/stack/SecondaryView.h>
+#include <corsika/stack/Stack.h>
+
+#include <testTestStack.h> // for testing: simple stack. This is a
+// test-build, and inluce file is obtained from CMAKE_CURRENT_SOURCE_DIR
+
+#include <boost/type_index.hpp>
+#include <type_traits>
+using boost::typeindex::type_id_with_cvr;
+
+#include <iomanip>
+#include <iostream>
+#include <vector>
+
+#define CATCH_CONFIG_MAIN // This tells Catch to provide a main() - only do this in one
+                          // cpp file
+#include <catch2/catch.hpp>
+
+using namespace corsika::stack;
+using namespace std;
+
+typedef Stack<TestStackData, TestParticleInterface> StackTest;
+typedef StackTest::ParticleType Particle;
+
+TEST_CASE("SecondaryStack", "[stack]") {
+
+  // helper function for sum over stack data
+  auto sum = [](const StackTest& stack) {
+    double v = 0;
+    for (const auto& p : stack) v += p.GetData();
+    return v;
+  };
+
+  SECTION("secondary view") {
+    StackTest s;
+    REQUIRE(s.GetSize() == 0);
+    s.AddParticle(std::tuple{9.9});
+    s.AddParticle(std::tuple{8.8});
+    const double sumS = 9.9 + 8.8;
+
+    auto particle = s.GetNextParticle();
+
+    typedef SecondaryView<TestStackData, TestParticleInterface> StackTestView;
+    StackTestView v(particle);
+    REQUIRE(v.GetSize() == 0);
+
+    {
+      auto proj = v.GetProjectile();
+      REQUIRE(proj.GetData() == particle.GetData());
+    }
+
+    v.AddSecondary(std::tuple{4.4});
+    v.AddSecondary(std::tuple{4.5});
+    v.AddSecondary(std::tuple{4.6});
+
+    REQUIRE(v.GetSize() == 3);
+    REQUIRE(s.GetSize() == 5);
+    REQUIRE(!v.IsEmpty());
+
+    auto sumView = [](const StackTestView& stack) {
+      double v = 0;
+      for (const auto& p : stack) { v += p.GetData(); }
+      return v;
+    };
+
+    REQUIRE(sum(s) == sumS + 4.4 + 4.5 + 4.6);
+    REQUIRE(sumView(v) == 4.4 + 4.5 + 4.6);
+
+    v.DeleteLast();
+    REQUIRE(v.GetSize() == 2);
+    REQUIRE(s.GetSize() == 4);
+
+    REQUIRE(sum(s) == sumS + 4.4 + 4.5);
+    REQUIRE(sumView(v) == 4.4 + 4.5);
+
+    auto pDel = v.GetNextParticle();
+    v.Delete(pDel);
+    REQUIRE(v.GetSize() == 1);
+    REQUIRE(s.GetSize() == 3);
+
+    REQUIRE(sum(s) == sumS + 4.4 + 4.5 - pDel.GetData());
+    REQUIRE(sumView(v) == 4.4 + 4.5 - pDel.GetData());
+
+    v.Delete(v.GetNextParticle());
+    REQUIRE(sum(s) == sumS);
+    REQUIRE(sumView(v) == 0);
+    REQUIRE(v.IsEmpty());
+
+    {
+      auto proj = v.GetProjectile();
+      REQUIRE(proj.GetData() == particle.GetData());
+    }
+  }
+}

Framework/StackInterface/testStackInterface.cc

@@ -9,20 +9,19 @@
  * the license.
  */
 
-#include <corsika/stack/SecondaryView.h>
 #include <corsika/stack/Stack.h>
 
+#include <testTestStack.h> // simple test-stack for testing. This is
+                           // for testing only: include from
+                           // CMAKE_CURRENT_SOURCE_DIR
+
 #include <boost/type_index.hpp>
 #include <type_traits>
 using boost::typeindex::type_id_with_cvr;
 
-template <typename T>
-struct ID {
-  typedef T type;
-};
-
 #include <iomanip>
 #include <iostream>
+#include <tuple>
 #include <vector>
 
 #define CATCH_CONFIG_MAIN // This tells Catch to provide a main() - only do this in one
@@ -32,73 +31,6 @@ struct ID {
 using namespace corsika::stack;
 using namespace std;
 
-// definition of stack-data object
-class TestStackData {
-
-public:
-  // these functions are needed for the Stack interface
-  void Init() {}
-  void Dump() const {}
-  void Clear() { fData.clear(); }
-  unsigned int GetSize() const { return fData.size(); }
-  unsigned int GetCapacity() const { return fData.size(); }
-  void Copy(const int i1, const int i2) { fData[i2] = fData[i1]; }
-  void Swap(const int i1, const int i2) {
-    double tmp0 = fData[i1];
-    fData[i1] = fData[i2];
-    fData[i2] = tmp0;
-  }
-
-  // custom data access function
-  void SetData(const int i, const double v) { fData[i] = v; }
-  double GetData(const int i) const { return fData[i]; }
-
-  // these functions are also needed by the Stack interface
-  void IncrementSize() { fData.push_back(0.); }
-  void DecrementSize() {
-    if (fData.size() > 0) { fData.pop_back(); }
-  }
-
-  // custom private data section
-private:
-  std::vector<double> fData;
-};
-
-// defintion of a stack-readout object, the iteractor dereference
-// operator will deliver access to these function
-template <typename StackIteratorInterface>
-class TestParticleInterface : public ParticleBase<StackIteratorInterface> {
-  using ParticleBase<StackIteratorInterface>::GetStack;
-  using ParticleBase<StackIteratorInterface>::GetStackData;
-  using ParticleBase<StackIteratorInterface>::GetIndex;
-  using ParticleBase<StackIteratorInterface>::GetIterator;
-
-public:
-  /* normally this function does not need to be specified here, it is
-   part of ParticleBase<StackIteratorInterface>.
-
-   However, when overloading this function again with a different
-   parameter set, as here, seems to require also the declaration of
-   the original function here...
-  */
-
-  // default version for particle-creation from input data
-  void SetParticleData(const double v) { SetData(v); }
-  void SetParticleData(TestParticleInterface<StackIteratorInterface>& /*parent*/,
-                       const double v) {
-    SetData(v);
-  }
-  /// alternative set-particle data for non-standard construction from different inputs
-  void SetParticleData(const double v, const double p) { SetData(v + p); }
-  void SetParticleData(TestParticleInterface<StackIteratorInterface>& /*parent*/,
-                       const double v, const double p) {
-    SetData(v + p);
-  }
-
-  void SetData(const double v) { GetStackData().SetData(GetIndex(), v); }
-  double GetData() const { return GetStackData().GetData(GetIndex()); }
-};
-
 typedef Stack<TestStackData, TestParticleInterface> StackTest;
 typedef StackTest::ParticleType Particle;
 
@@ -117,7 +49,7 @@ TEST_CASE("Stack", "[Stack]") {
     StackTest s;
     s.Init();
     s.Clear();
-    s.AddParticle(0.);
+    s.AddParticle(std::tuple{0.});
     s.Copy(s.cbegin(), s.begin());
     s.Swap(s.begin(), s.begin());
     REQUIRE(s.GetSize() == 1);
@@ -132,7 +64,7 @@ TEST_CASE("Stack", "[Stack]") {
   SECTION("write and read") {
 
     StackTest s;
-    s.AddParticle(9.9);
+    s.AddParticle(std::tuple{9.9});
     const double v = sum(s);
     REQUIRE(v == 9.9);
   }
@@ -141,7 +73,8 @@ TEST_CASE("Stack", "[Stack]") {
 
     StackTest s;
     REQUIRE(s.GetSize() == 0);
-    StackTest::StackIterator p = s.AddParticle(0.); // valid way to access particle data
+    StackTest::StackIterator p =
+        s.AddParticle(std::tuple{0.}); // valid way to access particle data
     p.SetData(9.9);
     REQUIRE(s.GetSize() == 1);
     s.Delete(p);
@@ -152,8 +85,8 @@ TEST_CASE("Stack", "[Stack]") {
 
     StackTest s;
     REQUIRE(s.GetSize() == 0);
-    auto p =
-        s.AddParticle(9.9); // also valid way to access particle data, identical to above
+    auto p = s.AddParticle(
+        std::tuple{9.9}); // also valid way to access particle data, identical to above
     REQUIRE(s.GetSize() == 1);
     p.Delete();
     REQUIRE(s.GetSize() == 0);
@@ -163,23 +96,23 @@ TEST_CASE("Stack", "[Stack]") {
 
     StackTest s;
     REQUIRE(s.GetSize() == 0);
-    auto iter = s.AddParticle(9.9);
+    auto iter = s.AddParticle(std::tuple{9.9});
     Particle& p = *iter; // also this is valid to access particle data
     REQUIRE(s.GetSize() == 1);
-    p.AddSecondary(4.4);
+    p.AddSecondary(std::tuple{4.4});
     REQUIRE(s.GetSize() == 2);
-    p.AddSecondary(3.3, 2.2);
+    /*p.AddSecondary(3.3, 2.2);
     REQUIRE(s.GetSize() == 3);
     double v = 0;
     for (auto& p : s) { v += p.GetData(); }
-    REQUIRE(v == 9.9 + 4.4 + 3.3 + 2.2);
+    REQUIRE(v == 9.9 + 4.4 + 3.3 + 2.2);*/
   }
 
   SECTION("get next particle") {
     StackTest s;
     REQUIRE(s.GetSize() == 0);
-    s.AddParticle(9.9);
-    s.AddParticle(8.8);
+    s.AddParticle(std::tuple{9.9});
+    s.AddParticle(std::tuple{8.8});
     auto particle = s.GetNextParticle(); // first particle
     REQUIRE(particle.GetData() == 8.8);
 
@@ -190,65 +123,4 @@ TEST_CASE("Stack", "[Stack]") {
 
     REQUIRE(s.GetSize() == 0);
   }
-
-  SECTION("secondary view") {
-    StackTest s;
-    REQUIRE(s.GetSize() == 0);
-    s.AddParticle(9.9);
-    s.AddParticle(8.8);
-    const double sumS = 9.9 + 8.8;
-
-    auto particle = s.GetNextParticle();
-
-    typedef SecondaryView<TestStackData, TestParticleInterface> StackTestView;
-    StackTestView v(particle);
-    REQUIRE(v.GetSize() == 0);
-
-    {
-      auto proj = v.GetProjectile();
-      REQUIRE(proj.GetData() == particle.GetData());
-    }
-
-    v.AddSecondary(4.4);
-    v.AddSecondary(4.5);
-    v.AddSecondary(4.6);
-
-    REQUIRE(v.GetSize() == 3);
-    REQUIRE(s.GetSize() == 5);
-    REQUIRE(!v.IsEmpty());
-
-    auto sumView = [](const StackTestView& stack) {
-      double v = 0;
-      for (const auto& p : stack) { v += p.GetData(); }
-      return v;
-    };
-
-    REQUIRE(sum(s) == sumS + 4.4 + 4.5 + 4.6);
-    REQUIRE(sumView(v) == 4.4 + 4.5 + 4.6);
-
-    v.DeleteLast();
-    REQUIRE(v.GetSize() == 2);
-    REQUIRE(s.GetSize() == 4);
-
-    REQUIRE(sum(s) == sumS + 4.4 + 4.5);
-    REQUIRE(sumView(v) == 4.4 + 4.5);
-
-    auto pDel = v.GetNextParticle();
-    v.Delete(pDel);
-    REQUIRE(v.GetSize() == 1);
-    REQUIRE(s.GetSize() == 3);
-
-    REQUIRE(sum(s) == sumS + 4.4 + 4.5 - pDel.GetData());
-    REQUIRE(sumView(v) == 4.4 + 4.5 - pDel.GetData());
-
-    v.Delete(v.GetNextParticle());
-    REQUIRE(sum(s) == sumS);
-    REQUIRE(sumView(v) == 0);
-    REQUIRE(v.IsEmpty());
-
-    {
-      auto proj = v.GetProjectile();
-      REQUIRE(proj.GetData() == particle.GetData());
-    }
-  }
 }

Framework/StackInterface/testTestStack.h

+#ifndef _include_corsika_stack_testTestStack_h_
+#define _include_corsika_stack_testTestStack_h_
+
+#include <corsika/stack/Stack.h>
+#include <tuple>
+#include <vector>
+
+/**
+ * definition of stack-data object for unit tests
+ *
+ * TestStackData contain only a single variable "Data" stored in fData
+ * with Get/SetData functions.
+ */
+class TestStackData {
+
+public:
+  // these functions are needed for the Stack interface
+  void Init() {}
+  void Clear() { fData.clear(); }
+  unsigned int GetSize() const { return fData.size(); }
+  unsigned int GetCapacity() const { return fData.size(); }
+  void Copy(const unsigned int i1, const unsigned int i2) { fData[i2] = fData[i1]; }
+  void Swap(const unsigned int i1, const unsigned int i2) {
+    double tmp0 = fData[i1];
+    fData[i1] = fData[i2];
+    fData[i2] = tmp0;
+  }
+
+  // custom data access function
+  void SetData(const unsigned int i, const double v) { fData[i] = v; }
+  double GetData(const unsigned int i) const { return fData[i]; }
+
+  // these functions are also needed by the Stack interface
+  void IncrementSize() { fData.push_back(0.); }
+  void DecrementSize() {
+    if (fData.size() > 0) { fData.pop_back(); }
+  }
+
+  // custom private data section
+private:
+  std::vector<double> fData;
+};
+
+/**
+ * From static_cast of a StackIterator over entries in the
+ * TestStackData class you get and object of type
+ * TestParticleInterface defined here
+ *
+ * It provides Get/Set methods to read and write data to the "Data"
+ * storage of TestStackData obtained via
+ * "StackIteratorInterface::GetStackData()", given the index of the
+ * iterator "StackIteratorInterface::GetIndex()"
+ *
+ */
+template <typename StackIteratorInterface>
+class TestParticleInterface
+    : public corsika::stack::ParticleBase<StackIteratorInterface> {
+
+public:
+  using corsika::stack::ParticleBase<StackIteratorInterface>::GetStackData;
+  using corsika::stack::ParticleBase<StackIteratorInterface>::GetIndex;
+
+  /*
+     The SetParticleData methods are called for creating new entries
+     on the stack. You can specifiy various parametric versions to
+     perform this task:
+  */
+
+  // default version for particle-creation from input data
+  void SetParticleData(const std::tuple<double> v) { SetData(std::get<0>(v)); }
+  void SetParticleData(TestParticleInterface<StackIteratorInterface>& /*parent*/,
+                       std::tuple<double> v) {
+    SetData(std::get<0>(v));
+  }
+  /// alternative set-particle data for non-standard construction from different inputs
+  /*
+  void SetParticleData(const double v, const double p) { SetData(v + p); }
+  void SetParticleData(TestParticleInterface<StackIteratorInterface>&,
+                       const double v, const double p) {
+    SetData(v + p);
+    }*/
+
+  // here are the fundamental methods for access to TestStackData data
+  void SetData(const double v) { GetStackData().SetData(GetIndex(), v); }
+  double GetData() const { return GetStackData().GetData(GetIndex()); }
+};
+
+#endif

nuc/Makefile_ifc

+FORTRAN=gfortran
+#FLAGS=-r8 -C -w95
+FLAGS=
+
+all : main_f77
+
+main04.o : main04.f
+	$(FORTRAN) $(FLAGS) -c $<
+
+main_f77 : main04.o
+	$(FORTRAN) $(FLAGS) -o $@ $^
+
+clean :
+	-rm -f *.o *.mod *~ main

nuc/main04.f

+      PROGRAM main
+      COMMON
+     *      /DATABA_SIZVEC/ NUMsphe, NUModd,  NUMdefo,
+     *                      NUMd3n,  NUMd3p,  NUMsupd,
+     *                      NUMmono, NUMdipo, NUMtwop,
+     *                      NUMdvpn, NUMterm, NUMqpSH
+C
+      CALL InitializeData()
+C
+      CALL GetData()
+C
+      CALL PrintData()
+C
+      STOP
+C
+      END
+C
+      INCLUDE 'input04.f'
+
+