
本教程介绍MATLAB中三个核心的科学可视化技术:等值线图、三维矢量场流场模拟以及曲面拟合与插值可视化。
文章目录
等值线图与等高线
等值线图是表示二维标量场的有效方法,广泛应用于地形图、温度分布、压力场等科学领域的可视化。
基本概念
- 等值线(contour):连接具有相同数值的点形成的曲线
- 等高线:特指地形高度的等值线
- 等值面填充:在等值线之间填充颜色以增强视觉效果
contour函数:线条等值线图
contour函数用于绘制等值线,语法格式:
contour(X, Y, Z, levels)
基础示例:二元函数的等值线
figure;
% 创建网格
[X, Y] = meshgrid(-3:0.1:3, -3:0.1:3);
% 定义二元函数:高斯函数
Z = exp(-(X.^2 + Y.^2));
% 绘制等值线图
figure;
contour(X, Y, Z, 15); % 15条等值线
colorbar;
xlabel('X');
ylabel('Y');
title('高斯函数等值线图');
grid on;
运行结果:

进阶应用:自定义等值线数值
figure;
% 指定特定等值线数值
levels = [0.1, 0.3, 0.5, 0.7, 0.9];
[C, h] = contour(X, Y, Z, levels);
clabel(C, h, 'FontSize', 12); % 添加数值标签
colorbar;
title('自定义等值线数值');
运行结果:

contourf函数:填充等值线图
contourf函数在等值线之间填充颜色,提供更直观的视觉效果:
figure;
% 创建更复杂的函数进行演示
[X, Y] = meshgrid(-2:0.05:2, -2:0.05:2);
Z = sin(X) .* cos(Y) + 0.1 * (X.^2 + Y.^2);
contourf(X, Y, Z, 20); % 20个填充级别
colorbar;
xlabel('X');
ylabel('Y');
title('填充等值线图:sin(x)cos(y) + 0.1(x² + y²)');
colormap(jet); % 设置颜色映射
运行结果:

地形图模拟示例
% 模拟山地地形
[X, Y] = meshgrid(-10:0.2:10, -10:0.2:10);
Z1 = 50 * exp(-((X-2).^2 + (Y-2).^2)/8); % 山峰1
Z2 = 30 * exp(-((X+3).^2 + (Y-1).^2)/12); % 山峰2
Z3 = 40 * exp(-(X.^2 + (Y+3).^2)/15); % 山峰3
Z = Z1 + Z2 + Z3 + 5; % 基础海拔
figure;
subplot(1,2,1);
contour(X, Y, Z, 15);
colorbar;
title('地形等高线图');
xlabel('X (km)');
ylabel('Y (km)');
grid on;
subplot(1,2,2);
contourf(X, Y, Z, 20);
colorbar;
title('地形填充图');
xlabel('X (km)');
ylabel('Y (km)');
运行结果:

三维矢量场与流场模拟
矢量场可视化是计算流体力学(CFD)中的重要技术,用于显示速度场、电磁场等矢量数据。
基本矢量场可视化
二维矢量场:quiver函数
% 创建二维网格
[X, Y] = meshgrid(-2:0.3:2, -2:0.3:2);
% 定义矢量场:旋转场
U = -Y; % X方向分量
V = X; % Y方向分量
figure;
quiver(X, Y, U, V, 'b');
axis equal;
xlabel('X');
ylabel('Y');
title('二维旋转矢量场');
grid on;
运行结果:

三维矢量场:quiver3函数
% 创建三维网格
[X, Y, Z] = meshgrid(-1:0.4:1, -1:0.4:1, -1:0.4:1);
% 定义三维矢量场:螺旋场
U = -Y;
V = X;
W = 0.5 * Z;
figure;
quiver3(X, Y, Z, U, V, W, 0.5);
xlabel('X');
ylabel('Y');
zlabel('Z');
title('三维螺旋矢量场');
axis equal;
grid on;
运行结果:

