喷雾机药罐液体纵向晃动等效力学模型建立与应用

陈 雨 路力权 张军昌 陈 军 胡宸玮 曹佳宇

(西北农林科技大学机械与电子工程学院， 陕西杨凌 712100)

0 引言

高地隙自走式喷雾机具有高效、环保等优势，得到了广泛应用[1-3]。但喷雾机运行工况复杂，地隙高，轮距相对于车身高度偏窄[4-5]，作业时药液易发生晃动。液体易于流动的特性和喷雾机运行状态的改变使得药液在外力作用下冲击药罐壁面，产生变化的动态负载，致使作用于整机的力矩增加，极大地降低了整机稳定性与作业质量。因此，研究药罐内液体晃动的动力学特性，对于抑制整机侧倾与俯仰运动，提高喷雾机稳定性与喷药效果具有实际意义[6]。

目前，车辆罐体内液体晃动动力学的研究方法主要有流体动力学法[7-8]、数值解析法[9]、等效力学模型法[10]等。其中，等效力学模型法是指通过力学模型近似描述液体晃动，将流体力学问题转化为机械运动问题。由于其更容易纳入到系统稳定性分析、控制器设计以及系统模拟中，因此得到了广泛的研究。例如KOLAEI等[11]将非满载罐体液体晃动等效力学模型与半挂液罐车五自由度动力学模型进行耦合分析，确定了液罐车侧倾稳定特性。NOKHBATOLFOGHAHAI等[12]、李金辉等[13]将液体晃动弹簧-质量模型整车多体动力学模型结合，分析了液体晃动对整车稳定性的影响。NICOLSEN等[14]和GROSSI等[15]通过集成基于拉格朗日连续介质的液体晃动模型和卡车多体系统模型，研究了液体晃动对罐车动力学特性的影响。郑雪莲等[16-17]、孙文财等[18]和赵伟强等[19]将液罐内液体晃动等效为单摆力学模型，分析不同工况条件下变充液比罐内液体晃动所带来的冲击效应对于罐式半挂车整车侧倾稳定性的影响，并以此提出了基于差动制动的液罐车防侧翻控制策略。

综上所述，液体晃动等效力学模型在罐式道路运输车辆的稳定性分析与控制研究中已广泛应用。与罐体车辆相比，喷雾机质心高，药罐以圆柱和椭圆柱为主，药液与底盘质量比大且作业时药液质量时变，对于果园喷雾机[20-22]及中小型植保机[23-25]，受限于罐体体积，药罐内未放置防波板，药液晃动对其整车稳定性和安全性影响更大。研究构建可准确描述喷雾机液体纵向晃动过程的等效力学模型，可探明整机平顺性及稳定性等影响因素，进而为开展喷雾机底盘结构改进、悬架设计、罐体形状和防晃装置优化、整机平顺稳定控制等提供依据。基于此，本文以喷雾机药液纵向晃动为例，基于液体晃动动力学特性及机械模型等效准则，研究建立药罐液体纵向晃动等效力学模型，并通过仿真与实验方法验证模型的准确性。在此基础上，建立1/2喷雾机四自由度垂向动力学模型，分析讨论液体晃动与药液时变对喷雾机行驶平顺性能的影响。

1 液体纵向晃动等效力学模型建立

1.1 液体纵向晃动特性分析

高地隙自走式喷雾机药罐形状多为圆柱体或椭圆柱体，其纵向截面为矩形。建立罐内液体晃动坐标系如图1所示。取自由液面中心为坐标原点o，x为喷雾机前进方向，z为垂直方向。图中，2a为罐体长度，H为自由液面高度，Ω为罐内流体域。

图1 罐体纵向晃动坐标系Fig.1 Coordinate system of tank longitudinal shaking

当外界激励较小时，罐内液体处于小幅晃动。在如图1所示的平面坐标系oxz中，设h(x,z,t)为自由液面波高函数，引入流体相对速度势函数Φ(x,z,t)，由于小幅晃动下流体可看作是无旋的，则有

(1)

液体粒子在罐体内表面处的法向速度满足

(2)

(3)

同时，自由液面上的液体粒子在z方向上的运动应满足

(4)

(5)

式中ρ——流体密度，kg/m3

