fsolve求解12阶行列式方程(超越方程),遇到错误。
이전 댓글 표시
E2=76e9;
rho2=10500;
%%%%
E1=160e9;
rho1=2300;
%%%%
h1=2*1e-6;
h2=5*1e-6;
h0=h2-h1;
a=150*1e-6;
b=350*1e-6;
l=500*1e-6;
z=10*1e-6;
%%%%
co=1; %缩小系数
syms omiga;
beta1=(10.15179*(co*omiga)^0.5) %beta1、2、3分别为以下行列式中的参数
beta2=(10.1533*(co*omiga)^0.5)
beta3=(10.15179*(co*omiga)^0.5)
for i=1:2000 %求解超越方程特征根,以下12阶行列式中omiga为待求解的未知数
root(i)=fsolve(@(omiga) det([0,1,0,1,0,0,0,0,0,0,0,0;
beta1,0,beta1,0,0,0,0,0,0,0,0,0;
0,0,0,0,0,0,0,0,sin(beta3*l),cos(beta3*l),sinh(beta3*l),cosh(beta3*l);
0,0,0,0,0,0,0,0,beta3*cos(beta3*l),-beta3*sin(beta3*l),beta3*cosh(beta3*l),beta3*sinh(beta3*l);
sin(beta1*a),cos(beta1*a),sinh(beta1*a),cosh(beta1*a),-sin(beta2*a),-cos(beta2*a),-sinh(beta2*a),-cosh(beta2*a),0,0,0,0;
0,0,0,0,sin(beta2*b),cos(beta2*b),sinh(beta2*b),cosh(beta2*b),-sin(beta3*b),-cos(beta3*b),-sinh(beta3*b),-cosh(beta3*b);
beta1*cos(beta1*a),-beta1*sin(beta1*a),beta1*cosh(beta1*a),beta1*sinh(beta1*a),-beta2*cos(beta2*a),beta2*sin(beta2*a),-beta2*cosh(beta2*a),-beta2*sinh(beta2*a),0,0,0,0;
0,0,0,0,beta2*cos(beta2*b),-beta2*sin(beta2*b),beta2*cosh(beta2*b),beta2*sinh(beta2*b),-beta3*cos(beta3*b),beta3*sin(beta3*b),-beta3*cosh(beta3*b),-beta3*sinh(beta3*b);
-((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta1)^2)*sin(beta1*a),-((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta1)^2)*cos(beta1*a),((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta1)^2)*sinh(beta1*a),((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta1)^2)*cosh(beta1*a),(2*E2*z*((h2^3-h1^3)/3))*((beta2)^2)*sin(beta2*a),(2*E2*z*((h2^3-h1^3)/3))*((beta2)^2)*cos(beta2*a),-(2*E2*z*((h2^3-h1^3)/3))*((beta2)^2)*sinh(beta2*a),-(2*E2*z*((h2^3-h1^3)/3))*((beta2)^2)*cosh(beta2*a),0,0,0,0;
0,0,0,0,-(2*E2*z*((h2^3-h1^3)/3))*((beta2)^2)*sin(beta2*b),-(2*E2*z*((h2^3-h1^3)/3))*((beta2)^2)*cos(beta2*b),(2*E2*z*((h2^3-h1^3)/3))*((beta2)^2)*sinh(beta2*b),(2*E2*z*((h2^3-h1^3)/3))*((beta2)^2)*cosh(beta2*b),((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta3)^2)*sin(beta3*b),((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta3)^2)*cos(beta3*b),-((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta3)^2)*sinh(beta3*b),-((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta3)^2)*cosh(beta3*b);
-((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta1)^3)*cos(beta1*a),((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta1)^3)*sin(beta1*a),((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta1)^3)*cosh(beta1*a),((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta1)^3)*sinh(beta1*a),(2*E2*z*((h2^3-h1^3)/3))*((beta2)^3)*cos(beta2*a),-(2*E2*z*((h2^3-h1^3)/3))*((beta2)^3)*sin(beta2*a),-(2*E2*z*((h2^3-h1^3)/3))*((beta2)^3)*cosh(beta2*a),-(2*E2*z*((h2^3-h1^3)/3))*((beta2)^3)*sinh(beta2*a),0,0,0,0;
0,0,0,0,-(2*E2*z*((h2^3-h1^3)/3))*((beta2)^3)*cos(beta2*b),(2*E2*z*((h2^3-h1^3)/3))*((beta2)^3)*sin(beta2*b),(2*E2*z*((h2^3-h1^3)/3))*((beta2)^3)*cosh(beta2*b),(2*E2*z*((h2^3-h1^3)/3))*((beta2)^3)*sinh(beta2*b),((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta3)^3)*cos(beta3*b),-((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta3)^3)*sin(beta3*b),-((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta3)^3)*cosh(beta3*b),-((2*z/3)*(E1*(h1)^3+E2*((h2)^3-(h1)^3)))*((beta3)^3)*sinh(beta3*b)]),i,optimset('display','off'));
end
ro1=root;
j_tem=1;
k_tem=1;
for i_tem=1:1999 %去重根
b_tem(i_tem)=ro1(i_tem+1);
end
b_tem(2000)=ro1(2000);
c_tem=b_tem-ro1;
c_tem_find_first= find( c_tem>0.1);
r_tem(j_tem)=ro1(c_tem_find_first(1));
r_tem(j_tem+1)=b_tem(c_tem_find_first(1));
j_tem=j_tem+2;
for i_c=c_tem_find_first(2):2000
if c_tem(i_c)>0.1
r_tem(j_tem)=b_tem(i_c);
j_tem=j_tem+1;
end
end
r_tem;
for i_re=1:j_tem-2
root_return(i_re)=r_tem(i_re+1);
end
root_return;
ri=root_return/co


채택된 답변
추가 답변 (0개)
카테고리
도움말 센터 및 File Exchange에서 Linear Algebra에 대해 자세히 알아보기
Community Treasure Hunt
Find the treasures in MATLAB Central and discover how the community can help you!
Start Hunting!