Тестирование алгоритма деления больших чисел на С++ с использованием Python C API

c9395c7a2f64b05a3502ecea55ceb6b3

Введение

Ранее был предложен некоторый Алгоритм деления 2W-битовых чисел с использованием операций над W-битовыми числами. Для тестирования использовались целые числа языка С++, что не позволяло проверять, например, 128-битные целые числа. Однако, в язык Python встроена поддержка целых чисел неограниченной ширины (Big Integer), а также имеется API для вызова методов Python из программ на языке С/С++. Это позволяет протестировать разные алгоритмы с числами, в том числе деление, используя в качестве результата строковое представление чисел.

В данной статье расписаны шаги для использования Python C API в программе на языке С++, а также показан пример вызова оператора деления двух целых чисел с возвратом результата в виде строки С. Использовалась следующая программная конфигурация:

  • Windows 11, WSL Ubuntu 22.04.3 LTS, GCC 12.3.0.

  • IDE Visual Studio Code актуальной версии.

  • Встроенный в Ubuntu Python 3.10.

Настройка среды

Предварительно требуется установить набор build-essential и средство разработчика python3-dev.

  1. Открываем консоль Ubuntu и создаем директорию с проектом:

mkdir division_test
  1. Переходим в созданный каталог и создаем точку входа в программу:

cd division_test
touch main.cpp
  1. Запускаем IDE Visual Studio Code. Устанавливаем расширения WSL и C/C++.

  2. Перезапускаем Visual Studio Code, и через IDE открываем директорию созданного проекта, division_test. В итоге, на панели EXPLORER должен отобразиться файл main.cpp, который надо открыть для редактирования.

  3. Копируем/создаем тестовый исходный код.

  4. Через меню IDE «Terminal/Configure Default Build Task…» создаем конфигурацию сборки tasks.json.

  5. В файл конфигурации tasks.json добавляем линковку библиотеки Python, а также дополнительные опции: стандарт С++20 и оптимизацию компиляции O3:

...
"args": [
    "-fdiagnostics-color=always",
    "-O3",
    "-std=c++20",
    "${file}",
    "-o",
    "${fileDirname}/${fileBasenameNoExtension}",
    "-lpython3.10"
],
...
  1. Активируем точку входа: файл main.cpp. Проверяем подключение заголовочного файла Python

#include 
  1. Открываем терминал «Terminal/New Terminal», активируем фокус (курсор) на терминале и нажимаем «Ctrl+Shift+B», не забыв переключиться на английский язык EN. Программа должна собраться. Запустить программу можно в терминале как обычно: ./main.

Использование API C Python

После настройки среды и сборки тестового проекта можно свободно пользоваться методами API C Python.

Для тестирования деления двух 128-битных чисел был создан вспомогательный класс PythonCaller.

Вспомогательный класс деления чисел на API C Python

class PythonCaller {
public:

    PythonCaller() {
        Py_Initialize();
        mMain = PyImport_AddModule("__main__");
        mGlobalDictionary = PyModule_GetDict(mMain);
        mLocalDictionary = PyDict_New();
        std::cout << "Python was initialized!\n";
    }

    ~PythonCaller() {
        Py_Finalize();
        std::cout << "~Python was finalized!\n";
    }

    /**
     * Делит два 128-битных числа.
     * @param X Делимое.
     * @param Y Делитель.
     * @return Частное от деления, объект Python.
    */
    PyObject* Divide(U128 X, U128 Y) const {
        const char* pythonScript = "quotient = nominator // denominator\n";
        PyDict_SetItemString(mLocalDictionary, "nominator", PyLong_FromString(X.value().c_str(), nullptr, 10));
        PyDict_SetItemString(mLocalDictionary, "denominator", PyLong_FromString(Y.value().c_str(), nullptr, 10));
        PyRun_String(pythonScript, Py_file_input, mGlobalDictionary, mLocalDictionary);
        return PyDict_GetItemString(mLocalDictionary, "quotient");
    }
    /**
     * Сравнивает два частных от деления.
     * @param quotient Частное от деления Python, объект Python.
     * @param reference Частное от деления С++, строка.
    */
    bool Compare(PyObject* quotient, const char* reference) const {
        PyObject* repr = PyObject_Repr(quotient);
        PyObject* str = PyUnicode_AsEncodedString(repr, "utf-8", "~E~");
        const char *bytes = PyBytes_AsString(str);
        const bool is_ok = strcmp(bytes, reference) == 0;
        if (!is_ok) {
            printf("Python: %s\n", bytes);
            printf("C++: %s\n", reference);
        }
        Py_XDECREF(repr);
        Py_XDECREF(str);
        return is_ok;
    }
private:

