1.材料定义

亚波长条件下周期性布置压电片，同时压电片外接电感作为分流电路

2.计算公式

Matlab实现

function [F,TL] = PZTMetamaterialPlate(rho,E,h,v,phi,lbx,lby,lpx,lpy,rhop,hp)
F=zeros(10000,1);
TL=zeros(10000,1);
%% 空气参数定义rho0 = 1.2;         %空气密度c0 = 340;           %空气声速
%% 平板参数定义D = E*h^3/(12*(1-v^2)); %板的弯曲刚度
%% 压电材料参数定义Ap = (lpx*lpy);L = 0.642;s11 = 16.5*10^(-12);s12 = -4.78*10^(-12);d31 = -2.74*10^(-10);e33 = 3.01*10^(-8);
%% 隔声系数计算
for f = 1:1:100000w = 2*pi*f;k0 = w/c0;s = 1i*w;Cp = Ap*e33/hp;Z = 1i*w*L;                                                                         %分流电路阻抗Ep = hp*(1+s*Z*Cp)/(hp*s11*(1+s*Z*Cp)-s*Z*d31^2*Ap);                                %压电片等效杨氏模量vp = -(s12*(1+s*Z*Cp)-s*Z*d31^2*Ap*hp^(-1))/(s11*(1+s*Z*Cp)-s*Z*d31^2*Ap*hp^(-1));  %压电片等效泊松比Daj = D+2*Ep*((h+2*hp)^3-h^3)/((24*(1-vp^2)));                                      %A区域内等效弯曲刚度Dbj = D;                                                                            %B区域内的等效弯曲刚度maj = rho*h+2*rhop*hp;                                                              %A区域内等效面质量mbj = rho*h;                                                                        %B区域内等效面质量alpha = (lpx*lpy)/(lbx*lby);                                                        %AB区域之间面积比meq = alpha*maj+(1-alpha)*mbj;                                                      %单胞元等效面质量Deq = Daj*Dbj/((1-alpha)*Daj+alpha*Dbj);                                            %单胞元等效弯曲刚度delta = Deq*(k0*sin(phi))^4-meq*w^2+2*1i*rho0*c0*w/cos(phi);t = (2*rho0*c0*w/(delta*cos(phi)))^2;F(f) = f;    TL(f) = 10*log10(1/t);
end
semilogx(F,TL,'k','linewidth',2);
end

3.参数定义

PZTMetamaterialPlate(2730,77.6e9,0.005,0.3519,pi/6,0.040,0.040,0.036,0.036,7500,0.001);

4.结果对比

5.参考文献

[1] Hao Zhang, Jihong Wen, Yong Xiao, Gang Wang, Xisen Wen,Sound transmission loss of metamaterial thin plates with periodic subwavelength arrays of shunted piezoelectric patches,Journal of Sound and Vibration, Volume 343,2015, Pages 104-120, ISSN 0022-460X

6.固支有限大情况拓展研究

不同于上述论文中的无限大平板采用的结构弯曲波法计算隔声量，固支条件下需要采用模态叠加法；本章在光板模态叠加法的基础上，增加了压电谐振等效介质方法来计算超材料板的等效动态弯曲刚度和等效质量密度。

等效介质法参考本文第二章；

模态叠加法参考：

固支边界有限大平板隔声理论及Matlab实现

在模态叠加法的计算过程中，将等效介质法得到的等效动态弯曲刚度 $D_{eq}$ 和等效质量密度 $\rho_{eq}$ 替换到动力学方程中进行计算即可；

$D_{eq}\nabla^4w+\rho_{eq}h\frac{\partial^2w}{\partial t^2}=jw\rho_0(\Phi_{in}-\Phi_{t})$

function RigidPlate()j = sqrt(-1);
%% 预分配内存F=zeros(10000,1);TL = zeros(10000,1);bar= waitbar(0,'Simulation in process');%显示运算进度
%% 空气参数定义rho0 = 1.21;%空气密度c0 = 343;%空气声速phi = 0;%入射角theta = 0;%方位角Pin = 1;%入射波幅值
%% 板参数定义rho = 2700;%板密度E = 77.6e9;%板杨氏模量v = 0.3519;%板泊松比h = 0.005;%板厚度Lx = 0.04;%板x方向长度Ly = 0.04;%板y方向长度
%% 压电材料参数定义lpx = 0.036;%压电片x方向长度lpy = 0.036;%压电片y方向长度hp = 0.001;%压电片厚度rhop = 7500;%压电片密度Ap = (lpx*lpy);%压电片面积L = 0.642;%电感x = (lpx*lpy)/(Lx*Ly);eta = 0.001;%压电系数s11 = 16.5*10^(-12);s12 = -4.78*10^(-12);d31 = -2.74*10^(-10);e33 = 3.01*10^(-8);
%% 平板参数定义D = E*(1+j*0.005)*h^3/(12*(1-v^2));                              %板的弯曲刚度%模态截断M = 20;                                                          N = 20;
%% 求解delta算子[delta11,delta12,delta13,delta21,delta22,delta23,delta24] = GetDelta(Lx,Ly,M,N);
%% 隔声系数计算fmin = 1;%计算下限fmax = 10000;%计算上限df = 1;for f = fmin:df:fmaxstr=['Simulation in process:' num2str(ceil((f-fmin)*100/(fmax-fmin))) '%'];waitbar((f-fmin)/(fmax-fmin),bar,str);%进度显示w = 2*pi*f;k0 = w/c0;kx = k0*sin(phi)*cos(theta);ky = k0*sin(phi)*sin(theta);kz = k0*cos(phi);% 等效动态弯曲刚度求解s = 1i*w;Cp = Ap*e33/hp;Z = 1i*w*L;%电感型电路Ep = hp*(1+s*Z*Cp)/(hp*s11*(1+s*Z*Cp)-s*Z*d31^2*Ap);vp = -(s12*(1+s*Z*Cp)-s*Z*d31^2*Ap*hp^(-1))/(s11*(1+s*Z*Cp)-s*Z*d31^2*Ap*hp^(-1));Daj = D+2*Ep*((h+2*hp)^3-h^3)/((24*(1-vp^2)));%单片压电是不乘以2，双片压电乘以2Dbj = D;Deq = Daj*Dbj/((1-x)*Daj+x*Dbj);rhoeq = (x*(Lx*Ly*h*rho)+(1-x)*(lpx*lpy*hp*rhop))/(Lx*Ly*h);%求解Mrs、FrsMrs = (4*Deq*(pi^4)*Lx*Ly).*(delta11+delta12+delta13)-(rhoeq*h*(w^2)+j*w*rho0*2*w/kz).*(delta21+delta22+delta23+delta24);[Frs,I] = GetF(M,N,Pin,kx,ky,Lx,Ly,w,rho0,j);abare = Mrs\Frs;F(f) = f;TL(f) = GetTL(abare,rho0,w,c0,Pin,kx,ky,kz,Lx,Ly,M,N,j);endsemilogx(F(fmin:df:fmax),TL(fmin:df:fmax),'k','linewidth',2);
end