シンプルな*ゲーム物理飛行機の作成

*-3次元で。



警告：さらなる推論は間違いかもしれません。私の経験は、フライトシミュレーターと理論力学のコースに限られています。 空力およびゲーム物理エンジンの知識は非常に不足しています。 注目を集める画像は、スクリーンショットではなく写真です。

「飛行機よりもシンプルなものは何ですか？」 揚力は速度の2乗に比例し、エンジンは前進します。すべてが単純です」-夏にそのような考えが浮かび、ゲームを書くために座った。 夏が過ぎ、いくつかの熊手が集められ、プロジェクトに追加する予定のリストが非常に増えました。

この記事では、スケジュールなどに気を取られることなく、物理的なコンポーネントについて書きます。 これが誰かの助けになることを願っています-このトピックに関する情報はインターネットにはあまりありません。

Scalaで動作する、あまりコードではない例があります（言語の本質を理解するために知る必要はありません）。 さらに、libgdxのベクトル、行列、および四元数のクラスを使用します。 実装の機能について質問がある場合は、エンジンコードは公開されています。 便宜上、ベクトルには+ =、* =、+、-のタイプのメソッドが追加されています。これは、scalaで可能であるためです。 メソッドを追加することも便利でした：=値による割り当て。

コード：

v:=(0, 1, 2) v2:=v 

同等の

 v.set(0,1,2) v2.set(v), 

四元数を除く：

set(x,y,z,w)しかし、 :=(w,x,y,z)

2番目のオプションは、私にとってはもっと身近なものであり、setは使用しません。 私はこれまでにコードをアップロードしません。小さなゴミがあります。それは時々認識を超えて書き換えられます。

それでは始めましょう。

モデル位置のクラスを書きましょう：

 class Position(val pos:Vector3f, val rot:Quaternion) 

ご想像のとおり、クラスには、モデルの重心の位置と世界に対する相対的な方向を記述する2つの最終フィールドがあります。
libgdxのマトリックスは、matrix.idt matrix.idt().translate(position.pos).rotate(position.rot)ように取得できますmatrix.idt().translate(position.pos).rotate(position.rot)

導関数のクラスを紹介します：

 class Derivative(val linear: Vector3f, val angular: Vector3f) 

このクラスは、速度（1次導関数）の記述だけでなく、加速度と力の記述にも便利です。

やめて クォータニオンはどのようにしてベクトルになりましたか？ 実際、角速度、加速度、モーメントは通常、ベクトルで記述されます。 角速度と向きの大きな違いは、向きが「ループ」していることです; 2パイの回転は「ゼロ」回転に相当します。 それどころか、角速度は原則として制限されない。

四元数qベクトルvの対数の操作を導入できます。これは回転軸の方向に向けられ、その長さはラジアン単位の回転角に等しくなります。 出展者は逆の操作です。 q == exp（v * 2）== exp（v）* exp（v）

マッピングベクトル->クォータニオンは一意ですが、その逆は一意ではありません。 dt中の角度alphaを通る回転は、角速度（alpha + 2 * pi * n）/ dtに対応できます。nは任意の整数です。

明らかに、角速度wでの時間dtに対して、回転q = exp（w * dt）です。

