Matlab 解非线性方程组2

非线性方程组—Matlab—fsolve—机械 CAD 论坛—huright
实例六：
对于求解非线性方程组一般用 fsolve 命令就可以了，但是对于方程组中某一系 数是变化的，该怎么求呢？ %定义方程组如下，其中 k 为变量 function F = myfun(x,k) H=0.32; Pc0=0.23;W=0.18; F=[Pc0+H*(1+1.5*(x(1)/W-1)-0.5*(x(1)/W-1)^3)-x(2); x(1)-k*sqrt(x(2))]; %求解过程 H=0.32; Pc0=0.23;W=0.18; x0 = [2*W; Pc0+2*H]; % 取初值 options = optimset('Display','off'); k=0:0.01:1; % 变量取值范围[0 1] for i=1:1:length(k) kk=k(i); x = fsolve(@(x) myfun(x,kk), x0, options);%求解非线性方程组 x1(i)=x(1); x2(i)=x(2); end plot(k,x1,'-b',k,x2,'-r'); xlabel('k') legend('x1','x2') （来自振动论坛，作者：studyboy）
（4）检验
fxy1=sin(xy(1)-xy(2));fxy2=cos(x y(1)+xy(2));disp([fxy1,fxy2]) 1.0e-006 *
---------------------------------------------------------------------
非线性方程组—Matlab—fsolve—机械 CAD 论坛—huolve—机械 CAD 论坛—huright
实例三：
求下列非线性方程组在(0.5,0.5) 附近的数值解。 (1) 建立函数文件 myfun.m。 function q=myfun(p) x=p(1); y=p(2); q(1)=x-0.6*sin(x)-0.3*cos(y); q(2)=y-0.6*cos(x)+0.3*sin(y); (2) 在给定的初值 x0=0.5,y0=0.5 下，调用 fsolve 函数求方程的根。 x=fsolve('myfun',[0.5,0.5]',optimset('Display','off')) x= 0.6354 0.3734
实例五： ①建立文件 myfxeq05.m function f=myfxeq05(xvect) x = xvect(1); y = xvect(2); z = xvect(3); f(1)=sin(x) + y^2 + log(z) - 7 f(2)=3*x + 2^y - z^3 + 1 f(3)=x + y + z – 5 ②解方程 clear all disp('Solve the following set of nonlinear algebraic equations:') disp(' sin(x) + y^2 + ln(z) - 7 = 0 ') disp(' 3*x + 2^y - z^3 + 1 = 0 ') disp(' x + y + z - 5 = 0 ') disp(' ') xguess=[1 1 1]'; xvect = fsolve('myfxeq05', xguess); x = xvect(1); y = xvect(2); z = xvect(3); disp('The roots from the default "fsolve" are: ') disp([' x = ', num2str(x) ]) disp([' y = ', num2str(y) ]) disp([' z = ', num2str(z) ]) % Repeat with a different set of options ------------------------------options(2) = 1.e-6; %Tolerance for x options(3) = 1.e-6; %Tolerance for f options(5) = 1; %Levenberg-Marquardt Method xvect = fsolve('myfxeq05', xguess, options); x = xvect(1); y = xvect(2); z = xvect(3); disp('The roots from "fsolve" with Levenberg-Marquardt are: ') disp([' x = ', num2str(x) ]) disp([' y = ', num2str(y) ]) disp([' z = ', num2str(z) ])
-1.5
-1
-0.5
0
0.5
1
1.5
2
图 5.6.3-1
两个二元函数 0 等位线的交点图
（2）从图形获取零点的初始近似值 在图 5.6.3-1 中，用 ginput 获取两个函 数 0 等位线（即三线组中间那条线）交点的 坐标。
[x0,y0]=ginput(2); %在图上取两个点的坐标 disp([x0,y0]) -0.7926 0.7926 -0.7843 0.7843
非线性方程组—Matlab—fsolve—机械 CAD 论坛—huright
实例四：
【例 5.6.3-1 】 求 解 二 元 函 数 方 程 组
2 1.5 1 0.5 0 -0.5 -1 -1.5 -2 -2
f 1 ( x, y ) sin( x y ) 0 的零点。 f 2 ( x, y ) cos( x y ) 0
（0）从三维坐标初步观察两函数图形相交 情况
x=-2:0.05:2;y=x;[X,Y]=meshgrid(x ,y); %产生 x-y 平面上网点坐标 F1=sin(X-Y);F2=cos(X+Y); F0=zeros(size(X)); surf(X,Y,F1), xlabel('x'),ylabel('y'), view([-31,62]),hold on, surf(X,Y,F2),surf(X,Y,F0), shading interp, hold off
图 5.6.3-0 两函数的三维相交图
0.7854
0.7854
（1）在某区域观察两函数 0 等位线的交点 情况
clear; x=-2:0.5:2;y=x;[X,Y]=meshgrid(x, y); %产生 x-y 平面上网点坐标 F1=sin(X-Y);F2=cos(X+Y); v=[-0.2, 0, 0.2]; %指定三个等位值， 是 为了更可靠地判断 0 等位线的存在。 contour(X,Y,F1,v) %画 F1 的三条等位线。 hold on,contour(X,Y,F2,v),hold off %画 F2 的三条等位线。
非线性方程组—Matlab—fsolve—机械 CAD 论坛—huright
非线性方程组求解-Matlab-fsolve
实例一： ①建立文件 fun.m： function y=fun(x) y=[x(1)-0.5*sin(x(1))-0.3*cos(x(2)), ... x(2) - 0.5*cos(x(1))+0.3*sin(x(2))]; ②>>clear;x0=[0.1,0.1];fsolve(@fun,x0,optimset('fsolve')) 注：...为续行符 m 文件必须以 function 为文件头， 调用符为@； 文件名必须与定义的函数名相同； fsolve（）主要求解复杂非线性方程和方程组，求解过程是一个逼近过程。
实例二： ①建立文件 fun.m
function F=myfun(x) F=[x(1)-3*x(2)-sin(x(1));2*x(1)+x(2)-cos(x(2))];
②然后在命令窗口求解：
>> x0=[0;0]; >> options=optimset('Display','iter'); >> [x,fv]=fsolve(@myfun,x0,options) %MATLAB 显示的优化过程 %设定求解初值 %设定优化条件 %优化求解 Norm of step
（ 3 ） 利 用 fsolve 求 精 确 解 。 以 求 （0.7926,7843）附近的解为例。 本例直接用字符串表达被解函数。注 意：在此，自变量必须写成 x(1), x(2)。假如 写成 xy(1), xy(2)，指令运行将出错。
fun='[sin(x(1)-x(2)),cos(x(1)+x( 2))]'; xy=fsolve(fun,[x0(2),y0(2)]) xy =
First-order Trust-region Iteration Func-count f(x) optimality radius 0 3 1 2 1 1 6 0.000423308 0.5 0.0617 1 2 9 5.17424e-010 0.00751433 4.55e-005 1.25 3 12 9.99174e-022 1.15212e-005 9.46e-011 1.25 Optimization terminated: first-order optimality is less than options.TolFun. x= 0.4966 0.0067 fv = 1.0e-010 * 0.3161 0.0018