💥💥💞💞欢迎来到本博客❤️❤️💥💥
🏆博主优势：🌞🌞🌞博客内容尽量做到思维缜密，逻辑清晰，为了方便读者。
⛳️座右铭：行百里者，半于九十。
📋📋📋本文目录如下：🎁🎁🎁
目录
⛳️赠与读者
👨💻做科研，涉及到一个深在的思想系统，需要科研者逻辑缜密，踏实认真，但是不能只是努力，很多时候借力比努力更重要，然后还要有仰望星空的创新点和启发点。当哲学课上老师问你什么是科学，什么是电的时候，不要觉得这些问题搞笑。哲学是科学之母，哲学就是追究终极问题，寻找那些不言自明只有小孩子会问的但是你却回答不出来的问题。建议读者按目录次序逐一浏览，免得骤然跌入幽暗的迷宫找不到来时的路，它不足为你揭示全部问题的答案，但若能让人胸中升起一朵朵疑云，也未尝不会酿成晚霞斑斓的别一番景致，万一它居然给你带来了一场精神世界的苦雨，那就借机洗刷一下原来存放在那儿的“躺平”上的尘埃吧。
或许，雨过云收，神驰的天地更清朗.......🔎🔎🔎
💥1 概述
文献来源：
本项目的目的是对薄壁铣削的颤振进行建模，并探索其借助MATLAB等计算机软件，研究其对表面粗糙度的影响。考虑到以下因素，此处对颤振进行了模拟，类似于弹簧-质量阻尼系统
两个自由度。阅读并参考了大量研究论文，以找出切割与表面粗糙度之间的关系参数和振动幅度，以及振动幅度与表面之间的关系粗糙度。我们在案例中选择的模型，通过MATLAB表示，采用了轴向切削深度和转速（单位为rpm）作为输入，并给出X方向和Y方向作为输出提供给用户。
2. 颤振
加工振动，技术上称为颤振，是指工件和刀具之间的相对运动工件和切削刀具。这些振动在加工表面上形成波浪状图案表面。这会影响许多典型的加工工艺，如铣削、车削和钻孔和一些非典型加工工艺，如研磨。颤痕是不规则的车轮在磨削过程中出现的不真实表面缺陷，或在转动时留下的常规标记由于加工振动，长工件在车床上发生弯曲。
3. 颤振建模
一般来说，通过求解频率和/或时间上的运动方程来模拟颤振时域，其中切削力和刀具尖端的一些模态动态参数使用。当铣削切削力显示为线性行为。这在传统铣削中非常典型。微细铣削是
其特征在于使用小直径的切削工具（例如，直径在10毫米至1毫米)和小进给率(每齿几微米)。
详细文章见第4部分。
在薄壁铣削过程中，chatter（颤振）振动是一个常见问题，它不仅影响工件的加工质量，特别是表面粗糙度，还可能导致刀具和机床的损坏，降低生产效率。以下是对薄壁铣削中chatter振动的建模及其对工件表面粗糙度影响的研究概述：
一、Chatter振动的建模
- 振动类型与成因
- 薄壁件铣削加工过程中产生的振动主要有两种：强迫振动和自激振动（颤振）。强迫振动主要由周期性铣削力引起，而颤振则是一个闭环系统，由切削过程中的动态不稳定性导致。
- 颤振的形成与多种因素有关，包括切削参数（如主轴转速、进给量、切削深度）、刀具几何参数、工件材料特性以及机床的动态特性等。
- 建模方法
- 动力学模型：基于切削力模型和工件-刀具系统的动态特性，建立薄壁件铣削的动力学模型。该模型通常考虑刀具与工件之间的相互作用力、系统的刚度和阻尼等因素。
- 数值仿真：利用有限元分析（FEA）或多体动力学仿真软件对铣削过程进行数值模拟，以预测颤振的发生和演变过程。
- 实验验证：通过设计实验，采集切削过程中的力信号、加速度信号等，验证模型的准确性和可靠性。
- 特征提取与识别
- 选取力信号和加速度信号作为监测对象，提取实时方差特征和小波能量特征等，用于构建颤振信号识别特征量。
- 设计基于机器学习（如支持向量机SVM）的颤振识别模型，对颤振状态进行分类和预测。
二、Chatter振动对工件表面粗糙度的影响
- 恶化表面质量
- Chatter振动会导致切削过程中的不稳定切削，使切削刃在工件表面产生不规则的运动轨迹，从而在工件表面留下明显的振纹和波纹，增加表面粗糙度。
- 振动还可能引起切削力的波动，进一步影响切削过程中的材料去除率和加工精度。
- 降低加工效率
- Chatter振动会加剧刀具和工件的磨损，缩短刀具寿命，增加更换刀具的频率和停机时间。
- 同时，由于表面质量的恶化，可能需要额外的加工步骤（如抛光）来修复工件表面，从而增加加工成本和时间。
- 控制措施
- 优化切削参数：合理选择主轴转速、进给量和切削深度等参数，以降低切削过程中的动态不稳定性。
- 改善刀具和夹具设计：采用合适的刀具几何形状和夹紧方式，提高系统的刚度和稳定性。
- 应用振动抑制技术：如采用主动振动控制、阻尼减振等措施来抑制颤振的发生。
综上所述，薄壁铣削中的chatter振动是一个复杂且重要的问题，其建模和控制对于提高工件表面质量和加工效率具有重要意义。通过深入研究颤振的成因和机理，建立准确的动力学模型，并采取相应的控制措施，可以有效地抑制颤振的发生，提高薄壁件的加工质量。
本文目标是编写MATLAB代码来模拟颤振现象。对于一个具备两个自由度的模型，这一现象可以用二阶延迟微分方程来表示。
本文的核心旨在通过MATLAB编程语言实现对复杂物理现象——颤振现象的高精度仿真与分析。具体而言，我们将聚焦于一个具有两个自由度的动力学系统模型，此模型被广泛认为是理解和研究颤振特性的理想框架。在此基础上，颤振现象的动态演化将借助一组精心设计的二阶延迟微分方程（DDEs）来精确描述，这类方程因其能有效捕捉到系统状态随时间变化及过去状态影响当前行为的特性，在模拟涉及时间滞后效应的物理问题中展现出显著优势。
本研究不仅限于单纯建立数学模型及编码实现，还将深入探讨如何优化算法参数以确保仿真结果的准确性和计算效率。我们将采用MATLAB强大的内置函数和工具箱，如dde23求解器，来高效解决所提出的延迟微分方程问题，并利用其丰富的可视化功能，直观展示颤振现象的动态过程，包括但不限于响应曲线、相平面图以及可能的频谱分析，以便更深入地洞察系统内部的振动模式与稳定性特征。
此外，为了全面评估所建模型的有效性与实用性，我们将对比模拟结果与已有的实验数据或理论预测，通过误差分析验证模型的准确度，并探讨模型参数变化对颤振特性的影响，为相关工程设计与控制策略提供理论依据与实践指导。总之，本研究工作旨在构建一个可靠、高效的MATLAB仿真平台，不仅促进对颤振机理的深入理解，也为未来在航空航天、机械工程、结构动力学等领域中预防与控制颤振问题提供强有力的工具支持。
📚2 运行结果
我们的结果似乎与理论部分得到的结果一致。振动分布图下图所示为3mm ADOC和2000 rpm。突然出现的山谷和山峰图中的曲线是由于两种不同振动模式的叠加。这些振动模式是
部分代码：
global omega_nx omega_ny xi_x xi_y Kt Kr kx Ky a N Z T theta time theta_entry theta_exit mx my fz f
global f1 c1 k1 f2 c2 k2 f1_sum f2_sum state
Kt = 600; % radial cutting coefficient (from data book) in MPa
Kr = 0.07; % tangential cutting coefficient (from data book)
xi_x = 0.035; % damping ratio in x dirn
xi_y = 0.035; % damping ratio in y dirn
omega_nx = 600; % natural frequency in Hz (x compnt)
omega_ny = 660; % natural frequency in Hz (y compnt)
kx = 5600*10^3; % stiffness in N/m
Ky = 5600*10^3; % stiffness in N/m
mx = 10; % modal mass
my = 10;
a = 3; % axial depth of cut, ADOC in mm
N = 2000; % rpm
Z = 4; % no. of teeth on cutter
T = 60/(N*Z); % tooth passing period
f = 1; % feed
fz = f/Z; % feed per tooth
lags = [T T]; % time lags for DDE equal to tooth passing period
tf = 3; % total time of simulation
theta = 0; % instantaneous theta initialized
theta_entry = 0 *pi/180; % in rads
theta_exit = 90 *pi/180; % in rads
x=0;
y=0;
state = 0; % represents time in steps of tooth passing time
time = linspace(0,tf,10000); % in secs
%% First calculation
f1_sum = 0;
f2_sum = 0;
fprintf('################ Simulation begins ################');
for tooth = 1:Z
theta = 2*pi*N*state/60 + Z*2*pi/N;
fprintf('< %d >-----',tooth);
engagement = g(theta);
f1 = ( Kt*Kr*a*(-sin(theta)) + Kt*a*(-cos(theta)) )*engagement;
f2 = ( Kt*Kr*a*(-cos(theta)) + Kt*a*(sin(theta)) )*engagement;
f1_sum = f1 + f1_sum;
f2_sum = f2 + f2_sum;
end
c1 = -2*xi_x*omega_nx;
k1 = -omega_nx^2;
c2 = -2*xi_y*omega_ny;
k2 = -omega_ny^2;
options = ddeset('Event',@ChtrEvents); % custom trigger event is set
state = state + T
sol = dde23(@ddefunc, lags, @yhist, [0 tf], options);
%% Simulation calculations till final time
while sol.x(end) < tf
state = state + T;
fprintf('\n______________Integration Restart at %5.6f_____________\n', sol.x(end));
fprintf('state value after this cycle = %f \n',state);
f1_sum = 0;
f2_sum = 0;
for tooth = 1:Z
theta = 2*pi*N*state/60 + (tooth-1)*2*pi/Z;
fprintf('< %d >-----',tooth);
engagement = g(theta);
f1 = ( Kt*Kr*a*(-sin(theta)) + Kt*a*(-cos(theta)) )*engagement;
f2 = ( Kt*Kr*a*(-cos(theta)) + Kt*a*(sin(theta)) )*engagement;
f1_sum = f1 + f1_sum;
f2_sum = f2 + f2_sum;
end
c1 = -2*xi_x*omega_nx;
k1 = -omega_nx^2;
c2 = -2*xi_y*omega_ny;
k2 = -omega_ny^2;
% calculating the solution
sol = dde23(@ddefunc, lags, sol, [sol.x(end) tf], options);
%y = deval(sol, t); % t is diif from time
end
fprintf('\n******************Simulation Done******************\n\n');
%% Roughness calculations
Ray = 0; % avg. roughness y dirn
for i = 1:length(sol.x)
itr = i-1;
if itr == 0
prevTime = 0;
else
prevTime = sol.x(1,itr);
end
Ray = Ray + abs(sol.y(1,i))*( sol.x(1,i) - prevTime )/tf ;
end
fprintf('Average Roughness in Y direction = %f mm\n',Ray*1000);
Rax = 0; % avg. roughness x dirn
for i = 1:length(sol.x)
itr = i-1;
if itr == 0
prevTime = 0;
else
prevTime = sol.x(1,itr);
end
Rax = Rax + abs(sol.y(3,i))*( sol.x(1,i) - prevTime )/tf ;
end
fprintf('Average Roughness in X direction = %f mm\n\n',Rax*1000);
Rqy = 0; % rms roughness y dirn
for i = 1:length(sol.x)
itr = i-1;
if itr == 0
prevTime = 0;
else
prevTime = sol.x(1,itr);
end
Rqy = Rqy + ( (sol.y(1,i))^2 )*( sol.x(1,i) - prevTime )/tf ;
end
Rqy = sqrt(Rqy);
fprintf('Root mean square Roughness in Y direction = %f mm\n',Rqy*1000);
Rqx = 0; % rms roughness x dirn
for i = 1:length(sol.x)
itr = i-1;
if itr == 0
prevTime = 0;
else
prevTime = sol.x(1,itr);
end
Rqx = Rqx + ( (sol.y(3,i))^2 )*( sol.x(1,i) - prevTime )/tf ;
end
Rqx = sqrt(Rqx);
fprintf('Root mean square Roughness in X direction = %f mm\n\n',Rqx*1000);
Rty = max(sol.y(1,:)) - min(sol.y(1,:));
Rtx = max(sol.y(3,:)) - min(sol.y(3,:));
🎉3 参考文献
文章中一些内容引自网络，会注明出处或引用为参考文献，难免有未尽之处，如有不妥，请随时联系删除。
[1]王学彬.微铣削表面粗糙度影响因素与刀具失效机制研究[D].哈尔滨工业大学,2015.
[2]丁悦,王焱,吴茂庆,等.铝锂合金高速铣削表面粗糙度试验研究[J].航空制造技术, 2015(6):92-94,99.
🌈4 Matlab代码、文档、PPT下载
资料获取，更多粉丝福利，MATLAB|Simulink|Python资源获取
扫一扫
536
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
