這篇文章給出不動點及牛頓法迭代是完整實作步驟，沒有考慮運行效率及記憶體損耗，僅供參考，
一、舉例計算
求下列方程的實根
(1) . x^2 - 3x + 2 - exp(x) = 0
(2). x^3 + 2x^2 +10x - 20 = 0
設計一種不動點及牛頓法迭代法，使迭代序列收斂，且計算到x(k) - x(k-1)< 10-8為止，結果輸出小數點 12 位及以上，通過迭代次數及精度等來比較方法的優劣，
二、不動點及牛頓法迭代實作
2.1 、實作環境
MATLAB 2018b
2.2 、實作原理
不動點迭代法的實作:通過合理構造j(x) 函式，滿足不動點迭代區域收斂的條件：即在真值 x* 的鄰域附近使得|j'(x)| < 1，在程式中首先需要大概計算下建構式的導數，看是否滿足條件，以其為后面迭代函式，而在不動點迭代函式MATLAB實作中，撰寫了一個函式myALG(x0, func, x_real, tol, MaxIter) ，輸入參數有初值、構造的j(x) 函式、經過vpasolve 計算得出的真實值、誤差容忍值及最大迭代次數引數，將這些引數傳輸至函式中后，函式通過定義中間變數來存盤迭代值，將初值送入至構造的迭代函式后的輸出賦值給中間變數，然后依據設定的誤差容忍程度來進行回圈，
而牛頓法迭代的實作需要重新構造j(x) 函式：j(x) = x - f(x)/f'(x) ，需要對所要計算的方程進行求導等進行構造，而在 MATLAB的實作也是和不動點迭代法類似，都是函式通過定義中間變數來存盤迭代值，將初值送入至構造的迭代函式后的輸出賦值給中間變數，然后依據設定的誤差容忍程度來進行回圈，兩種方法的具體實作見附錄代碼，
2.3 、資料測驗
方程 1 的實作： x^2 - 3x + 2 - exp(x) = 0
不動點迭代構造j(x) 函式：j(x) = x^2 / 3 + 2 / 3 - exp(x) / 3
牛頓法迭代構造j(x) 函式：j(x) = x - (x^2 - 3x + 2 - exp(x) ) /(2x - 3 - exp(x))初始化myALG(2.0,@ myFunc, real_value,1e-8,20) ，之后通過 MATLAB 進行測驗，經過vpasolve函式計算方程1，得到方程1解為x* =0.25753028543986076045536730493724，之后進行不動點及牛頓法迭代測驗，
不動點迭代得到如下結果:
牛頓法迭代得到如下結果:
方程 2 的實作： x3 + 2x2 +10x - 20 = 0
不動點迭代構造j(x) 函式：j(x) = x + (x^3 + 2x^2 +10x - 20) × (-1/ 20)
牛頓法迭代構造j(x) 函式：j(x) = x - (x^3 + 2x^2 +10x - 20) /(3x^2 + 4x +10)
初始化myALG(1.3@ myFunc, real _ value,1e - 8,15) ，之后通過 MATLAB 進行測驗， 經過 vpasolve 函式計算方程 2，得到方程 2有幾個不同解， 以解x* =1.3688081078213726352274143300213 為例來進行測驗，
不動點迭代得到如下結果:
牛頓法迭代得到如下結果:
2.4 、資料分析
通過構造了合適的不動點及牛頓法迭代，實作了方程 1 與方程 2 的求解，并使之與 xk - xk -1 < 10-8 及結果輸出小數點 12 位及以上等要求，經過測驗可以看出在一定條件下牛頓法收斂的速度要比不動點迭代法的收斂速度要快，甚至相差一定數量級，這與數值分析課本中對牛頓法及不動點迭代法收斂速度相對應，
附錄
不動點迭代代碼：
math_test_BDD.m
1 % syms x 2 % real_value=https://www.cnblogs.com/iFrank/p/vpasolve(x^2-3*x+2-exp(x) == 0, x) %不動點迭代方程式1 3 % cal_value=https://www.cnblogs.com/iFrank/p/myALG(2.0,@myFunc,real_value,1e-8,15); 4 % 5 % if abs(cal_value - real_value)< 1e-8 6 % fprintf('abs(cal_value - real_value)*10^8=%f\n',abs(cal_value - real_value)*1e8); 7 % end 8 % abs(cal_value - real_value) 9 10 11 syms x 12 real_value=https://www.cnblogs.com/iFrank/p/vpasolve(x^3 + 2*x^2 + 10*x - 20 == 0, x) %不動點迭代方程式2 13 real_value_r = real_value(1) 14 cal_value=https://www.cnblogs.com/iFrank/p/myALG_BDD(1.3,@myFunc_BDD,real_value_r,1e-8,15); 15 16 if abs(cal_value - real_value_r)< 1e-8 17 fprintf('abs(cal_value - real_value_r)*10^8=%f\n',abs(cal_value - real_value_r)*1e8); 18 end 19 abs(cal_value - real_value_r)
myALG_BDD.m
1 function [y] = myALG_BDD(x0,func,x_real,tol,MaxIter) 2 3 xn = x0; 4 fprintf('Iter 0: %16.14f\n',x0); 5 xnp1 = func(xn); 6 fprintf('Iter 1: %16.14f\n',xnp1); 7 %criterion = abs(xnp1-xn); 8 criterion = abs(xnp1-x_real); 9 xn = xnp1; 10 Iter = 1; 11 %while(criterion>tol) 12 while(abs(criterion)>tol) 13 xnp1 = func(xn); 14 %criterion = abs(xnp1-xn); 15 criterion = abs(xnp1-x_real); 16 xn = xnp1; 17 Iter = Iter + 1; 18 fprintf('Iter %2.0d: %16.14f\n',Iter,xnp1); 19 if Iter>=MaxIter 20 break; 21 end 22 end 23 y = xnp1; 24 end
myFunc_BDD.m
1 function [y] = myFunc_BDD(x) 2 3 %y = x^2-3*x+2-exp(x); 4 %y = x^2/3+2/3-exp(x)/3; %不動點迭代方程式1 5 6 %y = x^3 + 2*x^2 + 10*x - 20 y = -x^3/10 - x^2/5 + 2 7 y = x + (x^3 + 2*x^2 + 10*x - 20)*(-1/20); %不動點迭代方程式2 8 end
牛頓法迭代代碼：
math_test_ND.m
1 % syms x 2 % real_value=https://www.cnblogs.com/iFrank/p/vpasolve(x^2-3*x+2-exp(x) == 0, x) %不動點迭代方程式1 3 % cal_value=https://www.cnblogs.com/iFrank/p/myALG(2.0,@myFunc,real_value,1e-8,15); 4 % 5 % if abs(cal_value - real_value)< 1e-8 6 % fprintf('abs(cal_value - real_value)*10^8=%f\n',abs(cal_value - real_value)*1e8); 7 % end 8 % abs(cal_value - real_value) 9 10 11 % syms x 12 % real_value=https://www.cnblogs.com/iFrank/p/vpasolve(x^3 + 2*x^2 + 10*x - 20 == 0, x) %不動點迭代方程式2 13 % real_value_r = real_value(1) 14 % cal_value=https://www.cnblogs.com/iFrank/p/myALG(1.3,@myFunc,real_value_r,1e-8,15); 15 % 16 % if abs(cal_value - real_value_r)< 1e-8 17 % fprintf('abs(cal_value - real_value_r)*10^8=%f\n',abs(cal_value - real_value_r)*1e8); 18 % end 19 % abs(cal_value - real_value_r) 20 21 22 23 % syms x 24 % real_value=https://www.cnblogs.com/iFrank/p/vpasolve(x^2-3*x+2-exp(x) == 0, x) %牛頓法迭代1 25 % real_value_r = real_value 26 % cal_value=https://www.cnblogs.com/iFrank/p/myALG(2,@myFunc,real_value_r,1e-8,15); 27 % 28 % if abs(cal_value - real_value_r)< 1e-8 29 % fprintf('abs(cal_value - real_value_r)*10^8=%f\n',abs(cal_value - real_value_r)*1e8); 30 % end 31 % abs(cal_value - real_value_r) 32 33 34 syms x 35 real_value=https://www.cnblogs.com/iFrank/p/vpasolve(x^3 + 2*x^2 + 10*x - 20 == 0, x) %牛頓法迭代2 36 real_value_r = real_value(1) 37 cal_value=https://www.cnblogs.com/iFrank/p/myALG_ND(1.3,@myFunc_ND,real_value_r,1e-8,15); 38 39 if abs(cal_value - real_value_r)< 1e-8 40 fprintf('abs(cal_value - real_value_r)*10^8=%f\n',abs(cal_value - real_value_r)*1e8); 41 end 42 abs(cal_value - real_value_r)
myALG_ND.m
1 function [y] = myALG_ND(x0,func,x_real,tol,MaxIter) 2 3 xn = x0; 4 fprintf('Iter 0: %16.14f\n',x0); 5 xnp1 = func(xn); 6 fprintf('Iter 1: %16.14f\n',xnp1); 7 %criterion = abs(xnp1-xn); 8 criterion = abs(xnp1-x_real); 9 xn = xnp1; 10 Iter = 1; 11 %while(criterion>tol) 12 while(abs(criterion)>tol) 13 xnp1 = func(xn); 14 %criterion = abs(xnp1-xn); 15 criterion = abs(xnp1-x_real); 16 xn = xnp1; 17 Iter = Iter + 1; 18 fprintf('Iter %2.0d: %16.14f\n',Iter,xnp1); 19 if Iter>=MaxIter 20 break; 21 end 22 end 23 y = xnp1; 24 end
myFunc_ND.m
1 function [y] = myFunc_ND(x) 2 3 4 5 %y = x^2-3*x+2-exp(x); 6 %y = x^2/3+2/3-exp(x)/3; %不動點迭代方程式1 7 8 %y = x^3 + 2*x^2 + 10*x - 20 y = -x^3/10 - x^2/5 + 2 9 %y = x + (x^3 + 2*x^2 + 10*x - 20)*(-1/20); %不動點迭代方程2 10 11 %y = x-(x^2-3*x+2-exp(x))/(2*x-3-exp(x)); %牛頓法迭代1 12 13 y = x - (x^3 + 2*x^2 + 10*x - 20)/(3*x^2+4*x+10); %牛頓法迭代2 14 15 16 end
轉載請註明出處，本文鏈接：https://www.uj5u.com/qita/256220.html
標籤：其他
下一篇：Jacobi與SOR迭代法的實作與性能比較及均勻間距與Chebyshev插值的實作、性能分析及二者生成的插值誤差比較