I am making a function for a matrix class for finding the inverse matrix using gaussian elimination, The function works fine most of the time it would have random warnings that come in pairs that look like this

RuntimeWarning: invalid value encountered in double_scalars
RuntimeWarning: divide by zero encountered in double_scalars

The problem comes for the lines where I divide everything in the same row by a constant. However, the first part of the function already makes sure the all the pivot elements in the matrix would not be zero. I have also ran a lot of tests and the problem isn't coming from the pivot elements being zero, the pivot elements are the only elements used for division, so I don't understand where the error is coming from

def inverse(self):
        if self.row != self.col:
            raise ValueError("Inverse Matrices only possible with square matrices")
        det = self.determinant()
        if det == 0:
            raise ValueError("The provided matrix does not have an Inverse")

        # Makes sure that the pivot elements are not zero
        inverse = Matrix.identity(self.row)
        placeHolder = [[self.matrix[i, j] for j in range(self.col)] for i in range(self.row)]
        for i in range(self.col):
            if self.matrix[i, i] == 0:
                for j in range(self.row):
                    if self.matrix[j, i] != 0:
                        for k in range(self.col):
                            self.matrix[i, k]  = self.matrix[j, k]
                            inverse.matrix[i, k]  = inverse.matrix[j, k]
                        break

        # Row operations
        for i in range(self.row):
            constant = self.matrix[i, i]
            for j in range(self.col):
                self.matrix[i, j] /= constant
                inverse.matrix[i, j] /= constant

            row = i
            while row != 0:
                constant = self.matrix[row-1, i]
                for j in range(self.col):
                    self.matrix[row-1, j] -= self.matrix[i, j] * constant
                    inverse.matrix[row-1, j] -= inverse.matrix[i, j] * constant
                row -= 1

            row = i
            while row != self.row - 1:
                constant = self.matrix[row 1, i]
                for j in range(self.col):
                    self.matrix[row 1, j] -= self.matrix[i, j] * constant
                    inverse.matrix[row 1, j] -= inverse.matrix[i, j] * constant
                row  = 1

        self.set(placeHolder)
        return inverse

CodePudding user response：

While initially you check that the diagonal elements of a matrix are non-zero, these elements are modified as you perform row operations, and they may become zeros. A usual step in the row reduction of a matrix is to check that an element in a pivot position is non-zero before using it do modify other rows, since this cannot be determined beforehand.

As a side note, I don’t know how the function computing the determinant is implemented, but because floating point arithmetic is not exact, typically computations of a determinant can give a non-zero value even if matrix is singular. The same problems applies to row reduction of matrices: due to rounding errors these computations may give wrong results, with singular matrices appearing to be non-singular and vice versa. In numerical analysis there are a lot of considerations and methods that deal with such issues.