    PyObject* mMain = nullptr;
    PyObject* mGlobalDictionary = nullptr;
    PyObject* mLocalDictionary = nullptr;
};

В конструкторе класса инициализируется среда Python с вспомогательными словарями для ввода-вывода результатов счета.

Для установки длинного Python числа типа int из С++ используется API метод PyLong_FromString, принимающий строковое представление числа (const char*, первый аргумент) и основание системы счисления (третий аргумент). Скрипт Python создается в виде строки с непосредственной записью деления через переменные, имена которых совпадают с ключами обменного словаря mLocalDictionary. Поэтому результат деления (частное) будет доступен по соответствующему ключу автоматически через API метод PyDict_GetItemString в виде Python объекта.

Для преобразования Python объекта в строку последовательно используются три API метода:

Первый метод преобразует объект Python во внутреннее строковое представление Python, выдавая на выход измененный объект Python. Это аналог метода repr (.) Python.

Второй метод формирует байтовое представление строки, пригодное для передачи в некоторый канал связи. На выходе имеем объект Python.

Наконец, третий метод позволяет по готовому объекту Python получить указатель на байты искомой строки. В нашем случае анализируются числа, поэтому байты совпадают с символами: цифрами 0…9, знаком минус, поэтому эти байты можно напрямую сравнивать с байтами строкового представления других чисел — результатов счета некоторого алгоритма на С++.

Набор методов для тестирования приведен ниже.

Методы тестирования алгоритма деления 128-битных чисел

#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 
#include 

using u64 = uint64_t;
using ULOW = u64; // Тип половинок: старшей и младшей частей составного числа.

static_assert(CHAR_BIT == 8);

struct Quadrupole { // Структура для задания дроби (A*M + B) / (C*M + D).
                    // M - множитель системы счисления, 2^W, W - битовая ширина половинок.
    ULOW A;
    ULOW B;
    ULOW C;
    ULOW D;
};

struct Signess { // Структура для задания знаков двух чисел.
    bool s1;
    bool s2;
};

static auto const seed = std::random_device{}();

/***
 * Генератор случайных чисел.
 */
auto roll_ulow = [urbg = std::mt19937{seed},
                distr = std::uniform_int_distribution{}]() mutable -> ULOW {
    return distr(urbg);
};

auto roll_uint = [urbg = std::mt19937{seed},
                distr = std::uniform_int_distribution{}]() mutable -> uint {
    return distr(urbg);
};

auto roll_bool = [urbg = std::mt19937{seed},
                distr = std::uniform_int_distribution{}]() mutable -> bool {
    return distr(urbg) % 2;
};

/**
 * Тест деления двух 128-битных чисел.
 * @param z1 Делимое.
 * @param z2 Делитель.
 * @return Успех/неудача.
*/
bool test_div(U128 z1, U128 z2, PythonCaller& caller) {
    const U128 z3 = z1 / z2;
    PyObject* quotient = caller.Divide(z1, z2);
    return caller.Compare(quotient, z3.value().c_str());
}

/**
 * Полуслучайный тест деления 128-битных чисел используя ограниченный набор значений вблизи угловых и граничных.
 * @param N Количество внешних итераций.
*/
void test_division_semi_randomly(long long N) {
    if (N < 1) {
        std::cout << "Skipped!\n";
        return;
    }
    PythonCaller caller;
    const std::vector choice {0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 
                                    65535, 65534, 65533, 65532, 65531, 65530, 
                                    16384, 16383, 16382, 16385, 16386, 16387, 16388, 
                                    -1ull, -2ull, -3ull, -4ull, -5ull, -6ull, -7ull};
    auto make_test = [&caller](const Quadrupole q, const Signess s) -> bool {
        return test_div(U128{.mHigh = q.A, .mLow = q.B, .mSign = s.s1},
                        U128{.mHigh = q.C, .mLow = q.D, .mSign = s.s2},
                        caller);
    };
    auto get_quadrupole = [&choice]() -> Quadrupole {
        auto idx1 = roll_uint() % choice.size();
        auto idx2 = roll_uint() % choice.size();
        auto idx3 = roll_uint() % choice.size();
        auto idx4 = roll_uint() % choice.size();
        Quadrupole q {choice[idx1], choice[idx2], choice[idx3], choice[idx4]};
        return q;
    };
    long long counter = 0;
    long long external_iterations = 0;
    bool is_ok = true;
    while (external_iterations < N) {
        ++counter;
        const Quadrupole q = get_quadrupole();
        const Signess s{roll_bool(), roll_bool()};
        if (q.C == 0 && q.D == 0) {
            continue;
        }
        is_ok &= make_test(q, s);
        assert(is_ok);
        if (counter % (1ll << 20) == 0) {
            external_iterations++;
            std::cout << "... iterations: " << counter << ". External: " << 
                external_iterations << " from " << N << '\n';
        }
    }
}

