1
моя проблема связана с квантовой оптикой. Я пытаюсь определить функцию Вигнера для состояния Фока, которая включает интеграцию в комплексную плоскость. Мой код очень медленный, поэтому я хотел бы спросить, можете ли вы помочь мне улучшить скорость?Интеграция области в сложной плоскости с использованием python dblquad
from mpl_toolkits.mplot3d import Axes3D
from matplotlib import cm
import pylab
import numpy as np
import mpmath
import sys
mpmath.dps = 5
from sympy import laguerre_l
import scipy
from scipy.integrate import quad,dblquad
def complex_quadrature(func, a, b, **kwargs):
def real_func(x,y):
return scipy.real(func(complex(x,y)))
def imag_func(x,y):
return scipy.imag(func(complex(x,y)))
real_integral = dblquad(real_func, a, b, lambda x: a, lambda x: b, **kwargs)
imag_integral = dblquad(imag_func, a, b, lambda x: a, lambda x: b, **kwargs)
return complex(real_integral[0], imag_integral[0])
f = lambda z: 1/(np.pi**2)*complex_quadrature(lambda x:scipy.exp(x.conjugate()*z- x*z.conjugate())*complex(laguerre_l(1,0,abs(x)**2),0)*scipy.exp(-abs(x)**2), -5, 5)
fig = pylab.figure()
ax = Axes3D(fig)
X = np.arange(-4, 4, 0.1)
Y = np.arange(-4, 4, 0.1)
X, Y = np.meshgrid(X, Y)
xn, yn = X.shape
W = X*0
for xk in range(xn):
for yk in range(yn):
try:
z = complex(X[xk,yk],Y[xk,yk])
print (f(z))
w = abs(z)
if w != w:
raise ValueError
W[xk,yk] = w
except (ValueError, TypeError, ZeroDivisionError):
pass
ax.plot_surface(X, Y, W, rstride=1, cstride=1, cmap=cm.jet,linewidth=0,alpha=.9)
pylab.show()
Примечание: Большая часть этого была скопирована и изменена с StackOverflow, так что я не претендует на право собственности ...
Я предполагаю, что это будет намного проще и быстрее, по крайней мере, чтобы отобразить комплексная плоскость к двумерному вещественному пространству. – pedda
'if w! = W'? шутки в сторону? Также вставляемый код нуждается в улучшенном отступе – Tobia