Group Abstract

Message Boards

WOLFRAM COMMUNITY

13K Views

13 Replies

0 Total Likes

View groups...

Follow this post

Share this post:

GROUPS:

Mathematics Physics Calculus Wolfram Language Numerical Computation Symbolic Computations

Integration involving vectors?

Arny Toynbee

Posted 9 years ago

I am trying to integrate the following expression - where ${\bf p3}$, ${\bf q3}$, ${\bf p}$, ${\bf q}$ are vectors. The angle between each pair of vectors e.g. ${\bf q3}$, ${\bf p3}$ can't be ignored. $\mu$ is a scalar. $\int d{\bf p3} \int d{\bf q3} \int d{\bf q} \frac{1}{\sqrt{\mu ^2+\text{q3}^2} \sqrt{\mu ^2+({\bf p}-{\bf p3}-{\bf q3})^2} \sqrt{\mu ^2+(-{\bf p3}+{\bf q}-{\bf q3})^2} \left(\sqrt{\mu ^2+(-{\bf p3}+{\bf q}-{\bf q3})^2}+\sqrt{\mu ^2+{\bf q3}^2}\right)}$. Tried the following: q3v = Table[Subscript[q3, i], {i, 3}]; p3v = Table[Subscript[p3, i], {i, 3}]; qv = Table[Subscript[q, i], {i, 3}]; pv = Table[Subscript[p, i], {i, 3}]; and then a := 1/Sqrt[q3v^2 + \[Mu]^2]; b := 1/(Sqrt[q3v^2 + \[Mu]^2] + Sqrt[(-p3v + qv - q3v)^2 + \[Mu]^2]); c := 1/(Sqrt[(pv - p3v - q3v)^2 + \[Mu]^2] Sqrt[(-p3v + qv - q3v)^2 + \[Mu]^2]) The output of a.b.c gives {1/(Sqrt[\[Mu]^2 + (Subscript[p, 1] - Subscript[p3, 1] - Subscript[q3, 1])^2]Sqrt[\[Mu]^2 + (-Subscript[p3, 1] + Subscript[q, 1] - Subscript[q3, 1])^2]Sqrt[\[Mu]^2 + Subscript[q3, 1]^2]* (Sqrt[\[Mu]^2 + (-Subscript[p3, 1] + Subscript[q, 1] - Subscript[q3, 1])^2] + Sqrt[\[Mu]^2 + Subscript[q3, 1]^2])), 1/(Sqrt[\[Mu]^2 + (Subscript[p, 2] - Subscript[p3, 2] - Subscript[q3, 2])^2]Sqrt[\[Mu]^2 + (-Subscript[p3, 2] + Subscript[q, 2] - Subscript[q3, 2])^2]Sqrt[\[Mu]^2 + Subscript[q3, 2]^2]* (Sqrt[\[Mu]^2 + (-Subscript[p3, 2] + Subscript[q, 2] - Subscript[q3, 2])^2] + Sqrt[\[Mu]^2 + Subscript[q3, 2]^2])), 1/(Sqrt[\[Mu]^2 + (Subscript[p, 3] - Subscript[p3, 3] - Subscript[q3, 3])^2]Sqrt[\[Mu]^2 + (-Subscript[p3, 3] + Subscript[q, 3] - Subscript[q3, 3])^2]Sqrt[\[Mu]^2 + Subscript[q3, 3]^2]* (Sqrt[\[Mu]^2 + (-Subscript[p3, 3] + Subscript[q, 3] - Subscript[q3, 3])^2] + Sqrt[\[Mu]^2 + Subscript[q3, 3]^2]))} Now, let's try to evaluate the first of the triplet in the output of a.b.c Integrate[1/(Sqrt[\[Mu]^2 + (Subscript[p, 1] - Subscript[p3, 1] - Subscript[q3, 1])^2]Sqrt[\[Mu]^2 + (-Subscript[p3, 1] + Subscript[q, 1] - Subscript[q3, 1])^2]Sqrt[\[Mu]^2 + Subscript[q3, 1]^2]* (Sqrt[\[Mu]^2 + (-Subscript[p3, 1] + Subscript[q, 1] - Subscript[q3, 1])^2] + Sqrt[\[Mu]^2 + Subscript[q3, 1]^2])), {q3, -1, 1}] which Mathematica 11 will not evaluate. Is this approach correct?

I am trying to integrate the following expression - where ${\bf p3}$, ${\bf q3}$, ${\bf p}$, ${\bf q}$ are vectors. The angle between each pair of vectors e.g. ${\bf q3}$, ${\bf p3}$ can't be ignored. $\mu$ is a scalar.

