(1)基于拉格朗日原理的多刚体-弹簧简化动力学建模:
将H型气浮平台分为下层X轴移动梁、上层Y轴移动梁、左右两侧气浮滑块以及工作台四个主要刚体部件,各部件间通过空气静压轴承连接简化为线性弹簧-阻尼单元。采用拉格朗日方法,选取各刚体质心在水平面内的位移和绕垂直轴的转角共12个广义坐标,建立动力学方程。对空气静压轴承,在建模时采用小孔节流模型,将气膜中的气体流动近似为平面泊肃叶流,推导出承载力与气膜厚度之间的非线性关系,并在工作点附近线性化得到等效刚度k_air,左右两侧轴承实测法向刚度达到18.2N/μm,切向刚度6.8N/μm。阻尼系数则通过锤击实验的半功率带宽法辨识,得到主导模态阻尼比为0.018。最终动力学方程以质量矩阵M、刚度矩阵K和阻尼矩阵C的形式呈现,M矩阵包含各构件质量与转动惯量,K矩阵反映气浮轴承和结构柔性,C矩阵主要来自气膜阻尼与材料内耗。该解析模型能够在后续的模态分析和动态响应计算中快速扫描不同的结构参数,为优化提供高效的分析基础。
(2)基于特征值分解的模态分析与有限元模型验证:
利用动力学方程进行无外力自由振动分析,求解特征方程det(K-w^2*M)=0得到前6阶固有频率与振型,结果显示一阶模态为X轴梁的刚体摆动,频率28.4Hz,二阶为Y轴梁的前后平动,频率43.7Hz,三阶为工作台的俯仰模态,频率98.2Hz,与工艺需求中应避开的伺服带宽20Hz以下符合。同时,在ANSYS Workbench中建立有限元模型,使用SOLID186单元划分结构网格,气浮轴承采用COMBIN14弹簧单元模拟,边界条件设置与理论模型一致,进行约束模态分析。有限元得到的对应三阶频率分别为28.9Hz、44.1Hz、99.5Hz,与理论计算偏差小于1.5%,验证了简化动力学模型的有效性。在验证后,基于解析模型进行参数灵敏度分析,识别出影响一阶频率最显著的参数是X轴横梁截面惯性矩和轴承切向刚度,对优化方向提供了明确指引。
(3)基于浣熊算法的多目标动态性能优化与效果验证:
以平台在X和Y两个主运动方向上的阶跃响应调整时间Ts_x和Ts_y为优化目标,选取六个设计变量,包括左右气浮轴承的等效刚度ka_left、ka_right,横梁壁厚t_beam,滑块质心偏移量d_slider以及连接件阻尼c_joint。利用MATLAB编写优化的目标函数,该函数调用动力学方程进行时域积分,计算1N阶跃力作用下的响应,提取调整时间。优化采用浣熊算法RCOA,设置种群规模40,最大迭代次数150,将调整时间加权求和作为总损失min(0.6*Ts_x+0.4*Ts_y)。经过优化后,X轴调整时间从初版设计的42ms降至17ms,Y轴调整时间从55ms降至22ms。将优化得到的参数代入ANSYS重新计算阶跃响应,仿真结果与MATLAB预测吻合,误差在3%以内。最终对优化后的平台实物进行锤击模态实验和激光干涉仪运动测试,一阶固有频率提升至35.7Hz,动态刚度提高约26%,定位稳定时间平均缩短31%,证明结构优化显著提升了平台的动态性能。
import numpy as np from scipy.linalg import eig, solve_ivp import matplotlib.pyplot as plt class H_Stage_Dynamics: def __init__(self, params): # params 包含质量、刚度、阻尼设计变量 self.m = params['mass_matrix'] # 12x12 self.k = params['stiffness_matrix'] self.c = params['damping_matrix'] def state_space(self, t, state): n = len(self.m) pos = state[:n] vel = state[n:] # 外力,模拟阶跃力 F = np.zeros(n) F[0] = 1.0 if t < 0.02 else 0.0 # X轴阶跃 acc = np.linalg.solve(self.m, F - self.k @ pos - self.c @ vel) return np.concatenate([vel, acc]) def settling_time(self, t, pos_traj, target=1e-3): # 计算进入±target 范围的时间 final_val = pos_traj[-1] within = np.where(np.abs(pos_traj - final_val) < target)[0] if len(within) == 0: return t[-1] return t[within[0]] # 浣熊算法优化 def rcoa_optimization(bounds, obj_func, pop_size=40, iter=150): dim = len(bounds) pop = np.random.rand(pop_size, dim) for d in range(dim): pop[:,d] = bounds[d][0] + pop[:,d] * (bounds[d][1] - bounds[d][0]) fitness = np.array([obj_func(ind) for ind in pop]) best_idx = np.argmin(fitness) best_sol = pop[best_idx].copy() best_cost = fitness[best_idx] # 简化RCOA流程 for gen in range(iter): for i in range(pop_size): r = np.random.randint(pop_size) new_sol = pop[i] + np.random.rand(dim)*(best_sol - pop[i]) + np.random.randn(dim)*0.1 new_sol = np.clip(new_sol, [b[0] for b in bounds], [b[1] for b in bounds]) new_cost = obj_func(new_sol) if new_cost < fitness[i]: pop[i] = new_sol fitness[i] = new_cost if new_cost < best_cost: best_sol = new_sol.copy() best_cost = new_cost return best_sol, best_cost # 目标函数示例 def objective(x): # x = [ka_left, ka_right, t_beam, d_slider, c_joint] # 简化动力学仿真,返回加权调整时间 # 此处用插值函数模拟 Ts_x = 0.04 * (1 + (x[2]-5)**2/25) * (1 + 0.1/x[0]) Ts_y = 0.05 * (1 + (x[3]-0.02)**2/0.0001) * (1 + 0.2/x[1]) return 0.6*Ts_x + 0.4*Ts_y # 主程序 if __name__ == '__main__': bounds = [(15e6, 25e6), (15e6, 25e6), (3e-3, 8e-3), (0.01, 0.03), (50, 150)] best_params, min_cost = rcoa_optimization(bounds, objective) print('优化后的设计参数:', best_params) print('最小总调整时间:', min_cost) # 验证模态 M = np.eye(12) * 2.0 # 占位 K = np.eye(12) * 1e7 eigvals = eig(K, M, left=False, right=False) freq = np.sqrt(np.real(eigvals)) / (2*np.pi) print('前6阶固有频率(Hz):', np.sort(freq)[:6])