數學實驗第四章

1、MATLAB數學實驗第四章 函數和方程第四章 函數和方程4.1 預備知識：零點、極值和最小二乘法4.2 函數零點、極值和最小二乘擬合的MATLAB指令4.3 計算實驗：迭代法4.4 建模實驗：購房貸款的利率和最佳訂貨量 4.1 預備知識：零點非線性方程 f (x) = 0若對于數有f () = 0, 則稱為方程的解或根，也稱為函數f (x)的零點若f () = 0, f ()0 則稱為單根。若有k 1, f () = f () = = f (k-1)() = 0，但f (k)()0 , 稱為k重根非線性方程求解通常用數值方法求近似解，常見的有二分法、牛頓法等非線性方程（組）f (x) = 0，

2、 x=(x1, x2, , xn), f=(f1, f2, , fm) 4.1 預備知識：極值如果對于包含x=a的某個鄰域 ,有 f(a)f(x) (f(a)f(x)對任意x成立, 則稱a為f(x)的一個局部極小(大)值點。如果對任意xD，有f(a)f(x)（f(a)f(x)）成立，則稱a為f(x)在區(qū)域D上的一個全局極小(大)值點。設x為標量或向量，y=f(x)是xD上的標量值函數。4.1 預備知識：極值4.2 函數零點MATLAB指令多項式y(tǒng)=polyval(p,x) 求得多項式p在x處的值y，x可以是一個或多個點p3=conv(p1,p2) 返回多項式p1和p2的乘積p3,r=decon

3、v(p1,p2) p3返回多項式p1除以p2的商，r返回余項x=roots(p) 求得多項式p的所有復根.p=polyfit(x,y,k)用k次多項式擬合向量數據(x, y),返回多項式的降冪系數MATLAB中一個多項式用系數降冪排列向量來表示。例如：多項式x3+2x2-5，在MATLAB指令中表示為1 2 0 5。計算( x3+2x2-5) (x3+2x2-5)，并驗算。 p1=1 2 0 -5; p2=1 -1 2; p,r=deconv(p1,p2)p = 1 3r = 0 0 1 -11驗算 conv(p,p2)+r例1.求多項式x3 + 2 x2 - 5的根 p=1 2 0 -5;

4、x=roots(p) , polyval(p,x) x = -1.6209 + 1.1826i -1.6209 - 1.1826i 1.2419 ans = 1.0e-014 * 0.1776 - 0.5773i 0.1776 + 0.5773i -0.4441 例2.用2次多項式擬合下列數據. x 0.1 0.2 0.15 0 -0.2 0.3 y 0.95 0.84 0.86 1.06 1.50 0.72 解：clear; x=0.1,0.2,0.15,0,-0.2,0.3;y=0.95,0.84,0.86,1.06,1.50,0.72;p=polyfit(x,y,2) % 得到二次擬合多

5、項式% 畫擬合效果圖xi=-0.2:0.01:0.3;yi=polyval(p,xi);plot(x,y,o,xi,yi); Fun=inline(funstr,var) 定義一個Inline 函數,其中funstr是函數的表達式, var是變量名Fun=Mfun 定義一個函數句柄，這里Mfun是函數的M文件表達方式Fun=(var)funstr 定義匿名函數，其中var是變量名, funstr是函數的表達式 4.2 非線性函數的MATLAB表達 x=fzero(Fun, x0) 返回一元函數Fun的一 個零點,其中Fun為函數句柄、 內嵌函數或字符串表達方式。 x0為標量時, 返回函數在x0

6、附近的零點； x0為向量a, b時, 返回在a,b中的零點(要求在a, b的 函數值異號) 4.2 函數零點MATLAB指令 x,f,h=fsolve(Fun, x0) x返回一元或 多元函數Fun在x0附近的一個零點， 其中x0為迭代初值; f返回Fun在x的函數值, 應該接近0; h返回值如果大于0， 說明計算結果可靠， 否則計算結果不可靠。 例3 求函數y=xsin(x2-x-1)在(-2, -0.1)內的 零點 fun=inline(x*sin(x2-x-1),x) fun=inline(x*sin(x2-x-1),x) fzero(fun,-2 -0.1) %由于對參數x0用區(qū)間情形

7、，fzero要求區(qū)間兩端的函數值異號, 所以出現(xiàn)不能直接求解. %先作圖觀察一下fplot(fun,-2 0.1);grid on;%可見在x=-1.6和-0.6附近各有一個零點。我們分兩個小區(qū)間分別求解.fzero(fun,-2,-1.2), fzero(fun,-1.2,-0.1) %可以正確求解fzero(fun,-1.6), fzero(fun,-0.6) %參數x0也可以用一個點x,f,h=fsolve(fun,-1.6)x,f,h=fsolve(fun,-0.6) %也可以用fsolve求解例4 求方程組在原點附近的解xx(1)y x(2)解：若用函數句柄方式，先寫一個M函數%M函

8、數eg4_4fun.mfunction f=fun(x)f(1)=4*x(1)-x(2)+exp(x(1)/10-1;f(2)=-x(1)+4*x(2)+x(1)2/8;注意：x, y要合寫成向量變量xx,f,h=fsolve(eg4_4fun,0 0) % 0,0為初始值, 需m函數eg4_4fun.m%得到解為x=0.2326, y=0.0565, 兩個方程誤差分別為0.0908*10-6和0.1798*10-6，h0說明結果是可靠的。%也可以用下列更便捷的Inline函數或匿名函數方式求解。fun=inline(4*x(1)-x(2)+exp(x(1)/10-1,-x(1)+4*x(2)

9、+x(1)2/8,x); %Inline函數x,f,h=fsolve(fun,0 0) fun=(x)4*x(1)-x(2)+exp(x(1)/10-1,-x(1)+4*x(2)+x(1)2/8; %匿名函數 x,f,h=fsolve(fun,0 0)注1 fzero只能求零點附近變號的根，試以x=1.1為初值，用fzero 和fsolve 求解 (x-1)2=0，看看發(fā)生了什么？注2 fzero 和 fsolve 只能求實根，試用它們解x2+x+1=0，看看發(fā)生了什么？min(y) 返回向量y的最小值max(y) 返回向量y的最大值x,f=fminbnd(fun,a,b) x返回一元函數y=

10、f(x)在a,b內的 局部極小值點,f返回局部極小值 fun為函數句柄或inline。x,f=fminsearch(fun,x0) x返回多元函數y=f(x)在初始值x0 附近的局部極小值點,f返回局部極小值. x, x0均為向量。4.2 函數極值MATLAB指令例 5 .求二元函數f(x,y)= 5-x4-y4+4xy在原點附近的極大值。解：問題等價于求 f(x,y)的極小值 max fmin(-f) x x(1), y x(2)注：在使用fsolve, fminsearch等指令時， 多變量必須合寫成一個向量變量，如用x(1), x(2),。 fun=inline(x(1)4+x(2)4-

11、4*x(1)*x(2)-5); x,g=fminsearch(fun,0,0)4.2 最小二乘擬合MATLAB指令假設已知經驗公式y(tǒng)=f(c,x)(c和x均可為向量), 要求根據一批有誤差的數據(xi,yi), i=0,1,n, 確定參數c.這樣的問題稱為數據擬合。最小二乘法就是求c使得均方誤差最小化 Q(c)=當f關于c是線性函數,問題轉化為一個線性方程組求解，且其解存在唯一。 如果f關于c是非線性函數，問題轉化為函數極值問題c= lsqnonlin (Fun,c0) 使用迭代法搜索最 優(yōu)參數c. 其中Fun是以參數c(可以是向量) 為自變量的函數,表示誤差向量 y-f(c,x)(x, y為

12、數據)， c0為參數c的近似值,作為迭代初值c=lsqcurvefit(Fun2,c0, x, y) 從外部輸入數據， 這里Fun2為兩變量c和x的函數 f(c, x) 非線性最小二乘擬合考慮例4.2，若用lsqnonlin，先寫fitf.mfunction e=fitf(c)x=0.1 0.2 0.15 0 -0.2 0.3;y=0.95 0.84 0.86 1.06 1.50 0.72;e=y-(c(1)*x.2+c(2)*x+c(3);然后執(zhí)行：c=lsqnonlin(fitf,0,0,0)用lsqcurvefit更方便：fun2=inline(c(1)*x.2+c(2)*x+c(3),

13、c,x)x=0.1 0.2 0.15 0 -0.2 0.3;y=0.95 0.84 0.86 1.06 1.50 0.72;c=lsqcurvefit(fun2,0,0,0,x,y) 迭代法是從解的初始近似值x0(簡稱初值)開始，利用某種迭代格式x k+1 = g (x k ),求得一近似值序列x1, x2, , xk, xk+1, 逐步逼近于所求的解(稱為不動點)。最常用的迭代法是牛頓迭代法，其迭代格式為1 迭代法4.3 計算實驗：迭代法例6 求方程 x 2 - 3 x + e x = 2 的正根 (要求精度 = 10 -6)解 令f (x) = x 2 - 3 x + e x - 2, f

14、(0)=-1,當x 2, f (x) 0, f (x) 0即f (x)單調上升，所以根在0,2內。先用圖解法找初值，再用牛頓法程序newton.m求解。M函數newton.mfunction x=newton(fname,dfname,x0,e)%用途：牛頓迭代法解非線性方程f(x)=0%格式：x=nanewton(fname,dfname,x0,e) x返回數值解,% fname和dfname分別表示f(x)及其導函數f(x),% x0為迭代初值，e為精度要求(默認1e-4)。 if nargine&kclear; fun=inline(25.2*(1+r)360-(1+r)360- 1)/

15、r*0.1436 ,r)r=fzero(fun,0.0198/12);R=12*r得年利率為5.53%. (你知道最新利率嗎？)每次訂貨需要收取一定量的生產準備費。沒用完的配件，要在倉庫里儲存一段時間，為此要付出儲存費。若訂貨量很小，則需頻繁定貨，造成生產準備費的增加；反之，若訂貨量很大，定貨周期延長而使生產準備費減少但會造成儲存費的增加。如何確定合適的訂貨量？4.4 建模實驗：最佳訂貨量解 先作一些必要的假設將問題簡化1）汽車工廠對配件的日需求量是恒定的， 每日為r件；2）所訂配件按時一次性交貨， 生產準備費每次k1元；3）儲存費按當日實際儲存量計算， 儲存費每日每件k2元;4）你的工廠不允

16、許缺貨。設一次訂貨x件，則訂貨周期為 T= x/r, 第t天的儲存量為 q(t)= x-r t, 0tT第t天的儲存費為 k2q(t)一個周期的總儲存費為一個周期總費用 C(x) = k1+k2x2/(2r)優(yōu)化目標是使單位產品費用 f(x)=C(x)/x=k1/x+k2x/(2r) 達到最小由f(x)=0 ,即 -k1/x2+k2/(2r)=0 得（經濟批量訂貨公式）線性迭代要么收斂于它的不動點，要么趨于無窮大。不收斂的非線性迭代可能會趨于無窮大，也可能趨于一個周期解，但也有可能在一個有限區(qū)域內雜亂無章地游蕩，這類由確定性運動導致的貌似隨機的現(xiàn)象稱為混沌現(xiàn)象4.4 建模實驗：混沌*昆蟲數量的Logistic模型xk+1 = a x k (1 - x k), 0a4 xk表示第k代昆蟲數量(1表示理想資源環(huán)境最大可能昆蟲數量)。 a為資源系數0a4保證了xk在區(qū)間(0,1)上封閉。*平衡與穩(wěn)定若g () = ,稱為映射g(x)的不動點若對于不動點附近的初始值x0，迭代收斂于此不動點，稱此不動點是穩(wěn)定的當0a1, 在0，1內有一