xl2033 2022-08-16 11:23 采纳率: 57.1%
浏览 545
已结题

matlab fsolve求解非线性方程组为什么得到复数解

用matlab fsolve求解非线性方程组得到的是复数解,该方程组有实数解,而且发现求得的复数解的实部和答案基本吻合,这是为什么呢?下面是我建立的方程和fsolve求解的程序

function F = zw(X)
da= X(1);
dr= X(2);
ct= X(3);
dc=31;
Dw = 12.303;
A = 0.627; 
dm = 58;
a=35/180*pi;
Kn=3.0363*10^5;
ri=6.405;
Ri=dm/2+(ri-Dw/2)*cos(a);

    Fr=4500;Fa=3000;M=800000;
    p1 = 0:2*pi/13:24*pi/13;
    p2 = pi/13:2*pi/13:25*pi/13;
    F(1) = Fa-Kn*A^1.5*sum((sin(a) + da/A+Ri*(ct/A)*cos(p1)).*(((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*ct*cos(p1)/A).^2).^0.5 - 1).^1.5./((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*(ct/A)*cos(p1)).^2).^0.5)+Kn*A^1.5*sum((sin(a) - da/A-Ri*(ct/A)*cos(p2)).*(((sin(a) - da/A-Ri*(ct/A)*cos(p2)).^2 + (cos(a) + dr* cos(p2)/A-0.5*dc*ct*cos(p2)/A).^2).^0.5 - 1).^1.5./((sin(a) - da/A-Ri*(ct/A)*cos(p2)).^2 + (cos(a) + dr* cos(p2)/A-0.5*dc*(ct/A)*cos(p2)).^2).^0.5);
    F(2) = Fr-Kn*A^1.5*sum((cos(a) + (dr/A)*cos(p1)+0.5*dc*(ct/A)*cos(p1)).*cos(p1).*(((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*ct*cos(p1)/A).^2).^0.5 - 1).^1.5./((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*(ct/A)*cos(p1)).^2).^0.5)-Kn*A^1.5*sum((cos(a) + (dr/A)*cos(p2)-0.5*dc*(ct/A)*cos(p2)).*cos(p2).*(((sin(a) -da/A-Ri*(ct/A)*cos(p2)).^2 + (cos(a) + dr* cos(p2)/A-0.5*dc*ct*cos(p2)/A).^2).^0.5 - 1).^1.5./((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*(ct/A)*cos(p1)).^2).^0.5);
    F(3)=M-dm/2*Kn*A^1.5*sum((sin(a) + da/A+Ri*(ct/A)*cos(p1)).*cos(p1).*(((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*ct*cos(p1)/A).^2).^0.5 - 1).^1.5./((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*(ct/A)*cos(p1)).^2).^0.5)-dc/2*Kn*A^1.5*sum((cos(a) + (dr/A)*cos(p1)+0.5*dc*(ct/A)*cos(p1)).*cos(p1).*(((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*ct*cos(p1)/A).^2).^0.5 - 1).^1.5./((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*(ct/A)*cos(p1)).^2).^0.5)+dm/2*Kn*A^1.5*sum((sin(a) - da/A-Ri*(ct/A)*cos(p2)).*cos(p2).*(((sin(a) - da/A-Ri*(ct/A)*cos(p2)).^2 + (cos(a) + dr* cos(p2)/A-0.5*dc*ct*cos(p2)/A).^2).^0.5 - 1).^1.5./((sin(a) - da/A-Ri*(ct/A)*cos(p2)).^2 + (cos(a) + dr* cos(p2)/A-0.5*dc*(ct/A)*cos(p2)).^2).^0.5)+dc/2*Kn*A^1.5*sum((cos(a) + (dr/A)*cos(p2)-0.5*dc*(ct/A)*cos(p2)).*cos(p2).*(((sin(a) -da/A-Ri*(ct/A)*cos(p2)).^2 + (cos(a) + dr* cos(p2)/A-0.5*dc*ct*cos(p2)/A).^2).^0.5 - 1).^1.5./((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*(ct/A)*cos(p1)).^2).^0.5);
  
end

h=optimset;
h.Display='off';
h.TolX=1e-6;
[jg,fval,exitflag]=fsolve('zw',[0,0,0],h) %求解
  • 写回答

