refactor(linalg)!: GPU elementwise arithmetic typed std::visit dispatch (#1013)

[yingjerkao]

Problem

GPU elementwise arithmetic (Add/Sub/Mul/Div/Mod/Cpr, in-place and out-of-place)
still dispatched through per-dtype-pair function tables (cuAri_ii, the 121
cuArithmetic_internal_* dispatchers, and the per-dtype cu*_internal_* kernels)
and a parallel type_promote_gpu_t / cy_typeid_gpu_v promotion hierarchy — the
legacy pattern #650/#938/#1013 are migrating away from. The GPU path had also
drifted from the CPU semantics established in #941/#1056: integer division stayed
integer, in-place ops kept the LHS dtype instead of promoting, and non-contiguous
tensor *= tensor / tensor /= tensor was rejected outright (#988).

Fix

Convert the GPU arithmetic dispatch to the same std::visit design as the CPU
path (#941):

• Out-of-place — Sub/Mul/Div/Mod route through a shared
  cuArithmeticDispatchGPU<op> (cuArithmeticDispatch.cuh): std::visit over
  as_storage_variant() + type_promote_t, with the CUDA-native complex
  representation confined to the kernel boundary by a new to_cuda_t trait
  (cuTypeCvt.hpp). Cpr keeps its Bool-output kernels but drops the
  type_promote_gpu_t switch. Add follows in the removal commit (its dispatch
  shared a file with the legacy per-dtype defs).
• In-place — iAdd/iSub/iMul/iDiv use cuiArithmeticDispatch.cuh, mirroring
  the CPU DispatchInplaceArithmeticCPU: promote Lt's storage via
  storage_as_type_or_replace, then run the op at Lt's physical layout. Removes
  the tmpo clone dance and the narrow-LHS / zero-extent special-cases.
• Legacy removed — cuAri_ii, cuArithmetic_internal.{cu,hpp}, the four dead
  cu{Sub,Mul,Div,Mod}_internal.cu, the per-dtype declarations in the
  cu*_internal.hpp, and the obsolete bet.py generator. Net −35k lines.
  type_promote_gpu_t/cy_typeid_gpu_v/Type_list_gpu are kept — still used
  by Tensor.hpp's GPU storage variant.
• Build fix (pre-existing, first commit) — the CUDA-13 cccl/ include dir is
  now exposed to in-tree GPU consumers (gpu_test_main/benchmarks), which
  otherwise couldn't compile <cuda/std/complex> at all.

⚠️ Behavior changes (GPU now matches CPU)

These are marked ! on the commits and warrant a numerical-correctness look:

1. GPU true division — int / int → Double (was integer division).
2. GPU in-place promotion — int16 += int64 → the tensor becomes Int64;
   int64 /= int64 → Double.
3. Non-contiguous tensor *= tensor / tensor /= tensor now works in place
   (previously threw — the legacy cuMul/cuDiv kernels ignored the layout mappers,
   GPU: scalar in-place ops throw on non-contiguous tensors for *= and /=, and do a per-call H2D copy #988).

Testing

On an RTX 4070 Ti SUPER (CUDA 13):

• gpu_test_main — 68 Add tests + all-dtype-pair sweeps for Sub/Mul/Div/Mod/Cpr,
  in-place regression for all four ops, and the new
  InplacePromote_test.cpp (promotion + non-contiguous, validated against CPU) —
  all pass.
• compute-sanitizer memcheck — 0 errors on both in-place and out-of-place.
• CPU — debug build green, 211 CPU arithmetic/in-place/rank tests pass (CPU
  behavior unchanged).
• clang-format v14 clean.

Relates to #1013. Stacked on #1056 (base branch refactor/cpu-linalg-typed-dispatch).

🤖 Generated with Claude Code

[gemini-code-assist]

Code Review

This pull request refactors the GPU elementwise binary arithmetic and comparison dispatch to use std::visit over as_storage_variant(), aligning the GPU implementation with the CPU design and removing legacy per-dtype function tables. This change correctly handles true division on GPU and supports in-place operations on non-contiguous tensors. However, two critical bugs were identified in cuArithmeticDispatch.cuh and cuiArithmeticDispatch.cuh where device-allocated pointers are directly dereferenced and written to from host code during shape accumulator calculations. These values must be computed in host vectors first and then copied to the device using cudaMemcpy to prevent segmentation faults or undefined behavior.

Important

The consumer version of Gemini Code Assist on GitHub is being sunset. Starting June 18, 2026, new organization installations will be blocked, and all code review activity will officially cease on July 17, 2026.
For more details on the timeline and next steps, please review the Help Documentation.

[gemini-code-assist]

Directly dereferencing and writing to device pointers (m_accu_shape, m_old_accu_shapeL, m_old_accu_shapeR) from host code is a critical bug. In standard CUDA, device memory allocated via cuCalloc_gpu (which wraps cudaMalloc) is not host-writable and will cause a segmentation fault or undefined behavior.

To fix this, compute the values in host vectors first, and then copy them to the device pointers using cudaMemcpy. Additionally, since the entire buffers are overwritten immediately, you can use cuMalloc_gpu instead of cuCalloc_gpu to avoid redundant zero-initialization.

            cytnx_uint64 *m_accu_shape = reinterpret_cast<cytnx_uint64 *>(
              utils_internal::cuMalloc_gpu(shape.size() * sizeof(cytnx_uint64)));
            cytnx_uint64 *m_old_accu_shapeL = reinterpret_cast<cytnx_uint64 *>(
              utils_internal::cuMalloc_gpu(shape.size() * sizeof(cytnx_uint64)));
            cytnx_uint64 *m_old_accu_shapeR = reinterpret_cast<cytnx_uint64 *>(
              utils_internal::cuMalloc_gpu(shape.size() * sizeof(cytnx_uint64)));
            cytnx_uint64 *m_invmapper_L = reinterpret_cast<cytnx_uint64 *>(
              utils_internal::cuMalloc_gpu(invmapper_L.size() * sizeof(cytnx_uint64)));
            cytnx_uint64 *m_invmapper_R = reinterpret_cast<cytnx_uint64 *>(
              utils_internal::cuMalloc_gpu(invmapper_R.size() * sizeof(cytnx_uint64)));

            std::vector<cytnx_uint64> accu_shape(shape.size());
            std::vector<cytnx_uint64> old_accu_shapeL(shape.size());
            std::vector<cytnx_uint64> old_accu_shapeR(shape.size());

            cytnx_uint64 tmp1 = 1, tmp2 = 1, tmp3 = 1;
            for (cytnx_uint64 i = 0; i < shape.size(); i++) {
              accu_shape[shape.size() - 1 - i] = tmp1;
              tmp1 *= shape[shape.size() - 1 - i];

              old_accu_shapeL[shape.size() - 1 - i] = tmp2;
              tmp2 *= shape[invmapper_L[shape.size() - 1 - i]];

              old_accu_shapeR[shape.size() - 1 - i] = tmp3;
              tmp3 *= shape[invmapper_R[shape.size() - 1 - i]];
            }

            checkCudaErrors(cudaMemcpy(m_accu_shape, accu_shape.data(),
                                       sizeof(cytnx_uint64) * shape.size(),
                                       cudaMemcpyHostToDevice));
            checkCudaErrors(cudaMemcpy(m_old_accu_shapeL, old_accu_shapeL.data(),
                                       sizeof(cytnx_uint64) * shape.size(),
                                       cudaMemcpyHostToDevice));
            checkCudaErrors(cudaMemcpy(m_old_accu_shapeR, old_accu_shapeR.data(),
                                       sizeof(cytnx_uint64) * shape.size(),
                                       cudaMemcpyHostToDevice));
            checkCudaErrors(cudaMemcpy(m_invmapper_L, &invmapper_L[0],
                                       sizeof(cytnx_uint64) * invmapper_L.size(),
                                       cudaMemcpyHostToDevice));
            checkCudaErrors(cudaMemcpy(m_invmapper_R, &invmapper_R[0],
                                       sizeof(cytnx_uint64) * invmapper_R.size(),
                                       cudaMemcpyHostToDevice));

[yingjerkao]

Thanks, but this is a false positive — no host-write-to-device-memory happens here.

cuCalloc_gpu/cuMalloc_gpu allocate unified (managed) memory via cudaMallocManaged, not cudaMalloc — see src/backend/utils_internal_gpu/cuAlloc_gpu.cu:11 (cuCalloc_gpu) and :17 (cuMalloc_gpu). Managed memory is host-accessible, so writing m_accu_shape/m_old_accu_shapeL/m_old_accu_shapeR from host before the kernel launch is valid; the launch orders the host writes ahead of the device reads. (The invmapper buffers use cudaMemcpy only because they originate from a host std::vector, not because managed memory requires it.) This is the existing convention for the non-contiguous kernels — unchanged from the legacy per-dtype code this PR replaces.

Verified on an RTX 4070 Ti SUPER (sm_89): the non-contiguous GPU tests that exercise this exact path pass, and compute-sanitizer --tool memcheck reports 0 errors for both the out-of-place and in-place dispatch. Leaving as-is.

[gemini-code-assist]

Directly dereferencing and writing to device pointers (m_accu_shape, m_old_accu_shapeL, m_old_accu_shapeR) from host code is a critical bug. In standard CUDA, device memory allocated via cuCalloc_gpu (which wraps cudaMalloc) is not host-writable and will cause a segmentation fault or undefined behavior.

To fix this, compute the values in host vectors first, and then copy them to the device pointers using cudaMemcpy. Additionally, since the entire buffers are overwritten immediately, you can use cuMalloc_gpu instead of cuCalloc_gpu to avoid redundant zero-initialization.

          cytnx_uint64 *m_accu_shape = reinterpret_cast<cytnx_uint64 *>(
            utils_internal::cuMalloc_gpu(shape.size() * sizeof(cytnx_uint64)));
          cytnx_uint64 *m_old_accu_shapeL = reinterpret_cast<cytnx_uint64 *>(
            utils_internal::cuMalloc_gpu(shape.size() * sizeof(cytnx_uint64)));
          cytnx_uint64 *m_old_accu_shapeR = reinterpret_cast<cytnx_uint64 *>(
            utils_internal::cuMalloc_gpu(shape.size() * sizeof(cytnx_uint64)));
          cytnx_uint64 *m_invmapper_L = reinterpret_cast<cytnx_uint64 *>(
            utils_internal::cuMalloc_gpu(invmapper_L.size() * sizeof(cytnx_uint64)));
          cytnx_uint64 *m_invmapper_R = reinterpret_cast<cytnx_uint64 *>(
            utils_internal::cuMalloc_gpu(invmapper_R.size() * sizeof(cytnx_uint64)));

          std::vector<cytnx_uint64> accu_shape(shape.size());
          std::vector<cytnx_uint64> old_accu_shapeL(shape.size());
          std::vector<cytnx_uint64> old_accu_shapeR(shape.size());

          cytnx_uint64 tmp1 = 1, tmp2 = 1, tmp3 = 1;
          for (cytnx_uint64 i = 0; i < shape.size(); i++) {
            accu_shape[shape.size() - 1 - i] = tmp1;
            tmp1 *= shape[shape.size() - 1 - i];

            old_accu_shapeL[shape.size() - 1 - i] = tmp2;
            tmp2 *= shape[invmapper_L[shape.size() - 1 - i]];

            old_accu_shapeR[shape.size() - 1 - i] = tmp3;
            tmp3 *= shape[invmapper_R[shape.size() - 1 - i]];
          }

          checkCudaErrors(cudaMemcpy(m_accu_shape, accu_shape.data(),
                                     sizeof(cytnx_uint64) * shape.size(),
                                     cudaMemcpyHostToDevice));
          checkCudaErrors(cudaMemcpy(m_old_accu_shapeL, old_accu_shapeL.data(),
                                     sizeof(cytnx_uint64) * shape.size(),
                                     cudaMemcpyHostToDevice));
          checkCudaErrors(cudaMemcpy(m_old_accu_shapeR, old_accu_shapeR.data(),
                                     sizeof(cytnx_uint64) * shape.size(),
                                     cudaMemcpyHostToDevice));
          checkCudaErrors(cudaMemcpy(m_invmapper_L, &invmapper_L[0],
                                     sizeof(cytnx_uint64) * invmapper_L.size(),
                                     cudaMemcpyHostToDevice));
          checkCudaErrors(cudaMemcpy(m_invmapper_R, &invmapper_R[0],
                                     sizeof(cytnx_uint64) * invmapper_R.size(),
                                     cudaMemcpyHostToDevice));

[yingjerkao]

Thanks, but this is a false positive — no host-write-to-device-memory happens here.

cuCalloc_gpu/cuMalloc_gpu allocate unified (managed) memory via cudaMallocManaged, not cudaMalloc — see src/backend/utils_internal_gpu/cuAlloc_gpu.cu:11 (cuCalloc_gpu) and :17 (cuMalloc_gpu). Managed memory is host-accessible, so writing m_accu_shape/m_old_accu_shapeL/m_old_accu_shapeR from host before the kernel launch is valid; the launch orders the host writes ahead of the device reads. (The invmapper buffers use cudaMemcpy only because they originate from a host std::vector, not because managed memory requires it.) This is the existing convention for the non-contiguous kernels — unchanged from the legacy per-dtype code this PR replaces.

Verified on an RTX 4070 Ti SUPER (sm_89): the non-contiguous GPU tests that exercise this exact path pass, and compute-sanitizer --tool memcheck reports 0 errors for both the out-of-place and in-place dispatch. Leaving as-is.

[chatgpt-codex-connector]

💡 Codex Review

Here are some automated review suggestions for this pull request.

Reviewed commit: c13ebcbfba

ℹ️ About Codex in GitHub

Your team has set up Codex to review pull requests in this repo. Reviews are triggered when you

• Open a pull request for review
• Mark a draft as ready
• Comment "@codex review".

If Codex has suggestions, it will comment; otherwise it will react with 👍.

Codex can also answer questions or update the PR. Try commenting "@codex address that feedback".

[chatgpt-codex-connector]

Assert independent expected arithmetic results

CLAUDE.md's test guidance for arithmetic/dtype changes says tests should check independent expected values, but this new coverage uses the CPU in-place path as the oracle (ApplyInplace(op, cpu_l, cpu_r)) and only compares GPU against that result. If the shared typed promotion or non-contiguous indexing logic is wrong in the same way on both CPU and GPU, these tests still pass; please add a few literal expected dtype/value cases such as Int16 += Int64, Int64 /= Int64, and a non-contiguous tensor case.

Useful? React with 👍 / 👎.

[yingjerkao]

Good catch — addressed in 623c88b. Added InplacePromoteTest.gpu_inplace_literal_expected with hand-computed expected dtype and values, independent of any Cytnx execution path:

• Int16 += Int64 -> Int64 = [11,22,33]
• Int64 /= Int64 -> Double = [0.5,1.0,1.5] (true division)
• contiguous LHS + permuted RHS (Add), and permuted LHS + contiguous RHS (Sub) — so Lidx != Ridx is actually exercised and a wrong physical zip would fail.

Inputs are built with CPU arange then moved to the GPU, so the test data itself doesn't depend on GPU kernels.

[ianmccul]

#1003

[ianmccul]

This is a valuable correctness improvement: typed dispatch removes roughly 35,000 lines of generated dtype-pair machinery, fixes GPU true division and in-place promotion, and eliminates several unsafe pointer/dtype combinations. However, the non-contiguous path needs a more careful review than the current tests provide.

The new implementation handles non-contiguous operands by carrying forward the old permutation-decoding algorithm. For every logical output index, the kernel reconstructs a coordinate vector, applies invmapper_L and invmapper_R, and computes physical offsets. The out-of-place kernel writes out[idx]; the in-place kernel writes out[Lidx]. This is plausibly correct for Cytnx current permutation-only view model, but it is not the stride/layout design proposed in #1003.

The Gemini claim that the host writes themselves are invalid is a false positive: cuCalloc_gpu uses cudaMallocManaged, so those writes are legal. But legal managed-memory access does not make this a good implementation.

Findings

1. The new non-contiguous test barely exercises the mapper logic. Both operands use the same permutation:

   gpu_l = ...permute({1, 0});
   gpu_r = ...permute({1, 0});

   Consequently Lidx == Ridx, and even an incorrect physical-buffer zip can pass. The test then compares against the CPU path, which uses essentially the same mapper algorithm, so a shared indexing error can pass on both devices. Please add independently checked cases with literal expected values, including:

   • contiguous LHS and permuted RHS;
   • permuted LHS and contiguous RHS;
   • differently permuted operands with the same logical shape;
   • promoted in-place output with different LHS/RHS layouts;
   • at least subtraction and division, where an operand-order/indexing error is visible.

2. Dynamic shared-memory use imposes an undocumented rank-dependent launch limit. Each block requests

   512 * rank * sizeof(uint64_t) = 4096 * rank bytes
   

   for the per-thread coordinate vectors. Rank 13 requests 53,248 bytes, already above the common 48 KiB default dynamic-shared-memory limit. There is no rank check, launch-configuration adjustment, or opt-in shared-memory attribute, so sufficiently high-rank non-contiguous arithmetic can fail at kernel launch. At minimum this needs an explicit supported-rank check or a block size derived from rank; preferably the per-thread coordinate vector should not live in dynamic shared memory at all.

3. The old per-call metadata-allocation model remains. Every non-contiguous operation performs five managed-memory allocations, two copies, a kernel launch, and five immediate frees. This adds substantial launch overhead and prevents clean asynchronous lifetime management. Extents and strides are tiny launch metadata and should be preprocessed into a compact trivially-copyable layout passed as a kernel argument, as described in #1003.

4. Traversal remains duplicated. Nearly identical coordinate-decoding code now exists in cuArithmeticDispatch.cuh, cuiArithmeticDispatch.cuh, and cuCpr_dispatch.cu, in addition to the analogous CPU traversal introduced by #1056. The operation functor is shared within one GPU arithmetic family, but the actual elementwise traversal is not yet centralized across out-of-place, in-place, comparison, CPU, and GPU paths.

5. The generic out-of-place kernel is not uniformly reachable. The front ends for out-of-place GPU Mul, Div, Mod, and Cpr still reject non-contiguous tensor/tensor operands. Only Add and Sub, plus the newly enabled in-place paths, reach mapper-aware kernels. The PR description should not imply general non-contiguous GPU elementwise support.

I do not think the full #1003 architecture must be implemented in this PR before the typed-dispatch cleanup can land. The deletion of the generated tables and the dtype correctness fixes are independently worthwhile. However, the weak non-contiguous tests and the rank-dependent kernel-launch failure are correctness issues that should be addressed before approval.

After this lands, #1003 should be broadened explicitly to a shared CPU/GPU elementwise framework: common operation functors and dtype policies, common host-side stride/layout preprocessing, and backend-specific CPU/CUDA execution consuming the same compact layout. The mapper kernels in this PR should be treated as transitional code, not the final non-contiguous implementation.

Written by OpenAI Codex and posted at @ianmccul request.

[ianmccul]

One additional point about the in-place non-contiguous kernel: even if the current inverse-mapper representation is retained temporarily, there is no reason to map both operands.

For an in-place operation, the LHS is also the output, so its physical storage order should define the traversal:

const auto lhs_idx = physical_idx;
const auto rhs_idx = map_lhs_physical_to_rhs_physical(physical_idx);
out[lhs_idx] = op(lhs[lhs_idx], rhs[rhs_idx]);

The current kernel instead iterates a logical linear index, decodes logical coordinates, maps them to Lidx, maps them again to Ridx, and finally writes out[Lidx]. Calculating Lidx is redundant. Worse, this logical-order traversal can make both the LHS reads and output writes uncoalesced, even though those two buffers necessarily have the same layout.

Iterating LHS physical order guarantees coalesced LHS reads and output writes; only the RHS may require an index transformation. The two permutations can be composed during host-side preprocessing into a direct LHS-physical-to-RHS-physical mapping. If both operands have the same layout, that mapping is identity and the operation reduces to a normal physical-buffer zip.

This is independent of the broader stride-layout redesign in #1003. It is a straightforward simplification of the existing mapper approach and is another reason the test where both operands use the same permutation is insufficient: that case should not need mapping at all.

Written by OpenAI Codex and posted at @ianmccul request.

[ianmccul]

A performance concern follows directly from the identical-layout test case.

The front end chooses the simple physical-buffer kernel only when both tensors are individually contiguous:

if (Lt.is_contiguous() && Rt.is_contiguous()) {
  // physical zip
} else {
  // full mapper kernel
}

Two tensors carrying the same non-contiguous layout therefore take the most expensive mapper path, even though an in-place operation needs no mapping at all. If

Lt.shape() == Rt.shape() && Lt.invmapper() == Rt.invmapper()

then corresponding logical elements occupy corresponding physical offsets, so the correct implementation is simply:

for (std::size_t i = 0; i < Lt.size(); ++i) {
  lhs_data[i] = op(lhs_data[i], rhs_data[i]);
}

This remains true for non-commutative operations such as subtraction and division; it depends only on identical layouts, not on the operation.

The new test constructs both operands with the same permute({1, 0}), so it should exercise this physical-zip fast path. Instead it currently exercises the full path with five managed allocations, two copies, five frees, logical-index decoding, two mapper applications, integer division/modulo per axis, and rank-proportional dynamic shared memory. A kernel that ignored both mappers would pass that test, which is why it does not validate mapper correctness.

Before treating this as an acceptable non-contiguous implementation, I suggest adding benchmarks for at least:

1. contiguous a += b;
2. identically permuted a += b;
3. contiguous LHS with permuted RHS;
4. differently permuted operands with the same logical shape;
5. explicit contiguous() copies followed by the contiguous kernel, including the copy cost.

Run these over multiple tensor sizes and ranks. Case 2 should be close to case 1 because the required kernel is identical. With the current branch selection it instead pays the complete mapper overhead. For some sizes, even materializing contiguous copies may be faster than the nominal no-copy path.

The same-layout physical-zip check is a small improvement that can be made independently of #1003. The broader benchmark should then guide the stride-layout replacement rather than preserving this mapper path without performance evidence.

Written by OpenAI Codex and posted at @ianmccul request.

[yingjerkao]

CI is green ✅ — all checks pass on the current head: Formatting Check, the four BuildWheel jobs (ubuntu-24.04, ubuntu-24.04-arm, macos-14, macos-15-intel — each builds and runs ctest + pytest), check, and the docs build. ReleasePyPI/PublishNightlyAnaconda are release-only and skipped as expected.

Local validation on an RTX 4070 Ti SUPER (CUDA 13, sm_89) is also green: full GPU suite for all six ops (out-of-place + in-place, incl. promotion and non-contiguous coverage), compute-sanitizer --tool memcheck with 0 errors, and pytest pytests/ = 248 passed / 1 skipped.

Re: the two "critical" findings — both are false positives, addressed on the inline threads: cuCalloc_gpu/cuMalloc_gpu allocate managed memory (cudaMallocManaged, cuAlloc_gpu.cu:11,17), so the host writes to m_accu_shape/m_old_accu_shapeL/R are valid. That's confirmed by the non-contiguous tests passing and the clean memcheck run.

[ianmccul]

Managed memory makes these host writes legal, but it does not make this a good GPU metadata-transfer path.

After the CPU writes a cudaMallocManaged allocation, the first GPU access may fault and migrate/map the allocation at page granularity while the kernel is stalled. For these arrays the useful payload is only rank * sizeof(uint64_t) bytes, so the fixed fault and page-management overhead is disproportionately large. NVIDIA's own Unified Memory guidance recommends explicit prefetching or copying when the access pattern and destination are known in advance.

There is also a redundant initialization here. cuCalloc_gpu performs:

cudaMallocManaged(...);
cudaMemset(...);

but the CPU immediately overwrites every element of m_accu_shape, m_old_accu_shapeL, and m_old_accu_shapeR.

So yes: the Gemini report was wrong to classify the host writes as invalid, and passing tests plus a clean memcheck confirm legality. They do not establish that demand-paging five tiny managed allocations is efficient.

For these small, predictable arrays, an explicit transfer would generally be preferable to first-touch migration. However, five separate synchronous cudaMemcpy calls would still be a poor final design. The metadata should be packed into one compact trivially-copyable object or buffer and either:

1. passed directly as kernel launch parameters where practical;
2. copied once with cudaMemcpyAsync to one device allocation; or
3. at minimum, explicitly prefetched to the active GPU before launch.

The current non-contiguous path performs five managed allocations, two explicit copies, a kernel launch, and five immediate frees for metadata totaling only:

5 * rank * sizeof(uint64_t)

At rank 12 that is 480 bytes. Allocation, synchronization, page-fault, and free overhead will dominate moving the actual metadata. This is independent of whether the managed-memory access is technically valid.

Written by OpenAI Codex and posted at @ianmccul request.

[ianmccul]

The overhead of copying a few hundred bytes - or even a few thousand bytes - of data into the kernel as fixed-size parameters is essentially zero, it is immeasurably small in the kernel startup procedure. A suitable implementation would resemble

constexpr int max_rank = 32;

struct LayoutMetadata {
  std::array<std::uint64_t, max_rank> accu_shape;
  std::array<std::uint64_t, max_rank> lhs_strides;
  std::array<std::uint64_t, max_rank> rhs_strides;
  std::array<std::uint64_t, max_rank> lhs_mapper;
  std::array<std::uint64_t, max_rank> rhs_mapper;
  int rank;
};

__global__ void kernel(const __grid_constant__ LayoutMetadata layout, ...);

That is about 1.3 KiB for max_rank == 32. __grid_constant__ prevents the compiler from producing a private per-thread copy if the parameter’s address is taken. CUDA __grid_constant__ documentation (https://docs.nvidia.com/cuda/cuda-programming-guide/05-appendices/cpp-language-extensions.html)

NVIDIA documents a total parameter limit of 32,764 bytes on modern devices; older pre-Volta devices have a 4 KiB limit. CUDA parameter documentation (https://docs.nvidia.com/cuda/cuda-programming-guide/05-appendices/cpp-language-support.html)

set max_rank = 52 if you're paranoid. But since the dynamic memory used by the kernel currently imposes a maximum rank of 12, anything >= 12 is fine here.

[yingjerkao]

Thanks @ianmccul — the two "before approval" items are addressed and the architectural items are filed as follow-up #1069.

Finding 1 — weak non-contiguous tests (fixed, 623c88b). The old case used the same permutation on both operands (Lidx == Ridx), so a wrong physical zip would pass. Replaced with independent literal-value cases: Int16 += Int64 -> [11,22,33], Int64 /= Int64 -> [0.5,1.0,1.5], contiguous LHS + permuted RHS (Add), and permuted LHS + contiguous RHS (Sub). Values/dtypes are hand-computed, not compared against the CPU path.

Finding 2 (rank-dependent shared-memory launch limit) + the architecture items (per-call metadata model / __grid_constant__ layout, LHS-physical traversal, same-layout physical-zip fast path, uniform non-contiguous reachability for Mul/Div/Mod/Cpr, traversal de-duplication, and the benchmark matrix) are captured in #1069 as the #1003 broadening — since, as you noted, these mapper kernels are transitional and the full redesign needn't block the typed-dispatch cleanup. The rank limit is pre-existing (inherited verbatim from the legacy kernels), so it isn't a regression introduced here, but it's tracked to fix.

Separately, while writing the literal tests I hit what looks like a pre-existing, unrelated bug: GPU arange(start, end, step) with step != 1 seems to ignore the step — arange(10, 40, 10) on Device.cuda produced [10, 11, 12] instead of [10, 20, 30] (CPU is correct). I worked around it by building the test inputs with CPU arange + .to(cuda). Happy to file/investigate that separately.

[codecov]

Codecov Report

✅ All modified and coverable lines are covered by tests.
✅ Project coverage is 72.14%. Comparing base (2adb13b) to head (db84d56).
⚠️ Report is 6 commits behind head on master.
✅ All tests successful. No failed tests found.

Additional details and impacted files

@@            Coverage Diff             @@
##           master    #1067      +/-   ##
==========================================
+ Coverage   72.11%   72.14%   +0.03%     
==========================================
  Files         225      225              
  Lines       28207    28183      -24     
  Branches       71       71              
==========================================
- Hits        20341    20333       -8     
+ Misses       7845     7829      -16     
  Partials       21       21              

Flag	Coverage Δ
cpp	72.25% <100.00%> (+0.03%)	⬆️
python	64.13% <ø> (ø)

Flags with carried forward coverage won't be shown. Click here to find out more.

Components	Coverage Δ
C++ backend	71.60% <100.00%> (+0.03%)	⬆️
Python bindings	75.80% <ø> (ø)
Python package	64.13% <ø> (ø)

Files with missing lines	Coverage Δ
src/backend/linalg_internal_interface.cpp	100.00% <ø> (ø)
src/linalg/Div.cpp	86.91% <100.00%> (ø)
src/linalg/iAdd.cpp	95.23% <100.00%> (+13.75%)	⬆️
src/linalg/iDiv.cpp	95.23% <100.00%> (+13.75%)	⬆️
src/linalg/iMul.cpp	95.23% <100.00%> (+13.75%)	⬆️
src/linalg/iSub.cpp	95.23% <100.00%> (+13.75%)	⬆️

... and 1 file with indirect coverage changes

---

Continue to review full report in Codecov by Harness.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update 2adb13b...db84d56. Read the comment docs.

🚀 New features to boost your workflow:

• 📦 JS Bundle Analysis: Save yourself from yourself by tracking and limiting bundle sizes in JS merges.

[ianmccul]

Thanks, the literal expected-value tests in 623c88b2 address finding 1. Given that the rank-dependent shared-memory limit is inherited from the existing GPU mapper kernels, I am also fine with leaving that and the broader kernel/layout redesign to #1069.

I still see two consistency/build issues:

1. The new visitors state that a genuine Tensor RHS promotes the LHS storage, including a genuine rank-0 tensor. However, iAdd, iSub, iMul, and iDiv still reject every real-LHS/complex-RHS combination before reaching the visitor:

   cytnx_error_msg(!Type.is_complex(Lt.dtype()) && Type.is_complex(Rt.dtype()), ...);

   That guard made sense when the operation had to write into the existing real buffer. It no longer follows once storage_as_type_or_replace<TO> can install complex output storage. The deliberate weak-scalar behavior from fix(Tensor): scalar in-place arithmetic mutates storage in place #980 can still reject real op= complex scalar, but a genuine complex Tensor RHS should promote according to this PR's stated rule. This needs CPU and GPU tests.

2. The CUDA-13 CCCL build fix is only exported through:

   $<BUILD_INTERFACE:${_cytnx_cccl_dir}>

   An installed GPU Cytnx package still exports UNI_GPU, and its public Type.hpp includes <cuda/std/complex>, but the installed Cytnx::cytnx target receives no CCCL include directory. The new DownstreamFindPackage job uses the openblas-cpu preset, so it does not test this case. The install-side package config needs to recover/propagate the CCCL include path, or a CUDA-enabled install-tree consumer test should demonstrate that another exported CUDA target supplies it.

Apart from those, I am fine with the typed-dispatch cleanup landing before the single-kernel redesign in #1069.

Written by OpenAI Codex and posted at @ianmccul request.

[yingjerkao]

Thanks @ianmccul — item 1 is fixed in aab98bf3, with the two stale GPU tests it surfaced fixed in f2718151.

Item 1 — genuine complex Tensor RHS now promotes a real LHS in-place. The guard in iAdd/iSub/iMul/iDiv is narrowed so it fires only for a complex python weak scalar:

cytnx_error_msg(
  rhs_is_weak_scalar && !Type.is_complex(Lt.dtype()) && Type.is_complex(Rt.dtype()),
  "[iAdd] Cannot perform real += complex-scalar in-place: weak-scalar "
  "semantics preserve the real LHS dtype, so the complex scalar cannot be stored. "
  "Use a complex LHS, or a genuine complex tensor RHS to promote.%s", "\n");

A genuine complex Tensor RHS (including a rank-0 tensor) now promotes Lt's storage to complex via storage_as_type_or_replace<TO>, matching this PR's stated rule and the out-of-place operator. The deliberate #980 weak-scalar rejection is preserved. The guard stays device-independent on purpose: the GPU kernel's complex-into-real branch silently returns zero rather than throwing, so it relies on this host-side check — I noted that in the code comments and in cuiArithmeticDispatch.cuh.

Tests (independent, hand-computed expected values — not compared against a CPU path):

• CPU (TypedArithmeticCoverage_test.cpp): InplaceRealOpComplexTensorPromotes covers all four ops plus a genuine rank-0 complex RHS; InplaceRealOpComplexWeakScalarRejected covers the scalar operators, an explicit rhs_is_weak_scalar=true call (throws), and a genuine tensor RHS (promotes). Both new tests were confirmed to fail on the pre-fix guard.
• GPU (InplacePromote_test.cpp): new gpu_inplace_real_op_complex_weak_scalar_rejected, plus real←complex promotion added to the sweep and the literal-value cases.

Two stale GPU tests fixed (f2718151). Running the GPU suite locally surfaced two GPU-only tests in tests/gpu/Tensor_test.cpp that still asserted the pre-#1013 in-place behavior and had been red for the whole PR (CI has no GPU runner, so it never flagged them): GpuRankZeroTensorRhsInplacePreserveDtype (a rank-0 genuine tensor RHS now promotes) and GpuNoncontigTensorTensorMulDivThrows (non-contiguous *=//= now works). Updated both to the new behavior with independent expected values. The in-place/arithmetic GPU suites are green locally on the RTX 4070 Ti SUPER.

Item 2 — CCCL install-interface — not yet addressed. You're right that the fix is only exported through $<BUILD_INTERFACE:${_cytnx_cccl_dir}>, so an installed GPU package whose public Type.hpp pulls in <cuda/std/complex> gets no CCCL include dir, and DownstreamFindPackage (openblas-cpu) doesn't exercise it. I'll handle this in a follow-up commit on this branch — propagating the CCCL include path through the install interface — before this is ready to approve. Will ping you when it's up.

[yingjerkao]

I will merge this and move to #1003