压电式加速度传感器动态特性
压电式加速度传感器工作原理如下图所示：
将加速度传感器壳体位移设为x1(t)x_1(t)x1(t)，质块m运动设为x0(t)x_0(t)x0(t)，
质块相对于壳体的运动设为x01(t)x_{01}(t)x01(t)。
压电式加速度传感器的数学模型
依照传感器的力学模型，有达朗贝尔原理可以有以下运动微分方程：
md2x0dt2+c(dx0dt−dx1dt)+k(x0(t)−x1(t))m\frac{d^2x_0}{dt^2}+c(\frac{dx_0}{dt}-\frac{dx_1}{dt})+k(x_0(t)-x_1(t)) mdt2d2x0+c(dtdx0−dtdx1)+k(x0(t)−x1(t))
压电式传感器系统输入为壳体运动x1(t)x_1(t)x1(t)，输出为质块m与壳体的相对运动x01(t)x_{01}(t)x01(t)，此二者的关系为：
x01(t)=x0(t)−x1(t)x_{01}(t)=x_0(t)-x_1(t) x01(t)=x0(t)−x1(t)
于是运动微分方程可以写为：
md2x01dt2+c(dx01dt)+kx01(t)=−md2x1dt2m\frac{d^2x_01}{dt^2}+c(\frac{dx_{01}}{dt})+kx_{01}(t)=-m\frac{d^2x_1}{dt^2} mdt2d2x01+c(dtdx01)+kx01(t)=−mdt2d2x1
设被测振动为为谐波振动，即有：
x1(t)=X1sin(ωt)x_1(t)=X_1sin({\omega}t) x1(t)=X1sin(ωt)
于是运动微分方程：
md2x01dt2+c(dx01dt)+kx01(t)=mX1ω2sin(ωt)m\frac{d^2x_01}{dt^2}+c(\frac{dx_{01}}{dt})+kx_{01}(t)=mX_1{\omega}^2sin({\omega}t) mdt2d2x01+c(dtdx01)+kx01(t)=mX1ω2sin(ωt)
将其写为一般形式为：
d2x01dt2+2ξωn(dx01dt)+ω2x01(t)=X1ω2sin(ωt)\frac{d^2x_01}{dt^2}+2{\xi}{\omega_n}(\frac{dx_{01}}{dt})+{\omega^2}x_{01}(t)=X_1{\omega}^2sin({\omega}t) dt2d2x01+2ξωn(dtdx01)+ω2x01(t)=X1ω2sin(ωt)
次微分方程的解由通解和特解两部分组成，其中通解反映此系统的固有特性，特解反映其动态特性；
由于激励为一个正弦信号，因此其响应也为同频率的正弦信号，设其响应为：
x01(t)=X01sin(ωt+ϕ)x_{01}(t)=X_{01}\sin({\omega}t+\phi) x01(t)=X01sin(ωt+ϕ)
于是有：
dx01dt=ωX01cos(ωt+ϕ)=ωX01sin(ωt+ϕ+π2)dx012dt2=−ω2X01sin(ωt+ϕ)=dx012dt2=ω2X01sin(ωt+ϕ+π)\frac{dx_{01}}{dt}={\omega}X_{01}\cos{({\omega}t+\phi)}={\omega}X_{01}\sin{({\omega}t+\phi+\frac{\pi}{2})}\\ \frac{dx^2_{01}}{dt^2}=-{\omega}^2X_{01}\sin{({\omega}t+\phi)}=\frac{dx^2_{01}}{dt^2}={\omega}^2X_{01}\sin{({\omega}t+\phi+\pi)} dtdx01=ωX01cos(ωt+ϕ)=ωX01sin(ωt+ϕ+2π)dt2dx012=−ω2X01sin(ωt+ϕ)=dt2dx012=ω2X01sin(ωt+ϕ+π)
按照旋转向量计算，最终得出其响应为：
Ax(ω)=X01X1=(ω/ωn)[1−(ω/ωn)2]2+[2ξ(ω/ωn)]2ϕ(ω)=arctan(2ξ(ω/ωn)1−(ω/ωn)2)A_x(\omega)=\frac{X_01}{X_1}=\frac{({\omega}/{\omega}_n)}{\sqrt{[1-({\omega}/{\omega}_n)^2]^2+[2{\xi}({\omega}/{\omega}_n)]^2}}\\ {\phi}(\omega)=\arctan(\frac{2{\xi}({\omega}/{\omega}_n)}{1-({\omega}/{\omega}_n)^2}) Ax(ω)=X1X01=[1−(ω/ωn)2]2+[2ξ(ω/ωn)]2(ω/ωn)ϕ(ω)=arctan(1−(ω/ωn)22ξ(ω/ωn))
MATALB绘制幅频特性曲线与相频特性曲线
由于MATLAB中反三角函数arctanarctanarctan输出默认为[−π/2,π/2][-\pi/2,\pi/2][−π/2,π/2]，因此需要将自定义函数，将其输出改为[−π/2,π/2][-\pi/2,\pi/2][−π/2,π/2]，代码为：
function angle = atandpi(x)%将atand输出更改为0-pi% 此处显示详细说明y=atand(x);L=length(y);for i=1:Lif y(i)<0y(i)=y(i)+180;endendangle=y;end
然后绘制幅频与相频特性曲线，代码为：
%压电式振动传感器的动态特性clear;clc;r=0:0.01:10;xi=[0,0.3,0.4,0.5,0.6,0.7,1,2,5,10];%位移幅频特性Ax1=(r.^2)./sqrt(((1-(r.^2)).^2)+((2*(xi(1))*(r)).^2));Ax2=(r.^2)./sqrt(((1-(r.^2)).^2)+((2*(xi(2))*(r)).^2));Ax3=(r.^2)./sqrt(((1-(r.^2)).^2)+((2*(xi(3))*(r)).^2));Ax4=(r.^2)./sqrt(((1-(r.^2)).^2)+((2*(xi(4))*(r)).^2));Ax5=(r.^2)./sqrt(((1-(r.^2)).^2)+((2*(xi(5))*(r)).^2));Ax6=(r.^2)./sqrt(((1-(r.^2)).^2)+((2*(xi(6))*(r)).^2));Ax7=(r.^2)./sqrt(((1-(r.^2)).^2)+((2*(xi(7))*(r)).^2));Ax8=(r.^2)./sqrt(((1-(r.^2)).^2)+((2*(xi(8))*(r)).^2));Ax9=(r.^2)./sqrt(((1-(r.^2)).^2)+((2*(xi(9))*(r)).^2));Ax10=(r.^2)./sqrt(((1-(r.^2)).^2)+((2*(xi(10))*(r)).^2));subplot(2,1,1);plot(r,Ax1,r,Ax2,r,Ax3,r,Ax4,r,Ax5,r,Ax6,r,Ax7,r,Ax8,r,Ax9,r,Ax10);title('位移幅频特性曲线');xlim([0 5]);ylim([0.1 2.5]);grid on;%相频特性曲线phi1=atandpi((2*xi(1)*r)./(1-(r.^2)));phi2=atandpi((2*xi(2)*r)./(1-(r.^2)));phi3=atandpi((2*xi(3)*r)./(1-(r.^2)));phi4=atandpi((2*xi(4)*r)./(1-(r.^2)));phi5=atandpi((2*xi(5)*r)./(1-(r.^2)));phi6=atandpi((2*xi(6)*r)./(1-(r.^2)));phi7=atandpi((2*xi(7)*r)./(1-(r.^2)));phi8=atandpi((2*xi(8)*r)./(1-(r.^2)));phi9=atandpi((2*xi(9)*r)./(1-(r.^2)));phi10=atandpi((2*xi(10)*r)./(1-(r.^2)));subplot(2,1,2);plot(r,phi1,r,phi2,r,phi3,r,phi4,r,phi5,r,phi6,r,phi7,r,phi8,r,phi9,r,phi10);title('相频特性曲线');xlim([0 3])grid on;%速度幅频特性
运行结果如下：
