数値計算

粒子法での流体シミュレーションを少し並列化する

2008年の春ごろに物理シミュレーションにハマっていて、粒子法を使った流体計算をがんばっていたのです。 そのコードの並列化できる部分を並列化して、ちょっとだけ高速化しました。 こんな感じのシミュレーションです。これは10倍速にしています。 10倍速で…

このHMMがこんな出力をする確率は?

出力があるとして、HMMがその出力をする確率がどれだけあるか求めてみます。 前向きアルゴリズムといいますが、ほとんどビタビのアルゴリズムと同じプログラムになってます。 出力はこんな感じになりました。 出力:serrar rer te rr rrsertah reht 確率:0.00…

HMMの隠れたマルコフを覗いてみる

HMMではマルコフモデルが隠れていたのですが、こういう隠れているものを覗いてみたくなるのが人情というものです。 ということで、出力からマルコフ遷移を推定してみます。これにはビタビ(Viterbi)のアルゴリズムを使います。 ビタビのアルゴリズムは、要す…

隠れマルコフモデル(HMM)で遊ぶ

出力と状態が1対1に対応していなくて、状態を直接知ることができないようなものを隠れマルコフモデル(HMM:Hidden Markov Model)といいます。 で、その隠れマルコフモデルで遊ぶために、とりあえず状態遷移と出力文字候補を作ってみます。 この隠れマルコフモ…

マルコフモデルで遊ぶ

マルコフモデルというのは、要するに確率付きの状態遷移で状態ごとに出力があるものです。 このときの確率は、その状態だけではなく、一つ前の状態やもうひとつ前の状態に影響されたりもします。 で、これを言葉の単語の遷移でやったりすると、n-gramとかい…

積分でのsinとcosの関係をみてみる

微分と積分は逆の関係になります。 sinを微分するとcosになってたので、cosを積分するとsinになります。cosの微分は-sinだったので、-sinの積分がcos、つまりsinの積分は-cosになります。 このことがわかるプログラムを考えてみる。 dx間隔にデータを取ると…

微分でのsinとcosの関係

微分でのsinとcosの関係がわかるようなプログラムを考えてみる sinは微分するとcosに、cosは微分すると-sinになります。 f 微分 sin cos cos -sin で、これがどういうことかと。 sinは、0から始まって、だんだん増えていきます。で、その増え方はだんだん減…

ついでに、2階微分の差分式

上のふたつのテイラー展開を足すと となります。 ※f'''が消えたのでいままで省略されてたf''''をちゃんと書きました。 これを変形すると となります。 を誤差として切り捨てると、2階微分の差分式 が出てきます。 追記: これを使って、ラプラス方程式など2…

差分の数学的意味と誤差

数学では、なんとなく同じような値になるからといって、無条件に使えるような都合のいい話はゆるされません。 差分が微分の近似として使えるのにもちゃんとした理由があります。 ここで、テイラー展開を持ってきます ここで、x=x0+d、a=x0と置くと、次のよう…

微分は引き算

微分は引き算ってことがわかりやすくなるようなプログラムを考えてみる。 微分というのは、関数の傾きを求める操作です。 で、これを式として求めるときの規則として簡単なものに次のようなものがあります。 例えば次のような式を微分します。 そうすると、…

積分は足し算

積分は足し算だということを説明するためのプログラムを考えてみる。 aが加速度。dx,dyが速度。x,yが座標。 毎回x,yにdx,dyを足す。で、dyにaを足す。dxは一定。 そうすると、こんな動きに。 小学生のころ、わけわからず1、2、3、4と足していくとジャンプの…

サイン・コサインと円の関係

サインが縦、コサインが横 っていうのがわかりやすくなるようなプログラムを考えてみる //SinCos.java import java.awt.*; import java.awt.image.BufferedImage; import javax.swing.*; import static java.lang.Math.*; public class SinCos { public stat…

移流方程式の計算で数値計算の問題点をみる

数値計算でどんな微分方程式でもいい感じに計算できるかというと、そういうわけはもちろんなく。 という移流方程式をいくつかの計算をやってみて、理論どおりにならないことを見てみました。 左辺の第2項の偏微分について、いくつかのやりかたをしています…

ラプラシアンの計算をしてラプラス方程式を解く

温度の拡散などをあらわすラプラス方程式 を解いてみました。 ここでは、左右の辺が1℃、下の辺が0℃になるとき、板の温度分布がどうなるかという計算になっています。 ラプラシアンというのはのことです。 普通の偏微分であらわすと となります。 この形の偏…

ざっぱ〜んの処理を並列化する

