Is BJ maxBlockSize > 1 OOM on SYCL a known issue? Is there a fix or workaround planned for Ginkgo 2.0?

[heikogleu-dev]

Summary

Comprehensive testing of Ginkgo 1.10 SYCL preconditioners on Intel Arc Pro B70 Pro (Battlemage G31, 32 GB) reveals critical issues for CFD use cases.
Hardware (performing well — software is the bottleneck)

GPU: Intel Arc Pro B70 Pro (0xe223, 32 GB GDDR6, Ubuntu 26.04)
FP64: 1335 GFLOPS (93% of spec) ✅
VRAM: 530 GB/s sustained (87% of spec) ✅
oneAPI 2026.0.0, Compute Runtime 26.05.37020

Critical Issue 1: BJ maxBlockSize > 1 → gko::AllocationError (OOM)

Any maxBlockSize > 1 causes immediate OOM during generate(). Confirmed on: maxBlockSize = 2, 4, 8, 16, 32 — all fail identically. Only maxBlockSize = 1 (diagonal scaling) works.

Test configuration: 34M cells, 8 MPI ranks (4.25M cells/rank, ~270 MB/rank) dmesg: xe VM worker error: -12 (ENOMEM)

Suspected cause: SYCL BJ generate workspace is O(N × BS²) instead of O(N × BS). This makes maxBlockSize > 1 unusable on any reasonably-sized mesh.
Critical Issue 2: IC → gko::NotImplemented

sparselib_ic not implemented for SYCL backend:

dpcpp/factorization/ic_kernels.dp.cpp:21: sparselib_ic is not implemented

Critical Issue 3: ICT → DEVICE_LOST (GPU hardware hang)

ICT causes immediate GPU hang requiring driver reset.

sycl::_V1::exception: level_zero backend failed with error: 20
(UR_RESULT_ERROR_DEVICE_LOST)

Critical Issue 4: Multigrid → OOM in PGM coarsening

PGM coarsening crashes in generate_local() → Csr::clone() → raw_copy_to() with DEVICE_LOST (underlying ENOMEM). The one solve that completes diverges (Final residual > Initial residual).

Stack: gko::multigrid::Pgm::generate_local() → gko::matrix::Csr::apply_impl() → DpcppExecutor::raw_copy_to() → DEVICE_LOST
Impact on CFD

With only BJ maxBlockSize=1 available (point-Jacobi = diagonal scaling), the pressure solver never converges for 34M-cell meshes — always hits the 200-iteration cap. CPU GAMG (algebraic multigrid) converges in 5-10 iterations. GPU is 1.5× slower than CPU under fair comparison.
Full Documentation

https://github.com/heikogleu-dev/Openfoam13---GPU-Offloading-Intel-B70-Pro

Findings 02, 05, 08 contain stack traces and reproduction steps.

[greole]

Could you improve the description of the issue, please? Are you solving pressure and momentum on the device? How many devices are you using? Can you estimate the memory footprint per cell and rank?

[heikogleu-dev]

Update: Comprehensive testing complete — additional findings

Following up with results from exhaustive testing. All plausible
fvSolution parameters are now documented.

New Finding: GINKGO_JACOBI_FULL_OPTIMIZATIONS=ON does NOT fix BJ OOM

Rebuilt libOGL.so with:

• GINKGO_JACOBI_FULL_OPTIMIZATIONS=ON
• -fp-model=precise

Result: BS=2 and BS=4 still crash with gko::AllocationError identically.
The optimization flag does not address the underlying SYCL workspace
allocation bug. Confirmed on CR 26.05 post-reboot — not a stack issue.

New Finding: ISAI sparsityPower=3 → int32 overflow at 34M cells

With scaling=-1.0 (correct for SPD systems per OGL docs) and
sparsityPower=3:
SYCL error: Provided range and/or offset does not fit in int

The ISAI sparsity pattern generation uses int32 internally, which
overflows at ~4.25M cells/rank. sparsityPower=1 runs but diverges
(Final residual > Initial residual on first solve).

New Finding: IC NotImplemented with scaling=-1.0 — unchanged

