JASRACが怖くて日和りました。
Swiftの話だと思った?
残念、ジュリーでした。
前回(半年前だけどな)予告してた空気抵抗を考慮した自由落下を数値的に解くをやってみます。
まずは、数字見ても面白くもなんともないので、cube-600.htmlを加工して、平面の上に球を落とすことにします。
cube-600.htmlは「WebGLでいきますぜい」で紹介したページね。
↓これね
http://www.tetera.jp/xcc/ameba/cube-600.html
何やってるかや、使い方は以下を順に読みましょう。
altitude(高度)
という変数を用意して、このaltitudeに設定された高さから、高度0mに配置したcubeまでを1画面にとらえるようにします。
具体的には、カメラポジションとカメラの向きをaltitudeをもとに計算することで上のような画面になる。
var altitude = 828; // 初期高度 m ・・・ function threeStart() { ・・・ camera.position.set(0, altitude * 1.1, altitude * 1.1); camera.up.set(0, 1, 0); camera.lookAt( {x:0, y:altitude * 0.55, z:0 } );
ここで指定してる3次元座標の1はそのまま1mと考えてください。座標の単位をなんにするかは作る人の勝手です。m単位でもkm、mmなんでもいいけど、自分がどういった単位で座標を扱ってるかは、忘れないように。
物理ではMKS単位系で計算式を表現することが多いのでmにしました。
前回の抵抗を考慮した自由落下の式も単位はMKS単位系です。
↓こいつね
cubeの大きさやsphereの大きさも微妙に調整しますた。
ここら辺は、興味がある人は各自で
function arrengeObjects()
を見てください。特に重要な話じゃないです。
重要な点は、ページの表示から経過した時間によって、sphereを動かさなければいけないこと。
なので、時刻を記憶するための変数を用意。
var lastDate; // 最後の移動処理の時刻
こいつに、ページ表示時に呼ばれるthreeStart関数内で
lastDate = new Date();
とすることで、この処理が実行された時点の時刻を記憶させています。
そして繰り返し呼び出されるloop関数で
var now = new Date(); var deltasec = (now.getTime() - lastDate.getTime()) /1000;
として、記憶した時刻(lastDate)からの経過時間を秒単位で計算します。getTime()はミリ秒単位なので1000で割ってる。
これで経過時間が測定できたので、あとはこの経過時間を使ってオイラー法(前回説明した、移動距離を数値的に解くやり方の呼び名です)で位置を更新するわけですよ。そして最後に、次回のloop関数呼び出し時に、今の時刻からの差分を計算させるためにlastDateを更新する。
lastDate = now;
loop関数は、window.requestAnimationFrame(loop)によって、画面の更新間隔(大抵、毎秒60回)で呼び出されるんで、こんな面倒なことせずにdeltasecを1/60=0.17に固定してもいいんですけどね。実際、私のマシンではdeltasecは0.017秒くらいになってたし…
ただ、毎秒60回の呼び出しはあくまで努力目標なのと、必ずしも画面の更新間隔(リフレッシュレートと呼ぶのだよ)が60Hz(ヘルツと呼ぶのだよ、いや、わかっとるわ)とも限らないんで、今回は真面目に上記のような経過時間測定をしてます。
で、いよいよオイラー法。
オイラー法で位置を更新する部分はupdateTime関数というのを作ってまとめています。
当然、変数velocityとdisplacementは値を記憶するためにupdateTime関数外で定義ね。
var velocity = 0; // 初速度 m/s var displacement = 0; // 変位 m ・・・ function updateTime(deltasec) { ・・・ displacement += velocity * deltasec; sphere.position.set(0, altitude + displacement, 0); }
velocityとdisplacementがオイラー法で計算される速度と変位です。
altitudeの初期高度に現時刻での変位を加算して、現時刻で高度としてます。
それをsphereのy座標に設定して位置を移動させている。
で、ここで残念なお知らせ。
粘性抵抗側を求める式
6πrηv π=3.14 r=球体の半径(m)
しか見つからない。Lの人間の体積はだいたい
L = 0.06
みたいで、η(エータと読むのだよ)は1.8 x 10-5てのがわかったんだけどαが見つからない。
↓ここね
http://www.mterm-pro.com/machine-yougo/fluid-dynamics/water-air-bussei.html
しょうがないので、球体の体積が0.06になるrを計算してα=76と決めました。
もういっそのこと、粘性抵抗側は無かったことにしたいな。
まあ、気を取り直して慣性抵抗側。
こっちはスムーズに見つかって、βを形状を円柱とした1.2、ρ(ピーじゃないよローだよ)を1.2kg/m3と決定。
↓ここから
https://ja.wikipedia.org/wiki/抗力
進行方向断面積は適当に0.33m2、質量mは60kgとして、いざオイラー法!
function updateTime(deltasec) { var a = 76; // 粘性側係数 var b = 1.2; // 抵抗側係数(円柱) var L = 0.06; // 体積(人間、適当) m3 var n = 1.8 / 100000; // 粘度(空気) Pas var Lnv = L * n * velocity; var s = 0.33; // 進行方向断面積(人間、適当) m2 var p = 1.2; // 密度(空気) kg/m3 var spv2 = s / 2 * p * velocity * velocity; var m = 60; // 質量 kg velocity += - ((9.8 + (-a * Lnv - b * spv2) / m)) * deltasec; displacement += velocity * deltasec; sphere.position.set(0, altitude + displacement, 0); }
画面下にある現在の速度や経過時間、高度を更新する
function updateInfo(now)
に関しても説明はパスします。計算結果を画面に表示してるだけです。
「走れJavaScriptちゃん!」を読むとか、JavaScript getElementById innerHTMLあたりで検索すればすぐわかります。
あとthreeStart関数でやってる
light = new THREE.DirectionalLight(0xFFFFFF, 100); light.position.set( 0, 2000, 0 ); scene.add(light); ・・・ renderer.setClearColor(0xEEEEEE, 1.0); renderer.shadowMap.enabled = true; light.castShadow = true; light.shadow.camera.left = -altitude; light.shadow.camera.bottom = -altitude; light.shadow.camera.right = altitude; light.shadow.camera.top = altitude; light.shadow.camera.far = 3000; cube.receiveShadow = true; sphere.castShadow = true;