libstdc++
numeric
Go to the documentation of this file.
1 // <numeric> -*- C++ -*-
2 
3 // Copyright (C) 2001-2020 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 #include <ext/numeric_traits.h>
64 
65 #ifdef _GLIBCXX_PARALLEL
66 # include <parallel/numeric>
67 #endif
68 
69 /**
70  * @defgroup numerics Numerics
71  *
72  * Components for performing numeric operations. Includes support for
73  * complex number types, random number generation, numeric (n-at-a-time)
74  * arrays, generalized numeric algorithms, and mathematical special functions.
75  */
76 
77 #if __cplusplus >= 201402L
78 #include <type_traits>
79 #include <ext/numeric_traits.h>
80 
81 namespace std _GLIBCXX_VISIBILITY(default)
82 {
83 _GLIBCXX_BEGIN_NAMESPACE_VERSION
84 
85 namespace __detail
86 {
87  // Like std::abs, but supports unsigned types and returns the specified type,
88  // so |std::numeric_limits<_Tp>::min()| is OK if representable in _Res.
89  template<typename _Res, typename _Tp>
90  constexpr _Res
91  __abs_r(_Tp __val)
92  {
93  static_assert(sizeof(_Res) >= sizeof(_Tp),
94  "result type must be at least as wide as the input type");
95 
96  if (__val >= 0)
97  return __val;
98 #if defined _GLIBCXX_ASSERTIONS && defined _GLIBCXX_HAVE_BUILTIN_IS_CONSTANT_EVALUATED
99  if (!__builtin_is_constant_evaluated()) // overflow already detected in constexpr
100  __glibcxx_assert(__val != __gnu_cxx::__int_traits<_Res>::__min);
101 #endif
102  return -static_cast<_Res>(__val);
103  }
104 
105  template<typename> void __abs_r(bool) = delete;
106 
107  // GCD implementation
108  template<typename _Tp>
109  constexpr _Tp
110  __gcd(_Tp __m, _Tp __n)
111  {
112  static_assert(is_unsigned<_Tp>::value, "type must be unsigned");
113  return __m == 0 ? __n
114  : __n == 0 ? __m
115  : __detail::__gcd(__n, _Tp(__m % __n));
116  }
117 } // namespace __detail
118 
119 #if __cplusplus >= 201703L
120 
121 #define __cpp_lib_gcd_lcm 201606
122 // These were used in drafts of SD-6:
123 #define __cpp_lib_gcd 201606
124 #define __cpp_lib_lcm 201606
125 
126  /// Greatest common divisor
127  template<typename _Mn, typename _Nn>
128  constexpr common_type_t<_Mn, _Nn>
129  gcd(_Mn __m, _Nn __n) noexcept
130  {
131  static_assert(is_integral_v<_Mn> && is_integral_v<_Nn>,
132  "std::gcd arguments must be integers");
133  static_assert(_Mn(2) == 2 && _Nn(2) == 2,
134  "std::gcd arguments must not be bool");
135  using _Ct = common_type_t<_Mn, _Nn>;
136  const _Ct __m2 = __detail::__abs_r<_Ct>(__m);
137  const _Ct __n2 = __detail::__abs_r<_Ct>(__n);
138  return __detail::__gcd<make_unsigned_t<_Ct>>(__m2, __n2);
139  }
140 
141  /// Least common multiple
142  template<typename _Mn, typename _Nn>
143  constexpr common_type_t<_Mn, _Nn>
144  lcm(_Mn __m, _Nn __n) noexcept
145  {
146  static_assert(is_integral_v<_Mn> && is_integral_v<_Nn>,
147  "std::lcm arguments must be integers");
148  static_assert(_Mn(2) == 2 && _Nn(2) == 2,
149  "std::lcm arguments must not be bool");
150  using _Ct = common_type_t<_Mn, _Nn>;
151  const _Ct __m2 = __detail::__abs_r<_Ct>(__m);
152  const _Ct __n2 = __detail::__abs_r<_Ct>(__n);
153  if (__m2 == 0 || __n2 == 0)
154  return 0;
155  _Ct __r = __m2 / __detail::__gcd<make_unsigned_t<_Ct>>(__m2, __n2);
156 
157 #if defined _GLIBCXX_HAVE_BUILTIN_IS_CONSTANT_EVALUATED
158  if constexpr (is_signed_v<_Ct>)
159  if (__builtin_is_constant_evaluated())
160  return __r * __n2; // constant evaluation can detect overflow here.
161 #endif
162 
163  bool __overflow = __builtin_mul_overflow(__r, __n2, &__r);
164 #if defined _GLIBCXX_HAVE_BUILTIN_IS_CONSTANT_EVALUATED
165  if (__builtin_is_constant_evaluated())
166  if (__overflow)
167  _GLIBCXX_THROW_OR_ABORT("std::lcm result is out of range of type");
168 #endif
169  __glibcxx_assert(!__overflow);
170  return __r;
171  }
172 
173 #endif // C++17
174 
175 _GLIBCXX_END_NAMESPACE_VERSION
176 } // namespace std
177 
178 #endif // C++14
179 
180 #if __cplusplus > 201703L
181 #include <limits>
182 
183 namespace std _GLIBCXX_VISIBILITY(default)
184 {
185 _GLIBCXX_BEGIN_NAMESPACE_VERSION
186  // midpoint
187 # define __cpp_lib_interpolate 201902L
188 
189  template<typename _Tp>
190  constexpr
191  enable_if_t<__and_v<is_arithmetic<_Tp>, is_same<remove_cv_t<_Tp>, _Tp>,
192  __not_<is_same<_Tp, bool>>>,
193  _Tp>
194  midpoint(_Tp __a, _Tp __b) noexcept
195  {
196  if constexpr (is_integral_v<_Tp>)
197  {
198  using _Up = make_unsigned_t<_Tp>;
199 
200  int __k = 1;
201  _Up __m = __a;
202  _Up __M = __b;
203  if (__a > __b)
204  {
205  __k = -1;
206  __m = __b;
207  __M = __a;
208  }
209  return __a + __k * _Tp(_Up(__M - __m) / 2);
210  }
211  else // is_floating
212  {
213  constexpr _Tp __lo = numeric_limits<_Tp>::min() * 2;
214  constexpr _Tp __hi = numeric_limits<_Tp>::max() / 2;
215  const _Tp __abs_a = __a < 0 ? -__a : __a;
216  const _Tp __abs_b = __b < 0 ? -__b : __b;
217  if (__abs_a <= __hi && __abs_b <= __hi) [[likely]]
218  return (__a + __b) / 2; // always correctly rounded
219  if (__abs_a < __lo) // not safe to halve __a
220  return __a + __b/2;
221  if (__abs_b < __lo) // not safe to halve __b
222  return __a/2 + __b;
223  return __a/2 + __b/2; // otherwise correctly rounded
224  }
225  }
226 
227  template<typename _Tp>
228  constexpr enable_if_t<is_object_v<_Tp>, _Tp*>
229  midpoint(_Tp* __a, _Tp* __b) noexcept
230  {
231  static_assert( sizeof(_Tp) != 0, "type must be complete" );
232  return __a + (__b - __a) / 2;
233  }
234 _GLIBCXX_END_NAMESPACE_VERSION
235 } // namespace std
236 
237 #endif // C++20
238 
239 #if __cplusplus > 201402L
240 #include <bits/stl_function.h>
241 
242 namespace std _GLIBCXX_VISIBILITY(default)
243 {
244 _GLIBCXX_BEGIN_NAMESPACE_VERSION
245 
246 #if __cplusplus > 201703L
247 #define __cpp_lib_constexpr_numeric 201911L
248 #endif
249 
250  /// @addtogroup numeric_ops
251  /// @{
252 
253  /**
254  * @brief Calculate reduction of values in a range.
255  *
256  * @param __first Start of range.
257  * @param __last End of range.
258  * @param __init Starting value to add other values to.
259  * @param __binary_op A binary function object.
260  * @return The final sum.
261  *
262  * Reduce the values in the range `[first,last)` using a binary operation.
263  * The initial value is `init`. The values are not necessarily processed
264  * in order.
265  *
266  * This algorithm is similar to `std::accumulate` but is not required to
267  * perform the operations in order from first to last. For operations
268  * that are commutative and associative the result will be the same as
269  * for `std::accumulate`, but for other operations (such as floating point
270  * arithmetic) the result can be different.
271  */
272  template<typename _InputIterator, typename _Tp, typename _BinaryOperation>
273  _GLIBCXX20_CONSTEXPR
274  _Tp
275  reduce(_InputIterator __first, _InputIterator __last, _Tp __init,
276  _BinaryOperation __binary_op)
277  {
278  using __ref = typename iterator_traits<_InputIterator>::reference;
279  static_assert(is_invocable_r_v<_Tp, _BinaryOperation&, _Tp&, __ref>);
280  static_assert(is_invocable_r_v<_Tp, _BinaryOperation&, __ref, _Tp&>);
281  static_assert(is_invocable_r_v<_Tp, _BinaryOperation&, _Tp&, _Tp&>);
282  static_assert(is_invocable_r_v<_Tp, _BinaryOperation&, __ref, __ref>);
283  if constexpr (__is_random_access_iter<_InputIterator>::value)
284  {
285  while ((__last - __first) >= 4)
286  {
287  _Tp __v1 = __binary_op(__first[0], __first[1]);
288  _Tp __v2 = __binary_op(__first[2], __first[3]);
289  _Tp __v3 = __binary_op(__v1, __v2);
290  __init = __binary_op(__init, __v3);
291  __first += 4;
292  }
293  }
294  for (; __first != __last; ++__first)
295  __init = __binary_op(__init, *__first);
296  return __init;
297  }
298 
299  /**
300  * @brief Calculate reduction of values in a range.
301  *
302  * @param __first Start of range.
303  * @param __last End of range.
304  * @param __init Starting value to add other values to.
305  * @return The final sum.
306  *
307  * Reduce the values in the range `[first,last)` using addition.
308  * Equivalent to calling `std::reduce(first, last, init, std::plus<>())`.
309  */
310  template<typename _InputIterator, typename _Tp>
311  _GLIBCXX20_CONSTEXPR
312  inline _Tp
313  reduce(_InputIterator __first, _InputIterator __last, _Tp __init)
314  { return std::reduce(__first, __last, std::move(__init), plus<>()); }
315 
316  /**
317  * @brief Calculate reduction of values in a range.
318  *
319  * @param __first Start of range.
320  * @param __last End of range.
321  * @return The final sum.
322  *
323  * Reduce the values in the range `[first,last)` using addition, with
324  * an initial value of `T{}`, where `T` is the iterator's value type.
325  * Equivalent to calling `std::reduce(first, last, T{}, std::plus<>())`.
326  */
327  template<typename _InputIterator>
328  _GLIBCXX20_CONSTEXPR
329  inline typename iterator_traits<_InputIterator>::value_type
330  reduce(_InputIterator __first, _InputIterator __last)
331  {
332  using value_type = typename iterator_traits<_InputIterator>::value_type;
333  return std::reduce(__first, __last, value_type{}, plus<>());
334  }
335 
336  /**
337  * @brief Combine elements from two ranges and reduce
338  *
339  * @param __first1 Start of first range.
340  * @param __last1 End of first range.
341  * @param __first2 Start of second range.
342  * @param __init Starting value to add other values to.
343  * @param __binary_op1 The function used to perform reduction.
344  * @param __binary_op2 The function used to combine values from the ranges.
345  * @return The final sum.
346  *
347  * Call `binary_op2(first1[n],first2[n])` for each `n` in `[0,last1-first1)`
348  * and then use `binary_op1` to reduce the values returned by `binary_op2`
349  * to a single value of type `T`.
350  *
351  * The range beginning at `first2` must contain at least `last1-first1`
352  * elements.
353  */
354  template<typename _InputIterator1, typename _InputIterator2, typename _Tp,
355  typename _BinaryOperation1, typename _BinaryOperation2>
356  _GLIBCXX20_CONSTEXPR
357  _Tp
358  transform_reduce(_InputIterator1 __first1, _InputIterator1 __last1,
359  _InputIterator2 __first2, _Tp __init,
360  _BinaryOperation1 __binary_op1,
361  _BinaryOperation2 __binary_op2)
362  {
363  if constexpr (__and_v<__is_random_access_iter<_InputIterator1>,
364  __is_random_access_iter<_InputIterator2>>)
365  {
366  while ((__last1 - __first1) >= 4)
367  {
368  _Tp __v1 = __binary_op1(__binary_op2(__first1[0], __first2[0]),
369  __binary_op2(__first1[1], __first2[1]));
370  _Tp __v2 = __binary_op1(__binary_op2(__first1[2], __first2[2]),
371  __binary_op2(__first1[3], __first2[3]));
372  _Tp __v3 = __binary_op1(__v1, __v2);
373  __init = __binary_op1(__init, __v3);
374  __first1 += 4;
375  __first2 += 4;
376  }
377  }
378  for (; __first1 != __last1; ++__first1, (void) ++__first2)
379  __init = __binary_op1(__init, __binary_op2(*__first1, *__first2));
380  return __init;
381  }
382 
383  /**
384  * @brief Combine elements from two ranges and reduce
385  *
386  * @param __first1 Start of first range.
387  * @param __last1 End of first range.
388  * @param __first2 Start of second range.
389  * @param __init Starting value to add other values to.
390  * @return The final sum.
391  *
392  * Call `first1[n]*first2[n]` for each `n` in `[0,last1-first1)` and then
393  * use addition to sum those products to a single value of type `T`.
394  *
395  * The range beginning at `first2` must contain at least `last1-first1`
396  * elements.
397  */
398  template<typename _InputIterator1, typename _InputIterator2, typename _Tp>
399  _GLIBCXX20_CONSTEXPR
400  inline _Tp
401  transform_reduce(_InputIterator1 __first1, _InputIterator1 __last1,
402  _InputIterator2 __first2, _Tp __init)
403  {
404  return std::transform_reduce(__first1, __last1, __first2,
405  std::move(__init),
406  plus<>(), multiplies<>());
407  }
408 
409  /**
410  * @brief Transform the elements of a range and reduce
411  *
412  * @param __first Start of range.
413  * @param __last End of range.
414  * @param __init Starting value to add other values to.
415  * @param __binary_op The function used to perform reduction.
416  * @param __unary_op The function used to transform values from the range.
417  * @return The final sum.
418  *
419  * Call `unary_op(first[n])` for each `n` in `[0,last-first)` and then
420  * use `binary_op` to reduce the values returned by `unary_op`
421  * to a single value of type `T`.
422  */
423  template<typename _InputIterator, typename _Tp,
424  typename _BinaryOperation, typename _UnaryOperation>
425  _GLIBCXX20_CONSTEXPR
426  _Tp
427  transform_reduce(_InputIterator __first, _InputIterator __last, _Tp __init,
428  _BinaryOperation __binary_op, _UnaryOperation __unary_op)
429  {
430  if constexpr (__is_random_access_iter<_InputIterator>::value)
431  {
432  while ((__last - __first) >= 4)
433  {
434  _Tp __v1 = __binary_op(__unary_op(__first[0]),
435  __unary_op(__first[1]));
436  _Tp __v2 = __binary_op(__unary_op(__first[2]),
437  __unary_op(__first[3]));
438  _Tp __v3 = __binary_op(__v1, __v2);
439  __init = __binary_op(__init, __v3);
440  __first += 4;
441  }
442  }
443  for (; __first != __last; ++__first)
444  __init = __binary_op(__init, __unary_op(*__first));
445  return __init;
446  }
447 
448  /** @brief Output the cumulative sum of one range to a second range
449  *
450  * @param __first Start of input range.
451  * @param __last End of input range.
452  * @param __result Start of output range.
453  * @param __init Initial value.
454  * @param __binary_op Function to perform summation.
455  * @return The end of the output range.
456  *
457  * Write the cumulative sum (aka prefix sum, aka scan) of the input range
458  * to the output range. Each element of the output range contains the
459  * running total of all earlier elements (and the initial value),
460  * using `binary_op` for summation.
461  *
462  * This function generates an "exclusive" scan, meaning the Nth element
463  * of the output range is the sum of the first N-1 input elements,
464  * so the Nth input element is not included.
465  */
466  template<typename _InputIterator, typename _OutputIterator, typename _Tp,
467  typename _BinaryOperation>
468  _GLIBCXX20_CONSTEXPR
469  _OutputIterator
470  exclusive_scan(_InputIterator __first, _InputIterator __last,
471  _OutputIterator __result, _Tp __init,
472  _BinaryOperation __binary_op)
473  {
474  while (__first != __last)
475  {
476  auto __v = __init;
477  __init = __binary_op(__init, *__first);
478  ++__first;
479  *__result++ = std::move(__v);
480  }
481  return __result;
482  }
483 
484  /** @brief Output the cumulative sum of one range to a second range
485  *
486  * @param __first Start of input range.
487  * @param __last End of input range.
488  * @param __result Start of output range.
489  * @param __init Initial value.
490  * @return The end of the output range.
491  *
492  * Write the cumulative sum (aka prefix sum, aka scan) of the input range
493  * to the output range. Each element of the output range contains the
494  * running total of all earlier elements (and the initial value),
495  * using `std::plus<>` for summation.
496  *
497  * This function generates an "exclusive" scan, meaning the Nth element
498  * of the output range is the sum of the first N-1 input elements,
499  * so the Nth input element is not included.
500  */
501  template<typename _InputIterator, typename _OutputIterator, typename _Tp>
502  _GLIBCXX20_CONSTEXPR
503  inline _OutputIterator
504  exclusive_scan(_InputIterator __first, _InputIterator __last,
505  _OutputIterator __result, _Tp __init)
506  {
507  return std::exclusive_scan(__first, __last, __result, std::move(__init),
508  plus<>());
509  }
510 
511  /** @brief Output the cumulative sum of one range to a second range
512  *
513  * @param __first Start of input range.
514  * @param __last End of input range.
515  * @param __result Start of output range.
516  * @param __binary_op Function to perform summation.
517  * @param __init Initial value.
518  * @return The end of the output range.
519  *
520  * Write the cumulative sum (aka prefix sum, aka scan) of the input range
521  * to the output range. Each element of the output range contains the
522  * running total of all earlier elements (and the initial value),
523  * using `binary_op` for summation.
524  *
525  * This function generates an "inclusive" scan, meaning the Nth element
526  * of the output range is the sum of the first N input elements,
527  * so the Nth input element is included.
528  */
529  template<typename _InputIterator, typename _OutputIterator,
530  typename _BinaryOperation, typename _Tp>
531  _GLIBCXX20_CONSTEXPR
532  _OutputIterator
533  inclusive_scan(_InputIterator __first, _InputIterator __last,
534  _OutputIterator __result, _BinaryOperation __binary_op,
535  _Tp __init)
536  {
537  for (; __first != __last; ++__first)
538  *__result++ = __init = __binary_op(__init, *__first);
539  return __result;
540  }
541 
542  /** @brief Output the cumulative sum of one range to a second range
543  *
544  * @param __first Start of input range.
545  * @param __last End of input range.
546  * @param __result Start of output range.
547  * @param __binary_op Function to perform summation.
548  * @return The end of the output range.
549  *
550  * Write the cumulative sum (aka prefix sum, aka scan) of the input range
551  * to the output range. Each element of the output range contains the
552  * running total of all earlier elements, using `binary_op` for summation.
553  *
554  * This function generates an "inclusive" scan, meaning the Nth element
555  * of the output range is the sum of the first N input elements,
556  * so the Nth input element is included.
557  */
558  template<typename _InputIterator, typename _OutputIterator,
559  typename _BinaryOperation>
560  _GLIBCXX20_CONSTEXPR
561  _OutputIterator
562  inclusive_scan(_InputIterator __first, _InputIterator __last,
563  _OutputIterator __result, _BinaryOperation __binary_op)
564  {
565  if (__first != __last)
566  {
567  auto __init = *__first;
568  *__result++ = __init;
569  ++__first;
570  if (__first != __last)
571  __result = std::inclusive_scan(__first, __last, __result,
572  __binary_op, std::move(__init));
573  }
574  return __result;
575  }
576 
577  /** @brief Output the cumulative sum of one range to a second range
578  *
579  * @param __first Start of input range.
580  * @param __last End of input range.
581  * @param __result Start of output range.
582  * @return The end of the output range.
583  *
584  * Write the cumulative sum (aka prefix sum, aka scan) of the input range
585  * to the output range. Each element of the output range contains the
586  * running total of all earlier elements, using `std::plus<>` for summation.
587  *
588  * This function generates an "inclusive" scan, meaning the Nth element
589  * of the output range is the sum of the first N input elements,
590  * so the Nth input element is included.
591  */
592  template<typename _InputIterator, typename _OutputIterator>
593  _GLIBCXX20_CONSTEXPR
594  inline _OutputIterator
595  inclusive_scan(_InputIterator __first, _InputIterator __last,
596  _OutputIterator __result)
597  { return std::inclusive_scan(__first, __last, __result, plus<>()); }
598 
599  /** @brief Output the cumulative sum of one range to a second range
600  *
601  * @param __first Start of input range.
602  * @param __last End of input range.
603  * @param __result Start of output range.
604  * @param __init Initial value.
605  * @param __binary_op Function to perform summation.
606  * @param __unary_op Function to transform elements of the input range.
607  * @return The end of the output range.
608  *
609  * Write the cumulative sum (aka prefix sum, aka scan) of the input range
610  * to the output range. Each element of the output range contains the
611  * running total of all earlier elements (and the initial value),
612  * using `__unary_op` to transform the input elements
613  * and using `__binary_op` for summation.
614  *
615  * This function generates an "exclusive" scan, meaning the Nth element
616  * of the output range is the sum of the first N-1 input elements,
617  * so the Nth input element is not included.
618  */
619  template<typename _InputIterator, typename _OutputIterator, typename _Tp,
620  typename _BinaryOperation, typename _UnaryOperation>
621  _GLIBCXX20_CONSTEXPR
622  _OutputIterator
623  transform_exclusive_scan(_InputIterator __first, _InputIterator __last,
624  _OutputIterator __result, _Tp __init,
625  _BinaryOperation __binary_op,
626  _UnaryOperation __unary_op)
627  {
628  while (__first != __last)
629  {
630  auto __v = __init;
631  __init = __binary_op(__init, __unary_op(*__first));
632  ++__first;
633  *__result++ = std::move(__v);
634  }
635  return __result;
636  }
637 
638  /** @brief Output the cumulative sum of one range to a second range
639  *
640  * @param __first Start of input range.
641  * @param __last End of input range.
642  * @param __result Start of output range.
643  * @param __binary_op Function to perform summation.
644  * @param __unary_op Function to transform elements of the input range.
645  * @param __init Initial value.
646  * @return The end of the output range.
647  *
648  * Write the cumulative sum (aka prefix sum, aka scan) of the input range
649  * to the output range. Each element of the output range contains the
650  * running total of all earlier elements (and the initial value),
651  * using `__unary_op` to transform the input elements
652  * and using `__binary_op` for summation.
653  *
654  * This function generates an "inclusive" scan, meaning the Nth element
655  * of the output range is the sum of the first N input elements,
656  * so the Nth input element is included.
657  */
658  template<typename _InputIterator, typename _OutputIterator,
659  typename _BinaryOperation, typename _UnaryOperation, typename _Tp>
660  _GLIBCXX20_CONSTEXPR
661  _OutputIterator
662  transform_inclusive_scan(_InputIterator __first, _InputIterator __last,
663  _OutputIterator __result,
664  _BinaryOperation __binary_op,
665  _UnaryOperation __unary_op,
666  _Tp __init)
667  {
668  for (; __first != __last; ++__first)
669  *__result++ = __init = __binary_op(__init, __unary_op(*__first));
670  return __result;
671  }
672 
673  /** @brief Output the cumulative sum of one range to a second range
674  *
675  * @param __first Start of input range.
676  * @param __last End of input range.
677  * @param __result Start of output range.
678  * @param __binary_op Function to perform summation.
679  * @param __unary_op Function to transform elements of the input range.
680  * @return The end of the output range.
681  *
682  * Write the cumulative sum (aka prefix sum, aka scan) of the input range
683  * to the output range. Each element of the output range contains the
684  * running total of all earlier elements,
685  * using `__unary_op` to transform the input elements
686  * and using `__binary_op` for summation.
687  *
688  * This function generates an "inclusive" scan, meaning the Nth element
689  * of the output range is the sum of the first N input elements,
690  * so the Nth input element is included.
691  */
692  template<typename _InputIterator, typename _OutputIterator,
693  typename _BinaryOperation, typename _UnaryOperation>
694  _GLIBCXX20_CONSTEXPR
695  _OutputIterator
696  transform_inclusive_scan(_InputIterator __first, _InputIterator __last,
697  _OutputIterator __result,
698  _BinaryOperation __binary_op,
699  _UnaryOperation __unary_op)
700  {
701  if (__first != __last)
702  {
703  auto __init = __unary_op(*__first);
704  *__result++ = __init;
705  ++__first;
706  if (__first != __last)
707  __result = std::transform_inclusive_scan(__first, __last, __result,
708  __binary_op, __unary_op,
709  std::move(__init));
710  }
711  return __result;
712  }
713 
714  /// @} group numeric_ops
715 
716 _GLIBCXX_END_NAMESPACE_VERSION
717 } // namespace std
718 
719 // Parallel STL algorithms
720 # if _PSTL_EXECUTION_POLICIES_DEFINED
721 // If <execution> has already been included, pull in implementations
722 # include <pstl/glue_numeric_impl.h>
723 # else
724 // Otherwise just pull in forward declarations
725 # include <pstl/glue_numeric_defs.h>
726 # define _PSTL_NUMERIC_FORWARD_DECLARED 1
727 # endif
728 
729 // Feature test macro for parallel algorithms
730 # define __cpp_lib_parallel_algorithm 201603L
731 #endif // C++17
732 
733 #endif /* _GLIBCXX_NUMERIC */