Mathematicaで書いているポアソンソルバーの助けを探しています。コードは非常に長く、アレイが接続されていますが、詳細についてはhttp://pastebin.com/uSrSDcW6を参照してください。
ポアソンの式から導出された中央差分法を使用して、電荷密度が与えられた電圧を計算しています。電圧を計算した後、データセットの収束をテストします。収束しきい値を10^-1000+のオーダーに設定しています。フェイルセーフとして、何かがうまくいかない場合に備えて、10000回の反復後にキックアウトするようにループを設定しています。正気のためにループカウンターを設置しています。収束しきい値が10^-100に設定されている限り、プログラムは正常に実行されているようです。
私の質問はこれです：しきい値を更新しても、たとえば10 ^ -100、10 ^ -150の場合、計算は633回の反復後に停止し、ループから抜け出します。私はこれで助けをいただければ幸いです、私は完全に立ち往生しています。このフォーラムの誰にとっても説明できるコメントをプログラムに追加しました。繰り返しになりますが、この説明は限られていることを知っているので、完全なプログラムについては、添付のURLhttp ://pastebin.com/uSrSDcW6を参照してください。
*更新10/9/12 ***問題を16桁のマシン精度に切り分けました。マシンの最大精度が10^309になるまで開く必要があります。Mathematicaヘルプはこれを行う方法についてはまばらです。ex「N[MachinePrecision、50]」。これをプログラムのどこに設定して、すべての計算に適用しますか？それが助けになるなら、ここにループを貼り付けます
Vnew / Vold / RHOは10x10x34です行列イプシロンは定数です
（ConvergenceLoopをOに初期化します-これは、必要に応じてループからキックアウトするためのフェイルセーフとして機能します）
ConvergenceLoop = 0;
（収束をゼロに初期化します）
収束=0;
While [Convergence == 0 && ConvergenceLoop <10000、
（すべてのi、j、k要素を実行し、新しい電圧値を計算します）
Do [Vnew [[i]] [[j]] [[k]] =（1 /（2 / deltaX ^ 2 + 2 / deltaY ^ 2 + 2 / deltaZ ^ 2））*（（（Vold [[i + 1]] [[j]] [[k]] +
Vold[[i - 1]][[j]][[k]])/(deltaX^2)) + ((Vold[[i]][[j + 1]][[k]] +
Vold[[i]][[j - 1]][[k]])/(deltaY^2)) + ((Vold[[i]][[j]][[k + 1]] +
Vold[[i]][[j]][[k - 1]])/(deltaZ^2)) + ((Rho[[i]][[j]][[k]]/Epsilon))), {i, 2, 9}, {j, 2,9}, {k, 2, 33}];
（収束したと仮定して、テストが定義された収束しきい値を超える最初の値に達したときにループがトリガーされるようにします）
収束=1;
（これは収束テストです。ユーザー定義の収束しきい値）
Do[If[Vold[[i]][[j]][[k]] == 0, Null,
If[(Vnew[[i]][[j]][[k]] - Vold[[i]][[j]][[k]])/Vold[[i]][[j]][[k]] > .0000001, Convergence = 0;
(*This is purely diagnostic. I added a Tracker point to follow the convergence along.
任意の要素で定義されたユーザー*）
If[i == 5 && j == 5 && k == 10,
Print[ "Tracker Point" (Vnew[[i]][[j]][[k]] -
Vold[[i]][[j]][[k]])/Vold[[i]][[j]][[k]]], Null],Null]], {i, 2, 9}, {j, 2, 9}, {k, 2, 33}];
（VnewとVoldがゼロ以外になるまで、最初の反復を無視します）
If [ConvergenceLoop <2、Convergence = 0、Null];
（VlandをVnewで進化させる）
Vold = Vnew;
ConvergenceLoop ++;]
（将来の計画のためにSessionTimeを追加しました）
If [ConvergenceLoop == 10000、
Print ["Convergence Loop Limit Reached。"（SessionTime [] / 3600）]、
Print["収束ループ制限に達していません。"]];
（whileループから抜け出しました。つまり、データが収束したため、収束した値を出力します）
If[収束==1、
Print [ConvergenceLoop "おめでとうございます！" MatrixForm [Vnew]]、Print ["Did Not Converge！"]];