Метод Камера c++. Матрица 1000х1000. Floating point overflow

[cripple_nipple]

Написал метод Гаусса с выбором главного элемента, который считает определитель 1000х1000 (размерности от 2х2 до 7х7 считает правильно, надеюсь и 1000х1000 тоже, еще не проверял).

Теперь нужно сделать Крамера. За основу взял уже написанный мной метод Гаусса, добавив функцию вычисления определителя с подстановкой i-го столбца. Тут и появляется ошибка:

При этом подсвечивается

for (int i = j; i < n-1; i++)   
        {                        
            double coef = masTmp[i+1][j] / masTmp[j][j];  //<-вот эта строка
            for (int jj = j; jj < n; jj++)
            {
                double tmpM = coef * masTmp[j][jj];
                masTmp[i+1][jj] -= tmpM;
            }
        }

К слову, этот код также решает матрицы до размерности 7х7 включительно (другие не проверял), но 500х500 и 1000х1000 уже не получается.

Пробовал кучу всего, не могу решить.

[cripple_nipple]

Здесь по ссылке: https://vk.com/doc296627823_378907723 текстовые файлы с данными. Сначала нужно загружать файл с "m" в названии - это файл с матрицей, а затем с "v" - это столбец свободных членов. Очевидно, что файлы идут парами, то есть 1000m.txt + 1000v.txt, 7m.txt + 7v.txt и т.д. Вводить названия в консоль нужно с расширением

[challengerr]

Число с плавающей точкой переполняется. Длинную арифметику нужно дополнительно реализовывать, то есть вручную делить.

[cripple_nipple]

Сообщение от challengerr

Число с плавающей точкой переполняется. Длинную арифметику нужно дополнительно реализовывать, то есть вручную делить.

Гуглил уже это. Слишком напряжно. Есть другой способ?
Проблема еще в том, что даже алгоритмы не гуглятся. Я бы попробовал переписать прогу, но не знаю как.

[p51x]

Есть другой способ?

Перейти к целочисленной арифметике. Заюзать сопроцессор или спец. комманды, типа каких-то SSE*, AVX, ... но не уверен, что это будет проще, если вы не разобрались в длинной арифметике, материалов по которой навалом.

[krGorka]

double coef = masTmp[i+1][j] / masTmp[j][j];

скорее всего деление на 0. Обработать это событие заранее.

[cripple_nipple]

Сообщение от krGorka

double coef = masTmp[i+1][j] / masTmp[j][j];

скорее всего деление на 0. Обработать это событие заранее.

Как это событие можно обработать?

[cripple_nipple]

Подкиньте хотя бы алгоритмов нахождения определителя?