查看: 67|回复: 0|关注: 0

[未答复] fsolve求解线性方程组由于初值问题出现一些错误解

[复制链接]

新手

27 麦片

财富积分


050


0

主题

0

帖子

0

最佳答案
发表于 6 天前 | 显示全部楼层 |阅读模式
求解方程:其中方程phi,phi1,A和Aa未知。求解自己编的程序:

clear all %二自由度系统的解析解

clc

Omega=linspace(0,20,600)

allx = [];

zeta1=0.05;

zeta2=0.05;

p=0.05 %p=ma/m (质量比)

n=2;%kv2/kv 刚度比2

f1=0.005;

v=n/(n-(1+p))%kv1/kv 刚度比1

for i = 1:length(Omega)

    x0 = [0,0,0,0];

    x = fsolve(@(x)root2d(x, Omega(i)), x0);

    allx = [allx; x];

end
%
x1=abs((allx(:,3)));

x2=abs((allx(:,3)));

Z1=abs(allx(:,3));

plot(Omega,Z1)

Z2=abs(allx(:,4));

% TFl1=20*log(sqrt(Z1.^2+(2.*zeta1.*Omega.*Z1).^2)./f1);
%
% TFl1=20*log(sqrt(n.^2.*Z2.^2+(2.*zeta2.*Omega.*p).^2.*Z2.^2)./f1)


TFl1=20*log(sqrt(Z1.^2+(2.*zeta1.*Omega.*Z1).^2)./f1)

TFl2=20*log(sqrt(n.^2.*Z2.^2+(2.*p.*zeta2.*Omega).^2.*Z2.^2)./f1);


plot(Omega,TFl2(:,1),'or');


function F = root2d(x, Omega)

zeta1=0.05;

zeta2=0.05;

p=0.05;%p=0.001, p=0.005

n=2;

f1=0.005;

v=n/(n-(1+p))

F(1)=(v-Omega.^2).*x(3).*cos(x(1))-v.*x(4).*cos(x(2))-2.*zeta1.*x(3).*Omega.*sin(x(1))+2.*zeta1.*x(4).*Omega.*sin(x(2))-f1;% x(3)=A,x(4)=Aa,x(2)=phia.x(1)=phi

F(2)=(v-Omega.^2).*x(3).*sin(x(1))-v.*x(4).*sin(x(2))+2.*zeta1.*x(3).*Omega.*cos(x(1))-2.*zeta1.*x(4).*Omega.*cos(x(2));

F(3)=(v+n-p.*Omega.^2).*x(4).*cos(x(2))-x(3).*v.*cos(x(1))+2.*zeta1.*x(3).*Omega.*sin(x(1))-2.*(zeta1+p.*zeta2).*Omega.*x(4).*sin(x(2));

F(4)=(v+n-p.*Omega.^2).*x(4).*sin(x(2))-x(3).*v.*sin(x(1))-2.*zeta1.*x(3).*Omega.*cos(x(1))+2.*(zeta1+p.*zeta2).*Omega.*x(4).*cos(x(2));

end
结果
使用fsolve在第一个峰值出现根混乱的情况。跟实际情况不符。
问题:1.如何方程对初值敏感度比较高,fsolve不适用怎么求解。
           2.如何排除那些错误的解.

拓展:希望有求解过这种线性方程组的老铁能给个建议。本人论文只差这个地方,比较急。如果有懂求解非线性方程组的兄弟可以加个好友,解决这个问题可以挂您一个论文作者。有偿也可以。






您需要登录后才可以回帖 登录 | 注册

本版积分规则

关闭

站长推荐上一条 /3 下一条

快速回复 返回顶部 返回列表