Numerical Integration
- To: mathgroup at smc.vnet.net
- Subject: [mg67370] Numerical Integration
- From: Stefano Chesi <schesi at physics.purdue.edu>
- Date: Tue, 20 Jun 2006 02:15:09 -0400 (EDT)
- Sender: owner-wri-mathgroup at wolfram.com
Hi,
I have to calculate a double integral of a function (the definition is
long and I copy it
at the end of the mail). I can plot the result of the first integration:
INPUT:
FirstIntegration[kky_?NumberQ]:=
NIntegrate[int[kkx,kky],{kkx,-Sqrt[1-kky^2],Sqrt[1-kky^2]},MaxRecursion->20]
Plot[FirstIntegration[kky],{kky,-1,1}]
OUTPUT: Graphics
But the second integration doesn't work:
INPUT:
NIntegrate[FirstIntegration[kky],{kky,-1,1}]
OUTPUT: NIntegrate::inum: .... is not numerical at {kky} = {0.}
NIntegrate[FirstIntegration[kky],{kky,-1,1}]
The actual code is slightly more complicated, and I need to do the
integration
in two different steps (sometimes in the first integration I have poles
and I need
to take the principal value), in the way it is coded above.
Thank you very much,
Stefano
Here's the function definition. If I make it a simple function, e.g. by
adding " ; kky^2 "
before the last parenthesis "]", the second integration works. So, I
don't understand
what's the problem.
aa=0.1;
kx=0.2;
ky=0.1;
q=3;
int[kkx_,kky_]:=
Module[
{nn,n1,n2,n3,n4,c1,c2,c3,c4,s1,s2,s3,s4},
nn=Sqrt[(q+kx-kkx)^2+(ky-kky)^2];
n1=-Sqrt[kx^2+ky^2];
n2=Sqrt[(kx+q)^2+ky^2];
n3=-Sqrt[kkx^2+kky^2];
n4=Sqrt[(kkx-q)^2+kky^2];
c1=((kx+q)kx+ky ky)/n1/n2;
c2=(kkx(kx+q)+kky ky)/n2/n3;
c3=((kkx-q)kkx+kky kky)/n3/n4;
c4=(kx(kkx-q)+ky kky)/n4/n1;
s1=(ky kx-(kx+q) ky)/n1/n2;
s2=(kky(kx+q)-kkx ky)/n2/n3;
s3=(kky kkx-(kkx-q) kky)/n3/n4;
s4=(ky (kkx-q)-kx kky)/n4/n1;
(
(c1 c1 c3 c3/q- c1 c2 c3
c4/nn)/(q^2+kx*q-kkx*q-aa(n1+n2+n3+n4))+
(s1 s1 c3 c3/q+ s1 c2 c3 s4/nn)/(q^2+kx*q-kkx*q-
aa(-n1+n2+n3+n4))+
(s1 s1 c3 c3/q+s1 s2 c3 c4/nn)/(q^2+kx*q-kkx*q-
aa(+n1-n2+n3+n4))+
(c1 c1 s3 s3/q+ c1 s2 s3 c4/nn)/(q^2+kx*q-kkx*q-
aa(+n1+n2-n3+n4))+
(c1 c1 s3 s3/q+ c1 c2 s3 s4/nn)/(q^2+kx*q-kkx*q-
aa(+n1+n2+n3-n4))+
(c1 c1 c3 c3/q+ c1 s2 c3 s4/nn)/(q^2+kx*q-kkx*q-
aa(-n1-n2+n3+n4))+
(s1 s1 s3 s3/q- s1 s2 s3 s4/nn)/(q^2+kx*q-kkx*q-
aa(-n1+n2-n3+n4))+
(s1 s1 s3 s3/q+ s1 c2 s3 c4/nn)/(q^2+kx*q-kkx*q-
aa(-n1+n2+n3-n4))+
(s1 s1 s3 s3/q+ s1 c2 s3 c4/nn)/(q^2+kx*q-kkx*q-
aa(+n1-n2-n3+n4))+
(s1 s1 s3 s3/q- s1 s2 s3 s4/nn)/(q^2+kx*q-kkx*q-
aa(+n1-n2+n3-n4))+
(c1 c1 c3 c3/q+ c1 s2 c3 s4/nn)/(q^2+kx*q-kkx*q-
aa(+n1+n2-n3-n4))+
(s1 s1 c3 c3/q+ s1 c2 c3 s4/nn)/(q^2+kx*q-kkx*q-
aa(+n1-n2-n3-n4))+
(s1 s1 c3 c3/q+ s1 s2 c3 c4/nn)/(q^2+kx*q-kkx*q-
aa(-n1+n2-n3-n4))+
(c1 c1 s3 s3/q+ c1 s2 s3 c4/nn)/(q^2+kx*q-kkx*q-
aa(-n1-n2+n3-n4))+
(c1 c1 s3 s3/q+ c1 c2 s3 s4/nn)/(q^2+kx*q-kkx*q-
aa(-n1-n2-n3+n4))+
(c1 c1 c3 c3/q- c1 c2 c3 c4/nn)/(q^2+kx*q-kkx*q-
aa(-n1-n2-n3-n4))
)/4/Pi^2-
(4/q - 2/nn)/(q^2 + kx*q - kkx*q)/4/Pi^2
]