Scipy库中哪个函数可以求解一元函数的根?
本帖最后由 Rainyboy 于 2010-12-5 16:39 编辑如题,我找了半天,主要参考的是下面三个文档:
NumPy Reference
SciPy Reference Guide
NumPy User Guide
感觉整个Numpy库只是提供了数值计算的基础,最突出的特点就是扩展一系列比Python内建的列表类更好用,更接近数值计算需求的向量类,用它可以方便地声明各种各样的一维或多维数组。
SciPy吧倒是提供了很多工具,但是我找来找去就是找不到对一元函数求根的函数?类似于Matlab中的fzero.
给一段MATLAB中的fzero的用法吧,各位看看Python的这些库类中有没有类似的?没有的话估计我得自己写了……
p0=0.2969;
p1=-0.1260;
p2=-0.3516;
p3=0.2843;
p4=-0.1015;
%翼型内联函数
y_s = inline('p0*(k1*r+0.5)^0.5 + p1*(k1*r+0.5) + p2*(k1*r+0.5)^2 + p3*(k1*r+0.5)^3 + p4*(k1*r+0.5)^4 - k2*r','r',' p0' ,' p1', 'p2' ,'p3', 'p4','k1','k2');
k1 =1;
k2 =2;
r = fzero(y_s ,,[ ], p0,p1,p2,p3,p4,k1,k2);
本帖最后由 Rainyboy 于 2010-12-5 22:42 编辑
好吧……我还是找不到,自己写了一个牛顿法求根……请大家指教!
===================代码分割线======================
# -*- coding: cp936 -*-
#牛顿法求一元函数在制定初值附近的根
#范雨 2010/12/5
def FZeros(func,bound,paralist):
"""对指定的一元函数在指定的点求导数值
func 函数对象,其定义形式应为,(注意,该函数以一个 元组 为参数)
test((x,p1,p2,p3,p4,...))
bound 求根区间及初始值(x_min,x_max,x_init)
paralist (p1,p2,p3,p4,...)
"""
(x_min,x_max,x_init) = bound;
root_last = x_init;
root_now = x_init;
while True:
root_now = root_last - func((root_last,) + paralist)/diff(func,root_last,paralist);
if abs(root_last - root_now) <= 1e-5:
break;
if root_now >= x_min and root_now <= x_max:
root_last = root_now;
else:
raise NoRootException(bound,paralist);
return root_now;
def diff(func,x,paralist):
"""对指定的一元函数在指定的点求导数值
func 函数对象,其定义形式应为,(注意,该函数以一个 元组 为参数)
test((x,p1,p2,p3,p4,...))
x 自变量的值
paralist (p1,p2,p3,p4,...)
"""
dx = 1e-5;
return (func((x + dx,) + paralist) - func((x - dx,) + paralist))/(2*dx);
class NoRootException(Exception):
"""求根过程自定义的异常"""
def __init__(Self,bound,paralist):
Exception.__init__(Self);
Self.ebound = bound;
Self.eparalist = paralist;
if __name__ == '__main__':
def test((x,p1,p2,p3)):
y = p1*x**2 + p2*x + p3;
return y;
try:
print FZeros(test,(5,10,6),(1,-4,3));
except NoRootException ,e:
print "在所设置的区间" , e.ebound , "及参数表" , e.eparalist , "内没有根" ;
===================代码分割线======================
不清楚,似乎没有。。。 回复 3 # wqsong 的帖子
你是说循环试dx的大小么? 本帖最后由 wqsong 于 2010-12-5 23:44 编辑
回复 4 # Rainyboy 的帖子
嗯,调整步长,测试一下精度。。。
不知道python对浮点数的敏感度,也就是浮点数最小值。
设置一个步长在初始值和最小值之间缩放? 回复 5 # wqsong 的帖子
说得倒是,确实应该改改,万一人家被导的函数本身就是1e-5量级的怎么办呢…… Rainyboy 发表于 2010-12-5 23:31 static/image/common/back.gif
回复 5 # wqsong 的帖子
说得倒是,确实应该改改,万一人家被导的函数本身就是1e-5量级的怎么办呢……
嗯嗯,哈哈,本身数量级1e-5。。。 optimize库中的fsolve函数可以用来对非线性方程组进行求解。它的基本调用形式如下:
fsolve(func, x0)
func(x)是计算方程组误差的函数,它的参数x是一个矢量,表示方程组的各个未知数的一组可能解,func返回将x代入方程组之后得到的误差;x0为未知数矢量的初始值。如果要对如下方程组进行求解的话:
f1(u1,u2,u3) = 0
f2(u1,u2,u3) = 0
f3(u1,u2,u3) = 0
那么func可以如下定义:
def func(x):
u1,u2,u3 = x
return
下面是一个实际的例子,求解如下方程组的解:
5*x1 + 3 = 0
4*x0*x0 - 2*sin(x1*x2) = 0
x1*x2 - 1.5 = 0
程序如下: from scipy.optimize import fsolve
from math import sin,cos
def f(x):
x0 = float(x)
x1 = float(x)
x2 = float(x)
return [
5*x1+3,
4*x0*x0 - 2*sin(x1*x2),
x1*x2 - 1.5
]
result = fsolve(f, )
print result
print f(result)输出为:
[-0.70622057 -0.6 -2.5 ]
不好意思,没有看到你们在这里讨论。 回复 9 # smtmobly 的帖子
呵呵,谢谢!高手就是不一样啊,出手不凡,{:{36}:} http://hyry.dip.jp/pydoc
这儿有一本国人写的《用Python做数值计算》很不错,我发在论坛的很多代码是这里来的。大家可以看看,是一本Python方面很好的手册。使用与数值计算哈。 回复 8 # smtmobly 的帖子
出手不凡。。。
页:
[1]