数値安定性の問題をデバッグするための戦略？

Question

私はWilsonのスペクトル密度分解アルゴリズム[1]をPython用に実装しようとしています。このアルゴリズムは[QxQ]行列関数をその平方根に反復的に因数分解する（これはスペクトル密度行列のNewton-Raphson平方根ファインダの拡張の一種である）。

問題は、私の実装が45x45以下のサイズの行列にしか収束しないということです。したがって、20回の反復の後、マトリックス間の合計の二乗差は約2.45e-13である。しかし、私が46x46のサイズを入力すると、それは100回程度の反復まで収束しません。 47x47以上の場合、行列は収束しません。約100回反復する間に誤差は100と1000の間で変動し、その後非常に迅速に増加し始める。

このようなものをデバッグしようとするとどうなりますか？それが狂った具体的なポイントは何もないように見えますが、実際に計算を実際に試みるには行列が大きすぎます。誰もこのような奇妙な数値のバグを見つけるためのヒント/チュートリアル/ヒューリスティックを持っていますか？

私は前にこのようなものを扱ったことがありませんが、私はあなたのいくつかは持っている願っています...

おかげで、 - ダン

[1] G. T.ウィルソン。 "行列のスペクトル密度の因数分解"。 SIAM J. Appl。 Math（Vol 23、No. 4、1972年12月）

Answer 1

この質問はscipy-userメーリングリストで質問することをお勧めします。一般的に、リストに載っている人々は、数値計算に非常に熟練しているようで、実際に役立ちます。

そうでなければ、私はあなたはそれが数値精度/浮動小数点の丸めの問題だと思うなら、私は...任意のアイデアを持っていない怖い、あなたが試みることができる最初の事はアップfloat128にすべてdtypesバンプと見ています何か違いがあれば。

Answer 2

Interval arithmetic助けてもらえますが、性能が実際に興味のあるマトリクスサイズで意味のあるデバッグを可能にするのに十分なのかどうかはわかりません（2桁の減速度、 -WWは、SWが重い「間隔」の浮動小数点演算を支援し、どの間隔が広すぎるのか、いつ、どこで、どのような間隔で増加するかについてのチェックを追加しました。