MATLAB应用-求解非线性方程讲课讲稿.doc

MATLAB应用-求解非线性方程

精品文档
第7章 求解非线性方程
7.1 多项式运算在MATLAB中的实现
一、多项式的表达
n次多项式表达为：，是n+1项之和
在MATLAB中，n次多项式可以用n次多项式系数构成的长度为n+1的行向量表示
[a0, a1,……an-1,an]
二、多项式的加减运算
设有两个多项式和。它们的加减运算实际上就是它们的对应系数的加减运算。当它们的次数相同时，可以直接对多项式的系数向量进行加减运算。当它们的次数不同时，应该把次数低的多项式无高次项部分用0系数表示。
例2 计算
a=[1, -2, 5, 3]; b=[0 0, 6, -1]; c=a+b
例3 设，，求f(x)+g(x)
f=[3, -5, 2, -7, 5, 6]; g=[3, 5, -3]; g1=[0, 0, 0, g];%为了和f的次数找齐
f+g1, f-g1
三、多项式的乘法运算
conv(p1,p2)
例4 在上例中，求f(x)*g(x)
f=[3, -5, 2, -7, 5, 6]; g=[3, 5, -3];
conv(f, g)
四、多项式的除法运算
[Q, r]=deconv(p1, p2)
表示p1除以p2，给出商式Q(x)，余式r(x)。Q,和r仍为多项式系数向量
例4 在上例中，求f(x)/g(x)
f=[3, -5, 2, -7, 5, 6]; g=[3, 5, -3];
[Q, r]=deconv(f, g)
五、多项式的导函数
p=polyder(P)：求多项式P的导函数
p=polyder(P,Q)：求P·Q的导函数
[p,q]=polyder(P,Q)：求P/Q的导函数，导函数的分子存入p，分母存入q。
参数P,Q是多项式的向量表示，p,q也是多项式的向量表示。
例4 求有理分式的导函数
P=[3, 5, 0, -8, 1, -5]; %有理分式分子
Q=[10, 5, 0, 0, 6, 0, 0, 7, -1, 0, -100]; %有理分式分母
[p,q]=polyder(P,Q)
六、多项式求根
多项式求根就是求满足多项式p(x)＝0的x值。N次多项式应该有n个根。这些根可能是实根，也可能是若干对共轭复根。其调用格式是
x=roots(P)
其中P为多项式的系数向量，求得的根赋给向量x，即x(1),x(2),…,x(n)分别代表多项式的n个根。
该命令每次只能求一个一元多项式的根，该指令不能用于求方程组的解，必须把多项式方程变成Pn (x) = 0的形式；
例4 求方程的解。
首先将方程变成Pn (x) = 0的形式：
roots([1 -1 0 -1])

