torch.inverse based on cuSOLVER does not raise error for singular input

[IvanYashchuk]

🐛 Bug

The cuSOLVER path in torch.inverse seems not to raise an error for non-invertible input.

To Reproduce

Single matrix and 2xNxN batched input inversion is based on cuSOLVER

import torch
batch_dim = 2
A = torch.eye(3, 3, dtype=dtype, device=device)
A = A.reshape((1, 3, 3))
A = A.repeat(batch_dim, 1, 1)
A[0, -1, -1] = 0  # Now A[0] is singular
Ainv = torch.inverse(A) # Doesn't raise errors Ainv[0] contains NaNs

General batched input case is based on MAGMA

import torch
batch_dim = 3
A = torch.eye(3, 3, dtype=dtype, device=device)
A = A.reshape((1, 3, 3))
A = A.repeat(batch_dim, 1, 1)
A[0, -1, -1] = 0  # Now A[0] is singular
Ainv = torch.inverse(A) # Raises RuntimeError because input is singular

Expected behavior

torch.inverse should raise an error for singular input when input is single matrix or a batch of matrices with batch size less or equal than 2 (cuSOLVER's code path).

cc: @xwang233

cc @ngimel @vishwakftw @jianyuh @nikitaved @pearu @mruberry @heitorschueroff

[ezyang]

Might not be so easy to fix; singular test has to be done on CUDA

[ngimel]

@IvanYashchuk is it because info is not checked for cuSolver (that would be an easy fix then), or because cuSolver just returns bogus results?

[IvanYashchuk]

Taking a quick look at the code, info is being checked here.
cuSOLVER correctly returns info = 3 for this input (checked via CuPy).

[xwang233]

pytorch/aten/src/ATen/native/cuda/BatchLinearAlgebraLib.cu

Lines 28 to 35 in d75c402

	template <typename scalar_t>
	inline static void _apply_single_inverse_helper(scalar_t* self_ptr, scalar_t* self_inv_ptr, int* ipiv_ptr, int* info_ptr, int n) {
	// self_inv_ptr should already be an identity matrix
	auto handle = at::cuda::getCurrentCUDASolverDnHandle();
	at::cuda::solver::getrf<scalar_t>(handle, n, n, self_ptr, n, ipiv_ptr, info_ptr);
	at::cuda::solver::getrs<scalar_t>(handle, n, n, self_ptr, n, ipiv_ptr, self_inv_ptr, n, info_ptr);
	}

I found the problem. In your example, getrf returns info = 3, but it was later overwritten by getrs with info = 0. So we definitely need to check info returned by the first call as well.

I was thinking if we can do this asynchronously to avoid performance drop due to the device -> host memory copy.

[ngimel]

You can use 2-element infos to save both.