Home | History | Annotate | Line # | Download | only in std
numeric revision 1.1.1.9
      1 // <numeric> -*- C++ -*-
      2 
      3 // Copyright (C) 2001-2019 Free Software Foundation, Inc.
      4 //
      5 // This file is part of the GNU ISO C++ Library.  This library is free
      6 // software; you can redistribute it and/or modify it under the
      7 // terms of the GNU General Public License as published by the
      8 // Free Software Foundation; either version 3, or (at your option)
      9 // any later version.
     10 
     11 // This library is distributed in the hope that it will be useful,
     12 // but WITHOUT ANY WARRANTY; without even the implied warranty of
     13 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
     14 // GNU General Public License for more details.
     15 
     16 // Under Section 7 of GPL version 3, you are granted additional
     17 // permissions described in the GCC Runtime Library Exception, version
     18 // 3.1, as published by the Free Software Foundation.
     19 
     20 // You should have received a copy of the GNU General Public License and
     21 // a copy of the GCC Runtime Library Exception along with this program;
     22 // see the files COPYING3 and COPYING.RUNTIME respectively.  If not, see
     23 // <http://www.gnu.org/licenses/>.
     24 
     25 /*
     26  *
     27  * Copyright (c) 1994
     28  * Hewlett-Packard Company
     29  *
     30  * Permission to use, copy, modify, distribute and sell this software
     31  * and its documentation for any purpose is hereby granted without fee,
     32  * provided that the above copyright notice appear in all copies and
     33  * that both that copyright notice and this permission notice appear
     34  * in supporting documentation.  Hewlett-Packard Company makes no
     35  * representations about the suitability of this software for any
     36  * purpose.  It is provided "as is" without express or implied warranty.
     37  *
     38  *
     39  * Copyright (c) 1996,1997
     40  * Silicon Graphics Computer Systems, Inc.
     41  *
     42  * Permission to use, copy, modify, distribute and sell this software
     43  * and its documentation for any purpose is hereby granted without fee,
     44  * provided that the above copyright notice appear in all copies and
     45  * that both that copyright notice and this permission notice appear
     46  * in supporting documentation.  Silicon Graphics makes no
     47  * representations about the suitability of this software for any
     48  * purpose.  It is provided "as is" without express or implied warranty.
     49  */
     50 
     51 /** @file include/numeric
     52  *  This is a Standard C++ Library header.
     53  */
     54 
     55 #ifndef _GLIBCXX_NUMERIC
     56 #define _GLIBCXX_NUMERIC 1
     57 
     58 #pragma GCC system_header
     59 
     60 #include <bits/c++config.h>
     61 #include <bits/stl_iterator_base_types.h>
     62 #include <bits/stl_numeric.h>
     63 
     64 #ifdef _GLIBCXX_PARALLEL
     65 # include <parallel/numeric>
     66 #endif
     67 
     68 /**
     69  * @defgroup numerics Numerics
     70  *
     71  * Components for performing numeric operations. Includes support for
     72  * complex number types, random number generation, numeric (n-at-a-time)
     73  * arrays, generalized numeric algorithms, and mathematical special functions.
     74  */
     75 
     76 #if __cplusplus >= 201402L
     77 #include <type_traits>
     78 
     79 namespace std _GLIBCXX_VISIBILITY(default)
     80 {
     81 _GLIBCXX_BEGIN_NAMESPACE_VERSION
     82 
     83 namespace __detail
     84 {
     85   // std::abs is not constexpr and doesn't support unsigned integers.
     86   template<typename _Tp>
     87     constexpr
     88     enable_if_t<__and_<is_integral<_Tp>, is_signed<_Tp>>::value, _Tp>
     89     __abs_integral(_Tp __val)
     90     { return __val < 0 ? -__val : __val; }
     91 
     92   template<typename _Tp>
     93     constexpr
     94     enable_if_t<__and_<is_integral<_Tp>, is_unsigned<_Tp>>::value, _Tp>
     95     __abs_integral(_Tp __val)
     96     { return __val; }
     97 
     98   void __abs_integral(bool) = delete;
     99 
    100   template<typename _Mn, typename _Nn>
    101     constexpr common_type_t<_Mn, _Nn>
    102     __gcd(_Mn __m, _Nn __n)
    103     {
    104       return __m == 0 ? __detail::__abs_integral(__n)
    105 	: __n == 0 ? __detail::__abs_integral(__m)
    106 	: __detail::__gcd(__n, __m % __n);
    107     }
    108 
    109   /// Least common multiple
    110   template<typename _Mn, typename _Nn>
    111     constexpr common_type_t<_Mn, _Nn>
    112     __lcm(_Mn __m, _Nn __n)
    113     {
    114       return (__m != 0 && __n != 0)
    115 	? (__detail::__abs_integral(__m) / __detail::__gcd(__m, __n))
    116 	  * __detail::__abs_integral(__n)
    117 	: 0;
    118     }
    119 } // namespace __detail
    120 
    121 #if __cplusplus >= 201703L
    122 
    123 #define __cpp_lib_gcd_lcm 201606
    124 // These were used in drafts of SD-6:
    125 #define __cpp_lib_gcd 201606
    126 #define __cpp_lib_lcm 201606
    127 
    128   /// Greatest common divisor
    129   template<typename _Mn, typename _Nn>
    130     constexpr common_type_t<_Mn, _Nn>
    131     gcd(_Mn __m, _Nn __n)
    132     {
    133       static_assert(is_integral_v<_Mn>, "gcd arguments are integers");
    134       static_assert(is_integral_v<_Nn>, "gcd arguments are integers");
    135       static_assert(!is_same_v<remove_cv_t<_Mn>, bool>,
    136 		    "gcd arguments are not bools");
    137       static_assert(!is_same_v<remove_cv_t<_Nn>, bool>,
    138 		    "gcd arguments are not bools");
    139       return __detail::__gcd(__m, __n);
    140     }
    141 
    142   /// Least common multiple
    143   template<typename _Mn, typename _Nn>
    144     constexpr common_type_t<_Mn, _Nn>
    145     lcm(_Mn __m, _Nn __n)
    146     {
    147       static_assert(is_integral_v<_Mn>, "lcm arguments are integers");
    148       static_assert(is_integral_v<_Nn>, "lcm arguments are integers");
    149       static_assert(!is_same_v<remove_cv_t<_Mn>, bool>,
    150 		    "lcm arguments are not bools");
    151       static_assert(!is_same_v<remove_cv_t<_Nn>, bool>,
    152 		    "lcm arguments are not bools");
    153       return __detail::__lcm(__m, __n);
    154     }
    155 
    156 #endif // C++17
    157 
    158 _GLIBCXX_END_NAMESPACE_VERSION
    159 } // namespace std
    160 
    161 #endif // C++14
    162 
    163 #if __cplusplus > 201703L
    164 #include <limits>
    165 
    166 namespace std _GLIBCXX_VISIBILITY(default)
    167 {
    168 _GLIBCXX_BEGIN_NAMESPACE_VERSION
    169   // midpoint
    170 # define __cpp_lib_interpolate 201902L
    171 
    172   template<typename _Tp>
    173     constexpr
    174     enable_if_t<__and_v<is_arithmetic<_Tp>, is_same<remove_cv_t<_Tp>, _Tp>,
    175 			__not_<is_same<_Tp, bool>>>,
    176 		_Tp>
    177     midpoint(_Tp __a, _Tp __b) noexcept
    178     {
    179       if constexpr (is_integral_v<_Tp>)
    180 	{
    181 	  using _Up = make_unsigned_t<_Tp>;
    182 
    183 	  int __k = 1;
    184 	  _Up __m = __a;
    185 	  _Up __M = __b;
    186 	  if (__a > __b)
    187 	    {
    188 	      __k = -1;
    189 	      __m = __b;
    190 	      __M = __a;
    191 	    }
    192 	  return __a + __k * _Tp(_Up(__M - __m) / 2);
    193 	}
    194       else // is_floating
    195 	{
    196 	  constexpr _Tp __lo = numeric_limits<_Tp>::min() * 2;
    197 	  constexpr _Tp __hi = numeric_limits<_Tp>::max() / 2;
    198 	  const _Tp __abs_a = __a < 0 ? -__a : __a;
    199 	  const _Tp __abs_b = __b < 0 ? -__b : __b;
    200 	  if (__abs_a <= __hi && __abs_b <= __hi) [[likely]]
    201 	    return (__a + __b) / 2; // always correctly rounded
    202 	  if (__abs_a < __lo) // not safe to halve __a
    203 	    return __a + __b/2;
    204 	  if (__abs_b < __lo) // not safe to halve __b
    205 	    return __a/2 + __b;
    206 	  return __a/2 + __b/2;	    // otherwise correctly rounded
    207 	}
    208     }
    209 
    210   template<typename _Tp>
    211     constexpr
    212     enable_if_t<__and_v<is_object<_Tp>, bool_constant<sizeof(_Tp) != 0>>, _Tp*>
    213     midpoint(_Tp* __a, _Tp* __b) noexcept
    214     {
    215       return __a  + (__b - __a) / 2;
    216     }
    217 _GLIBCXX_END_NAMESPACE_VERSION
    218 } // namespace std
    219 
    220 #endif // C++20
    221 
    222 #if __cplusplus > 201402L
    223 #include <bits/stl_function.h>
    224 
    225 namespace std _GLIBCXX_VISIBILITY(default)
    226 {
    227 _GLIBCXX_BEGIN_NAMESPACE_VERSION
    228 
    229   /// @addtogroup numeric_ops
    230   /// @{
    231 
    232   /// @cond undocumented
    233   template<typename _It, typename _Traits = iterator_traits<_It>,
    234 	   typename _Cat = typename _Traits::iterator_category>
    235     using __is_random_access_iter
    236       = is_base_of<random_access_iterator_tag, _Cat>;
    237   /// @endcond
    238 
    239   /**
    240    *  @brief  Calculate reduction of values in a range.
    241    *
    242    *  @param  __first  Start of range.
    243    *  @param  __last  End of range.
    244    *  @param  __init  Starting value to add other values to.
    245    *  @param  __binary_op A binary function object.
    246    *  @return  The final sum.
    247    *
    248    *  Reduce the values in the range `[first,last)` using a binary operation.
    249    *  The initial value is `init`.  The values are not necessarily processed
    250    *  in order.
    251    *
    252    *  This algorithm is similar to `std::accumulate` but is not required to
    253    *  perform the operations in order from first to last. For operations
    254    *  that are commutative and associative the result will be the same as
    255    *  for `std::accumulate`, but for other operations (such as floating point
    256    *  arithmetic) the result can be different.
    257    */
    258   template<typename _InputIterator, typename _Tp, typename _BinaryOperation>
    259     _Tp
    260     reduce(_InputIterator __first, _InputIterator __last, _Tp __init,
    261 	   _BinaryOperation __binary_op)
    262     {
    263       using value_type = typename iterator_traits<_InputIterator>::value_type;
    264       static_assert(is_invocable_r_v<_Tp, _BinaryOperation&, _Tp&, _Tp&>);
    265       static_assert(is_convertible_v<value_type, _Tp>);
    266       if constexpr (__is_random_access_iter<_InputIterator>::value)
    267 	{
    268 	  while ((__last - __first) >= 4)
    269 	    {
    270 	      _Tp __v1 = __binary_op(__first[0], __first[1]);
    271 	      _Tp __v2 = __binary_op(__first[2], __first[3]);
    272 	      _Tp __v3 = __binary_op(__v1, __v2);
    273 	      __init = __binary_op(__init, __v3);
    274 	      __first += 4;
    275 	    }
    276 	}
    277       for (; __first != __last; ++__first)
    278 	__init = __binary_op(__init, *__first);
    279       return __init;
    280     }
    281 
    282  /**
    283    *  @brief  Calculate reduction of values in a range.
    284    *
    285    *  @param  __first  Start of range.
    286    *  @param  __last  End of range.
    287    *  @param  __init  Starting value to add other values to.
    288    *  @return  The final sum.
    289    *
    290    *  Reduce the values in the range `[first,last)` using addition.
    291    *  Equivalent to calling `std::reduce(first, last, init, std::plus<>())`.
    292    */
    293   template<typename _InputIterator, typename _Tp>
    294     inline _Tp
    295     reduce(_InputIterator __first, _InputIterator __last, _Tp __init)
    296     { return std::reduce(__first, __last, std::move(__init), plus<>()); }
    297 
    298   /**
    299    *  @brief  Calculate reduction of values in a range.
    300    *
    301    *  @param  __first  Start of range.
    302    *  @param  __last  End of range.
    303    *  @return  The final sum.
    304    *
    305    *  Reduce the values in the range `[first,last)` using addition, with
    306    *  an initial value of `T{}`, where `T` is the iterator's value type.
    307    *  Equivalent to calling `std::reduce(first, last, T{}, std::plus<>())`.
    308    */
    309   template<typename _InputIterator>
    310     inline typename iterator_traits<_InputIterator>::value_type
    311     reduce(_InputIterator __first, _InputIterator __last)
    312     {
    313       using value_type = typename iterator_traits<_InputIterator>::value_type;
    314       return std::reduce(__first, __last, value_type{}, plus<>());
    315     }
    316 
    317   /**
    318    *  @brief  Combine elements from two ranges and reduce
    319    *
    320    *  @param  __first1  Start of first range.
    321    *  @param  __last1  End of first range.
    322    *  @param  __first2  Start of second range.
    323    *  @param  __init  Starting value to add other values to.
    324    *  @param  __binary_op1 The function used to perform reduction.
    325    *  @param  __binary_op2 The function used to combine values from the ranges.
    326    *  @return  The final sum.
    327    *
    328    *  Call `binary_op2(first1[n],first2[n])` for each `n` in `[0,last1-first1)`
    329    *  and then use `binary_op1` to reduce the values returned by `binary_op2`
    330    *  to a single value of type `T`.
    331    *
    332    *  The range beginning at `first2` must contain at least `last1-first1`
    333    *  elements.
    334    */
    335   template<typename _InputIterator1, typename _InputIterator2, typename _Tp,
    336 	   typename _BinaryOperation1, typename _BinaryOperation2>
    337     _Tp
    338     transform_reduce(_InputIterator1 __first1, _InputIterator1 __last1,
    339 		     _InputIterator2 __first2, _Tp __init,
    340 		     _BinaryOperation1 __binary_op1,
    341 		     _BinaryOperation2 __binary_op2)
    342     {
    343       if constexpr (__and_v<__is_random_access_iter<_InputIterator1>,
    344 			    __is_random_access_iter<_InputIterator2>>)
    345 	{
    346 	  while ((__last1 - __first1) >= 4)
    347 	    {
    348 	      _Tp __v1 = __binary_op1(__binary_op2(__first1[0], __first2[0]),
    349 				      __binary_op2(__first1[1], __first2[1]));
    350 	      _Tp __v2 = __binary_op1(__binary_op2(__first1[2], __first2[2]),
    351 				      __binary_op2(__first1[3], __first2[3]));
    352 	      _Tp __v3 = __binary_op1(__v1, __v2);
    353 	      __init = __binary_op1(__init, __v3);
    354 	      __first1 += 4;
    355 	      __first2 += 4;
    356 	    }
    357 	}
    358       for (; __first1 != __last1; ++__first1, (void) ++__first2)
    359 	__init = __binary_op1(__init, __binary_op2(*__first1, *__first2));
    360       return __init;
    361     }
    362 
    363   /**
    364    *  @brief  Combine elements from two ranges and reduce
    365    *
    366    *  @param  __first1  Start of first range.
    367    *  @param  __last1  End of first range.
    368    *  @param  __first2  Start of second range.
    369    *  @param  __init  Starting value to add other values to.
    370    *  @return  The final sum.
    371    *
    372    *  Call `first1[n]*first2[n]` for each `n` in `[0,last1-first1)` and then
    373    *  use addition to sum those products to a single value of type `T`.
    374    *
    375    *  The range beginning at `first2` must contain at least `last1-first1`
    376    *  elements.
    377    */
    378   template<typename _InputIterator1, typename _InputIterator2, typename _Tp>
    379     inline _Tp
    380     transform_reduce(_InputIterator1 __first1, _InputIterator1 __last1,
    381 		     _InputIterator2 __first2, _Tp __init)
    382     {
    383       return std::transform_reduce(__first1, __last1, __first2,
    384 				   std::move(__init),
    385 				   plus<>(), multiplies<>());
    386     }
    387 
    388   /**
    389    *  @brief  Transform the elements of a range and reduce
    390    *
    391    *  @param  __first  Start of range.
    392    *  @param  __last  End of range.
    393    *  @param  __init  Starting value to add other values to.
    394    *  @param  __binary_op The function used to perform reduction.
    395    *  @param  __unary_op The function used to transform values from the range.
    396    *  @return  The final sum.
    397    *
    398    *  Call `unary_op(first[n])` for each `n` in `[0,last-first)` and then
    399    *  use `binary_op` to reduce the values returned by `unary_op`
    400    *  to a single value of type `T`.
    401    */
    402   template<typename _InputIterator, typename _Tp,
    403 	   typename _BinaryOperation, typename _UnaryOperation>
    404     _Tp
    405     transform_reduce(_InputIterator __first, _InputIterator __last, _Tp __init,
    406 		     _BinaryOperation __binary_op, _UnaryOperation __unary_op)
    407     {
    408       if constexpr (__is_random_access_iter<_InputIterator>::value)
    409 	{
    410 	  while ((__last - __first) >= 4)
    411 	    {
    412 	      _Tp __v1 = __binary_op(__unary_op(__first[0]),
    413 				     __unary_op(__first[1]));
    414 	      _Tp __v2 = __binary_op(__unary_op(__first[2]),
    415 				     __unary_op(__first[3]));
    416 	      _Tp __v3 = __binary_op(__v1, __v2);
    417 	      __init = __binary_op(__init, __v3);
    418 	      __first += 4;
    419 	    }
    420 	}
    421       for (; __first != __last; ++__first)
    422 	__init = __binary_op(__init, __unary_op(*__first));
    423       return __init;
    424     }
    425 
    426   /** @brief Output the cumulative sum of one range to a second range
    427    *
    428    *  @param __first  Start of input range.
    429    *  @param __last   End of input range.
    430    *  @param __result Start of output range.
    431    *  @param __init   Initial value.
    432    *  @param __binary_op Function to perform summation.
    433    *  @return The end of the output range.
    434    *
    435    *  Write the cumulative sum (aka prefix sum, aka scan) of the input range
    436    *  to the output range. Each element of the output range contains the
    437    *  running total of all earlier elements (and the initial value),
    438    *  using `binary_op` for summation.
    439    *
    440    *  This function generates an "exclusive" scan, meaning the Nth element
    441    *  of the output range is the sum of the first N-1 input elements,
    442    *  so the Nth input element is not included.
    443    */
    444   template<typename _InputIterator, typename _OutputIterator, typename _Tp,
    445 	   typename _BinaryOperation>
    446     _OutputIterator
    447     exclusive_scan(_InputIterator __first, _InputIterator __last,
    448 		   _OutputIterator __result, _Tp __init,
    449 		   _BinaryOperation __binary_op)
    450     {
    451       while (__first != __last)
    452 	{
    453 	  auto __v = __init;
    454 	  __init = __binary_op(__init, *__first);
    455 	  ++__first;
    456 	  *__result++ = std::move(__v);
    457 	}
    458       return __result;
    459     }
    460 
    461   /** @brief Output the cumulative sum of one range to a second range
    462    *
    463    *  @param __first  Start of input range.
    464    *  @param __last   End of input range.
    465    *  @param __result Start of output range.
    466    *  @param __init   Initial value.
    467    *  @return The end of the output range.
    468    *
    469    *  Write the cumulative sum (aka prefix sum, aka scan) of the input range
    470    *  to the output range. Each element of the output range contains the
    471    *  running total of all earlier elements (and the initial value),
    472    *  using `std::plus<>` for summation.
    473    *
    474    *  This function generates an "exclusive" scan, meaning the Nth element
    475    *  of the output range is the sum of the first N-1 input elements,
    476    *  so the Nth input element is not included.
    477    */
    478   template<typename _InputIterator, typename _OutputIterator, typename _Tp>
    479     inline _OutputIterator
    480     exclusive_scan(_InputIterator __first, _InputIterator __last,
    481 		   _OutputIterator __result, _Tp __init)
    482     {
    483       return std::exclusive_scan(__first, __last, __result, std::move(__init),
    484 				 plus<>());
    485     }
    486 
    487   /** @brief Output the cumulative sum of one range to a second range
    488    *
    489    *  @param __first  Start of input range.
    490    *  @param __last   End of input range.
    491    *  @param __result Start of output range.
    492    *  @param __binary_op Function to perform summation.
    493    *  @param __init   Initial value.
    494    *  @return The end of the output range.
    495    *
    496    *  Write the cumulative sum (aka prefix sum, aka scan) of the input range
    497    *  to the output range. Each element of the output range contains the
    498    *  running total of all earlier elements (and the initial value),
    499    *  using `binary_op` for summation.
    500    *
    501    *  This function generates an "inclusive" scan, meaning the Nth element
    502    *  of the output range is the sum of the first N input elements,
    503    *  so the Nth input element is included.
    504    */
    505   template<typename _InputIterator, typename _OutputIterator,
    506 	   typename _BinaryOperation, typename _Tp>
    507     _OutputIterator
    508     inclusive_scan(_InputIterator __first, _InputIterator __last,
    509 		   _OutputIterator __result, _BinaryOperation __binary_op,
    510 		   _Tp __init)
    511     {
    512       for (; __first != __last; ++__first)
    513 	*__result++ = __init = __binary_op(__init, *__first);
    514       return __result;
    515     }
    516 
    517   /** @brief Output the cumulative sum of one range to a second range
    518    *
    519    *  @param __first  Start of input range.
    520    *  @param __last   End of input range.
    521    *  @param __result Start of output range.
    522    *  @param __binary_op Function to perform summation.
    523    *  @return The end of the output range.
    524    *
    525    *  Write the cumulative sum (aka prefix sum, aka scan) of the input range
    526    *  to the output range. Each element of the output range contains the
    527    *  running total of all earlier elements, using `binary_op` for summation.
    528    *
    529    *  This function generates an "inclusive" scan, meaning the Nth element
    530    *  of the output range is the sum of the first N input elements,
    531    *  so the Nth input element is included.
    532    */
    533   template<typename _InputIterator, typename _OutputIterator,
    534 	   typename _BinaryOperation>
    535     _OutputIterator
    536     inclusive_scan(_InputIterator __first, _InputIterator __last,
    537 		   _OutputIterator __result, _BinaryOperation __binary_op)
    538     {
    539       if (__first != __last)
    540 	{
    541 	  auto __init = *__first;
    542 	  *__result++ = __init;
    543 	  ++__first;
    544 	  if (__first != __last)
    545 	    __result = std::inclusive_scan(__first, __last, __result,
    546 					   __binary_op, std::move(__init));
    547 	}
    548       return __result;
    549     }
    550 
    551   /** @brief Output the cumulative sum of one range to a second range
    552    *
    553    *  @param __first  Start of input range.
    554    *  @param __last   End of input range.
    555    *  @param __result Start of output range.
    556    *  @return The end of the output range.
    557    *
    558    *  Write the cumulative sum (aka prefix sum, aka scan) of the input range
    559    *  to the output range. Each element of the output range contains the
    560    *  running total of all earlier elements, using `std::plus<>` for summation.
    561    *
    562    *  This function generates an "inclusive" scan, meaning the Nth element
    563    *  of the output range is the sum of the first N input elements,
    564    *  so the Nth input element is included.
    565    */
    566   template<typename _InputIterator, typename _OutputIterator>
    567     inline _OutputIterator
    568     inclusive_scan(_InputIterator __first, _InputIterator __last,
    569 		   _OutputIterator __result)
    570     { return std::inclusive_scan(__first, __last, __result, plus<>()); }
    571 
    572   /** @brief Output the cumulative sum of one range to a second range
    573    *
    574    *  @param __first  Start of input range.
    575    *  @param __last   End of input range.
    576    *  @param __result Start of output range.
    577    *  @param __init   Initial value.
    578    *  @param __binary_op Function to perform summation.
    579    *  @param __unary_op Function to transform elements of the input range.
    580    *  @return The end of the output range.
    581    *
    582    *  Write the cumulative sum (aka prefix sum, aka scan) of the input range
    583    *  to the output range. Each element of the output range contains the
    584    *  running total of all earlier elements (and the initial value),
    585    *  using `__unary_op` to transform the input elements
    586    *  and using `__binary_op` for summation.
    587    *
    588    *  This function generates an "exclusive" scan, meaning the Nth element
    589    *  of the output range is the sum of the first N-1 input elements,
    590    *  so the Nth input element is not included.
    591    */
    592   template<typename _InputIterator, typename _OutputIterator, typename _Tp,
    593 	   typename _BinaryOperation, typename _UnaryOperation>
    594     _OutputIterator
    595     transform_exclusive_scan(_InputIterator __first, _InputIterator __last,
    596 			     _OutputIterator __result, _Tp __init,
    597 			     _BinaryOperation __binary_op,
    598 			     _UnaryOperation __unary_op)
    599     {
    600       while (__first != __last)
    601 	{
    602 	  auto __v = __init;
    603 	  __init = __binary_op(__init, __unary_op(*__first));
    604 	  ++__first;
    605 	  *__result++ = std::move(__v);
    606 	}
    607       return __result;
    608     }
    609 
    610   /** @brief Output the cumulative sum of one range to a second range
    611    *
    612    *  @param __first  Start of input range.
    613    *  @param __last   End of input range.
    614    *  @param __result Start of output range.
    615    *  @param __binary_op Function to perform summation.
    616    *  @param __unary_op Function to transform elements of the input range.
    617    *  @param __init   Initial value.
    618    *  @return The end of the output range.
    619    *
    620    *  Write the cumulative sum (aka prefix sum, aka scan) of the input range
    621    *  to the output range. Each element of the output range contains the
    622    *  running total of all earlier elements (and the initial value),
    623    *  using `__unary_op` to transform the input elements
    624    *  and using `__binary_op` for summation.
    625    *
    626    *  This function generates an "inclusive" scan, meaning the Nth element
    627    *  of the output range is the sum of the first N input elements,
    628    *  so the Nth input element is included.
    629    */
    630   template<typename _InputIterator, typename _OutputIterator,
    631 	   typename _BinaryOperation, typename _UnaryOperation, typename _Tp>
    632     _OutputIterator
    633     transform_inclusive_scan(_InputIterator __first, _InputIterator __last,
    634 			     _OutputIterator __result,
    635 			     _BinaryOperation __binary_op,
    636 			     _UnaryOperation __unary_op,
    637 			     _Tp __init)
    638     {
    639       for (; __first != __last; ++__first)
    640 	*__result++ = __init = __binary_op(__init, __unary_op(*__first));
    641       return __result;
    642     }
    643 
    644   /** @brief Output the cumulative sum of one range to a second range
    645    *
    646    *  @param __first  Start of input range.
    647    *  @param __last   End of input range.
    648    *  @param __result Start of output range.
    649    *  @param __binary_op Function to perform summation.
    650    *  @param __unary_op Function to transform elements of the input range.
    651    *  @return The end of the output range.
    652    *
    653    *  Write the cumulative sum (aka prefix sum, aka scan) of the input range
    654    *  to the output range. Each element of the output range contains the
    655    *  running total of all earlier elements,
    656    *  using `__unary_op` to transform the input elements
    657    *  and using `__binary_op` for summation.
    658    *
    659    *  This function generates an "inclusive" scan, meaning the Nth element
    660    *  of the output range is the sum of the first N input elements,
    661    *  so the Nth input element is included.
    662    */
    663   template<typename _InputIterator, typename _OutputIterator,
    664 	  typename _BinaryOperation, typename _UnaryOperation>
    665     _OutputIterator
    666     transform_inclusive_scan(_InputIterator __first, _InputIterator __last,
    667 			     _OutputIterator __result,
    668 			     _BinaryOperation __binary_op,
    669 			     _UnaryOperation __unary_op)
    670     {
    671       if (__first != __last)
    672 	{
    673 	  auto __init = __unary_op(*__first);
    674 	  *__result++ = __init;
    675 	  ++__first;
    676 	  if (__first != __last)
    677 	    __result = std::transform_inclusive_scan(__first, __last, __result,
    678 						     __binary_op, __unary_op,
    679 						     std::move(__init));
    680 	}
    681       return __result;
    682     }
    683 
    684   // @} group numeric_ops
    685 
    686 _GLIBCXX_END_NAMESPACE_VERSION
    687 } // namespace std
    688 
    689 // Parallel STL algorithms
    690 # if __PSTL_EXECUTION_POLICIES_DEFINED
    691 // If <execution> has already been included, pull in implementations
    692 #  include <pstl/glue_numeric_impl.h>
    693 # else
    694 // Otherwise just pull in forward declarations
    695 #  include <pstl/glue_numeric_defs.h>
    696 #  define __PSTL_NUMERIC_FORWARD_DECLARED 1
    697 # endif
    698 
    699 // Feature test macro for parallel algorithms
    700 # define __cpp_lib_parallel_algorithm 201603L
    701 #endif // C++17
    702 
    703 #endif /* _GLIBCXX_NUMERIC */
    704