$\int d{\bf p3} \int d{\bf q3} \int d{\bf q} \frac{1}{\sqrt{\mu ^2+\text{q3}^2} \sqrt{\mu ^2+({\bf p}-{\bf p3}-{\bf q3})^2} \sqrt{\mu ^2+(-{\bf p3}+{\bf q}-{\bf q3})^2} \left(\sqrt{\mu ^2+(-{\bf p3}+{\bf q}-{\bf q3})^2}+\sqrt{\mu ^2+{\bf q3}^2}\right)}$.

Tried the following:

q3v = Table[Subscript[q3, i], {i, 3}]; 
p3v = Table[Subscript[p3, i], {i, 3}];
qv = Table[Subscript[q, i], {i, 3}];
 pv = Table[Subscript[p, i], {i, 3}];

and then

a := 1/Sqrt[q3v^2 + \[Mu]^2]; 
b := 1/(Sqrt[q3v^2 + \[Mu]^2] + Sqrt[(-p3v + qv - q3v)^2 + \[Mu]^2]); 
c := 1/(Sqrt[(pv - p3v - q3v)^2 + \[Mu]^2] Sqrt[(-p3v + qv - 
        q3v)^2 + \[Mu]^2])

The output of a.b.c gives

{1/(Sqrt[\[Mu]^2 + (Subscript[p, 1] - Subscript[p3, 1] - Subscript[q3, 1])^2]*Sqrt[\[Mu]^2 + (-Subscript[p3, 1] + Subscript[q, 1] - Subscript[q3, 1])^2]*Sqrt[\[Mu]^2 + Subscript[q3, 1]^2]*
    (Sqrt[\[Mu]^2 + (-Subscript[p3, 1] + Subscript[q, 1] - Subscript[q3, 1])^2] + Sqrt[\[Mu]^2 + Subscript[q3, 1]^2])), 
  1/(Sqrt[\[Mu]^2 + (Subscript[p, 2] - Subscript[p3, 2] - Subscript[q3, 2])^2]*Sqrt[\[Mu]^2 + (-Subscript[p3, 2] + Subscript[q, 2] - Subscript[q3, 2])^2]*Sqrt[\[Mu]^2 + Subscript[q3, 2]^2]*
    (Sqrt[\[Mu]^2 + (-Subscript[p3, 2] + Subscript[q, 2] - Subscript[q3, 2])^2] + Sqrt[\[Mu]^2 + Subscript[q3, 2]^2])), 
  1/(Sqrt[\[Mu]^2 + (Subscript[p, 3] - Subscript[p3, 3] - Subscript[q3, 3])^2]*Sqrt[\[Mu]^2 + (-Subscript[p3, 3] + Subscript[q, 3] - Subscript[q3, 3])^2]*Sqrt[\[Mu]^2 + Subscript[q3, 3]^2]*
    (Sqrt[\[Mu]^2 + (-Subscript[p3, 3] + Subscript[q, 3] - Subscript[q3, 3])^2] + Sqrt[\[Mu]^2 + Subscript[q3, 3]^2]))}

Now, let's try to evaluate the first of the triplet in the output of a.b.c

Integrate[1/(Sqrt[\[Mu]^2 + (Subscript[p, 1] - Subscript[p3, 1] - Subscript[q3, 1])^2]*Sqrt[\[Mu]^2 + (-Subscript[p3, 1] + Subscript[q, 1] - Subscript[q3, 1])^2]*Sqrt[\[Mu]^2 + Subscript[q3, 1]^2]*
    (Sqrt[\[Mu]^2 + (-Subscript[p3, 1] + Subscript[q, 1] - Subscript[q3, 1])^2] + Sqrt[\[Mu]^2 + Subscript[q3, 1]^2])), {q3, -1, 1}]

which Mathematica 11 will not evaluate.

Is this approach correct?

POSTED BY: Arny Toynbee

13 Replies

Sort By:

Arny Toynbee

Posted 9 years ago

Thank you, Hans.

POSTED BY: Arny Toynbee

Hans Dolhaine

Hans Dolhaine, retired

Posted 9 years ago

jjm = jj /. m -> 2. implies that m is being set to 2. prior to performing the integration. Yes. To do the numerical integration m must have a definite (numerical) value. is it over all of p3v, q3v, pv, qv? Yes. All p3v, q3v, pv, qv? Look at intvariables2 to see the boundaries. See attached notebook. Attachments:

POSTED BY: Hans Dolhaine

Arny Toynbee

Posted 9 years ago

