杜小振, P. A. Mbango-Ngoma, 常 恒, 张 咪, 王 宇
(山东科技大学 机械电子工程学院,山东 青岛 266590)
压电振动能量收集获得广泛研究,拓展激振响应频率可有效提高压电发电效率和输出功率[1]。流致涡激振动能量采集[2]和优化流场钝体,使漩涡脱落频率与压电换能结构固有频率一致。钝体结构和流场结构设计可以有效提高压电发电输出能量。研究方法包括仿真建模、理论及实验等。Wu等[3]采用OpenFOAM软件建模方法模拟了雷诺数为500时与钝体刚度连接的涡激控制圆杆间隙比(钝体表面间隙与钝体直径比)及其对称布置分布角度等产生的涡激振动(VIV)效果。Hu G等[4]在风能采集器圆柱钝体上下两侧迎风面安装横截面形状为圆形、三角形和正方形涡激调节元件,在迎风面安装三角形结构呈120°角度时能有效地提高涡激横向风力。王世龙等[5]仿真分析振荡水柱气室排气风筒内涡激压电悬臂梁振动收集海洋能,能够将低频的波浪振动转换为高频风激压电悬臂梁振动,输出电压达到6.1 V,有效的提高了压电发电效率。Hu等[6]提出一种涡激振动压电风速传感器,采用有限元方法分析了功率谱与几何参数(压电元件尺寸和位置)之间的关系。Wang等[7]设计Y形高性能压电风能采集器,基于实验和LBM(Lattice-Boltzmann CFD method)计算流体动力学方法分析振幅和频率,证实了从涡激振动到驰振的转变。Dai等[8]实验研究钝体方向(底部、顶部、水平和垂直)与低风速和高风速(1~20 m/s)范围压电采集器的固有频率、阻尼、同步区等对输出功率的影响。Su等[9]提出采用U形梁、一对压电片和一个附着在梁中心的泡沫圆柱体组成的涡激双向压电能量采集器。宋汝君等[10]提出了一种压电悬臂梁和末端圆柱体组成的水流涡激振动复摆式压电俘能器,采用流-固-电耦合仿真方法,分析了较低水流流速下实现高频涡激共振特点。赵道利等[11]设计了一种适用于低流速水流中的悬臂梁式压电能量收集器,利用明渠流弯道水槽研究了五种截面形状钝体质量块结构的压电能量收集器,结果表明截面为三棱柱钝体质量块结构输出性能较好。Sun等[12]提出一种基于自调谐涡激振动压电能量收集系统,在悬臂梁上配置可滑动钝体,扩大了激振频率锁定范围。实验结果表明钝体高度0.12 m时,输出平均功率增加110%。Wang等[13]提出上下耦合式双晶压电悬臂梁结构,风洞模拟试验研究不同风荷载宽频发电获得最大输出电能为56.64 μW。王定标等[14]采用等效电路法对变三角截面驰振压电能量俘获器进行建模,借助风洞实验对钝体顶角分别为30°,60°,90°下的驰振压电俘能器的发电性能进行研究,实验结果表明钝体顶角为90°时发电性能最优。Hu等[15]分析了涡流脱落中压电能量采集区域的理论模型,通过数值计算和实验发现流体速度和钝体直径优化配合有助于提高系统发电性能。以上研究中,主要采用涡激尾流区域压电悬臂梁结构振动发电研究。涡激振动空气动力能量分布区域研究较少,多以尾流区压电悬臂梁发电结构为主,没有涉及涡激风洞管壁圆形压电发电研究。本文主要分析了钝体结构、尾流区域空气动力和管壁分布圆形压电阵子激振发电特性等,拓展了涡激压电发电有效利用模式。
涡激振动流场模型如图1。气流空气动力为流场动力源,流经钝体绕流产生交替脱落的漩涡,进而引起置于涡激振动尾流区流场中的压电薄膜在脱离漩涡交变压强差作用下激振形变发电。为了分析尾流区涡街脱落频率和漩涡强度,对比流场稳定性并探究压电薄膜的最佳取压位置,在钝体尾流区不同位置设置监测点,获取流场流动时压力周期性变化特征。本文采用结构化网格划分方法并在距钝体不同位置设置了5个监测点,如图1(a),B为钝体迎风面截面宽度(其中Bt、Bs、Br分别为三角形、方形、圆形钝体迎风截面边长或直径),W为流场宽度,涡街脱离区域流场监测点P0位于钝体右后方距离钝体右边界0.5B处,然后每隔0.5B的距离设置一个监测点,从左到右依次为P0~P4,提取监测点流场数据,采用快速傅里叶变换(FFT)求得流场中的涡街脱落频率和幅度,分析不同监测点位置的流场特征,确定压电换能能采集元件最佳位置。流场涡激压力幅值大小即涡街强弱,对压电薄膜发电性能影响较大。为此,对比分析了三角形、方形和圆形钝体尾迹漩涡脱落特性,钝体附近网格划分如图1(a),考虑到钝体周围流场变化较为剧烈,其周围网格做加密处理。通过二维流场仿真分析,得到圆形钝体涡激振动流场流线如图1(b),当流体绕流钝体形成周期性脱离漩涡,引起压力波动驱动压电薄膜振动发电。流场二维分析结果表明P1位置(距钝体尾缘1B)涡激产生气流压差最大。为了更准确地模拟压电薄膜在流场中的发电特性,以二维仿真分析为基础,建立三维流场模型,以三角形钝体为例的三维流场网格划分如图1(c),L、W、H分别为计算流体动力学(CFD)计算域长、宽、高,PVDF压电薄膜中心放置在距钝体尾缘1B的管壁上。计算域上下两壁面为无滑动壁面边界,左侧为空气速度入口边界,右侧为压力出口边界。
(1)
(2)
式中:N是按时间序列排列的结果数;ρ为空气密度;u为流体速度;B为迎风面截面宽度;fd(t)和fL(t)分别是波动的阻力和升力。
图1 涡激压电发电风能流场模型
图2 网格独立性验证
目前广泛采用k-ε湍流模型,标准k-ε模型计算精度,适合高雷诺数湍流,但仅适用于完全湍流流场模拟。为了提高精度和考虑涡流漩涡现象等,提出了各种改进的k-ε模型,RNGk-ε模型和Realizablek-ε是最具代表性的模型。与标准k-ε模型相比,RNGk-ε模型修正了湍动粘度并考虑了流体漩转情况,提高了计算精度,但不适用于近壁区和湍流雷诺数较低场合。Realizablek-ε模型为湍流黏性增加了一个公式,对平板和圆柱射流发散率计算结果更加精确,且对于漩转均匀剪切流、自由流、边界层流动和腔道流动等模拟效果好。据此,本文选择Realizablek-ε湍流模型进行流场仿真。本文的流场流体介质为空气,密度1.225 kg/m3,动力黏度1.789 4×105Pa·s。分别定义7个入口速度:2.6 m/s、4.8 m/s、7.5 m/s、9 m/s、11 m/s、14.7 m/s和18.8 m/s。为了确保计算过程稳定性,选择压力-速度耦合求解器,对应算法有SIMPLE、SIMPLEC和PISO三种,其中SIMPLE和SIMPLEC常用于计算定常流动,而且SIMPLEC算法稳定性好、求解速度快,适用于多种计算工况。PISO算法适用于非定常流计算和网格偏斜较严重的模型,但该算法计算量大,算时较长。本文选用SIMPLEC算法进行求解,依据本算例的网格划分形式及计算精度,选择二阶迎风离散格式模拟流场。
气流经相同钝体结构,虽然风速变化,在尾流区流场涡街脱落涡漩、压力及速度场的分布状态相似,但脱落频率与振幅存在差别,影响压电单元发电结构发电特性。以风速u=7.5 m/s为例对流场分析,空气流经三角形、方形和圆形钝体的流场压力云图和速度云图分别如图3(a)和(b)。压力云图显示高压区基本出现在迎风面区域,低压区则在下游流动时交替分布。当漩涡从钝体上脱落时,靠近管壁的区域会形成一个低压区,随着漩涡继续运动,低压区将向下游转移并被一个高涡所取代,而前一个高压区则被一个较低压力的漩涡所取代,中心压力降低,并沿径向递增。对比三角形、方形和圆形钝体的漩涡分布可知,圆形钝体低压区从钝体两侧开始脱落且尾流区的摆动幅度较小,漩涡主要聚集在钝体两侧;而方形钝体和三角形钝体的低压区从迎流面的两个尖点开始脱落且漩涡逐渐向管壁两侧扩展,摆动幅度较大。从中可以得出结论,两个尖点的存在增强了前缘位置的摩擦阻力和速度梯度,促使漩涡径向扩散。
图3 流经钝体流场压力与速度云图
尾流流场内涡街脱落对压电薄膜作用过程如图4,流场风速18.8 m/s时压电薄膜置于钝体尾流1B位置三维流场截面涡量云图。在t=0.05 s时,漩涡形成初期附着在三角形钝体两侧,未发生分离;当t=0.1 s时,漩涡1到达压电薄膜下表面,压迫薄膜向上变形;而当漩涡1从钝体脱离时漩涡2在钝体另一侧逐渐形成;当t=0.15 s时,漩涡2成型,漩涡1逐渐扩散,使薄膜下表面压力消散,压力消失,形变回复;当t=0.2 s时,漩涡3到达压电薄膜下表面,并在下表面形成压力差推动压电薄膜再次向上运动,当漩涡3向下扩散衰减后,钝体另一侧再次形成新的漩涡,如此往复。每一个漩涡的形成和消散都会经历成形-壮大-衰减过程。漩涡脱落分离过程致使压电薄膜下表面的压力差值不断改变,引起应力应变,利用压电圆片机电耦合效应转换涡激动能为电能。
图4 涡街漩涡与压电薄膜作用过程
Fig.4 Vortex and piezoelectric film interaction
压电晶圆换能结构在涡街尾流区受到波动气压作用,为了确定最优能量转换区域,需要进行涡街脱离流场分析。根据图5(a)可知,P1监测点压力时间变化曲线接近正弦曲线。随着监测点位置改变,流场中各点的压力随之变化,进而导致其他监测点压力波动为非正弦曲线状态。对压力时间曲线进行FFT分析,获得各监测点处漩涡脱落频率如图5(b)。对比五个监测点压力FFT变换曲线,发现P1监测点信号强度最大,且只存在一个主峰,说明发生体表面涡街分布均匀。而P0、P2~P4监测点存在多个峰值,主峰相对较弱,涡街脱离不稳定。根据涡街信号强度和稳定性分析,可以确定最佳压力监测点位置在P1(180,55)附近,将压电薄膜放置此区域范围内,系统能够获得较高能量采集性能。
为了进一步研究钝体结构对尾流区压力影响,在风速为7.5 m/s时,各监测点方形钝体和圆形钝体的压力波动曲线及FFT变换曲线分别如图5(a2)、(a3)和(b2)、(b3)。根据方形钝体压力云图可知,当风流经方形钝体时,P0、P1和P2监测点位置的漩涡比较稳定,因此在压力-时间曲线中表现为标准的正弦曲线,如图5(a2),而P3、P4监测点漩涡由于存在耗散,导致漩涡变得不稳定,压力-时间曲线表现为不规则正弦输出。观察各监测点的FFT变换曲线,发现P0~P4五个监测点压力幅值相差较小,P1监测点压强幅值略高于其它几个监测点的幅值,但由于方形钝体不同监测点的相位有所差异,FFT变换曲线存在不同程度的扰动,如图5(b2)。
五个监测点位置压强变化对比如图5(a)系列,与方形钝体相比,圆形钝体尾迹的压强幅度相差较大,但其幅值较小。主要由于钝体边界比较光滑,表面涡街分离不均衡导致反转速度梯度较小,压差较小。根据圆形钝体不同位置处FFT曲线可知,P0和P1监测点压力信号稳定,但P0监测点的幅值要明显要低于P1监测点,而P2~P4监测点波动不稳定且幅值较小。主要是P0监测点离钝体较近,涡街刚从钝体表面脱离,未得到完全发展,气流波动较小;P2监测点处的涡街漩涡虽得到了充分发展,但该区域压力较小,而对于P3和P4监测点,其距离钝体位置较远,形成的涡街不稳定导致能量发生分散,压强逐渐减弱。
涡街强度表征了由漩涡引起的流场振动特性。设计钝体迎风面宽度30 mm,以不同风速流经三角形钝体,各监测点的压力相差明显,为了确定最佳压力点位置,图6(a)给出了风速与各监测点涡街强度关系,数据取涡街强度最大峰值。随着监测点到钝体距离变化,涡街强度呈先增大后减小趋势,且风速越大,趋势越明显。同时,随着风速的增加,漩涡脱落强度明显增加。比较不同钝体在P1监测点的涡街强度幅值,如图6(b)。三种钝体产生的涡街强度相差较大,三角形钝体最高。由于三角形钝体的后缘边长较短甚至变为尖锐棱角,使钝体两侧分离的漩涡相遇时间变短,漩涡脱落速度加快,在每次漩涡脱落的过程中,都会伴随产生较大的压力脉动,导致柱体振动加强,幅值增大;对于圆形和方形钝体,其后缘边长导致两侧涡街相互作用时间被延长,振幅较小,信号强度弱。综上分析,气流流经三角形钝体获得的压力幅值最大。据此,三角形钝体涡激压电薄膜振动效果最优。
(a) 各监测点风速与涡街强度关系
压电薄膜能量采集性能与涡街频率密切相关,当薄膜固有频率与漩涡脱落频率相匹配时达到共振频率范围,可输出较大能量。漩涡脱落频率主要受风速和钝体尺寸影响,分析不同迎风面宽度和风速下漩涡脱落频率。三角形钝体尾流区获得涡街频率如图7。当迎风面宽度一定时,涡激频率随风速增加呈直线增大;随着迎风面宽度的增大,风速对频率的影响逐渐降低。当风速一定时,频率随迎风面宽度增加逐渐降低;随着风速增加,迎风面宽度对脱落频率影响增大。当迎风面宽度为30 mm,风速从2.6 m/s增加到18.8 m/s时,漩涡脱落频率由18.2 Hz增加至131.6 Hz,频率增加效果明显;风速增加,压电发电输出能量提高。
压电薄膜尺寸为φ35 mm×30 μm的PVDF压电薄膜周边固支,在流场风速18.8 m/s时,尾流漩涡激励作用产生应力分布和压电薄膜变形如图8。当压电薄膜发生弯曲变形时,压电薄膜的最大应力集中在边界位置处,且沿半径方向由内向外逐渐增大。涡激作用下不同压电薄膜厚度形变云图如图8(b)。
图7 漩涡脱落频率与风速和迎风面宽度关系
(a) 应力分布云图
图9(a)和(b)分别分析了压电薄膜厚度与形变量和开路电压等时间响应曲线。随着压电薄膜厚度的增加,薄膜总体变形量减小;当厚度由30 μm增加到100 μm时,最大形变量从0.989×10-3m逐渐减小到0.301×10-5m。压电薄膜的开路电压波形均呈正弦变化,输出电压幅值随形变量增加而增高。在钝体迎流面宽度30 mm,压电薄膜厚度30 μm,风速18.8 m/s时,压电振子输出电压达到8.07 V。图9(c)和(d)分别描述了风速和钝体直径改变时压电薄膜开路电压随时间响应曲线。随着风速的增加,压电薄膜开路电压显著增大;随着钝体直径的增加,压电薄膜的开路电压亦随之增加,而涡激频率随之减小,电压从8.07 V增加到8.97 V,对应的频率从131.6 Hz减小到78.9 Hz。
(a) 压电薄膜径向各点形变量
基于流体力学控制方程及钝体绕流相关理论,采用Fluent流体分析软件中的Realizable k-ε湍流模型仿真分析流致涡激压电风能采集系统。
(1) 流场绕流钝体二维流场分析结果获取涡街在流场中的压力分布特点、涡脱落规律。在钝体尾流区域设置监测点,研究了不同位置处的漩涡强度和稳定性,涡街强度和稳定性随监测点位置变化关系。对比不同钝体形状和风速下漩涡强度,结果表明:三角形钝体涡激强度最大,同时比较尾流区各监测点涡激强度,能够确定压电换能器最佳分布位置。
(2) 建立流致涡激压电风能采集系统的三维仿真模型,研究涡流激励发电性能。确定压电薄膜位移和开路电压随时间变化规律,结果表明,当压电薄膜的厚度为30 μm以及钝体直径为30 mm时,输出电压为正弦信号且最大电压为8.07 V。