Take the 2-minute tour ×
Mathematica Stack Exchange is a question and answer site for users of Mathematica. It's 100% free, no registration required.

I am trying to calculate a 2d Integral which has another 2d integral in the integrands. here is the setup:

alpha = 1; sigma1 = 4; sigma2 = 1; phi = 0.5;
a = 0.5*(Cos[phi]/sigma1)^2 + 0.5*(Sin[phi]/sigma2)^2
b = 0.25*Sin[2*phi]*(-1*sigma1^-2 + sigma2^-2)
c = 0.5*(Sin[phi]/sigma1)^2 + 0.5*(Cos[phi]/sigma2)^2
A = alpha/(2*Pi*sigma1*sigma2)

t[x1_, y1_, x_, y_] := A*Exp[-1*(a*(x - x1)^2 + b*(x - x1)*(y - y1) + c*(y - y1)^2)]

B2[x1_, y1_, x_, y_] := HeavisideLambda[x - x1, y - y1]

trans:= NIntegrate[B2[x1, y1, xz, yz]*
          NIntegrate[t[xz, yz, xp, yp]*(B2[x2, y2, xz, yz] - B2[x2, y2, xp, yp]),
             {xp, x2 - 10*sigma1, x2 + 10*sigma1},
             {yp, y2 - 10*sigma2, y2 + 10*sigma2}],
                     {xz, x1 - 1, x1 + 1}, {yz, y1 - 1, y1 + 1}]
trans[0,0,1,0]
trans[0,0,10,10]
trans[0,0,50,50]

I am trying to get numbers for trans[0,0,x,y] for different x,y, but I couldn't get any results yet. I tried analytic integration too (Integrate instead of NIntegrate) but it didn't help either.

I am quite new to Mathematica so I don't know what options in NIntegrate can help me here.

I have tried the same integral in Maple, there it works for some values of phi, and doesn't for others (so i'm sure the thing should some how work here too). for example, it works fine for phi=0, but for phi=0.5 it fails. For phi=Pi/4, it works but the results don't go to zero as x,y get far from zero (e.g. trans[0,0,50,50]), which is the expected behavior for this integral.

any ideas?

Thanks in advance

share|improve this question
    
It looks like all your 4-arguments dependencies are actually only on the differences; for instance, look at your B2. Reducing the number of dimensions might help. –  b.gatessucks Apr 17 at 18:48
1  
your definition trans:= does not create a function that takes four arguments. You should do trans[x1_,y1_..] (Its not immediately obvious what args should be) –  george2079 Apr 17 at 19:00
    
the other thing you need to take that inner integral and turn it into a function that takes only numeric arguments (eq inner[x2_?NumericQ, y2_?NumericQ, xz_?NumericQ, yz_?NumericQ] := NIntegrate[...] . Probably add ?NumericQ to the arguments of t and B2 as well –  george2079 Apr 17 at 19:15
    

Your Answer

 
discard

By posting your answer, you agree to the privacy policy and terms of service.

Browse other questions tagged or ask your own question.