DirectSolver yields incorrect result for random and close to random matrices
未关闭
还没有人认领这个 Issue。
- 主要语言
- Cython
- 星标
- 601
- 派生
- 46
- PR 合并指标
- 30 天内没有已合并 PR
描述
I was considering swapping cupy sparse solvers for nvmath Python/cuDSS DirectSolver, but I am unable to obtain correct results even for small problem instances.
Minimal example
In example01_cupy.py replace
n = 8
...
a += sp.diags([2.0] * n, format="csr", dtype="float64")
with
n = 100
...
# a += sp.diags([2.0] * n, format="csr", dtype="float64")
The same remains true, if I used small, but not tiny values for the diagonal like 0.01.
Example observation 1
||A||: 40.65133430135082
||b||: 14.142135623730951
det(A): 7.262125991606542e+28
02-10 17:35:07 userlogger INFO = SPECIFICATION PHASE =
02-10 17:35:07 userlogger INFO The LHS package is cupyx.
02-10 17:35:07 userlogger INFO The RHS package is cupy.
02-10 17:35:07 userlogger INFO The device_id=0, dtype = float64, index type = int32.
02-10 17:35:07 userlogger INFO The number of equations = 100.
02-10 17:35:07 userlogger INFO The operands' memory space is cuda, and the execution space is on device 0.
02-10 17:35:07 userlogger INFO The specified stream for the DirectSolver ctor is <cuda.core.experimental._stream.Stream object at 0x14b6709b87c0>.
02-10 17:35:07 userlogger INFO The library handle has been created: 94793135484384.
02-10 17:35:07 userlogger INFO The sparse direct solver operation has been created.
02-10 17:35:07 userlogger INFO Starting solver phase ANALYSIS...
02-10 17:35:07 userlogger INFO This call is non-blocking and will return immediately after the operation is launched on the device.
02-10 17:35:07 userlogger INFO Starting solver phase FACTORIZATION...
02-10 17:35:07 userlogger INFO This call is non-blocking and will return immediately after the operation is launched on the device.
02-10 17:35:07 userlogger INFO Starting solver phase SOLVE...
02-10 17:35:07 userlogger INFO This call is non-blocking and will return immediately after the operation is launched on the device.
02-10 17:35:07 userlogger INFO The DirectSolver object's resources have been released.
CuDSS:
||x||: 29391823863.993862
||Ax - b||: 92575181911.9934
CuPy:
||x||: 10.647411227999223
||Ax - b||: 4.289313285599112e-14
NumPy:
||x||: 10.647411227999223
||Ax - b||: 4.29904389535207e-14
Note that I computed the determinant to check whether A is invertible.
Example observation 2
Vice versa I made the observation that DirectSolver sometimes returns not nan even though the determinant is 0. Note that here b = cp.ones((n, 1), order="F") was only a vector
||A||: 5.6922584694712395
||b||: 10.0
det(A): 0.0
02-10 17:30:48 userlogger INFO = SPECIFICATION PHASE =
02-10 17:30:48 userlogger INFO The LHS package is cupyx.
02-10 17:30:48 userlogger INFO The RHS package is cupy.
02-10 17:30:48 userlogger INFO The device_id=0, dtype = float64, index type = int32.
02-10 17:30:48 userlogger INFO The number of equations = 100.
02-10 17:30:48 userlogger INFO The operands' memory space is cuda, and the execution space is on device 0.
02-10 17:30:48 userlogger INFO The specified stream for the DirectSolver ctor is <cuda.core.experimental._stream.Stream object at 0x14c6b6504580>.
02-10 17:30:48 userlogger INFO The library handle has been created: 94064601263200.
02-10 17:30:48 userlogger INFO The sparse direct solver operation has been created.
02-10 17:30:48 userlogger INFO Starting solver phase ANALYSIS...
02-10 17:30:48 userlogger INFO This call is non-blocking and will return immediately after the operation is launched on the device.
02-10 17:30:48 userlogger INFO Starting solver phase FACTORIZATION...
02-10 17:30:48 userlogger INFO This call is non-blocking and will return immediately after the operation is launched on the device.
02-10 17:30:48 userlogger INFO Starting solver phase SOLVE...
02-10 17:30:48 userlogger INFO This call is non-blocking and will return immediately after the operation is launched on the device.
02-10 17:30:48 userlogger INFO The DirectSolver object's resources have been released.
CuDSS:
||x||: 1.5379474527090538e+102
||Ax - b||: 1.5379474527090537e+89
CuPy:
||x||: nan
||Ax - b||: nan
NumPy:
||x||: nan
||Ax - b||: nan
贡献指南
从这里开始
- 先读完整个 Issue,再读项目的贡献指南。
- 在 Issue 下留言说明你要接手 —— 这能避免两个人做同样的事。
- Fork 仓库,在一个分支上完成修改。
- 提交 Pull Request,并在描述里引用这个 Issue 编号。
调研方向
从 example01_cupy.py 开始,通过改变 n 和对角线值来复现报告中的案例。将 DirectSolver 的结果与 CuPy 和 NumPy 的输出进行比较,包括 issue 中显示的范数和残差。完成的标准是确定并修正非奇异矩阵和奇异矩阵产生错误结果的原因。
由索引模型根据 Issue 内容生成。
评估
- 技术栈
- python
- 领域
- backend
- Issue 类型
- 缺陷
- 难度
- 4/5
- 预计耗时
- 3-5 天
- 活跃度
- 停滞
- 描述清晰度
- 基本清楚
- 新手友好度
- 35/100