コード：

 class QuaternionS extends Quaternion{ ... def :=(log: Vector3): Unit = { val len = log.len if (len > 0.0001f) { val a = len / 2 val m = Math.sin(a).toFloat / len this :=(Math.cos(a).toFloat, log.x * m, log.y * m, log.z * m) } else this :=(1, 0, 0, 0) } def log(result: Vector3): Unit = { val xyz = Math.sqrt(x*x + y*y + z*z) //      if (xyz > 0.0001f) { val m = (Math.acos(w).toFloat * 2f) / xyz result.set(m * x, m * y, m * z) } else { result.set(0f, 0f, 0f) } } } 

抽象的観点からの飛行機の飛行とは何ですか？ 二次微分方程式系を解くことにより！ 航空機の状態を設定します。

 trait TimeStamp{ def time: Long //   } trait VisibleState extends TimeStamp{ def position: Position def speed: Derivative } 

航空機の実際の状態のクラスは、モデルの機能とその物理的な開発の程度に依存しますが、このインターフェイスを実装します。これは、画面に航空機を描画するために必要です。

タスクは2つの独立した部分に分けることができます

1）所定の状態で航空機に作用する力の計算
2）ディファーの数値解

最も興味深いので、最初の部分は記事の最後に残します。

グラフィック部分については、インターフェイスを記述します。

 abstract class Avatar(val player: Player, val model: Airplane) { def getState(time: Long): VisibleState } 

実装はしませんが、本質は単純です。時間t1、t2、t3など、t（n + 1）> t（n）の状態のシーケンスが内部に格納されます。 必要に応じて状態が完成します; getStateメソッドでは、最も近い2つが補間されます。 したがって、たとえば、物理を毎秒10回考慮し、同時に60 fpsで滑らかな動きを観察することができます。

次のインターフェイスを作成しましょう。

 trait Airplane{ def mass: Float def I: Matrix4 private val tempM = new Matrix4() private val tempQ = new QuaternionS() def getInertiaTensor(orientation: Quaternion): Matrix4 = { tempM.set(orientation) tempM.mul(I) tempQ := orientation tempQ.conjugate() tempM.rotate(tempQ) tempM } def getForceAndMoment(state: State, input: PlaneInput): Derivative } 

角運動量はL = I * wであり、Lとw（角速度）はベクトルとして変換されます。 したがって、変換L '= qL、w' = qwでは、次のようになります。
L '= I' * w '
qL = I '* qw
L = q ^（-1）* I '* q * w

I = q ^（-1）* I '* q、またはI' = q * I * q ^（-1）を取得します。

変換w '= position.rot * wは、角速度をローカル座標系からグローバル座標系に変換します。

getForceAndMomentメソッドについては後述し、航空機に作用する力とトルクを計算します。

加速度で移動および回転するモデルの動きを正確に計算する方法がよくわからないため、20msの固定ステップで最も簡単な方法を選択しました。

 class StateGenerator { //      ,           // , ,     private val inversed = new Matrix4() private val omega = new Vector3f() private val iOm = new Vector3f() private val angularAcceleration = new Vector3f() private val inv = new QuaternionS() def nextState(now: State, timeStep: Long, model: Airplane, planeInput: PlaneInput): State = { assert(timeStep > 0) val dt = timeStep * 0.001f //   val next = new State(now.time + timeStep) next.position := now.position next.speed := now.speed //      val forces = model.getForceAndMoment(now, planeInput) //linear next.speed.linear += forces.linear * (dt / model.mass) next.position.pos += (now.speed.linear + next.speed.linear) * (0.5f * dt) //angular val I = model.getInertiaTensor(now.position.rot) inversed.set(I).inv() omega := now.speed.angular iOm := omega iOm.mul(I) angularAcceleration.setCross(iOm, omega) angularAcceleration += forces.angular angularAcceleration.mul(inversed) next.speed.angular.madd(angularAcceleration, dt) next.speed.angular *= 1f - dt * 0.1f //.    .     - 10 . val angSp = (next.speed.angular + now.speed.angular) / 2 next.position.rot := new QuaternionS(angSp * dt) * now.position.rot next //  Scala    return } } 

回転の詳細については、この記事habrahabr.ru/post/264099を参照してください 。 正直なところ、私はテンソルのファンではなく、角加速度を得るためにそこからベクトル形式の式を取りました。 計算は、世界の座標系で行われます。 ところで、外力がオフになったとき、私はジャニベコフ効果に非常に似た動きを観察することができました。

飛行機に作用する力

飛行機を制御する必要があります：

 trait PlaneInput extends TimeStamp { def yaw: Float // right is positive def pitch: Float // up is positive def roll: Float // clockwise is positive def engineTrust: Float } 

値は、間隔[-1、1]、エンジン推力0〜1に切り捨てられます。

さて、最も重要な部分に移りましょう-私たちは力を見つけます。 ここでは、Terra incognitaについて少し説明しますが、空気力学に関する私の知識は非常に表面的なものです。 エラーや不正確な可能性があります。

最初に頭に浮かぶのはリフトです。 ゴミを作成しないために、インターネット上で、攻撃角度に応じた揚力係数のグラフを含む航空プロファイルのディレクトリが見つかりました。 本質は非常に単純であることが判明しました-y（揚力係数）は臨界角までかなり直線的に増加し、約1に達し、その後流れは翼を打ち破り、揚力は減少し始めます。 また、抽象翼の係数グラフは英語版ウィキペディアで見つけることができます：

ここに私を待っているレーキがあります-ドラッグに関する別のウィキペディアの記事を読むと、何らかの誘導ドラッグがあることがわかります。 驚いたことに、揚力は速度の方向に垂直な方向にあり、翼の表面に垂直ではないと考えられます（私が思ったように）。 翼の上下の空気圧の差は依然として翼の表面に垂直な力をもたらすため、この力の動きとは反対の方向への投影はゼロではありません。 私が正しく理解していれば、これは誘導力です。
しかし、ありません。 以下のコメントを参照してください。 さらにテキストとコードは変更されません。

揚力が航空機の基準枠内で上方に向けられていると仮定した場合、誘導力は必要ないように思われます-すでに考慮されています。 軸の向きはopenGLと同じです：

牛-右
オイ-アップ
オズ-戻る

 //   def cy(angle: Float) = { val pi = 3.1415928535f if (angle > maxAngle) (pi - angle) / (pi - maxAngle) else if (angle < -maxAngle) -(pi + angle) / (pi - maxAngle) else angle / maxAngle } def liftingForce(speed2: Float, angle: Float, airDensity: Float): Float = cy(angle) * airDensity * speed2 * 0.5f * planesS //      . val forceLocal = new Vector3f val speedLocal = state.position.rot^(-1) * state.speed.linear val sp2 = state.speed.linear.length2 val airDensity = WorldParams.atmosphere.aitDensity(state.position.pos) val attackAngle = -Math.asin(speedLocal.y / (speedLocal.length + 0.1f)).toFloat val steeringAngle = Math.asin(speedLocal.x / (speedLocal.length + 0.1f)).toFloat //       ,      ( 0.1 /)   ,   ""  .  -    . forceLocal.y += liftingForce(sp2, attackAngle, airDensity) forceLocal.x -= steeringForce(sp2, steeringAngle, airDensity) forceLocal += speedLocal/(0.1f + speedLocal.length) * -dragForce(sp2, attackAngle, steeringAngle, airDensity) 

揚力に加えて、空気抵抗力dragForceと、飛行機が少し横に飛んでいる場合に発生する力dragForceが必要になります。

空力に関する十分な知識がありません。 主な目標は、式の単純さと、可能であれば、攻撃と滑空の飛行角度に対する航空機の適切な動作です。 0.5f-式の除数2の結果。 0.1f-フィット係数の結果。

 private def steeringForce(speed2: Float, angle: Float, airDensity: Float): Float = angle * airDensity * speed2 * 0.5f * planeSVert private def dragForce(speed2: Float, attackAngle: Float, steeringAngle: Float, airDensity: Float): Float = { speed2 * (0.5f * airDensity * crossSectionalArea + Math.abs(attackAngle) * 0.1f * planesS + Math.abs(steeringAngle) * 0.1f * planeSVert) } 

モータートラクションを追加

このモデルは可能な限りシンプルです。プロペラのステップは不要で、エンジンがすべての力を航空機の加速に費やします。 慣性モーメントにボーナスはありません。回転数を変更してもトルクはありません。 ただし、ターンもありません。 電力=電力*速度。 飛行機がロケットのように離陸できないように、最大​​力を制限します（最小速度制限を使用）。

 val speed = Math.max(minSpeed, -speedLocal.z) forceLocal.z -= input.engineTrust * maxPower / speed 

運営管理

興味深い点があります-ネジによって分散された空気は尾の制御面に直接行き、原則として、飛行機は滑走路上であっても尾によってわずかに制御されます。 さらに、空気抵抗は、プロペラの反対の回転方向に飛行機を回転させようとしています。 そして、ヒープに-エンジンには慣性モーメントがあり、回転速度を増減しながら、航空機は少し「ねじれ」ます。 私はこれをすべて無視しました...

翼に関しては、おなじみの要素が速度と空気密度の2乗に現れます。

 val mul = spLocal.z * spLocal.z * airDensity result.angular.x = (input.pitch + 0.3f) * mul * 2f result.angular.y = -input.yaw * mul * 1f result.angular.z = -input.roll * mul * 5 

ピッチに対称性はありません。プレーン（およびパイロット）は、負のオーバーロードよりも正のオーバーロードにはるかに耐えることができます。 さらに、航空機自体は安定しており（Su-47でない場合）、「前方位置」の位置に戻る傾向があります。

 result.angular.x -= attackAngle * airDensity * mul * 5f result.angular.y -= steeringAngle * airDensity * mul * 5f 

何も忘れていませんか？

行動がより面白くなる別の力があります。 航空機を正面または背面から見ると、翼がラテン文字Vでわずかに上方に曲がっていることがわかります。これは、飛行の安定性への懸念によって決まります。 飛行機が真っ直ぐに飛んでいないが、わずかに横に動いている場合、左右の持ち上げ力が異なり、回転し始めます。

 result.angular.z += forceLocal.y * steeringAngle * 1f 

forceLocal.y-揚力

回転に摩擦を追加します。

私の美意識が抗議する何かが起こりましたが、そうでなければモデルは非常に複雑になります。 強度を追加する前に、私はまだそれを正当化しようとします。 たとえば、左へのロールによってまっすぐな飛行面が回転すると、左翼の迎え角が上がり、右翼が逆になり、この効果により回転が遅くなります。 他の軸に沿っておそらく似たようなものがあります（StateGeneratorクラスでは、計算回路の安定性のために回転中の摩擦はほとんどありませんでしたが、ここでは、平面が振り子のようにならないようにします）：

 result.angular.mulAdd(spAngLocal, -airDensity * 5000f) 

グローバル座標系に変換します：

 result.angular.mul(state.position.rot) forceLocal.mul(state.position.rot) result.linear := forceLocal 

ご注意

単位系-メートル、キログラム、秒。 「調整係数」は理由のために与えられます-私はI-16パラメータのためにそれらを選択しようとしました。 質量1400、出力750hp、または（750 * 735.5）ワット。 （概算によると）慣性モーメントはOX、OYに沿って5000で、OZに沿ってはるかに小さくなります（航空機の胴体に集中する主質量など、かなり長い）。
Imp5は、より正確なデータを報告しました。それぞれ、軸に沿った慣性モーメント2440、5520、3080が 「前方」、「上方」、「右側」です。

この物理モデルは、航空機の回転を考慮に入れていないため、コルク抜きに陥ることはありません。 将来的には、各翼にいくつかのポイントを取り、各ポイントについて、攻撃角度と空気に対する速度を個別に計算する予定です。 ウィングピースのパラメーターを変更するようにテールとエルロンを制御します。 おそらく、航空機の回転は、空気抵抗のために正直に消えます。

航空機の強度と動きを計算するコードは、いつでもより深刻なものに置き換えることができます。

PS第二次世界大戦の航空機に関するシミュレーターの国内の開発者に帽子を脱ぎます。そこから、昔から航空への興味が始まりました。

物理学を自分で書いて、彼らがどのようなタイタニックな仕事をしたかを理解するのは価値がありました。 たとえば、縦方向に配置されたエンジンの回転は、一方向に曲がると航空機の機首が上下に移動するという事実につながります。 他の多くの人と同じように、この効果はありません。 一方では、些細なことですが、そのような些細なことから、各モデルに固有の動作が形成されます。