Целочисленные вычисления
А что это я всё о политике, да о политике... Можно поговорить и о куда более приятных вещах.
Недавно моей жене понадобилось ну чисто для домашнего хозяйства посчитать определители некоего семейства целочисленных матриц. Я, стремясь быть паинькой, не торопясь запрограмировал это дело. Определитель я решил считать примерно как в школе учили, приводя матрицу к верхнетреугольной, но, чтобы ограничиться целочисленными вычислениями и не лезть во всякие double'ы, занулять элементы матрицы решил не вычисляя нецелый коэффициент, а, по алгоритму Евклида попытаться занулить либо текущий диагональный элемент, либо поддиагональный, уж как получится, если занулится диагональный, то переставить две строчки. Для тех, кто не понял -- псевдокод.
Процедура Занулить поддиагональную часть столбца i
Для всех j > i
Дело нехитрое, написал эти 20 строчек, отладил на небольших матрицах запустил... Вот блин. Переполнение int'а. Невлезаем в 32 бита на матрице 8 на 8.
Удивился. Почесал репу. Заменил 32-битное целое 64-битным. Запустил -- переполнился на матрице 11 на 11.
Ё-моё. В исходной матрице все элементы по модулю не больше 4, в основном она из нулей и единиц состоит. Какого хрена????
Подумал. Применил небольшую магию -- перед занулением поддиагональных элементов переставляю две строки матрицы так, чтобы на диагональ встал самый маленький по модулю элемент. Полгечало. Переполнения нет, всё посчиталось. Но вопрос остался.
Какого хрена оно переполняется? Откуда такой бешеный рост значений???
Недавно моей жене понадобилось ну чисто для домашнего хозяйства посчитать определители некоего семейства целочисленных матриц. Я, стремясь быть паинькой, не торопясь запрограмировал это дело. Определитель я решил считать примерно как в школе учили, приводя матрицу к верхнетреугольной, но, чтобы ограничиться целочисленными вычислениями и не лезть во всякие double'ы, занулять элементы матрицы решил не вычисляя нецелый коэффициент, а, по алгоритму Евклида попытаться занулить либо текущий диагональный элемент, либо поддиагональный, уж как получится, если занулится диагональный, то переставить две строчки. Для тех, кто не понял -- псевдокод.
Процедура Занулить поддиагональную часть столбца i
Для всех j > i
- Пока matrix[i][i] != 0 И matrix[i][j] != 0
- если abs( matrix[i][i] ) > abs( matrix[i][j] )
- Из строки i вычесть строку j умноженную на ( matrix[i][i] / matrix[i][j] )
- иначе
- Из строки j вычесть строку i умноженную на ( matrix[i][j] / matrix[i][i] ) конец цикла
Дело нехитрое, написал эти 20 строчек, отладил на небольших матрицах запустил... Вот блин. Переполнение int'а. Невлезаем в 32 бита на матрице 8 на 8.
Удивился. Почесал репу. Заменил 32-битное целое 64-битным. Запустил -- переполнился на матрице 11 на 11.
Ё-моё. В исходной матрице все элементы по модулю не больше 4, в основном она из нулей и единиц состоит. Какого хрена????
Подумал. Применил небольшую магию -- перед занулением поддиагональных элементов переставляю две строки матрицы так, чтобы на диагональ встал самый маленький по модулю элемент. Полгечало. Переполнения нет, всё посчиталось. Но вопрос остался.
Какого хрена оно переполняется? Откуда такой бешеный рост значений???