流场模拟案例
案例1:圆柱体绕流
% 圆柱体绕流模拟(势流近似)
[X, Y] = meshgrid(-3:0.15:3, -2:0.15:2);
% 圆柱体参数
cylinder_radius = 0.5;
U_inf = 1; % 来流速度
% 计算距离
R = sqrt(X.^2 + Y.^2);
theta = atan2(Y, X);
% 势流解(无粘性近似)
% 避免圆柱体内部计算
mask = R > cylinder_radius;
U = zeros(size(X));
V = zeros(size(Y));
% 只在圆柱体外部计算速度场
U(mask) = U_inf * (1 - (cylinder_radius^2) ./ (R(mask).^2)) .* (1 + cos(2*theta(mask)));
V(mask) = -U_inf * (cylinder_radius^2) ./ (R(mask).^2) .* sin(2*theta(mask));
% 绘制流场
figure;
% 绘制流线
streamslice(X, Y, U, V);
hold on;
% 绘制圆柱体
theta_circle = 0:0.1:2*pi;
x_circle = cylinder_radius * cos(theta_circle);
y_circle = cylinder_radius * sin(theta_circle);
fill(x_circle, y_circle, 'k');
axis equal;
xlim([-3 3]);
ylim([-2 2]);
xlabel('X');
ylabel('Y');
title('圆柱体绕流流线图');
colorbar;
运行结果:

案例2:涡旋流场
% 双涡旋相互作用
[X, Y] = meshgrid(-4:0.2:4, -4:0.2:4);
% 涡旋参数
gamma1 = 10; % 涡旋1强度
gamma2 = -8; % 涡旋2强度(反向)
x1 = -1; y1 = 0; % 涡旋1位置
x2 = 1; y2 = 0; % 涡旋2位置
% 计算距离
R1 = sqrt((X-x1).^2 + (Y-y1).^2);
R2 = sqrt((X-x2).^2 + (Y-y2).^2);
% 避免奇点
R1(R1 < 0.1) = 0.1;
R2(R2 < 0.1) = 0.1;
% 计算速度场
U1 = -gamma1 * (Y-y1) ./ (2*pi * R1.^2);
V1 = gamma1 * (X-x1) ./ (2*pi * R1.^2);
U2 = -gamma2 * (Y-y2) ./ (2*pi * R2.^2);
V2 = gamma2 * (X-x2) ./ (2*pi * R2.^2);
U = U1 + U2;
V = V1 + V2;
figure;
% 绘制速度大小的填充图
speed = sqrt(U.^2 + V.^2);
contourf(X, Y, speed, 20);
colorbar;
hold on;
% 叠加流线
streamslice(X, Y, U, V, 2);
xlabel('X');
ylabel('Y');
title('双涡旋相互作用流场');
axis equal;
运行结果:

粒子轨迹追踪
% 在涡旋场中追踪粒子轨迹
% 使用前面定义的涡旋流场
% 初始粒子位置
x0 = [-2, -1.5, 0, 1.5, 2];
y0 = [1, -1, 1.5, 1, -0.5];
figure;
contourf(X, Y, speed, 20);
colorbar;
hold on;
% 为每个粒子计算轨迹
dt = 0.05;
t_final = 10;
time = 0:dt:t_final;
colors = lines(length(x0));
for i = 1:length(x0)
x_traj = zeros(size(time));
y_traj = zeros(size(time));
x_traj(1) = x0(i);
y_traj(1) = y0(i);
for j = 2:length(time)
% 插值获取当前位置的速度
u_interp = interp2(X, Y, U, x_traj(j-1), y_traj(j-1), 'linear', 0);
v_interp = interp2(X, Y, V, x_traj(j-1), y_traj(j-1), 'linear', 0);
% 更新位置(欧拉法)
x_traj(j) = x_traj(j-1) + u_interp * dt;
y_traj(j) = y_traj(j-1) + v_interp * dt;
end
% 绘制轨迹
plot(x_traj, y_traj, 'Color', colors(i,:), 'LineWidth', 2);
plot(x0(i), y0(i), 'o', 'Color', colors(i,:), 'MarkerSize', 8, 'MarkerFaceColor', colors(i,:));
end
xlabel('X');
ylabel('Y');
title('粒子在涡旋场中的轨迹');
axis equal;
xlim([-4 4]);
ylim([-4 4]);
运行结果:

