I have the following function in python:
def cramer(coef_matrix, free_coef_matrix):
    num_rows, num_cols = coef_matrix.shape
    if num_rows != num_cols:
        raise Exception("Cannot apply cramer's rule to a non-square matrix")
    determinant = np.linalg.det(coef_matrix)
    if determinant == 0:
        raise Exception("Cannot apply cramer's rule to a matrix with determinant 0")
    solution = num_rows*[0]
    for i in range(num_rows):
        coef_matrix_copy = coef_matrix.copy()
        coef_matrix_copy[:, i] = free_coef_matrix
        current_det = np.linalg.det(coef_matrix_copy)
        solution[i] = current_det / determinant
    return solution
But the problem is when I run the function on the input:
coef_matrix:
[2, 1, 1],
[1, -1, -1]
[1, 2, 1]
free_coef_matrix:
[3, 0, 0]
It returns [1.0, -2.0000000000000004, 3.0] instead of [1.0, -2.0, 3.0].
I know that there might be problems related to floating point, but I want to know if there might be a solution?
