python可以计算复杂积分吗_python – Scipy：加快2D复数积分的计算

使用Cython可以获得10倍的速度,如下所示：

In [87]: %timeit cythonmodule.doit(lam=lam,y0=y0,zxp=zxp,z=z,k=k,ra=ra)

1 loops,best of 3: 501 ms per loop

In [85]: %timeit doit()

1 loops,best of 3: 4.97 s per loop

这可能还不够,坏消息是这可能是

相当接近(可能是最多2个因素)到C / Fortran速度

—如果使用相同的算法进行自适应集成. (scipy.integrate.quad

本身已经在Fortran了.)

要进一步,你需要考虑不同的方法来做

积分.这需要一些思考 – 不能提供太多的东西

我的头顶上现在.

或者,您可以减少积分的容差

被评估.

# Do in Python

#

# >>> import pyximport; pyximport.install(reload_support=True)

# >>> import cythonmodule

cimport numpy as np

cimport cython

cdef extern from "complex.h":

double complex csqrt(double complex z) nogil

double complex cexp(double complex z) nogil

double creal(double complex z) nogil

double cimag(double complex z) nogil

from libc.math cimport sqrt

from scipy.integrate import dblquad

cdef class Params:

cdef public double lam,y0,k,ra

def __init__(self,lam,ra):

self.lam = lam

self.y0 = y0

self.k = k

self.zxp = zxp

self.z = z

self.ra = ra

@cython.cdivision(True)

def integrand_real(double x,double y,Params p):

R1 = sqrt(x**2 + (y-p.y0)**2 + p.z**2)

R2 = sqrt(x**2 + y**2 + p.zxp**2)

return creal(cexp(1j*p.k*(R1-R2)) * (-1j*p.z/p.lam/R2/R1**2) * (1+1j/p.k/R1))

@cython.cdivision(True)

def integrand_imag(double x,Params p):

R1 = sqrt(x**2 + (y-p.y0)**2 + p.z**2)

R2 = sqrt(x**2 + y**2 + p.zxp**2)

return cimag(cexp(1j*p.k*(R1-R2)) * (-1j*p.z/p.lam/R2/R1**2) * (1+1j/p.k/R1))

def ymax(double x,Params p):

return sqrt(p.ra**2 + x**2)

def doit(lam,ra):

p = Params(lam=lam,ra=ra)

rr,err = dblquad(integrand_real,lambda x: -ymax(x,p),lambda x: ymax(x,args=(p,))

ri,err = dblquad(integrand_imag,))

return rr + 1j*ri