Thank you. Sorry for the following really dumb question - jjm = jj /. m -> 2. implies that m is being set to 2. prior to performing the integration. What is the @@ intvariables2 in this case, is it over all of p3v, q3v, pv, qv?

POSTED BY: Arny Toynbee

Hans Dolhaine

Hans Dolhaine, retired

Posted 9 years ago

At last I noticed that ther is a typo in the definition of the integrand. It should read (at least I think so) a := 1/Sqrt[m^2 + q3v.q3v] b := 1/Sqrt[m^2 + (p3v + q3v - pv).(p3v + q3v - pv)] c := 1/Sqrt[m^2 + (p3v + q3v - qv).(p3v + q3v - qv)] d := 1/(Sqrt[m^2 + (p3v + q3v - qv).(p3v + q3v - qv)] + Sqrt[m^2 + q3v.q3v]) With this In[9]:= jj = a b c d // FullSimplify Out[9]= 1/(Sqrt[ m^2 + q3[1]^2 + q3[2]^2 + q3[3]^2] \[Sqrt](m^2 + (-p[1] + p3[1] + q3[1])^2 + (-p[2] + p3[2] + q3[2])^2 + (-p[3] + p3[3] + q3[3])^2) Sqrt[m^2 + (p3[1] - q[1] + q3[1])^2 + (p3[2] - q[2] + q3[2])^2+ (p3[3] - q[3] + q3[3])^2] (Sqrt[ m^2 + q3[1]^2 + q3[2]^2 + q3[3]^2] + \[Sqrt](m^2 + (p3[1] - q[1] + q3[1])^2 + (p3[2] - q[2] + q3[2])^2 + (p3[3] - q[3] + q3[3])^2))) and the denominator should never vanish. Then jjm = jj /. m -> 2. and In[15]:= NIntegrate[jjm, Evaluate[Sequence @@ intvariables2], MaxPoints -> 100] During evaluation of In[15]:= NIntegrate::maxp: The integral failed to converge after 13227 integrand evaluations. NIntegrate obtained 60.37299593570131` and 8.705681178065333` for the integral and error estimates. >> Out[15]= 60.373 The Integration is reasonably fast. Without the MaxPoint-Option it takes quite a time, so I aborted it. But it could be an idea to let it run and compare the result with that with the MaxPoint-Option.

At last I noticed that ther is a typo in the definition of the integrand. It should read (at least I think so)

a := 1/Sqrt[m^2 + q3v.q3v]
b := 1/Sqrt[m^2 + (p3v + q3v - pv).(p3v + q3v - pv)]
c := 1/Sqrt[m^2 + (p3v + q3v - qv).(p3v + q3v - qv)]
d := 1/(Sqrt[m^2 + (p3v + q3v - qv).(p3v + q3v - qv)] + Sqrt[m^2 + q3v.q3v])

With this

In[9]:= jj = a b c d // FullSimplify

Out[9]= 1/(Sqrt[
   m^2 + q3[1]^2 + q3[2]^2 +  q3[3]^2] \[Sqrt](m^2 + (-p[1] + p3[1] + q3[1])^2 + 
  (-p[2] + p3[2] + q3[2])^2 + (-p[3] + p3[3] + q3[3])^2) Sqrt[m^2 + (p3[1] - q[1] + q3[1])^2 + (p3[2] - q[2] + q3[2])^2+
   (p3[3] - q[3] + q3[3])^2] (Sqrt[ m^2 + q3[1]^2 + q3[2]^2 + q3[3]^2] + \[Sqrt](m^2 + (p3[1] - q[1] + q3[1])^2 +
   (p3[2] - q[2] + q3[2])^2 + (p3[3] - q[3] + q3[3])^2)))

and the denominator should never vanish.

Then

jjm = jj /. m -> 2.

and

In[15]:= NIntegrate[jjm, Evaluate[Sequence @@ intvariables2], 
 MaxPoints -> 100]

During evaluation of In[15]:= NIntegrate::maxp: The integral failed to converge after 13227 integrand evaluations. NIntegrate obtained 60.37299593570131` and 8.705681178065333` for the integral and error estimates. >>

Out[15]= 60.373

The Integration is reasonably fast. Without the MaxPoint-Option it takes quite a time, so I aborted it. But it could be an idea to let it run and compare the result with that with the MaxPoint-Option.

POSTED BY: Hans Dolhaine

Hans Dolhaine

Hans Dolhaine, retired

Posted 9 years ago

