Description
Describe the issue:
I expect numpy.array_api.linalg.cholesky(a, upper=True) to return the same array as scipy.linalg.cholesky(a, lower=False), the upper Cholesky decomposition of a. Instead, it returns the conjugate.
There is a possibility that this is a bug on scipy:main, but I assume that it is here instead. I can't check with np.linalg.cholesky since that only provides lower Cholesky decompositions.
xp.linalg.cholesky: https://data-apis.org/array-api/latest/extensions/generated/array_api.linalg.cholesky.html#array_api.linalg.cholesky
scipy.linalg.cholesky: https://docs.scipy.org/doc/scipy/reference/generated/scipy.linalg.cholesky.html#scipy.linalg.cholesky
Maybe a conjugation should be included on this line?:
numpy/numpy/array_api/linalg.py
Line 56 in 2ce972e
Reproduce the code example:
In [1]: import numpy as np
In [2]: import numpy.array_api as npx
<ipython-input-2-f2c2051749b9>:1: UserWarning: The numpy.array_api submodule is still experimental. See NEP 47.
import numpy.array_api as npx
In [3]: from scipy.linalg import cholesky
In [4]: a = np.asarray([[ 10. +0.j, 13. +9.j, 15. -5.j],
...: [ 13. -9.j, 33. +0.j, 33.-10.j],
...: [ 15. +5.j, 33.+10.j, 143. +0.j]])
In [5]: cholesky(a, lower=False)
Out[5]:
array([[3.16227766+0.j , 4.11096096+2.84604989j,
4.74341649-1.58113883j],
[0. +0.j , 2.82842712+0.j ,
6.36396103+3.53553391j],
[0. +0.j , 0. +0.j ,
8.06225775+0.j ]])
In [6]: b = npx.asarray([[ 10. +0.j, 13. +9.j, 15. -5.j],
...: [ 13. -9.j, 33. +0.j, 33.-10.j],
...: [ 15. +5.j, 33.+10.j, 143. +0.j]])
In [7]: npx.linalg.cholesky(b, upper=True)
Out[7]:
Array([[3.16227766+0.j , 4.11096096-2.84604989j,
4.74341649+1.58113883j],
[0. +0.j , 2.82842712+0.j ,
6.36396103-3.53553391j],
[0. +0.j , 0. +0.j ,
8.06225775+0.j ]], dtype=complex128)Error message:
No response
Runtime information:
import sys, numpy; print(numpy.__version__); print(sys.version):
1.25.0
3.10.12 | packaged by conda-forge | (main, Jun 23 2023, 22:40:32) [GCC 12.3.0]
print(numpy.show_runtime()):
[{'numpy_version': '1.25.0',
'python': '3.10.12 | packaged by conda-forge | (main, Jun 23 2023, 22:40:32) '
'[GCC 12.3.0]',
'uname': uname_result(system='Linux', node='lucas-pc', release='6.2.0-26-generic', version='#26~22.04.1-Ubuntu SMP PREEMPT_DYNAMIC Thu Jul 13 16:27:29 UTC 2', machine='x86_64')},
{'simd_extensions': {'baseline': ['SSE', 'SSE2', 'SSE3'],
'found': ['SSSE3',
'SSE41',
'POPCNT',
'SSE42',
'AVX',
'F16C',
'FMA3',
'AVX2'],
'not_found': ['AVX512F',
'AVX512CD',
'AVX512_SKX',
'AVX512_CLX',
'AVX512_CNL',
'AVX512_ICL',
'AVX512_SPR']}},
{'filepath': '/home/lucas/mambaforge/envs/scipy-dev/lib/libmkl_rt.so.2',
'internal_api': 'mkl',
'num_threads': 6,
'prefix': 'libmkl_rt',
'threading_layer': 'intel',
'user_api': 'blas',
'version': '2022.2-Product'},
{'filepath': '/home/lucas/mambaforge/envs/scipy-dev/lib/libomp.so',
'internal_api': 'openmp',
'num_threads': 12,
'prefix': 'libomp',
'user_api': 'openmp',
'version': None},
{'architecture': 'Zen',
'filepath': '/home/lucas/mambaforge/envs/scipy-dev/lib/libopenblasp-r0.3.23.so',
'internal_api': 'openblas',
'num_threads': 12,
'prefix': 'libopenblas',
'threading_layer': 'pthreads',
'user_api': 'blas',
'version': '0.3.23'}]
None
Context for the issue:
I am working on array API support in scipy.linalg. This is failing tests over there for numpy.array_api arrays.