【问题标题】:Using nquad for a double integral使用 nquad 进行二重积分
【发布时间】:2016-05-16 15:39:28
【问题描述】:

这里有问题。到目前为止,这是我的代码:

from scipy import integrate
import math
import numpy as np

a = 0.250
s02 = 214.0
a_s = 0.0163

def integrand(r, R, s02, a_s, a):
        return 2.0 * r * (r/a)**(-0.1) * (1.0 + (r**2/a**2))**(-2.45)\\
*(math.sqrt(r**2 - R**2))**(-1.0) * (a_s/(1 + (R-0.0283)**2/a_s**2 ))

def bounds_R(s02, a_s, a):
        return [0, np.inf]
def bounds_r(R, s02, a_s, a):
        return [R, np.inf]

result = integrate.nquad(integrand, [bounds_r(R, s02, a_s, a), bounds_R(s02, a_s, a)])

a、s02 和 a_s 是常量。我需要对 r 执行第一个积分,然后对 R 执行第二个积分。我认为的问题是 R 出现在第一个积分的限制中(称为 Abel 变换)。尝试了一些东西,每次都得到一个错误,即边界函数中的参数太少或太少。

请帮忙!

【问题讨论】:

  • 我需要计算 Abel 变换,您可以使用 github.com/PyAbel/PyAbel 中已经实现的算法之一,这将比使用 quad 计算效率更高。

标签: python numpy scipy integration


【解决方案1】:

如果你写integrate.nquad(integrand, [bounds_r(R, s02, a_s, a), bounds_R(s02, a_s, a)]),python 期望你影响R 的值。但你没有,因为集成是通过 R 进行的。

这个语法应该可以工作:

result = integrate.nquad(integrand, [bounds_r, bounds_R], args=(s02,a_s,a))

看看documentation of integrate.nquad中的第二个例子。

【讨论】:

  • 谢谢,这最终奏效了,但是当我实际上想将一个值传递给 R 以获得第二个积分时,似乎出现了问题。所以如果我在边界的定义中为 R 输入一个数值,它就可以了,但是我需要在循环中为不同的 R 值执行这个积分,除了重新定义我的边界之外,我想不出一种方法来做到这一点循环中的每一步。
  • 当然可以通过将此值作为被积函数中的参数之一来解决。谢谢!
  • 是的,您需要在被积函数和边界函数中定义一个附加参数,然后在 nquad 的第三个参数 (args=...) 中添加您的迭代值。如果对您有帮助,请支持我的回答;)
  • 看来我还不够先进,无法投票:(
猜你喜欢
  • 2016-12-13
  • 2012-05-07
  • 2020-05-01
  • 1970-01-01
  • 2015-09-03
  • 2018-10-29
  • 2021-05-21
  • 2011-01-27
  • 2019-05-07
相关资源
最近更新 更多