日本語プログラミング言語なでしこで数値計算(1) まずは何はともあれ4段4次のルンゲクッタ法でローレンツ方程式を計算してみる。
最近プログラム記事のネタ切れ感が半端ない、、、Visual C#を使い始めたがまあ仕事で使っている部分もあるのであまりかけない。
流行りのRUSTなどより、使っている人がとても少なそうな言語でいろいろ遊ぼうと思ったり。
そこで、日本語でプログラミングできる「なでしこ」の存在を思い出した。
おお、いつの間にかブラウザベースでプログラミングができるようになっているのか。これはちょうどいい。
しかも保存して公開できるし。誰もやってなさそうな数値計算ネタをやっていこう。
まずは何はともあれローレンツ方程式だ。
リンクはこちら:
https://n3s.nadesi.com/index.php?page=1875&action=show
こんな感じで日本語でプログラムして、
実行すると、、、おお!おなじみのグラフが描けた。
普通に何でもプログラムできそうだ(まあJavascriptに変換されているそうなのでJavascriptで出来ることは何でもできそう)。
ちょっと遊んで行ってみよう(続く)。マンデルブロ集合とかSIRもでるとかどうだろう。
テキストでのプログラム:
# ローレンツ方程式
# dx/dt に-σ(x-y)
# dy/dt に-xz + rz -y
# dz/dt にxy - bz
# を4段4次のルンゲクッタ法で計算する。
全描画クリア。
# 画面設定
画面幅=描画中キャンバスの「width」をDOM属性取得。
画面高=描画中キャンバスの「height」をDOM属性取得。
「24px sans-serif」に描画フォント設定。
[25,20]に「ローレンツ方程式(ルンゲクッタ法で計算)」を文字描画。
x最大値に25を代入する。
x最小値に-25を代入する。
y最大値に50を代入する。
y最小値に-30を代入する。
x最小値と0からx最大値と0まで黒色で直線プロット。
0とy最小値から0とy最大値まで黒色で直線プロット。
時間に0を代入する。
時間ステップに0.01を代入する。
#初期値
X0に0.1を代入する。
Y0に0.1を代入する。
Z0に0.1を代入する。
10000回繰り返す
#ルンゲクッタ法のメイン計算
KX1にX0とY0とZ0の関数fXを代入する。
KY1にX0とY0とZ0の関数fYを代入する。
KZ1にX0とY0とZ0の関数fZを代入する。
KX2に(X0+時間ステップ×KX1/2)と(Y0+時間ステップ×KY1/2)と(Z0+時間ステップ×KZ1/2)の関数fXを代入する。
KY2に(X0+時間ステップ×KX1/2)と(Y0+時間ステップ×KY1/2)と(Z0+時間ステップ×KZ1/2)の関数fYを代入する。
KZ2に(X0+時間ステップ×KX1/2)と(Y0+時間ステップ×KY1/2)と(Z0+時間ステップ×KZ1/2)の関数fZを代入する。
KX3に(X0+時間ステップ×KX2/2)と(Y0+時間ステップ×KY2/2)と(Z0+時間ステップ×KZ2/2)の関数fXを代入する。
KY3に(X0+時間ステップ×KX2/2)と(Y0+時間ステップ×KY2/2)と(Z0+時間ステップ×KZ2/2)の関数fYを代入する。
KZ3に(X0+時間ステップ×KX2/2)と(Y0+時間ステップ×KY2/2)と(Z0+時間ステップ×KZ2/2)の関数fZを代入する。
KX4に(X0+時間ステップ×KX3)と(Y0+時間ステップ×KY3)と(Z0+時間ステップ×KZ3)の関数fXを代入する。
KY4に(X0+時間ステップ×KX3)と(Y0+時間ステップ×KY3)と(Z0+時間ステップ×KZ3)の関数fYを代入する。
KZ4に(X0+時間ステップ×KX3)と(Y0+時間ステップ×KY3)と(Z0+時間ステップ×KZ3)の関数fZを代入する。
X1 にX0 + 時間ステップ×(KX1+2*KX2+2*KX3+KX4)/6を代入する。
Y1 にY0 + 時間ステップ×(KY1+2*KY2+2*KY3+KY4)/6を代入する。
Z1 にZ0 + 時間ステップ×(KZ1+2*KZ2+2*KZ3+KZ4)/6を代入する。
X0とY0からX1とY1まで赤色で直線プロット。
X0とZ0からX1とZ1まで青色で直線プロット。
X0=X1
Y0=Y1
Z0=Z1
ここまで。
#ローレンツ方程式
●(XとYとZの)関数fXとは
-10*(X-Y)を戻すこと
ここまで
●(XとYとZの)関数fYとは
-Y-X*Z+26*Xを戻すこと
ここまで
●(XとYとZの)関数fZとは
X*Y - 8*Z/3を戻すこと
ここまで
●(xの)x座標変換とは
それ=画面幅 * (x - x最小値)/(x最大値 - x最小値)
ここまで
●(yの)y座標変換とは
それ=画面高 * (y - y最大値)/(y最小値 - y最大値)
ここまで
●(x1とy1からx2とy2までcで)直線プロットとは
cに線色設定
[x1のx座標変換, y1のy座標変換]から[x2のx座標変換, y2のy座標変換]まで線描画
ここまで
« 「クジラアタマの王様」(伊坂幸太郎さん)を読んだ。面白かった!最初挿絵?と思っていたイラストが文と並んでお話を形成していることに気付いてから特に!ハシビロコウの絵が特に好き!コロナ前に書かれた作品なのにそれを予期したかのような展開も。 | トップページ | 高周波(RF・マイクロ波・ミリ波・5G)関連ニュース(10/15) IEEE Microwave Magazineはマイクロ波磁気工学特集。ランダウ・リフシッツ・ギルバート(LLG)方程式!Microwave JournalではiPhone14が衛星通信に使うGlobalstarの記事!アナデバとキーサイトがフェイズドアレイアンテナで協業、エリクソンはEバンド(70/80GHz)が大事とレポート、など。 »
「パソコン・インターネット」カテゴリの記事
「学問・資格」カテゴリの記事
「日記・コラム・つぶやき」カテゴリの記事
« 「クジラアタマの王様」(伊坂幸太郎さん)を読んだ。面白かった!最初挿絵?と思っていたイラストが文と並んでお話を形成していることに気付いてから特に!ハシビロコウの絵が特に好き!コロナ前に書かれた作品なのにそれを予期したかのような展開も。 | トップページ | 高周波(RF・マイクロ波・ミリ波・5G)関連ニュース(10/15) IEEE Microwave Magazineはマイクロ波磁気工学特集。ランダウ・リフシッツ・ギルバート(LLG)方程式!Microwave JournalではiPhone14が衛星通信に使うGlobalstarの記事!アナデバとキーサイトがフェイズドアレイアンテナで協業、エリクソンはEバンド(70/80GHz)が大事とレポート、など。 »




コメント