【发布时间】：2015-04-12 18:20:09
【问题描述】：
我正在尝试运行此代码，但标题中出现错误。我查找了有关 tplquad 的文档和示例，但我无法理解我的问题。非常感谢您！
这是我的代码：
from numpy import *
from pylab import *
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from math import *
from scipy.integrate import quad,dblquad,tplquad
from scipy.integrate import nquad
fig_width = 6.
fig_height = fig_width*0.75
fig_size = [fig_width,fig_height]
params = {'backend': 'TkAgg',
'axes.labelsize': 30,
'text.fontsize': 20,
'title.fontsize': 20,
'legend.fontsize': 20,
'xtick.labelsize': 20,
'ytick.labelsize': 20,
'text.usetex': False,
'font.family': 'sans-serif',
'figure.figsize': fig_size}
rcParams.update(params)
pi=3.14
pt_T=3.
#T=0.47
thetaP= -pi
precision=5
y=0
M_T=linspace(1.,7.,precision)
integral1d=[0]*precision
#chi now is def with a plus instead of the minus in the article
def chi(thetap1,p1,thetaP,T,M_T):
return abs((2*p1*T*sqrt(pt_T**2+(M_T**2+pt_T**2)*sinh(y)**2)*sin(thetaP)*sin(thetap1))**2 - (2*p1*T*(sqrt(M_T**2+pt_T**2)*cosh(y)- sqrt(pt_T**2+(M_T**2+pt_T**2)*sinh(y)**2)*cos(thetaP) *cos(thetap1) )-(T**2)* M_T**2)**2)+1
def p1max(thetaP, thetap1,T,M_T):
return (M_T**2)*T/(2*(sqrt(M_T**2+pt_T**2)*cosh(y)- sqrt(pt_T**2+(M_T**2+pt_T**2))*sinh(y)**2*cos(thetaP-thetap1)))-0.1
def p1min(thetaP, thetap1,T,M_T):
#NOT SURE ABOUT THE T AT DENOMINATOR
return (M_T**2)*T/(2*(sqrt(M_T**2+pt_T**2)*cosh(y)- sqrt(pt_T**2+(M_T**2+pt_T**2))*sinh(y)**2*cos(thetaP+thetap1))) +0.1
def integral(thetaP,T,M_T):
area =dblquad(lambda p1, thetap1: 5*(1/(18*pi**5))*sin(thetap1)*(p1/(sqrt(chi(thetap1,p1,thetaP,T,M_T))))*(1/(exp(p1/T) + 1))*(1/(exp((sqrt(M_T**2 + pt_T**2)*cosh(y) - p1/T) +1))) , -pi+0.1, -0.1, lambda p1: p1min(thetaP, p1,T,M_T), lambda p1: p1max(thetaP,p1,T,M_T)) #CHANGE 1., lambda p1:10.)
return area[0]
def integrand(M_T, p1,thetap1,T):
return pt_T*T*2*pi*5*(1/(18*pi**5))*sin(thetap1)*(p1/(sqrt(chi(thetap1,p1,thetaP,T,M_T))))*(1/(exp(p1/T) + 1))*(1/(exp((sqrt(M_T**2 + pt_T**2)*cosh(y) - p1/T) +1)))
def formula151(M_T):
area =tplquad(lambda p1, thetap1,T: pt_T*T*2*pi*5*(1/(18*pi**5))*sin(thetap1)*(p1/(sqrt(chi(thetap1,p1,thetaP,T,M_T))))*(1/(exp(p1/T) + 1))*(1/(exp((sqrt(M_T**2 + pt_T**2)*cosh(y) - p1/T) +1))) ,0.333, 20./3,lambda thetap1: -pi+0.1, -0.1, lambda thetap1, p1: p1min(thetaP, p1,T,M_T),lambda thetap1,p1: p1max(thetaP,p1,T,M_T) )
return area[0]
#solving the integral
for ind in range(0, precision):
integral1d[ind]=formula151( M_T[ind])
print integral1d[ind]
plot(M_T,integral1d)
xlabel('M/T')
ylabel('prod rate')
title('thetaP =-3.12')
plt.yscale('log')
#plt.xscale('log')
show()
错误来自第 57 行，其中使用了 tplquad，完整的回溯是
---------------------------------------------------------------------------
TypeError Traceback (most recent call last)
/usr/lib/python2.7/dist-packages/IPython/utils/py3compat.pyc in execfile(fname, *where)
202 else:
203 filename = fname
--> 204 __builtin__.execfile(filename, *where)
/home/chiara/Scrivania/formula15a.py in <module>()
61 #solving the integral
62 for ind in range(0, precision):
---> 63 integral1d[ind]=formula151( M_T[ind])
64 print integral1d[ind]
65
/home/chiara/Scrivania/formula15a.py in formula151(M_T)
55
56 def formula151(M_T):
---> 57 area =tplquad(lambda p1, thetap1,T: pt_T*T*2*pi*5*(1/(18*pi**5))*sin(thetap1)*(p1/(sqrt(chi(thetap1,p1,thetaP,T,M_T))))*(1/(exp(p1/T) + 1))*(1/(exp((sqrt(M_T**2 + pt_T**2)*cosh(y) - p1/T) +1))) ,0.333, 20./3,lambda thetap1: -pi+0.1, -0.1, lambda thetap1, p1: p1min(thetaP, p1,T,M_T),lambda thetap1,p1: p1max(thetaP,p1,T,M_T) )
58 return area[0]
59
/usr/lib/python2.7/dist-packages/scipy/integrate/quadpack.pyc in tplquad(func, a, b, gfun, hfun, qfun, rfun, args, epsabs, epsrel)
498
499 """
--> 500 return dblquad(_infunc2,a,b,gfun,hfun,(func,qfun,rfun,args),epsabs=epsabs,epsrel=epsrel)
501
502
/usr/lib/python2.7/dist-packages/scipy/integrate/quadpack.pyc in dblquad(func, a, b, gfun, hfun, args, epsabs, epsrel)
433
434 """
--> 435 return quad(_infunc,a,b,(func,gfun,hfun,args),epsabs=epsabs,epsrel=epsrel)
436
437
/usr/lib/python2.7/dist-packages/scipy/integrate/quadpack.pyc in quad(func, a, b, args, full_output, epsabs, epsrel, limit, points, weight, wvar, wopts, maxp1, limlst)
252 args = (args,)
253 if (weight is None):
--> 254 retval = _quad(func,a,b,args,full_output,epsabs,epsrel,limit,points)
255 else:
256 retval = _quad_weight(func,a,b,args,full_output,epsabs,epsrel,limlst,limit,maxp1,weight,wvar,wopts)
/usr/lib/python2.7/dist-packages/scipy/integrate/quadpack.pyc in _quad(func, a, b, args, full_output, epsabs, epsrel, limit, points)
317 if points is None:
318 if infbounds == 0:
--> 319 return _quadpack._qagse(func,a,b,args,full_output,epsabs,epsrel,limit)
320 else:
321 return _quadpack._qagie(func,bound,infbounds,args,full_output,epsabs,epsrel,limit)
/usr/lib/python2.7/dist-packages/scipy/integrate/quadpack.pyc in _infunc(x, func, gfun, hfun, more_args)
379 def _infunc(x,func,gfun,hfun,more_args):
380 a = gfun(x)
--> 381 b = hfun(x)
382 myargs = (x,) + more_args
383 return quad(func,a,b,args=myargs)[0]
TypeError: 'float' object is not callable
注意：函数“integral”和“integrand”已定义但最终没有使用......我只是把它们留在那里
【问题讨论】：
-
哪一行包含错误？你能显示完整的
Traceback -
对不起...我用错误行和回溯编辑了问题！
-
@MarkDickinson 我刚刚尝试了您的建议，但随后出现此错误：TypeError: tplquad() 需要至少 7 个参数（给定 6 个）
-
您可能不应该将所有这些回调和巨大的表达式塞进一行。尝试将这些巨大的表达式分解成有意义的部分，并为这些部分定义临时变量；它会让事情更容易阅读和调试。
-
你真的想让
pi成为3.14吗？您已经导入了 3 次更好的近似值！ （来自math、numpy和pylab中的每一个，尽管pi的pylab对象与numpy的对象相同。）