Matlab中的牛顿拉夫森法？

Question

Newtons-Raphsons 方法在 Mathematica 中很容易实现，但在 Matlab 中似乎有点困难。我不知道是否可以将函数传递给函数以及如何将导数用作函数。

newtonRaphson[f_, n_, guess_] := 
 If[n == 0, guess, newtonRaphson[f, n - 1, guess - f[guess]/f'[guess]]]
newtonRaphsonOptimize[f_, n_, guess_] := 
 If[n == 0, guess, 
  newtonRaphsonOptimize[f, n - 1, guess - f'[guess]/f''[guess]]]

似乎您既无法派生函数句柄，也无法派生文件中定义的函数，但我可能是错的。

Answer 1

您可以使用这样的实现：

function x = newton(f,dfdx,x0,tolerance)
err = Inf;
x = x0;
while abs(err) > tolerance
   xPrev = x;
   x = xPrev - f(xPrev)./dfdx(xPrev);
   % stop criterion: (f(x) - 0) < tolerance
   err = f(x); % 
   % stop criterion: change of x < tolerance
   % err = x - xPrev;
end

并将函数及其派生函数的函数句柄传递给它。该导数可以通过一些不同的方法获得：手动微分、符号微分或自动微分。您还可以以数字方式执行微分，但这既慢又需要您使用修改后的实现。所以我假设您已经以任何合适的方式计算了导数。然后就可以调用代码了：

f = @(x)((x-4).^2-4);
dfdx = @(x)(2.*(x-4));
x0 = 1;
xRoot = newton(@f,@dfdx,x0,1e-10);

Answer 2

没有办法以代数方式获取函数句柄或 m 文件中定义的函数的导数。您必须通过在多个点评估函数并逼近导数来以数值方式来完成此操作。

您可能想要做的是符号方程的微分 ，并且您需要符号数学工具箱。下面是使用 Newton-Raphson 方法 查找根的示例：

>> syms x            %# Create a symbolic variable x
>> f = (x-4)^2-4;    %# Create a function of x to find a root of
>> xRoot = 1;        %# Initial guess for the root
>> g = x-f/diff(f);  %# Create a Newton-Raphson approximation function
>> xRoot = subs(g,'x',xRoot)  %# Evaluate the function at the initial guess

xRoot =

    1.8333

>> xRoot = subs(g,'x',xRoot)  %# Evaluate the function at the refined guess

xRoot =

    1.9936

>> xRoot = subs(g,'x',xRoot)  %# Evaluate the function at the refined guess

xRoot =

    2.0000

您可以看到经过几次迭代后，xRoot 的值就接近真根的值（2）。您还可以将函数求值放在 while 循环中，并设置一个条件来检查每个新猜测与前一个猜测之间有多大差异，当差异足够小时停止（即已找到根）：

xRoot = 1;                     %# Initial guess
xNew = subs(g,'x',xRoot);      %# Refined guess
while abs(xNew-xRoot) > 1e-10  %# Loop while they differ by more than 1e-10
  xRoot = xNew;                %# Update the old guess
  xNew = subs(g,'x',xRoot);    %# Update the new guess
end
xRoot = xNew;                  %# Update the final value for the root

Answer 3

% Friday June 07 by Ehsan Behnam.
% b) Newton's method implemented in MATLAB.
% INPUT:1) "fx" is the equation string of the interest. The user 
% may input any string but it should be constructable as a "sym" object. 
% 2) x0 is the initial point.
% 3) intrvl is the interval of interest to find the roots.
% returns "rt" a vector containing all of the roots for eq = 0
% on the given interval and also the number of iterations to
% find these roots. This may be useful to find out the convergence rate
% and to compare with other methods (e.g. Bisection method).
%
function [rt iter_arr] = newton_raphson(fx, x, intrvl)
n_seeds = 10; %number of initial guesses!
x0 = linspace(intrvl(1), intrvl(2), n_seeds);
rt = zeros(1, n_seeds);

% An array that keeps the number of required iterations.
iter_arr = zeros(1, n_seeds);
n_rt = 0;

% Since sometimes we may not converge "max_iter" is set.
max_iter = 100;

% A threshold for distinguishing roots coming from different seeds. 
thresh = 0.001;

for i = 1:length(x0)
    iter = 0;
    eq = sym(fx);
    max_error = 10^(-12);
    df = diff(eq);
    err = Inf;
    x_this = x0(i);
    while (abs(err) > max_error)
        iter = iter + 1;
        x_prev = x_this;

        % Iterative process for solving the equation.
        x_this = x_prev - subs(fx, x, x_prev) / subs(df, x, x_prev);
        err = subs(fx, x, x_this);
        if (iter >= max_iter)
            break;
        end
    end
    if (abs(err) < max_error)
        % Many guesses will result in the same root.
        % So we check if the found root is new
        isNew = true;
        if (x_this >= intrvl(1) && x_this <= intrvl(2))
            for j = 1:n_rt
                if (abs(x_this - rt(j)) < thresh)
                    isNew = false;
                    break;
                end
            end
            if (isNew)
                n_rt = n_rt + 1;
                rt(n_rt) = x_this;
                iter_arr(n_rt) = iter;
            end
        end
    end        
end
rt(n_rt + 1:end) = [];
iter_arr(n_rt + 1:end) = [];