final var threshold = EPSILON * Math.max(rows, cols) * s[0];
            final var invW = new Matrix(cols, cols);
            invW.initialize(0.0);
            double value;
            for (var n = 0; n < cols; n++) {
                value = s[n];
                if (value > threshold) {
                    invW.setElementAt(n, n, 1.0 / value);
                }
            }
            // V * invW * U', where U' is transposed of U
            u.transpose();
            return v.multiplyAndReturnNew(invW.multiplyAndReturnNew(u));
        } catch (final DecomposerException e) {
            throw e;
        } catch (final Exception e) {
            throw new DecomposerException(e);
        }
    }

    /**
     * Computes Moore-Penrose pseudo-inverse of provided array considering it as
     * a column matrix.
     *
     * @param array Input array
     * @return More-Penrose pseudo-inverse
     * @throws DecomposerException Exception thrown if matrix cannot be
     *                             inverted, usually because matrix contains numerically unstable values
     *                             such as NaN or inf.
     * @see #pseudoInverse(Matrix)
     */
    public static Matrix pseudoInverse(final double[] array) throws DecomposerException {

}

    /**
     * Solves a linear system of equations of the following form:
     * A * X = B, where A is the input matrix provided for QR decomposition,
     * X is the solution to the system of equations, and B is the parameters
     * vector/matrix.
     * Note: This method can be reused for different b vector/matrices without
     * having to recompute QR decomposition on the same input matrix.
     * Note: Provided b matrix must have the same number of rows as provided
     * input matrix A, otherwise a WrongSizeException will be raised.
     * Note: Provided input matrix "A" must have at least as many rows as columns,
     * otherwise a WrongSizeException will be raised as well. For input matrices
     * having a higher number of rows than columns, the system of equations will
     * be overdetermined and the least squares solution will be found.
     * Note: If provided input matrix A is rank deficient, a
     * RankDeficientMatrixException will be thrown.
     * Note: In order to execute this method, a QR decomposition must be
     * available, otherwise a NotAvailableException will be raised. In order to
     * avoid this exception call decompose() method first.
     *
     * @param b      Parameters matrix that determine a linear system of equations.
     *               Provided matrix must have the same number of rows as provided input
     *               matrix for QR decomposition. Besides, each column on parameters matrix
     *               will represent a new system of equations, whose solution will be returned
     *               on appropriate column as an output of this method.
     * @param result Matrix containing solution of linear system of equations on
     *               each column for each column of provided matrix of parameters b. Provided
     *               matrix will be resized if needed
     * @throws NotAvailableException        Exception thrown if attempting to call this
     *                                      method before computing QR decomposition. To avoid this exception call
     *                                      decompose() method first.
     * @throws WrongSizeException           Exception thrown if attempting to call this
     *                                      method using an input matrix with fewer rows than columns; or if provided
     *                                      parameters matrix (b) does not have the same number of rows as input
     *                                      matrix being QR decomposed.
     * @throws RankDeficientMatrixException Exception thrown if provided input
     *                                      matrix to be QR decomposed is rank deficient. In this case linear system
     *                                      of equations cannot be solved.
     */
    public void solve(final Matrix b, final Matrix result) throws NotAvailableException, WrongSizeException,
            RankDeficientMatrixException {
        solve(b, DEFAULT_ROUND_ERROR, result);
    }

    /**
     * Solves a linear system of equations of the following form:
     * A * X = B, where A is the input matrix provided for QR decomposition,
     * X is the solution to the system of equations, and B is the parameters
     * vector/matrix.
     * Note: This method can be reused for different b vector/matrices without
     * having to recompute QR decomposition on the same input matrix.
     * Note: Provided b matrix must have the same number of rows as provided
     * input matrix A, otherwise a WrongSizeException will be raised.
     * Note: Provided input matrix "A" must have at least as many rows as columns,
     * otherwise a WrongSizeException will be raised as well. For input matrices
     * having a higher number of rows than columns, the system of equations will
     * be overdetermined and the least squares solution will be found.
     * Note: If provided input matrix A is rank deficient, a
     * RankDeficientMatrixException will be thrown.
     * Note: In order to execute this method, a QR decomposition must be
     * available, otherwise a NotAvailableException will be raised. In order to
     * avoid this exception call decompose() method first.
     *
     * @param b             Parameters matrix that determine a linear system of equations.
     *                      Provided matrix must have the same number of rows as provided input
     *                      matrix for QR decomposition. Besides, each column on parameters matrix
     *                      will represent a new system of equations, whose solution will be returned
     *                      on appropriate column as an output of this method.
     * @param roundingError threshold to determine whether matrix b has full rank or not.
     *                      By default, this is typically a tiny value close to zero.
     * @param result        Matrix containing solution of linear system of equations on
     *                      each column for each column of provided matrix of parameters b. Provided
     *                      matrix will be resized if needed
     * @throws NotAvailableException        Exception thrown if attempting to call this
     *                                      method before computing QR decomposition. To avoid this exception call
     *                                      decompose() method first.
     * @throws WrongSizeException           Exception thrown if attempting to call this
     *                                      method using an input matrix with fewer rows than columns; or if provided
     *                                      parameters matrix (b) does not have the same number of rows as input
     *                                      matrix being QR decomposed.
     * @throws RankDeficientMatrixException Exception thrown if provided input
     *                                      matrix to be QR decomposed is rank deficient. In this case linear system
     *                                      of equations cannot be solved.
     * @throws IllegalArgumentException     if provided rounding error is negative.
     */
    public void solve(final Matrix b, final double roundingError, final Matrix result)
            throws NotAvailableException, WrongSizeException, RankDeficientMatrixException {

        if (!isDecompositionAvailable()) {
            throw new NotAvailableException();
        }

        final var rows = inputMatrix.getRows();
        final var columns = inputMatrix.getColumns();
        final var rowsB = b.getRows();
        final var colsB = b.getColumns();
        double s;