1
$\begingroup$

I am attempting to create a python function to assist in calculating the following numerical integration of the Huygens Fresnel integral in the form of def fresnel(x_,y_,x,y,z) where x_ and y_ are $x'$ and $y'$ coordinates.

$$E(x,y,z)=\frac{1}{i\lambda}\iint_{-\infty}^{+\infty}E(x',y',0)\frac{ze^{ikr}}{r}dx'dy'$$

where

$$r=\sqrt{(x-x')^2+(y-y')^2+z^2}$$

$k=\frac{2\pi}{\lambda}$is the wavenumber, and $\lambda$ is the wavelength of the lightsource, and it is assumed that light field is uniform, $E(x',y',0)=1$

The few examples I have seen online use

scipy.integrate.dblquad(func, a, b, gfun, hfun, 
args=(), epsabs=1.49e-08, epsrel=1.49e-08)

to perform the double integration, however, I've been struggling to implement the gfun and hfun callables of the dblquad function as I'm still quite new to this.

$\endgroup$

1 Answer 1

2
$\begingroup$

The limits of integration on dblquad can be functions (callable type) or floats. In your case you need the limits to be infinite.

Following is an example integrating a Gaussian function over $\mathbb{R}^2$, the result should be $\pi$.

from numpy import exp, pi, inf
from scipy.integrate import dblquad


def gaussian(x, y):
    return exp(-x**2-y**2)


inte = dblquad(gaussian, -inf, inf, -inf, inf)

print(inte)
print(pi)

And the result is

(3.141592653589777, 2.5173086244657047e-08)
3.141592653589793
$\endgroup$
0

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.