Должны ли матричные элементы, уступающие машинной точности, обрезаться до нуля?

У меня есть код в сложной арифметике с двойной точностью, который формирует большие разреженные матрицы (благодаря PETSc) для решения задач, где обычно требуется высокая точность (результаты сходятся, скажем, до 7/8 цифр)

Полученная линейная система Ax=b решается с помощью параллельного LU-разложения. Я стремлюсь решать большие трехмерные задачи, чтобы главный размер матрицы мог достигать нескольких десятков миллионов.

Матрица содержит элементы, полученные в результате умножения весов метода конечных разностей на коэффициенты различных решаемых физических уравнений и метрические члены.

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

  • «Очень маленькие» записи могут загрязнить процесс инверсии LU и сделать его нестабильным при увеличении заполнения.

  • Если в процессе используются элементы матрицы, близкие к машинной точности, я бы предположил, что результат любой арифметической операции с этими элементами можно считать «ненадежным»? Если да, то я не понимаю, почему они должны быть сохранены.

Я знаю, что машинная точность - это другое понятие, чем наименьшее представимое положительное число. Поэтому возможно, что мои рассуждения концептуально неверны, и я хотел бы получить ваше мнение или исправления. Спасибо!


person circuitbreaker    schedule 25.07.2018    source источник
comment
Не мной (не знаю об этом достаточно), но это может быть отклонено или закрыто, потому что это, кажется, вопрос мнения, и это не по теме здесь.   -  person Rudy Velthuis    schedule 25.07.2018
comment
Абсолютный «машинный эпсилон», как правило, не является правильным порогом для использования. Мы не знаем, каковы значимые значения в вашей матрице: 1, 1e-20 или 1e20, поэтому мы не знаем, каковы они по отношению к машинному эпсилону. Более целесообразно установить порог, равный около величине значимых значений в ваших вычислениях, умноженной на машинный эпсилон. То, что означает «вокруг», зависит от контекста. Любая отдельная операция имеет максимальную ошибку ½ эпсилона, но объединение нескольких операций может привести к большим ошибкам…   -  person Eric Postpischil    schedule 25.07.2018
comment
… Если ваша инверсия LU численно нестабильна, то любые изменения плохи. Любые малые значения, возникшие из-за ошибок округления, тем не менее являются полезной информацией, которая является лучшим результатом, который может дать арифметика с плавающей запятой. Другими словами, их не следует считать ненадежными. Низкая информация лучше, чем информация. Их отбрасывание приведет к отбрасыванию информации, которая, вероятно, отдалит вычисленные результаты от правильных ответов, а нестабильные процессы составят ошибки. Таким образом, реальный вопрос заключается в том, могут ли ваши вычисления хорошо переносить небольшие изменения, и будут ли в результате...   -  person Eric Postpischil    schedule 25.07.2018
comment
…увеличение производительности стоит снижения качества результата. Это то, что требует больше информации, чтобы ответить.   -  person Eric Postpischil    schedule 25.07.2018
comment
@EricPostpischil: подумайте о том, чтобы сделать это ответом, даже если вопрос, возможно, не по теме.   -  person President James K. Polk    schedule 26.07.2018
comment
@JamesKPolk: Это немного выходит за рамки моей обычной практики, поэтому я не готов писать авторитетно, по крайней мере, до тех пор, пока вопрос не будет доработан, и у других будет возможность ответить.   -  person Eric Postpischil    schedule 26.07.2018


Ответы (1)


Ответ на этот вопрос «это зависит», как уже предлагали другие люди. Однако можно точно сказать, от чего зависит ответ.

Предположим, вы хотите решить линейную систему,

Ax = b,

где A и b — матрица и вектор правильной размерности. Теперь вы решаете возмущенную линейную систему

(A + E)  \шляпа х = b ,

где матрица возмущения E удовлетворяет условию

\| E \|_\infty \le \alpha \| A \|_\infty .

Если \alpha \kappa_\infty  (A) \le \tfrac12, то

\frac{\  | x - \hat x \|_\infty}{\| x \|_\infty} \le 4 \alpha\,\kappa_\infty(A) ,

где http://latex.codecogs.com/png.latex?%5Ckappa_%5Cinfty  %28A%29 — это состояние A, измеренное в норме суммы строк.

(Этот результат можно найти, например, в книге Голуба ван Лоана Матричные вычисления, глава 3.)

Другими словами, если вы измените элементы матрицы, результат изменится на сумму возмущений по строке (\| \cdot \|_\infty - норма суммы строк), умноженная на кратное условие матрицы.

Следовательно, возмущение матрицы разрешено, если суммарное возмущение по строке относительно нормы суммы строк A все еще мало при умножении на условие матрицы.

Этот результат справедлив, конечно, только для точных вычислений. Однако, как мы видели, для того, чтобы иметь возможность возмущать матрицу, число обусловленности матрицы должно быть достаточно малым. Таким образом, в этом случае, если возмущение мало, число обусловленности возмущенной матрицы также будет мало, и в этом случае LU-разложение не представляет проблем.

person H. Rittich    schedule 26.07.2018