Feature/integration part6

cmake/dependencies.cmake

@@ -11,6 +11,8 @@ FetchContent_Declare( eigen
     GIT_REPOSITORY  https://gitlab.com/libeigen/eigen.git
     GIT_TAG         d0bfdc1658ca0b4c659fd3702c351d2c2cdc876c # 3.4.1 branch on July 26, 2023
     )
+# the next line prevents eigen from populating the .cmake/packages folder
+set( CMAKE_EXPORT_NO_PACKAGE_REGISTRY ON )
 
 #######################################################################
 # Load dependencies

cmake/unit_testing.cmake

@@ -39,6 +39,9 @@ add_subdirectory( src/scion/integration/LinearLinear/test )
 add_subdirectory( src/scion/integration/LinearLogarithmic/test )
 add_subdirectory( src/scion/integration/LogarithmicLinear/test )
 add_subdirectory( src/scion/integration/LogarithmicLogarithmic/test )
+add_subdirectory( src/scion/integration/FirstMomentHistogram/test )
+add_subdirectory( src/scion/integration/FirstMomentLinearLinear/test )
+add_subdirectory( src/scion/integration/FirstMomentLinearLogarithmic/test )
 add_subdirectory( src/scion/integration/GaussLegendre/2/test )
 add_subdirectory( src/scion/integration/GaussLegendre/4/test )
 add_subdirectory( src/scion/integration/GaussLegendre/6/test )

src/scion/integration.hpp

@@ -4,5 +4,6 @@
 #include "scion/integration/LinearLogarithmic.hpp"
 #include "scion/integration/LogarithmicLinear.hpp"
 #include "scion/integration/LogarithmicLogarithmic.hpp"
+#include "scion/integration/FirstMomentLinearLinear.hpp"
 #include "scion/integration/GaussLegendre.hpp"
 #include "scion/integration/GaussLobatto.hpp"

src/scion/integration/FirstMomentHistogram.hpp

@@ -0,0 +1,54 @@
+#ifndef NJOY_SCION_INTEGRATION_FIRSTMOMENTHISTOGRAM
+#define NJOY_SCION_INTEGRATION_FIRSTMOMENTHISTOGRAM
+
+// system includes
+
+// other includes
+#include "scion/integration/IntegratorBase.hpp"
+
+namespace njoy {
+namespace scion {
+namespace integration {
+
+  /**
+   *  @class
+   *  @brief First raw moment of a histogram panel (y is constant in x)
+   *
+   *  The first raw moment or mean is defined as the integral of x * f(x)
+   */
+  class FirstMomentHistogram : public IntegratorBase< FirstMomentHistogram > {
+
+    /* friend declarations */
+    friend class IntegratorBase< FirstMomentHistogram >;
+
+    /* interface implementation functions */
+
+    /**
+     *  @brief Perform first raw moment integration of a histogram panel (y is constant in x)
+     *
+     *  @param[in] xLeft    the left value on the x interval
+     *  @param[in] xRight   the right value on the x interval
+     *  @param[in] yLeft    the left value on the y interval
+     *  @param[in] yRight   the right value on the y interval
+     */
+    template < typename X, typename Y,
+               typename I = decltype( std::declval< X >() * std::declval< X >() * std::declval< Y >() ) >
+    I integrate( const X& xLeft, const X& xRight,
+                 const Y& yLeft, const Y& ) const noexcept {
+
+      return 0.5 * yLeft * ( xRight - xLeft ) * ( xRight + xLeft );
+    }
+
+  public:
+
+    using IntegratorBase::operator();
+  };
+
+  // integration function
+  static constexpr FirstMomentHistogram firstMomentHistogram;
+
+} // integration namespace
+} // scion namespace
+} // njoy namespace
+
+#endif

src/scion/integration/FirstMomentHistogram/test/CMakeLists.txt

@@ -0,0 +1 @@
+add_cpp_test( integration.FirstMomentHistogram FirstMomentHistogram.test.cpp )

src/scion/integration/FirstMomentHistogram/test/FirstMomentHistogram.test.cpp

