Loading [MathJax]/jax/output/HTML-CSS/jax.js
Ask Your Question
1

2D Gaussian integrals

asked 1 year ago

Unwise gravatar image

updated 1 year ago

FrédéricC gravatar image

Hello. I have the following Gaussian function G(r)=1πR0er22R20,R0>0 of a 2D vector r=(x,y) as well as a 'transfer' function Tk(r,r,d)=eikd+ik2d|rr|2id(2π/k),d,k>0. My objective is to performe the symbolic integration G(r1)=d2r0G(r0)Tk(r) in Sage. First of all: this integral is actually extremely simple to compute by hand, by applying the standard Gaussian formula. However, I'm using this example to learn how to do symbolic integration in Sage.

My first question: is there some particular syntax to do this 2D integral directly, or is it always necessary to do the integrals in dx and dy separately? Following this route, after expanding the modulus squared we get G(r1)=1piR0eikdeik2d(x21+y21)id(2π/k)dx0 e(12R20il2d)x20ikdx0x1dy1 e(12R20il2d)y20ikdy0y1. Since the two integrals are analogous, I focus on the first. The code I've written to try and solve it is the following:

#define all variables
var('x_0, x_1, y_0, y_1, d, k, R_0')
assume(d>0, k>0, R_0>0)

#define constant prefactors
a = 1/(2*R_0^2) - i*k/(2*d)
b = i*k*x_1 / d

#define integrand
I(x_0) = exp(-a*x_0^2 + b*x_0)

#perform Gaussian integral
result = I.integral(x_0, -oo, oo)
show(result)

However, this only returns several warnings of the type

 Warning, need to choose a branch for the root of a polynomial with parameters. This might be wrong.

on top of the error Error trying to find limit of -sqrt(pi)*(-2*i)/(2*i*sageVARR_0^4*sageVARd*sageVARk/(-2*sageVARR_0^2*sageVARd^2+2*sqrt(sageVARR_0^8*sageVARd^2*sageVARk^2+sageVARR_0^4*sageVARd^4))+1)/sqrt(-sageVARR_0^2*sageVARd^2+sageVARR_0^2*sageVARd*sqrt(sageVARR_0^4*sageVARk^2+sageVARd^2))*sageVARR_0^2*sageVARd*exp((-i)*sageVARR_0^2*sageVARk^2*sageVARx_1^2/(2*sageVARR_0^2*sageVARd*sageVARk+2*i*sageVARd^2))/2*erf(-1/(-2*i)*(2*i*sageVARR_0^4*sageVARd*sageVARk/(-2*sageVARR_0^2*sageVARd^2+2*sqrt(sageVARR_0^8*sageVARd^2*sageVARk^2+sageVARR_0^4*sageVARd^4))+1)*sqrt(-sageVARR_0^2*sageVARd^2+sageVARR_0^2*sageVARd*sqrt(sageVARR_0^4*sageVARk^2+sageVARd^2))/sageVARR_0^2/sageVARd*(sageVARx_0+i*sageVARk*sageVARx_1/sageVARd/(i*sageVARR_0^2*sageVARk-sageVARd)*2*sageVARR_0^2*sageVARd/2)).

With a finite choice of extremal points, the integration is performed correctly. I assume the problem here is that Sage has trouble finding the limits of a complex functions, but how do I circumvent the issue?

Preview: (hide)

1 Answer

Sort by » oldest newest most voted
0

answered 1 year ago

achrzesz gravatar image

Your sage code can be replaced by:

var('a b')
assume(a>0)
I0=integrate(exp((-a*x^2+b*x)),x,-oo,oo)
var('x_0, x_1, y_0, y_1, d, k, R_0')
I1=I0.subs({a:1/(2*R_0^2) - i*k/(2*d),b:i*k*x_1 / d})
%display latex
I1
Preview: (hide)
link

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Add Answer

Question Tools

1 follower

Stats

Asked: 1 year ago

Seen: 294 times

Last updated: May 30 '23