如何在CFD速度场contourf插值中屏蔽固体区域零值

问题背景 #

从Ansys Fluent导出包含x、y、z（非规则网格节点坐标）和u（流向速度幅值）的CSV数据，需绘制x-z平面的contourf等高线图，要求屏蔽u=0的固体区域，生成带不规则边界的等高线图（匹配Ansys CFD-Post的输出效果）。

数据参数与预处理代码 #

### 几何参数（单位：m）
L = 0.3048 
hwall = 0.0095
Lin = L - hwall * 2
H = 0.073 
Lpad = Lin / 10
Hpad = H / 10

### 读取数据并截取目标区域
import pandas as pd
df = pd.read_csv('test2.csv', names=['x', 'y', 'z', 'u'])
dfsub = df.loc[(df.x > -Lin/2 - Lpad) & (df.x < Lin/2 + Lpad) & (df.z < H + Hpad)]

已尝试的无效方案 #

• 将u=0的数值替换为NaN：触发ValueError: z array must not contain non-finite values within the triangulation错误
• 直接删除u=0的数据行：插值时自动填充固体区域，无法保留不规则边界

可行解决方法 #

方案1：三角剖分+掩码过滤固体区域 #

核心逻辑：仅对流体区域的点做Delaunay三角剖分，通过掩码排除覆盖固体区域的三角形

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.tri import Triangulation

# 分离流体区域（用极小值替代0，规避浮点误差）
df_fluid = dfsub[dfsub.u > 1e-6]
x_fluid = df_fluid.x.values
z_fluid = df_fluid.z.values
u_fluid = df_fluid.u.values

# 生成三角剖分
tri = Triangulation(x_fluid, z_fluid)

# 定义固体区域判断函数（需匹配实际几何边界）
def is_solid_region(x_c, z_c):
    # 示例：根据给定几何参数定义固体范围
    return (x_c < -Lin/2) | (x_c > Lin/2) | (z_c > H)

# 计算每个三角形的质心，生成掩码
x_centroid = np.mean(x_fluid[tri.triangles], axis=1)
z_centroid = np.mean(z_fluid[tri.triangles], axis=1)
tri_mask = is_solid_region(x_centroid, z_centroid)
tri.set_mask(tri_mask)

# 绘制等高线图
plt.figure(figsize=(10, 6))
contour = plt.tricontourf(tri, u_fluid, levels=20, cmap='viridis')
plt.colorbar(contour, label='流向速度u (m/s)')
plt.xlabel('x (m)')
plt.ylabel('z (m)')
plt.title('x-z平面流向速度等高线（屏蔽固体区域）')
plt.gca().set_aspect('equal')
plt.show()

方案2：规则网格插值+手动掩码标记 #

核心逻辑：先对目标区域生成规则网格，用流体区域数据插值后，手动将固体区域设为NaN

import numpy as np
import matplotlib.pyplot as plt
from scipy.interpolate import griddata

# 分离流体区域
df_fluid = dfsub[dfsub.u > 1e-6]

# 生成目标规则网格
x_grid = np.linspace(dfsub.x.min(), dfsub.x.max(), 500)
z_grid = np.linspace(dfsub.z.min(), dfsub.z.max(), 500)
X, Z = np.meshgrid(x_grid, z_grid)

# 用流体区域数据完成插值
u_grid = griddata((df_fluid.x, df_fluid.z), df_fluid.u, (X, Z), method='cubic')

# 手动标记固体区域为NaN（匹配实际几何边界）
solid_mask = (X < -Lin/2) | (X > Lin/2) | (Z > H)
u_grid[solid_mask] = np.nan

# 绘制等高线图
plt.figure(figsize=(10, 6))
contour = plt.contourf(X, Z, u_grid, levels=20, cmap='viridis')
plt.colorbar(contour, label='流向速度u (m/s)')
plt.xlabel('x (m)')
plt.ylabel('z (m)')
plt.title('x-z平面流向速度等高线（屏蔽固体区域）')
plt.gca().set_aspect('equal')
plt.show()

关键注意事项 #

• 避免直接用u == 0判断，改用u > 1e-6过滤浮点误差导致的极小值
• 固体区域的掩码定义必须与实际几何边界完全匹配，否则会出现边界错误
• 三角剖分方案更适配非规则网格的原始分布，规则网格插值方案适合生成规整的可视化结果

内容的提问来源于stack exchange，提问作者lsmontal