超越方程的牛顿法
Newton method for transcendental equation
超越方程:
tan(x)/x + b = 0,其中 b 是任何实数。
我需要引入n并给我这个方程的n个解
我的代码(Python):
from math import tan, cos, pi, sqrt, sin,exp
import numpy as np
from matplotlib.figure import Figure
import matplotlib.pyplot as plt
def f(x,b):
return tan(x)/x + b
def f1(x,b):
return (x/(cos(x)*cos(x)) - tan(x))/(x**2)
e = 0.00000000001
def newtons_method(x0, f, f1, e):
x0 = float(x0)
while True:
x1 = x0 - (f(x0,b) / f1(x0,b))
if abs(x1 - x0) < e:
return x1
x0 = x1
result = []
n = int(input("Input n: "))
b = float(input("Input b: "))
for i in range(2,4*n,1):
result.append(newtons_method(i, f , f1, e))
lambda_result = sorted(list(set(result)))
print(len(lambda_result))
我用 1.The 步更改初始近似值,根以 ~pi 为周期重复,因此第二个参数为 4*n。我通过集合排除重复的解决方案。如果 n 是 50 那么他只能找到 18 个解。需要修复什么才能使其正常工作?请帮助我。
在您的交叉post https://math.stackexchange.com/q/2942400/115115中,Yves Daoust 强烈建议您将牛顿法基于函数
f(x)=sin(x) + b*x*cos(x)
或
f(x)=sin(x)/x + b*cos(x)
因为这些函数没有极点或其他奇点(如果 x 不接近 0)。
至少对于大 c
,解接近初始值 (i+0.5)*pi for i in range(n)
。结果不需要排序或缩短结果。
人们可以在余弦的根处使用正弦项具有交替符号。这使得应用 括号法 的完美情况如 regula falsi (illinois), Dekker's fzeroin, Brent's method,...这些方法几乎与牛顿法一样快,并且保证在区间内找到根。
唯一复杂的是区间 (0,pi/2)
,因为 b<-1
将有一个非零根。必须从 x=0
中删除寻根过程,这对于接近 -1
.
的 b
来说非常重要
牛顿法 只有当初始点离根足够近时,才能可靠地使根为零。如果点离得更远,接近函数的极值,则切线的根可能离得很远。因此,要成功应用牛顿法,需要从一开始就找到好的根近似值。为此,可以使用全局收敛的定点迭代或所考虑函数的结构简单的近似值。
使用收缩不动点迭代:k*pi
附近的解也是等价方程x+arctan(b*x)=k*pi
的根。这给出了近似解x=g(k*pi)=k*pi-arctan(b*k*pi)
。由于弧切线即使对于小 k
也相当平坦,这给出了一个很好的近似值。
如果b<-1
k=0
有一个正根,即在区间(0,pi/2)
内。之前的方法在这种情况下不起作用,因为弧切线在此区间内有一个围绕 1
的斜率。根是由于方程的高阶非线性项,因此需要方程的一种等效形式的非线性近似。近似 tan(x)=x/(1-(2*x/pi)^2)
在 +-pi/2
处给出相同的极点并且两者之间足够接近。将此近似值代入给定方程式并求解得到初始根近似值 x=pi/2*sqrt(1+1/b)
.
在实现中,必须移动 b<-1
的根集以包含附加的第一个解决方案。
from math import tan, cos, pi, sqrt, sin, atan
def f(x,b):
return sin(x)/x+b*cos(x)
def f1(x,b):
return cos(x)/x-(b+x**-2)*sin(x)
e = 1e-12
def newtons_method(x0, f, f1, e):
x0 = float(x0)
while True:
x1 = x0 - (f(x0,b) / f1(x0,b))
if abs(x1 - x0) < e:
return x1
x0 = x1
result = []
n = int(input("Input n: "))
b = float(input("Input b: "))
for i in range(n):
k=i;
if b >= -1: k=k+1
x0 = pi/2*sqrt(1+1/b) if k==0 else k*pi-atan(b*k*pi)
result.append(newtons_method(x0, f , f1, e))
lambda_result = sorted(list(set(result)))
print(len(result), len(lambda_result))
超越方程: tan(x)/x + b = 0,其中 b 是任何实数。 我需要引入n并给我这个方程的n个解
我的代码(Python):
from math import tan, cos, pi, sqrt, sin,exp
import numpy as np
from matplotlib.figure import Figure
import matplotlib.pyplot as plt
def f(x,b):
return tan(x)/x + b
def f1(x,b):
return (x/(cos(x)*cos(x)) - tan(x))/(x**2)
e = 0.00000000001
def newtons_method(x0, f, f1, e):
x0 = float(x0)
while True:
x1 = x0 - (f(x0,b) / f1(x0,b))
if abs(x1 - x0) < e:
return x1
x0 = x1
result = []
n = int(input("Input n: "))
b = float(input("Input b: "))
for i in range(2,4*n,1):
result.append(newtons_method(i, f , f1, e))
lambda_result = sorted(list(set(result)))
print(len(lambda_result))
我用 1.The 步更改初始近似值,根以 ~pi 为周期重复,因此第二个参数为 4*n。我通过集合排除重复的解决方案。如果 n 是 50 那么他只能找到 18 个解。需要修复什么才能使其正常工作?请帮助我。
在您的交叉post https://math.stackexchange.com/q/2942400/115115中,Yves Daoust 强烈建议您将牛顿法基于函数
f(x)=sin(x) + b*x*cos(x)
或
f(x)=sin(x)/x + b*cos(x)
因为这些函数没有极点或其他奇点(如果 x 不接近 0)。
至少对于大 c
,解接近初始值 (i+0.5)*pi for i in range(n)
。结果不需要排序或缩短结果。
人们可以在余弦的根处使用正弦项具有交替符号。这使得应用 括号法 的完美情况如 regula falsi (illinois), Dekker's fzeroin, Brent's method,...这些方法几乎与牛顿法一样快,并且保证在区间内找到根。
唯一复杂的是区间 (0,pi/2)
,因为 b<-1
将有一个非零根。必须从 x=0
中删除寻根过程,这对于接近 -1
.
b
来说非常重要
牛顿法 只有当初始点离根足够近时,才能可靠地使根为零。如果点离得更远,接近函数的极值,则切线的根可能离得很远。因此,要成功应用牛顿法,需要从一开始就找到好的根近似值。为此,可以使用全局收敛的定点迭代或所考虑函数的结构简单的近似值。
使用收缩不动点迭代:
k*pi
附近的解也是等价方程x+arctan(b*x)=k*pi
的根。这给出了近似解x=g(k*pi)=k*pi-arctan(b*k*pi)
。由于弧切线即使对于小k
也相当平坦,这给出了一个很好的近似值。如果
b<-1
k=0
有一个正根,即在区间(0,pi/2)
内。之前的方法在这种情况下不起作用,因为弧切线在此区间内有一个围绕1
的斜率。根是由于方程的高阶非线性项,因此需要方程的一种等效形式的非线性近似。近似tan(x)=x/(1-(2*x/pi)^2)
在+-pi/2
处给出相同的极点并且两者之间足够接近。将此近似值代入给定方程式并求解得到初始根近似值x=pi/2*sqrt(1+1/b)
.
在实现中,必须移动 b<-1
的根集以包含附加的第一个解决方案。
from math import tan, cos, pi, sqrt, sin, atan
def f(x,b):
return sin(x)/x+b*cos(x)
def f1(x,b):
return cos(x)/x-(b+x**-2)*sin(x)
e = 1e-12
def newtons_method(x0, f, f1, e):
x0 = float(x0)
while True:
x1 = x0 - (f(x0,b) / f1(x0,b))
if abs(x1 - x0) < e:
return x1
x0 = x1
result = []
n = int(input("Input n: "))
b = float(input("Input b: "))
for i in range(n):
k=i;
if b >= -1: k=k+1
x0 = pi/2*sqrt(1+1/b) if k==0 else k*pi-atan(b*k*pi)
result.append(newtons_method(x0, f , f1, e))
lambda_result = sorted(list(set(result)))
print(len(result), len(lambda_result))