まとめn回の計算に1回だけ結果を取り出す方法VBAでルンゲ・クッタ法のような数値計算を行う中で、刻み幅を細かくしていくと、その分だけ計算結果の数が多くなります。0から10秒まで計算したいとき、刻み幅を0.1にすれば100個、0.01にすれば1000個、0.001にすれば1万個の計算結果を得ます。ストーク刻み幅をな~細かくしていくと、それだけ精度の良い計算ができるし、イマドキのパソコンだと10万回や100万回くらい簡単にループできる。10万個のプロットでグラフを書くのは、なんつーか、だるい。100万個とかありえへん。シママあ~わかるな~その気持ち。100万個のプロットとか描いたことないけど。グラフ重いし、刻み幅を変更したらいちいちグラフの仕様を変更しなくちゃいけないのがメンドクサイ。ストークこういう作業をいちいちやるのも嫌だけん、何とかならん?シママなるわよ!要は、n回に1回だけ出力するプログラムを書けばいいんでしょ?ワタシが書いたやつ見せてあげる!というわけで、本記事では、VBAにおいて、n回に1回だけ出力するプログラムの備忘録を残したいと思います。実装例シママ今回は、振動の運動方程式を例題にしてみるね。この状態でワークシートを作って・・・シママプログラムを実行すると・・・http://shimaphoto03.com/wp-content/uploads/2018/01/rk.mp4シママこんな風に結果が出るの!スゴイでしょ?ストークこれは・・・。50秒間積分するために刻み幅を0.01秒にしているから5000回ループしている。10回に1回しか出力しないように設定しているから、グラフのプロット数は500個だけなんやな。シママそういうこと!ソースコードシママソースコードは以下の通り。このコードは前回の振動計算の記事と同じものよ。使うと、様々な自然現象を記述することができます。その常微分方程式の数値解法としてルンゲクッタ(Runge-kutta)法があります。*/ソースコード解説シママ太字の部分見てみて!刻み幅hに対して、新たに「h2」という変数を用意して、h2回に1回結果を出力するようにしたの。h2ってなんだ?シママ6は位置合わせ。計算結果が6行目から出てるでしょ?それで、i/h2っていうがポイント!ストークえっと、h2=10なら、iが0のとき、jは6,iが1のとき、jは6.1iが2のとき、jは6.2・・・iが10のとき、jは7。あ、これでi=10の時の結果を7行目に、i=20の時の結果を8行目に・・・ってことか。シママそうそう!そして、jって整数型(integer)だから結局i=1~9のとき6か7になるんだけど、i=1~9のときに出力しないでね!って言うのが次のif文。確かにi=1~9のとき、このifはFALSEになるたいね。i=0,10,20・・・のとき、iがh2で割り切れてTrueになって、i=0のときj=6になって6行目、i=10のときj=7になって7行目、i=20のときj=8になって8行目に出力されるのか~。ストークあ~このちょっと、こそばゆいところに手が届く感じええなぁ~。ラクできるとこはラクしなきゃね~。マクロってそういうものだからね~。ストークこれなら、新しい方程式を100万回くらい計算させて、10万回計算させた結果と照合して問題なければ、「その方程式に関しては刻み幅を大きくしてもよい」という判断ができるから、ムダがないな。非線形性の大きい問題にもチャレンジできる。100万回計算させている間にタバコ吸いに行ける(※)からさらにムダがない。8GBで100万回計算させるのに必要な時間は3分4秒でした。シママひせんけい?ストーク2次関数とかexpとか非線形やろ?直線でないってことたい。4次関数とかタチ悪いたいね。※ちなみに、生物生存闘争モデル(例:ウサギとキツネの個体数の時間変化)のLotka-Volterra方程式も非線形の微分方程式です。右辺のxyの部分が2次になっています。解曲線も極大を超えたら急激に減少するような挙動を示します。TS38数値計算による微分方程式解法の基本~runge-kutta法~ルンゲ・クッタ(runge-kutta)法は、コンピュータを使った微分方程式解法の基本として工学部の基礎教育に登場します。中でも4次のルンゲ・クッタ法はアルゴリズムが簡単な上に高い精度...シママそっか・・・。4次関数は精度が荒いと2次関数と間違える可能性もあるよね・・。いう危険も孕んでいるから、数値計算結果の過信はよくないんやが・・・。100万回の計算結果と10万回の計算結果が同じなら、計算の信頼性が強くなるたいね。ストーク例えばな、今回、n回に1回出力するプログラムが出来たから、100万秒の現象を刻み幅1秒で計算して1万回に1回出力したときのグラフと100万秒の現象を刻み幅10秒で計算して1,000回に1回出力したときのグラフを比較して一致すれば刻み幅10秒で信頼たる計算が出来とると判断できるやろ?一応言うが、支配方程式が間違っとるっていうのはナシやで。両方ともプロットの数は100個だもんね。グラフのフォーマットも参照範囲も変えなくていいよね。(そうか・・・プログラミングを自然現象の解析に使っている人はそういう事を考えるのか・・・・。これは助かるよ!ありがとう!シママいえいえ。いつも色々教えてもらってるから、こういうことくらいは・・・ね。ストーク専門分野を異にするダチってええなぁ。まとめRunge-kutta法で微分方程式を解いていて、n回に1回出力できるようにならないかな~と思って、試行錯誤した結果をまとめました。VBAではCなどとは違って、一工夫が必要で、Cと同じようにMODだけ使っていると、こうなります。結果を見れば「そりゃ、そうか。」って感じなのですけどね。VBAはあくまでツールで、それ専門じゃないので、またVBAでプログラムを組むのは数か月後ということは、あり得ることです。いった意味で、同じところで苦戦しないためにも備忘録を残しました。コメントをキャンセルメールアドレスが公開されることはありません。