The OGL README correctly notes: "IC/ISAI require positive diagonal —
scale pressure equation by -1.0". Applied scaling=-1.0, retested IC:
still gko::NotImplemented. The scaling is correct but cannot overcome
the missing SYCL implementation.

New Finding: Kernel-launch latency is 5.6 µs — NOT ~100 µs

Direct measurement on B70 Pro / Level Zero:
Trivial kernel + sync: 5.6 µs/launch (CUDA-par)
Trivial kernel async: 1.5 µs/launch

This means kernel launches are NOT the CFD bottleneck. Our initial
documentation was wrong.

Profiling: GPU idle 66% of wall clock

Direct measurement using xe driver gtidle/idle_residency_ms:
Wall clock (3 timesteps): 211.4 s
GPU active: 71.9 s (34%)
GPU idle: 139.5 s (66%)

When active: average 2528 MHz (90% of 2800 spec) — GPU compute is fine.
Idle time = MPI Allreduce wait + forced host-buffer copies
(forceHostBuffer=true is mandatory — no GPU-aware MPI for xe).

forceHostBuffer=false crashes immediately in
gko::mpi::all_to_all_v with device pointer.

Discrepancy question: Tsai 2023 paper vs Ginkgo 1.10

Tsai et al. 2023 (Wiley CCPE) explicitly states:

"We have ported ParILU, ParIC, ParILUT, ParICT, ISAI to the DPC++ backend"

Ginkgo 1.10 release notes list new IC/ILU for:
"CUDA, HIP, OpenMP, and Reference backends" — SYCL is absent.

Our measurements confirm: IC → NotImplemented on SYCL in Ginkgo 1.10.

Question: Is the SYCL IC/ILU work from Tsai 2023 in a separate
branch that did not make it into 1.10? What is the roadmap for SYCL
IC/ILU in Ginkgo 2.0?

Summary: All options exhausted in Ginkgo 1.10 SYCL

Preconditioner	Result
BJ maxBlockSize=1	✅ Works — always hits 200-iter cap, never converges
BJ maxBlockSize=2,4,8,16	❌ gko::AllocationError (JACOBI_OPT does not help)
ISAI sparsityPower=1, scaling=-1.0	❌ Diverges
ISAI sparsityPower=3, scaling=-1.0	❌ int32 overflow
IC, scaling=-1.0	❌ NotImplemented
ICT	❌ DEVICE_LOST
Multigrid	❌ OOM + diverges

Full documentation with stack traces, profiling data, and benchmarks:
https://github.com/heikogleu-dev/Openfoam13---GPU-Offloading-Intel-B70-Pro

[greole]

Regarding the OOM error for BJ with maxBlockSize>1 and Multigrid, I think this is likely due to the limited memory available. A typical value to estimate the memory footprint for icoFoam is about 500bytes per cell, which gives you approx 13.6GB in total. However, on the GPU side the memory footprint strongly depends on several factors, 1. Are you offloading all transport equations, i.e. pressure and momentum or just pressure 2. which preconditioners do you use? Additionally, we keep all the matrices (for U and p) on the GPU to avoid repeated copies to the device.
As a rough estimate:
A single linear system on the GPU takes approximately (Ncells * NNZperRow * (TScalar + 2TLabelIdx + TReorder) + NCells*2*TScalar) where Tscalar is the size of a single scalar ie 8byte for DP and TLabel is the size of a Label, typically 4byte, TReorder is the stored reorder map also 4byte. This puts the memory footprint of a single (unpreconditioned) linear system at approximately 5.3GB; a preconditioner with the same storage requirement (scalar Jacobi requires less, multigrid a lot more) as the system matrix increases the footprint to 7.2GB. If you offload pressure and momentum, you end up with 28.8GB.

[heikogleu-dev]

Thank you — this is really helpful context on the memory footprint.

Your calculation makes the Multigrid OOM immediately clear. A question
on a key assumption: in our setup, only p uses the SYCL executor.
U, k, and ω are solved on CPU (PBiCGStab + DILU). Does OGL still
transfer and keep all field matrices (U, k, ω) on the GPU device,
even when only p has executor sycl?

