简介:这份PDF文档为2022年全国大学生数学建模竞赛本科组B题的完整研究论文,面向备赛数模竞赛的高校学生、指导教师及关注无人机纯方位无源定位的读者。全文围绕圆形与锥形编队下无人机位置偏差调整展开,涵盖问题重述、模型假设、符号说明,以及基于二维直角坐标系和余弦定理的被动接收信号无人机定位模型;对额外发射信号无人机数量、三圆共点判定、步长遍历位置调整和锥形编队扩展方案均有推导与讨论,并给出MATLAB编程实现思路。资源包内共1个PDF文件,大小约2.63MB,便于直接阅读和打印研习。目前已有8041人学习下载,可作为赛题论文写作范本、建模思路梳理与MATLAB复现的参考材料,尤其适合需要理解方向信息解算、误差阈值设置与编队迭代定位流程的读者。
1. 从“看不见的偏差”切入:这套纯方位无源定位方案到底解决什么
编队飞行里最让人后背发凉的不是某架无人机飞偏,而是它飞偏之后仍然以为自己在队形里。2022年数学建模竞赛本科组B题把这个问题摊开:十架无人机做圆形编队,靠纯方位无源定位,只用接收到的方向角把位置找回来,再给出调整方案;换成锥形编队后,同一套几何思路还要继续复用。做这道题的人通常不是缺公式,而是缺一条从坐标基线、定位模型、发射机数量到位置调整和验收的完整链路。
这套方案适合三类人:正在做数学建模课设、需要交完整论文和 MATLAB 程序的人;做无人机编队仿真、想复现纯方位定位流程的人;以及要把“角度反解坐标”讲成可验收模块的从业者。它不追求把每架机都写成复杂控制对象,而是先把几何关系锁死,再用枚举和迭代把偏差一点点压下去。
2. 编队几何与总体方案:圆形十机、锥形十五机的坐标基线怎么搭
纯方位无源定位的第一步不是写代码,而是把“理想位置”定义清楚。没有稳定的坐标基线,后面所有角度都会变成玄学。圆形编队和锥形编队看起来是两种队形,底层做法一致:先确定无偏差发射机的几何位置,再让接收机根据方向角反解自身位置,最后用一个小步长循环逼近理想点。整套 MATLAB 程序按“定位模型、额外发射机数量、圆形调整、锥形调整”四块拆,输入输出边界清楚,后面改参数才不至于一改全崩。
2.1 圆形编队:十个机位、三个已知点和四十度夹角
圆形编队由十架无人机组成,0号机固定在圆心,1到9号机均匀分布在半径为 r 的圆周上。相邻两架圆周机之间的圆心角是 40 度,所以 1 号机放在 x 轴正方向后,2 号机是 40 度,3 号机是 80 度,一直到 9 号机 320 度。问题一第三小问给的是 r=100 m,这个半径不是装饰参数,后面弦长、余弦定理和步长遍历都要用它。
| 编号 | 极角 | 直角坐标(r=100 m) |
|---|---|---|
| 0 | 圆心 | (0, 0) |
| 1 | 0° | (100.0, 0.0) |
| 2 | 40° | (76.6, 64.3) |
| 3 | 80° | (17.4, 98.5) |
| 4 | 120° | (-50.0, 86.6) |
| 5 | 160° | (-94.0, 34.2) |
| 6 | 200° | (-94.0, -34.2) |
| 7 | 240° | (-50.0, -86.6) |
| 8 | 280° | (17.4, -98.5) |
| 9 | 320° | (76.6, -64.3) |
这张表在程序里通常不手写死,而是用循环生成:theta=(i-1)*40,x=r*cosd(theta),y=r*sind(theta)。手写表容易在第 5、6 号机的负号上翻车,尤其是复制粘贴后只改一半。坐标基线一旦固定,问题一第一小问的“两架圆周机加圆心发射、另一架接收”就有了统一参照:0号机是原点,i、j 是圆周上的两个定点,k 是待定位的接收机。
2.2 锥形编队:十五机等边三角形网格与两机无偏差假设
锥形编队不是把圆形队形简单压扁,而是由等边三角形网格铺出来的。直线上相邻两架无人机间距相等,相邻三架构成等边三角形。按题目给的锥形队形,可以整理出十五个机位的理想坐标,整体呈锥形展开:有的机位在轴线上,有的在两侧斜线上,横向和纵向间距都由同一个基准长度 d 推出来。写程序时不要直接把十五个坐标硬编码成一张表,最好用“行、列、基准间距”生成,再给每个机位编号,这样第二轮、第三轮换发射机组时不会把编号和坐标弄混。
锥形编队比圆形多了一个现实约束:如果只有 0 架或 1 架无人机位置无偏差,根本没法给其他机提供足够的方向信息,定位方程欠定。所以通常假设 1 号和 5 号机位置无偏差,并且二者距离保持不变,把它们当作初始基准。接下来先定位 2 号机:让 1、2 号发射,5 号接收,得到一个方向信息;再让 2、5 号发射,1 号接收,又得到一个方向信息。三角形 1-2-5 里两个角已知,第三个角自然确定,2 号机的实际位置就被闭合出来。这个“两角定第三点”的思路,是锥形编队能继续往下调度的关键。
2.3 MATLAB 程序分工与数据流
四个附录程序不要写成四份互不相干的脚本。比较稳的做法是抽一个公共函数:输入发射机编号、接收机编号、半径和角度,输出理想坐标、实际坐标或三个模型角。定位模型负责“角度到坐标”,额外发射机数量程序负责“枚举表加平方和最小”,圆形调整程序负责“0.01 米步长遍历”,锥形调整程序负责“两角闭合加三轮调度”。数据流是单向的:几何基线生成理想坐标,定位函数反解实际坐标,误差函数算偏差,调整循环改坐标,调度表决定下一轮谁发射。
function [p_ideal, p_actual] = locate_one(k, i, j, r, alpha) % k: 接收机编号;i,j: 发射机编号;r: 编队半径 % alpha: 三个方位角,单位度,顺序为 [角Oki, 角Okj, 角jki] theta = (0:9)' * 40; % 圆周机极角 P = [0 0; r*cosd(theta(1:9)), r*sind(theta(1:9))]; p_ideal = P(k, :); % 该机理想位置 % 下面进入余弦定理反解,实际位置由 alpha 和 i,j 几何决定 % 具体方程组见第3章,这里保留统一入口,避免四个脚本各写一套 end这段函数签名把“编号、半径、角度”固定成参数,后面换圆形或锥形时只换坐标生成部分,不动调用方式。参数说明要写进注释:k是接收机编号,i、j是发射机编号,r是圆形编队半径,alpha是三个方向角。很多队伍把角度顺序写反,导致alpha1和alpha2对调,最后定位结果左右对称,看起来像“模型错了”,其实是输入顺序错了。
3. 被动接收机定位模型:三组方位角如何反解实际位置与理想位置
问题一第一小问有两个容易混淆的目标:一个是根据接收到的方向信息确定接收机的实际位置,另一个是确定它应当到达的准确位置。实际位置用余弦定理反解,准确位置用三条圆弧轨迹的交点。两者都建立在同一个二维直角坐标系上,但前者吃的是“实测角”,后者吃的是“无偏差发射机几何”。把这两件事分开写函数,后面调整模块才不会把“当前在哪”和“该去哪”混成一个变量。
3.1 实际位置:余弦定理反解三条边
设 0 号机在圆心 O,i、j 是两架圆周发射机,k 是接收机。已知三个方向角:α1 是接收机看 O 与 i 的夹角,α2 是看 O 与 j 的夹角,α3 是看 j 与 i 的夹角。未知量是三条距离:k0=Ok,a2=kj,a3=ki。圆周上两发射机之间的距离 a1=ij 可由半径和夹角直接算:
a1 = 2*r*sind(abs(i-j)*20); % 圆周相邻40度,半角20度,|i-j|个间隔然后在三角形 O-k-i、O-k-j、i-k-j 里分别用余弦定理,得到三个方程:
syms k0 a2 a3 eq1 = (k0^2 + a3^2 - r^2)/(2*k0*a3) - cosd(alpha(1)); eq2 = (k0^2 + a2^2 - r^2)/(2*k0*a2) - cosd(alpha(2)); eq3 = (a2^2 + a3^2 - a1^2)/(2*a2*a3) - cosd(alpha(3)); S = vpasolve([eq1 eq2 eq3], [k0 a2 a3], [0 2*r; 0 2*r; 0 2*r]);这三个方程的几何含义很直白:第一个方程把 α1 锁在三角形 O-k-i 里,第二个把 α2 锁在 O-k-j 里,第三个把 α3 锁在 i-k-j 里。解出 k0、a2、a3 后,再用“以 O 为圆心 k0、以 i 为圆心 a3”的两个圆求交点,会得到上下两个候选点,用 α2 或 α3 的符号筛掉另一个。参数上要注意:r用米,角度统一用度,cosd和sind不要混成cos、sin;vpasolve的搜索区间要覆盖 0 到 2r,否则可能只返回一个不合理的负根。
实际坐标解出来后,和理想坐标p_ideal一比,偏差向量就有了。这个偏差不是直接拿来做控制,而是先记录误差,再交给调整循环。很多题解在这里直接开始“向北还是向南”,没有把实际位置函数独立出来,导致后面圆形和锥形各写一套反解,代码越写越像黑匣子。
3.2 理想位置:三条圆弧轨迹与三圆共点
理想位置不是靠“猜”出来的。平面几何里有一个基本事实:一个动点到两个定点的夹角是定值,动点轨迹是由上下对称两段圆弧合成的封闭曲线。夹角是钝角时轨迹像“眼睛”,夹角是锐角时轨迹像“葫芦”。当有三架位置无偏差的发射机时,接收机会得到三条这样的轨迹;三条轨迹的公共交点,就是它应当到达的准确位置。所谓“三圆共点”,本质就是三个圆弧约束同时成立。
程序里不需要把“眼睛”和“葫芦”画出来再目测交点,直接用两组定点加一个角度构造轨迹方程,联立求交即可。设两组定点为 (O,i) 和 (O,j),轨迹上的点 k 满足 ∠Ok i=α1、∠Ok j=α2,再把第三组 (i,j) 的 ∠jki=α3 加进来做筛选。实际写的时候,可以先用前两组轨迹求出候选交点,再用第三组角度验算,误差最小的候选点就是理想位置。这样比三组方程直接硬解更稳,也方便排查“交点跑到圆外”的问题。
这里要区分“实际位置”和“理想位置”的数据来源:实际位置由接收到的三个实测角反解,理想位置由无偏差发射机的几何关系确定。题目里说“定位无人机应当到达的准确位置”,说的就是后者。调整模块比较的是实际坐标和理想坐标,不是拿实测角和模型角直接相减就完事。角度残差可以用来判断匹配,但最终停止条件要落到坐标距离上。
3.3 额外一架就够:56 种组合与平方和最小
问题一第二小问规定 0 号和 1 号发射信号,问还需要额外增加几架发射机,才能让其他接收机准确定位。先看只增加 0 架:除 0、1 外,某架接收机只能收到一个方向信息,它的可能位置形成一条“眼睛”或“葫芦”轨迹,无法唯一确定。只增加 1 架后,接收机能收到三个方向信息,轨迹有三条,三圆共点,位置就能唯一确定。增加 2 架、3 架当然更冗余,但题目问的是“需要额外增加几架”,最小可行解就是 1 架。
枚举空间是 8×7=56 种:除 0、1 外有 8 架可能的接收机,每架接收机再从剩下的 7 架里选一架作为额外发射机。对每一种组合,用无偏差几何算出三个模型角,和实测三个角做平方和:
best = inf; best_j = -1; for k = 2:9 for j = 2:9 if j == k, continue; end alpha_model = angle_table(j, k, r); % 三个模型角 cost = sum((beta(:) - alpha_model(:)).^2); if cost < best best = cost; best_j = j; end end endbeta是实测三个方向角,alpha_model是查表得到的模型角,cost是三个角残差的平方和。遍历完取最小,对应的best_j就是需要额外增加的发射机编号。注意表 2 的横纵轴有方向:横向通常是接收信号无人机,纵向是发射信号无人机,转置后索引会整体错位。程序里最好用结构体或二维数组加显式行列标签,不要只写alpha(i,j)然后靠记忆判断谁是谁。
举个检查例子:当圆形编队半径取 10,接收机为 8 号时,枚举应能锁定额外发射机为 9 号;接收机为 3 号时,应能锁定额外发射机为 5 号。这两个例子不是用来背答案的,而是用来验证枚举表的行列方向和平方和计算有没有写反。只要有一个例子对不上,先查角度顺序,再查表转置,最后查编号是否从 0 开始。
4. 圆形与锥形编队的位置调整:步长遍历、延长线盲区与三轮调度
定位模型给出“在哪”和“该去哪”,调整模型负责“怎么过去”。圆形编队和锥形编队的调整逻辑有共同点:都用小步长向理想位置逼近,都把“调整时间不计”当作假设,所以可以认为无人机一旦算出目标方向,就能瞬间完成微调。区别在于锥形编队存在“延长线盲区”,某些机位在发射机三角形边的延长线上时只能收到一个有效方向信息,必须靠换发射机组来绕开。
4.1 圆形编队:0.01 米步长与极坐标更新
圆形编队里 0 号和 1 号已经处于理想位置,飞行过程中保持二者相对位置关系不变。先定位 2 号机:0、1 发射,2 接收,用第 3 章的模型算出 2 号机实际位置,再让它向理想位置调整。调整时设定步长 0.01 m,向四周全方位遍历,每走一步算一次与理想位置的距离,当距离不超过 0.01 m 时,就近似认为到达理想位置。2 号机调好后加入“无偏差元组”,接着用 0、1、2 定位 3 号,再调 3 号,如此逐个推进。
step = 0.01; tol = 0.01; maxIter = 20000; for iter = 1:maxIter if norm(p_actual - p_ideal) <= tol, break; end bestErr = inf; bestPhi = 0; for phi = 0:359 p_try = p_actual + step * [cosd(phi), sind(phi)]; err = norm(p_try - p_ideal); if err < bestErr bestErr = err; bestPhi = phi; end end p_actual = p_actual + step * [cosd(bestPhi), sind(bestPhi)]; errLog(iter) = bestErr; end这段循环里step是每次移动的步长,tol是停止容差,phi是 0 到 359 度的方向遍历。实际写的时候可以把“向北/向南”两类方向合并成 360 度全遍历,避免漏掉斜向调整。errLog用来记录每轮误差,答辩时能直接画收敛曲线。注意步长不要大于容差,否则会在目标点两侧来回跳;如果初始偏差很大,可以先粗步长走几轮,再切到 0.01 m 精调。
4.2 锥形编队:两个夹角闭合第三点与延长线盲区
锥形编队先假设 1 号和 5 号位置无偏差,且二者距离保持不变。定位 2 号机时,让 1、2 发射,5 接收,得到一个方向信息;再让 2、5 发射,1 接收,又得到一个方向信息。三角形 1-2-5 中两个角已知,第三个角由内角和确定,2 号机的实际位置就被闭合出来。接着用和圆形编队一样的 0.01 m 步长遍历,把 2 号机调到理想位置。此时 1、2、5 都是无偏差机,可以作为下一轮发射机组。
锥形编队的坑在“延长线”。当某架待定位机正好位于三架发射机组成的三角形的边的延长线上时,它只能收到一个方向信息,无法准确定位。第一轮之后,4、7、9、11、13、14 号机就属于这种情况,第二轮无法直接定位;第二轮能定位的是 3、6、8、10、12、15 号机。把这些机调好后,第三轮选 2、6、8 作为发射机,剩下的未定位机不再处于这三架机三角形的延长线上,于是都能定位并调整。
| 轮次 | 无偏差/发射机组 | 可定位并调整 | 暂不能定位 |
|---|---|---|---|
| 第一轮 | 1、5 | 2 | 4、7、9、11、13、14 等 |
| 第二轮 | 1、2、5 | 3、6、8、10、12、15 | 4、7、9、11、13、14 |
| 第三轮 | 2、6、8 | 4、7、9、11、13、14 | 无 |
这张调度表是锥形编队的核心交付物之一。它说明“定位模型能用”不等于“调度方案能跑通”:如果第二轮硬用 1、2、5 去定位延长线上的机位,程序会一直得到不唯一解,看起来像方程无解,其实是发射机组选错了。
4.3 调整循环的停止条件与误差记录
调整循环的停止条件要统一成坐标距离,不要用“迭代了多少次”代替。圆形编队用norm(p_actual-p_ideal)<=0.01,锥形编队同样用这个标准。每次调整后把误差写进日志,最后统计最大误差、平均误差和收敛轮数。如果某架机误差降不下去,先看它是不是在延长线上,再看步长是否大于容差,最后看理想坐标有没有生成错。
err = norm(p_actual - p_ideal); fprintf('编号=%d, 轮次=%d, 误差=%.4f m\n', k, iter, err); if err <= 0.01 status(k) = "已调整"; else status(k) = "未收敛"; endstatus数组按编号记录每架机的调整状态,后面生成调度表时直接查状态,不用靠人工记“谁调过了”。误差打印保留四位小数,方便看出 0.01 m 附近的抖动。参数上,maxIter要给够,但不要无限循环;如果超过上限还没收敛,应该输出当前误差和候选方向,而不是静默退出。
5. 避坑与常见问题排查:坐标、角度、枚举、延长线四类翻车点
这类题的代码量不算大,翻车点却非常集中:角度单位、枚举表方向、延长线盲区、步长收敛。每个问题单看都像小毛病,叠在一起就会让结果“看起来差不多,实际全错”。下面四类是我在复现时优先排查的点,每条按现象、原因、解决写,照着查能省很多时间。
5.1 角度单位和象限:cosd 与 acosd 混用会翻车
现象是定位结果左右对称,或者明明输入角在 0 到 180 度之间,反解坐标却跑到圆外。原因通常是cosd、sind和acos、asin混用,或者acosd只返回 0 到 180 度,丢了上下半平面的符号。解决方法是全篇统一用度:cosd、sind、acosd、asind配套使用;反解出两个候选交点时,用第三组角度或叉积符号筛掉另一个。写代码前先手算一个 1 号机在 x 轴正方向、2 号机在 40 度的例子,确认符号方向再往下跑。
5.2 枚举表转置与编号:56 种组合对不上
现象是每架接收机都匹配到同一个额外发射机,或者表 2 里明明有的组合程序里找不到。原因是枚举表的横纵轴定义和循环里的i,j写反了,或者编号从 1 开始、坐标数组却从 1 开始存 0 号机。解决方法是先把表头打印出来:横向是接收机编号,纵向是发射机编号;循环里用if j==k continue排除自己;匹配前把beta和alpha_model都整理成 1×3 行向量,避免beta(:)和alpha_model(:)方向不一致。
5.3 延长线盲区:锥形第二轮漏定位
现象是锥形编队第二轮里 4、7、9、11、13、14 号机反复算不出唯一位置,程序报无解或给出两个对称点。原因是这些机位在 1、2、5 号机三角形的边的延长线上,只能收到一个有效方向信息,方程欠定。解决方法是换发射机组,第三轮用 2、6、8 号机发射,把延长线上的机位重新纳入可定位范围。调度表要提前写进程序,不要等循环失败再手动改编号。
5.4 步长与收敛:0.01 米不是万能
现象是误差在 0.01 到 0.03 m 之间来回跳,或者初始偏差大时循环几十轮都到不了目标点。原因是步长大于容差,或者一开始就用 0.01 m 走远距离。解决方法是先粗后细:初始偏差大时用 0.1 m 甚至 0.5 m 走几轮,进入 0.05 m 范围后再切 0.01 m;容差设成 0.01 m 时,步长不要超过 0.005 m 到 0.01 m。每次循环记录误差,误差连续多轮不下降就停下来检查理想坐标和方向角符号。
6. 验收与答辩:用“可复现证据链”讲清定位精度与调度轮次
验收这份方案,我一般不看“图好不好看”,只看证据链能不能复现。先把圆形编队的理想坐标表打印出来,核对 1 到 9 号机的极角是不是 0、40、80 一直到 320 度,半径是不是 100 m;再固定一组发射机编号和三个方向角,跑定位函数,看反解坐标和理想坐标的误差是否落在 0.01 m 以内。第二步跑 56 种枚举,检查每架接收机的最小平方和是否唯一,额外发射机编号是否和表 2 对得上。第三步从有偏差的初始位置跑圆形调整循环,统计每架机最终误差和收敛轮数。第四步按三轮调度跑锥形编队,确认第三轮后 4、7、9、11、13、14 号机都能定位并调整。
| 验收项 | 检查方法 | 通过标准 |
|---|---|---|
| 定位模型 | 固定 r=100,输入三组模型角反解坐标 | 与理想坐标误差≤0.01 m |
| 额外发射机数量 | 遍历 56 种组合,看最小 cost | 每架接收机有唯一最优解 |
| 圆形调整 | 从偏差初始位置跑步长循环 | 所有机误差≤0.01 m |
| 锥形调度 | 按三轮发射机组跑 | 第三轮后无未定位机 |
| 代码复现 | 清空变量,按附录顺序运行 | 无未定义变量、无硬编码编号 |
答辩时不要从“无人机很复杂”讲起,先讲三个假设:飞行无外界干扰、接收方向信息准确、调整时间不计。然后画出圆形编队的坐标系,指出 0 号机在圆心、1 到 9 号机每 40 度一个。接着现场跑一个输入,把三个方向角输进去,展示实际坐标、理想坐标和误差。评委如果追问“为什么额外一架就够”,就讲单一方位角的眼睛/葫芦轨迹不唯一,增加一架后三条轨迹共点。追问锥形为什么第二轮漏机,就讲延长线盲区,再给三轮调度表。最后把误差日志和收敛曲线放出来,比任何口头解释都硬。
从那以后我每次带这类定位题,都强制先跑一遍“坐标表—单点反解—枚举表—调整日志”四步,不跳过任何一步。看起来慢,其实比后面满屏找符号错误快得多。希望帮到你。
本文还有配套的精品资源,点击获取