第三节:用NumPy与Pandas处理力学仿真数据
——从力学公式到可训练的数据集
本节摘要
第二节中,我们完成了Python开发环境搭建,并绘制了基本的应力—应变曲线。本节将进一步学习NumPy和Pandas在力学仿真中的实际用法。
我们将以均匀直杆轴向拉伸为例,批量计算钢、铝合金和钛合金在不同载荷下的应力、应变与位移,并将结果整理成结构化数据表。随后,通过数据检查、分组统计、图形对比和CSV导出,建立一条完整的力学数据处理流程。
一、本节学习目标
完成本节后,你将能够:
使用NumPy数组保存和计算多组力学工况;
理解向量化计算与普通循环的区别;
使用Pandas建立力学仿真数据表;
检查数据中的缺失值、重复值和物理异常;
按材料类型进行分组统计;
绘制载荷—位移曲线和材料对比图;
将计算结果导出为机器学习可读取的CSV数据集。
二、为什么要学习力学数据处理?
有限元软件最终会产生大量数据,例如:
不同节点的位移;
不同积分点的应力和应变;
不同载荷步的响应;
不同材料参数下的计算结果;
不同裂纹长度对应的J积分;
不同工况下的最大应力和最大位移。
如果只计算一个模型,我们可以手动查看结果。
但如果需要计算100种材料参数、50种载荷和20种几何尺寸,就可能产生:
组仿真结果。
这时,逐个打开结果文件已经不现实。我们需要让Python自动完成:
图3-1 力学仿真数据处理流程
这套流程也是后续建立机器学习代理模型的数据基础。
三、工程案例:三种材料的轴向拉伸
假设有一根均匀直杆,其参数如下:
| 参数 | 符号 | 数值 |
|---|---|---|
| 杆长 | L | 1000 mm |
| 截面积 | A | 100 mm^2 |
| 载荷 | F | 2000~10000 N |
| 钢的弹性模量 | E_s | 210000 MPa |
| 铝合金弹性模量 | E_a | 70000 MPa |
| 钛合金弹性模量 | E_t | 110000 MPa |
本案例统一采用N、mm和MPa单位制。
因为:
所以在公式中可以直接使用上述数据,不需要额外进行单位换算。
1. 应力计算公式
均匀直杆的轴向应力为:
2. 应变计算公式
在线弹性范围内,根据胡克定律:
3. 端部位移公式
均匀直杆的轴向位移为:
也可以写成:
我们的任务是计算三种材料在五种载荷下的力学响应,共得到:
组数据。
四、使用NumPy进行批量计算
1. 导入计算库
import numpy as np np.set_printoptions(precision=4, suppress=True) print("NumPy版本:", np.__version__)预期输出
版本号会随安装环境略有不同,例如:
NumPy版本: 2.2.6只要程序能够正常显示版本号,就说明NumPy已经安装成功。
2. 建立载荷数组
loads = np.array( [2000, 4000, 6000, 8000, 10000], dtype=float ) print("载荷数组:", loads) print("工况数量:", loads.size)预期输出
载荷数组: [ 2000. 4000. 6000. 8000. 10000.] 工况数量: 5这里的loads不是一个普通数字,而是包含五个载荷工况的NumPy数组。
3. 批量计算钢杆的应力、应变和位移
length = 1000.0 area = 100.0 elastic_modulus = 210000.0 stress = loads / area strain = stress / elastic_modulus displacement = strain * length print("应力(MPa):", stress) print("应变:", strain) print("位移(mm):", displacement)预期输出
应力(MPa): [ 20. 40. 60. 80. 100.] 应变: [0.0001 0.0002 0.0003 0.0004 0.0005] 位移(mm): [0.0952 0.1905 0.2857 0.3810 0.4762]为了方便显示,程序保留了四位小数,因此应变显示值有所简化。以最大载荷为例,其准确结果为:
0.00047619
对应端部位移为:
0.47619 mm
五、什么是向量化计算?
如果不使用NumPy,我们可能需要编写循环:
stress_list = [] for load in loads: value = load / area stress_list.append(value) print(stress_list)预期输出为:
[20.0, 40.0, 60.0, 80.0, 100.0]使用NumPy后,只需要一行代码:
stress = loads / area这叫作向量化计算。
它的特点是:
代码更加简洁;
更接近力学公式;
处理大量数据时速度更快;
可以减少循环和下标错误。
可以把NumPy数组想象成“一批同时参与计算的工况”。
六、使用Pandas建立力学数据表
NumPy擅长数值计算,但不适合为每一列数据添加明确的工程含义。Pandas可以将计算结果整理成带有行、列标签的数据表。
1. 批量计算三种材料
import numpy as np import pandas as pd length = 1000.0 area = 100.0 loads = np.array( [2000, 4000, 6000, 8000, 10000], dtype=float ) materials = { "钢": 210000.0, "铝合金": 70000.0, "钛合金": 110000.0 } records = [] for material, elastic_modulus in materials.items(): stress = loads / area strain = stress / elastic_modulus displacement = strain * length for load, s, e, u in zip( loads, stress, strain, displacement ): records.append({ "材料": material, "弹性模量_MPa": elastic_modulus, "载荷_N": load, "应力_MPa": s, "应变": e, "端部位移_mm": u }) data = pd.DataFrame(records) print(data.head().to_string(index=False)) print("\n数据形状:", data.shape)预期输出
材料 弹性模量_MPa 载荷_N 应力_MPa 应变 端部位移_mm 钢 210000.0 2000.0 20.0 0.000095 0.095238 钢 210000.0 4000.0 40.0 0.000190 0.190476 钢 210000.0 6000.0 60.0 0.000286 0.285714 钢 210000.0 8000.0 80.0 0.000381 0.380952 钢 210000.0 10000.0 100.0 0.000476 0.476190 数据形状: (15, 6)(15, 6)表示这张数据表共有:
15行数据,即15组材料—载荷组合;
6列变量,即材料、弹性模量、载荷、应力、应变和位移。
七、理解DataFrame的结构
在机器学习中,一般把数据分为输入和输出。
本案例中的输入特征可以写成:
输出响应可以写成:
也就是说,模型需要学习下面的映射:
目前使用的是理论公式。后续课程中,我们会让机器学习模型从数据中自动学习这种关系。
八、进行数据质量检查
在真实的Abaqus数据中,可能出现:
某个作业没有收敛,结果为空;
同一工况被重复计算;
单位使用错误;
位移或应力超出合理范围;
ODB结果提取不完整。
因此,建模前必须先检查数据。
1. 检查缺失值
print(data.isna().sum())预期输出
材料 0 弹性模量_MPa 0 载荷_N 0 应力_MPa 0 应变 0 端部位移_mm 0 dtype: int64所有列均为0,表示没有缺失值。
2. 检查重复数据
duplicate_count = data.duplicated().sum() print("重复行数量:", duplicate_count)预期输出
重复行数量: 03. 检查物理异常
对于当前的拉伸问题,弹性模量、载荷和位移都应为正值。
invalid_data = data[ (data["弹性模量_MPa"] <= 0) | (data["载荷_N"] <= 0) | (data["端部位移_mm"] < 0) ] print("物理异常数据数量:", len(invalid_data))预期输出
物理异常数据数量: 0这里的判断并不是单纯检查数据格式,而是把力学知识加入数据检查过程。
这正是“物理感知数据处理”的基本思想。
九、按照材料进行分组统计
现在,我们想知道在给定载荷范围内,哪种材料的变形最小。
summary = ( data.groupby("材料", as_index=False) .agg( 最大位移_mm=("端部位移_mm", "max"), 平均位移_mm=("端部位移_mm", "mean") ) .sort_values("最大位移_mm") ) print( summary.to_string( index=False, float_format=lambda x: f"{x:.6f}" ) )预期输出
材料 最大位移_mm 平均位移_mm 钢 0.476190 0.285714 钛合金 0.909091 0.545455 铝合金 1.428571 0.857143从结果可以看出:
在几何和载荷相同的情况下,弹性模量越大,位移越小。
因为:
当F、L和A保持不变时:
这说明数据统计结果与力学理论一致。
十、绘制三种材料的载荷—位移曲线
import matplotlib.pyplot as plt plt.rcParams["font.sans-serif"] = [ "Microsoft YaHei", "SimHei" ] plt.rcParams["axes.unicode_minus"] = False colors = { "钢": "#1769aa", "铝合金": "#f28e2b", "钛合金": "#00a7a0" } for material in materials: subset = data[data["材料"] == material] plt.plot( subset["载荷_N"] / 1000, subset["端部位移_mm"], marker="o", linewidth=2.5, label=material, color=colors[material] ) plt.xlabel("载荷 F(kN)") plt.ylabel("端部位移 u(mm)") plt.title("相同几何与载荷下的材料位移对比") plt.grid(alpha=0.3) plt.legend() plt.tight_layout() plt.show()预期绘图结果
程序会生成三条近似直线:
钢的曲线斜率最小;
铝合金的曲线斜率最大;
钛合金位于两者之间;
三条曲线均随载荷线性增加。
图3-2 三种材料的载荷—位移曲线
注意:曲线更陡并不意味着材料“更强”,而只表示在当前线弹性模型下,它的轴向刚度较低、变形更大。
材料强度和材料刚度不是同一个概念:
弹性模量主要反映刚度;
屈服强度反映开始发生塑性变形的能力;
极限强度反映材料能够承受的最大应力水平。
十一、绘制最大载荷下的位移对比图
max_load_data = data[ data["载荷_N"] == data["载荷_N"].max() ] plt.bar( max_load_data["材料"], max_load_data["端部位移_mm"], color=[ colors[name] for name in max_load_data["材料"] ] ) plt.xlabel("材料") plt.ylabel("端部位移 u(mm)") plt.title("10 kN载荷下的端部位移") plt.grid(axis="y", alpha=0.3) plt.tight_layout() plt.show()预期绘图结果
在10 kN载荷下:
| 材料 | 端部位移 |
| 钢 | 0.4762 mm |
| 铝合金 | 1.4286 mm |
| 钛合金 | 0.9091 mm |
图3-3 最大载荷下的材料位移对比
十二、将结果保存为CSV文件
整理完成的数据可以保存为CSV文件,供后续机器学习模型读取。
from pathlib import Path output_dir = Path("outputs") output_dir.mkdir(exist_ok=True) file_path = output_dir / "bar_simulation_data.csv" data.to_csv( file_path, index=False, encoding="utf-8-sig" ) print("数据已保存到:", file_path.resolve())预期输出
数据已保存到: D:\你的项目目录\outputs\bar_simulation_data.csv具体路径会根据项目保存位置有所不同。
encoding="utf-8-sig"可以减少CSV文件在中文版Excel中打开时出现中文乱码的概率。
十三、重新读取并验证数据
保存文件后,不要立即认为任务已经完成。一个可靠的数据流程还应重新读取文件并检查数据是否完整。
loaded_data = pd.read_csv(file_path) print("读取后的数据形状:", loaded_data.shape) print("列名:", loaded_data.columns.tolist())预期输出
读取后的数据形状: (15, 6) 列名: ['材料', '弹性模量_MPa', '载荷_N', '应力_MPa', '应变', '端部位移_mm']如果保存前后都是15行、6列,并且列名一致,就说明数据导出基本正确。
十四、从这个案例可以发现什么?
1. 应力为什么与材料无关?
当前模型是截面积相同的均匀直杆,并且受到相同的轴向载荷,因此:
应力只由载荷和截面积决定,与弹性模量无关。
所以三种材料在同一载荷下的应力相同。
2. 位移为什么与材料有关?
位移公式为:
弹性模量$E$越大,杆件越不容易变形,因此位移越小。
3. 为什么三条曲线都是直线?
因为当前使用的是线弹性模型:
并且$E$保持不变,所以载荷与位移呈线性关系。
如果材料进入塑性阶段,或者存在接触、损伤、大变形等非线性因素,曲线就可能不再是直线。
十五、常见错误与排查方法
错误一:找不到Pandas或Matplotlib
ModuleNotFoundError: No module named 'pandas'这通常表示当前Python环境没有安装对应的库,或者VS Code选择了错误的解释器。
可以先检查:
python --version python -m pip --version然后安装:
python -m pip install numpy pandas matplotlib错误二:图中的中文显示为方框
可以在绘图前加入:
plt.rcParams["font.sans-serif"] = [ "Microsoft YaHei", "SimHei" ] plt.rcParams["axes.unicode_minus"] = False错误三:位移结果相差1000倍
这通常是单位混用造成的。例如:
长度使用m;
截面积使用
;
弹性模量使用Pa;
载荷使用N。
进行力学计算前,应先统一单位体系。
错误四:CSV中的中文乱码
保存时使用:
encoding="utf-8-sig"十六、本节练习
基础练习
将载荷数组修改为:
loads = np.arange(1000, 11000, 1000)观察数据量和曲线发生了什么变化。
进阶练习
增加一种材料:
materials["铜合金"] = 120000.0重新生成数据表和对比曲线。
思考练习
如果保持材料和载荷不变,将截面积从:
修改为:
请先不要运行程序,预测应力和位移将如何变化,再用代码验证。
理论上:
当截面积减小一半时,应力和位移都将增大为原来的两倍。
十七、本节小结
本节完成了一个小型但完整的力学数据处理项目:
需要记住的不是某一行代码,而是这条完整的数据工作流。
在后面的有限元项目中,理论公式会被Abaqus或FEniCS求解器取代,但数据处理思想不会改变:
当这样的数据积累到一定规模后,就可以训练机器学习模型,让模型快速预测位移、应力、应变甚至J积分。
下一节将进入线性回归,学习如何让机器从这些数据中自动寻找输入参数与力学响应之间的数学关系。