/**
 * Случайный тест деления 128-битных чисел.
 * @param N Количество внешних итераций.
*/
void test_division_randomly(long long N) {
    if (N < 1) {
        std::cout << "Skipped!\n";
        return;
    }
    PythonCaller caller;
    auto make_test = [&caller](const Quadrupole q, const Signess s) -> bool {
        return test_div(U128{.mHigh = q.A, .mLow = q.B, .mSign = s.s1},
                        U128{.mHigh = q.C, .mLow = q.D, .mSign = s.s2},
                        caller);
    };
    auto get_quadrupole = []() -> Quadrupole {
        Quadrupole q {roll_ulow(), roll_ulow(), roll_ulow(), roll_ulow()};
        return q;
    };
    long long counter = 0;
    long long external_iterations = 0;
    bool is_ok = true;
    while (external_iterations < N) {
        ++counter;
        const Quadrupole q = get_quadrupole();
        const Signess s{roll_bool(), roll_bool()};
        if (q.C == 0 && q.D == 0) {
            continue;
        }
        is_ok &= make_test(q, s);
        assert(is_ok);
        if (counter % (1ll << 20) == 0) {
            external_iterations++;
            std::cout << "... iterations: " << counter << ". External: " << 
                external_iterations << " from " << N << '\n';
        }
    }
}


int main(int argc, char* argv[]) {
    long long N = 10;
    if (argc > 1) {
        N = atoi(argv[1]);
        std::cout << "You set the number of external iterations N: " << N << '\n';
    }
    std::cout << "Run semi-random test...\n";
    test_division_semi_randomly(N);
    std::cout << "Ok\n";

    std::cout << "Run random test...\n";
    test_division_randomly(N);
    std::cout << "Ok\n";
    return 0;
}

Здесь используются тесты двух типов:

  1. Полу-случайный тест,

  2. Случайный тест.

Полу-случайный тест проходится по граничным и угловым случаям проверяемого алгоритма и завершается, что позволяет относительно быстро выявить изъяны алгоритма. Случайный же тест позволяет просмотреть в том числе и некоторые отличные от граничных случаев.

Реализация 128-битных чисел U128 приведена ниже.

Структура U128

static constexpr char DIGITS[10]{'0', '1', '2', '3', '4',
                                 '5', '6', '7', '8', '9'};

