在python中集成一个返回相同大小数组而不是元组的数组？

Question

对于冗长的标题，我深表歉意。我已经在一个项目上工作了一段时间，并且我在代码中的某个部分陷入了困境。我会尽力做到彻底。

我有一个大小为 167197 的质量 M 数组。

## Non-constant
M = data['m200'] # kg // Mass of dark matter haloes

R = [] # Km // Radius of sphere 
for masses in M:
    R.append(((3*masses)/(RHO_C*4*(3.14))**(1.0/3.0)))

我有拟合函数，k 的独立值是我的问题的一部分。 k 是我的代码中定义的值。

def T(k): # Fitting Function // Assuming a lambdaCDM model
    q = k/((OMEGA_M)*H**2)*((T_CMB)/27)**2
    L = np.log(euler+1.84*q)
    C = 14.4 + 325/(1+60.5*q**1.11)
    return L/(L+C*q**2)

##############################################################################

def P(k): # Linear Power Spectrum
    A = 0.75 # LambdaCDM Power Normalization
    n = 0.95 # current constraints from WMAP+LSS
    return A*k**n*T(k)**2

* 针对实际问题 *

我有一个 W(kR) 的傅里叶变换

def W(R):# Fourier Transfrom of Top Hat function
    return (3*(np.sin(k*R)-(k*R)*np.cos(k*R)))/(k*R)**(3)
W_a = [] 
for radii in R:
     W_a.append(W(radii))

在这种情况下，我将 R 视为独立值而不是 kR 组合

打印 W_a 的长度给了我与 mu numpy 数组完全相同的大小，所以一切都很好。

此函数将与 sigma 的此函数中包含的积分一起发挥作用

def sigma(R): # Mass Varience
    k1 = lambda k: k**2*P(k)*W(R)**2
    norm1 = 1/(2*np.pi**2)
    return (integrate.quad(k1, 0, np.Inf))
sigma_a = []
for radii in R:
    sigma_a.extend(sigma(radii))

当然，积分会创建一个元组。但是对于 R 中的每个值。我想创建一个列表或数组。因此，当使用.extend() 时，我的数组的长度现在翻了一番，现在的长度为 334394。

如何将其校正到积分计算 W(kR) 中的每个 R 返回相同大小的数组 167197 的位置？

Answer 1

首先只是一个 Python 注释：

R = [] # Km // Radius of sphere 
for masses in M:
    R.append(((3*masses)/(RHO_C*4*(3.14))**(1.0/3.0)))

可以表示为：

R = [((3*masses)/(RHO_C*4*(3.14))**(1.0/3.0)) for masses in M]

在：

return (integrate.quad(k1, 0, np.Inf))

() 的外部集合没有区别。

return integrate.quad(k1, 0, np.Inf)

应该返回相同的东西。

现在翻倍从何而来？在 quad 文档中，我们看到它返回 2 个值，即整数和错误项。这在某些示例中显示为元组，但在其他示例中也未打包：

y, err = integrate.quad(f, 0, 1, args=(3,))

如果你只想要整数而不是错误，你可以索引，integate...()[0]。

sigma_a = []
for radii in R:
    sigma_a.append(sigma(radii)[0])

或

sigma_a = [sigma(radii)[0] for radii in R]

或

def sigma1(R): # Mass Varience
    k1 = lambda k: k**2*P(k)*W(R)**2
    norm1 = 1/(2*np.pi**2)
    y, err = integrate.quad(k1, 0, np.Inf)
    return y # return just the integral

sigma_a = [sigma1(radii) for radii in R]

如果您想同时收集 y 和 err，但在单独的列表中，请使用 zip* 重新打包它们（类似于 numpy 转置）。

ll = [sigma(radii) for radii in R]
# [(y0,err0),(y1,err1), ...]
ys, errs = zip(*ll)