设函数 
     
      
       
       
         f 
        
       
         ( 
        
       
         x 
        
       
         ) 
        
       
      
        f(x) 
       
      
    f(x),在 
     
      
       
       
         [ 
        
       
         a 
        
       
         , 
        
       
         b 
        
       
         ] 
        
       
      
        [a,b] 
       
      
    [a,b]上二阶连续可微且有唯一的最小值点 
     
      
       
        
        
          x 
         
        
          0 
         
        
       
      
        x_0 
       
      
    x0。由于 
     
      
       
       
         f 
        
       
         ( 
        
       
         x 
        
       
         ) 
        
       
      
        f(x) 
       
      
    f(x)是 
     
      
       
       
         [ 
        
       
         a 
        
       
         , 
        
       
         b 
        
       
         ] 
        
       
      
        [a,b] 
       
      
    [a,b]上的单峰函数,故 
     
      
       
        
        
          f 
         
         
         
           ′ 
          
         
           ′ 
          
         
        
       
         ( 
        
       
         x 
        
       
         ) 
        
       
         > 
        
       
         0 
        
       
      
        f''(x)>0 
       
      
    f′′(x)>0, 
     
      
       
       
         x 
        
       
         ∈ 
        
       
         ( 
        
       
         a 
        
       
         , 
        
       
         b 
        
       
         ) 
        
       
      
        x\in(a,b) 
       
      
    x∈(a,b)。对给定 
     
      
       
        
        
          x 
         
        
          k 
         
        
       
         ∈ 
        
       
         [ 
        
       
         a 
        
       
         , 
        
       
         b 
        
       
         ] 
        
       
      
        x_k\in[a,b] 
       
      
    xk∈[a,b],函数 
     
      
       
       
         f 
        
       
         ( 
        
       
         x 
        
       
         ) 
        
       
      
        f(x) 
       
      
    f(x)在 
     
      
       
        
        
          x 
         
        
          k 
         
        
       
      
        x_k 
       
      
    xk的近旁近似为
  
      
       
        
         
         
           q 
          
         
           k 
          
         
        
          ( 
         
        
          x 
         
        
          ) 
         
        
          = 
         
        
          f 
         
        
          ( 
         
         
         
           x 
          
         
           k 
          
         
        
          ) 
         
        
          + 
         
         
         
           f 
          
         
           ′ 
          
         
        
          ( 
         
         
         
           x 
          
         
           k 
          
         
        
          ) 
         
        
          ( 
         
        
          x 
         
        
          − 
         
         
         
           x 
          
         
           k 
          
         
        
          ) 
         
        
          + 
         
         
         
           1 
          
         
           2 
          
         
        
          f 
         
        
          " 
         
        
          ( 
         
         
         
           x 
          
         
           k 
          
         
        
          ) 
         
        
          ( 
         
        
          x 
         
        
          − 
         
         
         
           x 
          
         
           k 
          
         
         
         
           ) 
          
         
           2 
          
         
        
          . 
         
        
       
         q_k(x)=f(x_k)+f'(x_k)(x-x_k)+\frac{1}{2}f"(x_k)(x-x_k)^2. 
        
       
     qk(x)=f(xk)+f′(xk)(x−xk)+21f"(xk)(x−xk)2.
 由于 
     
      
       
        
        
          q 
         
        
          k 
         
        
          ′ 
         
        
       
         ( 
        
        
        
          x 
         
        
          k 
         
        
       
         ) 
        
       
         = 
        
        
        
          f 
         
        
          ′ 
         
        
       
         ( 
        
        
        
          x 
         
        
          k 
         
        
       
         ) 
        
       
      
        q'_k(x_k)=f'(x_k) 
       
      
    qk′(xk)=f′(xk), 
     
      
       
        
        
          q 
         
        
          k 
         
         
         
           ′ 
          
         
           ′ 
          
         
        
       
         ( 
        
        
        
          x 
         
        
          k 
         
        
       
         ) 
        
       
         = 
        
        
        
          f 
         
         
         
           ′ 
          
         
           ′ 
          
         
        
       
         ( 
        
        
        
          x 
         
        
          k 
         
        
       
         ) 
        
       
         > 
        
       
         0 
        
       
      
        q''_k(x_k)=f''(x_k)>0 
       
      
    qk′′(xk)=f′′(xk)>0,故 
     
      
       
        
        
          x 
         
        
          k 
         
        
       
         − 
        
        
         
          
          
            f 
           
          
            ′ 
           
          
         
           ( 
          
          
          
            x 
           
          
            k 
           
          
         
           ) 
          
         
         
          
          
            f 
           
           
           
             ′ 
            
           
             ′ 
            
           
          
         
           ( 
          
          
          
            x 
           
          
            k 
           
          
         
           ) 
          
         
        
       
      
        x_k-\frac{f'(x_k)}{f''(x_k)} 
       
      
    xk−f′′(xk)f′(xk)是 
     
      
       
        
        
          q 
         
        
          k 
         
        
       
         ( 
        
       
         x 
        
       
         ) 
        
       
      
        q_k(x) 
       
      
    qk(x)的驻点,也是唯一的极小值点。
 
 我们用以下策略计算 
     
      
       
       
         f 
        
       
         ( 
        
       
         x 
        
       
         ) 
        
       
      
        f(x) 
       
      
    f(x)的最优解 
     
      
       
        
        
          x 
         
        
          0 
         
        
       
      
        x_0 
       
      
    x0的搜索序列:取 
     
      
       
        
        
          x 
         
        
          1 
         
        
       
         ∈ 
        
       
         [ 
        
       
         a 
        
       
         , 
        
       
         b 
        
       
         ] 
        
       
      
        x_1\in[a,b] 
       
      
    x1∈[a,b]充分接近 
     
      
       
        
        
          x 
         
        
          0 
         
        
       
      
        x_0 
       
      
    x0,从 
     
      
       
       
         k 
        
       
         = 
        
       
         1 
        
       
      
        k=1 
       
      
    k=1开始作迭代,若 
     
      
       
        
        
          f 
         
        
          ′ 
         
        
       
         ( 
        
        
        
          x 
         
        
          k 
         
        
       
         ) 
        
       
         = 
        
       
         0 
        
       
      
        f'(x_k)=0 
       
      
    f′(xk)=0,由 
     
      
       
       
         f 
        
       
         ( 
        
       
         x 
        
       
         ) 
        
       
      
        f(x) 
       
      
    f(x)在 
     
      
       
       
         [ 
        
       
         a 
        
       
         , 
        
       
         b 
        
       
         ] 
        
       
      
        [a,b] 
       
      
    [a,b]上的可微性及单峰性知,找到最优解 
     
      
       
        
        
          x 
         
        
          0 
         
        
       
         = 
        
        
        
          x 
         
        
          k 
         
        
       
      
        x_0=x_k 
       
      
    x0=xk,如上图(a)所示。否则,我们用 
     
      
       
        
        
          q 
         
        
          k 
         
        
       
         ( 
        
       
         x 
        
       
         ) 
        
       
      
        q_k(x) 
       
      
    qk(x)的唯一极小值点 
     
      
       
        
        
          x 
         
        
          k 
         
        
       
         − 
        
        
         
          
          
            f 
           
          
            ′ 
           
          
         
           ( 
          
          
          
            x 
           
          
            k 
           
          
         
           ) 
          
         
         
          
          
            f 
           
           
           
             ′ 
            
           
             ′ 
            
           
          
         
           ( 
          
          
          
            x 
           
          
            k 
           
          
         
           ) 
          
         
        
       
      
        x_k-\frac{f'(x_k)}{f''(x_k)} 
       
      
    xk−f′′(xk)f′(xk)作为第 
     
      
       
       
         k 
        
       
         + 
        
       
         1 
        
       
      
        k+1 
       
      
    k+1个迭代点 
     
      
       
        
        
          x 
         
         
         
           k 
          
         
           + 
          
         
           1 
          
         
        
       
      
        x_{k+1} 
       
      
    xk+1,即
  
      
       
        
         
         
           x 
          
          
          
            k 
           
          
            + 
           
          
            1 
           
          
         
        
          = 
         
         
         
           x 
          
         
           k 
          
         
        
          − 
         
         
          
           
           
             f 
            
           
             ′ 
            
           
          
            ( 
           
           
           
             x 
            
           
             k 
            
           
          
            ) 
           
          
          
           
           
             f 
            
            
            
              ′ 
             
            
              ′ 
             
            
           
          
            ( 
           
           
           
             x 
            
           
             k 
            
           
          
            ) 
           
          
         
        
          . 
         
        
       
         x_{k+1}=x_k-\frac{f'(x_k)}{f''(x_k)}. 
        
       
     xk+1=xk−f′′(xk)f′(xk).
 无论 
     
      
       
        
        
          x 
         
        
          k 
         
        
       
      
        x_k 
       
      
    xk是处于 
     
      
       
       
         f 
        
       
         ( 
        
       
         x 
        
       
         ) 
        
       
      
        f(x) 
       
      
    f(x)的上升区间如上图(b)或处于 
     
      
       
       
         f 
        
       
         ( 
        
       
         x 
        
       
         ) 
        
       
      
        f(x) 
       
      
    f(x)的下降区间如上图(c)所示, 
     
      
       
        
        
          x 
         
         
         
           k 
          
         
           + 
          
         
           1 
          
         
        
       
      
        x_{k+1} 
       
      
    xk+1都有望比 
     
      
       
        
        
          x 
         
        
          k 
         
        
       
      
        x_k 
       
      
    xk更接近 
     
      
       
       
         f 
        
       
         ( 
        
       
         x 
        
       
         ) 
        
       
      
        f(x) 
       
      
    f(x)的最小值点 
     
      
       
        
        
          x 
         
        
          0 
         
        
       
      
        x_0 
       
      
    x0。循环往复,直至遇到 
     
      
       
        
        
          f 
         
        
          ′ 
         
        
       
         ( 
        
        
        
          x 
         
        
          k 
         
        
       
         ) 
        
       
         = 
        
       
         0 
        
       
      
        f'(x_k)=0 
       
      
    f′(xk)=0。
 按此策略的搜索算法称为牛顿方法,代码实现牛顿算法。
from scipy.optimize import OptimizeResult
def newton(fun, x1, gtol, **options):
    xk=x1
    f1,f2=der_scalar(fun,xk)
    k=1
    while abs(f1)>=gtol:
        k+=1
        xk-=f1/f2
        f1,f2=der_scalar(fun,xk)
    bestx=xk
    besty=fun(bestx)
    return OptimizeResult(fun=besty, x=bestx, nit=k)
程序的第2~12行定义实现牛顿算法的函数newton。参数fun表示函数 
     
      
       
       
         f 
        
       
         ( 
        
       
         x 
        
       
         ) 
        
       
      
        f(x) 
       
      
    f(x),x1表示初始点 
     
      
       
        
        
          x 
         
        
          1 
         
        
       
      
        x_1 
       
      
    x1,gtol表示容错误差 
     
      
       
       
         ε 
        
       
      
        \varepsilon 
       
      
    ε,options用来使minimize_scalar将x1和gtol等实际参数传递给newton。第3行将迭代点xk初始化为x1。第4行调用导数计算函数der_scalar详见博文《一元函数导数的数值计算》)计算 
     
      
       
       
         f 
        
       
         ( 
        
       
         x 
        
       
         ) 
        
       
      
        f(x) 
       
      
    f(x)的一、二阶导数置于f1和f2。第5行将迭代次数置为1。第6~9行的while循环执行迭代计算:第7行迭代次数k自增1,第8行计算迭代点 
     
      
       
        
        
          x 
         
        
          k 
         
        
       
         − 
        
        
         
          
          
            f 
           
          
            ′ 
           
          
         
           ( 
          
          
          
            x 
           
          
            k 
           
          
         
           ) 
          
         
         
          
          
            f 
           
           
           
             ′ 
            
           
             ′ 
            
           
          
         
           ( 
          
          
          
            x 
           
          
            k 
           
          
         
           ) 
          
         
        
       
      
        x_k-\frac{f'(x_k)}{f''(x_k)} 
       
      
    xk−f′′(xk)f′(xk)更新xk。第9行再次调用der_scalar计算更新后迭代点处的一、二阶导数。循环往复,直至 
     
      
       
       
         ∣ 
        
        
        
          f 
         
        
          ′ 
         
        
       
         ( 
        
        
        
          x 
         
        
          k 
         
        
       
         ) 
        
       
         ∣ 
        
       
         < 
        
       
         ε 
        
       
      
        |f'(x_k)|<\varepsilon 
       
      
    ∣f′(xk)∣<ε。第10、11行分别计算最优解近似值bestx及最优解处函数近似besty。值12行返回用计算结果besty(最优解处函数值 
     
      
       
       
         f 
        
       
         ( 
        
        
        
          x 
         
        
          0 
         
        
       
         ) 
        
       
      
        f(x_0) 
       
      
    f(x0)的近似值)、bestx(最优解 
     
      
       
        
        
          x 
         
        
          0 
         
        
       
      
        x_0 
       
      
    x0的近似值)和k(迭代次数)构建OptimizeResult类型对象作为返回值。
 例1 用上列程序定义的newton函数,计算函数 
     
      
       
       
         f 
        
       
         ( 
        
       
         x 
        
       
         ) 
        
       
         = 
        
        
        
          x 
         
        
          2 
         
        
       
         + 
        
       
         4 
        
       
         cos 
        
       
          
        
       
         x 
        
       
      
        f(x)=x^2+4\cos{x} 
       
      
    f(x)=x2+4cosx在 
     
      
       
        
        
          x 
         
        
          1 
         
        
       
         = 
        
       
         1.5 
        
       
      
        x_1=1.5 
       
      
    x1=1.5近旁的极小值点。
 解:下列代码完成计算。
import numpy as np                                                  #导入numpy
from scipy.optimize import minimize_scalar                          #导入minimize_scalar
f=lambda x:x**2+4*np.cos(x)                                         #设置目标函数
minimize_scalar(f,method=newton,options={'x1':1.5,'eps':1.48e-8})   #计算最优解
利用代码中的注释信息不难理解本程序。运行程序,输出
 fun: 2.316808419788213
 nit: 7
   x: -1.8954942647118507
读者可将此运行结果与博文《一元函数搜索算法——二分法》中例1对同一函数使用二分法计算的结果相比。在相同的容错误差水平下,牛顿法的效率(仅用7次迭代)显然高于二分法(用27次迭代)。



















