diff --git a/bench/mo_exact_decimal.cpp b/bench/mo_exact_decimal.cpp new file mode 100644 index 0000000000..0554d3ccba --- /dev/null +++ b/bench/mo_exact_decimal.cpp @@ -0,0 +1,120 @@ +/* Copyright 2026 Sirius Contributors. SPDX-License-Identifier: Apache-2.0 */ +#include "numeric/exact_decimal_gpu.hpp" + +#include +#include +#include + +#include +#include +#include + +#include +#include +#include + +using namespace sirius::mo_decimal; +namespace { +constexpr cudf::size_type rows = 1 << 18; + +std::unique_ptr input(decimal_type type, + rmm::cuda_stream_view stream, + rmm::device_async_resource_ref mr) +{ + auto result = make_decimal_column(type, rows, cudf::mask_state::ALL_VALID, stream, mr); + auto view = result->mutable_view(); + std::vector ones(rows, 1); + if (type.bits == 64) { + CUDF_CUDA_TRY(cudaMemcpyAsync(view.data(), + ones.data(), + rows * sizeof(uint64_t), + cudaMemcpyHostToDevice, + stream.value())); + } else { + for (int child = 0; child < 3; ++child) + CUDF_CUDA_TRY(cudaMemsetAsync( + view.child(child).data(), 0, rows * sizeof(uint64_t), stream.value())); + CUDF_CUDA_TRY(cudaMemcpyAsync(view.child(3).data(), + ones.data(), + rows * sizeof(uint64_t), + cudaMemcpyHostToDevice, + stream.value())); + } + stream.synchronize(); + return result; +} + +template +void measure(char const* name, rmm::cuda_stream_view stream, Evaluate evaluate) +{ + // Include allocations, masks, error storage and required stream quiescence. + // Exclude input upload, warm-up and output destruction for every case. + for (int repetition = -2; repetition < 7; ++repetition) { + auto start = std::chrono::steady_clock::now(); + auto result = evaluate(); + stream.synchronize(); + auto ms = + std::chrono::duration(std::chrono::steady_clock::now() - start).count(); + if (repetition >= 0) + std::cout << name << ',' << rows << ',' << repetition << ',' << ms << ',' + << rows * 1000.0 / ms << '\n'; + } +} +} // namespace + +int main() +{ + rmm::cuda_stream stream; + // Match the pooled allocation model of a task resource, and use the same + // bounded caller-owned pool for cuDF controls and checked operations. + rmm::mr::pool_memory_resource pool( + cuda::mr::any_resource{rmm::mr::cuda_memory_resource{}}, + 64u << 20, + 128u << 20); + auto mr = rmm::to_device_async_resource_ref_checked(&pool); + decimal_type narrow{64, 18, 0}, wide{256, 65, 0}; + auto a = input(narrow, stream.view(), mr); + auto b = input(narrow, stream.view(), mr); + auto x = input(wide, stream.view(), mr); + auto y = input(wide, stream.view(), mr); + std::cout << "case,rows,repetition,wall_ms,rows_per_second\n"; + for (auto operation : {decimal_op::add, decimal_op::multiply}) { + auto out = operation == decimal_op::add ? narrow : decimal_type{128, 38, 0}; + auto cudf_op = + operation == decimal_op::add ? cudf::binary_operator::ADD : cudf::binary_operator::MUL; + measure(operation == decimal_op::add ? "cudf_add64" : "cudf_widen_multiply64_to128", + stream.view(), + [&] { + auto output_type = cudf::data_type{ + operation == decimal_op::add ? cudf::type_id::DECIMAL64 : cudf::type_id::DECIMAL128, + 0}; + // cuDF arithmetic runs at the operand width. Widen both operands for + // an equivalent multiply control; changing only the output type wraps. + std::unique_ptr left, right; + if (operation == decimal_op::multiply) { + left = cudf::cast(a->view(), output_type, stream.view(), mr); + right = cudf::cast(b->view(), output_type, stream.view(), mr); + } + auto result = cudf::binary_operation(left ? left->view() : a->view(), + right ? right->view() : b->view(), + cudf_op, + output_type, + stream.view(), + mr); + stream.synchronize(); + return result; + }); + measure(operation == decimal_op::add ? "checked_add64" : "checked_multiply64_to128", + stream.view(), + [&] { + return evaluate_decimal_columns( + operation, a->view(), narrow, b->view(), narrow, out, nullptr, stream.view(), mr); + }); + } + for (auto operation : {decimal_op::add, decimal_op::divide}) + measure( + operation == decimal_op::add ? "checked_add256" : "checked_divide256", stream.view(), [&] { + return evaluate_decimal_columns( + operation, x->view(), wide, y->view(), wide, wide, nullptr, stream.view(), mr); + }); +} diff --git a/cmake/sirius-embedding-tests.cmake b/cmake/sirius-embedding-tests.cmake index 50120436e5..3e8bae102a 100644 --- a/cmake/sirius-embedding-tests.cmake +++ b/cmake/sirius-embedding-tests.cmake @@ -1,5 +1,20 @@ # Embedding control tests need neither a GPU nor a DuckDB runtime. Keep them # independently runnable from the linked ABI and GPU integration tests. +add_executable( + sirius_exact_decimal_unittest test/cpp/exec/execution_evidence_unittest.cpp + test/cpp/embedding/test_exact_decimal.cpp) +target_compile_features(sirius_exact_decimal_unittest PRIVATE cxx_std_20) +target_include_directories( + sirius_exact_decimal_unittest PRIVATE "${CMAKE_SOURCE_DIR}/third_party/catch" + "${CMAKE_CURRENT_SOURCE_DIR}/src") + +add_executable(sirius_exact_decimal_benchmark EXCLUDE_FROM_ALL + bench/mo_exact_decimal.cpp) +target_compile_features(sirius_exact_decimal_benchmark PRIVATE cxx_std_20) +target_include_directories(sirius_exact_decimal_benchmark + PRIVATE "${CMAKE_CURRENT_SOURCE_DIR}/src") +target_link_libraries(sirius_exact_decimal_benchmark PRIVATE sirius_extension) + add_executable( sirius_native_control_unittest test/cpp/exec/execution_evidence_unittest.cpp @@ -52,6 +67,7 @@ add_executable( test/cpp/embedding/test_native_gpu.cpp test/cpp/embedding/test_native_admission.cpp test/cpp/embedding/test_native_result_codec.cpp + test/cpp/embedding/test_exact_decimal_gpu.cpp test/cpp/embedding/test_tae_gpu.cpp src/embedding/c_api.cpp) target_compile_features(sirius_native_gpu_unittest PRIVATE cxx_std_20) diff --git a/docs/mo-exact-decimal-kernels.csv b/docs/mo-exact-decimal-kernels.csv new file mode 100644 index 0000000000..3f2fb7ec39 --- /dev/null +++ b/docs/mo-exact-decimal-kernels.csv @@ -0,0 +1,43 @@ +case,rows,repetition,wall_ms,rows_per_second +cudf_add64,262144,0,0.20445,1.28219e+09 +cudf_add64,262144,1,0.189379,1.38423e+09 +cudf_add64,262144,2,0.197139,1.32974e+09 +cudf_add64,262144,3,0.195489,1.34097e+09 +cudf_add64,262144,4,0.19587,1.33836e+09 +cudf_add64,262144,5,0.195839,1.33857e+09 +cudf_add64,262144,6,0.204669,1.28082e+09 +checked_add64,262144,0,0.02729,9.60586e+09 +checked_add64,262144,1,0.02719,9.64119e+09 +checked_add64,262144,2,0.02758,9.50486e+09 +checked_add64,262144,3,0.02637,9.94099e+09 +checked_add64,262144,4,0.033949,7.7217e+09 +checked_add64,262144,5,0.03364,7.79263e+09 +checked_add64,262144,6,0.03487,7.51775e+09 +cudf_widen_multiply64_to128,262144,0,0.502128,5.22066e+08 +cudf_widen_multiply64_to128,262144,1,0.495698,5.28838e+08 +cudf_widen_multiply64_to128,262144,2,0.499908,5.24384e+08 +cudf_widen_multiply64_to128,262144,3,0.497848,5.26554e+08 +cudf_widen_multiply64_to128,262144,4,0.499428,5.24888e+08 +cudf_widen_multiply64_to128,262144,5,0.499538,5.24773e+08 +cudf_widen_multiply64_to128,262144,6,0.507368,5.16674e+08 +checked_multiply64_to128,262144,0,0.505538,5.18545e+08 +checked_multiply64_to128,262144,1,0.500578,5.23683e+08 +checked_multiply64_to128,262144,2,0.500228,5.24049e+08 +checked_multiply64_to128,262144,3,0.495228,5.2934e+08 +checked_multiply64_to128,262144,4,0.512808,5.11193e+08 +checked_multiply64_to128,262144,5,0.498808,5.25541e+08 +checked_multiply64_to128,262144,6,0.499539,5.24772e+08 +checked_add256,262144,0,0.263789,9.93764e+08 +checked_add256,262144,1,0.262309,9.99371e+08 +checked_add256,262144,2,0.258849,1.01273e+09 +checked_add256,262144,3,0.262029,1.00044e+09 +checked_add256,262144,4,0.262509,9.9861e+08 +checked_add256,262144,5,0.262309,9.99371e+08 +checked_add256,262144,6,0.263589,9.94518e+08 +checked_divide256,262144,0,4.67064,5.61259e+07 +checked_divide256,262144,1,4.66343,5.62127e+07 +checked_divide256,262144,2,4.66613,5.61802e+07 +checked_divide256,262144,3,4.66244,5.62246e+07 +checked_divide256,262144,4,4.65689,5.62916e+07 +checked_divide256,262144,5,4.65525,5.63115e+07 +checked_divide256,262144,6,4.65287,5.63403e+07 diff --git a/docs/mo-exact-decimal-primitives.md b/docs/mo-exact-decimal-primitives.md new file mode 100644 index 0000000000..123cfe6c31 --- /dev/null +++ b/docs/mo-exact-decimal-primitives.md @@ -0,0 +1,118 @@ +# MO exact-decimal primitives + +This is the first native implementation increment for +[MO #28968](https://github.com/matrixorigin/matrixone/issues/28968), following +the approved numeric design revision +`42a89f09a1d168d02b9583cb3ea7b4de6dbb5634` in MO #29449. + +It preserves the existing Decimal64/128 carriers and input/result ownership. +It does not register the MO exact-decimal Substrait family, advertise its +capability, change ordinary DuckDB decimal execution, or admit new MO queries. +The importer is pinned to merged duckdb-substrait #3; its scoped handler is +installed by the following native integration increment. + +## Representation and arithmetic + +`sirius::mo_decimal::decimal_type` stores physical coefficient width, precision +and scale independently. The native wire remains little-endian signed +two's-complement coefficients: 8, 16 and 32 bytes, with separate validity. +Working Decimal256 precision reaches 76; public result precision stops at 65. +In particular, Decimal256(15,2) remains physically 32 bytes. + +Decimal64/128 use existing cuDF fixed-point columns. Decimal256 uses a STRUCT +of signed high64 followed by three unsigned limbs in descending significance; +validity belongs to the parent, not to individual limbs. Input decoding and +result interleaving reuse the existing admitted storage and scratch slab. + +The checked host/device core supports arithmetic, comparisons, proven-safe +negation and checked normal rescaling. It uses fixed 512-bit scratch, never +per-row allocation, floating point or numeric strings. Multiplication rescales +and rounds before checking the final physical and declared domains. Division +consumes the supplied result scale; it does not derive MO types or session +settings. Scale reduction and division round half away from zero; DIV and +modulo truncate. NULL, an inactive execution mask and SELECT zero divisors +produce NULL without arithmetic failure. + +Scalar failure is `decimal_error::out_of_range`; checked cast failure is +`decimal_error::invalid_input`. These are internal errors here, not new public +C statuses. The native integration increment maps the operation-owned status +without parsing messages or inferring from the enclosing aggregate. + +## Ownership and use + +`evaluate_decimal_columns` takes the caller's task stream and reservation-aware +RMM resource. Operand and mask owners must outlive that stream's work. Result +values and the bounded per-row error vector are owned together. Failure first +proves stream quiescence; an unprovable owner is retained and reported through +the existing fatal GPU path. No device-wide synchronization or independent +allocator/pool is installed. The pinned cuDF API requires an explicit result +null count, computed on the supplied stream. + +## Local validation + +Use the frozen Sirius MO Pixi environment in the existing checkout: + +```sh +pixi run --frozen -e mo mo-configure-embedding-sdk +pixi run --frozen -e mo cmake --build build/mo --target \ + sirius_exact_decimal_unittest sirius_native_control_unittest \ + sirius_native_gpu_unittest -j8 +pixi run --frozen -e mo build/mo/extension/sirius/sirius_exact_decimal_unittest +pixi run --frozen -e mo build/mo/extension/sirius/sirius_native_control_unittest +pixi run --frozen -e mo build/mo/extension/sirius/sirius_native_gpu_unittest \ + '[exact_decimal_gpu],[native_gpu],[native_result_gpu]' +``` + +The CPU suite includes independent integer/rational boundary vectors, not +expectations computed by the native implementation. GPU tests check those same +independent full-width results, masks, 64/128-bit controls, sliced wide output, +and constant/NULL/retry input decoding. The existing native suite retains its +1/2/4-worker checks. Small native tests are not SF10/all-22 correctness or +performance evidence. Full capability, MO lowering and public/resource/timing +gates remain required before default cutover. + +## Kernel measurements + +The standalone benchmark uses 262,144 rows, two warm-ups and seven measured +repetitions per operation. One caller-owned RMM pool starts at 64 MiB and is +capped at 128 MiB. Both the cuDF controls and checked kernels use that pool +and one task stream. Timing includes output/error/mask allocation and stream +quiescence, and excludes input upload, warm-ups and output destruction. + +```sh +pixi run --frozen -e mo cmake --build build/mo \ + --target sirius_exact_decimal_benchmark -j8 +pixi run --frozen -e mo build/mo/extension/sirius/sirius_exact_decimal_benchmark \ + > decimal-kernels.csv +``` + +The cuDF multiply control widens both operands to Decimal128 first. Merely +requesting a Decimal128 output still computes at Decimal64 width in the +pinned cuDF and can wrap; the signed physical-domain test uses an independent +host `__int128` oracle. The add control measures cuDF's unchecked arithmetic +on the same inputs. Decimal256 has no equivalent cuDF decimal carrier. + +The CSV retains every sample. These kernel measurements establish the cost +of these primitives; all-22 numeric query timing and the rollout regression +gate belong to the native plan integration and public campaign. + +### Local snapshot + +[Raw repetitions](mo-exact-decimal-kernels.csv), measured on 2026-10-08 from +this increment based on Sirius `e2e2f08f9fd1eaa1253297493df0ab11a6078664`: +RTX 3070, NVIDIA driver 615.71.09, CUDA 13.3.73, GCC 14.4.0, frozen Pixi `mo` +environment, DuckDB `069cc9f9b5be802405797faecc284961b07c70ef`, importer +`dd4cab14b82754ca919633913436cd530f496e00`. Pixi lock SHA-256: +`ed808ec91232d916b49e11abb14aba0dcdd82fa13d908342bc970ce584ce4fef`. + +| Operation | Median wall time (ms) | +| --- | ---: | +| cuDF add64, unchecked | 0.195870 | +| Checked add64 | 0.027580 | +| cuDF widened multiply64-to128 | 0.499538 | +| Checked multiply64-to128 | 0.500228 | +| Checked add256 | 0.262309 | +| Checked divide256 | 4.662440 | + +The checked multiply differs by approximately 0.14% in this sample. Wide +operations report absolute throughput because cuDF has no equivalent carrier. diff --git a/src/CMakeLists.txt b/src/CMakeLists.txt index 375e04b2ea..39ce239035 100644 --- a/src/CMakeLists.txt +++ b/src/CMakeLists.txt @@ -27,6 +27,7 @@ add_subdirectory(io) add_subdirectory(late_mat) add_subdirectory(log) add_subdirectory(memory) +add_subdirectory(numeric) add_subdirectory(op) add_subdirectory(parallel) add_subdirectory(pipeline) diff --git a/src/cuda/scan/native_mo_input.cu b/src/cuda/scan/native_mo_input.cu index 114b568cbb..a6d8003c3f 100644 --- a/src/cuda/scan/native_mo_input.cu +++ b/src/cuda/scan/native_mo_input.cu @@ -1,5 +1,6 @@ /* Copyright 2026 Sirius Contributors. SPDX-License-Identifier: Apache-2.0 */ #include "embedding/native_gpu.hpp" +#include "numeric/exact_decimal_gpu.hpp" #include "pipeline/gpu_stream_quiescence_error.hpp" #include @@ -87,6 +88,27 @@ __global__ void scatter(device_slice const* slices, int32_t const* offsets, uint dst[b] = src[b]; } } +__global__ void decode_wide(device_slice const* slices, + uint64_t* high, + uint64_t* mid_high, + uint64_t* mid_low, + uint64_t* low, + uint32_t* validity) +{ + auto const& d = slices[blockIdx.y]; + for (uint32_t r = blockIdx.x * blockDim.x + threadIdx.x; r < d.rows; + r += blockDim.x * gridDim.x) { + auto row = d.begin + r, out = d.output + r; + bool null = is_null(d, row); + if (null) atomicAnd(validity + out / 32, ~(1u << (out % 32))); + auto source_row = d.column.vector_class == SIRIUS_VECTOR_CONSTANT ? 0 : row; + auto source = null ? nullptr : d.base + d.column.data_offset + uint64_t(source_row) * 32; + uint64_t* limbs[]{low, mid_low, mid_high, high}; + for (int limb = 0; limb < 4; ++limb) + limbs[limb][out] = + null ? 0 : uint64_t(u32(source + limb * 8)) | (uint64_t(u32(source + limb * 8 + 4)) << 32); + } +} cudf::data_type type(sirius_input_column const& c) { using id = cudf::type_id; @@ -104,6 +126,7 @@ cudf::data_type type(sirius_input_column const& c) case 31: return cudf::data_type{id::FLOAT64}; case 32: return cudf::data_type{id::DECIMAL64, -c.scale}; case 33: return cudf::data_type{id::DECIMAL128, -c.scale}; + case 34: return cudf::data_type{id::STRUCT}; case 50: return cudf::data_type{id::TIMESTAMP_DAYS}; case 52: return cudf::data_type{id::TIMESTAMP_MICROSECONDS}; default: return cudf::data_type{id::STRING}; @@ -126,8 +149,17 @@ std::unique_ptr convert_native_input(input_unit const& unit, if (!unit.rows) { if (schema.empty()) owner->columns.push_back(cudf::make_empty_column(cudf::data_type{cudf::type_id::INT8})); - for (auto const& c : schema) - owner->columns.push_back(cudf::make_empty_column(type(c))); + for (auto const& c : schema) { + if (c.oid == 34) + owner->columns.push_back(mo_decimal::make_decimal_column( + {256, static_cast(c.width), static_cast(c.scale)}, + 0, + cudf::mask_state::UNALLOCATED, + stream, + mr)); + else + owner->columns.push_back(cudf::make_empty_column(type(c))); + } return std::make_unique(std::move(owner->columns)); } for (auto const& slice : unit.slices) { @@ -167,7 +199,24 @@ std::unique_ptr convert_native_input(input_unit const& unit, desc.data(), host.data(), desc.size(), cudaMemcpyHostToDevice, stream.value())); auto d = static_cast(desc.data()); dim3 grid(128, host.size()); - if (!input_string_type(schema[c].oid)) { + if (schema[c].oid == 34) { + auto col = mo_decimal::make_decimal_column( + {256, static_cast(schema[c].width), static_cast(schema[c].scale)}, + unit.rows, + cudf::mask_state::ALL_VALID, + stream, + mr); + owner->columns.push_back(std::move(col)); + auto v = owner->columns.back()->mutable_view(); + decode_wide<<>>(d, + v.child(0).data(), + v.child(1).data(), + v.child(2).data(), + v.child(3).data(), + v.null_mask()); + CUDF_CUDA_TRY(cudaGetLastError()); + owner->columns.back()->set_null_count(unit.nulls[c]); + } else if (!input_string_type(schema[c].oid)) { auto col = cudf::make_fixed_width_column( type(schema[c]), unit.rows, cudf::mask_state::ALL_VALID, stream, mr); owner->columns.push_back(std::move(col)); diff --git a/src/embedding/control.cpp b/src/embedding/control.cpp index 32f69cc588..32721bb770 100644 --- a/src/embedding/control.cpp +++ b/src/embedding/control.cpp @@ -126,6 +126,8 @@ owned_column copy_column(sirius_column const& c, std::size_t& bytes) if (c.reserved || c.nullable > 1) throw failure(SIRIUS_INVALID_ARGUMENT, "invalid column"); sirius_input_column scalar{c.oid, c.width, c.scale, c.nullable}; validate_input_schema({&scalar, 1}); + if (c.oid == 34 && c.width > 65) + throw failure(SIRIUS_INVALID_ARGUMENT, "public decimal result precision exceeds 65"); owned_column out{c.oid, c.width, c.scale, diff --git a/src/embedding/input.cpp b/src/embedding/input.cpp index 55f5bbaeb1..69bd3a5074 100644 --- a/src/embedding/input.cpp +++ b/src/embedding/input.cpp @@ -63,6 +63,7 @@ std::size_t input_element_size(uint32_t oid) case 32: case 52: return 8; case 33: return 16; + case 34: return 32; default: if (input_string_type(oid)) return 24; } @@ -74,10 +75,13 @@ void validate_input_schema(std::span columns) for (auto const& c : columns) { (void)input_element_size(c.oid); require(c.nullable <= 1, "invalid native column nullability"); - if (c.oid == 32 || c.oid == 33) - require( - c.width >= 1 && c.width <= (c.oid == 32 ? 18 : 38) && c.scale >= 0 && c.scale <= c.width, - "invalid native decimal precision or scale"); + if (c.oid == 32 || c.oid == 33 || c.oid == 34) + require(c.width >= 1 && + c.width <= (c.oid == 32 ? 18 + : c.oid == 33 ? 38 + : 76) && + c.scale >= 0 && c.scale <= c.width, + "invalid native decimal precision or scale"); } } input_batch::~input_batch() diff --git a/src/embedding/result_codec.cpp b/src/embedding/result_codec.cpp index a1bad41b57..549584982f 100644 --- a/src/embedding/result_codec.cpp +++ b/src/embedding/result_codec.cpp @@ -1,6 +1,7 @@ /* Copyright 2026 Sirius Contributors. SPDX-License-Identifier: Apache-2.0 */ #include "embedding/result_codec.hpp" +#include "numeric/exact_decimal_gpu.hpp" #include "pipeline/gpu_stream_quiescence_error.hpp" #include @@ -72,6 +73,7 @@ cudf::data_type expected(owned_column const& c) case 31: return cudf::data_type{id::FLOAT64}; case 32: return cudf::data_type{id::DECIMAL64, -c.scale}; case 33: return cudf::data_type{id::DECIMAL128, -c.scale}; + case 34: return cudf::data_type{id::STRUCT}; case 50: return cudf::data_type{id::TIMESTAMP_DAYS}; case 52: return cudf::data_type{id::TIMESTAMP_MICROSECONDS}; default: return cudf::data_type{id::STRING}; @@ -115,7 +117,11 @@ result_slice_layout size_native_result(cudf::table_view table, if (table.num_columns() != static_cast(schema.size())) throw failure(SIRIUS_EXECUTION_FAILED, "native result column count mismatch"); for (std::size_t c = 0; c < schema.size(); ++c) - if (table.column(c).type() != expected(schema[c])) + if (table.column(c).type() != expected(schema[c]) || + (schema[c].oid == 34 && + !mo_decimal::decimal_column_matches( + table.column(c), + {256, static_cast(schema[c].width), static_cast(schema[c].scale)}))) throw failure(SIRIUS_UNSUPPORTED, "native result physical type does not match MO contract"); auto size = [&](uint32_t rows) { std::size_t position{}; @@ -151,7 +157,39 @@ void encode_native_result(cudf::table_view table, auto& vector = output.columns[c]; vector = column_layout(col, schema[c], begin, layout.rows, position, scratch, stream); auto width = input_element_size(schema[c].oid); - if (!input_string_type(schema[c].oid)) { + if (schema[c].oid == 34) { + // Reuse the admitted scratch slab; never allocate a second complete + // result while interleaving the four fixed-width device children. + if (scratch.storage->size() < 32) + throw failure(SIRIUS_RESOURCE_EXHAUSTED, "wide result scratch too small"); + auto chunk = std::min(1024, scratch.storage->size() / 32); + for (uint32_t first = 0; first < layout.rows; first += chunk) { + auto count = std::min(chunk, layout.rows - first); + for (int limb = 0; limb < 4; ++limb) { + auto child = col.child(3 - limb); + download(*scratch.storage, + limb * count * 8, + child.head() + child.offset() + col.offset() + begin + first, + count * 8, + stream); + } + auto destination = vector.data_offset + first * 32ULL; + // Visit output blocks once per chunk. A per-row write would rescan + // every output block for each coefficient in a large result slice. + output.storage->visit([&](std::size_t base, std::span block) { + auto low = std::max(base, destination); + auto high = std::min(base + block.size(), destination + count * 32); + for (auto offset = low; offset < high;) { + auto relative = offset - destination; + auto row = relative / 32, limb = relative % 32 / 8, byte = relative % 8; + auto length = std::min(8 - byte, high - offset); + scratch.storage->read(limb * count * 8 + row * 8 + byte, + block.subspan(offset - base, length)); + offset += length; + } + }); + } + } else if (!input_string_type(schema[c].oid)) { download(*output.storage, vector.data_offset, col.head() + (static_cast(col.offset()) + begin) * width, diff --git a/src/numeric/CMakeLists.txt b/src/numeric/CMakeLists.txt new file mode 100644 index 0000000000..c88cba83f6 --- /dev/null +++ b/src/numeric/CMakeLists.txt @@ -0,0 +1 @@ +target_sources(sirius_objects PRIVATE exact_decimal_gpu.cu) diff --git a/src/numeric/exact_decimal.hpp b/src/numeric/exact_decimal.hpp new file mode 100644 index 0000000000..9e31c3a2c2 --- /dev/null +++ b/src/numeric/exact_decimal.hpp @@ -0,0 +1,399 @@ +/* Copyright 2026 Sirius Contributors. SPDX-License-Identifier: Apache-2.0 */ +#pragma once + +#include + +#ifdef __CUDACC__ +#define SIRIUS_DECIMAL_HD __host__ __device__ +#else +#define SIRIUS_DECIMAL_HD +#endif + +namespace sirius::mo_decimal { +// MO physical width is independent of SQL precision. This descriptor does not +// change ordinary DuckDB/cuDF decimal type inference. +struct decimal_type { + uint16_t bits; + uint8_t precision, scale; + SIRIUS_DECIMAL_HD bool valid() const + { + auto limit = bits == 64 ? 18 : bits == 128 ? 38 : bits == 256 ? 76 : 0; + return precision > 0 && precision <= limit && scale <= precision; + } + SIRIUS_DECIMAL_HD uint32_t bytes() const { return bits / 8; } +}; + +struct coefficient { + uint32_t words[8]{}; // Little-endian, signed two's complement, sign extended. + SIRIUS_DECIMAL_HD bool negative() const { return (words[7] >> 31) != 0; } +}; +enum class decimal_error : uint8_t { none, out_of_range, invalid_input }; +enum class decimal_op : uint8_t { + add, + subtract, + multiply, + divide, + integer_divide, + modulo, + negate, + cast, + equal, + not_equal, + less, + less_equal, + greater, + greater_equal +}; +struct decimal_result { + coefficient value{}; + decimal_error error{decimal_error::none}; + bool valid{false}; +}; + +// Same-scale 64-bit operations need no wide scratch. Unsigned arithmetic makes +// physical overflow checks defined, including signed-minimum negation. +SIRIUS_DECIMAL_HD inline decimal_result evaluate_same_scale64(decimal_op op, + uint64_t left, + uint64_t right, + uint64_t precision_bound, + bool left_valid, + bool right_valid, + bool active) +{ + bool unary = op == decimal_op::cast || op == decimal_op::negate; + if (!active || !left_valid || (!unary && !right_valid)) return {}; + bool ln = (left >> 63) != 0, rn = (right >> 63) != 0; + uint64_t value = left; + auto error = op == decimal_op::cast ? decimal_error::invalid_input : decimal_error::out_of_range; + if (op >= decimal_op::equal && op <= decimal_op::greater_equal) { + bool less = ln != rn ? ln : left < right; + bool answer = op == decimal_op::equal ? left == right + : op == decimal_op::not_equal ? left != right + : op == decimal_op::less ? less + : op == decimal_op::less_equal ? less || left == right + : op == decimal_op::greater ? !less && left != right + : !less; + decimal_result result; + result.valid = true; + result.value.words[0] = answer; + return result; + } + if (op == decimal_op::add) { + value = left + right; + if (ln == rn && bool(value >> 63) != ln) return {{}, error, false}; + } else if (op == decimal_op::subtract) { + value = left - right; + if (ln != rn && bool(value >> 63) != ln) return {{}, error, false}; + } else if (op == decimal_op::negate) { + if (left == (uint64_t(1) << 63)) return {{}, error, false}; + value = ~left + 1; + } else if (op != decimal_op::cast) + return {{}, decimal_error::invalid_input, false}; + auto negative = (value >> 63) != 0; + auto magnitude = negative ? ~value + 1 : value; + if (magnitude >= precision_bound) return {{}, error, false}; + decimal_result result; + result.valid = true; + result.value.words[0] = value; + result.value.words[1] = value >> 32; + for (int i = 2; i < 8; ++i) + result.value.words[i] = negative ? UINT32_MAX : 0; + return result; +} + +SIRIUS_DECIMAL_HD inline coefficient load_coefficient(uint8_t const* data, uint32_t bytes) +{ + coefficient result; + auto extension = (data[bytes - 1] & 128) ? UINT32_MAX : 0; + for (uint32_t i = 0; i < 8; ++i) + result.words[i] = extension; + for (uint32_t i = 0; i < bytes / 4; ++i) { + auto p = data + i * 4; + result.words[i] = + uint32_t(p[0]) | uint32_t(p[1]) << 8 | uint32_t(p[2]) << 16 | uint32_t(p[3]) << 24; + } + return result; +} +SIRIUS_DECIMAL_HD inline void store_coefficient(coefficient const& value, + uint8_t* data, + uint32_t bytes) +{ + for (uint32_t i = 0; i < bytes; ++i) + data[i] = value.words[i / 4] >> (i % 4 * 8); +} + +namespace decimal_detail { +// Fixed scratch only: no host/device heap, floating point or numeric strings. +struct magnitude { + uint32_t words[16]{}; + SIRIUS_DECIMAL_HD bool zero() const + { + for (auto w : words) + if (w) return false; + return true; + } +}; +SIRIUS_DECIMAL_HD inline int compare(magnitude const& a, magnitude const& b) +{ + for (int i = 15; i >= 0; --i) { + if (a.words[i] < b.words[i]) return -1; + if (a.words[i] > b.words[i]) return 1; + } + return 0; +} +SIRIUS_DECIMAL_HD inline bool add(magnitude& a, magnitude const& b) +{ + uint64_t carry = 0; + for (int i = 0; i < 16; ++i) { + auto sum = uint64_t(a.words[i]) + b.words[i] + carry; + a.words[i] = sum; + carry = sum >> 32; + } + return carry != 0; +} +SIRIUS_DECIMAL_HD inline void subtract(magnitude& a, magnitude const& b) +{ + uint64_t borrow = 0; + for (int i = 0; i < 16; ++i) { + auto rhs = uint64_t(b.words[i]) + borrow; + auto lhs = uint64_t(a.words[i]); + a.words[i] = lhs - rhs; + borrow = lhs < rhs; + } +} +SIRIUS_DECIMAL_HD inline bool multiply_small(magnitude& a, uint32_t factor) +{ + uint64_t carry = 0; + for (int i = 0; i < 16; ++i) { + auto value = uint64_t(a.words[i]) * factor + carry; + a.words[i] = value; + carry = value >> 32; + } + return carry != 0; +} +SIRIUS_DECIMAL_HD inline bool scale_up(magnitude& a, uint32_t digits) +{ + for (uint32_t i = 0; i < digits; ++i) + if (multiply_small(a, 10)) return false; + return true; +} +SIRIUS_DECIMAL_HD inline magnitude power10(uint32_t digits) +{ + magnitude result; + result.words[0] = 1; + scale_up(result, digits); // Callers bound digits to at most 152. + return result; +} +SIRIUS_DECIMAL_HD inline magnitude absolute(coefficient const& c) +{ + magnitude result; + uint64_t carry = c.negative() ? 1 : 0; + for (int i = 0; i < 8; ++i) { + auto word = c.negative() ? ~c.words[i] : c.words[i]; + auto value = uint64_t(word) + carry; + result.words[i] = value; + carry = value >> 32; + } + return result; +} +SIRIUS_DECIMAL_HD inline magnitude multiply(magnitude const& a, magnitude const& b) +{ + magnitude result; + for (int i = 0; i < 8; ++i) { + uint64_t carry = 0; + for (int j = 0; j < 8; ++j) { + auto value = uint64_t(a.words[i]) * b.words[j] + result.words[i + j] + carry; + result.words[i + j] = value; + carry = value >> 32; + } + result.words[i + 8] = carry; + } + return result; +} +SIRIUS_DECIMAL_HD inline bool shift_left(magnitude& a, uint32_t bit) +{ + uint32_t carry = bit; + for (int i = 0; i < 16; ++i) { + auto next = a.words[i] >> 31; + a.words[i] = (a.words[i] << 1) | carry; + carry = next; + } + return carry != 0; +} +SIRIUS_DECIMAL_HD inline void divide(magnitude const& numerator, + magnitude const& denominator, + magnitude& quotient, + magnitude& remainder) +{ + quotient = {}; + remainder = {}; + for (int bit = 511; bit >= 0; --bit) { + auto carry = shift_left(remainder, (numerator.words[bit / 32] >> (bit % 32)) & 1); + if (carry || compare(remainder, denominator) >= 0) { + subtract(remainder, denominator); + quotient.words[bit / 32] |= uint32_t(1) << (bit % 32); + } + } +} +SIRIUS_DECIMAL_HD inline void round(magnitude& quotient, + magnitude const& remainder, + magnitude const& denominator) +{ + magnitude half; + for (int i = 0; i < 16; ++i) + half.words[i] = (denominator.words[i] >> 1) | (i == 15 ? 0 : denominator.words[i + 1] << 31); + magnitude one; + one.words[0] = 1; + if (denominator.words[0] & 1) add(half, one); + if (compare(remainder, half) >= 0) add(quotient, one); +} +SIRIUS_DECIMAL_HD inline bool rescale(magnitude& value, int from, int to, bool rounding) +{ + if (to >= from) return scale_up(value, to - from); + auto denominator = power10(from - to); + magnitude quotient, remainder; + divide(value, denominator, quotient, remainder); + if (rounding) round(quotient, remainder, denominator); + value = quotient; + return true; +} +SIRIUS_DECIMAL_HD inline decimal_result finish(magnitude const& value, + bool negative, + decimal_type type, + decimal_error error, + bool precision_check = true) +{ + magnitude physical; + physical.words[(type.bits - 1) / 32] = uint32_t(1) << ((type.bits - 1) % 32); + auto comparison = compare(value, physical); + if (comparison > 0 || (comparison == 0 && !negative) || + (precision_check && compare(value, power10(type.precision)) >= 0)) + return {{}, error, false}; + decimal_result result; + result.valid = true; + uint64_t carry = negative && !value.zero() ? 1 : 0; + for (int i = 0; i < 8; ++i) { + auto word = negative && !value.zero() ? ~value.words[i] : value.words[i]; + auto v = uint64_t(word) + carry; + result.value.words[i] = v; + carry = v >> 32; + } + return result; +} +} // namespace decimal_detail + +// Descriptor validation is repeated by the host launch wrapper. This guard +// also keeps direct host/device consumers total for malformed descriptors. +SIRIUS_DECIMAL_HD inline decimal_result evaluate_decimal(decimal_op op, + coefficient const& left, + decimal_type left_type, + bool left_valid, + coefficient const& right, + decimal_type right_type, + bool right_valid, + decimal_type output_type, + bool active = true) +{ + using namespace decimal_detail; + auto error = op == decimal_op::cast ? decimal_error::invalid_input : decimal_error::out_of_range; + if (!left_type.valid() || !right_type.valid() || !output_type.valid() || + static_cast(op) > static_cast(decimal_op::greater_equal)) + return {{}, decimal_error::invalid_input, false}; + bool unary = op == decimal_op::cast || op == decimal_op::negate; + if (!active || !left_valid || (!unary && !right_valid)) return {}; + if (left_type.bits == 64 && right_type.bits == 64 && output_type.bits == 64 && + left_type.scale == right_type.scale && left_type.scale == output_type.scale && + (op == decimal_op::add || op == decimal_op::subtract || op == decimal_op::negate || + op == decimal_op::cast || op >= decimal_op::equal)) { + uint64_t bound = 1; + for (int i = 0; i < output_type.precision; ++i) + bound *= 10; + return evaluate_same_scale64(op, + uint64_t(left.words[0]) | (uint64_t(left.words[1]) << 32), + uint64_t(right.words[0]) | (uint64_t(right.words[1]) << 32), + bound, + left_valid, + right_valid, + active); + } + auto a = absolute(left), b = absolute(right); + bool negative = left.negative(); + int scale = left_type.scale; + switch (op) { + case decimal_op::cast: break; + case decimal_op::negate: negative = !negative; break; + case decimal_op::multiply: + a = multiply(a, b); + scale += right_type.scale; + negative = left.negative() != right.negative(); + break; + case decimal_op::divide: + case decimal_op::integer_divide: { + if (b.zero()) return {}; + auto target_scale = op == decimal_op::integer_divide ? 0 : output_type.scale; + int exponent = target_scale + int(right_type.scale) - int(left_type.scale); + if (exponent >= 0) { + // If this overflows 512 bits, division by any 256-bit coefficient + // cannot produce a representable signed 256-bit result. + if (!scale_up(a, exponent)) return {{}, error, false}; + } else if (!scale_up(b, -exponent)) + return {{}, error, false}; + magnitude quotient, remainder; + divide(a, b, quotient, remainder); + if (op == decimal_op::divide) round(quotient, remainder, b); + negative = left.negative() != right.negative(); + if (op == decimal_op::integer_divide) + return finish(quotient, negative, {64, 18, 0}, error, false); + return finish(quotient, negative, output_type, error); + } + default: { + scale = left_type.scale > right_type.scale ? left_type.scale : right_type.scale; + if (!scale_up(a, scale - left_type.scale) || !scale_up(b, scale - right_type.scale)) + return {{}, error, false}; + if (op == decimal_op::modulo) { + if (b.zero()) return {}; + magnitude quotient, remainder; + divide(a, b, quotient, remainder); + a = remainder; + break; + } + auto rhs_negative = right.negative() != (op == decimal_op::subtract); + auto cmp = compare(a, b); + if (op >= decimal_op::equal) { + if (a.zero() && b.zero()) + cmp = 0; + else if (left.negative() != right.negative()) + cmp = left.negative() ? -1 : 1; + else if (left.negative()) + cmp = -cmp; + bool answer = op == decimal_op::equal ? cmp == 0 + : op == decimal_op::not_equal ? cmp != 0 + : op == decimal_op::less ? cmp < 0 + : op == decimal_op::less_equal ? cmp <= 0 + : op == decimal_op::greater ? cmp > 0 + : cmp >= 0; + decimal_result result; + result.valid = true; + result.value.words[0] = answer; + return result; + } + if (op != decimal_op::add && op != decimal_op::subtract) + return {{}, decimal_error::invalid_input, false}; + if (left.negative() == rhs_negative) { + if (add(a, b)) return {{}, error, false}; + } else if (cmp >= 0) + subtract(a, b); + else { + subtract(b, a); + a = b; + negative = rhs_negative; + } + break; + } + } + if (!rescale(a, scale, output_type.scale, op != decimal_op::modulo)) return {{}, error, false}; + return finish(a, negative, output_type, error); +} +} // namespace sirius::mo_decimal + +#undef SIRIUS_DECIMAL_HD diff --git a/src/numeric/exact_decimal_gpu.cu b/src/numeric/exact_decimal_gpu.cu new file mode 100644 index 0000000000..f636f9b483 --- /dev/null +++ b/src/numeric/exact_decimal_gpu.cu @@ -0,0 +1,319 @@ +/* Copyright 2026 Sirius Contributors. SPDX-License-Identifier: Apache-2.0 */ +#include "numeric/exact_decimal_gpu.hpp" +#include "pipeline/gpu_stream_quiescence_error.hpp" + +#include +#include +#include +#include + +#include +#include + +namespace sirius::mo_decimal { +namespace { +struct operand { + uint8_t const* data{}; + uint64_t const* limbs[4]{}; + cudf::bitmask_type const* mask{}; + cudf::size_type offset{}; + uint32_t bytes{}; + __device__ bool valid(cudf::size_type row) const + { + auto bit = offset + row; + return !mask || ((mask[bit / 32] >> (bit % 32)) & 1); + } + __device__ coefficient value(cudf::size_type row) const + { + if (bytes != 32) return load_coefficient(data + uint64_t(row + offset) * bytes, bytes); + coefficient result; + for (int i = 0; i < 4; ++i) { + auto value = limbs[3 - i][row + offset]; + result.words[2 * i] = value; + result.words[2 * i + 1] = value >> 32; + } + return result; + } +}; +operand view(cudf::column_view const& column, decimal_type type) +{ + operand result; + result.bytes = type.bytes(); + result.offset = column.offset(); + result.mask = column.null_mask(); + if (type.bits != 256) + result.data = column.head(); + else + for (int i = 0; i < 4; ++i) + result.limbs[i] = column.child(i).head() + column.child(i).offset(); + return result; +} +__global__ void evaluate_fast64(decimal_op op, + operand left, + operand right, + cudf::size_type rows, + uint64_t bound, + uint8_t* output, + cudf::bitmask_type* validity, + uint8_t* errors, + uint32_t* valid_count, + bool const* active) +{ + for (int64_t base = blockIdx.x * blockDim.x; base < rows; base += blockDim.x * gridDim.x) { + auto row = base + threadIdx.x; + decimal_result result; + if (row < rows) { + bool enabled = !active || active[row]; + bool lv = enabled && left.valid(row), rv = enabled && right.valid(row); + uint64_t a = lv ? reinterpret_cast(left.data)[row + left.offset] : 0; + uint64_t b = rv ? reinterpret_cast(right.data)[row + right.offset] : 0; + result = evaluate_same_scale64(op, a, b, bound, lv, rv, enabled); + errors[row] = static_cast(result.error); + if (op >= decimal_op::equal) + output[row] = result.value.words[0]; + else + reinterpret_cast(output)[row] = + uint64_t(result.value.words[0]) | (uint64_t(result.value.words[1]) << 32); + } + auto bits = __ballot_sync(0xffffffff, row < rows && result.valid); + if ((threadIdx.x % 32) == 0 && row < rows) { + validity[row / 32] = bits; + atomicAdd(valid_count, static_cast(__popc(bits))); + } + } +} +__global__ void evaluate(decimal_op op, + operand left, + decimal_type left_type, + operand right, + decimal_type right_type, + decimal_type output_type, + cudf::size_type rows, + uint8_t* output, + uint64_t* high, + uint64_t* mid_high, + uint64_t* mid_low, + uint64_t* low, + cudf::bitmask_type* validity, + uint8_t* errors, + bool const* active) +{ + for (int64_t row = blockIdx.x * blockDim.x + threadIdx.x; row < rows; + row += blockDim.x * gridDim.x) { + bool enabled = !active || active[row]; + bool lv = enabled && left.valid(row), rv = enabled && right.valid(row); + coefficient a{}, b{}; + if (lv) a = left.value(row); + if (rv) b = right.value(row); + auto result = evaluate_decimal(op, a, left_type, lv, b, right_type, rv, output_type, enabled); + errors[row] = static_cast(result.error); + if (result.valid) atomicOr(validity + row / 32, uint32_t(1) << (row % 32)); + if (op >= decimal_op::equal) + output[row] = result.value.words[0]; + else if (op == decimal_op::integer_divide) + store_coefficient(result.value, output + uint64_t(row) * 8, 8); + else if (output_type.bits != 256) + store_coefficient( + result.value, output + uint64_t(row) * output_type.bytes(), output_type.bytes()); + else { + uint64_t* limbs[]{low, mid_low, mid_high, high}; + for (int i = 0; i < 4; ++i) + limbs[i][row] = + uint64_t(result.value.words[i * 2]) | (uint64_t(result.value.words[i * 2 + 1]) << 32); + } + } +} +} // namespace + +std::unique_ptr make_decimal_column(decimal_type type, + cudf::size_type rows, + cudf::mask_state mask, + rmm::cuda_stream_view stream, + rmm::device_async_resource_ref mr) +{ + if (!type.valid() || rows < 0 || mask == cudf::mask_state::UNINITIALIZED) + throw std::invalid_argument("invalid MO exact decimal descriptor or mask state"); + if (type.bits != 256) + return cudf::make_fixed_width_column( + cudf::data_type{type.bits == 64 ? cudf::type_id::DECIMAL64 : cudf::type_id::DECIMAL128, + -static_cast(type.scale)}, + rows, + mask, + stream, + mr); + std::vector> children; + for (int i = 0; i < 4; ++i) + children.push_back(cudf::make_numeric_column( + cudf::data_type{i == 0 ? cudf::type_id::INT64 : cudf::type_id::UINT64}, + rows, + cudf::mask_state::UNALLOCATED, + stream, + mr)); + auto validity = cudf::create_null_mask(rows, mask, stream, mr); + auto nulls = mask == cudf::mask_state::ALL_NULL ? rows : 0; + return cudf::create_structs_hierarchy( + rows, std::move(children), nulls, std::move(validity), stream, mr); +} +bool decimal_column_matches(cudf::column_view const& column, decimal_type type) +{ + if (!type.valid()) return false; + if (type.bits != 256) + return column.type() == + cudf::data_type{type.bits == 64 ? cudf::type_id::DECIMAL64 : cudf::type_id::DECIMAL128, + -static_cast(type.scale)}; + if (column.type().id() != cudf::type_id::STRUCT || column.num_children() != 4) return false; + for (int i = 0; i < 4; ++i) { + auto child = column.child(i); + if (child.type().id() != (i == 0 ? cudf::type_id::INT64 : cudf::type_id::UINT64) || + child.nullable() || child.size() < column.offset() + column.size()) + return false; + } + return true; +} +decimal_column_result evaluate_decimal_columns(decimal_op op, + cudf::column_view left, + decimal_type left_type, + cudf::column_view right, + decimal_type right_type, + decimal_type output_type, + bool const* active, + rmm::cuda_stream_view stream, + rmm::device_async_resource_ref mr) +{ + if (!decimal_column_matches(left, left_type) || !decimal_column_matches(right, right_type) || + !output_type.valid() || left.size() != right.size() || + static_cast(op) > static_cast(decimal_op::greater_equal)) + throw std::invalid_argument("invalid MO exact decimal operands"); + auto rows = left.size(); + // These signatures contain the entire signed 64-bit physical input domain, + // not merely sampled values. No rounding or physical/precision overflow is + // possible, so retain the existing cuDF implementation rather than emulate it. + bool safe_add = (op == decimal_op::add || op == decimal_op::subtract) && left_type.bits == 64 && + right_type.bits == 64 && output_type.bits == 128 && output_type.precision >= 20 && + left_type.scale == right_type.scale && left_type.scale == output_type.scale; + bool safe_multiply = op == decimal_op::multiply && left_type.bits == 64 && + right_type.bits == 64 && output_type.bits == 128 && + output_type.precision == 38 && + output_type.scale == left_type.scale + right_type.scale; + if (!active && (safe_add || safe_multiply)) { + struct widened_owner { + decimal_column_result result; + std::unique_ptr left, right; + }; + auto owner = std::make_unique( + widened_owner{{nullptr, rmm::device_uvector(rows, stream, mr)}, nullptr, nullptr}); + try { + auto operation = op == decimal_op::add ? cudf::binary_operator::ADD + : op == decimal_op::subtract ? cudf::binary_operator::SUB + : cudf::binary_operator::MUL; + owner->left = cudf::cast( + left, + cudf::data_type{cudf::type_id::DECIMAL128, -static_cast(left_type.scale)}, + stream, + mr); + owner->right = cudf::cast( + right, + cudf::data_type{cudf::type_id::DECIMAL128, -static_cast(right_type.scale)}, + stream, + mr); + owner->result.values = cudf::binary_operation( + owner->left->view(), + owner->right->view(), + operation, + cudf::data_type{cudf::type_id::DECIMAL128, -static_cast(output_type.scale)}, + stream, + mr); + if (rows) + CUDF_CUDA_TRY(cudaMemsetAsync(owner->result.errors.data(), 0, rows, stream.value())); + // cuDF may enqueue readers of the widened temporaries. Do not let their + // owner disappear before those reads and result writes are quiescent. + CUDF_CUDA_TRY(cudaStreamSynchronize(stream.value())); + return std::move(owner->result); + } catch (...) { + auto error = std::current_exception(); + if (cudaStreamSynchronize(stream.value()) != cudaSuccess) { + (void)owner.release(); + throw pipeline::gpu_stream_quiescence_error( + "exact decimal fast path could not prove quiescence"); + } + std::rethrow_exception(error); + } + } + std::unique_ptr output; + if (op >= decimal_op::equal || op == decimal_op::integer_divide) + output = cudf::make_numeric_column( + cudf::data_type{op >= decimal_op::equal ? cudf::type_id::BOOL8 : cudf::type_id::INT64}, + rows, + cudf::mask_state::ALL_NULL, + stream, + mr); + else + output = make_decimal_column(output_type, rows, cudf::mask_state::ALL_NULL, stream, mr); + struct kernel_owner { + decimal_column_result result; + rmm::device_uvector count; + }; + auto owner = std::make_unique( + kernel_owner{{std::move(output), rmm::device_uvector(rows, stream, mr)}, + rmm::device_uvector(1, stream, mr)}); + try { + if (rows) { + auto ov = owner->result.values->mutable_view(); + bool wide = ov.type().id() == cudf::type_id::STRUCT; + bool fast64 = left_type.bits == 64 && right_type.bits == 64 && output_type.bits == 64 && + left_type.scale == right_type.scale && left_type.scale == output_type.scale && + (op == decimal_op::add || op == decimal_op::subtract || + op == decimal_op::negate || op == decimal_op::cast || op >= decimal_op::equal); + if (fast64) { + CUDF_CUDA_TRY(cudaMemsetAsync(owner->count.data(), 0, sizeof(uint32_t), stream.value())); + uint64_t bound = 1; + for (int i = 0; i < output_type.precision; ++i) + bound *= 10; + evaluate_fast64<<<128, 256, 0, stream.value()>>>(op, + view(left, left_type), + view(right, right_type), + rows, + bound, + ov.data(), + ov.null_mask(), + owner->result.errors.data(), + owner->count.data(), + active); + } else + evaluate<<<128, 256, 0, stream.value()>>>(op, + view(left, left_type), + left_type, + view(right, right_type), + right_type, + output_type, + rows, + wide ? nullptr : ov.data(), + wide ? ov.child(0).data() : nullptr, + wide ? ov.child(1).data() : nullptr, + wide ? ov.child(2).data() : nullptr, + wide ? ov.child(3).data() : nullptr, + ov.null_mask(), + owner->result.errors.data(), + active); + CUDF_CUDA_TRY(cudaGetLastError()); + if (fast64) { + uint32_t count; + CUDF_CUDA_TRY(cudaMemcpyAsync( + &count, owner->count.data(), sizeof(count), cudaMemcpyDeviceToHost, stream.value())); + CUDF_CUDA_TRY(cudaStreamSynchronize(stream.value())); + owner->result.values->set_null_count(rows - count); + } else + owner->result.values->set_null_count(cudf::null_count(ov.null_mask(), 0, rows, stream)); + } + return std::move(owner->result); + } catch (...) { + auto failure = std::current_exception(); + if (cudaStreamSynchronize(stream.value()) != cudaSuccess) { + (void)owner.release(); + throw pipeline::gpu_stream_quiescence_error( + "exact decimal kernel could not prove quiescence"); + } + std::rethrow_exception(failure); + } +} +} // namespace sirius::mo_decimal diff --git a/src/numeric/exact_decimal_gpu.hpp b/src/numeric/exact_decimal_gpu.hpp new file mode 100644 index 0000000000..80bb3769b3 --- /dev/null +++ b/src/numeric/exact_decimal_gpu.hpp @@ -0,0 +1,38 @@ +/* Copyright 2026 Sirius Contributors. SPDX-License-Identifier: Apache-2.0 */ +#pragma once +#include "numeric/exact_decimal.hpp" + +#include +#include +#include + +#include + +#include + +namespace sirius::mo_decimal { +// Decimal64/128 retain cuDF's existing carriers. Decimal256 uses the numeric +// ordering [signed high64, unsigned next64, unsigned next64, unsigned low64]. +std::unique_ptr make_decimal_column(decimal_type type, + cudf::size_type rows, + cudf::mask_state mask, + rmm::cuda_stream_view stream, + rmm::device_async_resource_ref mr); +bool decimal_column_matches(cudf::column_view const& column, decimal_type type); + +struct decimal_column_result { + std::unique_ptr values; + rmm::device_uvector errors; +}; +// Caller retains operands/mask and the returned buffers until its task stream +// is quiescent. All allocations use that stream's reservation-aware resource. +decimal_column_result evaluate_decimal_columns(decimal_op op, + cudf::column_view left, + decimal_type left_type, + cudf::column_view right, + decimal_type right_type, + decimal_type output_type, + bool const* active_rows, + rmm::cuda_stream_view stream, + rmm::device_async_resource_ref mr); +} // namespace sirius::mo_decimal diff --git a/substrait b/substrait index b131fd6659..dd4cab14b8 160000 --- a/substrait +++ b/substrait @@ -1 +1 @@ -Subproject commit b131fd6659fe5fa0c6e39c52059e29ac903b2d59 +Subproject commit dd4cab14b82754ca919633913436cd530f496e00 diff --git a/test/cpp/embedding/exact_decimal_wide_cases.hpp b/test/cpp/embedding/exact_decimal_wide_cases.hpp new file mode 100644 index 0000000000..fa77b1a5ac --- /dev/null +++ b/test/cpp/embedding/exact_decimal_wide_cases.hpp @@ -0,0 +1,783 @@ +/* Copyright 2026 Sirius Contributors. SPDX-License-Identifier: Apache-2.0 */ +// Generated with an independent Python integer/rational oracle; no native +// implementation is used to compute these expected coefficients. +struct wide_case { + decimal_op op; + char const* a; + uint8_t sa; + char const* b; + uint8_t sb; + uint8_t so; + char const* expected; + bool overflow; +}; +static constexpr wide_case wide_cases[]{ + {decimal_op::add, + "-7481905979689374862490126991986085102344607008417686552012812223529", + 34, + "-18090579963602448256349102859279172589208048057396334649679720663", + 55, + 7, + "-7481905979689374862490145082566048704793", + false}, + {decimal_op::subtract, + "43906875401718534596961826478519352409127109597548431056984155153232", + 1, + "-40792453652211172261959882904744283005917976142778270196936455287", + 65, + 5, + "", + true}, + {decimal_op::multiply, + "36354419922597055300759184901284868200286529704023033979366398816632", + 51, + "24212700195838636108555084364432828062549200431806332300085855784", + 31, + 16, + "", + true}, + {decimal_op::divide, + "18241025075303801359299646781236858613374033036953916164633170313584", + 51, + "31532464949533435334935616580309565774604327479014944317663772565", + 38, + 16, + "578484", + false}, + {decimal_op::modulo, + "-83560928004337159522895524415359983899694695261859092459529625366696", + 39, + "-12436612052996461982061751959033724199425939681177274601703138480", + 9, + 16, + "-835609280043371595228955244153599838996946952", + false}, + {decimal_op::add, + "45517215243963473200301786005384648457213564347127805099730635646623", + 41, + "79320958666170412280115184041223940151211217493370165031134513229", + 33, + 9, + "79321413838322851914847187059083993997696", + false}, + {decimal_op::subtract, + "57487653057518032540796079806280277135522531160563731943581546454309", + 59, + "58265383250333113258773120738744693833228050455071776578147805706", + 4, + 4, + "-58265383250333113258773120738744693833228050455071770829382499954", + false}, + {decimal_op::multiply, + "-56590798613317731982372132141549151062594200830634165288028208786704", + 33, + "2368721266459173022227231024880987478075056816235598394414731537", + 47, + 16, + "", + true}, + {decimal_op::divide, + "-73511782348633980095238618012534901670403878071657348715837925192065", + 35, + "-31595571065939402633434114833087490256408708044729856922310810031", + 7, + 1, + "0", + false}, + {decimal_op::modulo, + "47449588332921981644397541675218544538024411417500251526105508090889", + 52, + "-81626323729248385395059507364479197536756748671414628265996370232", + 47, + 11, + "474495883329219816443975416", + false}, + {decimal_op::add, + "49321720532957072007674649903205040551436170333697150454777895403450", + 19, + "79217218933285288286127375781794427873730347749291949263329225754", + 62, + 20, + "", + true}, + {decimal_op::subtract, + "34351874742536557559938925025237157400856817225761550418564433606445", + 32, + "-92767296834264129690938139223448532051070699536387716588738424347", + 15, + 15, + "92767296834264473209685564589024131440320951907961725156910681963", + false}, + {decimal_op::multiply, + "-82218873649850732452002383199172643235078742441848958781711279027542", + 32, + "-61833369839738304888696379719115521722334300002416606516658149028", + 64, + 13, + "50838700221979347232084974288418132562187100916654", + false}, + {decimal_op::divide, + "41033035752691861092169834920174950108209436040127811410353689141547", + 60, + "4056895594565733856430283365983712423618811636521773549774954542", + 11, + 25, + "0", + false}, + {decimal_op::modulo, + "29766948079690554957758475699764350090776118761929488278716729719861", + 41, + "92994954102497624620249448668351723738862840538623968894186122742", + 69, + 19, + "17568247353069", + false}, + {decimal_op::add, + "-6646167618007906481047688537967352098712044874412753814953364779795", + 1, + "-56319694931092815309305644225219620774461772307248408972705964970", + 66, + 13, + "", + true}, + {decimal_op::subtract, + "-76688162943424109969346014126230980072698546765002039661559289482261", + 21, + "70156095919622732096232672851123974250531927614420096985268786423", + 5, + 24, + "", + true}, + {decimal_op::multiply, + "-12313563951373437425514652488408507392655848374793628539025025033118", + 69, + "-83959552346227197780652118761866872317936615422264857403624030181", + 42, + 4, + "10338413171439543325231229", + false}, + {decimal_op::divide, + "77064271301066435037924886313349375305925928505642732446185937992739", + 54, + "30406601727782775428016462362960007250095208194675680199392002945", + 56, + 17, + "25344585360439059751750", + false}, + {decimal_op::modulo, + "47133515847494873413305163120289842835403451552528436970112572430408", + 6, + "55617286653107075897849449067194834890868246508085518234483083027", + 35, + 8, + "1720736603735804625518602673766944581", + false}, + {decimal_op::add, + "-27553014280157562546075144307131744987008950172934356410140997918856", + 14, + "10198890656897175427168996621209070618826301827904373709539354580", + 5, + 29, + "", + true}, + {decimal_op::subtract, + "49643732802384828577082811938853296174023616314108332442158131477794", + 31, + "16216277604330379540008770975330709279340310473316035321568100893", + 36, + 28, + "49643732640222052533779016538765586420716523520705227708997778262", + false}, + {decimal_op::multiply, + "9242598567971653689309514688657133749216693352133238408460035745888", + 11, + "-92764005818737212231115528347209401055426360046426503996816309494", + 38, + 16, + "", + true}, + {decimal_op::divide, + "79402539588926734555621589200183634923498483671433602726149291090474", + 49, + "-55611139563582772909343330447378850791409804312938033345083607355", + 59, + 9, + "-14278171641878001807355", + false}, + {decimal_op::modulo, + "4467359765444642000246478872806763023149342653222471686797536395235", + 53, + "-55140143426754081313115580484516970830740413636827083270541525633", + 9, + 13, + "446735976544464200024647887", + false}, + {decimal_op::add, + "-62532871468114014368965251934357933873041719379698296053895670367620", + 8, + "-57977065487320986289861316474104018050647343669031923332264540770", + 26, + 9, + "", + true}, + {decimal_op::subtract, + "-89868445128071931797907101224423825256735670963671554256753833206603", + 30, + "84233921243271531439947067364049507124173028507827402583808326224", + 13, + 14, + "", + true}, + {decimal_op::multiply, + "4583752179308648586799784846042204789592020746973383275466462152161", + 0, + "52542953259315895665592332960721947200061122656513381923963122768", + 51, + 16, + "", + true}, + {decimal_op::divide, + "21754392811944196540484950901756720376220272599837720276164334160769", + 20, + "-86188760572739037015731626582092557913251094863425131316481138527", + 10, + 0, + "0", + false}, + {decimal_op::modulo, + "43638326046124535096926703332841072158583521743428030117187298792917", + 10, + "35097747944787662895859400985289728803650792043333951534669136391", + 4, + 9, + "", + true}, + {decimal_op::add, + "48781448168678828792476406505295227714069453823719230568733913737113", + 33, + "67161479142937478968140284431664728394714155837667591597489128511", + 63, + 17, + "4878144816867882879247640650529529487554859676119820", + false}, + {decimal_op::subtract, + "-6560386506164668004737780034282060761861440577904747964249430340843", + 4, + "-17719674935886243353080286315246890686351792909439534013164570249", + 26, + 29, + "", + true}, + {decimal_op::multiply, + "44677772555189417860452086999137562636035232152763423507199519452207", + 64, + "50004011810945784093091447642540091422030203529106278434593408060", + 15, + 19, + "", + true}, + {decimal_op::divide, + "-43538825580732423711693172142903991300520710406748132924322703076813", + 35, + "-50107779658446339494281380341839924927464647579523773961870254419", + 30, + 25, + "86890350914587710504699", + false}, + {decimal_op::modulo, + "-57190523796342753444249716381387935284193192288695592276240162525148", + 57, + "-78486415589302233594884492762501316785818526346874997432837592062", + 25, + 23, + "-5719052379634275344424971638138793", + false}, + {decimal_op::add, + "-80834409030848307436393210019819153960196379159461041331989082451588", + 42, + "63949870220517844500938025886148393104950727611857089228480391347", + 55, + 14, + "-8083440903084830104140618796803470386639", + false}, + {decimal_op::subtract, + "-40938093545632709075864830134685746175551648149822986570498824713644", + 15, + "-75128935223636747817814531086521650518149556771295117197205239039", + 66, + 13, + "", + true}, + {decimal_op::multiply, + "76399162111648096761324458233371366488578961581508871977109779443789", + 69, + "40537815269532036442645929627567027313645451330913038272190989992", + 50, + 16, + "309705512042902163807701901412", + false}, + {decimal_op::divide, + "-67652072890002630757476596629621675576086227278874315857634785574820", + 62, + "2234226560345107644033773365169304945212223328297922455201081027", + 52, + 22, + "-30279862432372489", + false}, + {decimal_op::modulo, + "86039386734972433594798493627353169421408220902982662203585748263938", + 3, + "33289369160061679013643649290705173758508381506067401712086933882", + 45, + 29, + "2892374593361874656974966417141499850896667439108", + false}, + {decimal_op::add, + "-95691307932891187750780067352718531937708850061351970194887843139510", + 31, + "36465351443536655860796990798929588346718623853617659419850047328", + 47, + 22, + "-95691307932891187747133532208364866351629150981459011360216", + false}, + {decimal_op::subtract, + "53109111594343446362564029593228905004462610662065864406254172869291", + 8, + "-16161447588363623694234246132221164318151851454832781739995354218", + 0, + 29, + "", + true}, + {decimal_op::multiply, + "-29373322844011887745758601923747493825595883611663794360430885836052", + 14, + "667277503926172257715430273358634674838030964568080254964747540", + 21, + 29, + "", + true}, + {decimal_op::divide, + "-85508083527099839239263035166948993620206264140824438545192672096128", + 54, + "-56618430043507036291919528971184669283308826761991174083549585962", + 22, + 1, + "0", + false}, + {decimal_op::modulo, + "-88946757424837344243657319544248322360749316484074328720471343053196", + 26, + "-71355746322721020171155970917635748469919043842551213984287476843", + 12, + 25, + "", + true}, + {decimal_op::add, + "3201431874656376589168394707261202188525971549757058271744817314702", + 32, + "-49653782316895305336375863927389839682804832602464034632563211596", + 10, + 30, + "", + true}, + {decimal_op::subtract, + "-22544886520506045163032628459119142475641908011243673404511820977411", + 63, + "20602562485593010600679518720007463374244698497298859311084787288", + 5, + 14, + "", + true}, + {decimal_op::multiply, + "-85351670091996523467496536701893868091552573553618269007380246452499", + 65, + "48469002426700998759205220604429913207936569570879716396090367881", + 45, + 29, + "-4136910304811962553960238469659103452792705033257317", + false}, + {decimal_op::divide, + "90955162854384512890554059365307616635522956627525886907649648204843", + 53, + "-55989928166705099371167935576676265224642939600676989574709228905", + 56, + 1, + "-16244915", + false}, + {decimal_op::modulo, + "94829335522448787709297090606459738874085864317190389685287049183913", + 23, + "54834281820832589312763344638667960536312010467599073979703166405", + 53, + 15, + "285203838500889899614743117", + false}, + {decimal_op::add, + "29671945580959755056135721606272448238402604674491588839687081369381", + 41, + "-9764411330112823042690235398109225273938722560816965125577075526", + 16, + 30, + "", + true}, + {decimal_op::subtract, + "57654565860977595223158323575708546981534271798470716833731472568936", + 41, + "-88298782072246210910451884374373471966938321214939541855359194278", + 45, + 29, + "57654574690855802447779414620896984418881468492302838328", + false}, + {decimal_op::multiply, + "13843945728122959590400061545725091094221822573651536299619825612227", + 69, + "14793377422572250198271925680858116759208013594869460103139242757", + 58, + 8, + "2047987141737", + false}, + {decimal_op::divide, + "-18816268399970172886553852366605808052280963673099778352127787513034", + 24, + "-38424403187607020429085414123373986877780446472448362899612165439", + 65, + 27, + "", + true}, + {decimal_op::modulo, + "37527054384815833754940132503670872167220525248775345506372558302722", + 65, + "-93200998685286181640742146635292431056998067832549694326209693628", + 43, + 7, + "3752705438", + false}, + {decimal_op::add, + "-45902675558909750521081495575406074199017090519844716425337054104843", + 8, + "-61818068607829728158563559810436909571668703559606335150338440399", + 41, + 23, + "", + true}, + {decimal_op::subtract, + "45447403516764477220497742524099798780854129943206986224754458716225", + 4, + "-81428544310898493137041864935961256388383793822854638306341678470", + 42, + 6, + "", + true}, + {decimal_op::multiply, + "30104695795624922600328777794177915989961945278394783003638856212235", + 35, + "-31984180757593199183319916186553390553387218621651260528021406500", + 7, + 13, + "", + true}, + {decimal_op::divide, + "-56611066261256018216247316010106281566251802150124646590832192824607", + 46, + "-75700975960476369249169538914102228448915330811038896344979846357", + 57, + 24, + "74782478750092692310535344203151626985", + false}, + {decimal_op::modulo, + "77241280487362246010398828532927486124919843242046705132791520122649", + 56, + "-70002141497770118678666718232282517565010607262820258860283817906", + 39, + 2, + "77241280487362", + false}, + {decimal_op::add, + "91685001081989133734539387604051855170689021635877529048128718945726", + 6, + "-91527678171078254678706948522067056550292039427295200931243967812", + 46, + 20, + "", + true}, + {decimal_op::subtract, + "-32485598978725423927065335311464603529763046045442314528048545031612", + 60, + "-14600367052685151800666967741791826449767497227527422355080929419", + 54, + 20, + "1456788145370642637673990240648", + false}, + {decimal_op::multiply, + "-66535584858275731124313220582321998310682035607273936746757770844670", + 16, + "-93519592238177723369910564801609344584251943082625173391614076850", + 54, + 24, + "", + true}, + {decimal_op::divide, + "-93248998929224503895486708797074204975844787192449872246606281945304", + 6, + "94619130281437803591390985143149551867975689371804746860624558057", + 30, + 19, + "-9855195101863868073207450885420092862911802041", + false}, + {decimal_op::modulo, + "-23508830058089209777688028499382590518400019018305811491715975571247", + 2, + "-72301543272633525975725649494809673899103703878785074871256155909", + 21, + 13, + "-676756512104023782873133132794930730110143423008178244677", + false}, + {decimal_op::add, + "21695269739260146911103096741960583000095297416675957291118859863485", + 44, + "59336321680604279579995705866433333898669537357388931850494503259", + 30, + 9, + "59336321680821232277388307335544364866089143", + false}, + {decimal_op::subtract, + "68550912592089252102295348446217274063372108190264796909694757475482", + 30, + "-33552523974612157511083324391839586091659743949573444993545610877", + 61, + 14, + "6855091259208925210229534844621727741862450565148055", + false}, + {decimal_op::multiply, + "-45401877812543875327519031917628748661221015370546936409958575135240", + 31, + "35974848592716408986324436107069216483576049737946224971911259955", + 25, + 5, + "", + true}, + {decimal_op::divide, + "-31702934359710622794455420459297055034151263470199963485471171936741", + 26, + "-94599633276565895648136479038750963878160906882156573972393589280", + 28, + 25, + "335127455167038523287957237059", + false}, + {decimal_op::modulo, + "-76317216301960949082501339569419379367110816773170516601985326983546", + 30, + "32499501680586744552813138176334408539202049471332639287345491366", + 16, + 8, + "-7631721630196094908250133956941937936711081677", + false}, + {decimal_op::add, + "97326373482047965509120234622202228182580945736045793147938780478857", + 58, + "-32015417709426732154819937570296661505111621151501214757130169692", + 59, + 20, + "973231719402770228359047526284", + false}, + {decimal_op::subtract, + "25089679972506708545676831501255473372441199354436860283788955660909", + 59, + "-32279545762468728971169744388257977377608737953938564614558797227", + 57, + 9, + "283176345487535814", + false}, + {decimal_op::multiply, + "-56327812931026929462885997786464193282169413576863751347739245824904", + 47, + "71616107834204353281281241469775414834731037412291453131051086032", + 38, + 26, + "", + true}, + {decimal_op::divide, + "-82207807432672694303015622295724905291969934678273729487725218497990", + 27, + "63299234371615072141850666654947681070827392543707732524380366597", + 12, + 13, + "-13", + false}, + {decimal_op::modulo, + "57777136675011706665787252982491847550025107315963279086939347555183", + 61, + "-94522045618341660772124896559945140710183111428151016095037440834", + 18, + 2, + "577771366", + false}, + {decimal_op::add, + "87749741336498637671355273480251458815828512085725075992945471961504", + 16, + "94055293893281701898247862647550480752769605258639162675723731122", + 50, + 15, + "", + true}, + {decimal_op::subtract, + "-65431424333147490037509632273200203463776346631550119331618081585460", + 0, + "41130341761653272597034803201840340425173532365393746215659583777", + 33, + 21, + "", + true}, + {decimal_op::multiply, + "61943808792623274099160258958191877988905608291591001276427913761844", + 25, + "30093210910203633992193529841316090827702300437125832325928250541", + 2, + 12, + "", + true}, + {decimal_op::divide, + "-93203608596794982347575921238101222097399800198728085536586801650787", + 8, + "-14702574708971255601293087803703493679627503060445067817602516411", + 0, + 11, + "6339271", + false}, + {decimal_op::modulo, + "8328464685863787954281851516535787563067727468857232472972311222044", + 55, + "-12905787694824881758493769143618174332702554645151882577274426096", + 9, + 29, + "83284646858637879542818515165357875630677", + false}, + {decimal_op::add, + "-5378596046242743298358019890348219773914647562530499629946199836342", + 53, + "-26587200000644441883427602299966187118927150380652588923459356320", + 22, + 4, + "-26587200000644441883427602300504046723551424710", + false}, + {decimal_op::subtract, + "-39043750391651823468181717985011529778205632876042340727779866887564", + 24, + "9763836384081754881438300750940765675621232967826039599421815439", + 41, + 14, + "-3904375039165182346827935634885234732702001588355174838454", + false}, + {decimal_op::multiply, + "-57942317984997089054444561367039214180803752724482737805934081000924", + 39, + "-56160588090317144676961521868164513515363862806940190511207102226", + 48, + 9, + "3254074653353596415955462275178848268722352641742388912", + false}, + {decimal_op::divide, + "-3149460375497038165505632863797591832061945529233231324369887456917", + 41, + "15806432848631760288234558598178429107549147611017432368489344765", + 42, + 2, + "-199252", + false}, + {decimal_op::modulo, + "12584951688411334248524493169238444899865284719332849563640581203873", + 15, + "34866428325844646458507259737817050857451467497360734085546263837", + 58, + 15, + "1707364813568782039446", + false}, + {decimal_op::add, + "47815834734894925924171607242036396249758236215413448610936728733880", + 50, + "-9159781376552138843403314616416336599226911001388631613117685506", + 60, + 22, + "4781583473489400994603395202815205591830", + false}, + {decimal_op::subtract, + "-91180102197969433614185045890629960805003806905814869809878403441952", + 65, + "92790345173570783320383717092711836189064853875606257444803436378", + 11, + 25, + "", + true}, + {decimal_op::multiply, + "-20892175077832515040274173931336558602871729160716211580608627084368", + 24, + "-59670916225989460565336396491384812805256075781211053642259904560", + 26, + 30, + "", + true}, + {decimal_op::divide, + "84150380736383242613740178394054037593049354450112564202423480018541", + 48, + "-30799397332983254652478358510843950143821159655627925430222719582", + 11, + 21, + "0", + false}, + {decimal_op::modulo, + "-90800147423919206486487571192349973037832533516190924784329302451538", + 64, + "27355500375419605662287703019429954052411053091570626888030958739", + 41, + 17, + "-908001474239192064864", + false}, + {decimal_op::add, + "51004494118968476605381921673756358327458483425549941007995429422430", + 46, + "-22602226491494325593029500867237804937847369030150426289177107064", + 60, + 29, + "510044941189684763793596567588131023971634747531719", + false}, + {decimal_op::subtract, + "-77116105843252114087905047622457197709974437253849684618873015072922", + 55, + "78115334356130919264390054998282314699645111856370044932725996149", + 67, + 4, + "-77116105843252192", + false}, + {decimal_op::multiply, + "39046526482868385270201290403762423719559746523962443741018282763175", + 53, + "-53380067106443357269987925204191229161391868875068248695986233502", + 59, + 23, + "-20843062039290321264731625741090142707627919", + false}, + {decimal_op::divide, + "-49736178715357074322666105677424578596691514442675328225075828425395", + 61, + "-54546013932461772015515216757571526519206048618438162862178125060", + 3, + 14, + "0", + false}, + {decimal_op::modulo, + "-98671719838755245581902603601295466253981291125323236318890420060556", + 6, + "68751117389653248691808481281802631661403138191847297410384736313", + 41, + 3, + "-670344989637517760786971088", + false}, + {decimal_op::add, + "81781201245570203618395086061072833681512346314201813921786580822665", + 37, + "6985232314992467960536624235678561780905700323743844072529508595", + 18, + 13, + "69852323149924761386567487926989236204143064310272122237641", + false}, +}; diff --git a/test/cpp/embedding/test_exact_decimal.cpp b/test/cpp/embedding/test_exact_decimal.cpp new file mode 100644 index 0000000000..ba64a7609e --- /dev/null +++ b/test/cpp/embedding/test_exact_decimal.cpp @@ -0,0 +1,219 @@ +/* Copyright 2026 Sirius Contributors. SPDX-License-Identifier: Apache-2.0 */ +#include "numeric/exact_decimal.hpp" + +#include + +#include +#include +#include + +using namespace sirius::mo_decimal; +#include "exact_decimal_wide_cases.hpp" +namespace { +// Test-only independent decimal-string fixture builder, not a production codec +// or an arithmetic oracle that calls the implementation under test. +coefficient number(std::string digits) +{ + coefficient result; + bool negative = digits[0] == '-'; + for (auto digit : digits.substr(negative ? 1 : 0)) { + uint64_t carry = digit - '0'; + for (auto& word : result.words) { + auto next = uint64_t(word) * 10 + carry; + word = next; + carry = next >> 32; + } + REQUIRE(carry == 0); + } + if (negative) { + uint64_t carry = 1; + for (auto& word : result.words) { + auto next = uint64_t(~word) + carry; + word = next; + carry = next >> 32; + } + } + return result; +} +void expect(decimal_result const& result, std::string const& value) +{ + REQUIRE(result.error == decimal_error::none); + REQUIRE(result.valid); + auto expected = number(value); + for (int i = 0; i < 8; ++i) + CHECK(result.value.words[i] == expected.words[i]); +} +decimal_result run(decimal_op op, + std::string const& a, + decimal_type at, + std::string const& b, + decimal_type bt, + decimal_type out) +{ + return evaluate_decimal(op, number(a), at, true, number(b), bt, true, out); +} +} // namespace + +TEST_CASE("MO exact descriptor and coefficient bytes preserve physical width", "[exact_decimal]") +{ + CHECK((decimal_type{256, 15, 2}.valid())); + CHECK_FALSE((decimal_type{128, 39, 0}.valid())); + CHECK_FALSE((decimal_type{256, 77, 0}.valid())); + CHECK_FALSE((decimal_type{64, 15, 16}.valid())); + CHECK_FALSE((decimal_type{32, 9, 0}.valid())); + for (auto width : {8u, 16u, 32u}) { + for (auto n : {"0", "1", "-1", "9007199254740993", "-123456789"}) { + auto original = number(n); + std::array bytes{}; + store_coefficient(original, bytes.data(), width); + auto observed = load_coefficient(bytes.data(), width); + for (int i = 0; i < 8; ++i) + CHECK(observed.words[i] == original.words[i]); + } + } +} + +TEST_CASE("MO exact scalars widen and round at the declared boundary", "[exact_decimal]") +{ + auto max38 = std::string(38, '9'); + expect(run(decimal_op::add, max38, {128, 38, 0}, "1", {64, 1, 0}, {256, 39, 0}), + "1" + std::string(38, '0')); + expect(run(decimal_op::subtract, "-12", {64, 2, 0}, "5", {64, 1, 0}, {64, 3, 0}), "-17"); + expect(run(decimal_op::add, "100", {64, 3, 2}, "-2", {64, 1, 0}, {64, 3, 2}), "-100"); + for (auto sign : {std::string(), std::string("-")}) { + expect(run(decimal_op::cast, sign + "125", {64, 3, 2}, "0", {64, 1, 0}, {64, 2, 1}), + sign + "13"); + expect(run(decimal_op::divide, sign + "125", {64, 3, 2}, "1", {64, 1, 0}, {128, 38, 1}), + sign + "13"); + expect(run(decimal_op::integer_divide, sign + "125", {64, 3, 2}, "1", {64, 1, 0}, {64, 18, 0}), + sign + "1"); + expect(run(decimal_op::modulo, sign + "125", {64, 3, 2}, "10", {64, 2, 1}, {64, 3, 2}), + sign + "25"); + } + // Raw product 10^80 is wider than 256 bits; its once-scaled 10^50 result fits. + expect(run(decimal_op::multiply, + "1" + std::string(40, '0'), + {256, 41, 30}, + "1" + std::string(40, '0'), + {256, 41, 30}, + {256, 65, 30}), + "1" + std::string(50, '0')); + for (auto scale : {2u, 6u, 12u, 30u}) + expect(run(decimal_op::divide, + "100", + {64, 10, 2}, + "300", + {64, 10, 2}, + {256, 65, static_cast(scale)}), + std::string(scale, '3')); + expect(run(decimal_op::divide, + "1" + std::string(36, '0'), + {128, 38, 37}, + "2", + {64, 1, 0}, + {128, 38, 30}), + "5" + std::string(28, '0')); +} + +TEST_CASE("MO exact errors NULLs masks and comparison remain distinct", "[exact_decimal]") +{ + auto maximum = std::string(65, '9'); + expect(run(decimal_op::add, maximum, {256, 65, 0}, "0", {64, 1, 0}, {256, 65, 0}), maximum); + CHECK(run(decimal_op::add, maximum, {256, 65, 0}, "1", {64, 1, 0}, {256, 65, 0}).error == + decimal_error::out_of_range); + CHECK(run(decimal_op::cast, "100", {64, 3, 0}, "0", {64, 1, 0}, {64, 2, 0}).error == + decimal_error::invalid_input); + auto zero = number("0"), one = number("1"); + CHECK(evaluate_decimal( + static_cast(255), one, {64, 1, 0}, true, one, {64, 1, 0}, true, {64, 1, 0}) + .error == decimal_error::invalid_input); + for (auto op : {decimal_op::divide, decimal_op::integer_divide, decimal_op::modulo}) { + auto r = evaluate_decimal(op, one, {64, 1, 0}, true, zero, {64, 1, 0}, true, {64, 18, 0}); + CHECK_FALSE(r.valid); + CHECK(r.error == decimal_error::none); + } + auto r = evaluate_decimal(decimal_op::add, + number(maximum), + {256, 65, 0}, + true, + one, + {64, 1, 0}, + true, + {256, 65, 0}, + false); + CHECK_FALSE(r.valid); + CHECK(r.error == decimal_error::none); + r = evaluate_decimal(decimal_op::multiply, + number(maximum), + {256, 65, 0}, + false, + number(maximum), + {256, 65, 0}, + true, + {256, 65, 0}); + CHECK_FALSE(r.valid); + CHECK(r.error == decimal_error::none); + expect(run(decimal_op::equal, "100", {64, 3, 2}, "10000", {128, 5, 4}, {64, 1, 0}), "1"); + expect(run(decimal_op::less, + "-9007199254740993", + {64, 18, 0}, + "-9007199254740992", + {64, 18, 0}, + {64, 1, 0}), + "1"); + expect(run(decimal_op::equal, + "9007199254740993", + {64, 18, 0}, + "9007199254740992", + {64, 18, 0}, + {64, 1, 0}), + "0"); + expect(run(decimal_op::integer_divide, + "-9223372036854775808", + {128, 38, 0}, + "1", + {64, 1, 0}, + {64, 18, 0}), + "-9223372036854775808"); + CHECK( + run( + decimal_op::integer_divide, "9223372036854775808", {128, 38, 0}, "1", {64, 1, 0}, {64, 18, 0}) + .error == decimal_error::out_of_range); +} + +TEST_CASE("MO exact small-domain operations match independent integer arithmetic", + "[exact_decimal]") +{ + std::mt19937 random(28968); + for (int i = 0; i < 256; ++i) { + int64_t a = int64_t(random() % 2000001) - 1000000; + int64_t b = int64_t(random() % 2000000) - 1000000; + if (!b) b = 1; + auto av = std::to_string(a), bv = std::to_string(b); + expect(run(decimal_op::add, av, {64, 18, 0}, bv, {64, 18, 0}, {128, 38, 0}), + std::to_string(a + b)); + expect(run(decimal_op::subtract, av, {64, 18, 0}, bv, {64, 18, 0}, {128, 38, 0}), + std::to_string(a - b)); + expect(run(decimal_op::multiply, av, {64, 18, 0}, bv, {64, 18, 0}, {128, 38, 0}), + std::to_string(a * b)); + expect(run(decimal_op::integer_divide, av, {64, 18, 0}, bv, {64, 18, 0}, {64, 18, 0}), + std::to_string(a / b)); + expect(run(decimal_op::modulo, av, {64, 18, 0}, bv, {64, 18, 0}, {128, 38, 0}), + std::to_string(a % b)); + } +} + +TEST_CASE("MO exact wide operations match an independent integer rational oracle", + "[exact_decimal]") +{ + for (auto const& test : wide_cases) { + INFO("op=" << static_cast(test.op) << " a=" << test.a << " b=" << test.b); + auto result = + run(test.op, test.a, {256, 76, test.sa}, test.b, {256, 76, test.sb}, {256, 65, test.so}); + if (test.overflow) { + CHECK(result.error == decimal_error::out_of_range); + CHECK_FALSE(result.valid); + } else + expect(result, test.expected); + } +} diff --git a/test/cpp/embedding/test_exact_decimal_gpu.cpp b/test/cpp/embedding/test_exact_decimal_gpu.cpp new file mode 100644 index 0000000000..3ec4014d72 --- /dev/null +++ b/test/cpp/embedding/test_exact_decimal_gpu.cpp @@ -0,0 +1,442 @@ +/* Copyright 2026 Sirius Contributors. SPDX-License-Identifier: Apache-2.0 */ +#include "embedding/result_codec.hpp" +#include "numeric/exact_decimal_gpu.hpp" + +#include + +#include + +#include + +#include + +#include +#include +#include + +using namespace sirius::mo_decimal; +#include "exact_decimal_wide_cases.hpp" +namespace { +coefficient decimal_fixture(char const* text) +{ + coefficient result; + bool negative = *text == '-'; + if (negative) ++text; + for (; *text; ++text) { + uint64_t carry = *text - '0'; + for (auto& word : result.words) { + auto value = uint64_t(word) * 10 + carry; + word = value; + carry = value >> 32; + } + REQUIRE(carry == 0); + } + if (negative) { + uint64_t carry = 1; + for (auto& word : result.words) { + auto value = uint64_t(~word) + carry; + word = value; + carry = value >> 32; + } + } + return result; +} +coefficient small(int64_t value) +{ + coefficient result; + auto bits = static_cast(value); + result.words[0] = bits; + result.words[1] = bits >> 32; + for (int i = 2; i < 8; ++i) + result.words[i] = value < 0 ? UINT32_MAX : 0; + return result; +} +std::unique_ptr input_column(decimal_type type, + std::vector const& values, + uint32_t validity, + rmm::cuda_stream_view stream) +{ + auto column = make_decimal_column(type, + values.size(), + cudf::mask_state::ALL_VALID, + stream, + cudf::get_current_device_resource_ref()); + auto view = column->mutable_view(); + if (type.bits != 256) { + std::vector raw(values.size() * type.bytes()); + for (std::size_t i = 0; i < values.size(); ++i) + store_coefficient(values[i], raw.data() + i * type.bytes(), type.bytes()); + REQUIRE( + cudaMemcpyAsync( + view.data(), raw.data(), raw.size(), cudaMemcpyHostToDevice, stream.value()) == + cudaSuccess); + stream.synchronize(); + } else { + std::vector lane(values.size()); + for (int i = 0; i < 4; ++i) { + for (std::size_t row = 0; row < values.size(); ++row) + lane[row] = uint64_t(values[row].words[(3 - i) * 2]) | + uint64_t(values[row].words[(3 - i) * 2 + 1]) << 32; + REQUIRE(cudaMemcpyAsync(view.child(i).data(), + lane.data(), + lane.size() * 8, + cudaMemcpyHostToDevice, + stream.value()) == cudaSuccess); + stream.synchronize(); + } + } + REQUIRE(cudaMemcpyAsync(view.null_mask(), &validity, 4, cudaMemcpyHostToDevice, stream.value()) == + cudaSuccess); + stream.synchronize(); + column->set_null_count(cudf::null_count(view.null_mask(), 0, values.size(), stream)); + return column; +} +std::vector download(cudf::column_view column, + uint32_t bytes, + rmm::cuda_stream_view stream) +{ + std::vector result(column.size()); + if (column.type().id() == cudf::type_id::STRUCT) { + std::vector lane(column.size()); + for (int i = 0; i < 4; ++i) { + REQUIRE(cudaMemcpyAsync(lane.data(), + column.child(i).head() + column.offset(), + lane.size() * 8, + cudaMemcpyDeviceToHost, + stream.value()) == cudaSuccess); + stream.synchronize(); + for (std::size_t row = 0; row < result.size(); ++row) { + result[row].words[(3 - i) * 2] = lane[row]; + result[row].words[(3 - i) * 2 + 1] = lane[row] >> 32; + } + } + } else { + std::vector raw(result.size() * bytes); + REQUIRE( + cudaMemcpyAsync( + raw.data(), column.head(), raw.size(), cudaMemcpyDeviceToHost, stream.value()) == + cudaSuccess); + stream.synchronize(); + for (std::size_t row = 0; row < result.size(); ++row) { + if (bytes == 1) + result[row].words[0] = raw[row]; + else + result[row] = load_coefficient(raw.data() + row * bytes, bytes); + } + } + return result; +} +struct storage : sirius::embedding::input_storage { + std::vector data; + explicit storage(std::size_t size) : data(size) {} + std::size_t size() const override { return data.size(); } + void visit(std::function)> const& visitor) override + { + for (std::size_t at = 0; at < data.size(); at += 13) + visitor(at, std::span(data).subspan(at, std::min(13, data.size() - at))); + } +}; +struct pool : sirius::embedding::input_pool { + std::size_t rounded(std::size_t size) const override { return size; } + std::unique_ptr allocate(std::size_t size) override + { + return std::make_unique(size); + } +}; +} // namespace + +TEST_CASE("MO checked decimal GPU kernels preserve masks and all coefficient widths", + "[exact_decimal_gpu]") +{ + rmm::cuda_stream stream; + std::vector left{small(125), small(-125), small(9007199254740993LL), small(0)}; + std::vector right{small(100), small(100), small(1), small(0)}; + bool active[]{true, true, false, true}; + rmm::device_buffer mask(active, sizeof(active), stream.view()); + stream.synchronize(); + for (auto width : {64u, 128u, 256u}) { + decimal_type type{static_cast(width), + static_cast(width == 64 ? 18 + : width == 128 ? 38 + : 65), + 2}; + auto a = input_column(type, left, 0b0111, stream.view()); + auto b = input_column(type, right, 0b1111, stream.view()); + for (auto op : {decimal_op::add, + decimal_op::subtract, + decimal_op::multiply, + decimal_op::divide, + decimal_op::integer_divide, + decimal_op::modulo, + decimal_op::negate, + decimal_op::cast, + decimal_op::equal, + decimal_op::not_equal, + decimal_op::less, + decimal_op::less_equal, + decimal_op::greater, + decimal_op::greater_equal}) { + INFO("width=" << width << " op=" << static_cast(op)); + auto result = evaluate_decimal_columns(op, + a->view(), + type, + b->view(), + type, + type, + static_cast(mask.data()), + stream.view(), + cudf::get_current_device_resource_ref()); + auto bytes = op >= decimal_op::equal ? 1 + : op == decimal_op::integer_divide ? 8 + : type.bytes(); + auto values = download(result.values->view(), bytes, stream.view()); + std::array errors{}; + uint32_t validity{}; + REQUIRE(cudaMemcpyAsync( + errors.data(), result.errors.data(), 4, cudaMemcpyDeviceToHost, stream.value()) == + cudaSuccess); + REQUIRE(cudaMemcpyAsync(&validity, + result.values->view().null_mask(), + 4, + cudaMemcpyDeviceToHost, + stream.value()) == cudaSuccess); + stream.synchronize(); + for (uint32_t row = 0; row < 4; ++row) { + auto expected = evaluate_decimal( + op, left[row], type, row != 3, right[row], type, true, type, active[row]); + CHECK(errors[row] == static_cast(expected.error)); + CHECK(bool((validity >> row) & 1) == expected.valid); + if (expected.valid) + for (int limb = 0; limb < 8; ++limb) + CHECK(values[row].words[limb] == expected.value.words[limb]); + } + CHECK(result.values->null_count() == 2); + } + // The operation owns the failure class even when the two operations + // overflow the same declared coefficient domain. Inactive and NULL rows + // must retain neither a value nor an arithmetic failure. + for (auto op : {decimal_op::add, decimal_op::cast}) { + auto result = evaluate_decimal_columns(op, + a->view(), + type, + a->view(), + type, + {type.bits, 2, 2}, + static_cast(mask.data()), + stream.view(), + cudf::get_current_device_resource_ref()); + std::array errors{}; + REQUIRE(cudaMemcpyAsync(errors.data(), + result.errors.data(), + errors.size(), + cudaMemcpyDeviceToHost, + stream.value()) == cudaSuccess); + stream.synchronize(); + auto expected = + op == decimal_op::cast ? decimal_error::invalid_input : decimal_error::out_of_range; + CHECK(errors[0] == static_cast(expected)); + CHECK(errors[1] == static_cast(expected)); + CHECK(errors[2] == 0); + CHECK(errors[3] == 0); + CHECK(result.values->null_count() == 4); + } + // A sliced parent carries both a data offset and a validity offset; + // Decimal256 children remain in their original physical allocation. + auto sliced_a = cudf::slice(a->view(), {2, 4}); + auto sliced_b = cudf::slice(b->view(), {2, 4}); + auto sliced = evaluate_decimal_columns(decimal_op::cast, + sliced_a.front(), + type, + sliced_b.front(), + type, + type, + nullptr, + stream.view(), + cudf::get_current_device_resource_ref()); + auto observed = download(sliced.values->view(), type.bytes(), stream.view()); + CHECK(sliced.values->null_count() == 1); + for (int word = 0; word < 8; ++word) + CHECK(observed[0].words[word] == left[2].words[word]); + } +} + +TEST_CASE("MO Decimal64 validity spans complete warps and a partial tail", "[exact_decimal_gpu]") +{ + rmm::cuda_stream stream; + decimal_type type{64, 18, 0}; + std::vector ones(65, small(1)); + auto a = input_column(type, ones, UINT32_MAX, stream.view()); + auto b = input_column(type, ones, UINT32_MAX, stream.view()); + uint32_t validity[]{UINT32_MAX ^ (uint32_t(1) << 31), UINT32_MAX ^ 1u, 1}; + REQUIRE(cudaMemcpyAsync(a->mutable_view().null_mask(), + validity, + sizeof(validity), + cudaMemcpyHostToDevice, + stream.value()) == cudaSuccess); + std::array active; + active.fill(true); + active[63] = false; + rmm::device_buffer mask(active.data(), active.size() * sizeof(bool), stream.view()); + stream.synchronize(); + a->set_null_count(2); + auto result = evaluate_decimal_columns(decimal_op::add, + a->view(), + type, + b->view(), + type, + type, + static_cast(mask.data()), + stream.view(), + cudf::get_current_device_resource_ref()); + auto values = download(result.values->view(), 8, stream.view()); + std::array errors{}; + REQUIRE( + cudaMemcpyAsync( + errors.data(), result.errors.data(), errors.size(), cudaMemcpyDeviceToHost, stream.value()) == + cudaSuccess); + REQUIRE(cudaMemcpyAsync(validity, + result.values->view().null_mask(), + sizeof(validity), + cudaMemcpyDeviceToHost, + stream.value()) == cudaSuccess); + stream.synchronize(); + CHECK(result.values->null_count() == 3); + for (uint32_t row = 0; row < ones.size(); ++row) { + auto valid = row != 31 && row != 32 && row != 63; + CHECK(bool((validity[row / 32] >> (row % 32)) & 1) == valid); + CHECK(errors[row] == 0); + CHECK(values[row].words[0] == (valid ? 2 : 0)); + } +} + +TEST_CASE("MO Decimal256 result codec honors slices and bounded scratch", "[exact_decimal_gpu]") +{ + using namespace sirius::embedding; + rmm::cuda_stream stream; + std::vector values{small(17), small(-17), small(9007199254740993LL), small(99)}; + auto column = input_column({256, 15, 2}, values, 0b1011, stream.view()); + auto view = column->view(); + cudf::column_view sliced(view.type(), + 3, + nullptr, + view.null_mask(), + 1, + 1, + {view.child(0), view.child(1), view.child(2), view.child(3)}); + CHECK(decimal_column_matches(sliced, {256, 15, 2})); + std::vector schema{{34, 15, 2, true, "wide"}}; + auto result = std::make_shared(); + result->activate(std::make_shared()); + std::shared_ptr scratch; + // Two rows fit the scratch slab: the three-row slice crosses its chunk + // boundary. The mock storage also splits rows and limbs into 13-byte blocks. + REQUIRE(result->try_allocate(64, 0, scratch) == SIRIUS_OK); + for (uint32_t begin : {0, 1}) { + auto layout = + size_native_result(cudf::table_view({sliced}), begin, schema, *scratch, stream.view(), 4096); + REQUIRE(layout.rows == 3 - begin); + std::shared_ptr batch; + REQUIRE(result->try_allocate(layout.bytes, 1, batch) == SIRIUS_OK); + encode_native_result( + cudf::table_view({sliced}), begin, schema, layout, *batch, *scratch, stream.view()); + std::array encoded{}; + batch->storage->read(batch->columns[0].data_offset, + {reinterpret_cast(encoded.data()), layout.rows * 32}); + for (uint32_t row = 0; row < layout.rows; ++row) { + auto observed = load_coefficient(encoded.data() + row * 32, 32); + auto expected = row + begin == 1 ? coefficient{} : values[row + begin + 1]; + for (int word = 0; word < 8; ++word) + CHECK(observed.words[word] == expected.words[word]); + } + uint64_t nulls{}; + batch->storage->read(batch->columns[0].null_offset, + {reinterpret_cast(&nulls), sizeof(nulls)}); + CHECK(nulls == (begin == 0 ? 2 : 1)); + } + result->cancel(); +} + +TEST_CASE("MO full-width GPU arithmetic matches independent wide rational vectors", + "[exact_decimal_gpu]") +{ + rmm::cuda_stream stream; + bool enabled[]{true, false}; + rmm::device_buffer active(enabled, sizeof(enabled), stream.view()); + stream.synchronize(); + for (auto const& test : wide_cases) { + INFO("op=" << static_cast(test.op) << " a=" << test.a << " b=" << test.b); + auto av = decimal_fixture(test.a), bv = decimal_fixture(test.b); + auto a = input_column({256, 76, test.sa}, {av, av}, 3, stream.view()); + auto b = input_column({256, 76, test.sb}, {bv, bv}, 3, stream.view()); + auto result = evaluate_decimal_columns(test.op, + a->view(), + {256, 76, test.sa}, + b->view(), + {256, 76, test.sb}, + {256, 65, test.so}, + static_cast(active.data()), + stream.view(), + cudf::get_current_device_resource_ref()); + auto observed = download(result.values->view(), 32, stream.view()); + uint8_t errors[2]{}; + uint32_t validity{}; + REQUIRE( + cudaMemcpyAsync(errors, result.errors.data(), 2, cudaMemcpyDeviceToHost, stream.value()) == + cudaSuccess); + REQUIRE( + cudaMemcpyAsync( + &validity, result.values->view().null_mask(), 4, cudaMemcpyDeviceToHost, stream.value()) == + cudaSuccess); + stream.synchronize(); + CHECK(errors[0] == + static_cast(test.overflow ? decimal_error::out_of_range : decimal_error::none)); + CHECK(bool(validity & 1) == !test.overflow); + CHECK(errors[1] == 0); + CHECK((validity & 2) == 0); + if (!test.overflow) { + auto expected = decimal_fixture(test.expected); + for (int i = 0; i < 8; ++i) + CHECK(observed[0].words[i] == expected.words[i]); + } + } +} + +TEST_CASE("MO proven 64-to-128 cuDF paths cover the full signed physical input domain", + "[exact_decimal_gpu]") +{ + rmm::cuda_stream stream; + std::vector av{small(INT64_MAX), small(INT64_MIN), small(INT64_MAX)}; + std::vector bv{small(INT64_MAX), small(-1), small(INT64_MIN)}; + auto a = input_column({64, 18, 0}, av, 7, stream.view()); + auto b = input_column({64, 18, 0}, bv, 7, stream.view()); + for (auto op : {decimal_op::add, decimal_op::subtract, decimal_op::multiply}) { + auto result = evaluate_decimal_columns(op, + a->view(), + {64, 18, 0}, + b->view(), + {64, 18, 0}, + {128, 38, 0}, + nullptr, + stream.view(), + cudf::get_current_device_resource_ref()); + auto observed = download(result.values->view(), 16, stream.view()); + uint8_t errors[3]{}; + REQUIRE( + cudaMemcpyAsync(errors, result.errors.data(), 3, cudaMemcpyDeviceToHost, stream.value()) == + cudaSuccess); + stream.synchronize(); + for (int row = 0; row < 3; ++row) { + int64_t left[] = {INT64_MAX, INT64_MIN, INT64_MAX}; + int64_t right[] = {INT64_MAX, -1, INT64_MIN}; + __int128 lhs = left[row], rhs = right[row]; + __int128 value = op == decimal_op::add ? lhs + rhs + : op == decimal_op::subtract ? lhs - rhs + : lhs * rhs; + auto expected = load_coefficient(reinterpret_cast(&value), 16); + CHECK(errors[row] == 0); + for (int word = 0; word < 8; ++word) { + CHECK(observed[row].words[word] == expected.words[word]); + } + } + } +} diff --git a/test/cpp/embedding/test_native_control.cpp b/test/cpp/embedding/test_native_control.cpp index 10b288c3c7..d2770b525d 100644 --- a/test/cpp/embedding/test_native_control.cpp +++ b/test/cpp/embedding/test_native_control.cpp @@ -3,6 +3,7 @@ */ #include "embedding/buffer_budget.hpp" #include "embedding/control.hpp" +#include "embedding/input.hpp" #include "pipeline/completion_handler.hpp" #include @@ -215,6 +216,29 @@ TEST_CASE("native query identity is an owned opaque byte string", "[native_contr CHECK(q->contract->query_id == expected); } +TEST_CASE("native Decimal256 public and working precision have separate bounds", "[native_control]") +{ + fixture f; + auto q = f.create(); + char id = 'q'; + sirius_column output{34, 76, 2, 1, "wide", 4, 0}; + sirius_query_contract contract{ + sizeof(contract), SIRIUS_ABI_VERSION, 0, 0, &id, 1, {0}, &output, 1}; + REQUIRE_THROWS_AS(f.control.bind_query(q, contract), failure); + CHECK_FALSE(q->contract); + output.width = 65; + REQUIRE_NOTHROW(f.control.bind_query(q, contract)); + CHECK(q->contract->outputs[0].oid == 34); + CHECK(q->contract->outputs[0].width == 65); + sirius_input_column working{34, 76, 76, 1}; + REQUIRE_NOTHROW(validate_input_schema({&working, 1})); + working.width = 77; + REQUIRE_THROWS_AS(validate_input_schema({&working, 1}), failure); + working = {34, 15, 2, 0}; + REQUIRE_NOTHROW(validate_input_schema({&working, 1})); + CHECK(input_element_size(34) == 32); +} + TEST_CASE("native query identity retains its nonempty bounded contract", "[native_control]") { fixture f; diff --git a/test/cpp/embedding/test_native_gpu.cpp b/test/cpp/embedding/test_native_gpu.cpp index fb7022d707..3006b77fcb 100644 --- a/test/cpp/embedding/test_native_gpu.cpp +++ b/test/cpp/embedding/test_native_gpu.cpp @@ -423,11 +423,11 @@ TEST_CASE("native GPU scalar decoding preserves bits nulls constants and epochs" auto* gpu = manager.get_memory_spaces_for_tier(cucascade::memory::Tier::GPU).front(); auto pool = make_native_input_pool(*host); rmm::cuda_stream stream; - for (uint32_t oid : {10, 20, 21, 22, 23, 25, 26, 27, 28, 30, 31, 32, 33, 50, 52}) { + for (uint32_t oid : {10, 20, 21, 22, 23, 25, 26, 27, 28, 30, 31, 32, 33, 34, 50, 52}) { for (uint32_t cls : {SIRIUS_VECTOR_FLAT, SIRIUS_VECTOR_CONSTANT, SIRIUS_VECTOR_NULL}) { INFO("oid=" << oid << " vector_class=" << cls); auto width = input_element_size(oid); - sirius_input_column schema{oid, oid == 32 ? 18 : 38, 2, 1}; + sirius_input_column schema{oid, oid == 32 ? 18 : oid == 34 ? 65 : 38, 2, 1}; auto input = std::make_shared( 1, std::vector{schema}, std::stop_token{}, clock::now() + 30s); input->activate(pool); @@ -439,7 +439,12 @@ TEST_CASE("native GPU scalar decoding preserves bits nulls constants and epochs" if (oid == 50) lo += 719162; if (oid == 52) lo += 62135596800000000ULL; std::memcpy(raw.data() + row * width, &lo, std::min(width, 8)); - if (width == 16) std::memcpy(raw.data() + row * width + 8, &hi, 8); + if (width >= 16) std::memcpy(raw.data() + row * width + 8, &hi, 8); + if (width == 32) { + uint64_t mid = 0xfedcba987654321ULL, high = row % 2 ? UINT64_MAX : 7; + std::memcpy(raw.data() + row * width + 16, &mid, 8); + std::memcpy(raw.data() + row * width + 24, &high, 8); + } } if (cls == SIRIUS_VECTOR_FLAT) raw[physical * width + 7] = std::byte{128}; // row 63 NULL, omitted row-64 word is zero @@ -458,11 +463,24 @@ TEST_CASE("native GPU scalar decoding preserves bits nulls constants and epochs" auto table = convert_native_input(*unit, input->schema, *gpu, stream.view()); std::vector values(65 * width); uint32_t mask[3]{}; - cudaMemcpyAsync(values.data(), - table->view().column(0).head(), - values.size(), - cudaMemcpyDeviceToHost, - stream.value()); + if (oid == 34) { + std::vector lane(65); + for (int limb = 0; limb < 4; ++limb) { + cudaMemcpyAsync(lane.data(), + table->view().column(0).child(3 - limb).head(), + lane.size() * 8, + cudaMemcpyDeviceToHost, + stream.value()); + stream.synchronize(); + for (uint32_t row = 0; row < 65; ++row) + std::memcpy(values.data() + row * 32 + limb * 8, &lane[row], 8); + } + } else + cudaMemcpyAsync(values.data(), + table->view().column(0).head(), + values.size(), + cudaMemcpyDeviceToHost, + stream.value()); cudaMemcpyAsync(mask, table->view().column(0).null_mask(), sizeof(mask),