p——液体压力，N

g——重力加速度，m/s2

去掉式(5)非线性项，可得

(6)

则罐内液体产生的总压力可写为

(7)

在自由表面上，液体满足的动力学条件为

(8)

式(1)～(4)、(8)构成了罐内液体受迫晃动的完整动力学方程。

在罐内液体小幅晃动时，其一阶模态起主导作用。基于此，在对液体晃动的控制方程求解时不考虑其多阶模态的叠加，仅考虑一阶模态，此时可将液体的速度势函数Φ(x,z,t)展开为

(9)

显然，式(9)满足式(1)～(3)，将其代入式(4)、(8)中得

(10)

(11)

(12)

再将式(12)回代至式(11)中得

(13)

式中ω1——罐内液体一阶晃动自然频率

将式(13)右端的-x展开为三角级数，得到

(14)

由Duhamel积分得到t=0条件下的方程(14)的解为

(15)

由式(15)可得液体的液动压力为

(16)

(17)

(18)

1.2 等效力学模型构建与参数求解

目前，使用较多的等效模型为单摆模型和弹簧-质量-阻尼模型[27]。这两种模型都可以较为准确地描述罐内液体往复冲击的非线性运动，而后者在低频小幅的振动下对液体运动的描述更为精确。因此，基于液体纵向晃动特性分析，建立如图2所示等效力学模型。图中，G为液体质心，hg为质心高度；m0为固定质量块质量，h0为固定质量块高度；mk、kk、ck分别为等效质量块质量、等效弹簧刚度以及等效阻尼系数，hk为等效质量块高度，φ为车身俯仰角。

图2 纵向弹簧-质量-阻尼等效模型Fig.2 Longitudinal spring-mass-damp equivalent model

(19)

进而，根据系统动力学理论，该系统对罐壁在水平方向上产生的力与力矩为[20]

(20)

(21)

式中ωk——k阶液体晃动固有频率

根据前文对液体晃动过程的分析，在求解模型参数时同样只分析其一阶振动模态。根据力学相似原则，罐内实际液体冲击产生的力与力矩与等效模型产生的力与力矩相等，进而根据方程两边时间函数的系数相等可求得等效模型各参数，结果为

(22)

式中M——罐内液体总质量，kg

m1——等效移动质量块质量，kg

(23)

式中h1——等效移动质量块高度，m

(24)

(25)

(26)

本研究中喷雾机药罐截面半径600 mm，长1 000 mm。 为分析不同充液比下液体晃动冲击作用，根据其几何尺寸及式(22)～(26)，以0.1为变化步长，分别求得充液比0.1～0.9时纵向等效力学模型各参数如表1、2所示。

表1 纵向等效模型等效质量和高度Tab.1 Equivalent mass and height of longitudinal equivalent model

表2 纵向等效模型等效刚度和等效阻尼系数Tab.2 Equivalent damping and equivalent stiffness of longitudinal equivalent model

1.3 等效力学数值解析模型建立

纵向等效模型在求解过程中采用数值解析法，通过力学等价原则求解得到模型参数。为后续进一步验证模型的准确性，需结合求得的参数及液体晃动力矩表达式建立等效力学数值解析模型。

根据式(19)～(21)求得在药罐沿纵向晃动时药液产生的作用力矩为

(27)

n——液体晃动模态阶数，取1

结合式(27)在Matlab/Simulink平台中搭建数值仿真模型，如图3所示。

图3 Simulink等效力学模型系统Fig.3 Simulink equivalent mechanics model system

2 等效模型有效性验证

将等效数值解析模型在Matlab/Simulink中的仿真结果与液体晃动在Fluent中的仿真结果进行对比，以验证所建等效模型的准确性。

2.1 Fluent仿真模型建立与验证

在Fluent流体模拟软件中对罐内药液晃动产生的冲击作用进行模拟。按照喷雾机药罐实际尺寸建立药罐模型并标定罐内液体区域如图4所示。设置充液比变化范围为0.1～0.9，变化步长为0.1。

图4 Fluent仿真液罐模型Fig.4 Fluent simulation liquid tank model