曲面拟合与插值可视化
在实验数据处理和数值分析中,经常需要从散乱的数据点重构连续的曲面。
griddata函数:散点插值
griddata函数可以将散乱分布的数据插值到规则网格上。
基本语法:
Zi = griddata(x, y, z, Xi, Yi, method)
基础示例:随机数据插值
% 生成随机散点数据
n = 50;
x = 4 * rand(n, 1) - 2; % [-2, 2]范围内随机点
y = 4 * rand(n, 1) - 2;
z = sin(x) .* cos(y) + 0.2 * randn(n, 1); % 添加噪声
% 创建规则网格用于插值
[Xi, Yi] = meshgrid(-2:0.1:2, -2:0.1:2);
% 不同插值方法比较
methods = {'linear', 'cubic', 'nearest'};
figure;
for i = 1:3
subplot(2, 3, i);
Zi = griddata(x, y, z, Xi, Yi, methods{i});
surf(Xi, Yi, Zi);
hold on;
plot3(x, y, z, 'ro', 'MarkerSize', 6, 'MarkerFaceColor', 'r');
title([methods{i} ' 插值']);
xlabel('X'); ylabel('Y'); zlabel('Z');
subplot(2, 3, i+3);
contourf(Xi, Yi, Zi, 15);
hold on;
plot(x, y, 'ro', 'MarkerSize', 4, 'MarkerFaceColor', 'r');
colorbar;
title([methods{i} ' 插值等值线']);
xlabel('X'); ylabel('Y');
axis equal;
end
运行结果:

实际应用:地形数据插值
% 模拟GPS测量的地形高程数据
n_points = 100;
x_measure = 20 * rand(n_points, 1); % [0, 20] km
y_measure = 15 * rand(n_points, 1); % [0, 15] km
% 模拟地形:多个山丘
z_measure = 200 + 100 * sin(x_measure/5) .* cos(y_measure/4) + ...
50 * exp(-((x_measure-10).^2 + (y_measure-7).^2)/20) + ...
30 * randn(n_points, 1); % 测量误差
% 创建高分辨率网格
[Xi, Yi] = meshgrid(0:0.5:20, 0:0.5:15);
% 三次插值
Zi = griddata(x_measure, y_measure, z_measure, Xi, Yi, 'cubic');
figure;
subplot(1, 2, 1);
scatter3(x_measure, y_measure, z_measure, 50, z_measure, 'filled');
colorbar;
xlabel('X (km)'); ylabel('Y (km)'); zlabel('高程 (m)');
title('GPS测量点');
subplot(1, 2, 2);
surf(Xi, Yi, Zi);
colorbar;
xlabel('X (km)'); ylabel('Y (km)'); zlabel('高程 (m)');
title('插值重构地形');
shading interp;
运行结果:

scatteredInterpolant类:高效插值
scatteredInterpolant是MATLAB中更现代和高效的散点插值工具,特别适合大数据集。
基本用法:
% 使用scatteredInterpolant进行插值
n = 200;
x = 6 * rand(n, 1) - 3;
y = 6 * rand(n, 1) - 3;
z = peaks(x, y) + 0.1 * randn(n, 1); % 使用peaks函数加噪声
% 创建插值对象
F = scatteredInterpolant(x, y, z, 'linear', 'nearest');
% 生成查询网格
[Xq, Yq] = meshgrid(-3:0.1:3, -3:0.1:3);
Zq = F(Xq, Yq);
figure;
subplot(1, 2, 1);
scatter(x, y, 50, z, 'filled');
colorbar;
xlabel('X'); ylabel('Y');
title('原始散点数据');
axis equal;
subplot(1, 2, 2);
contourf(Xq, Yq, Zq, 20);
hold on;
plot(x, y, 'k.', 'MarkerSize', 4);
colorbar;
xlabel('X'); ylabel('Y');
title('scatteredInterpolant插值结果');
axis equal;
运行结果:

例程文件
上述例程汇总的mlx文件下载链接:
https://download.youkuaiyun.com/download/callmeup/91914095
下载后运行即可得到上面的结果。
如需帮助,或有导航、定位滤波相关的代码定制需求,请点击下方卡片联系作者

被折叠的 条评论
为什么被折叠?



