1. 從“角度”切入飛行姿態調整的數學本質最近在準備一個數學建模競賽題目是關于飛機飛行的。很多同學一看到“調整飛行角度”這個表述第一反應可能就是去翻空氣動力學教材找升力公式、阻力公式然后開始列微分方程。這當然沒錯但往往容易陷入復雜的物理細節而忽略了數學建模的核心——用合適的數學模型抓住問題的關鍵矛盾并給出可執行的解決方案。“調整飛行角度使飛機順利飛行”這個命題聽起來像是一個控制問題或者優化問題。它的核心是什么我認為是在動態變化的環境約束下如氣流、任務要求尋找一個最優或可行的飛行姿態序列使得飛機能夠穩定、高效、安全地從A點到達B點。這里的“順利”可以量化成多種目標可能是能耗最低、時間最短、顛簸最小或者是滿足特定航跡要求。那么“飛行角度”具體指什么在航空領域這通常指的是三個核心歐拉角俯仰角Pitch、滾轉角Roll、偏航角Yaw。俯仰角控制飛機爬升或下降滾轉角控制飛機轉向傾斜偏航角控制機頭左右指向。調整它們就相當于駕駛員在操縱桿和方向舵上的輸入。我們的數學模型就是要代替駕駛員計算出在什么時間點給出多大的“角度調整指令”。所以這個建模問題可以剝離掉復雜的飛機外形和發動機細節抽象為一個剛體運動學與動力學控制問題。我們需要建立描述飛機質心移動和姿態變化的微分方程運動學以及描述外力重力、氣動力、推力如何影響這種變化的方程動力學。然后“調整角度”就轉化為對這些微分方程施加控制輸入即改變歐拉角或其變化率并求解控制律的問題。2. 模型基石如何用數學描述飛機的“一舉一動”要調整角度首先得知道角度是如何影響飛機運動的。這就必須建立一個簡化的飛機運動模型。對于數學建模競賽而言使用完整的六自由度6-DOF模型可能過于復雜我們可以從更基礎的質點模型和剛體二自由度模型入手這已經能解決很多“順利飛行”的問題了。2.1 核心坐標系與歐拉角首先必須定義清楚坐標系這是所有計算的基準。地面坐標系慣性系固定在地面用于描述飛機的位置X, Y, Z和航跡。機體坐標系固定在飛機上原點在質心。X軸指向機頭Y軸指向右翼Z軸垂直向下按右手定則。我們關心的俯仰角θ、滾轉角φ、偏航角ψ就是描述機體坐標系相對于地面坐標系的姿態。歐拉角變換是有順序的通常按照“偏航-俯仰-滾轉”Z-Y-X的順序進行。這意味著我們先讓飛機繞地面系的Z軸轉ψ角偏航再繞新Y軸轉θ角俯仰最后繞新X軸轉φ角滾轉。這個順序很重要因為歐拉角變換不可交換。2.2 質點運動學模型先把飛機看成一個點如果我們只關心飛機飛行的路徑而不關心其姿態細節可以將其視為一個有質量的點。此時“調整角度”簡化為調整其速度矢量的方向。模型非常簡單位置更新[X_{k1}, Y_{k1}, Z_{k1}]^T [X_k, Y_k, Z_k]^T V_k * Δt * [cos(γ)cos(χ), cos(γ)sin(χ), -sin(γ)]^T這里V是空速相對于空氣的速度大小γ是航跡傾斜角爬升角χ是航跡方位角。通過調整γ和χ我們就間接調整了飛行路徑。這個模型常用于航路規劃和初步的軌跡優化。2.3 剛體二自由度模型引入姿態動力學要真正理解“調整俯仰角/滾轉角”就必須引入轉動動力學。一個非常經典且實用的簡化是縱向模型和橫側向模型的解耦分析。縱向運動模型主要在垂直平面關注俯仰角θ和空速V、高度H的關系。核心方程涉及升力、阻力、重力和推力在縱向的平衡。一個典型的線性化小擾動縱向方程可以寫成狀態空間形式? A*x B*u其中狀態變量x可能包含空速變化量u 攻角變化量α 俯仰角速率q 俯仰角θ。輸入u就是升降舵偏角δ_e它是駕駛員拉桿操縱直接導致俯仰角變化。在這個模型里我們的目標就是設計控制律通過δ_e來控制θ進而控制爬升率和高度。橫側向運動模型主要在水平平面關注滾轉角φ和偏航角ψ。核心是協調轉彎。要使飛機順利轉彎而不側滑需要滿足“協調轉彎”條件滾轉角φ產生的升力水平分量提供向心力同時用方向舵抵消偏航阻尼。轉彎半徑R與滾轉角φ、空速V的關系可近似為R ≈ V^2 / (g * tan(φ))。這里調整滾轉角φ就直接決定了轉彎的急緩。注意在實際建模中我們很少直接控制歐拉角θ和φ因為它們是“姿態”。我們控制的是舵面升降舵、副翼、方向舵。舵面偏轉產生氣動力矩力矩引起角加速度角加速度積分得到角速度角速度積分才得到歐拉角。這是一個二階關系。因此更精確的模型會以角速度p, q, r作為狀態變量。2.4 力與力矩的來源氣動系數是關鍵上述動力學方程中的力升力L、阻力D和力矩俯仰力矩M、滾轉力矩L_m、偏航力矩N如何計算這離不開氣動系數。它們是攻角α、側滑角β、舵面偏角、以及無量綱角速度如p*b/(2V),q*c/(2V)的函數。通常采用線性化模型C_L C_{L0} C_{Lα}*α C_{Lδ_e}*δ_e ...C_m C_{m0} C_{mα}*α C_{mδ_e}*δ_e C_{mq}*(q*c/(2V)) ...這些系數C_{Lα},C_{mα}等是飛機的固有屬性可以從風洞實驗或CFD計算中獲得在建模比賽中往往作為已知參數給出。我們的調整策略最終都體現在如何影響這些氣動系數上。3. “順利飛行”的目標函數與約束條件建模定義了模型接下來要定義什么是“順利”。這需要我們將模糊的“順利”轉化為精確的數學目標函數和約束條件。3.1 常見的目標函數我們想優化什么最短時間J ∫ dt從初始時刻到終端時刻的積分最小。適用于緊急任務或競速。最省燃料J ∫ ?_fuel dt燃料消耗率對時間的積分最小。這是民航領域最核心的經濟性指標。燃料消耗率通常與發動機推力或功率有關。最舒適顛簸最小J ∫ (a_z^2) dt將垂直加速度的平方對時間積分最小化。這模擬了乘客感受到的顛簸程度。軌跡跟蹤誤差最小J ∫ [(X-X_{ref})^2 (Y-Y_{ref})^2 (Z-Z_{ref})^2] dt使實際航跡盡可能貼近一條給定的理想航跡。控制能量最小J ∫ (δ_e^2 δ_a^2 δ_r^2) dt最小化舵面偏轉的幅度這有助于減少舵機磨損和能耗。在實際建模中可能需要組合多個目標例如“在保證不超過最大過載的前提下盡可能省油地飛抵目標”。3.2 必須遵守的約束條件飛行的“紅線”狀態約束過載約束法向過載n_z必須保持在安全范圍內如民航客機通常為-1g到2.5gn_z L / (mg)。大過載不僅乘客不適還會損傷飛機結構。姿態角約束滾轉角φ通常有絕對值限制如±30°以防失速或失控。俯仰角θ也有安全包線。空速約束速度必須介于失速速度V_stall和最大允許速度V_max之間。高度約束必須高于地面障礙物同時可能低于某個空域上限。控制輸入約束舵面偏角限制|δ_e| ≤ δ_{e,max}|δ_a| ≤ δ_{a,max}|δ_r| ≤ δ_{r,max}。舵片只能偏轉那么多。舵面偏轉速率限制舵面偏轉的速度也是有限的|δ?| ≤ δ?_max。這限制了姿態變化的快慢。終端約束在最終時刻t_f飛機必須到達指定位置(X_f, Y_f, Z_f)并且可能要求速度矢量方向與跑道對齊著陸場景或者速度、高度達到特定值。3.3 動態環境約束風與紊流的影響“順利飛行”幾乎總是在有風的環境中進行的。風尤其是突風是破壞飛行順利的主要因素。我們需要將風矢量W [W_x, W_y, W_z]^T納入模型。地速與空速飛機相對于地面的速度地速V_g等于空速V_a加上風速WV_g V_a W。導航和軌跡跟蹤是基于地速的而氣動力升力、阻力的計算依賴于空速。這是一個關鍵點你為了對抗側風保持航跡就需要讓機頭偏向上風方向形成偏航角這會產生側滑為了補償側滑又可能需要壓一點坡度滾轉角。這個調整過程本身就是一個經典的“航跡修正”控制問題。紊流建模對于高保真度的“順利”建模可能需要考慮紊流。紊流可以建模為在平均風基礎上的隨機擾動例如使用Dryden或Von Karman紊流譜模型生成隨機的時間序列疊加到風速上。這會使得控制問題變成一個隨機最優控制或魯棒控制問題。4. 求解策略從經典控制到最優控制模型和目標都有了怎么解如何算出那一系列“該調整的角度”這里有幾個層次的解決方法。4.1 層次一經典PID控制與自動駕駛儀模擬這是最直觀、最工程化的方法。我們為每一個需要控制的“環”設計一個PID控制器。高度保持環外環。輸入是高度誤差輸出是目標俯仰角。內環是俯仰角保持環輸入是俯仰角誤差輸出是升降舵指令。航向保持環外環。輸入是航向誤差輸出是目標滾轉角。內環是滾轉角保持環輸入是滾轉角誤差輸出是副翼指令。空速保持環通過油門控制。這樣當我們想讓飛機爬升到新高度時就給高度環一個目標值它會產生一個目標俯仰角比如5°上仰俯仰角環會驅動升降舵使飛機抬頭直到實際俯仰角達到5°。在這個過程中空速可能會變化空速環會自動調整油門。這種層級控制結構就是現代自動駕駛儀的核心邏輯。在建模中我們可以用Simulink等工具搭建這樣的控制框圖通過調整PID參數比例、積分、微分系數來獲得“順利”的響應曲線快速、無超調、穩態誤差小。4.2 層次二基于狀態空間的線性二次型調節器LQR當我們的模型是線性狀態空間形式? Ax Bu且目標函數是狀態和控制輸入的二次型J ∫ (x^T Q x u^T R u) dt時LQR可以提供最優的全狀態反饋控制律u -K x。矩陣K可以通過求解Riccati方程得到。這非常適合“穩定飛行”場景。例如我們設計一個LQR控制器讓飛機在平飛狀態某個平衡點附近能夠自動抵抗小的擾動如陣風迅速恢復到原狀態。Q和R矩陣的權重選擇決定了控制器是更激進地消除狀態偏差大Q還是更溫柔地使用舵面大R。通過調整Q和R我們可以在“響應速度”和“乘坐舒適度/能耗”之間取得平衡。4.3 層次三軌跡優化與最優控制這是最通用、最強大的方法直接求解我們之前定義的最優控制問題。問題可以表述為 尋找控制輸入歷史u(t)和對應的狀態軌跡x(t)在滿足動力學方程? f(x, u, t)和路徑約束g(x, u, t) ≤ 0的前提下最小化目標函數J φ(x(t_f)) ∫ L(x, u, t) dt。求解方法主要有兩種直接法將連續時間問題離散化。把整個飛行時間分成N段控制輸入u在每個時間段內假設為常數或簡單函數如線性狀態x通過數值積分如龍格-庫塔法向前傳播。這樣最優控制問題被轉化為一個大規模的非線性規劃問題NLP可以使用IPOPT、SNOPT等求解器求解。工具推薦MATLAB的fmincon優化函數或者更專業的工具包如GPOPS-II、CasADi。這種方法直觀能直接處理各種復雜約束。間接法龐特里亞金極小值原理引入協態變量拉格朗日乘子將最優控制問題轉化為一個兩點邊值問題。這需要推導一階最優性條件哈密頓系統求解通常更困難但對理解問題本質很有幫助。在建模中的應用假設題目要求飛機從機場A爬升并巡航到機場B要求總油耗最小。我們可以將整個過程離散化為多個階段爬升、巡航、下降建立包含發動機油耗模型、氣動模型的完整動力學方程以油門和舵面為控制量以位置、速度、姿態為狀態量以過載、空速、姿態角為約束以總油耗為目標用直接法進行求解。最終得到的結果就是一套最優的飛行剖面高度-距離曲線和對應的控制歷史油門、舵面隨時間變化曲線從中我們可以解讀出每個時刻“應該調整的角度”是多少。4.4 層次四模型預測控制MPC——應對不確定性的利器MPC結合了最優控制和反饋控制的優點。它在每個控制周期內根據當前飛機狀態預測未來一段時間預測時域的狀態演化。求解一個有限時域的最優控制問題通常用直接法得到未來一段控制時域內的最優控制序列。只實施該序列的第一個控制指令。到下一個周期重復以上步驟。MPC的強大之處在于其處理約束和應對擾動的能力。風擾動了飛機沒關系下一個控制周期MPC會根據新的狀態重新規劃自動調整角度來抵消風的影響。它實時地、滾動地求解“調整角度”的問題是實現復雜環境下“順利飛行”的先進策略。計算復雜度是其主要挑戰但隨著硬件發展其在航空中的應用越來越廣。5. 建模實戰一個簡化案例的完整推演讓我們用一個極度簡化的二維平面飛行案例把上述概念串起來。任務飛機在恒定高度飛行初始位置(0,0)初始航向朝東0度。目標飛到(10000, 0)點并讓機頭重新朝東。要求飛行時間盡量短但轉彎時滾轉角不能超過30度且法向過載不能超過1.5g。假設空速V恒定。步驟1建立模型我們使用質點模型并忽略高度變化。狀態變量為x東向位置y北向位置χ航跡方位角。控制變量為φ滾轉角假設協調轉彎側滑角為0。 動力學方程? V * cos(χ)? V * sin(χ)χ? (g * tan(φ)) / V由協調轉彎公式R V^2/(g*tanφ)和χ? V/R推導而來步驟2定義目標與約束目標函數最小化終端時間t_f。這是一個最小時間問題。 路徑約束|φ| ≤ φ_max 30°|n_z| ≤ 1.5。在協調轉彎中法向過載n_z 1/cos(φ)。因此|φ| ≤ arccos(1/1.5) ≈ 48.2°。顯然30°的滾轉角約束更嚴格成為主導約束。 終端約束x(t_f)10000,y(t_f)0,χ(t_f)0。步驟3問題分析與求解這是一個典型的Bang-Bang控制問題。為了最快到達目標點飛機應該以最大允許滾轉角±30°進行轉彎。觀察任務需要先轉彎調整航向角χ再直飛最后可能再微調航向。 實際上由于起點航向0°已經指向目標點方向似乎不需要轉彎但目標點(10000,0)在正東起點(0,0)航向0°這已經是直線。這個例子太簡單。我們修改一下讓目標點在(5000, 5000)即東北方向45°。那么最優策略很可能是以最大正滾轉角右壓坡度導致左轉彎轉彎直到航向角χ指向目標點方向45°然后改平φ0直線飛行在接近目標點時可能需要一個反向的轉彎來精確對準最終航向0°。這個“接近目標點”的時刻需要優化。步驟4數值求解直接法思路將時間t離散為N個點t_0, t_1, ..., t_Nt_N t_f。控制量φ在每個區間為常數設為優化變量φ_0, φ_1, ..., φ_{N-1}。狀態量x, y, χ通過數值積分如歐拉法從初始狀態開始傳播χ_{k1} χ_k (g * tan(φ_k) / V) * Δtx_{k1} x_k V * cos(χ_k) * Δty_{k1} y_k V * sin(χ_k) * Δt優化變量φ_0, ..., φ_{N-1} 以及總時間t_f或等價的Δt。目標最小化t_f。約束路徑約束|φ_k| ≤ 30°。終端約束x_N ≈ 5000,y_N ≈ 5000,χ_N ≈ 0可以設為軟約束即允許小誤差在目標函數中加懲罰項。使用非線性規劃求解器如fmincon求解。最終得到的φ_k序列就是“調整飛行角度”的最優策略。步驟5結果解讀求解器會給出一個最優的φ(t)曲線。我們很可能會看到一條先為30°后為0°最后為-30°的曲線。對應的航跡是一條Dubins路徑直線-圓弧-直線的變體。這就是數學建模給出的“順利飛行”此處為最快飛行的角度調整方案。6. 從模型到代碼實現與驗證的關鍵環節建好模型、設計好求解策略后必須通過編程實現和數值仿真來驗證方案的有效性。這里有幾個容易踩坑的地方。6.1 數值積分器的選擇在仿真動力學方程? f(x, u)時積分步長和方法至關重要。歐拉法最簡單x_{k1} x_k f(x_k, u_k)*dt。但精度低穩定性差對于剛性問題或大步長容易發散。僅適用于最簡單的教學模型和極小的dt。龍格-庫塔法RK4最常用的中級方法。精度為四階在大多數飛行仿真中表現良好。MATLAB中的ode45變步長RK是首選。在自行編寫時RK4是一個可靠的基準。對于剛性問題如果模型包含快慢相差很大的動態如剛體姿態動力學與發動機動力學耦合可能需要使用剛性求解器如MATLAB的ode15s。實操心得始終用ode45開始你的仿真。如果發現計算異常緩慢或報錯“剛度檢測”再換用ode15s。仿真步長dt應小于系統最快動態時間常數的1/10。例如俯仰角速率回路可能非常快時間常數0.1秒那么dt最好小于0.01秒。6.2 單位制的統一與量綱檢查這是新手最容易出錯的地方錯誤會導致結果完全失真甚至數值溢出。堅持使用國際單位制長度用米(m)速度用米/秒(m/s)角度用弧度(rad)質量用千克(kg)力用牛頓(N)。特別注意很多氣動數據表提供的角度可能是度(°)在代入三角函數sin(),cos(),tan()前務必轉換為弧度。在MATLAB中sind(),cosd()等函數可以直接接受度但為了統一我強烈建議在代碼開始就把所有輸入轉換為弧度。進行量綱檢查在寫出每一個方程后手動檢查等式兩邊的量綱是否一致。例如χ? g*tan(φ)/V左邊是弧度/秒(rad/s)右邊是(m/s2) / (m/s) 1/s量綱一致弧度是無量綱的。這一步能避免很多低級錯誤。6.3 控制器設計與調參的工程經驗無論是PID還是LQR調參都有技巧。PID調參遵循“先P后I最后D”的原則。比例P逐漸增大P直到系統出現等幅振蕩。此時增益記為K_u振蕩周期記為T_u。積分I根據齊格勒-尼科爾斯法則對于PI控制器K_p 0.45*K_u,T_i 0.85*T_uK_i K_p / T_i。然后微調。微分DD項能抑制超調但會放大噪聲。在實際飛行控制中對姿態角速率如q的微分反饋比直接對姿態角微分更常用這相當于增加了阻尼。LQR權重選擇Q和R矩陣通常取對角陣。一個實用的方法是Bryson法則將每個狀態變量和控制量的最大允許偏差的平方的倒數作為對應Q或R對角線元素的初始值。例如如果希望俯仰角θ的偏差不要超過0.1 rad那么Q_θ可以設為1/(0.1^2) 100。然后在此基礎上進行縮放和微調。6.4 可視化讓結果自己說話仿真結果必須通過圖表清晰呈現。必須有的圖三維航跡圖顯示飛機在空間中的飛行路徑。時間序列圖將關鍵狀態量高度、空速、俯仰角、滾轉角和控制量舵面偏角、油門隨時間的變化畫出來。相平面圖對于分析穩定性非常有用例如畫出俯仰角θ與俯仰角速率q的關系。進階分析蒙特卡洛仿真為了驗證控制系統的魯棒性可以在關鍵參數如氣動系數、風速上添加隨機擾動進行上百次仿真統計性能指標如著陸點散布、最大過載的分布情況。頻域分析對線性化模型進行波特圖、奈奎斯特圖分析評估系統的穩定裕度幅值裕度、相位裕度。7. 總結與延伸思考通過上面的討論我們可以看到“調整飛行角度使飛機順利飛行”這個看似具體的問題背后是一套完整的數學建模與控制理論體系。從最基本的運動學關系到復雜的非線性最優控制解決問題的工具層次豐富。在數學建模競賽中關鍵在于根據題目提供的具體信息數據粒度、問題焦點選擇合適的模型復雜度和求解方法。我個人在類似建模中的體會是不要一開始就追求最復雜的模型。從一個能反映核心物理的簡化模型如質點模型協調轉彎假設入手把優化問題的框架搭起來得到初步結果。這個結果往往已經能揭示最優策略的形態比如Bang-Bang控制。然后再考慮增加模型保真度如加入俯仰動力學、發動機模型看最優策略如何變化。這種由簡入繁的迭代過程既能保證項目進度又能加深對問題本質的理解。最后這個問題的外延非常廣。所謂的“順利”在今天可以引申為很多智能化的要求。例如在密集空域中的沖突解脫本質上就是多架飛機各自調整飛行角度航向、高度在滿足安全間隔約束下的協同優化問題。再比如考慮突風擾動的魯棒軌跡規劃或者考慮發動機故障的應急著陸軌跡生成都是這一核心問題的變體和深化。掌握了這套從建模到求解再到驗證的完整流程就擁有了分析和解決一大類運動體控制與優化問題的鑰匙。