私はPython用のウィルソンのスペクトル密度因数分解アルゴリズム[1]の実装を書き込もうとしています。このアルゴリズムは、[QxQ]行列関数を反復的に平方根に因数分解します(これは、スペクトル密度行列のニュートンラプソン平方根ファインダーの拡張のようなものです)。
問題は、私の実装がサイズ45x45以下の行列に対してのみ収束することです。したがって、20回の反復後、行列間の二乗差の合計は約2.45e-13になります。ただし、サイズ46x46の入力を行うと、100回程度の反復まで収束しません。47x47以上の場合、行列は収束しません。エラーは約100回の反復で100から1000の間で変動し、その後非常に急速に大きくなり始めます。
このようなものをデバッグするにはどうすればよいですか?気が狂うような特定のポイントはないようで、行列が大きすぎて実際に手動で計算を試みることはできません。誰かがこのような奇妙な数値のバグを見つけるためのヒント/チュートリアル/ヒューリスティックを持っていますか?
私はこれまでこのようなことを扱ったことがありませんが、あなたの何人かが持っていることを願っています...
ありがとう、-ダン
[1]GTウィルソン。「マトリックススペクトル密度の因数分解」。SIAMJ.Appl。数学(第23巻、第4号、1972年12月)