At last I noticed that ther is a typo in the definition of the integrand. It should read (at least I think so) a := 1/Sqrt[m^2 + q3v.q3v] b := 1/Sqrt[m^2 + (p3v + q3v - pv).(p3v + q3v - pv)] c := 1/Sqrt[m^2 + (p3v + q3v - qv).(p3v + q3v - qv)] d := 1/(Sqrt[m^2 + (p3v + q3v - qv).(p3v + q3v - qv)] + Sqrt[m^2 + q3v.q3v]) With this In[9]:= jj = a b c d // FullSimplify Out[9]= 1/(Sqrt[ m^2 + q3[1]^2 + q3[2]^2 + q3[3]^2] \[Sqrt](m^2 + (-p[1] + p3[1] + q3[1])^2 + (-p[2] + p3[2] + q3[2])^2 + (-p[3] + p3[3] + q3[3])^2) Sqrt[m^2 + (p3[1] - q[1] + q3[1])^2 + (p3[2] - q[2] + q3[2])^2+ (p3[3] - q[3] + q3[3])^2] (Sqrt[ m^2 + q3[1]^2 + q3[2]^2 + q3[3]^2] + \[Sqrt](m^2 + (p3[1] - q[1] + q3[1])^2 + (p3[2] - q[2] + q3[2])^2 + (p3[3] - q[3] + q3[3])^2))) and the denominator should never vanish. Then jjm = jj /. m -> 2. and In[15]:= NIntegrate[jjm, Evaluate[Sequence @@ intvariables2], MaxPoints -> 100] During evaluation of In[15]:= NIntegrate::maxp: The integral failed to converge after 13227 integrand evaluations. NIntegrate obtained 60.37299593570131` and 8.705681178065333` for the integral and error estimates. >> Out[15]= 60.373 The Integration is reasonably fast. Without the MaxPoint-Option it takes quite a time, so I aborted it. But it could be an idea to let it run and compare the result with that with the MaxPoint-Option.

At last I noticed that ther is a typo in the definition of the integrand. It should read (at least I think so)

a := 1/Sqrt[m^2 + q3v.q3v]
b := 1/Sqrt[m^2 + (p3v + q3v - pv).(p3v + q3v - pv)]
c := 1/Sqrt[m^2 + (p3v + q3v - qv).(p3v + q3v - qv)]
d := 1/(Sqrt[m^2 + (p3v + q3v - qv).(p3v + q3v - qv)] + Sqrt[m^2 + q3v.q3v])

With this

In[9]:= jj = a b c d // FullSimplify

Out[9]= 1/(Sqrt[
   m^2 + q3[1]^2 + q3[2]^2 +  q3[3]^2] \[Sqrt](m^2 + (-p[1] + p3[1] + q3[1])^2 + 
  (-p[2] + p3[2] + q3[2])^2 + (-p[3] + p3[3] + q3[3])^2) Sqrt[m^2 + (p3[1] - q[1] + q3[1])^2 + (p3[2] - q[2] + q3[2])^2+
   (p3[3] - q[3] + q3[3])^2] (Sqrt[ m^2 + q3[1]^2 + q3[2]^2 + q3[3]^2] + \[Sqrt](m^2 + (p3[1] - q[1] + q3[1])^2 +
   (p3[2] - q[2] + q3[2])^2 + (p3[3] - q[3] + q3[3])^2)))

and the denominator should never vanish.

Then

jjm = jj /. m -> 2.

and

In[15]:= NIntegrate[jjm, Evaluate[Sequence @@ intvariables2], 
 MaxPoints -> 100]

During evaluation of In[15]:= NIntegrate::maxp: The integral failed to converge after 13227 integrand evaluations. NIntegrate obtained 60.37299593570131` and 8.705681178065333` for the integral and error estimates. >>

Out[15]= 60.373

POSTED BY: Hans Dolhaine

Hans Dolhaine

Hans Dolhaine, retired

Posted 9 years ago

I don't know. I got ( integrating for all variables from -1 to 1 ) NIntegrate::eincr: The global error of the strategy GlobalAdaptive has increased more than 2000 times. The global error is expected to decrease monotonically after a number of integrand evaluations. Suspect one of the following: the working precision is insufficient for the specified precision goal; the integrand is highly oscillatory or it is not a (piecewise) smooth function; or the true value of the integral is 0. Increasing the value of the GlobalAdaptive option MaxErrorIncreases might lead to a convergent numerical integration. NIntegrate obtained 70.76236375923493` and 0.010212021932754334` for the integral and error estimates. >> With the above code I did (# -> 0 & /@ intvariables) jjm /. Drop[(# -> 0 & /@ intvariables), {1}] jjm /. Drop[(# -> 0 & /@ intvariables), {5}] jjm /. Drop[(# -> 0 & /@ intvariables), {12}] and Plot[jjm /. Drop[(# -> 0 & /@ intvariables), {1}], {q3[1], -1, 1}] Plot[jjm /. Drop[(# -> 0 & /@ intvariables), {5}], {p3[2], -1, 1}] Plot[jjm /. Drop[(# -> 0 & /@ intvariables), {12}], {p[3], -1, 1}] to plot the integrand for some examples with all but one variable = zero. Seems to show a reasonable behavior. then I tried NIntegrate[jj /. m -> N[2., 30], Evaluate[Sequence @@ intvariables2], WorkingPrecision -> 30] but this is still running (I didn't notice but I think for more than an hour)

