4
$\begingroup$

I've spent the last week searching for a way to numerically calculate the pole contributions to this multi-dimensional complex contour integral while avoiding the singularity at z=0: $$ \oint_{C_1(0)}\oint_{C_p}\frac{z_{1}z_{2}}{z_{2}^{-1}p+(1-p)z_{1}-1} \, e^{\frac{1}{z_{1}} + z_{1} + \frac{1}{z_{2}} + z_{2}} \quad \text d z_1 \, \text d z_2 $$ where $p \in (0,1)$, $C_1(0)$ is a circle of radius 1 centered at 0 and $C_p$ is a contour surrounding the poles only.

I used Finding residues of multi-dimensional complex functions to compute this integral. When I give the whole integral to Mathematica, I don't have the contribution of the poles because the contour cannot be deformed due to the singularity. So I'm trying to find a way to exclude the singularity of the contour, maybe with a small contour around $z = 0$ using Piecewise or something similar to get only the contribution of the poles.

Code

NContourIntegrate[f_, par : (z_ -> g_), {t_, a_, b_}] := 
  NIntegrate[Evaluate[D[g, t] (f /. par) /. t -> t1], {t1, a, b}]  

Clear[Pinfz1];  
Pinfz1[p_?NumericQ, z2_?NumericQ] := 
  NContourIntegrate[1/(2*I*Pi)^2*z1*z2/(p + (1 - p)*z1*z2-z2)*E^(1/z1 + z1 + 1/z2 + z2), 
   z1 -> Exp[I t], 
   {t, 0, 2*Pi}]  

Clear[Pinfz];  
Pinfz[p_?NumericQ] := NContourIntegrate[Pinfz1[p, z2], z2 -> Exp[I t], {t, 0, 2*Pi }]  
$\endgroup$
  • 1
    $\begingroup$ Where does this computation arise? I'm wondering if a double contour integration is even what's actually wanted. $\endgroup$ – Daniel Lichtblau Jun 18 '13 at 16:03
  • $\begingroup$ This come from a physical situation (particle diffusion on lattice), so there is actually 2 integrals $\endgroup$ – Jprog Jun 18 '13 at 22:54
3
$\begingroup$

I'm not familiar with the physics so I cannot say whether integrating over this set of 2 real dimension, in C^2, is what is wanted. I think the code below will cover the product space of the contours that are requested.

ii = z1*z2/(p/z2 + (1 - p) z1 - 1)*Exp[1/z1 + z1 + 1/z2 + z2];
i1 = (ii /. {z1 -> Exp[I*t1], z2 -> Exp[I*t2]})*I*Exp[I*t1]*I*
   Exp[I*t2];
f[pval_] := 
 Quiet[NIntegrate[
   Evaluate[i1 /. p -> pval], {t1, 0, 2*Pi}, {t2, 0, 2*Pi}]]
f[1/2]
f[1/3]

(* Out[114]= 89.3519154317 - 5.40968238685*10^-8 I

Out[115]= 62.3532282953 + 0.000124769719229 I *)

Also:

f[.1]

(* Out[116]= 34.4365223939 + 0.000414467508089 I *)

f[.9]

(* Out[117]= 170.924622109 + 0.00264992865898 I *)

Offhand I cannot say much about how accurate these might be. When I did them at higher working precision I got consistent results, but that does not prove anything.

$\endgroup$

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.