使用 Python 建模牛顿-拉夫森方法
在本教程中,我将向您展示如何借助称为牛顿-拉夫森方法的数值方法来评估多项式或超越方程的根。这是一种迭代方法,我们从一个初始猜测(自变量)开始,然后根据猜测评估 𝑥 的新值。并且该过程持续进行,直到达到收敛。该方法通过下图所示的图表进行解释。
基于 $x_{g}$,评估函数 $(f^{'} \left ( x_{g} \right ))$ 的值。然后在该点绘制一条切线,使其与 𝑥 轴相交于 $x_{n}$。现在我们有两个点 $(x_{g},f\left ( x_{g} \right ))$ 和 $(x_{n} ,0)$。通过这些点经过的直线的斜率可以写成公式 1 所示。
$$\mathrm{f^{'} \left ( x_{g} \right )=\frac{0-f\left ( x_{g} \right )}{x_{n}-x_{g}}}$$
因此,$x_{n}$ 可以计算为 -
$$\mathrm{x_{n}-x_{_{g}}=-\frac{f\left ( x_{g} \right )}{f^{'}\left ( x_{g} \right )}}$$
$$\mathrm{x_{n}=x_{_{g}}-\frac{f\left ( x_{g} \right )}{f^{'}\left ( x_{g} \right )}}$$
现在,这个新的 x 值将作为下一步的猜测。基于这个新的猜测,再次评估函数的值,再次评估斜率,并重复该过程,直到获得根 $(i.e\left | x_{g} -x_{n}\right |<10^{-5})$。
该方法速度很快,但一次只能给你一个根。要获得另一个根,您必须从另一个猜测开始并再次重复该过程。
牛顿-拉夫森方法的 Python 实现
假设我们想找到方程 $x^{2}+4x+3=0$ 的根。牛顿-拉夫森方法的 Python 实现如下 -
导入包 -
from pylab import *
只使用了一个模块,即 pylab,因为它包含 numpy。因此,无需单独导入它。
形成多项式及其导数函数,即 𝑓(𝑥) 和 𝑓'(x)。
f=lambda x: x**2+4*x+3 dfdx=lambda x: 2*x+4
我使用了 'lambda',因为函数中只有一条语句。如果需要,也可以使用 'def' 方法。
使用 "linspace" 函数为 "x" 创建一个数组。
# Array of x x=linspace(-15,10,50)
现在,此步骤是可选的。考虑适当的域绘制函数。我还将向您展示如何绘制切线,以及解决方案如何收敛。因此,如果您对视觉效果感兴趣,则可以遵循此步骤。
# Plotting the function
figure(1,dpi=150)
plot(x,f(x))
plot([-15,10],[0,0],'k--')
xlabel('x')
ylabel('f(x)')
假设 𝑥 的初始猜测以开始第一次迭代。还将误差 $(\left | x_{g}-x_{n} \right |)$ 设置为大于收敛标准的值。在本文中,我将收敛标准设置为 $<10^{-5}$,但您可以根据所需的精度级别进行设置。并将循环计数器设置为 1。
# Initial Guess xg=10 # Setting initial error and loop counter error=1 count=1
在 "for" 循环中,使用上述收敛标准求解公式 (2)。此外,绘制误差和切线。切线绘制在名为 figure(1) 的绘图中,误差在 figure(2) 中。此外,还将 $x_{g}$ 和 $f\left ( x_{g} \right )$ 制成表格以显示不同的值。
# For printing x_g and f(x_g) at different steps
print(f"{'xg':^15}{'f(xg)':^15}")
print("===========================")
# Starting iterations
while error>1.E-5:
# Solving Eq. 1
xn=xg-f(xg)/dfdx(xg)
# Printing x_g and f(x_g)
print(f'{round(xg,5):^15}{round(f(xg),5):^15}')
# Plotting tangents
figure(1,dpi=300)
plot([xg,xn],[f(xg),0])
plot([xn,xn],[0,f(xn)],'--',label=f'xg={xg}')
legend(bbox_to_anchor=(1.01, 1),loc='upper left', borderaxespad=0)
# Evaluating error and plotting it
error=abs(xn-xg)
figure(2,dpi=300)
semilogy(count,error,'ko')
xlabel('Iteration count')
ylabel('Error')
# Setting up new value as guess for next step
xg=xn
# Incrementing the loop counter
count=count+1
# printing the final answer
print("===========================")
print("\nRoot ",round(xn,5))
show()
收敛后,打印根。并显示绘图。
在上述情况下,我将初始猜测设为 10。因此,程序输出将如下所示 -
xg f(xg)
======================================
10 143
4.04167 35.50174
1.10359 8.63228
-0.2871 1.93403
-0.85165 0.31871
-0.99042 0.01926
-0.99995 9e-05
-1.0 0.0
========================================
Root -1.0
误差图如下所示 -
带有切线的函数图在下图中显示。
因此,对应于 $x_{g}=10$,根为 -1。对于第二个根,我们必须更改猜测,将其设为 -10。然后程序输出将如下所示 -
xg f(xg)
===========================
-10 63
-6.0625 15.50391
-4.15433 3.64112
-3.30925 0.71415
-3.03652 0.07438
-3.00064 0.00129
-3.0 0.0
===========================
Root -3.0
现在,误差图将如下所示 -
函数图将如下所示 -
因此,对应于 $x_{g}=-10$,根为 -3。
完整 Python 代码
完整代码如下 -
# Importing module
from pylab import *
# Funciton for f(x) and f'(x)
f = lambda x: x ** 2 + 4 * x + 3
dfdx = lambda x: 2 * x + 4
# Array of x
x = linspace(-15, 10, 50)
# Plotting the function
figure(1, figsize=(7.20, 4.00))
plot(x, f(x))
plot([-15, 10], [0, 0], 'k--')
xlabel('x')
ylabel('f(x)')
# Initial Guess
xg = 10
# Setting initial error and loop counter
error = 1
count = 1
# For printing x_g and f(x_g) at different steps
print(f"{'xg':^15}{'f(xg)':^15}")
print("===========================")
# Starting iterations
while error > 1.E-5:
# Solving Eq. 1
xn = xg - f(xg) / dfdx(xg)
# Printing x_g and f(x_g)
print(f'{round(xg, 5):^15}{round(f(xg), 5):^15}')
# Plotting tangents
figure(1)
plot([xg, xn], [f(xg), 0])
plot([xn, xn], [0, f(xn)], '--', label=f'xg={xg}')
legend(bbox_to_anchor=(0.4, 1.1), loc='upper left', borderaxespad=0)
# Evaluating error and plotting it
error = abs(xn - xg)
figure(2, figsize=(7.20, 4.00))
semilogy(count, error, 'ko')
xlabel('Iteration count')
ylabel('Error')
# Settingup new value as guess for next step
xg = xn
# Incrementing the loop counter
count = count + 1
# printing the final answer
print("===========================")
print("\nRoot ", round(xn, 5))
show()
您可以将代码直接复制到 Jupyter Notebook 中并运行它。
对于您选择的多项式,您可以更改上述代码中所示的函数和导数多项式,并根据您的猜测值,您将获得输出。例如,如果您想找到 𝑥3−sin2(𝑥)−𝑥=0 的根,则在上述代码中,函数及其导数将更改为 -
# Function for f(x) and f'(x) f=lambda x: x**3-(sin(x))**2-x dfdx=lambda x: 3*x**2-2*sin(x)*cos(x)-1
然后对于 1 的猜测值,程序输出将为 -
xg f(xg)
===========================
1 -0.70807
1.64919 1.84246
1.39734 0.36083
1.31747 0.0321
1.30884 0.00036
1.30874 0.0
===========================
Root 1.30874
并且,函数图将如下所示 -
结论
在本教程中,您学习了如何使用牛顿-拉夫森方法求解方程的根。您还学习了如何在 pyplots 中绘制切线并显示根的收敛。
数据结构
网络
RDBMS
操作系统
Java
iOS
HTML
CSS
Android
Python
C 编程
C++
C#
MongoDB
MySQL
Javascript
PHP