在Python中计算具有可变边界的定积分

Question

我目前一直在计算一个定积分。我需要得到一些数值输出，但是我不能。我做了两次尝试来解这个积分。我如何运行下面这些方法之一?

计算积分的第一种方法是

import numpy as np
from sympy import *
Nt = 17                                         
alpha = .99
t = np.linspace(0, .85, Nt)     
s = np.linspace(0, .85, Nt) 
for k in range(Nt): 
for i in range(k): 
Int = integrate( (t[k] - s) ** - int(alpha), (s, t[i], t[i + 1] )) 
print(Int) 

错误:

ValueError: Invalid limits given: ((array([0.      , 0.053125, 0.10625 , 0.159375, 0.2125  , 0.265625,
0.31875 , 0.371875, 0.425   , 0.478125, 0.53125 , 0.584375,
0.6375  , 0.690625, 0.74375 , 0.796875, 0.85    ]), 0.0, 0.053125),)

我的第二个方法是

import numpy as np
from sympy import *
from numba import jit, prange
Nt = 17                                        
alpha = .99
t = np.linspace(0, .85, Nt)     
s = np.linspace(0, .85, Nt) 
@jit(nopython=True)
def Int(alpha):
Int = 0
for k in prange(Nt): 
for i in prange(k): 
Int = Int + integrate( (t[k] - s) ** - int(alpha), (s, t[i], t[i + 1] )) 
print(Int)   
return Int

没有任何输出。

编辑:对于我的第一种方法，我做了一个小的改变，得到了一些数值结果。我对结果不确定。如有任何意见，我们仍欢迎。

Nt = 17                                        
alpha = .99
t = np.linspace(0, .85, Nt)     
s = symbols('s')
for k in range(Nt): 
for i in range(k): 
Int = integrate((t[k] - s) ** - int(alpha), (s, t[i], t[i + 1] )) 
print(Int) 

输出:

0.0531250000000000
0.0531250000000000
0.0531250000000000
0.0531250000000000...

Answer 1

如果你使用SymPy，你可以得到一个符号表达式:

>>> var('s t e ti tj');i=integrate(1/(t-s)**e,(t,ti,tj));i
Piecewise(
((t - ti)**(1 - e)/(1 - e) - (t - tj)**(1 - e)/(1 - e),
(e > -oo) & (e < oo) & Ne(e, 1)), 
(log(t - ti) - log(t - tj), 
True))

对于e = -int(.99) = 0——你确定这是你想要的吗?——这是

>>> i0 = i.subs(e,-int(.99)); i0
-ti + tj

在这种情况下，积分不依赖于t，只依赖于极限的差。正如奥斯卡所说，你现在可以代入兴趣值。或者，知道它只是点的差，你可以直接计算这些差:

>>> ts = [.85/16*_ for _ in range(17)]
>>> i0.subs({ti:ts[0], tj:ts[1]})
0.053125
>>> ts[1] - ts[0]
0.053125

这和你得到的一样。您现在(希望)明白了为什么每次都得到相同的值。