Метод Гаусса для вычисления определителя

int gauss_method(int razm, double **matrix, double **&matrix_obr)
{
    int i, j, k;
    //создание единичной матрицы 
    for (i = 0; i<razm; i++)
    {
        for (j = 0; j<razm; j++)matrix_obr[i][j] = 0;
        matrix_obr[i][i] = 1;
    }
    //приводим к верхнетреугольному виду
    double a, b, det=1;
    for (i = 0; i<razm; i++)
    {
        a = matrix[i][i];
        for (j = i + 1; j<razm; j++)
        {
            b = matrix[j][i];
            for (k = 0; k<razm; k++)
            {
                matrix[j][k] = matrix[i][k] * b - matrix[j][k] * a;
                matrix_obr[j][k] = matrix_obr[i][k] * b - matrix_obr[j][k] * a;
            }
        }
    }
    for (i = 0; i < razm; i++)
        det *= matrix[i][i];
 
    //print_matrix(razm, razm, matrix);
    cout << "Определитель = " << det << endl;
    
    //обратный ход нахождения элементов обратной матрицы
    double sum;
    for (i = 0; i<razm; i++)
    {
        for (j = razm - 1; j >= 0; j--)
        {
            sum = 0;
            for (k = razm - 1; k>j; k--)
                sum += matrix[j][k] * matrix_obr[k][i];
            matrix_obr[j][i] = (matrix_obr[j][i] - sum) / matrix[j][j];
        }
    }
    return 1;
}

matrix[j][k] -= matrix[i][k] * b / a;