@@ -0,0 +1,57 @@
+// include Catch2
+#include <catch2/catch_test_macros.hpp>
+#include <catch2/matchers/catch_matchers_floating_point.hpp>
+using Catch::Matchers::WithinRel;
+
+// what we are testing
+#include "scion/integration/FirstMomentHistogram.hpp"
+
+// other includes
+
+// convenience typedefs
+using namespace njoy::scion;
+
+SCENARIO( "FirstMomentHistogram" ) {
+
+  GIVEN( "FirstMomentHistogram integration object" ) {
+
+    WHEN( "integrating an interval" ) {
+
+      integration::FirstMomentHistogram integrator{};
+
+      THEN( "the integration is performed correctly" ) {
+
+        double xLeft = 1.0;
+        double xRight = 2.0;
+        double yLeft = 1.0;
+        double yRight = 4.0;
+
+        // both y values are the same
+        CHECK_THAT( 1.5, WithinRel( integrator( xLeft, xRight, yLeft, yLeft ) ) );
+
+        // the y values are different
+        CHECK_THAT( 1.5, WithinRel( integrator( xLeft, xRight, yLeft, yRight ) ) );
+      } // THEN
+    } // WHEN
+  } // GIVEN
+
+  GIVEN( "linlin integration function" ) {
+
+    WHEN( "integrating an interval" ) {
+
+      THEN( "the integration is performed correctly" ) {
+
+        double xLeft = 1.0;
+        double xRight = 2.0;
+        double yLeft = 1.0;
+        double yRight = 4.0;
+
+        // both y values are the same
+        CHECK_THAT( 1.5, WithinRel( integration::firstMomentHistogram( xLeft, xRight, yLeft, yLeft ) ) );
+
+        // the y values are different
+        CHECK_THAT( 1.5, WithinRel( integration::firstMomentHistogram( xLeft, xRight, yLeft, yRight ) ) );
+      } // THEN
+    } // WHEN
+  } // GIVEN
+} // SCENARIO

src/scion/integration/FirstMomentLinearLinear.hpp

@@ -0,0 +1,58 @@
+#ifndef NJOY_SCION_INTEGRATION_FIRSTMOMENTLINEARLINEAR
+#define NJOY_SCION_INTEGRATION_FIRSTMOMENTLINEARLINEAR
+
+// system includes
+
+// other includes
+#include "scion/integration/IntegratorBase.hpp"
+
+namespace njoy {
+namespace scion {
+namespace integration {
+
+  /**
+   *  @class
+   *  @brief First raw moment of a linear-linear panel (y is linear in x)
+   *
+   *  The first raw moment or mean is defined as the integral of x * f(x)
+   */
+  class FirstMomentLinearLinear : public IntegratorBase< FirstMomentLinearLinear > {
+
+    /* friend declarations */
+    friend class IntegratorBase< FirstMomentLinearLinear >;
+
+    /* interface implementation functions */
+
+    /**
+     *  @brief Perform first raw moment integration of a linear-linear panel (y is linear in x)
+     *
+     *  @param[in] xLeft    the left value on the x interval
+     *  @param[in] xRight   the right value on the x interval
+     *  @param[in] yLeft    the left value on the y interval
+     *  @param[in] yRight   the right value on the y interval
+     */
+    template < typename X, typename Y,
+               typename I = decltype( std::declval< X >() * std::declval< X >() * std::declval< Y >() ) >
+    I integrate( const X& xLeft, const X& xRight,
+                 const Y& yLeft, const Y& yRight ) const noexcept {
+
+      auto delta = ( xRight - xLeft );
+      auto slope = ( yRight - yLeft ) / delta / 3.;
+      auto constant = 0.5 * ( xRight * yLeft - xLeft * yRight ) / delta;
+      return xRight * xRight * ( slope * xRight + constant )
+             - xLeft * xLeft * ( slope * xLeft + constant );
+    }
+
+  public:
+
+    using IntegratorBase::operator();
+  };
+
+  // integration function
+  static constexpr FirstMomentLinearLinear firstMomentLinLin;
+
+} // integration namespace
+} // scion namespace
+} // njoy namespace
+
+#endif

src/scion/integration/FirstMomentLinearLinear/test/CMakeLists.txt

@@ -0,0 +1 @@
+add_cpp_test( integration.FirstMomentLinearLinear FirstMomentLinearLinear.test.cpp )

src/scion/integration/FirstMomentLinearLinear/test/FirstMomentLinearLinear.test.cpp

