\mathbb{K}^{m \times k}, S \in \mathbb{R}^k, V \in \mathbb{K}^{k
      \times n}}
where \operatorname{diag}(S) \in \mathbb{K}^{k \times k}. In this
   case, U and V also have orthonormal columns.
Supports input of float, double, cfloat and cdouble dtypes. Also
   supports batches of matrices, and if "A" is a batch of matrices
   then the output has the same batch dimensions.
The returned decomposition is a named tuple (U, S, Vh) which
   corresponds to U, S, V^{\text{H}} above.
The singular values are returned in descending order.
The parameter "full_matrices" chooses between the full (default)
   and reduced SVD.
The "driver" kwarg may be used in CUDA with a cuSOLVER backend to
   choose the algorithm used to compute the SVD. The choice of a
   driver is a trade-off between accuracy and speed.


If "A" is well-conditioned (its condition number is not too
     large), or you do not mind some precision loss.

For a general matrix: 'gesvdj' (Jacobi method)