根据液体受迫晃动的作用特点设置仿真求解器为基于压力的瞬态流动求解器，罐内同时存在水与空气的情况下，设置其物理模型为VOF多相流模型，湍流模型为标准k-ε模型[22-23]。通过用户自定义函数UDF(Userdefined function)实现对药罐初始激励的加载[24]。仿真时间步长0.01 s，仿真步数500步。

设计并制作了液体晃动实验台架以验证Fluent仿真的准确性，如图5所示，在该平台上可对按照喷雾机药罐等比例制作的小型液罐进行一定频率和幅值的简谐激励。

图5 模拟液体晃动实验台架Fig.5 Experimental bench for simulating liquid shaking1.液罐 2.滑块导轨 3.液罐承托板及固定托 4.高速摄像机 5.控制器显示单元 6.调速电机 7.摄影灯 8.曲柄偏心轮

由于液罐模型为圆柱形，难以直接测得液晃实验中罐内液体对罐壁产生的冲击力，因此，本实验通过在液体晃动过程中监测液面状态的变化，将其与Fluent仿真中的实时液面状态进行对比，以验证其仿真结果的有效性。液面监测采用OLYMPUS公司i-SPEED TR高速摄像机，分辨率为1 280像素×1 024像素，最高帧率达10 000 f/s，适用于绝大多数测试，同时该摄像机配合独特的控制器显示元件(Control and display unit，CDU)来使用。

与仿真实验中初始条件给予的速度激励不同，在液晃模拟实验中给予罐体的是位移激励

Sx=0.02πsin(8πt)

(28)

根据式(28)，调节电机转速及偏心轮偏心距。将获得的高速摄像结果与Fluent实时液面状态对比，图6为充液比α=0.5条件下，不同时刻对比结果，左侧为实验结果，右侧为仿真结果。

图6 不同时刻实验结果与仿真结果对比Fig.6 Comparison of experimental results and simulation results

由图6可知，在相同时刻下，Fluent仿真结果与模拟液晃实验结果高度相似。为量化两者液面匹配程度，借用Matlab中的grabit插件分别提取一个运动周期内不同时刻的液面曲线，如图7所示。

图7 实验与仿真结果液面曲线Fig.7 Comparison of liquid level curves between experimental and simulation results

从图7中可以看出，在同一时刻下两者的液面曲线在变化趋势上基本一致，但在液面上水波高度的峰值存在略微差异。而从直观曲线中难以评价两种结果的液面曲线是否具有较高的拟合度。因此，在评价Fluent结果与高速摄像结果曲线的拟合度时，用correl函数来表示，即

(29)

根据式(29)计算两种结果下液面曲线的相关系数，结果如表3所示。

表3 不同充液比下各时刻的相关系数Tab.3 Correlation coefficient results at each moment under different liquid filling conditions

结合correl函数定义，两组数据相关系数越接近1时，其相关性越强，且呈现正相关。从表3可知，Fluent仿真液面结果与高速摄像结果两种曲线的相关系数在不同充液比时，不同时刻下基本处于0.85以上，说明Fluent仿真能够在一定误差范围内反映液体晃动的真实状态，可进一步提高等效模型验证的合理性。同时相关系数计算结果中亦存在小于0.6，甚至小于0.2的情形，由于进行液晃实验时，液体晃动产生的惯性力带给电机一定的运动阻力，因此，实际的晃动激励可能与期望的原始激励存在一定偏差，造成实验的运动周期相对Fluent仿真存在一定滞后。此外，由于高速摄像过程中液罐处于运动状态，因此存在拍摄角度差异，同时在进行液面曲线提取过程中亦存在一定误差。

2.2 等效模型准确性验证

根据图3数值解析仿真模型在Matlab/Simulink中进行仿真，其中各模块初始值设置为表1、2中的等效模型参数。选择Fixed-step固定步长模式的ode5求解器，设置仿真步长为0.001 s，仿真步数为5 000步。通过step模块与sinewave模块给予系统一定初始激励，在此激励的作用下将力学系统产生的作用力矩结果记录并导出。同时，使用Fluent流体模拟软件对罐内药液晃动产生的力矩进行仿真模拟。

为提高仿真实验验证的精确性，结合喷雾机实际作业工况，选择将阶跃和简谐两种初始激励条件下的仿真结果进行对比。其中，阶跃激励为