とりあえず並列化について。 マルチコアなら並列化しないと、CPUがもったいないです。 こんな感じで、内側のループが独立しているところについて、Executorで並列処理するといいかも。 int procs = ManagementFactory.getOperatingSystemMXBean().getAvailab…

爆発しなくなった、ざっぱ〜んのソース

とりあえず貼り付けときますね。 「時間感覚」のところと「粒子が近づきすぎてないか」のところが追加した処理。

ざっぱ〜ん、完成した!

粒子が近づきすぎないようにする処理を入れたら、爆発しなくなりました。 でも、処理が重くて1秒分のシミュレーションに1時間かかるし、前回のようなFlashでもスムーズには作れないし、どうしたものかと思ってたら時間が過ぎてしまいました。 ってことで、JM…

粒子法で、ざっぱ〜ん!

時間間隔がまずかったみたいで、ちゃんとざっぱ〜んとなりました。 最後爆発しちゃいますが。 あと、実際には粒子の速さによって計算する時間間隔を調整しないといけなかったり、粒子が近づき過ぎないようにする処理が必要になるのですが、まぁ適当に動いた…

その後のさっぱ〜ん

必要な計算はすべて行ってみたんですが、爆発しちゃいます。 圧力計算がおかしいくさい

爆発した!

笑える。 http://rps13.180.jp/parc.swf

ざっぱ〜ん、とやりたい

こいつがざっぱ〜んと動くには、どうしたらいいんやろか?

弦のシミュレーションをJavaScriptで

昨日の、弦の動きのシミュレーションをJavaScriptでやってみました。手軽に見れます。 Œ·‚Ì“®‚«‚̃Vƒ~ƒ ƒŒ[ƒVƒ‡ƒ“

弦の動きのアニメーション

昨日のはほぼ本の通りのプログラムだったんで、それはそれであれなので、ちゃんと自力で組んでみました。 ということで、アニメーションさせてみました。すごく、なまめかしい動きをします。なかなか面白い動きなので、ぜひ自分でも動かしてみてください。こ…

波動方程式を解いてみるよ

機械学習は飽きてきたので、今日は波動方程式を解いてみた。 弦をひっぱって手を離したらどういう運動をするかというシミュレーション。右下方向が時間のたつ向き つっても、ほぼ、本に載ってる通りなんだけど。初期値が違うくらい。 ちなみに数値計算の参考…

最急降下法で極小値を求める

関数の極小値を求めるために最急降下法を使ってみます。 y=f(x)の極小値を求めるとして、まず適当な値aをとります。f(x)のx=aでの微分f'(a)をとったとき、これが正であればaの値を減らし、f'(a)が負であればaの値を増やすと極小値に近づくはずというのが、最…

ルンゲクッタ法で2階連立微分方程式を解く

今回は、空気抵抗があるときの物体の軌跡を描いてみました。下のプログラムでは、実際にはアニメーションさせてます。 ※追記 2023/4/30 動画を投稿しました空気抵抗がある中でモノを投げるシミュレーションhttps://t.co/zd6SCX8qE1 pic.twitter.com/IM2HbKom…

ルンゲクッタ法で1階微分方程式を解く

微分方程式を解くルンゲクッタ法を試してみます。 今回は という微分方程式を解いて、空気抵抗がある自由落下の速度を求めてみます。g=9.8、k=0.01としてます。 空気抵抗により、速度が一定になることがわかります。 プログラム中、funcメソッドが求める微分…

ガウス・ザイデル法で解ける連立方程式を生成

対角要素の絶対値が大きい方程式だと収束するということなので、そういうデータをでっちあげる。 public class LinearSystem { private static int getRandom(int n){ int i = (int)(Math.random() * n + 1) * ((Math.random() < .5) ? 1 : -1); return i; }…

連立一次方程式を解く

ガウス・ザイデル法という解き方。 だんだん解が収束していきます。 4w +1x +2y= 21 2w +5x -1y= 14 1w +3x -7y=-24 w(0)= 5.3 x(0)= 0.7 y(0)= 4.5 w(1)= 2.8 x(1)= 2.6 y(1)= 4.9 w(2)= 2.1 x(2)= 2.9 y(2)= 5.0 w(3)= 2.0 x(3)= 3.0 y(3)= 5.0 w(4)= 2.0 …

回帰直線を引く

基本的な処理として。 public class Regression { public static void main(String[] args){ //データ作成 List<Point2D.Double> points = new ArrayList<Point2D.Double>(); double a = 10; double b = .7; int n = 100; Random r = new Random(); for(int i = 0; i < n; ++i){ double x = r</point2d.double></point2d.double>…