// High/Low структура 128-битного числа со знаком и флагом переполнения.
// Для иллюстрации алгоритма деления двух U128 чисел реализованы основные
// арифметические операторы, кроме умножения двух U128 чисел.
struct U128 {
    // Битовая полуширина половинок.
    static constexpr int mHalfWidth = (sizeof(ULOW) * CHAR_BIT) / 2;
    // Наибольшее значение половинок, M-1.
    static constexpr ULOW mMaxULOW = ULOW(-1);
    ULOW mHigh = 0;
    ULOW mLow = 0;
    bool mSign = 0;
    bool mOverflow = 0;

bool is_zero() const {
    return (mLow | mHigh | mOverflow) == 0;
}

bool is_negative() const {
    return mSign;
}

bool is_non_negative() const {
    return !mSign && !mOverflow;
}

bool is_nonzero_negative() const {
    return (mLow | mHigh) && mSign && !mOverflow;
}

bool is_overflow() const {
    return mOverflow;
}

constexpr static U128 get_zero() {
    return U128 {.mHigh = 0, .mLow = 0, .mSign = 0, .mOverflow = 0};
}

constexpr static U128 get_unit() {
    return U128 {.mHigh = 0, .mLow = 1, .mSign = 0, .mOverflow = 0};
}

constexpr static U128 get_unit_neg() {
    return U128 {.mHigh = 0, .mLow = 1, .mSign = 1, .mOverflow = 0};
}

U128 operator+(U128 rhs) const {
    U128 result{};
    U128 X = *this;
    if (X.is_negative() && !rhs.is_negative()) {
        X.mSign = 0;
        result = rhs - X;
        return result;
    }
    if (!X.is_negative() && rhs.is_negative()) {
        rhs.mSign = 0;
        result = X - rhs;
        return result;
    }
    result.mLow = X.mLow + rhs.mLow;
    const ULOW c1 = result.mLow < std::min(X.mLow, rhs.mLow);
    result.mHigh = X.mHigh + rhs.mHigh;
    const int c2 = result.mHigh < std::min(X.mHigh, rhs.mHigh);
    ULOW tmp = result.mHigh;
    result.mHigh = tmp + c1;
    const int c3 = result.mHigh < std::min(tmp, c1);
    result.mOverflow = c2 || c3;
    if (X.mSign && rhs.mSign) {
        result.mSign = 1;
    }
    return result;
}

U128& operator+=(U128 other) {
    *this = *this + other;
    return *this;
}

U128 operator-(U128 rhs) const {
    U128 result{};
    U128 X = *this;
    if (X.is_negative() && !rhs.is_negative()) {
        rhs.mSign = 1;
        result = rhs + X;
        return result;
    }
    if (!X.is_negative() && rhs.is_negative()) {
        rhs.mSign = 0;
        result = X + rhs;
        return result;
    }
    if (X.is_negative() && rhs.is_negative()) {
        rhs.mSign = 0;
        X.mSign = 0;
        result = rhs - X;
        return result;
    }
    if (X.is_zero()) {
        result = rhs;
        result.mSign = rhs.mSign ^ 1;
        return result;
    }
    result.mLow = X.mLow - rhs.mLow;
    result.mHigh = X.mHigh - rhs.mHigh;
    const bool borrow = X.mLow < rhs.mLow;
    const bool hasUnit = X.mHigh > rhs.mHigh;
    if (borrow && hasUnit) {
        result.mHigh -= ULOW(1);
    }
    if (borrow && !hasUnit) {
        result = rhs - X;
        result.mSign ^= 1;
        return result;
    }
    if (!borrow && X.mHigh < rhs.mHigh) {
        result.mHigh = -result.mHigh - ULOW(result.mLow != 0);
        result.mLow = -result.mLow;
        result.mSign = 1;
    }
    return result;
}

U128& operator-=(U128 other) {
    *this = *this - other;
    return *this;
}

U128 mult64(ULOW x, ULOW y) const {
    constexpr ULOW MASK = (ULOW(1) << mHalfWidth) - 1;
    const ULOW x_low = x & MASK;
    const ULOW y_low = y & MASK;
    const ULOW x_high = x >> mHalfWidth;
    const ULOW y_high = y >> mHalfWidth;
    const ULOW t1 = x_low * y_low;
    const ULOW t = t1 >> mHalfWidth;
    const ULOW t21 = x_low * y_high;
    const ULOW q = t21 >> mHalfWidth;
    const ULOW p = t21 & MASK;
    const ULOW t22 = x_high * y_low;
    const ULOW s = t22 >> mHalfWidth;
    const ULOW r = t22 & MASK;
    const ULOW t3 = x_high * y_high;
    U128 result{};
    result.mLow = t1;
    const ULOW div = (q + s) + ((p + r + t) >> mHalfWidth);
    const ULOW mod = (t21 << mHalfWidth) + (t22 << mHalfWidth);
    result.mLow += mod;
    result.mHigh += div;
    result.mHigh += t3;
    result.mOverflow = result.mHigh < t3 ? 1 : 0;
    return result;
}

U128 operator*(ULOW rhs) const {
    U128 result = mult64(mLow, rhs);
    U128 tmp = mult64(mHigh, rhs);
    tmp.mHigh = tmp.mLow;
    tmp.mLow = 0;
    result += tmp;
    result.mSign = !result.is_zero() ? this->mSign : 0;
    return result;
}

U128 div10() const { // Специальный метод деления на 10 для строкового представления числа.
    U128 X = *this;
    const int sign = X.mSign;
    X.mSign = 0;
    ULOW Q = X.mHigh / 10;
    ULOW R = X.mHigh % 10;
    ULOW N = R * (mMaxULOW / 10) + (X.mLow / 10);
    U128 result{.mHigh = Q, .mLow = N};
    U128 E = X - result * 10;
    while (E.mHigh != 0 || E.mLow >= 10) {
        Q = E.mHigh / 10;
        R = E.mHigh % 10;
        N = R * (mMaxULOW / 10) + (E.mLow / 10);
        const U128 tmp {.mHigh = Q, .mLow = N};
        result += tmp;
        E -= tmp * 10;
    }
    result.mSign = sign;
    return result;
}

U128 operator/(ULOW y) const {
    assert(y != 0);
    const U128 X = *this;
    ULOW Q = X.mHigh / y;
    ULOW R = X.mHigh % y;
    ULOW N = R * (mMaxULOW / y) + (X.mLow / y);
    U128 result{.mHigh = Q, .mLow = N, .mSign = X.mSign};
    U128 E = X - result * y; // Ошибка от деления: остаток от деления.
    while (1) {
        Q = E.mHigh / y;
        R = E.mHigh % y;
        N = R * (mMaxULOW / y) + (E.mLow / y);
        const U128 tmp {.mHigh = Q, .mLow = N, .mSign = E.mSign};
        if (tmp.is_zero()) {
            break;
        }
        result += tmp;
        E -= tmp * y;
    }
    if (E.is_nonzero_negative()) {
        result -= get_unit();
        E += U128{.mHigh = 0, .mLow = y};
    }
    return result;
}

U128& operator/=(ULOW y) {
    *this = *this / y;
    return *this;
}

U128 operator/(const U128 other) const {
    U128 X = *this;
    U128 Y = other;
    constexpr U128 ZERO = get_zero();
    constexpr U128 UNIT = get_unit();
    constexpr U128 UNIT_NEG = get_unit_neg();
    if (Y.mHigh == 0) {
        X.mSign ^= Y.mSign;
        U128 result = X / Y.mLow;
        return result;
    }
    const bool make_sign_inverse = X.mSign != Y.mSign;
    X.mSign = make_sign_inverse;
    Y.mSign = 0;
    const ULOW Q = X.mHigh / Y.mHigh;
    const ULOW R = X.mHigh % Y.mHigh;
    const ULOW Delta = mMaxULOW - Y.mLow;
    const U128 DeltaQ = mult64(Delta, Q);
    U128 W1 = U128{.mHigh = R, .mLow = 0} - U128{.mHigh = Q, .mLow = 0};
    W1 = W1 + DeltaQ;
    const ULOW C1 = (Y.mHigh < mMaxULOW) ? Y.mHigh + ULOW(1) : mMaxULOW;
    const ULOW W2 = mMaxULOW - Delta / C1;
    U128 Quotient = W1 / W2;
    Quotient = Quotient / C1;
    U128 result = U128{.mHigh = 0, .mLow = Q} + Quotient;
    assert(result.mHigh == 0);
    result.mSign ^= make_sign_inverse;
    U128 N = Y * result.mLow;
    N.mSign ^= make_sign_inverse;
    assert(!N.is_overflow());
    U128 Error = X - N;
    U128 More = Error - Y;
    bool do_inc = More.is_non_negative();
    bool do_dec = Error.is_nonzero_negative();
    while (do_dec || do_inc) {
        result += (do_inc ? UNIT : (do_dec ? UNIT_NEG : ZERO));
        if (do_dec) {
            Error += Y;
        }
        if (do_inc) {
            Error -= Y;
        }
        More = Error - Y;
        do_inc = More.is_non_negative();
        do_dec = Error.is_nonzero_negative();
    }
    return result;
}

/**
 * Возвращает строковое представление числа.
 */
std::string value() const {
    std::string result{};
    if (this->is_overflow()) {
        result = "Overflow";
        return result;
    }
    U128 X = *this;
    constexpr int multiplier_mod10 = mMaxULOW % 10 + 1;
    while (!X.is_zero()) {
        const int d =
            ((X.mLow % 10) + multiplier_mod10 * (X.mHigh % 10)) % 10;
        result.push_back(DIGITS[d]);
        X = X.div10();
    }
    if (this->is_negative() && !this->is_zero()) {
        result.push_back('-');
    }
    std::reverse(result.begin(), result.end());
    return result.length() != 0 ? result : "0";
}
};

Выводы

Представлены С++ тесты деления 128-битных чисел, использующие Python C API интерфейс для получения эталонных значений.

Приведены шаги настройки среды и проекта с использованием WSL Ubuntu 22.04 и Visual Studio Code.

© Habrahabr.ru