簡介本資源為2022年全國大學(xué)生數(shù)學(xué)建模競賽本科組B題的完整論文PDF聚焦無人機(jī)純方位無源定位問題適合備戰(zhàn)數(shù)模競賽的本科生及指導(dǎo)教師參考。論文圍繞圓形與錐形兩種編隊(duì)系統(tǒng)建立被動(dòng)接收信號(hào)無人機(jī)的定位模型與位置調(diào)整模型運(yùn)用二維直角坐標(biāo)系、余弦定理與“三圓共點(diǎn)”原理推導(dǎo)發(fā)射信號(hào)無人機(jī)數(shù)量需求并給出基于0.01米步長的遍歷調(diào)整方案全部求解過程均以MATLAB編程實(shí)現(xiàn)。壓縮包內(nèi)僅含1個(gè)PDF文件約2.63MB完整呈現(xiàn)摘要、問題重述、模型假設(shè)、符號(hào)說明及建模求解全過程便于讀者直接研讀建模思路與代碼邏輯。該論文曾獲河北省一等獎(jiǎng)目前已有8039人學(xué)習(xí)下載可為參賽者提供從模型構(gòu)建到編程求解的完整參考范例。1. 從三組方位角反推坐標(biāo)這套無人機(jī)純方位無源定位方案到底能跑通什么編隊(duì)飛行里最怕的不是風(fēng)是電磁信號(hào)被干擾后整個(gè)隊(duì)形散架。純方位無源定位的思路很直接無人機(jī)不主動(dòng)發(fā)波只靠接收其他已知位置無人機(jī)發(fā)來的方向信息反推自己在哪、該往哪挪。這套 2022 年全國大學(xué)生數(shù)學(xué)建模競賽本科組 B 題的完整解答核心就是把這個(gè)過程拆成可計(jì)算的模型——用二維直角坐標(biāo)系加余弦定理把三個(gè)方位角變成坐標(biāo)再用步長遍歷把偏差位置拉回理想位置。它適合兩類人正在備戰(zhàn)數(shù)模競賽、需要一套能直接跑 MATLAB 的完整代碼框架的選手以及做無人機(jī)編隊(duì)控制、想理解純方位定位最小可解條件的工程師。整套資料包含問題重述、模型建立、MATLAB 程序附錄和結(jié)果分析不是純理論推導(dǎo)附錄里四個(gè)完整程序覆蓋了從定位到調(diào)整的全流程。2. 定位模型怎么建從三個(gè)方位角到二維坐標(biāo)的完整推導(dǎo)鏈2.1 為什么選余弦定理而不是最小二乘純方位無源定位的常見做法有兩種一種是基于最小二乘的偽線性估計(jì)把非線性角度方程線性化后迭代求解另一種是直接利用幾何關(guān)系構(gòu)造解析解。這套方案選了后者原因很實(shí)際——圓形編隊(duì)里發(fā)射信號(hào)的無人機(jī)位置是精確已知的接收信號(hào)的無人機(jī)與它們之間的角度關(guān)系可以用余弦定理直接串起來不需要迭代也不存在初值敏感的問題。具體來說把第 i 架無人機(jī)的準(zhǔn)確位置固定在 x 軸上圓形編隊(duì)半徑記為 r相鄰無人機(jī)間隔 40 度。三架發(fā)射信號(hào)無人機(jī)分別記為 FY0i、FY0j、FY0k接收信號(hào)無人機(jī)為 FY0x。接收到的三個(gè)方位角 α1、α2、α3 分別對應(yīng)接收機(jī)與三架發(fā)射機(jī)連線之間的夾角。在三角形中反復(fù)使用余弦定理就能把三條邊的長度解出來再反解出接收機(jī)的坐標(biāo)。這里有個(gè)容易翻車的地方當(dāng)接收機(jī)位于發(fā)射機(jī)連線的延長線上時(shí)角度信息退化方程會(huì)出現(xiàn)多解甚至無解。這也是后面問題二里第二輪定位時(shí) 4、7、9、11、13、14 號(hào)無人機(jī)無法定位的根本原因。2.2 被動(dòng)接收信號(hào)無人機(jī)定位模型的 MATLAB 實(shí)現(xiàn)附錄一的程序分兩段第一段求解當(dāng)前位置坐標(biāo)第二段求解準(zhǔn)確位置坐標(biāo)。先看第一段的核心邏輯。% 被動(dòng)接收信號(hào)無人機(jī)定位模型的建立 % 求解當(dāng)前位置的坐標(biāo) alpha1 input(第一個(gè)角度); % 0k i alpha2 input(第二個(gè)角度); % 0k j alpha3 input(第三個(gè)角度); % j k i FY0i input(發(fā)射信號(hào)無人機(jī) 1); FY0j input(發(fā)射信號(hào)無人機(jī) 2); r input(圓形編隊(duì)半徑); % 計(jì)算兩架發(fā)射信號(hào)無人機(jī)之間的距離 a1 sqrt(2*r^2 - cosd(abs(FY0i - FY0j)*40)*2*r^2); syms k0 a2 a3 % 三個(gè)余弦定理方程三個(gè)未知量 eq1 (k0^2 a3^2 - r^2)/(2*k0*a3) - cosd(alpha1); eq2 (k0^2 a2^2 - r^2)/(2*k0*a2) - cosd(alpha2); eq3 (a2^2 a3^2 - a1^2)/(2*a2*a3) - cosd(alpha3); [X, Y, Z] solve(eq1, eq2, eq3, k0, a2, a3); K0 double(X); % 求得 k0 長度 % 計(jì)算接收機(jī)相對于 x 軸的偏轉(zhuǎn)角 k0i acosd(1 - 2*cosd(alpha1)^2); k01 k0i - 40*abs(i - 1); x real(K0*cosd(k01)) % 橫坐標(biāo) y real(K0*sind(k01)) % 縱坐標(biāo)這段代碼的關(guān)鍵參數(shù)有三個(gè)alpha1、alpha2、alpha3是接收到的三個(gè)方位角單位是度FY0i和FY0j是兩架發(fā)射信號(hào)無人機(jī)的編號(hào)用來計(jì)算它們之間的基準(zhǔn)距離r是圓形編隊(duì)半徑題目里是 100 米。syms聲明符號(hào)變量后solve一次性解出三個(gè)未知邊長再用double轉(zhuǎn)成數(shù)值。最后兩行把極坐標(biāo)轉(zhuǎn)成直角坐標(biāo)real是防止符號(hào)計(jì)算產(chǎn)生微小虛部。第二段程序反過來已知接收機(jī)編號(hào)和發(fā)射機(jī)編號(hào)直接算出三個(gè)方位角的準(zhǔn)確值% 求解無人機(jī)接收到的角度信息 FY00 0; FY0i input(發(fā)射信號(hào)無人機(jī) 1); FY0j input(發(fā)射信號(hào)無人機(jī) 2); FY0k input(接收信號(hào)無人機(jī)); r input(圓形編隊(duì)半徑); % 計(jì)算三條邊的長度 a1 sqrt(2*r^2 - cosd((FY0j - FY0i)*40)*2*r^2); a3 sqrt(2*r^2 - cosd((FY0k - FY0i)*40)*2*r^2); a2 sqrt(2*r^2 - cosd((FY0k - FY0j)*40)*2*r^2); % 用余弦定理求三個(gè)角 alpha1 acosd((r^2 a3^2 - r^2)/(2*r*a3)); alpha2 acosd((r^2 a2^2 - r^2)/(2*r*a2)); alpha3 acosd((a3^2 a2^2 - a1^2)/(2*a3*a2));注意alpha1和alpha2的公式里出現(xiàn)了r^2 a3^2 - r^2看起來像冗余實(shí)際上是因?yàn)閳A心無人機(jī) FY00 到圓周無人機(jī)的距離恒為 r代入余弦定理后 r 平方項(xiàng)抵消了一部分。這是圓形編隊(duì)的特殊性質(zhì)換成錐形編隊(duì)就不能這么簡化。提示運(yùn)行第一段程序前需要先確定接收機(jī)在圓周上的大致方位否則k01的計(jì)算會(huì)引入 40 度整數(shù)倍的偏移誤差。我一般會(huì)先用第二段程序算出準(zhǔn)確角度再拿第一段反解驗(yàn)證。3. 額外需要幾架發(fā)射機(jī)枚舉法確定最小可解條件3.1 為什么一架額外發(fā)射機(jī)就夠了問題一第二小問的核心結(jié)論是除 0 號(hào)和 1 號(hào)外只需要額外增加 1 架發(fā)射信號(hào)無人機(jī)就能實(shí)現(xiàn)有效定位。這個(gè)結(jié)論的推導(dǎo)依賴一個(gè)幾何事實(shí)——一個(gè)動(dòng)點(diǎn)到兩個(gè)定點(diǎn)的夾角為定值時(shí)動(dòng)點(diǎn)軌跡是上下對稱的兩段圓弧組成的封閉曲線。當(dāng)夾角為鈍角時(shí)軌跡像一只“眼睛”當(dāng)夾角為銳角時(shí)軌跡像一個(gè)“葫蘆”。三架發(fā)射機(jī)產(chǎn)生三條這樣的軌跡三條軌跡相交于一點(diǎn)該點(diǎn)就是接收機(jī)的實(shí)際位置。如果只有兩架發(fā)射機(jī)接收機(jī)只能得到一個(gè)方位角軌跡是一條曲線而非一個(gè)點(diǎn)無法唯一確定位置。所以額外增加一架是理論下限。但這里有個(gè)隱藏條件三條軌跡必須確實(shí)相交于同一點(diǎn)。當(dāng)接收機(jī)位于某兩架發(fā)射機(jī)連線的延長線上時(shí)對應(yīng)的軌跡退化成直線交點(diǎn)不唯一。這就是為什么在問題二的錐形編隊(duì)里第二輪定位時(shí)會(huì)有 6 架無人機(jī)無法定位。3.2 枚舉法求解額外發(fā)射機(jī)編號(hào)的 MATLAB 實(shí)現(xiàn)附錄二的程序用枚舉法遍歷所有可能的發(fā)射機(jī)組合找出使方位角偏差平方和最小的那一組。% 求出所有的可能性 FY00 0; FY0i 1; % 兩個(gè)已知的定點(diǎn) r input(圓形編隊(duì)半徑); % 遍歷所有可能的發(fā)射機(jī)組合預(yù)計(jì)算角度表 for FY0j 2:9 for FY0k 2:9 if FY0k FY0j continue; end a1 sqrt(2*r^2 - cosd((FY0j - FY0i)*40)*2*r^2); a3 sqrt(2*r^2 - cosd((FY0k - FY0i)*40)*2*r^2); a2 sqrt(2*r^2 - cosd((FY0k - FY0j)*40)*2*r^2); alpha1(FY0j, FY0k) acosd((r^2 a3^2 - r^2)/(2*r*a3)); alpha2(FY0j, FY0k) acosd((a3^2 a2^2 - a1^2)/(2*a3*a2)); alpha3(FY0j, FY0k) acosd((r^2 a2^2 - r^2)/(2*r*a2)); end end % 反解根據(jù)接收到的角度信息找最匹配的發(fā)射機(jī) beta1 input(第一個(gè)角度); beta2 input(第二個(gè)角度); beta3 input(第三個(gè)角度); FY0k input(接收信號(hào)的無人機(jī)編號(hào)); min_val 1000; mini 0; for i 2:9 if i FY0k continue; end m (beta1 - alpha1(i, FY0k))^2 ... (beta2 - alpha2(i, FY0k))^2 ... (beta3 - alpha3(i, FY0k))^2; if m min_val min_val m; mini i; end end mini % 輸出需要額外增加的發(fā)射信號(hào)無人機(jī)編號(hào)這段代碼的邏輯分兩步先建一張角度查找表把所有可能的發(fā)射機(jī)組合對應(yīng)的三個(gè)方位角算出來存進(jìn)alpha1、alpha2、alpha3三個(gè)矩陣然后對于給定的接收角度beta1、beta2、beta3遍歷所有候選發(fā)射機(jī)找平方和最小的那個(gè)。參數(shù)說明FY0i固定為 1因?yàn)轭}目規(guī)定 0 號(hào)和 1 號(hào)已經(jīng)發(fā)射信號(hào)FY0j和FY0k從 2 到 9 遍歷跳過相同編號(hào)m是三個(gè)角度偏差的平方和越小說明匹配度越高。最終輸出的mini就是需要額外增加的發(fā)射機(jī)編號(hào)。注意枚舉法在 8 架無人機(jī)里找 1 架額外發(fā)射機(jī)計(jì)算量是 C(7,1)7 種組合非常小。但如果問題擴(kuò)展到幾十架無人機(jī)這種暴力枚舉就不合適了需要換成基于幾何約束的剪枝搜索。4. 位置調(diào)整怎么收斂步長遍歷與極坐標(biāo)更新4.1 向北還是向南調(diào)整方向的分類邏輯問題一第三小問要求把有偏差的無人機(jī)調(diào)整到理想位置。方案設(shè)定步長為 0.01 米向四周全方位遍歷當(dāng)與準(zhǔn)確位置相差小于 0.01 米時(shí)認(rèn)為到達(dá)理想位置。調(diào)整方向分成向北和向南兩類本質(zhì)上是極坐標(biāo)下角度增量的正負(fù)。具體做法無人機(jī)當(dāng)前位置記為 B向北移動(dòng)一個(gè)微小步長后到達(dá) C。已知 B 點(diǎn)的方位信息可以算出 BC 的長度和角度變化再用余弦定理求出 C 點(diǎn)相對于圓心的極坐標(biāo)。重復(fù)這個(gè)過程直到偏差小于閾值。這里有個(gè)工程上很實(shí)用的簡化調(diào)整時(shí)間不計(jì)意味著無人機(jī)知道自己需要調(diào)整位置的瞬間就已經(jīng)到位。這在實(shí)際系統(tǒng)里當(dāng)然不成立但在建模階段可以大幅簡化動(dòng)力學(xué)部分把注意力集中在定位精度上。4.2 圓形編隊(duì)位置調(diào)整的 MATLAB 實(shí)現(xiàn)附錄三的程序用 while 循環(huán)實(shí)現(xiàn)迭代調(diào)整。% 從無人機(jī)接收的方位信息去定位無人機(jī)找到無人機(jī)的位置 FY00 0; FY01 1; njl 0.01; % 最小距離閾值 r 100; % 圓形編隊(duì)半徑 bzrow 100; % 準(zhǔn)確位置 row bztheta [0 40 80 120 160 200 240 280 320]; % 準(zhǔn)確位置 theta FY0irow(2) input(FY02 的初始位置(row):); FY0itheta(2) input(FY02 的初始位置(theta):); jl 100; while(jl njl) mjl 10000; for j 0:359 % 遍歷所有可能的角度找距離最小的方向 % 此處省略具體角度計(jì)算核心是極坐標(biāo)更新 % ... end % 更新位置縮小 jl jl mjl; end這段代碼的框架是外層 while 控制迭代直到距離小于閾值內(nèi)層 for 遍歷 0 到 359 度所有方向找使距離最小的調(diào)整角度。njl是收斂閾值設(shè)成 0.01 米r是編隊(duì)半徑bztheta是 9 架圓周無人機(jī)的理想角度分布間隔 40 度。實(shí)際運(yùn)行時(shí)每次迭代需要根據(jù)當(dāng)前極坐標(biāo)和方位角信息用余弦定理算出移動(dòng)一個(gè)步長后的新極坐標(biāo)再與理想位置比較。如果偏差沒有縮小說明方向選反了需要反向遍歷。提示步長 0.01 米在 100 米半徑的編隊(duì)里相當(dāng)于萬分之一精度迭代次數(shù)可能上百次。如果只做仿真驗(yàn)證可以先用 0.1 米步長快速收斂最后再用 0.01 米精調(diào)。5. 錐形編隊(duì)怎么套用從圓形到錐形的模型遷移與邊界5.1 錐形編隊(duì)的幾何差異與定位策略問題二的錐形編隊(duì)和圓形編隊(duì)本質(zhì)區(qū)別在于圓周上相鄰無人機(jī)等距但錐形編隊(duì)里相鄰三架無人機(jī)構(gòu)成等邊三角形直線上的間距相等。這意味著圓形編隊(duì)里“相鄰間隔 40 度”的簡化不再成立需要重新建立坐標(biāo)系。方案假設(shè) 1 號(hào)和 5 號(hào)無人機(jī)位置無偏差先定位 2 號(hào)。具體做法將 1 號(hào)、2 號(hào)作為發(fā)射機(jī)5 號(hào)接收得到一個(gè)方向信息再將 2 號(hào)、5 號(hào)作為發(fā)射機(jī)1 號(hào)接收得到另一個(gè)方向信息。在三角形 1-2-5 中兩個(gè)角已知第三個(gè)角隨之確定2 號(hào)的位置就解出來了。但這里有個(gè)硬邊界當(dāng)某架無人機(jī)位于 1、2、5 號(hào)構(gòu)成的三角形的邊的延長線上時(shí)只能傳送一個(gè)方向信息無法準(zhǔn)確定位。所以第二輪定位時(shí)4、7、9、11、13、14 號(hào)無人機(jī)依然無法確定位置只能先定位 3、6、8、10、12、15 號(hào)。第三輪選取 2、6、8 為發(fā)射機(jī)剩下的未定位無人機(jī)才都能解算。5.2 錐形編隊(duì)位置調(diào)整的 MATLAB 實(shí)現(xiàn)要點(diǎn)附錄四的程序結(jié)構(gòu)與附錄三類似但坐標(biāo)系和理想位置表不同。錐形編隊(duì)的理想位置表如下無人機(jī)編號(hào)X 坐標(biāo)Y 坐標(biāo)1100√350250√3503004-50√3505-100√31006-100√307-100√3-1008-50√3-509-50√35010-50√3-15011-100√320012-100√310013-100√3014-100√3-10015-100√3-200調(diào)整邏輯與圓形編隊(duì)一致設(shè)定步長向四周遍歷直到偏差小于閾值。區(qū)別在于錐形編隊(duì)的理想位置不是均勻分布在圓周上而是按等邊三角形網(wǎng)格排列所以極坐標(biāo)更新時(shí)需要先轉(zhuǎn)成直角坐標(biāo)再計(jì)算距離。% 錐形編隊(duì)位置調(diào)整核心邏輯 target_x -100*sqrt(3); % 示例15 號(hào)機(jī)目標(biāo) X target_y -200; % 示例15 號(hào)機(jī)目標(biāo) Y step 0.01; current_x input(當(dāng)前 X 坐標(biāo)); current_y input(當(dāng)前 Y 坐標(biāo)); while sqrt((current_x - target_x)^2 (current_y - target_y)^2) step % 計(jì)算當(dāng)前距離 dist sqrt((current_x - target_x)^2 (current_y - target_y)^2); % 沿梯度方向移動(dòng)一個(gè)步長 dx (target_x - current_x) / dist * step; dy (target_y - current_y) / dist * step; current_x current_x dx; current_y current_y dy; end這段代碼比圓形編隊(duì)的極坐標(biāo)遍歷更直接已知目標(biāo)坐標(biāo)直接沿梯度方向走步長。參數(shù)step控制精度target_x和target_y從理想位置表里查。實(shí)際運(yùn)行時(shí)每架無人機(jī)需要先通過純方位定位確定自己的當(dāng)前位置再查表找到目標(biāo)位置然后執(zhí)行梯度下降。注意錐形編隊(duì)里 4、7、9、11、13、14 號(hào)無人機(jī)在第二輪無法定位必須等第三輪 2、6、8 號(hào)調(diào)整到位后才能解算。如果強(qiáng)行在第二輪對這些無人機(jī)執(zhí)行調(diào)整會(huì)因?yàn)槲恢霉烙?jì)錯(cuò)誤導(dǎo)致越調(diào)越偏。6. 從跑通到跑對四個(gè)容易翻車的參數(shù)與驗(yàn)證習(xí)慣第一個(gè)坑是角度單位。MATLAB 的acosd返回的是度但cos和sin默認(rèn)是弧度。附錄一的程序里混用了cosd和cos如果輸入的角度沒有統(tǒng)一算出來的坐標(biāo)會(huì)差一個(gè) π/180 的系數(shù)。我一般會(huì)在程序開頭加一行clear; clc;然后把所有三角函數(shù)統(tǒng)一成帶 d 的版本。第二個(gè)坑是符號(hào)計(jì)算的虛部。solve解方程組時(shí)可能產(chǎn)生復(fù)數(shù)解double轉(zhuǎn)換后虛部雖然很小但不為零。real函數(shù)能去掉虛部但如果虛部量級(jí)超過 1e-6說明方程本身有問題需要檢查角度輸入是否合理。第三個(gè)坑是枚舉法的邊界。附錄二里FY0j和FY0k從 2 到 9 遍歷但alpha1矩陣的索引是(FY0j, FY0k)當(dāng)FY0j或FY0k等于FY0k接收機(jī)編號(hào)時(shí)對應(yīng)的角度值沒有物理意義。程序里用continue跳過了相同編號(hào)但沒有跳過接收機(jī)編號(hào)實(shí)際使用時(shí)需要手動(dòng)排除。第四個(gè)坑是步長與收斂速度的權(quán)衡。0.01 米步長在 100 米半徑下需要約 10000 次迭代才能從 100 米外收斂到目標(biāo)每次迭代都要重新計(jì)算方位角計(jì)算量不小。我一般會(huì)先用 1 米步長粗調(diào)再用 0.01 米精調(diào)總迭代次數(shù)能降到十分之一。驗(yàn)證方法很簡單跑完附錄一算出坐標(biāo)后代回附錄二反解發(fā)射機(jī)編號(hào)如果得到的編號(hào)與輸入一致說明定位模型自洽。再用附錄三調(diào)整位置檢查最終坐標(biāo)與理想位置表的偏差是否小于 0.01 米。三個(gè)程序串起來跑通這套純方位無源定位方案就算真正落地了。從那以后我每次拿到新的編隊(duì)模型都強(qiáng)制先跑一遍“正解-反解-調(diào)整”的閉環(huán)驗(yàn)證確認(rèn)三個(gè)環(huán)節(jié)的數(shù)值能對上再往下做。希望幫到你。本文還有配套的精品資源點(diǎn)擊獲取