2条回答 默认 最新

  • 斗迷飞鸟 2022-08-19 02:49
    关注

    (1)问题分析
    出现复数解应该是因为你的“zw”函数中有1.5次方和0.5次方,而fsolve搜索解的时候对应幂的底数可能是负的,这样对应的1.5次方和0.5次方就会出现复数,进而导致目标函数值为复数,最终导致得到的解为复数。
    (2)问题的可能改进代码
    考虑在zw中将最后的函数值取模,这样可以得到实数解

    function F = zw_real(X)
    da= X(1);
    dr= X(2);
    ct= X(3);
    dc=31;
    Dw = 12.303;
    A = 0.627;
    dm = 58;
    a=35/180*pi;
    Kn=3.0363*10^5;
    ri=6.405;
    Ri=dm/2+(ri-Dw/2)*cos(a);
    
    Fr=4500;Fa=3000;M=800000;
    p1 = 0:2*pi/13:24*pi/13;
    p2 = pi/13:2*pi/13:25*pi/13;
    F(1) = Fa-Kn*A^1.5*sum((sin(a) + da/A+Ri*(ct/A)*cos(p1)).*(((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*ct*cos(p1)/A).^2).^0.5 - 1).^1.5./((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*(ct/A)*cos(p1)).^2).^0.5)+Kn*A^1.5*sum((sin(a) - da/A-Ri*(ct/A)*cos(p2)).*(((sin(a) - da/A-Ri*(ct/A)*cos(p2)).^2 + (cos(a) + dr* cos(p2)/A-0.5*dc*ct*cos(p2)/A).^2).^0.5 - 1).^1.5./((sin(a) - da/A-Ri*(ct/A)*cos(p2)).^2 + (cos(a) + dr* cos(p2)/A-0.5*dc*(ct/A)*cos(p2)).^2).^0.5);
    F(2) = Fr-Kn*A^1.5*sum((cos(a) + (dr/A)*cos(p1)+0.5*dc*(ct/A)*cos(p1)).*cos(p1).*(((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*ct*cos(p1)/A).^2).^0.5 - 1).^1.5./((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*(ct/A)*cos(p1)).^2).^0.5)-Kn*A^1.5*sum((cos(a) + (dr/A)*cos(p2)-0.5*dc*(ct/A)*cos(p2)).*cos(p2).*(((sin(a) -da/A-Ri*(ct/A)*cos(p2)).^2 + (cos(a) + dr* cos(p2)/A-0.5*dc*ct*cos(p2)/A).^2).^0.5 - 1).^1.5./((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*(ct/A)*cos(p1)).^2).^0.5);
    F(3)=M-dm/2*Kn*A^1.5*sum((sin(a) + da/A+Ri*(ct/A)*cos(p1)).*cos(p1).*(((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*ct*cos(p1)/A).^2).^0.5 - 1).^1.5./((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*(ct/A)*cos(p1)).^2).^0.5)-dc/2*Kn*A^1.5*sum((cos(a) + (dr/A)*cos(p1)+0.5*dc*(ct/A)*cos(p1)).*cos(p1).*(((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*ct*cos(p1)/A).^2).^0.5 - 1).^1.5./((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*(ct/A)*cos(p1)).^2).^0.5)+dm/2*Kn*A^1.5*sum((sin(a) - da/A-Ri*(ct/A)*cos(p2)).*cos(p2).*(((sin(a) - da/A-Ri*(ct/A)*cos(p2)).^2 + (cos(a) + dr* cos(p2)/A-0.5*dc*ct*cos(p2)/A).^2).^0.5 - 1).^1.5./((sin(a) - da/A-Ri*(ct/A)*cos(p2)).^2 + (cos(a) + dr* cos(p2)/A-0.5*dc*(ct/A)*cos(p2)).^2).^0.5)+dc/2*Kn*A^1.5*sum((cos(a) + (dr/A)*cos(p2)-0.5*dc*(ct/A)*cos(p2)).*cos(p2).*(((sin(a) -da/A-Ri*(ct/A)*cos(p2)).^2 + (cos(a) + dr* cos(p2)/A-0.5*dc*ct*cos(p2)/A).^2).^0.5 - 1).^1.5./((sin(a) + da/A+Ri*(ct/A)*cos(p1)).^2 + (cos(a) + dr* cos(p1)/A+0.5*dc*(ct/A)*cos(p1)).^2).^0.5);
    F=abs(F);  %给函数返回值取模
    end
    
    clear;clc;
    h=optimset;
    h.Display='off';
    h.TolX=1e-6;
    format long;
    [jg,fval,exitflag]=fsolve('zw_real',[0,0,0],h) %求解
    format short;
    

    (3)代码运行结果截图
    这是将zw的返回值取模得到的结果,不知道能不能对上你的理论解,欢迎讨论

    img

    本回答被题主选为最佳回答 , 对您是否有帮助呢?
    评论 编辑记录
查看更多回答(1条)

报告相同问题?

问题事件

  • 系统已结题 8月27日
  • 已采纳回答 8月19日
  • 创建了问题 8月16日

悬赏问题

  • ¥15 phython如何实现以下功能?查找同一用户名的消费金额合并—
  • ¥15 孟德尔随机化怎样画共定位分析图
  • ¥18 模拟电路问题解答有偿速度
  • ¥15 CST仿真别人的模型结果仿真结果S参数完全不对
  • ¥15 误删注册表文件致win10无法开启
  • ¥15 请问在阿里云服务器中怎么利用数据库制作网站
  • ¥60 ESP32怎么烧录自启动程序
  • ¥50 html2canvas超出滚动条不显示
  • ¥15 java业务性能问题求解(sql,业务设计相关)
  • ¥15 52810 尾椎c三个a 写蓝牙地址