If yes, that would put us at approximately:
p (5.3 GB) + U_x/y/z (3 × 5.3 = 15.9 GB) + k (5.3 GB) + ω (5.3 GB)
≈ 31.8 GB total → exceeds 27.9 GB usable VRAM on B70 Pro

This would explain both Multigrid OOM and possibly BJ maxBlockSize>1
OOM as genuine memory pressure rather than SYCL allocation bugs.

However, with only the p-matrix kept on device (5.3 GB), BJ(2) would
only add ~0.5 GB workspace — well within limits. So the OOM at BS=2
is puzzling unless all matrices are indeed kept on device.

Two follow-up questions:

1. Is there a way to tell OGL to only keep the p-matrix on GPU
   device memory — since we are only solving p on GPU?
   This would reduce the footprint from ~32 GB to ~6 GB for our
   case and potentially make Multigrid and BJ(maxBlockSize>1) viable.

2. If Multigrid requires "a lot more" than the system matrix —
   can you give a rough estimate? At 27.9 GB usable VRAM we're
   curious whether GPU Multigrid would ever fit on a 32 GB card
   at this problem size, even with only the p-field on device.

Hardware: B70 Pro (BMG-G31), 32 GB GDDR6 ECC, 27.9 GB usable VRAM.

[heikogleu-dev]

Thank you for the memory analysis — we instrumented the run with direct
VRAM measurement via /sys/kernel/debug/dri/0/tile0/vram_mm (xe driver
debugfs) sampled at 5 Hz to test the hypothesis empirically.

Results (B70 Pro, 27.9 GB usable VRAM)

Configuration	Peak VRAM	Δ over idle	Outcome
BJ maxBlockSize=1	9.38 GB	8.19 GB	✅ runs (200-iter cap)
BJ maxBlockSize=2	8.41 GB	7.22 GB	❌ SIGABRT

So BJ(2) crashes with less VRAM than BJ(1), and 19.9 GB headroom
is unused at the time of crash. This confirms your "tight coupling"
hypothesis at the level of which fields are on GPU — only p is, not
all 5 fields (which would be ~30 GB and exceed available VRAM). The 8.2
GB Δ for BJ(1) matches your 5.3 GB matrix estimate plus Krylov vectors
and Ginkgo workspace.

But the BJ(2) crash is NOT memory pressure

Looking into the crash more carefully, the actual error message is:
terminate called after throwing an instance of 'gko::AllocationError'
what(): /opt/ogl-src/build/_deps/ginkgo-src/dpcpp/base/executor.dp.cpp:104:
DPC++: failed to allocate memory block of 18446744073709551615B

18446744073709551615 = 2^64 − 1 = (size_t)-1. Stack frame above:
gko::kernels::dpcpp::jacobi::find_blocks<double, int>
→ gko::array::array(executor, count = ~0ULL)
→ DpcppExecutor::raw_alloc(18446744073709551615)

So find_blocks is computing a negative block count and casting to
size_t before the array allocation. This is a deterministic integer
arithmetic bug in dpcpp/preconditioner/jacobi/ for any
maxBlockSize > 1 on this mesh size, not an OOM.

We rebuilt with GINKGO_JACOBI_FULL_OPTIMIZATIONS=ON to verify — same
crash, same underflow value. That flag affects the CUDA adaptive-precision
path, not the SYCL block-counting kernel where the underflow originates.

Implications

This means the available 19.9 GB VRAM headroom is real — there is room
for stronger preconditioners (Multigrid PGM, ParIc factorization, larger
BJ blocks) if the SYCL implementation bugs are addressed. Hardware is
not the constraint here, despite the earlier OOM appearance being
misleading.

Full VRAM traces, logs, and analysis preserved at:
https://github.com/heikogleu-dev/Openfoam13---GPU-Offloading-Intel-B70-Pro/blob/main/profiling/vram_analysis.md

Updated finding with correct root cause:
https://github.com/heikogleu-dev/Openfoam13---GPU-Offloading-Intel-B70-Pro/blob/main/findings/02_bj_blocksize_int_underflow.md

Happy to provide additional measurements or test other configurations
on B70 Pro if useful for tracking down the find_blocks underflow.