[Issue] Device dotc uses unconjugated cuBLAS dotu and breaks complex GMRES orthogonalization
Nobody has claimed this yet.
- Dominant language
- Cuda
- Stars
- 692
- Forks
- 197
- PR merge metrics
- No merged PRs in 30d
Description
Describe the issue
AMGX v2.5.0 appears to use an unconjugated cuBLAS dot product in the device implementation of dotc, which causes incorrect behavior for complex-valued GMRES.
In src/blas.cu, the device implementation of dotc(...) calls:
Cublas::dot(size, a_raw, 1, b_raw, 1, &result);
For complex double precision, Cublas::dot(...) is eventually mapped to cublasZdotu, which computes an unconjugated dot product.
However, the expected semantics of dotc(x, y) for complex-valued vectors are:
$$
x^H y
$$
i.e. the first argument should be conjugated.
AMGX already provides Cublas::dotc(...), which maps to cublasZdotc.
There is also a related issue in the GMRES Arnoldi orthogonalization in src/solvers/gmres_solver.cu. The current code computes:
m_H(k, i) = dot(A,
m_V_vectors[i + 1],
m_V_vectors[k]);
After restoring the expected dotc(x,y) = x^H y semantics, this computes
$$
w^H v_k
$$
while Arnoldi / Modified Gram-Schmidt requires
$$
v_k^H w.
$$
Therefore the two arguments need to be swapped for the complex case.
With the original AMGX v2.5.0 implementation, the minimal 3x3 complex system below does not converge after 500 GMRES iterations.
With the following two changes:
- Cublas::dot(size, a_raw, 1, b_raw, 1, &result);
+ Cublas::dotc(size, a_raw, 1, b_raw, 1, &result);
in the device dotc implementations, and
- m_H(k, i) = dot(A, m_V_vectors[i + 1], m_V_vectors[k]);
+ m_H(k, i) = dot(A, m_V_vectors[k], m_V_vectors[i + 1]);
in GMRES, the same 3x3 system converges in 3 iterations to approximately machine precision.
Observed result
| Original AMGX v2.5.0 | Patched | |
|---|---|---|
| Matrix size | 3x3 complex | 3x3 complex |
| Initial residual | 1.270827e+01 |
1.270827e+01 |
| Total iterations | 500 |
3 |
| Final residual | 2.955383e-05 |
4.981581e-16 |
| Status | Did not converge | Converged |
The original solver terminates with:
Total Iterations: 500
Final Residual: 2.955383e-05
Total Reduction in Residual: 2.325560e-06
***Solver Did Not Converge
The patched solver gives:
Ini 1.270827e+01
0 2.975743e+00
1 1.027928e+00
2 4.981581e-16
Total Iterations: 3
Final Residual: 4.981581e-16
Total Reduction in Residual: 3.919954e-17
Because this is only a 3x3 system and the GMRES restart length is larger than the system dimension, the patched result is also consistent with the expected behavior of full GMRES for this small problem.
Environment information
- OS:
Ubuntu 24.04.4 LTS - CUDA runtime:
CUDA 13.0 - MPI version: N/A (
CMAKE_NO_MPI=ON) - AMGX version:
v2.5.0 - AMGX commit hash:
cc1cebdbb32b14d33762d4ddabcb2e23c1669f47 - NVIDIA driver:
591.86 - CUDA driver API reported by AMGX:
13.1 - NVIDIA GPU:
NVIDIA GeForce RTX 5070 Ti - Build type:
Release - AMGX mode:
dZZI - Relevant environment variables: none required for reproduction other than the library path used to load the corresponding AMGX build.
AMGX startup information:
AMGX version 2.5.0
Built on Aug 25 2026
Compiled with CUDA Runtime 13.0, using CUDA driver 13.1
AMGX solver configuration
The following configuration was used:
{
"config_version": 2,
"determinism_flag": 1,
"exception_handling": 1,
"solver": {
"scope": "main",
"solver": "GMRES",
"preconditioner": {
"scope": "noprec",
"solver": "NOSOLVER"
},
"max_iters": 500,
"gmres_n_restart": 10,
"tolerance": 1e-8,
"norm": "L2",
"convergence": "RELATIVE_INI",
"monitor_residual": 1,
"store_res_history": 1,
"print_solve_stats": 1,
"obtain_timings": 1
}
}
No preconditioner is used so that the reproducer isolates the behavior of GMRES itself.
Matrix Data
The issue can be reproduced using the following 3x3 complex-valued Matrix Market system:
%%MatrixMarket matrix coordinate complex general
%%NVAMG 1 1 rhs solution
3 3 7
1 1 4.0 0.0
1 2 1.0 1.0
2 1 1.0 -1.0
2 2 3.0 0.0
2 3 1.0 2.0
3 2 1.0 -2.0
3 3 5.0 0.0
3
7.0 5.0
6.0 -4.5
-5.0 -2.5
3
0.0 0.0
0.0 0.0
0.0 0.0
The matrix is
$$
A=
\begin{bmatrix}
4 & 1+i & 0\
1-i & 3 & 1+2i\
0 & 1-2i & 5
\end{bmatrix}.
$$
The right-hand side is
$$
b=
\begin{bmatrix}
7+5i\
6-4.5i\
-5-2.5i
\end{bmatrix}.
$$
A complete reproducer containing the matrix, solver configuration, original log, patched log, environment information, and patch is attached to this issue.
Reproduction steps
AMGX was built without MPI.
Example configuration:
mkdir build
cd build
cmake \
-DCMAKE_BUILD_TYPE=Release \
-DCMAKE_NO_MPI=ON \
-DCMAKE_CUDA_ARCHITECTURES=120 \
..
make -j8
Run the reproducer using the device complex-double mode:
LD_LIBRARY_PATH="$PWD/build" \
./build/examples/amgx_capi \
-mode dZZI \
-m A.mtx \
-c gmres_complex.json
Original AMGX v2.5.0
The original implementation runs all 500 iterations and does not satisfy the requested tolerance:
Total Iterations: 500
Final Residual: 2.955383e-05
***Solver Did Not Converge
Patched version
Changing the device dotc calls from:
Cublas::dot(...)
to:
Cublas::dotc(...)
and changing the GMRES projection from:
dot(A, m_V_vectors[i + 1], m_V_vectors[k])
to:
dot(A, m_V_vectors[k], m_V_vectors[i + 1])
produces:
Total Iterations: 3
Final Residual: 4.981581e-16
No other solver changes were made for this comparison.
Additional context
For complex vectors, the distinction between the unconjugated and conjugated dot products is fundamental.
For example, for
$$
x=[1,i]^T,
$$
the unconjugated product gives
$$
x^T x = 1+i^2=0,
$$
while the Hermitian inner product gives
$$
x^H x=|1|^2+|i|^2=2.
$$
Therefore an unconjugated dot product cannot be used as the complex inner product required by Arnoldi orthogonalization.
For real-valued inputs this problem is hidden because
$$
x^T y=y^T x,
$$
whereas for complex-valued inputs
$$
x^H y
$$
and
$$
y^H x
$$
are generally different.
The minimal reproducer shows a significant behavioral difference:
Original v2.5.0:
500 iterations
final residual = 2.955383e-05
did not converge
Patched:
3 iterations
final residual = 4.981581e-16
converged
I can also prepare a pull request with the two changes above and a regression test if this diagnosis is confirmed.
Contributor guide
No contributing guide indexed for this repository
First steps
- Read the whole issue, then the project's contributing guide.
- Comment on the issue to say you are picking it up — it saves two people doing the same work.
- Fork the repository and make your change on a branch.
- Open a pull request that references the issue number.
Research direction
Inspect the device dotc implementation in src/blas.cu and the Arnoldi projection in src/solvers/gmres_solver.cu, then build the provided 3x3 complex reproducer with the documented CMake and amgx_capi commands. Verify that complex GMRES converges in about three iterations with the expected residual, and add regression coverage if the project’s test structure supports it.
Written by the indexing model from the issue text.
Assessment
- Tech stack
- cpp
- Domain
- hpc
- Issue type
- Bug
- Difficulty
- 3/5
- Estimated time
- 1-2 days
- Activity status
- Active
- Clarity
- Clearly specified
- Newbie friendliness
- 72/100