首页 文章

如何在scipy中创建数学表达式？

提问于
浏览 次
3

我一直在为数学项目探索scipy和scipy的核心包 .

我需要对某些方程进行微积分运算......所以对于学习scipy我决定测试一个简单的方程（PDF of a normal random variable） . 我需要在微积分操作期间保持常量......而不是为它赋值 .

我能够使用sympy成功创建它 . 这是代码：

from sympy import *

x = Symbol('x')
mu = Symbol('mu')
sigma = Symbol('sigma')

def normpdfeqn():
    y = exp(-(((x-mu)**2)) / (2*(sigma**2))) / (sqrt(2*pi*(sigma**2)))
    return y

print(integrate(normpdfeqn(), (x)))

并得到一个适当的输出：

sigma*erf(sqrt(2)*(-mu + x)/(2*sigma))/(2*sqrt(sigma**2))

然后我试着用scipy做同样的事情 .

我一直在阅读http://docs.scipy.org/doc/scipy/reference/tutorial/integrate.html，但无法弄清楚如何为它创建等式 . 这是我到目前为止所尝试的（它与上面的代码几乎相同）：

from sympy import exp
from sympy import sqrt
from sympy import pi
from scipy.integrate import quad
from sympy import Symbol
x = Symbol('x')
mu = Symbol('mu')
sigma = Symbol('sigma')

def integrand():
    y = exp(-((x-mu)**2) / (2*(sigma**2))) / (sqrt(2*pi*(sigma**2)))
    return y


I = quad(integrand(), 0, 1,)
print(I)

代码可能远未完成，我不知道如何使其工作 .

如果我一直在处理上面展示的方程式，我是否应该费心去学习scipy整合？还是继续使用sympy和numpy？

python python-2.7 python-3.x scipy sympy

1 回答

  • 1

    这个简化的代码可能是迈向目标的一步：

    它需要一个函数，这里是一个二次方来保持事物尽可能简单，并评估其封闭形式积分（使用sympy）

    然后它使用scipy在数值上评估0和3之间的积分值 .

    并且使用相同的限制与闭合形式积分 .

    这是一个原型，但似乎做你需要的 .

    from sympy import *
from scipy.integrate import quad

x = Symbol('x')

def func():
    y = x**2 + 2*x +3
    return y

def integrand(expr):
    y = lambda x: eval(expr)
    return y 


y = repr(func())             # a str representation of the original function

res = integrate(func(), (x)) # a sympy representation of the closed form integral 
                             # of the original function

I = quad(integrand(y), 0, 3) # a scipy numerical integration of the original function

print(y)        # sympy repr of original function
print(I)        # scipy numerical integration of the original function
print()
print(res)      # sympy repr of closed form integral of original function
print(integrand(repr(res))(3) - integrand(repr(res))(0)) # closed form integer evaluation

    Resulting output:

    x**2 + 2*x + 3
(27.0, 2.9976021664879227e-13)

x**3/3 + x**2 + 3*x
27.0
    回复于

相关问题