简介：基于MATLAB的船舶运动仿真资源，面向船舶运动控制、建模与仿真的学习者及工程师，可用于船舶操纵性、推进控制等场景的仿真分析，覆盖船舶运动模型构建、控制器设计、仿真结果可视化等完整环节。压缩包共12个文件，包含5个可运行的m源文件、3个asv自动备份文件、2份docx设计说明文档、1份txt全局变量说明以及1张运行结果图，整体大小仅278KB，轻量易部署；其中m文件实现核心仿真算法，docx文档讲解程序设计思路与运行步骤，txt说明全局变量定义，jpg直观展示运行效果。已有347人学习下载，口碑与实用性得到初步验证。基于MATLAB 2019b开发的仿真代码可直接运行，用户解压后将文件放入当前目录并调用主函数即可获得完整船舶运动仿真结果，无需额外配置环境；代码通过舵、螺旋桨等控制模块模拟运动响应，结构清晰、便于理解，适合作为课程设计或毕业设计的基础参考，也支持在此基础上进行二次开发与功能扩展，是学习船舶运动仿真的实用入门资料。
1. 船舶运动仿真为什么先确认坐标系而不是先跑代码
一艘 30 米工作船在 Matlab 里做回转仿真，从压舵到画出第一个圆，通常只要几十秒；同样的机动在海面上要花十几分钟。仿真跑得快是优势，但也容易掩盖问题。很多下载了「船舶运动仿真 Matlab 源码」的同行，第一反应是解压、双击 main.m、盯着轨迹看结果；真正该先确认的是坐标系、单位制和状态向量顺序。这个标题背后是一套完整需求：把船体运动数学模型、推进器与舵的输入、以及可视化环境拼成一个可复用的仿真原型，支撑航迹控制、无人艇路径跟踪和船用设备半实物仿真。对算法工程师来说，这套模型还能用来生成训练数据。最终要的不是“能动的模型”，而是“能验证对错的模型”。
2. 三自由度船体运动模型：坐标系、耦合项与参数表的取舍
2.1 固定坐标系与随船坐标系的换算
船舶运动仿真里最常出现的一个隐蔽错误，是把船体坐标系下的速度直接当成了轨迹速度。船体坐标系固连在船上，x 轴指向船艏、y 轴指向右舷，z 轴向下；而轨迹图关心的是固定坐标系（北东坐标）下的 X、Y 位置。两者之间靠一个绕 z 轴的旋转矩阵转换，偏航角记为 ψ，单位必须是弧度：
$$ \begin{bmatrix} \dot{x}_E \ \dot{y}_E \end{bmatrix}
\begin{bmatrix} \cos\psi & -\sin\psi \ \sin\psi & \cos\psi \end{bmatrix} \begin{bmatrix} u \ v \end{bmatrix} $$
这里 u 是纵荡速度（前进），v 是横荡速度（横向漂移）。很多人把 v 的符号写反，或者把 ψ 当成角度制传入 sin/cos，导致轨迹呈现出违背物理的螺旋外扩。对于只做水平面机动、船体横摇纵摇影响不大的场景，这个三自由度模型已经够用；只有研究耐波性或波浪载荷时，才需要把横摇、纵摇和垂荡补成六自由度。
2.1.1 Matlab 里坐标转换的一段最小代码
function XY = body2earth(uv, psi) % uv: 2x1 向量，[u; v]，船体坐标系下的速度 % psi: 偏航角，弧度 R = [cos(psi) -sin(psi); sin(psi) cos(psi)]; XY = R * uv; end调用前先统一单位。常见做法是把外界输入的航向角统一转成弧度再参与运算，而不是在公式里临时乘 pi/180。一旦仿真发散，第一步就应该检查所有传给三角函数的角度量是否已经转弧度。
2.2 动力学方程里的耦合项不能随便删
三自由度动力学模型通常写成“广义质量 × 加速度 = 水动力 + 控制力”的形式。设状态为[u, v, r]，其中 r 是艏摇角速度，广义质量分别为：
- m_x = m - X_dot_u，纵荡方向质量加附加质量
- m_y = m - Y_dot_v，横荡方向质量加附加质量
- m_r = I_zz - N_dot_r，艏摇转动惯量加附加惯量
典型的非线性方程如下：
$$ m_x \dot{u} = m_y v r - X_u u - X_{u|u|} u|u| + \tau_u $$
$$ m_y \dot{v} = -m_x u r - Y_v v - Y_{v|v|} v|v| + \tau_v $$
$$ m_r \dot{r} = (m_x - m_y) u v - N_r r - N_{r|r|} r|r| + \tau_r $$
其中 τ_u 是螺旋桨推力，τ_v 和 τ_r 来自舵力，写成：
$$ \tau_v = Y_\delta \delta, \quad \tau_r = N_\delta \delta $$
m_y v r和-m_x u r是科氏力和向心力耦合项。直航时 v=0、r=0，如果删掉耦合项，模型就退化成三条互相独立的微分方程；此时压舵虽然能产生力矩，但因为缺少“前进速度 + 艏摇角速度 → 横荡速度”的耦合路径，回转试验里那种“边转边横漂”的物理过程就完全复现不出来。很多源码包把耦合项放在注释里，移植时顺手删掉，结果模型只剩一个空壳。
2.3 船体参数表与估算方法
给一艘 30 米级工作船做仿真，参数可以按半经验方式估算。附加质量通常取排水质量的 10%～20%，转动惯量比质量多一到两个量级，水动力阻尼则按船型估算。以下是一组能跑出合理回转轨迹的参数参考：
| 参数 | 含义 | 参考值 | 单位 |
|---|---|---|---|
| m | 排水质量 | 1.2e6 | kg |
| m_x | 纵荡广义质量 | 1.5e6 | kg |
| m_y | 横荡广义质量 | 1.8e6 | kg |
| m_r | 艏摇广义惯量 | 1.4e8 | kg·m² |
| X_u | 纵荡线性阻尼 | 2.0e4 | N·s/m |
| X_uabs | 纵荡二次阻尼 | 1.2e4 | N·s²/m² |
| Y_v | 横荡线性阻尼 | 5.0e4 | N·s/m |
| Y_vabs | 横荡二次阻尼 | 3.0e5 | N·s²/m² |
| N_r | 艏摇线性阻尼 | 4.0e7 | N·m·s |
| N_rabs | 艏摇二次阻尼 | 4.0e8 | N·m·s² |
| Y_delta | 舵力系数 | 3.5e5 | N/rad |
| N_delta | 舵力矩系数 | 9.0e6 | N·m/rad |
拿到任何一份源码，先找参数表对应的变量名。不同源码对附加质量的处理方式不一样：有的把X_udot单独建模，有的像这里直接合并成 m_x。参数对不上，仿真结果就会南辕北辙。
3. 用 Matlab 脚本还原船舶运动仿真源码的最小可跑边界
3.1 状态向量、输入接口与 ODE45 的接法
常见做法是把状态定义成长度 6 的向量：前三个是船体坐标系下的 u、v、r，后三个是固定坐标系下的 x_E、y_E、ψ。这样运动学方程和动力学方程放在同一个导数函数里，ODE45 一次搞定。以下是一份可运行的模型代码：
function dx = shipODE(t, x, p) % 状态向量: x = [u; v; r; xE; yE; psi] u = x(1); v = x(2); r = x(3); psi = x(6); % 控制输入从参数结构体读取 tauU = p.tauU; % 螺旋桨推力(N) delta = p.delta; % 舵角(rad)，负值表示右舵 % 动力学：广义质量 x 加速度 = 水动力 + 控制力 u_dot = (p.m_y * v * r - p.X_u * u - p.X_uabs * u * abs(u) + tauU) / p.m_x; v_dot = (-p.m_x * u * r - p.Y_v * v - p.Y_vabs * v * abs(v) + p.Y_delta * delta) / p.m_y; r_dot = ((p.m_x - p.m_y) * u * v - p.N_r * r - p.N_rabs * r * abs(r) + p.N_delta * delta) / p.m_r; % 运动学：船体速度转到固定坐标系 xE_dot = u * cos(psi) - v * sin(psi); yE_dot = u * sin(psi) + v * cos(psi); psi_dot = r; dx = [u_dot; v_dot; r_dot; xE_dot; yE_dot; psi_dot]; end这个函数里，p 是参数结构体，所有模型参数都在外部定义一次，避免在模型内部写死数字。好处是后面做参数扫描、做控制器联合仿真时，只需要改 p 的字段，不需要动模型本体。ODE45 对刚性问题会自适应缩小步长，但是如果参数量纲不一致，u、r 之间会出现量级差异，照样会发散。
3.2 从参数表到可运行模型的初始化脚本
主脚本负责初始化参数、设定初始状态、调用积分器并画图：
% 主脚本：初始化参数并运行回转仿真 p.m_x = 1.5e6; p.m_y = 1.8e6; p.m_r = 1.4e8; p.X_u = 2.0e4; p.X_uabs = 1.2e4; p.Y_v = 5.0e4; p.Y_vabs = 3.0e5; p.N_r = 4.0e7; p.N_rabs = 4.0e8; p.Y_delta = 3.5e5; p.N_delta = 9.0e6; p.tauU = 1.2e5; % 恒定推力 p.delta = -10 * pi / 180; % 右舵 10 度 x0 = [2; 0; 0; 0; 0; 0]; % 初始航速 2 m/s，位置在原点 tSpan = [0 600]; opt = odeset('RelTol', 1e-6, 'AbsTol', 1e-6); [t, X] = ode45(@(t, x) shipODE(t, x, p), tSpan, x0, opt); plot(X(:,4), X(:,5), 'LineWidth', 1.5); axis equal; grid on; xlabel('x_E (m)'); ylabel('y_E (m)'); title('船舶回转仿真：右舵 10 度');这里X(:,4)和X(:,5)是固定坐标系下的轨迹。如果跑出螺旋状的轨迹，且半径随时间越来越大，说明模型里的阻尼项符号可能写错，或者atan2处理航向角时有跳变。仿真结尾还可以用min、max计算轨迹的范围，粗略评估回转直径是否在合理区间。
3.3 结果怎么判断“跑对了”
固定舵角下的船舶运动仿真，应该在几十秒后形成稳定的回转圆，而不是越转越大或者直接发散。判断标准有三条：一是 u 最终收敛在一个非零值附近，不会持续加速到离谱量级；二是 r 收敛到常数附近，轨迹呈现等间距圆环；三是轨迹方向与舵角方向一致，右舵对应右转。这些判断不需要精确数据，一眼能看出来的轨迹形状，通常比计算精度更早暴露问题。
4. 从脚本到 Simulink：把源码移植成框图时要盯住的 4 个点
4.1 为什么脚本能跑，框图却发散的常见原因
脚本能跑的模型搬进 Simulink 后发散，最常见原因是积分器初值没有同步。脚本里的x0明确写了初值，Simulink 的 Integrator 模块默认初值为 0，如果忘记把x0填进去，初始速度是 0，螺旋桨推力一加上去，船体瞬间获得巨大加速度，后续轨迹完全失真。第二个原因是代数环：MATLAB Function 模块如果同时把输入和输出接到同一个信号线上，会产生代数环，求解器每步都要迭代解隐式方程，仿真速度骤降甚至直接报错。第三是固定步长设置过大，像 0.1 秒的步长对于艏摇角速度 0.05 rad/s 的模型可能足够，但碰到舵角阶跃瞬间，瞬态量级可能会要求步长降到 0.01 秒以下。
4.2 MATLAB Function 模块与积分器的搭法
常见的 Simulink 搭建方式是：一个 MATLAB Function 模块放状态导数，输出接 Integrator，积分器输出再接回 Function 模块的输入，同时用 Mux 把舵角和推力接到第二个输入端口。Function 内部代码可以直接复用脚本模型，但要注意 Simulink 的 MATLAB Function 模块对结构体的支持有限，参数通常直接写在函数内部：
function xd = fcn(x, u) % x: 6x1 状态向量; u: [delta; tauU] delta = u(1); tauU = u(2); mx = 1.5e6; my = 1.8e6; mr = 1.4e8; Xu = 2.0e4; Xuabs = 1.2e4; Yv = 5.0e4; Yvabs = 3.0e5; Nr = 4.0e7; Nrabs = 4.0e8; Yd = 3.5e5; Nd = 9.0e6; xd = zeros(6, 1); xd(1) = (my * x(2) * x(3) - Xu * x(1) - Xuabs * x(1) * abs(x(1)) + tauU) / mx; xd(2) = (-mx * x(1) * x(3) - Yv * x(2) - Yvabs * x(2) * abs(x(2)) + Yd * delta) / my; xd(3) = ((mx - my) * x(1) * x(2) - Nr * x(3) - Nrabs * x(3) * abs(x(3)) + Nd * delta) / mr; xd(4) = x(1) * cos(x(6)) - x(2) * sin(x(6)); xd(5) = x(1) * sin(x(6)) + x(2) * cos(x(6)); xd(6) = x(3); endSimulink 模块表如下：
| 模块 | 作用 | 关键设置 |
|---|---|---|
| MATLAB Function | 状态导数计算 | 输入 x(6x1)、u(2x1)，内部直接写参数 |
| Integrator | 连续状态积分 | 初始条件填x0，与脚本一致 |
| Mux | 合并 delta 与 tauU | 输入为 2，输出 2x1 |
| To Workspace | 保存仿真结果 | 变量名设为 simX，输出格式选 Array |
代码生成时，如果参数需要频繁调试，可以把参数放到 MATLAB Function 模块的 Parameter 端口，或者用 Data Store Memory，但这类做法在快速原型阶段反而增加工作量大，直接写在函数里最方便调参。
4.3 二维和三维可视化的低成本方案
Simulink 跑完后，可视化不一定非要加 3D Animation 工具箱。二维轨迹图加一个随航向旋转的船形 patch 对象，已经够验证控制算法。常见做法是在仿真结束后用脚本画图：
% 用 patch 画一个简易船形轮廓，随航向角旋转 hull = [-0.6 -0.4 0.5 0.7 0.5 -0.4; -0.2 -0.25 0 0.2 -0.25 -0.2]; for k = 1:100:length(X) psi = X(k, 6); R = [cos(psi) -sin(psi); sin(psi) cos(psi)]; hullRot = R * hull + X(k, 4:5)'; patch(hullRot(1,:), hullRot(2,:), 'b'); hold on; end如果一定要三维效果，Matlab 2023a 之后的 3D Animation 工具箱支持将轨迹发送到虚拟场景，但需要额外授权。多数船舶运动仿真项目用二维轨迹加船形轮廓就足够表达问题。
5. 用回转试验验证船舶运动仿真模型是否正确
5.1 回转试验的判别标准
船舶操纵性验证的核心试验是回转试验：给定一个固定舵角，让船从直航进入稳定回转，测量稳态回转直径 D 与船长 L 的比值。货船通常在 2.5~4 倍船长之间，工作船机动性更好，D/L 可能低到 1.5~2.5。如果仿真出来的 D/L 大于 8 或者小于 1，模型参数大概率有问题。跑多组舵角对比是检查模型一致性的最简单方式：
deltaList = [-35 -10 10 35] * pi / 180; for k = 1:length(deltaList) p.delta = deltaList(k); [t, X] = ode45(@(t, x) shipODE(t, x, p), [0 600], x0, opt); plot(X(:,4), X(:,5), 'LineWidth', 1.2); hold on; end legend('-35 deg', '-10 deg', '10 deg', '35 deg');5.2 针对仿真发散的一组排查顺序
模型跑出发散轨迹时，按以下顺序检查：先看角度单位，舵角、航向角是否全为弧度；再看阻尼符号，所有线性阻尼项前的符号应该取负，二次阻尼项u*abs(u)形式写错会让能量持续注入系统；然后看初值，Simulink 积分器初值是否与脚本一致；最后检查步长，固定步长仿真时把步长缩到 0.01 秒试一次，如果轨迹发生明显变化，说明原步长已经不满足数值稳定性要求。一个实用技巧是监控 u 的数值，如果出现持续增大到几十米每秒，说明推力远大于阻尼，属于参数量级问题。
5.3 把仿真输出导出 CSV 做时间轴回放验证
验证模型可信度的一条实用路径：将状态与输入一起导出，生成带时间戳的 CSV 文件，再用 Python 或 Excel 复绘轨迹。具体做法是用writematrix把[t, X, deltaSeq, tauSeq]写成 CSV，回放时检查舵角变化时刻与轨迹曲率变化时刻是否对齐。这一步能发现数据记录时的索引错位问题，也能让不熟悉 Matlab 的同事直接看到仿真输入输出关系。
本文还有配套的精品资源，点击获取
