0
$\begingroup$

Consider a data file distr.dat. It is a tabulated grid in the form x1,x2, func[x1,x2]. I interpolate it in two different ways:

distrdata = 
  Import[FileNameJoin[{NotebookDirectory[], "distr.dat"}], "Table"];
x1minmax = MinMax[distrdata[[All, 1]]]
x2minmax = MinMax[distrdata[[All, 2]]]
int1[x1_, x2_] = 
  10^(Interpolation[{Log10[#[[1]]], Log10[#[[2]]], 
         Log10[#[[3]] + 10^-90]} & /@ distrdata, 
      InterpolationOrder -> 1][Log10[x1], Log10[x2]]);
int2[x1_, x2_] = 
  Interpolation[{#[[1]], #[[2]], #[[3]]} & /@ distrdata, 
    InterpolationOrder -> 1][x1, x2];

The interpolations are expected to differ in the domains where the function changes quickly - for large values of x1, x2.

Next, I want to compare the integrals

NIntegrate[
 int1[x1, x2], {x1, x1minmax[[1]], x1minmax[[2]]}, {x2, x2minmax[[1]],
   x2minmax[[2]]}, Method -> "InterpolationPointsSubdivision"]
NIntegrate[
 int2[x1, x2], {x1, x1minmax[[1]], x1minmax[[2]]}, {x2, x2minmax[[1]],
   x2minmax[[2]]}, Method -> "InterpolationPointsSubdivision"]

I get two completely different results:

0.838699

7.06604*10^-6

Next, if I take into account that the int1, int2 quickly drop for x1 > 0.2, I get

NIntegrate[
 int1[x1, x2], {x1, x1minmax[[1]], 0.2}, {x2, x2minmax[[1]],
   x2minmax[[2]]}, Method -> "InterpolationPointsSubdivision"]
NIntegrate[
 int2[x1, x2], {x1, x1minmax[[1]], 0.2}, {x2, x2minmax[[1]],
   x2minmax[[2]]}, Method -> "InterpolationPointsSubdivision"]

0.83866

1.04657

Questions:

  1. If integrating int2 over all domain of the definition of x1 (i.e., from x1minmax[[1]] to x1minmax[[2]]), what is the reason for the method "InterpolationPointsSubdivision" to behave like "MonteCarlo", i.e., return small results if the integrand is non-zero only within some narrow domain?

  2. Can the difference between the values of the integrals when integrating over x1 up to 0.2 be attributed solely to the way of the interpolation, i.e., adjusting the value of the functions between the data points?

To me, the answer to question 2 is not obvious. Indeed, let us integrate over x2 only and look at the x1 integrand:

tabb = ParallelTable[{10^x1, 
    Quiet[NIntegrate[
      int1[10^x1, x2], {x2, x2minmax[[1]], x2minmax[[2]]}, 
      Method -> "InterpolationPointsSubdivision"]], 
    Quiet[NIntegrate[
      int2[10^x1, x2], {x2, x2minmax[[1]], x2minmax[[2]]}, 
      Method -> "InterpolationPointsSubdivision"]]}, {x1, 
    Log10[x1minmax[[1]]], Log10[0.2], 0.1}];
ListLogLogPlot[{{#[[1]], #[[2]]/#[[3]]} & /@ tabb}, Joined -> True, 
 Frame -> True, ImageSize -> Large]

The ratio of the integrals of int1, int2 departs from 1 everywhere and does not jump except for the domain where the function quickly drops, which I find strange.

enter image description here

$\endgroup$

1 Answer 1

1
$\begingroup$

I doubt very much about the sense of your attempt "interpolate a distribution".

With a little help NIntegrate evaluates similar results. You only have to extend the integration range.

Therfore I assume that integrand takes significant values in the range x1<.04and 20<x2<100

NIntegrate[int1[x1, x2], {x1, x1minmax[[1]], 0.04, x1minmax[[2]]}
,{x2,x2minmax[[1]], 20, 100, x2minmax[[2]]} , PrecisionGoal -> 3, AccuracyGoal -> 6, Method -> "InterpolationPointsSubdivision"  ]
(*0.838117*)

NIntegrate[int2[x1, x2], {x1, x1minmax[[1]], .04, x1minmax[[2]]}
, {x2,x2minmax[[1]], 20, 100, x2minmax[[2]]} , PrecisionGoal -> 3,AccuracyGoal -> 6, Method -> "InterpolationPointsSubdivision" ]
(*1.04585*)
$\endgroup$
2
  • $\begingroup$ Thanks! But what do you mean by having doubts? $\endgroup$ Commented Mar 28, 2023 at 6:37
  • 1
    $\begingroup$ A distribution isn't plotable nor interpolable, though let's call your function narrow banded or something like that.. $\endgroup$ Commented Mar 28, 2023 at 6:42

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.