最优化方法Python计算:一元函数导数计算
定义1 给定连续函数 f ( x ) f(x) f(x), x ∈ Ω ⊆ R x\\in\\Omega\\subseteq\\text{ℝ} x∈Ω⊆R。设 x , x 1 ∈ Ω x,x_1\\in\\Omega x,x1∈Ω,
Δ x = x − x 1 \\Delta x=x-x_1 Δx=x−x1
称为变量 x x x的差分。此时, x = x 1 + Δ x x=x_1+\\Delta x x=x1+Δx。称
f ( x 1 + Δ x ) − f ( x 1 ) Δ x \\frac{f(x_1+\\Delta x)-f(x_1)}{\\Delta x} Δxf(x1+Δx)−f(x1)
为 f ( x ) f(x) f(x)在 x 1 x_1 x1处的前向差分。称
f ( x 1 ) − f ( x 1 − Δ x ) Δ x \\frac{f(x_1)-f(x_1-\\Delta x)}{\\Delta x} Δxf(x1)−f(x1−Δx)
为 f ( x ) f(x) f(x)在 x 1 x_1 x1处的后向差分。称
f ( x 1 + Δ x / 2 ) − f ( x 1 − Δ x / 2 ) Δ x \\frac{f(x_1+\\Delta x/2)-f(x_1-\\Delta x/2)}{\\Delta x} Δxf(x1+Δx/2)−f(x1−Δx/2)
为 f ( x ) f(x) f(x)在 x 1 x_1 x1处的中心差分。
若函数 f ( x ) f(x) f(x)在 x 1 x_1 x1处可微,则 f ( x ) f(x) f(x)在 x 1 x_1 x1处的导数为
f ′ ( x 1 ) = lim Δ x → 0 f ( x ) − f ( x 1 ) Δ x = lim Δ x → 0 f ( x 1 + Δ x ) − f ( x 1 ) Δ x f'(x_1)=\\lim_{\\Delta x\\rightarrow0}\\frac{f(x)-f(x_1)}{\\Delta x}=\\lim_{\\Delta x\\rightarrow0}\\frac{f(x_1+\\Delta x)-f(x_1)}{\\Delta x} f′(x1)=Δx→0limΔxf(x)−f(x1)=Δx→0limΔxf(x1+Δx)−f(x1)
即 f ( x ) f(x) f(x)在 x 1 x_1 x1处的导数 f ′ ( x 1 ) f'(x_1) f′(x1)是 f ( x ) f(x) f(x)在 x 1 x_1 x1处的前向差分 f ( x 1 + Δ x ) − f ( x 1 ) Δ x \\frac{f(x_1+\\Delta x)-f(x_1)}{\\Delta x} Δxf(x1+Δx)−f(x1)当 x x x的差分 Δ x \\Delta x Δx趋于0时的极限。因此,当 ∣ Δ x ∣ |\\Delta x| ∣Δx∣充分小时,
f ′ ( x 1 ) ≈ f ( x 1 ) − f ( x 1 − Δ x ) Δ x . f'(x_1)\\approx \\frac{f(x_1)-f(x_1-\\Delta x)}{\\Delta x}. f′(x1)≈Δxf(x1)−f(x1−Δx).
在 f ( x ) f(x) f(x)在 x 1 x_1 x1处的后向差分 f ( x 1 ) − f ( x 1 − Δ x ) Δ x \\frac{f(x_1)-f(x_1-\\Delta x)}{\\Delta x} Δxf(x1)−f(x1−Δx)中对变量差分 Δ x \\Delta x Δx作变换 Δ x ′ = − Δ x \\Delta x'=-\\Delta x Δx′=−Δx,得
f ( x 1 ) − f ( x 1 − Δ x ) Δ x = f ( x 1 ) − f ( x 1 + Δ ′ x ) − Δ x ′ = f ( x 1 + Δ x ′ ) − f ( x 1 ) Δ x ′ \\frac{f(x_1)-f(x_1-\\Delta x)}{\\Delta x}=\\frac{f(x_1)-f(x_1+\\Delta'x)}{-\\Delta x'}=\\frac{f(x_1+\\Delta x')-f(x_1)}{\\Delta x'} Δxf(x1)−f(x1−Δx)=−Δx′f(x1)−f(x1+Δ′x)=Δx′f(x1+Δx′)−f(x1)
由于 Δ x ′ → 0 \\Delta x'\\rightarrow0 Δx′→0当且进当 Δ x → 0 \\Delta x\\rightarrow0 Δx→0,故当 ∣ Δ x ∣ |\\Delta x| ∣Δx∣充分小时,有
f ′ ( x 1 ) ≈ f ( x 1 ) − f ( x 1 − Δ x ) Δ x . f'(x_1)\\approx\\frac{f(x_1)-f(x_1-\\Delta x)}{\\Delta x}. f′(x1)≈Δxf(x1)−f(x1−Δx).
对于 f ( x ) f(x) f(x)在 x 1 x_1 x1处的中心差分
f ( x 1 + Δ x / 2 ) − f ( x 1 − Δ x / 2 ) Δ x = f ( x 1 + Δ x / 2 ) − f ( x 1 ) + f ( x 1 ) − f ( x 1 − Δ x / 2 ) Δ x = 1 2 f ( x 1 + Δ x / 2 ) − f ( x 1 ) + f ( x 1 ) − f ( x 1 − Δ x / 2 ) Δ x / 2 = 1 2 [ f ( x 1 + Δ x / 2 ) − f ( x 1 ) Δ x / 2 + f ( x 1 ) − f ( x 1 − Δ x / 2 ) Δ x / 2 ] \\frac{f(x_1+\\Delta x/2)-f(x_1-\\Delta x/2)}{\\Delta x}=\\frac{f(x_1+\\Delta x/2)-f(x_1)+f(x_1)-f(x_1-\\Delta x/2)}{\\Delta x}\\\\ =\\frac{1}{2}\\frac{f(x_1+\\Delta x/2)-f(x_1)+f(x_1)-f(x_1-\\Delta x/2)}{\\Delta x/2}\\\\ =\\frac{1}{2}\\left[\\frac{f(x_1+\\Delta x/2)-f(x_1)}{\\Delta x/2}+\\frac{f(x_1)-f(x_1-\\Delta x/2)}{\\Delta x/2}\\right] Δxf(x1+Δx/2)−f(x1−Δx/2)=Δxf(x1+Δx/2)−f(x1)+f(x1)−f(x1−Δx/2)=21Δx/2f(x1+Δx/2)−f(x1)+f(x1)−f(x1−Δx/2)=21[Δx/2f(x1+Δx/2)−f(x1)+Δx/2f(x1)−f(x1−Δx/2)]
所以
f ′ ( x 1 ) ≈ f ( x 1 + Δ x / 2 ) − f ( x 1 − Δ x / 2 ) Δ x . f'(x_1)\\approx\\frac{f(x_1+\\Delta x/2)-f(x_1-\\Delta x/2)}{\\Delta x}. f′(x1)≈Δxf(x1+Δx/2)−f(x1−Δx/2).
令
g ( x ) = f ( x + Δ x / 2 ) − f ( x − Δ x / 2 ) Δ x g(x)=\\frac{f(x+\\Delta x/2)-f(x-\\Delta x/2)}{\\Delta x} g(x)=Δxf(x+Δx/2)−f(x−Δx/2)
称其为 f ( x ) f(x) f(x)的广义导数。对 f ( x ) f(x) f(x)的广义导数 g ( x ) g(x) g(x)计算其在 x 1 x_1 x1处的中心差分,可得 f ( x ) f(x) f(x)在 x 1 x_1 x1处的二阶导数近似计算公式
f ′ ′ ( x 1 ) ≈ f ( x 1 + Δ x ) + f ( x 1 − Δ x ) − 2 f ( x 1 ) Δ x 2 . f''(x_1)\\approx\\frac{f(x_1+\\Delta x)+f(x_1-\\Delta x)-2f(x_1)}{\\Delta x^2}. f′′(x1)≈Δx2f(x1+Δx)+f(x1−Δx)−2f(x1).
利用 f ′ ( x 1 ) f'(x_1) f′(x1)和 f " ( x 1 ) f"(x_1) f"(x1)的近似公式,实现一元函数二阶可微的一、二阶导数计算的Python函数定义如下
def der_scalar(f,x1): #f表示函数,x1表示自变量值dx=1e-7 #设置Δx=0.0000001f1=(f(x1+dx/2)-f(x1-dx/2))/dx #计算一阶导数f2=(f(x1+dx)+f(x1-dx)-2*f(x1))/dx**2 #计算二阶导数return f1,f2
程序中第1~5定义函数der_scalar。两个参数之一f表示二阶可微函数 f ( x ) f(x) f(x),x1表示自变量值 x 1 x_1 x1。第2行设置自变量增量 Δ x = 1 0 − 7 \\Delta x=10^{-7} Δx=10−7为dx。第3行按 f ′ ( x 1 ) f'(x_1) f′(x1)的中心差分近似式计算 x 1 x_1 x1处的一阶导数 f ′ ( x 1 ) f'(x_1) f′(x1)为f1。第4行按广义导数的中心差分计算 x 1 x_1 x1处的二阶导数 f ′ ′ ( x 1 ) f''(x_1) f′′(x1)为f2。
例1 考虑函数 f ( x ) = sin x 2 f(x)=\\sin{x^2} f(x)=sinx2,其一阶导数 f ′ ( x ) = 2 x cos x 2 f'(x)=2x\\cos{x^2} f′(x)=2xcosx2,二阶导数 f ′ ′ ( x ) = 2 ( cos x 2 − 2 x 2 sin x 2 ) f''(x)=2(\\cos{x^2}-2x^2\\sin{x^2}) f′′(x)=2(cosx2−2x2sinx2)。故对 x 1 = π 3 x_1=\\frac{\\pi}{3} x1=3π可算得 f ′ ( x 1 ) = 0.9563 f'(x_1)=0.9563 f′(x1)=0.9563, f ′ ′ ( x 1 ) = − 2.9893 f''(x_1)=-2.9893 f′′(x1)=−2.9893。下列代码比较直接由导函数解析式计算和调用上列程序定义的der_scalar函数计算的结果。
import numpy as np #导入numpy
f=lambda x:np.sin(x**2) #设置函数f(x)
f1=lambda x:2*x*np.cos(x**2) #设置一阶导数f'(x)
f2=lambda x:2*(np.cos(x**2)-2*np.sin(x**2)*x**2) #设置二阶导数f"(x)
x1=np.pi/3 #设置自变量x1
print((f1(x1),f2(x1))) #解析计算
print(der_scalar(f,x1)) #数值计算
程序的2、3、4行分别设置 f ( x ) f(x) f(x)、 f ′ ( x ) f'(x) f′(x)和 f ′ ′ ( x ) f''(x) f′′(x)的解析式。第5行设置自变量的值 x 1 = π 3 x_1=\\frac{\\pi}{3} x1=3π。第6行输出解析计算结果,第7行调用der_scalar函数数值地计算 f ′ ( x 1 ) f'(x_1) f′(x1)和 f ′ ′ ( x 1 ) f''(x_1) f′′(x1)。运行程序,输出
(0.9563079109495481, -2.9893240816612874)
(0.9563079098976839, -3.0000693287648676)
输出的第1行为解析计算的结果,第2行为数值计算结果。比较两者可见der_scalar函数计算的一阶导数结果精度是很高的,二阶导数的计算精度稍差。这是因为,这是对“广义导数”进行差分计算的结果。换句话说,在计算二阶导数时,累积了两次计算误差,产生了误差的“放大”效应。尽管如此,der_scalar函数在大多数场合都能满足应用的要求。