外科手术器械设计 超声刀设计思路和原理 - 04

根据上文牛顿第二定律的谐振频率模型可以得到：

又有，在无阻尼和无外力的情况下，可简化为：

根据简谐振动的情况，里面的t可以用x代替；根据二阶齐次微分方程的通解可以求得振速，和应力的关系：

syms k1 k2 k3 k4 k5 k6; 
syms l1 l2 l3 l4 l5 l6;
syms A1 A2 A3 A4 A5 A6;
syms B1 B2 B3 B4 B5 B6;
syms Z1 Z2 Z3 Z4 Z5 Z6;

syms x1 x2 x3 x4 x5 x6; % variables

syms vb vf; %input velocity and output velocity

v1 = A1*sin(k1*x1)+B1*cos(k1*x1); %velocity equation
v2 = A2*sin(k2*x2)+B2*cos(k2*x2);
v3 = A3*sin(k3*x3)+B3*cos(k3*x3);
v4 = A4*sin(k4*x4)+B4*cos(k4*x4);
v5 = A5*sin(k5*x5)+B5*cos(k5*x5);
v6 = A6*sin(k6*x6)+B6*cos(k6*x6);

f1 = -1i*Z1*(A1*cos(k1*x1)-B1*sin(k1*x1)); %force equation
f2 = -1i*Z2*(A2*cos(k2*x2)-B2*sin(k2*x2));
f3 = -1i*Z3*(A3*cos(k3*x3)-B3*sin(k3*x3));
f4 = -1i*Z4*(A4*cos(k4*x4)-B4*sin(k4*x4));
f5 = -1i*Z5*(A5*cos(k5*x5)-B5*sin(k5*x5));
f6 = -1i*Z6*(A6*cos(k6*x6)-B6*sin(k6*x6));

因为力和速度是连续的，所以可以有下面边界条件：

tempA3 = solve(subs(f3,x3,0),A3);
tempB3 = solve(subs(v3,x3,0)==vb, B3);
tempB2 = solve(subs(v2,x2,0)==subs(subs(subs(v3,x3,l3),A3,tempA3),B3,tempB3),B2);
tempA2 = solve(subs(f2,x2,0)==subs(subs(subs(f3,x3,l3),A3,tempA3),B3,tempB3),A2);
tempB1 = solve(subs(v1,x1,0)==subs(subs(subs(v2,x2,l2),A2,tempA2),B2,tempB2),B1);
tempA1 = solve(subs(subs(subs(v1,x1,l1),B1,tempB1),B2,tempB2)==0,A1);
tempB4 = 0;
tempB5 = vf*cos(k5*l5);
tempA4 = vf*cos(k5*l5)/sin(k4*l4);
tempA5 = vf*sin(k5*l5);

%
coefficient(1,1) = -(cos(k1*l1)*(vb*cos(k2*l2)*cos(k3*l3) - (Z3*vb*sin(k2*l2)*sin(k3*l3))/Z2))/sin(k1*l1);
coefficient(1,2) = -(Z3*vb*sin(k3*l3))/Z2;
coefficient(1,3) = 0;
coefficient(1,4) = (vf*cos(k5*l5))/sin(k4*l4);
coefficient(1,5) = vf*sin(k5*l5);
coefficient(2,1) = vb*cos(k2*l2)*cos(k3*l3) - (Z3*vb*sin(k2*l2)*sin(k3*l3))/Z2;
coefficient(2,2) = vb*cos(k3*l3);
coefficient(2,3) = vb;
coefficient(2,4) = 0;
coefficient(2,5) = vf*cos(k5*l5);
%

再根据边界条件中多余的两个等式可求得节点平面左边和右边的频率方程 

FreqEquationLeft = Z2*(vb*cos(k3*l3)*sin(k2*l2) + (Z3*vb*cos(k2*l2)*sin(k3*l3))/Z2) - (Z1*cos(k1*l1)*(vb*cos(k2*l2)*cos(k3*l3) - (Z3*vb*sin(k2*l2)*sin(k3*l3))/Z2))/sin(k1*l1);

%if l1=0 the equation can simplify:
FreqEquationLeft = (Z3*tan(k2*l2)*tan(k3*l3))/Z2 == 1
%

FreqEquationRight = (Z4*vf*cos(k4*l4)*cos(k5*l5))/sin(k4*l4) - Z5*vf*sin(k5*l5);

 至此我们可以进行简单的计算了：

%input parameters%
tempZ3 = 7900*sqrt(19.5E+10/7900)*1E-6; %SS
tempZ2 = 7750*sqrt(8E+10/7750)*1E-6; %PZT
tempk2 = 2*pi*55E+3/sqrt(8E+10/7750);
templ2 = 0.0025*4;
tempk3 = 2*pi*55E+3/sqrt(19.5E+10/7900);

%=>
templ3 = atan(tempZ2/(tempZ3*tan(tempk2*templ2)))/tempk3
%0.047%
%

%input parameters%
ratio = vf/vb == -(tempZ2/tempZ5)*(cos(k4l4)/sin(k2l2))*sqrt((1+(tempZ5/tempZ4)^2*tan(k4*l4)^2)/(1+(tempZ2/tempZ3)^2*atan(k2*l2)^2));
templ4 = 0.015;
tempZ4 = 2790*(sqrt(7.15E+10/2790))*1E-6;
tempZ5 = 2790*(sqrt(7.15E+10/2790))*1E-6
R4/R5 = 2.4; %the radius reduction

%=>
templ5 = atan(tempZ4/(tempZ5*tan(tempk4*templ4)))/tempk5
%0.009%
vf/vb
%10.7%
%