0
$\begingroup$

I am trying to do the following integration

$$\int_{ps_1}^{ps_2} \int_{pt_1}^{pt_2} \frac{e^{-k\,r(s,t)}}{r(s,t)}\, ds \, dt$$

Where $r(s,t)$ is the distance between points in two distinct line segment ($ps_1$ to $ps_2$ and $pt_1$ to $pt_2$), $k$ is a complex constant which $Re(k) < 0$.

I thought the easier way to go would be using a multiplying factor (to determine the point along the line) as integrand such that

$$u_s = 0 \to ps_1$$

and

$$u_s = 1 \to ps_2$$

The implemented integration being

$$\int_{0}^{1} \int_{0}^{1} \frac{e^{-k\,r(u_s,u_t)}}{r(u_s,u_t)}\, du_s \, du_t$$

A colleague of mine working on the same problem got different results from what I got. Comparing our implementations, I say they are mathematically identical...

Function definitions

integral[snd_, rcv_, k_] := Block[
   {ei, ef, ri, rf, Ls, Lr, int, r, pe, pr},
   ei = snd[[1]];
   ef = snd[[2]];
   ri = rcv[[1]];
   rf = rcv[[2]];
   Ls = Norm[ef - ei];
   Lr = Norm[rf - ri];
   int = NIntegrate[
     pe = (ef - ei) ts/Ls + ei;
     pr = (rf - ri) tr/Lr + ri;
     r = Norm[pe - pr];
     Exp[-k r]/r,
     {tr, 0., Lr}, {ts, 0., Ls},
     WorkingPrecision -> 10, MaxRecursion -> 200
     ];
   Return[int]
   ];

integral2[snd_, rcv_, k_] := Block[
   {ei, ef, ri, rf, Ls, Lr, int, r, pe, pr},
   ei = snd[[1]];
   ef = snd[[2]];
   ri = rcv[[1]];
   rf = rcv[[2]];
   Ls = Norm[ef - ei];
   Lr = Norm[rf - ri];
   int = NIntegrate[
     pe = (ef - ei) ts + ei;
     pr = (rf - ri) tr + ri;
     r = Norm[pe - pr];
     Exp[-k r]/r,
     {tr, 0., 1.}, {ts, 0., 1.},
     WorkingPrecision -> 10, MaxRecursion -> 200
     ];
   Return[int]
   ];

Testing

h = 0.3;
L = 1.0;
snd = {{0., 0., h}, {L, 0., h}};
rcv = {{0., L, -h}, {L, 0., -h}};
integral[snd, rcv, -1]
integral2[snd, rcv, -1]

Results in

3.949134609
2.792459862

Where is the mistake, is it implementation-wise or is the assumption I've made (mathematical identity) wrong? Or is this something related to numerical conditioning?

Any help is apprecciated.

$\endgroup$
2
$\begingroup$

By the transformation formula for integrals, the second integral should read

Ls Lr NIntegrate[
  pe = (ef - ei) ts + ei;
  pr = (rf - ri) tr + ri;
  r = Norm[pe - pr];
  Exp[-k r]/r, 
  {tr, 0., 1.}, {ts, 0., 1.}]
$\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.