29.4 Numerikus algoritmusok és kiértékelési sorrend

A <numeric> fejléc algoritmusai első ránézésre egyszerű összegzések, mégis a legtöbb meglepetés nem a műveletben, hanem a kiértékelési sorrendben és az accumulator típusában rejlik. Ez a rész a Rendezés, keresés és comparator szerződés (29.3) comparator-szerződésére épül: ott a rendezés és keresés előfeltételeit rögzítettük; itt azt vizsgáljuk, milyen szerződést ír elő az accumulate, a reduce, a prefixösszegek és az inner_product. A lebegőpontos hibák hátterét a Lebegőpontos reprezentáció és numerikus hibák (4.3), a signed túlcsordulásét a Promóciók, signedness és túlcsordulás (4.2) tárgyalta. A cél nem a numerikus analízis teljes kurzusa, hanem hogy a standard algoritmus kiválasztásakor a sorrendi és típusbeli következményeket tudatosan kezeld.

29.4.1 Az accumulator típusa és az accumulate szerződése

Az std::accumulate C++98 óta a <numeric>-ban található. A legegyszerűbb alak:

template<class InputIt, class T>
T accumulate(InputIt first, InputIt last, T init);

A hatás: init = init + *first, majd init = init + *(first+1), és így tovább, balról jobbra, a megadott sorrendben. A T a kezdeti érték típusa, az accumulator típusa és a visszatérési típus. A *it minden lépésben T-vé konvertálódik. Ha a tartomány üres, az eredmény init. A bináris műveletet elfogadó változat:

template<class InputIt, class T, class BinaryOp>
T accumulate(InputIt first, InputIt last, T init, BinaryOp op);

Itt init = op(init, *it), szintén balról jobbra. Az op nem módosíthatja a bemeneti tartományt, és nem hívhat olyan műveletet, amely a sorrendet megváltoztatná; a szerződés a balról jobbra haladás.

Az accumulator típusa nem a vektor értéktípusa, hanem a T. Ez a leggyakoribb hiba forrása:

