当前位置: 动力学知识库 > 问答 > 编程问答 >

numpy - Underflow while using dblquad in Python

问题描述:

I have looked for an answer to this problem, but no questions address the issue of underflow occurrence during the calculation of integrals. I am especially interested in the dblquad function, contained in the NumPy module, in Python.

The problem, essentially is the following: I am trying to calculate a definite double integral, of a rather complicated function:

from math import *

import numpy as np

import scipy.special as special

import scipy.integrate as integrate

#Specific values as example for certain parameters of the general problem

rg=1.5e6

H=5e8*rg

lowerlimit=8.53E+06

theta=0.4

vi=10**14

#A(v,z) is merely a function that appears in the integrand

def A(v,z):

return sqrt(v*z)

#This function is defined in terms of an integral of the modified Bessel function of the 2nd kind. It also appears in the integrand

def F(x):

i=integrate.quad(lambda xi: special.kv(5/3.,xi),x,np.inf)

return x*i[0]

#This is the integrand

def integrand(v,z,x):

value=F(x)*x**(-3/2.)*(A(v,z)**2/x-1/2.)*exp(-A(v,z)/(theta*sqrt(x)))

return value

#Calculation of the integral

def integral(v):

i=integrate.dblquad(lambda x,z: integrand(v,z,x),lowerlimit,H,lambda z:0.,lambda z:np.inf)

return i[0]

#Output

print 'integral = %.2E'%integral(vi)

The terminal output of the above is

integral = 0.00E+00

I know this cannot be true. I do know that the integral must have a finite, non-zero value.

I am not a programmer, and I have not dealt with this issue before. I cannot figure out a way to surpass this problem. If anyone has had experience with similar issues, your help would be much appreciated.

分享给朋友:
您可能感兴趣的文章:
随机阅读: