计算数学
非线性方程数值解法
Newton-Raphson迭代法
Newton-Rahpson迭代方法是解非线性方程比较著名而且也比较有效的方法之一。 如果初值比较接近根,收敛速度是很快的。Newton-Rahpson迭代法也是工程上广泛采用的方法。 Newton迭代法有比较直观的几何意义。
函数方程\(f(x) = 0\)的根是曲线\(y = f(x)\) 与\(x\)轴的交点的横坐标。 通过当前的点做曲线的切线(与一次导数有关),切线方程为 \[y = f({x_k}) + f'({x_k})(x - {x_k})\]
在上式中令\(y=0\)得到切线与\(x\)轴交点的横坐标 \[x = {x_k} - \frac{{f({x_k})}}{{f'({x_k})}}\] 于是我们可以构造迭代公式 \[{x_{k + 1}} = {x_k} - \frac{{f({x_k})}}{{f'({x_k})}}\] 如此反复迭代便可以比较快的得到方程的根,所以Newton法又叫切线法。
关于Newton迭代法如果函数 \(f(x)\)有二阶以上连续倒数,\(p\)是方程\(f(x) = 0\)的单根, 则当\(x\) 充分接近\(p\)时,Newton迭代法收敛,而且最少二阶收敛。
[例] 用Newton-Raphson方法计算方程 \[f(x) = {x^3} + 2{x^2} + 10x - 20 = 0\] 在1.5附近的根。下面给出go语言源代码。
// newton
package main
import (
"fmt"
"math"
)
func main() {
/*------------------------------------------------------
! Author : Song Yezhi
! verison : 2020-5-18 10:29
! -----------------------------------------------------
! Input Parameters :
!
! Output Parameters :
!
------------------------------------------------------*/
var x0 float64 = 0.0
x, fx := newton(x0,myfx,dmyfx)
//把函数作为参数传入
fmt.Printf("the result: \n")
fmt.Printf("x= %12.7f, f(x)= %12.7f \n", x, fx)
}
func myfx(x float64) float64 {
/*------------------------------------------------------
! Author : Song Yezhi
! verison : 2020-5-18 10:28
!
------------------------------------------------------*/
fx := x*x*x + 2.0*x*x + 10.0*x - 20.0
return fx
}
func dmyfx(x float64) float64 {
/*------------------------------------------------------
! Author : Song Yezhi
! verison : 2020-5-18 10:29
!
------------------------------------------------------*/
df := 3.0*x*x + 4.0*x + 10.0
return df
}
func newton(x0 float64,funcX func(float64) float64,
dfuncX func(float64) float64) (x1, fx float64) {
/*------------------------------------------------------
! Author : Song Yezhi
! verison : 2020-5-18 10:30
! 牛顿法计算方程根
! -----------------------------------------------------
! Input Parameters :
! x0-----初值
! funcX------目标函数
! dfuncX-----导函数
! Output Parameters :
!
-------------------------------------------------------*/
var imax int = 20
var tol float64 = 1e-8
for i := 0; i < imax; i++ {
x1 = x0 - funcX(x0)/dfuncX(x0)
fx = funcX(x0)
dfx := dfuncX(x0)
x1 = x0 - fx/dfx
dx := math.Abs(x1 - x0)
if dx < tol {
break
}
x0 = x1
fmt.Printf("i= %4d x= %12.7f f(x)=%12.7f \n", i, x1, fx)
}
return x1,fx
}