std::vector<int> v{1'000'000'000, 1'000'000'000, 1'000'000'000};
auto a = std::accumulate(v.begin(), v.end(), 0);    // int accumulator
auto b = std::accumulate(v.begin(), v.end(), 0LL);  // long long accumulator

Az a kiszámítása során a második lépésben 2'000'000'000 még ábrázolható int-ben, a harmadik viszont már nem: signed overflow, ami a C++ szabvány szerint nem definiált viselkedés. A b esetén az összeadás long long-ban történik, mert a T long long, és a *it int értéke az összeadás előtt long long-ra konvertálódik. A 0LL tehát nem kozmetika: az accumulator típusát választja meg. Ugyanígy 0.0 és 0.0f más lebegőpontos pontosságot ad. Az accumulate determinisztikus: ugyanaz a bemenet és ugyanaz a T ugyanazt az eredményt adja, feltéve, hogy a művelet determinisztikus. Ez teszi alkalmassá arra, amikor a reprodukálhatóság követelmény.

29.4.2 reduce, transform_reduce és a sorrend lazítása

A C++17 óta a <numeric> olyan algoritmusokat is tartalmaz, amelyek a sorrendet nem rögzítik. A reduce legegyszerűbb alakja:

template<class InputIt, class T>
T reduce(InputIt first, InputIt last, T init);

A reduce a std::plus<>-t használja, és a részeredményeket tetszőleges sorrendben és csoportosításban vonhatja össze. A standard előfeltételezi, hogy a binary_op asszociatív és kommutatív; ha nem az, az eredmény nem specifikált. A T itt is az accumulator típusa. A reduce nem garantál balról jobbra sorrendet, ezért lebegőpontoson eltérhet az accumulate-tól. A transform_reduce alakja:

template<class InputIt, class T, class BinaryOp, class UnaryOp>
T transform_reduce(InputIt first, InputIt last, T init,
                   BinaryOp binary_op, UnaryOp unary_op);

Az unary_op-ot minden elemre alkalmazza, majd a binary_op-pal redukál. A sorrend itt is lazított. Az unary_op és a binary_op ne legyen mellékhatásos; párhuzamos policy esetén a mellékhatás adatversenyhez vezethet. A std::execution::par nem automatikus gyorsítás, és a policy a sorrendet tovább lazíthatja; a részleteket az Execution policy és párhuzamos algoritmusok (29.5) tárgyalja.

Fontos terminológiai különbség: a reduce nem „nem definiált”, hanem a sorrendre nézve nem specifikált. A nem asszociatív művelet előfeltétel-sértés; a lebegőpontos plus pedig tipikusan nem asszociatív. Ezért a reduce lebegőpontoson csak akkor használható, ha a nem specifikált sorrendet és az ebből következő eltérést tudatosan elfogadjuk.

29.4.3 Prefixösszegek: partial_sum, inclusive_scan, exclusive_scan

template<class InputIt, class OutputIt>
OutputIt partial_sum(InputIt first, InputIt last, OutputIt d_first);

Hatás: balról jobbra halad, és minden lépésben a kimenetre írja az addigi összeget. A részösszegek típusát nem a kimenet, hanem a bemenet értéktípusa adja ([C++20] [partial.sum]): az akkumulátor InputIt value_type-jú változó (acc = *first, majd acc = std::move(acc) + *i), és csak az így képzett értéket írja a kimenetre. Egy int bemenet tehát int-ben összegződik akkor is, ha a kimenet long long. Széles akkumulációhoz a bemenetet kell előbb szélesíteni, vagy olyan algoritmust választani, amelynél a részeredmény típusa igazoltan a széles típus, például inclusive_scan-t long long kezdőértékkel:

// partial_sum_accumulator.cpp — C++20: az akkumulátor a BEMENET típusa
#include <climits>
#include <cstdio>
#include <functional>
#include <numeric>
#include <vector>

int main(int argc, char**) {
    const std::vector<int> in{INT_MAX, 1};
    std::vector<long long> out(in.size());

    if (argc > 1) {   // negatív ág: a részösszeg int-ben képződik
        std::partial_sum(in.begin(), in.end(), out.begin());
    } else {          // javított: a scan long long kezdőértékkel akkumulál
        std::inclusive_scan(in.begin(), in.end(), out.begin(),
                            std::plus<>{}, 0LL);
    }
    std::printf("%lld %lld\n", out[0], out[1]);
}

Buildjegyzőkönyv (lab-gcc profil, -fsanitize=undefined -fno-sanitize-recover=undefined). Argumentum nélkül (javított ág) a kimenet 2147483647 2147483648. Argumentummal (partial_sum ág) az UBSan signed integer overflow: 2147483647 + 1 cannot be represented in type 'int' hibát jelez, a long long cél ellenére.

template<class InputIt, class OutputIt>
OutputIt inclusive_scan(InputIt first, InputIt last, OutputIt d_first);

template<class InputIt, class OutputIt, class BinaryOp>
OutputIt inclusive_scan(InputIt first, InputIt last, OutputIt d_first,
                        BinaryOp op);

template<class InputIt, class OutputIt, class BinaryOp, class T>
OutputIt inclusive_scan(InputIt first, InputIt last, OutputIt d_first,
                        BinaryOp op, T init);

template<class InputIt, class OutputIt, class T>
OutputIt exclusive_scan(InputIt first, InputIt last, OutputIt d_first, T init);

Az inclusive_scan negyedik paramétere tehát a bináris művelet, nem a kezdőérték; kezdőértéket csak az ötparaméteres alak fogad, a művelet után. Az inclusive_scan(a.begin(), a.end(), out.begin(), 0) ezért nem fordul le: a 0-t műveletként próbálná hívni. Az exclusive_scan-nél fordított a helyzet, ott a negyedik paraméter a kötelező kezdőérték.

Az inclusive_scan a d_first[i]-be az első i+1 elem összegét írja (az ötparaméteres alakban az init-tel együtt); az exclusive_scan a d_first[i]-be az init és az első i elem összegét, tehát az első kimenet init. A scan-ek sorrendje nem specifikált, és a műveletnek asszociatívnak kell lennie. Kezdőérték megadásakor a T az init típusa, és a részeredmények T-ben képződnek; kezdőérték nélkül a bemenet értéktípusában. A partial_sum és az inclusive_scan közötti különbség tehát nem csak az init: az előbbi determinisztikus, az utóbbi párhuzamosítható és nem determinisztikus. Egész számokon, ahol az összeadás asszociatív, a két eredmény megegyezik; lebegőpontoson nem feltétlenül. A bemenet és a kimenet átfedése a standard szövegében ellenőrizendő; a biztonságos minta külön kimeneti buffer.

29.4.4 inner_product, overflow és a kezdeti érték csapdája

template<class InputIt1, class InputIt2, class T>
T inner_product(InputIt1 first1, InputIt1 last1, InputIt2 first2, T init);

Hatás: init = init + (*first1 * *first2), balról jobbra. A T a kezdeti érték és az accumulator típusa. A csapda: a *first1 * *first2 szorzás a szokásos aritmetikai konverzió szerint történik, még mielőtt az init-hez adódna. Ha mindkét elem int, a szorzat int, és már a szorzat is túlcsordulhat. A 0LL kezdeti érték önmagában nem véd:

std::vector<int> a{100'000}, b{100'000};
auto rossz = std::inner_product(a.begin(), a.end(), b.begin(), 0LL);
// a * b == 10'000'000'000, de a szorzat int-ben túlcsordul

auto jo = std::inner_product(a.begin(), a.end(), b.begin(), 0LL,
                             std::plus<long long>{},
                             std::multiplies<long long>{});

Ekkor a szorzás is long long-ban történik. Ugyanez transform_reduce-szal:

auto jo2 = std::transform_reduce(a.begin(), a.end(), b.begin(), 0LL,
                                 std::plus<long long>{},
                                 [](int x, int y) {
                                     return static_cast<long long>(x) * y;
                                 });

A transform_reduce C++17-es; az inner_product bináris op-os változatában az első művelet az összeadás, a második a szorzás. Az inner_product sorrendje balról jobbra, determinisztikus; a transform_reduce sorrendje nem specifikált.

29.4.5 Lebegőpontos asszociativitás és reprodukálhatóság

A lebegőpontos összeadás nem asszociatív. A laborprofilban rögzített x86-64 System V ABI és IEEE 754 binary32 feltételezés mellett a float 24 bites szignifikáns pontossággal rendelkezik. A

2^{24} = 16\,777\,216

és a

2^{25}

között a szomszédos ábrázolható értékek különbsége 2, ezért a

16\,777\,217

nem ábrázolható; a 16777216.0f + 1.0f eredménye a legközelebbi ábrázolható érték, azaz 16777216.0f. Ebből a levezetésből:

A két csoportosítás különbözik. Ezért a reduce és a scan-ek, amelyek a csoportosítást megváltoztathatják, lebegőpontoson nem adnak garantáltan ugyanolyan eredményt, mint az accumulate. A std::numeric_limits<float>::is_iec559 értéke a laborban ellenőrizendő; a C++ szabvány nem garantálja az IEEE 754-et, ezért a fenti példa a laborprofilra érvényes, nem minden platformra.

Reprodukálhatóság: ha a követelmény az, hogy ugyanaz a bemenet minden futásban bitre ugyanazt az eredményt adja, akkor az accumulate vagy a partial_sum a választás, vagy egy rögzített blokkolású saját összegzés. A reduce és a scan-ek ezt nem garantálják; a par policy tovább növeli a változékonyságot. A long double nem feltétlenül hordozhatóbb: x86-64-en 80 bites, AArch64-en a fordító és a libc döntése szerint eltérhet. A Kahan-összegzés csökkenti a kerekítési hibát, de nem teszi asszociatívvá a műveletet, és nem ad bitre reprodukálható eredményt minden sorrendben.

Az alábbi példa a fő kísérleteket gyűjti össze. A kód a laborban fordítandó és futtatandó; a kimenet a mérési jegyzőkönyvbe kerül, nem tekintjük előre rögzített eredménynek.

// numlab.cpp — C++20, normál labor: UBSan alatt is hibamentes
#include <numeric>
#include <vector>
#include <iostream>
#include <limits>
#include <functional>

int main() {
    // 1. accumulator típusa
    std::vector<int> v{1'000'000'000, 1'000'000'000, 1'000'000'000};
    // auto bad = std::accumulate(v.begin(), v.end(), 0); // signed overflow: UB
    auto good = std::accumulate(v.begin(), v.end(), 0LL);
    std::cout << "accumulate 0LL: " << good << '\n';

    // 2. lebegőpontos sorrend
    std::vector<float> f{16777216.0f, 1.0f, 1.0f};
    auto acc = std::accumulate(f.begin(), f.end(), 0.0f);
    auto red = std::reduce(f.begin(), f.end(), 0.0f);
    std::cout << "accumulate float: " << acc << '\n';
    std::cout << "reduce float (unspecified order): " << red << '\n';

    // 3. prefixösszegek
    std::vector<int> a{1, 2, 3, 4};
    std::vector<int> ps(a.size()), inc(a.size()), exc(a.size());
    std::partial_sum(a.begin(), a.end(), ps.begin());
    std::inclusive_scan(a.begin(), a.end(), inc.begin());   // nincs init
    std::exclusive_scan(a.begin(), a.end(), exc.begin(), 0);
    std::cout << "partial_sum: ";
    for (int x : ps) std::cout << x << ' ';
    std::cout << "\ninclusive_scan: ";
    for (int x : inc) std::cout << x << ' ';
    std::cout << "\nexclusive_scan: ";
    for (int x : exc) std::cout << x << ' ';
    std::cout << '\n';

    // 4. inner_product overflow
    std::vector<int> x{100'000}, y{100'000};
    // A szorzás is long long-ban történik; a 0LL kezdőérték önmagában
    // nem elég (lásd a külön negatív próbát: inner_product_overflow.cpp).
    auto good_ip = std::inner_product(x.begin(), x.end(), y.begin(), 0LL,
                                      std::plus<long long>{},
                                      std::multiplies<long long>{});
    std::cout << "inner_product long long: " << good_ip << '\n';

    // 5. numeric_limits
    std::cout << "int max: " << std::numeric_limits<int>::max() << '\n';
    std::cout << "float is_iec559: "
              << std::numeric_limits<float>::is_iec559 << '\n';
}

g++ -std=c++20 -O2 -Wall -Wextra -pedantic numlab.cpp -o numlab
./numlab

A -O0/-O2 önmagában nem változtathatja meg a szabványos lebegőpontos eredményt; a -ffast-math viszont igen, ezért a reprodukálhatósági kísérletben ne használd. A normál labor nem tartalmaz végrehajtott nem definiált viselkedést, ezért UBSan alatt is hibamentesen fut. Az inner_product túlcsordulását külön negatív forrás mutatja: a 0LL kezdőérték csak az összeadás akkumulátorát szélesíti, a szorzás int * int marad.

// inner_product_overflow.cpp — NEGATÍV próba: szándékos signed overflow
#include <iostream>
#include <numeric>
#include <vector>

int main() {
    std::vector<int> x{100'000}, y{100'000};
    // A 0LL csak az összeadás akkumulátora; a szorzás int * int marad,
    // és 10^10 nem fér el int-ben: nem definiált viselkedés.
    auto bad_ip = std::inner_product(x.begin(), x.end(), y.begin(), 0LL);
    std::cout << bad_ip << '\n';
}

g++ -std=c++20 -O1 -g -fsanitize=undefined -fno-sanitize-recover=undefined \
    numlab.cpp -o numlab_ubsan
./numlab_ubsan                       # hibamentes
g++ -std=c++20 -O1 -g -fsanitize=undefined -fno-sanitize-recover=undefined \
    inner_product_overflow.cpp -o ip_overflow
./ip_overflow 2>&1 | grep 'signed integer overflow'   # várt diagnosztika

Buildjegyzőkönyv (lab-gcc profil). A numlab.cpp figyelmeztetés nélkül fordult, kimenete: accumulate 0LL: 3000000000, accumulate float: 1.67772e+07, reduce float (unspecified order): 1.67772e+07, partial_sum: 1 3 6 10, inclusive_scan: 1 3 6 10, exclusive_scan: 0 1 3 6, inner_product long long: 10000000000, int max: 2147483647, float is_iec559: 1. UBSan alatt hibamentes. A negatív próba a signed integer overflow: 100000 * 100000 cannot be represented in type 'int' diagnosztikát adta. A javítás előtti változat (inclusive_scan(..., 0)) le sem fordult: '__binary_op' cannot be used as a function.

Kísérlet	Algoritmus	Accumulator	Sorrend	Reprodukálható?
A	`accumulate`	`float`	balról jobbra	igen
B	`reduce`	`float`	nem specifikált	nem
C	`accumulate`	`double`	balról jobbra	igen
D	`reduce`	`double`	nem specifikált	nem
E	`partial_sum`	`int`	balról jobbra	igen
F	`inclusive_scan`	`int`	nem specifikált	igen (asszociatív)

29.4.6 numeric_limits, gcd/lcm és ellenőrzött bemenet

A std::numeric_limits a típus tulajdonságait adja: max(), lowest(), epsilon(), digits10, max_digits10, is_integer, is_signed, is_iec559. A max_digits10 a round-triphez szükséges decimális számjegyek száma; a digits10 a garantáltan pontos decimális számjegyek száma. A lowest() a legkisebb véges érték, amely negatív is lehet; a min() a legkisebb pozitív normalizált érték. Az epsilon() az 1.0 és a következő ábrázolható érték különbsége.

A C++17-es std::gcd és std::lcm a <numeric>-ban található. Ha a std::lcm(a,b) eredménye nem ábrázolható a common_type-ban, a viselkedés nem definiált. A gcd/lcm a |a| és |b| értékekkel dolgozik; ha a vagy b INT_MIN, a |INT_MIN| nem ábrázolható int-ben, ezért már a bemenet is problémás. Ellenőrzött bemenetnél a std::from_chars használható egész számokra C++17 óta: nem dob, nem allokál, és nem függ a locale-tól. A lebegőpontos from_chars támogatása a standard library verziójától függ; a std::stod/std::strtod fallback lehet, de az allokálhat és locale-függő. A from_chars eredménye std::errc, amelyet explicit ellenőrizni kell. A feature-test makrókról a Feature-test, standard library és build capability (33.5) ad további támpontot.

template<class T>
bool lcm_ok(T a, T b) {
    if (a == 0 || b == 0) return true;
    const T g = std::gcd(a, b);
    const T q = a / g;
    // |q * b| <= max?
    return std::abs(q) <= std::numeric_limits<T>::max() / std::abs(b);
}

Ez a vázlat a T előjeles esetére; a std::abs(INT_MIN) csapdáját külön kell kezelni. A __builtin_mul_overflow GCC/Clang kiterjesztés, nem hordozható; a laborban jelölni kell.

29.4.7 Tipikus hibák és javításuk

29.4.8 Feladatok

29.4.9 Megoldások

double sum_ordered(const std::vector<double>& v) {
    return std::accumulate(v.begin(), v.end(), 0.0);
}

double sum_unordered(const std::vector<double>& v) {
    return std::reduce(v.begin(), v.end(), 0.0);
}

Az accumulate balról jobbra halad, ezért ugyanarra a bemenetre ugyanazt a bitmintát adja, feltéve, hogy a fordítás nem lazítja a lebegőpontos szabályokat. A reduce sorrendje nem specifikált, és a double összeadás nem asszociatív, ezért az eredmény futásonként vagy implementációnként eltérhet. Build: g++ -std=c++20 -O2 -Wall -Wextra -pedantic main.cpp -o sumlab.

auto jo = std::inner_product(a.begin(), a.end(), b.begin(), 0LL,
                             std::plus<long long>{},
                             std::multiplies<long long>{});

Ezzel a szorzás és az összegzés is long long-ban történik. Ellenőrzési feltétel: a szorzatok és a részösszegek a long long tartományában maradnak; ezt a bemenet mérete és a std::numeric_limits<long long>::max() alapján kell igazolni.

29.4.10 Összefoglaló

A numerikus algoritmusok kiválasztása három kérdésre vezethető vissza: milyen típusú az accumulator, milyen sorrendi szerződést vállal az algoritmus, és milyen tartományban marad a művelet. Az accumulate és a partial_sum balról jobbra determinisztikus; a reduce, a scan-ek és a transform_reduce sorrendje nem specifikált, ezért lebegőpontoson csak tudatos kompromisszummal használhatók. Az inner_product szorzattípusa külön csapda. A numeric_limits, a gcd/lcm és a from_chars az ellenőrzött bemenet eszközei. A laborban a mérési jegyzőkönyv rögzítse a fordítót, az optimalizálást és a kapott bitmintát.

29.4.1 Az accumulator típusa és az `accumulate` szerződése

29.4.2 `reduce`, `transform_reduce` és a sorrend lazítása

29.4.3 Prefixösszegek: `partial_sum`, `inclusive_scan`, `exclusive_scan`

29.4.4 `inner_product`, overflow és a kezdeti érték csapdája

29.4.5 Lebegőpontos asszociativitás és reprodukálhatóság

29.4.6 `numeric_limits`, `gcd`/`lcm` és ellenőrzött bemenet

29.4.7 Tipikus hibák és javításuk

29.4.8 Feladatok

29.4.9 Megoldások

29.4.10 Összefoglaló