(30)

简谐激励为

vx=0.16πcos(8πt)

(31)

图8为阶跃激励下仿真结果对比，图9为简谐激励下仿真结果对比。

图8 阶跃激励下仿真结果对比Fig.8 Comparison of results under step excitation

图9 简谐激励下仿真结果对比Fig.9 Comparison of results under simple excitation

由图8、9可知，在不同初始激励条件下Matlab/Simulink仿真中等效模型数值分析与Fluent流体仿真的药液晃动作用力矩随时间变化规律基本一致，说明所建立的纵向等效力学模型的准确性与合理性。通过分析图中的曲线误差可以得出，在药罐近乎空载即充液比为0.1时，两者误差略大。由于药罐接近空载时，罐内液体晃动的非线性特性较强，两者仿真结果存在一定的误差。此外，在阶跃激励条件下，最大误差出现在激励结束的时刻，由于仿真实验中阶跃激励施加时间短暂，此时系统处于非稳定状态，出现响应滞后的情形，增大了误差。

3 等效模型应用

3.1 1/2喷雾机垂向动力学模型建立

基于上述建立的等效力学模型与车辆动力学模型耦合，建立综合药罐液体晃动和不考虑药罐液体晃动的1/2喷雾机垂向四自由度动力学模型，通过仿真对比，探明罐内药液晃动对喷雾机行驶过程中整机行驶平顺性的影响。

以“地面-轮胎-悬架-车身”构成的系统为对象，建立不考虑药罐液体晃动的1/2喷雾机四自由度垂向动力学模型，如图10所示。传统模型建立过程中，假设喷雾机车身与药罐为刚体，不考虑罐体内液体晃动和车辆姿态引起的轮胎载荷转移的影响。图中zb为质心垂向位移；zF、zR为前后簧载质量的垂向位移；ztF、ztR为前后轮胎垂向位移；zwF、zwR为路面位移激励；msc为半车簧载质量，此时的簧载质量可看作是将罐内液体及车架视为同一刚体的简化结果；Isc为车身转动惯量；ksF、csF、ksR、csR分别为前后悬架等效刚度与等效阻尼系数；ktF、ksR为车轮等效刚度；lf、lr为车身质心至前后轴在地面的投影距离。

图10 不含液体晃动的喷雾机半车动力学模型Fig.10 Dynamics model of sprayer half-car without liquid sloshing

结合图10，根据车身质心处的垂向位移zb和俯仰角φ，动力学方程可表示为

(32)

(33)

(34)

(35)

当俯仰角φ较小时，近似有

zF=zb+lfφ

(36)

zR=zb-lrφ

(37)

结合式(27)、(32)～(37)，建立综合液体纵向晃动的1/2喷雾机非线性垂向动力学模型，如图11所示。

图11 含液体晃动的喷雾机半车动力学模型Fig.11 Dynamic model of sprayer half-car with liquid sloshing

图11中msc1为罐内液体质量，msc2为车架质量，且有msc1+msc2=msc。此时，由于要分析液体晃动的影响，故将液体与车架质量分开考虑。动力学方程表示为

(38)

(39)

(40)

3.2 Simulink动力学仿真

为探究罐内液体晃动对喷雾机行驶稳定性的影响，结合式(32)～(40)在Matlab/Simulink中搭建动力学仿真模型，如图12所示。

图12 Simulink半车动力学仿真模型Fig.12 Simulink half-car dynamics simulation model

仿真参数根据本团队已设计的小型喷雾机设置，半车车架质量msc2为10.3 kg，前后非簧载质量mtF、mtR为15.8 kg，悬架等效刚度ksF、ksR为16 000 N/m，悬架等效阻尼系数csF、csR为500 N·s/m，车轮等效刚度ktF、ktR为20 000 N/m，前车轴与质心投影距离lf为0.75 m，后车轴与质心投影距离lr为0.5 m，喷雾机行驶速度v为1.4 m/s，其中，仿真过程中随着罐内充液比的不同，半车簧载质量msc1、液体质量msc2、车身转动惯量Isc如表4所示。

表4 喷雾机仿真参数Tab.4 Simulation parameters of sprayer