I don't know.

I got ( integrating for all variables from -1 to 1 )

NIntegrate::eincr: The global error of the strategy GlobalAdaptive has increased more than 2000 times. The global error is expected to decrease monotonically after a number of integrand evaluations. Suspect one of the following: the working precision is insufficient for the specified precision goal; the integrand is highly oscillatory or it is not a (piecewise) smooth function; or the true value of the integral is 0. Increasing the value of the GlobalAdaptive option MaxErrorIncreases might lead to a convergent numerical integration. NIntegrate obtained 70.76236375923493` and 0.010212021932754334` for the integral and error estimates. >>

With the above code I did

(# -> 0 & /@ intvariables)
jjm /. Drop[(# -> 0 & /@ intvariables), {1}]
jjm /. Drop[(# -> 0 & /@ intvariables), {5}]
jjm /. Drop[(# -> 0 & /@ intvariables), {12}]

and

Plot[jjm /. Drop[(# -> 0 & /@ intvariables), {1}], {q3[1], -1, 1}]
Plot[jjm /. Drop[(# -> 0 & /@ intvariables), {5}], {p3[2], -1, 1}]
Plot[jjm /. Drop[(# -> 0 & /@ intvariables), {12}], {p[3], -1, 1}]

to plot the integrand for some examples with all but one variable = zero. Seems to show a reasonable behavior.

then I tried

NIntegrate[jj /. m -> N[2., 30], Evaluate[Sequence @@ intvariables2], WorkingPrecision -> 30]

but this is still running (I didn't notice but I think for more than an hour)

POSTED BY: Hans Dolhaine

Arny Toynbee

Posted 9 years ago

Thank you for the corrections, and the caution on using subscripts. I have suffered for a long time with subscripts, without realizing this potential bug. NIntegrate between -1, 1 produces the following output:\ NIntegrate::slwcon: Numerical integration converging too slowly; suspect one of the following: singularity, value of the integration is 0, highly oscillatory integrand, or WorkingPrecision too small. \ and after 10+ minutes\ NIntegrate::eincr: The global error of the strategy GlobalAdaptive has increased more than 2000 times. The global error is expected to decrease monotonically after a number of integrand evaluations. Suspect one of the following: the working precision is insufficient for the specified precision goal; the integrand is highly oscillatory or it is not a (piecewise) smooth function; or the true value of the integral is 0. Increasing the value of the GlobalAdaptive option MaxErrorIncreases might lead to a convergent numerical integration. NIntegrate obtained 441.218 -211.434 I and 1135.6427949097094` for the integral and error estimates.\ 441.218293633602 - 211.43356703452153*I Is there a way in Mathematica to analyze the integrand behaviour in a way that can catch any ill-behaved regions? Thanks

POSTED BY: Arny Toynbee

Hans Dolhaine

Hans Dolhaine, retired

Posted 9 years ago

POSTED BY: Hans Dolhaine

Frank Kampas

Frank Kampas, Physicist at Large Consulting

Posted 9 years ago

suggest you try a simpler version of the problem which will be easier to troubleshoot

POSTED BY: Frank Kampas

Arny Toynbee

Posted 9 years ago

NIntegrate doesn't seem to work either.

POSTED BY: Arny Toynbee

Frank Kampas

Frank Kampas, Physicist at Large Consulting

Posted 9 years ago

You could do numerical integration for various values of lambda.

POSTED BY: Frank Kampas

Arny Toynbee

Posted 9 years ago

Maple 17 appears to give an output when integrated wrt q3 for example, but the answer is very (several pages print) long. Is there a way you would suggest to integrate numerically, say between 0 and \lambda ? Thanks.

POSTED BY: Arny Toynbee

Frank Kampas

Frank Kampas, Physicist at Large Consulting

Posted 9 years ago

Assuming your code is correct, it's probably too complicated to be integrated symbolically.

POSTED BY: Frank Kampas

Reply to this discussion

Reply Preview

Attachments

Remove Add a file to this post

Follow this discussion

or Discard

Feedback