pyqcu-cpp-qcu

Files (2)

• AGENTS.md
• SKILL.md

SKILL.md (first lines)

# cpp/cuda/qcu

Primary C++ CUDA backend for PyQCU. Hand-tuned CUDA kernels with MPI halo exchange for Wilson/Clover Dirac operators, BiStabCG/CG solvers, multigrid, and gauge field generation.

Production execution topology is fixed:

```text
one Python process = one MPI rank = one GPU/DCU/GPGPU
```

Cross-device solves use MPI domain decomposition and C++ halo/global-reduction
paths. The Python-level `pyqcu.cuda._multi_gpu.MultiGpuMultigrid` driver is a
test-only, explicitly constructed exception and must not be used as the default
template for new CUDA programs.

Normative project rules:
[`docs/cuda_execution_topology.md`](../../docs/cuda_execution_topology.md).

## Build

```bash
source ./env.sh       # CUDA toolkit paths, MPI, etc.
bash ./make.sh        # symlinks CMakeLists-nv.txt → CMakeLists.txt, then cmake + make
```

Output: `libqcu.so` — dynamically linked library loaded by the Cython bridge.

CMake 默认自动使用 `env.sh` 探测到的 GPU compute capability（cubin + PTX），
不需要为单一 GPU 型号锁死架构；V100 生产基准仍可显式设置：

```bash
QCU_CUDA_ARCHITECTURES='70-real;70-virtual' bash ./build.sh
cuobjdump --list-elf cpp/cuda/qcu/libqcu.so | head
```

`cuobjdump` 必须显示 `sm_70.cubin`，而不是仅显示 `sm_60.cubin`。不要把
sm_60 JIT 与 sm_70 native 的 solve 数字混在同一速度比中；至少记录
`QCU_CUDA_ARCHITECTURES`、libqcu SHA256 和实际 `cuobjdump` 架构。当前
V100 side-only 复测中 sm70 大 c64 为 `0.657074 s`、大 c128 为
`2.062569 s`；随后 full formal 独立复测大 c64 为 `0.655570 s`。
sm70 相对默认 sm60 JIT 改善有限，但必须保留架构 provenance。
benchmark worker 已自动写入 `provenance.pyqcu_library`；不要依赖口头
记录或只看 `LD_LIBRARY_PATH`。

## Source Organization

```
include/          — 24 header files (templated C++ with CUDA kernels inline)
├── define.h      — Parameter index constants (must mirror pyqcu/cuda/define.py)
├── include.h     — 系统/第三方头聚合（cuda/cublas/curand/mpi/chrono/...）
├── lattice_complex.h   — Complex number arithmetic (operator*= was fixed for overwrite bug)
├── lattice_set.h       — Lattice geometry, grid layout, site indexing
├── lattice_cuda.h      — CUDA utility functions (stream management, etc.)
├── lattice_mpi.h       — MPI halo exchange helpers (blocking Sendrecv)
├── qcu.h               — Top-level include aggregator
├── dslash.h            — Dslash dispatch
├── wilson_dslash.h     — Wilson dslash kernel
├── clover_dslash.h     — Clover dslash entry
├── lattice_wilson_dslash.h   — Wilson dslash implementation
├── lattice_clover_dslash.h   — Clover dslash implementation
├── bistabcg.h          — BiCGStab algorithm (GPU kernels)
├── cg.h                — Conjugate gradient algorithm
├── lattice_wilson_bistabcg.h — Wilson BiStabCG solver
├── lattice_wilson_cg.h        — Wilson CG solver
├── lattice_clover_bistabcg.h  — Clover BiStabCG solver
├── lattice_sap.h              — SAP（Schwarz 交替过程）mask/update 内核，5-stream MG 的平滑器
├── multigrid.h                — MG restrict/prolong/coarse-dslash
├── lattice_multigrid.h        — MG implementation
├── lattice_clover_multigrid.h — Clover multigrid solver（5745 行，5-stream 架构）
├── laplacian.h / lattice_laplacian.h — Laplacian operator
└── gauss_gauge.h              — Gaussian gauge field generation

src/              — 32 个 .cu 文件：#include 对应头并做模板实例化与启动封装
python/
└── pyqcu.h       — C API declarations (extern "C"，37 个入口); must match pyqcu/cuda/qcu/qcu_api.pxd
```

## Parameter Protocol

Parameters are passed from Python as flat arrays. Index constants in `include/define.h` must stay in sync with `pyqcu/cuda/define.py`:

- **`params`** (CPU int32[58]): lattice dims, grid sizes, data types, iteration counts, plan selection, MG level configs and controls. Slots 54–57 are `_MG_USE_GCR_`, `_MG_USE_DEFLATE_`, `_MG_MU_PRE_`, and `_MG_USE_INIT_GUESS_`; slot 57 is the 0/1 warm-start flag and the ABI remains length 58.
- **`argv`** (real[8]): mass, atol, sigma, MG tolerances, csw；c64 为 float32，
  c128 为 float64。csw 位于槽 7，Clover 系数为 `-kappa*csw/8`。
- **`set_ptrs`** (CPU int64[100]): LatticeSet/scratch pointers. Per strict transition `t`, slots `60+4*t+0..3` are blocked `V`, optional raw `Y`, `Yhat`, and onsite pair `(X,X^-1)`; slot `80` is the persistent hierarchy handle.

`_SET_PLAN_` (params[16]) selects the kernel plan:
- `-2` = Laplacian, `-1` = Gauss gauge, `0` = Wilson dslash, `1` = BiStabCG/CG, `2` = Clover dslash

## Strict QUDA-style MultiGrid

The strict path is parallel to the retained legacy MG implementation. Select it with `hierarchy_mode="strict"` or `StrictMultigrid`; `setup_operator="schur"` belongs to the legacy compact odd-Schur path and is not a synonym for Strict. The currently supported all-level mode is `coarse_grid_solution_type=matpc` with `smoother_solve_type=direct_pc`. At level `l`, split the full operator as

\[
D_l=X_l+H_l,\qquad \widehat D_l=X_l^{-1}D_l,\qquad
D_{l+1}=R_l\widehat D_lP_l.
\]

The full/raw coarse onsite block and directional links are dense
\(E\times E\) matrices, not \(SU(3)\) links:
\[
X_C(X)=\sum_{x\in\mathcal A_X}V_X^\dagger C V_X,
\]
\[
Y^{\rm f}_{\mu}(X)=-\kappa\sum_{x\in\mathcal A_X}
(AV_X)^\dagger U_\mu V_{X+\hat\mu},
\]
\[
Y^{\rm b}_{\mu}(X-\hat\mu)=-\kappa\sum_{x\in\mathcal A_X}
V_X^\dagger U^\dagger_\mu AV_{X-\hat\mu},
\]
where \(AV=C^{-1}V\) for Clover-PC and \(AV=V\) for the raw path.
For a Wilson-like operator the onsite contribution also includes the
identity \(I_E\); a coarse-to-coarse construction projects the current
\(X_\ell\), not the fine Clover tensor. In the current Strict left-
preconditioned path, \(AV\) is constructed first and the coarse local
operator is based on \(X_f^{-1}D_f\); therefore \(V^\dagger CV\) describes
the full/raw Clover construction and must not be stated as an unconditional
formula for the Strict PC path.
The raw coarse action is

From ZhangXin8069/PyQCU, MIT.

More from ZhangXin8069/PyQCU