trapz

auto eve::trapz = functor<trapz_t>

inlineconstexpr

Header file

#include <eve/module/core.hpp>

Callable Signatures

namespace eve
{
   // Regular overloads
   constexpr auto trapz(eve::non_empty_product_type auto const& x,
                        eve::non_empty_product_type auto const& y)                   noexcept; // 1
   constexpr auto trapz(floating_value auto ... ys)                                  noexcept; // 2
   constexpr auto trapz(eve::non_empty_product_type auto const& y)                   noexcept; // 2
   constexpr auto trapz(floating_value h,
                        eve::non_empty_product_type auto const& y)                   noexcept; // 2
   constexpr auto trapz(eve::invocable f, floating_value auto ... xs)                noexcept; // 3
   constexpr auto trapz(eve::invocable f, eve::non_empty_product_type auto const& x) noexcept; // 3
 
   // Semantic options
   constexpr auto trapz[widen](/*any of the above overloads*/)                       noexcept; // 5
   constexpr auto trapz[kahan](/*any of the above overloads*/)                       noexcept; // 6
}

Parameters

• xs ... : real values ordered arguments in increasing order.
• x, : tuple of real values ordered arguments in increasing order.
• y : tuple of real arguments
• h : floating value scaling the abscissas
• f : invocable

Return value

1. Computes elementwise the integral of the piecewise linear function defined by \(f(x_i) = y_i\) (trapezoidal rule)
2. the missing x parameter is assumed equal to be an arithmetic progression of common difference h (1 if h is omitted).
3. the missing y parameter is assumed equal to the call of f applied to the x elements
4. the computation is made with the upgraded types
5. Internal summation use kahan algorithm for better accuracy

Note

definition of f or of y ?

• If f is a parameter the y values are defined by \(y_i = \mathrm{f}(x_i)\) and x or xsmust be defined
• If y (or ys...) is a parameter f is defined by \(f(x_i) = y_i\). and the values of \(x_i\) are not needed as they are equally spaced with space h (defaulted to one).
• If x (or xs ...) is a parameter its must be sorted in increasing order.

External references

• Wikipedia trapezoidal rule

  Example

// revision 0
#include <eve/module/core.hpp>
#include <iostream>
#include <iomanip>
 
 
int main()
{
 
  {
    kumi::tuple tt{0.01,1.01,2.01,3.01,-1.01,-2.01,-3.01,-4.01};
    kumi::tuple x{1.0,2.0,3.0,4.0,5.0,6.0,7.0,8.0};
    std::cout << " < - x                                                 = " << x << std::endl;
    std::cout << " < - tt                                                = " << tt << std::endl;
    std::cout << " -> trapz(0.01,1.01,2.01,3.01,-1.01,-2.01,-3.01,-4.01) = " << eve::trapz(0.01,1.01,2.01,3.01,-1.01,-2.01,-3.01,-4.01) << std::endl;
    std::cout << " -> trapz(tt)                                          = " << eve::trapz(tt) << std::endl;
    std::cout << " -> trapz(x, tt)                                       = " << eve::trapz(x, tt) << std::endl;
    using w_t = eve::wide<double>;
    auto f = [](auto p){return w_t([p](auto q){return p*p+q*3; }); };
    auto g = [](auto p){return 2*(p+1); };
    w_t h([](auto p){return p+1; });
    kumi::tuple tt1(kumi::generate<8>(f));
    kumi::tuple xx(kumi::generate<8>(g));
    std::cout << " <- xx                                                 = " << xx << std::endl;
    std::cout << " <- tt1                                                = " << tt1 << std::endl;
    std::cout << " <- h                                                  = " << h   << std::endl;
    std::cout << " -> trapz(tt1)                                         = " << eve::trapz(tt1) << std::endl;
    std::cout << " -> trapz(tt1)                                         = " << eve::trapz(tt1) << std::endl;
    std::cout << " -> trapz(2.0, tt1)                                    = " << eve::trapz(2.0, tt1) << std::endl;
    std::cout << " -> trapz(h, tt1)                                      = " << eve::trapz(h, tt1) << std::endl;
    std::cout << " -> trapz(x   tt1)                                     = " << eve::trapz(x, tt1) << std::endl;
    std::cout << " -> trapz[widen](x   tt1)                              = " << eve::trapz[eve::widen](x, tt1) << std::endl;
  }
 
}