仿真中以生成的随机路面为初始激励信号。结合喷雾机实际作业路况，选择仿真的路面激励模型为E级路面[28-29]。采用一阶滤波带白噪声法建立路面激励时域模型

(41)

式中w(t)——路面位移输入，m

q(t)——限带白噪声

v——喷雾机行驶车速，m/s

nd0——路面空间截止频率,Hz

考虑前后车轮的路面激励存在一定的时间延迟，结合前后轴距及喷雾机行驶速度，生成车轮随机路面激励信号如图13所示。

图13 前后车轮随机路面激励Fig.13 Random road input for front and rear wheels

将喷雾机车身加速度与俯仰角速度作为信号输出，并以此作为评价指标，对比有无液体晃动时的仿真结果，如图14所示。

图14 有无液体晃动时仿真结果对比Fig.14 Comparison of simulation results

由图14可知，当喷雾机在随机农田路面行驶时，由于药罐内液体晃动力作用，其车身加速度与俯仰角速度增大，降低了喷雾机行驶的平顺性。说明喷雾机行驶作业时，若将药罐药液视为刚体，忽略液体晃动的冲击作用时，难以准确分析喷雾机整机运动状态及车身姿态的变化。因此，本文所建立的液体纵向晃动等效力学模型可准确分析喷雾机整机行驶平顺性、稳定性。

为进一步分析喷雾机药液质量时变下液体晃动对整机性能影响，本文对不同充液比条件下的喷雾机垂向动力学特性进行对比分析。图15所示为充液比0.1、0.5、0.9下，喷雾机车身加速度与俯仰角速度的变化情况，并将不同充液比下两评价指标的最大值与均方根进行对比分析，如表5所示。

图15 不同充液比下仿真结果对比Fig.15 Comparison of simulation results at different filling ratios

由图15可知，随着药液质量的增大，液体晃动对喷雾机行驶作业的影响逐渐增大，且整车在路面激励下响应的滞后效应越明显。由表5可知，在相同质量下，若忽略液体晃动影响，将其视为簧载车身质量的一部分时，整车的俯仰角与车身加速度显著降低。对比是否考虑液晃影响下的数据结果，充液质量越大，两者的结果差值越大，进一步说明了液体晃动对整车行驶性能的影响。此外，当充液比为0.8时，整车的俯仰角与车身加速度达到最大值，此时的喷雾机行驶平顺性较差。由于空载与满载(即充液比为0和1)时，不存在罐内液体晃动作用，因此，表5中上述两种特殊情况下是否考虑液体晃动的仿真结果相同。通过以上分析，本文所建立的药罐液体晃动等效模型，在揭示液体晃动作用对喷雾机行驶稳定性影响的同时，也能反映药液时变因素对其所造成的影响，为进一步研究喷雾机不同作业工况和药液质量时变下底盘结构优化及平顺稳定分析提供实施途径和参考方案。

表5 不同充液比仿真结果对比Tab.5 Comparison of simulation results of each filling ratio

4 结论

(1)建立了喷雾机药罐内液体纵向晃动的弹簧-质量-阻尼等效力学模型，并对模型参数进行了求解。应用Fluent建立了罐内液体晃动仿真模型，设计并进行了模拟液体晃动实验，实验结果表明，实际晃动液面曲线与Fluent仿真液面曲线相关系数大于0.85，说明了基于Fluent的液体晃动仿真模型的准确性。

(2)应用Matlab/Simulink软件，建立等效力学数值解析模型，结合Fluent流体仿真结果，对比分析了充液比为0.1、0.5、0.9时液体纵向晃动作用在容器壁的力矩变化规律。对比结果表明，所建立的等效力学模型力矩变化与Fluent中的力矩变化过程一致，验证了所建立的液体纵向晃动等效力学模型的正确性。

(3)建立了综合药罐液体晃动和不考虑药罐液体晃动的半车喷雾机四自由度垂向动力学模型，通过仿真对比分析表明，药罐内液体晃动使车身加速度与俯仰角速度增大，降低了喷雾机行驶的平顺性。通过分析不同充液比条件下的喷雾机垂向动力学特性表明，当充液比为0.8时，整车行驶平顺性较差，在进行转场或喷药作业时应尽量避免充液比为0.8的情况。