例5 求多项式x4+8x3-10的根。
A=[1,8,0,0,-10];
x=roots(A)
若已知多项式的全部根，则可以用poly函数建立起该多项式，其调用格式为：
P=poly(x)
若x为具有n个元素的向量，则poly(x)建立以x为其根的多项式，且将该多项式的系数赋给向量P。
例6 已知 f(x)=3x5+4x3-5x2-7.2x+5
(1) 计算f(x)=0 的全部根。
(2) 由方程f(x)=0的根构造一个多项式g(x)，并与f(x)进行对比。
P=[3,0,4,-5,-7.2,5];
X=roots(P) %求方程f(x)=0的根
G=poly(X) %求多项式g(x)
将这个结果乘以3，就与f(x)一致
7.2 求解非线性方程f ( x ) = 0
方程求根的一般形式是求下列方程的根：
f ( x ) = 0 (l)
实际上，就是寻找使函数 f ( x）等于零的变量x，所以求方程（l）的根，也叫求函数 f ( x）的零点。如果变量x是列阵，则方程（l）就代表方程组。
当方程（l）中的函数 f (x）是有限个指数、对数、三角、反三角或幂函数的组合时，则方程（l）被称为超越方程，例如 e-x - sin（πx / 2 ) ＋lnx = 0 就是超越方程。
当方程（l）中的函数f（x）是多项式时，即 f（x）＝Pn（x）= anxn + an-1xn ＋ … ＋ alx + a0，则方程（l）就成为下面的多项式方程，也称代数方程：
Pn（x）= anxn + an-1xn ＋ … ＋ alx + a0 = 0 ( 2 )
Pn（x）的最高次数n等于2、3时，用代数方法可以求出方程（2）的解析解，但是，当n ≥ 5时，伽罗瓦（Galois）定理已经证明它是没有代数求根方法的。至于超越方程，通常很难求出其解析解。所以，方程（l）的求解经常使用作图法或数值法，而计算机的发展和普及又为这些方法提供了广阔的发展前景，使之成为科学和工程中最实用的方法之一。
本章首先介绍求解 f ( x ) = 0 的 MATLAB 符号法指令，然后介绍求方程数值解的基本原理，最后再介绍求解 f ( x ) = 0 的 MATLAB 数值法指令。
一、符号方程求解
在MATLAB中，求解用符号表达式表示的代数方程可由函数solve实现，其调用格式为：
solve(s)：求解符号表达式s的代数方程，求解变量为默认变量。
当方程右端为0时，方程可以不标出等号和0，仅标出方程的左端。
solve(s,v)：求解符号表达式s的代数方程，求解变量为v。
solve(s1,s2,…,sn,v1,v2,…,vn)：求解符号表达式s1,s2,…,sn组成的代数方程组，求解变量分别v1,v2,…,vn。
例1． 解下列方程。
1．
x= solve('1/(x+2)+4*x/(x^2-4)=1+2/(x-2)', 'x')
2．
f=sym('x-(x^3-4*x-7)^(1/3)=1')
x= solve(f)
3．
x= solve('2*sin(3*x-pi/4)=1')
4．
x= solve('x+x*exp(x)-10', 'x') %仅标出方程的左端
二、求方程f ( x ) = 0数值解的基本方法
并非所有的方程 f ( x ) = 0 都能求出精确解或解析解，不存在这种解的方程就需要用数值解法求出近似解，有几种常见的数值解法基本原理：二分法。
1 求实根的二分法原理
设方程 f (x) =0中的函数 f ( x）为实函数，且满足：
① 函数 f (x）在[ a , b]上单调、连续；
② 方程 f (x) = 0 在（a , b）内只有一个实根 x*。
则求方程 f (x) = 0 的根，就是在（a, b）内找出使f (x）为零的点x*：f (x*) = 0 ，即求函数 f ( x ) 的零点。因为 f (x）单调连续，由连续函数的性质可知，若任意两点aj，bjÎ[ a , b] ，而且满足条件 f (aj) f (bj) < 0 ，则闭区间[aj , bj] 上必然存在方程的根x*，即 x*Î[aj , bj]。
据此原理提出求实根的二分法如下图所示，

图1 方程求根二分法原理示意图
先用中点将区间[a, b]平分为两个子区间 (a,b1）和（b1, b)，方程的根必然在子区间两端点上函数值之积小于零的那一半中，即不在(a ,b1）内，就在（b1 ,b )内，除非 f(b1) = 0 ，于是寻根的范围缩小了一半。图1中的根x*在区间中点左侧，即 x*Î(a , bl）。再将新的含根区间( a , b1）分成两半，重复上述步骤确定出更新的含根子区间。如此重复n次，设含根区间缩小为（an, bn）则方程的根x*Î（an, bn）, 这一系列含根的子区间满足：
( a , b ) D É ( al , bl ) É ( a2 , b2 ) É … É ( a0， b0）É …
由于含根区间范围每次减半，子区间的宽度为 (n = 1,2,….)，显然当n→¥时，(bn一an）→0，即子区间收敛于一点x*，这个点就是方程的根。若n为有限整数，取最后一个子区间的中点 作为方程根的近似值，它满足 f ( xn）≈0 ，于是有：

这就是近似值xn的绝对误差限。假定预先要求的误差为e，由便可以求出满足误差要求的最小等分次数n。
下面是二分法的程序
function [c,err,yc] =bisect (f,a,b,delta)
%Input - f is the function input as a string ‘f’
% - a and b are the left and right end points
%. - delta is the tolerance
%Output - c is the zero
% - yc=f(c)
% - err is the error estimate for c
ya=feval (f,a);
yb=feval (f,b);
if ya*yb>0, break, end %表示无解，结束
maxl=l+round( (log (b-a) -log (delta))/log (2)); %从误差表达式得到最小等分次数n
for k=1:max1
c=(a+b)/2; %取区间中点
yc=feval (f,c);
if yc==0
a=c;
b=c; %这时解已经找到
elseif yb*yc>0
b=c; %区间减半
yb=yc;
else
a=c;
ya=yc;
end
if b-a < delta, break, end
end
c=(a+b)/2;
err=abs(b-a);
yc=feval (f, c)

2 迭代法
迭代法是计算数学中的一种重要方法，用途很广，求解线性方程组和矩阵特征值时也要用到它。这里结合非线性方程的迭代法求解，介绍一下它的基本原理。
迭代法基本原理
迭代法的基本原理就是构造一个迭代公式，反复用它得出一个逐次逼近方程根的数列，数列中每个元素都是方程根的近似值，只是精度不同。
迭代法求解方程
f ( x ) = 0 （1）
时，先把方程等价地变换成形式
f ( x ) = x－g(x) = 0 ， （2）
移项得出：
x = g(x) ( 3 )
若函数g (x）连续，则称（3）为迭代函数。用它构造出迭代公式：
xk+1= g ( xk） , k = 0 , l , 2 ， … ( 4 )
从初始值 x0出发，便可得出迭代序列：
{ x k ｝= x0, x1, x2,….xk,….. （ 5 ）
如果迭代序列（5 ）收敛，且收敛于x*，则由式（4）有：

可见 x*便是方程（l）的根。
迭代法几何意义：
如下图所示，解方程f ( x ) = 0可以等价地变换成求解 x = g ( x )，

图 4 - 2 方程求根迭代法原理示意图
在几何上，就等价求曲线y＝x 和y＝g ( x ）交点P*的坐标 x*。求迭代序列（5) ，就等于从图中x0点出发，由函数y＝g ( x0）得出y＝P0，代入函数y＝x中得出Q1，再把Q1的x坐标 x1代入方程y= g ( x ）得出P1，如此继续下去，便可在曲线y＝g ( x )上得到一系列的点P0，P1， … ，Pk， … ，这些点的x坐标便是迭代数列 xl , x2 ， … ， xk， … ，它趋向于方程（l ) 的根 x*，数列的元素就是方程根的近似值。数列的收敛就等价于曲线y＝x 和y＝g ( x ）能够相交于一点。
迭代公式收敛定理
要想用迭代法求出方程根的近似值，迭代序列（4 - 5）必须收敛。下面的定理给出了迭代法的收敛条件，同时也给出了迭代公式的误差。
收敛定理：方程 x = g ( x ）在（ a , b ）内有根 x*，如果：
① 当xÎ[a,b]时，g( x）Î[a,b];
② g ( x）可导，且存在正数 q < 1，使得对于任意xÎ[a,b]都有|g’( x )|£ q < 1，则有以下结论。
① 方程x = g ( x）在（a , b）内有唯一的根x*。
② 迭代公式 xk+1 = g ( xk）对（a , b）内任意初始近似根 x0均收敛于x*。
③ 近似根 xk的误差估计公式为：

(4 - 6)

3 切线法
切线法就是从函数曲线上的一点出发，不断用曲线的切线代替曲线，求得收敛于根的数列。
切线法原理：
解非线性方程f( x ) = 0 的切线法也称牛顿法，它是把方程线性化的一种近似方法，用函数 f (x）的切线代替曲线产生一个收敛于方程根的迭代序列，从而得到方程的近似根。
把函数 f ( x）在某一初始值 x 。点附近展开成泰勒级数：
(4 - 7)
取其线性部分，近似地代替函数f (x）可得方程的近似式：

设，解该近似方程可得：

把函数 f ( x ）在xl点附近展开成泰勒级数，取其线性部分替代函数f (x)，设，得：

如此继续做下去，就可以得到牛顿迭代公式：
，k = 0 , l , 2 ， … ( 8 )
由式（8）得出的迭代序列 xl , x2 ， … ，xk … ，在一定的条件下收敛于方程的根x*。
2 ．几何意义
图 4 一 3 方程求根切线法原理示意图

选取初值x0后，过点作曲线的切线，其方程为。设切线与X釉的交点为x1，则，再过作切线，与x轴的交点为， 如此不断作切线，求与x轴的交点，便可得出的一系列的交点x1，x2，…，xk，…，它们逐渐逼近方程的根x*。
3 ．切线法的收敛性
理论可以证明，在有根区间[a, b]上，如果、连续且不变号，则只要选取的初始近似根x0满足 f ( x 。），切线法必定收敛。它的收敛速度经推导可得出：
( 9)
是个常数，式(9)表明用牛顿迭代公式在某次算得的误差，与上次误差的平方成正比，可见牛顿迭代公式的收敛速度很快。
4 . 2 . 4 割线法（弦截法）
应用切线法的牛顿迭代公式时，每次都得计算导数，若将该导数用差商代替，就成为割线法（有时称快速弦截法）的迭代公式：
，k = 0 , l , 2 ， … ( 4 一 10 )
割线法的几何意义也很明显。如图 所示，

图 4 一 4 方程求根割线法原理示意图

过点（x0，f( x0））和（x1，f (x1)）作函数y＝f(x）曲线的割线，交X轴于点x2，再过点（x1，f (x1)）和（x2，f (x2)）作曲线的割线，交X轴于点x3，……一直做下去，则得到一系列割线与 X 轴的交点，这些交点序列将趋于方程的根 x*。
非线性方程的数值解法还有许多，这里仅介绍了几种基本方法的原理。
二分法简单方便，但收敛速度慢；
迭代法虽然收敛速度稍微快点，但需要判断能否收敛；
只要初值选取得当，切线法具有恒收敛且收敛速度快的优点，但需要求出函数的导数；
弦截法不需要求导数，特别是前面介绍的快速弦截法，收敛速度很快，但是需要知道两个近似的初始根值才能作出弦，要求的初始条件较多。
这些方法各有千秋，需根据具体情况选用。
三、方程f(x) = 0数值解的MATLAB实现
MATLAB中求方程数值解的办法很多，有的是专用指令，有的是根据方程性质而借用其他专用指令求得的。
4 . 3 . 2 求函数零点指令 fzero
求解方程f ( x ) = 0的实数根也就是求函数f ( x）的零点。MATLAB中设有求函数f (x）零点的指令fzero，可用它来求方程的实数根。该指令的使用格式为：
fzero (fun, x0, options)
① 输入参数fun为函数f (x）的字符表达式、内联函数名或M函数文件名
② 输入参数x0为函数某个零点的大概位置（不要取零）或存在的区间［xi，xj］，要求函数f (x）在x0点左右变号，即f (xi)f (xj) < 0。
③ 输入参数options可有多种选择，若用optimset ('disp', 'iter'）代替options 时，将输出寻找零点的中间数据。
④ 该指令无论对多项式函数还是超越函数都可以使用，但是每次只能求出函数的一个零点，因此在使用前需摸清函数零点数目和存在的大体范围。为此，一般先用绘图指令plot, fplot或ezplot 画出函数f (x）的曲线，从图上估计出函数零点的位置。
例 4 一 5 求方程 x2 + 4sinx) = 25 的实数根（-2π＜x < 2π）。
解：
一、fun为函数f (x）的字符表达式
（l）首先要确定方程实数根存在的大致范围。为此，先将方程变成标准形式f(x) = x2 + 4sin(x) - 25 = 0。作f(x)的曲线图：
x=-2*pi:0.1:2*pi;
f=x.^2+4*sin(x)-25;
plot(x,f);grid on;
从曲线上可以看出，函数的零点大约在x1 ≈ - 4和x2 ≈ 5附近。
(2）直接使用指令 fzero 求出方程在x1 ≈ - 4时的根。
x1= fzero ('x^2+4*sin(x)-25',-4)
若键入：fzero ('x^2+4*sin(x)-25',-4, optimset('disp', 'iter'))，将显示迭代过程。
中间数据表明，求根过程中不断缩小探测范围，最后得出- 4附近满足精度的近似根。
(3）求x2 ≈ 5的根：
x2= fzero ('x^2+4*sin(x)-25',5)
二、fun为函数f (x）的M函数文件名
将方程x2 + 4sin(x) = 25编成M函数文件（实用中在函数较为复杂、而又多次重复调用时，才这样做），用fzero求解。
(1)在M文件编辑调试窗中键入：
function yy = li4_5(x)
yy= x^2+4*sin(x)-25;
以li4_5为文件名存盘，退出编辑调试窗，回到指令窗。
(2）确定根的大体位置;
(3) 在指令窗中键入下述指令可求出 - 4 附近的根：
x1= fzero (' li4_5',-4)
键入下述指令可求出5附近的根：
x2= fzero (' li4_5',5)
三、fun为函数f (x）的内联函数名
内联函数是MATLAB提供的一个对象(Object)。它的性状表现和函数文件一样，但内联函数的创建比较容易。
inline('CE')
'CE'是字符串，CE为不包含赋值符号“=”的表达式。
上式把串表达式转化为输入宗量自动生成的内联函数。
上述调用格式将自动地对CE进行辫识，把CE中由字母／数字组成的连续字符认做变量，除“预定义变量名（如 i , j , pi ) ”和“常用函数名（如 sin )”以外的由字母／数字组成的连续字符将被认做变量。但注意：若连续字符后紧接“左圆括号”，那么将不被当作输入宗量。如 x ( 1 ) ，就不会认做输入宗量处理。
inline('CE',arg1,arg2,…)
上述调用格式把串表达式转化为arg1,arg2等指定输入宗量的内联函数；这种调用格式是创建内联函数的最稳妥、可靠途径。输入宗量字符可表达得更自如。
将函数f (x）写成内联函数的形式：
f = inline ('x^2+4*sin(x)-25')
这时内联函数名为f
分别求x1 ≈ - 4和x2 ≈ 5时的根：
x1= fzero (f,-4), x2= fzero (f,5)
例4-6 求f(x)=x-10x+2=0在x0=0.5附近的根。
从f(x)的曲线看(x=-2.5:0.01:0.5; fx=x-10.^x+2; plot(x,fx))，

曲线的零点有两个，一个在x=-2附近，另一个在x=0.5附近
(1) 建立函数文件funx.m。(function [输出变量列表]=函数名(输入变量列表))
function fx=funx(x)
fx=x-10.^x+2;
(2) 调用fzero函数求根。
z=fzero('funx',0.5)
例2 求在和作为迭代初值时的零点。
从f(x)的曲线看(x=-7:0.01:2; f= x-1./x+5; plot(x,f))，

曲线的零点有两个，一个在x=-5附近，另一个在x=1附近
(1) 建立函数文件fz.m。
function f=fz(x)
f= x-1./x+5;
(2) 调用fz.函数求根。
fzero('fz',-5) %以-5作为迭代初值
fzero('fz',1) %以1作为迭代初值
7. 3 求解非线性方程组数值解的迭代法
一、符号方程组求解
在MATLAB中，求解用符号表达式表示的方程组仍然可由函数solve实现，其调用格式与解用符号表达式表示的方程一样。
例 解下列方程组。
1．
[x y]=solve('1/x^3+1/y^3=28', '1/x+1/y=4', 'x,y')
2．
[u v]=solve('u^3+v^3=98', 'u+v=2', 'u,v')
3．
[x y]=solve('x+y=98', 'x^(1/3)+y^(1/3)=2', 'x,y')
回车后出现下面的提示
Warning: Explicit solution could not be found.
> In D:\MATLAB6p5\toolbox\symbolic\solve.m at line 136
如果做代换：，方程3就变成方程2，就可解
这个问题说明，符号求解并不是万能的。
如果用MATLAB得出无解或未找到所期望的解时，应该用其它方法试探求解。
4．
[x y]=solve('x^2+y^2=5', '2*x^2-3*x*y-2*y^2') %变量由默认规则确定
二、求解非线性方程组的基本方法
对于非线性方程组（以二元方程组为例，其他可以类推）
( 11)
的数值解求法，跟一元非线性方程的切线法（牛顿法）雷同，也是把非线性函数线性化，近似替代原方程得出数值解，所以也叫作牛顿迭代法。
假设方程组（11）的初始估计值为（x0，y0），可以把方程组（11）中的两个函数f1(x，y)和f2(x，y)在（x0，y0）处用二元泰勒级数展开，只取线性部分，移项得出：
(12)
若系数矩阵行列式，则方程组（12）的解为：
，
方程组（11）中的两个函数f1(x，y)和f2(x，y)在(x1，y1)处，再用二元泰勒级数展开，只取线性部分，… …如此继续替代下去，直到方程组的根达到所要求的精度，就完成了方程组的求解。
求解方程组（11）还有许多其他办法，如“最速下降法”，它是利用方程组（11）构成所谓模函数，通过求模函数极小值的方法得到方程组的数值解，诸如此类在此不再一一列举。
三、求方程组数值解的指令
fsolve是用最小二乘法求解非线性方程组 F (X）＝ 0 的指令，变量X可以是向量或矩阵，方程组可以由代数方程或者超越方程构成。它的使用格式为：
fsolve ( 'fun' , X0 , OPTIONS )
① 参数 fun 是编辑并存盘的M函数文件的名称，可以用＠代替单引号对它进行标识。 M函数文件主要内容是方程 F ( X ) = 0 中的函数 F (X) ，即方程左边的函数。

② 参数X0是向量或矩阵，为探索方程组解的起始点。求解将从X0出发，逐渐趋向，最终得到满足精度要求、最接近X0的近似根X* : F ( X*）≈ 0 。由于X0是向量或矩阵，无法用画图方法进行估计，实际问题中常常是根据专业知识、物理意义等进行估计。
③ 该指令输出一个与X0同维的向量或矩阵，为方程组的近似数值解。
④ 参数OPTIONS为设置选项，用它可以设置过程显示与否、误差、算法… …，具体内容可用 help查阅。通常可以省略该项内容。
例 求方程组在x0＝1 ，y0＝1，z0＝1附近的数值解。
解
（l）在文本编辑调试窗中编辑 M函数文件。首先将方程组变换成 F ( X ) 0 的形式， x, y, z看成向量X的三个分量。
function ms=li4_7(X)
ms(1) = X(1).^2 - 10* X(1) + X(2).^2 + X(3) +7;
ms(2) = X(1).*X(2).^2 - 2* X(3);
ms(3) = X(1).^2 + X(2).^2 - 3* X(2) + X(3).^2;
这里，X (l）＝x , X (2）＝y , X (3）＝z，输出参量 ms 也有三个分量
用“li4_7”为M函数文件名存盘。
(2）在指令窗中键入：
fsolve (' li4_7',[1 1 1])
若键入：
x = fsolve (@li4_7,[1 1 1], optimset('Display', 'iter'))
则得出求解过程。
该方程也可以用 MATLAB 的符号指令solve求解，但结果非常冗长。
例 求解方程组，在 x0＝0.1，y0＝0.1和z0＝ - 0.1附近的数值解。
解：
首先将方程组变换成fj (x, y, z）＝f (X）＝0 ( j＝1 , 2 , 3）的形式，设X为一个三维向量，令X (l）＝x , X (2）＝y , X (3）＝z ，则三维向量yy3＝f (X）＝fj (x, y, z），然后编程计算。
(l) 在文本编辑调试窗中编辑M函数文件
function yy3=li4_8(X)
yy3(1) = 3*X(1) - cos(X(2)*X(3)) - 0.5;
yy3(2) = 2*X(1)^2 - 81*(X(2) + 0.1)^2 +sin(X(3)) + 0.06;
yy3(3) = exp(-X(1)*X(2)) + 20*X(3) + 10*pi/3 -1;
(2）在指令窗中键入：
fsolve('li4_8',[0.1 0.1 -0.1])
这个方程组用符号指令solve无法得出最终结果。
例1 求解方程组在(1,1,1)附近的解，并对结果进行验证。
(1) 建立函数文件myfun.m。
function Fmyfun(X) %X、F是性质相同的向量
x=X(1);
y=X(2);
z=X(3);
F(1)= sin(x)+y+z^2.*exp(x);
F(2)=x+y+z;
F(3)=x.*y.*z;
(2) 在给定的初值x0=1,y0=1,z0=1下，调用fsolve函数求方程的根。
X=fsolve('myfun',[1,1,1],optimset('Display','off'))
(3) 将求得的解代回原方程，可以检验结果是否正确，命令如下：
q=myfun(X)
例2 求非线性方程组在(0.5,0.5) 附近的数值解。
(1) 建立函数文件myfun1.m。
function F=myfun1(X)
x=X(1);
y=X(2);
F(1)=x-0.6*sin(x)-0.3*cos(y);
F(2)=y-0.6*cos(x)+0.3*sin(y);
function F=myfun1(X)
x=X(1);
y=X(2);
F(1)=x-0.6*sin(x)-0.3*cos(y);
F(2)=y-0.6*cos(x)+0.3*sin(y);

(2) 在给定的初值x0=0.5,y0=0.5下，调用fsolve函数求方程的根。
X=fsolve('myfun1',[0.5,0.5]',optimset('Display','off'))
将求得的解代回原方程，可以检验结果是否正确，命令如下：
q=myfun1(X)

收集于网络，如有侵权请联系管理员删除