@@ -0,0 +1,57 @@
+// include Catch2
+#include <catch2/catch_test_macros.hpp>
+#include <catch2/matchers/catch_matchers_floating_point.hpp>
+using Catch::Matchers::WithinRel;
+
+// what we are testing
+#include "scion/integration/FirstMomentLinearLinear.hpp"
+
+// other includes
+
+// convenience typedefs
+using namespace njoy::scion;
+
+SCENARIO( "FirstMomentLinearLinear" ) {
+
+  GIVEN( "FirstMomentLinearLinear integration object" ) {
+
+    WHEN( "integrating an interval" ) {
+
+      integration::FirstMomentLinearLinear integrator{};
+
+      THEN( "the integration is performed correctly" ) {
+
+        double xLeft = 1.0;
+        double xRight = 2.0;
+        double yLeft = 1.0;
+        double yRight = 4.0;
+
+        // both y values are the same
+        CHECK_THAT( 1.5, WithinRel( integrator( xLeft, xRight, yLeft, yLeft ) ) );
+
+        // the y values are different
+        CHECK_THAT( 4.0, WithinRel( integrator( xLeft, xRight, yLeft, yRight ) ) );
+      } // THEN
+    } // WHEN
+  } // GIVEN
+
+  GIVEN( "linlin integration function" ) {
+
+    WHEN( "integrating an interval" ) {
+
+      THEN( "the integration is performed correctly" ) {
+
+        double xLeft = 1.0;
+        double xRight = 2.0;
+        double yLeft = 1.0;
+        double yRight = 4.0;
+
+        // both y values are the same
+        CHECK_THAT( 1.5, WithinRel( integration::firstMomentLinLin( xLeft, xRight, yLeft, yLeft ) ) );
+
+        // the y values are different
+        CHECK_THAT( 4.0, WithinRel( integration::firstMomentLinLin( xLeft, xRight, yLeft, yRight ) ) );
+      } // THEN
+    } // WHEN
+  } // GIVEN
+} // SCENARIO

src/scion/integration/FirstMomentLinearLogarithmic.hpp

@@ -0,0 +1,72 @@
+#ifndef NJOY_SCION_INTEGRATION_FIRSTMOMENTLINEARLOGARITHMIC
+#define NJOY_SCION_INTEGRATION_FIRSTMOMENTLINEARLOGARITHMIC
+
+// system includes
+#include <cmath>
+
+// other includes
+#include "scion/integration/IntegratorBase.hpp"
+
+namespace njoy {
+namespace scion {
+namespace integration {
+
+  /**
+   *  @class
+   *  @brief First raw moment of a linear-logarithmic panel (y is linear in ln(x))
+   *
+   *  The first raw moment or mean is defined as the integral of x * f(x)
+   *
+   *  The first raw moment of a panel is calculated analytically using the following
+   *  primitive function:
+   *    int x ln(x) = x**2 / 2 ( ln(x) - 1 / 2 )
+   *
+   *  Since y is linear in ln(x), the function to be integrated is:
+   *    y = a x ln(x/x1) + y1 x with a = (y2 - y1) / ln(x2/x1)
+   *
+   *  The integral over the panel is then given by:
+   *    int[x1,x2] y dx = int[x1,x2] a x ln(x/x1) dx + y1 (x2**2 - x1**2) / 2
+   *
+   *  By substituting z = x/x1 in the first term, we obtain:
+   *    int[x1,x2] y dx = int[1,x2/x1] a x1**2 z ln(z) dz + y1 (x2**2 - x1**2) / 2
+   *  which simplifies to (using the primitive):
+   *    int[x1,x2] y dx = a / 2 (x2**2 (ln(x2/x1) - 1 / 2) + x1**2 / 2) + y1 (x2**2 - x1**2) / 2
+   */
+  class FirstMomentLinearLogarithmic : public IntegratorBase< FirstMomentLinearLogarithmic > {
+
+    /* friend declarations */
+    friend class IntegratorBase< FirstMomentLinearLogarithmic >;
+
+    /* interface implementation functions */
+
+    /**
+     *  @brief Perform first raw moment integration of a linear-logarithmic panel (y is linear in ln(x))
+     *
+     *  @param[in] xLeft    the left value on the x interval
+     *  @param[in] xRight   the right value on the x interval
+     *  @param[in] yLeft    the left value on the y interval
+     *  @param[in] yRight   the right value on the y interval
+     */
+    template < typename X, typename Y,
+               typename I = decltype( std::declval< X >() * std::declval< X >() * std::declval< Y >() ) >
+    I integrate( const X& xLeft, const X& xRight,
+                 const Y& yLeft, const Y& yRight ) const noexcept {
+
+      return 0.5 * ( ( yRight - yLeft ) / std::log( xRight / xLeft )
+             * ( xRight * xRight * ( std::log( xRight / xLeft ) - 0.5 ) + 0.5 * xLeft * xLeft )
+             + yLeft * ( xRight * xRight - xLeft * xLeft ) );
+    }
+
+  public:
+
+    using IntegratorBase::operator();
+  };
+
+  // integration function
+  static constexpr FirstMomentLinearLogarithmic firstMomentLinLog;
+
+} // integration namespace
+} // scion namespace
+} // njoy namespace
+
+#endif

src/scion/integration/FirstMomentLinearLogarithmic/test/CMakeLists.txt

@@ -0,0 +1 @@
+add_cpp_test( integration.FirstMomentLinearLogarithmic FirstMomentLinearLogarithmic.test.cpp )