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:
If integrating
int2
over all domain of the definition ofx1
(i.e., fromx1minmax[[1]]
tox1minmax[[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?Can the difference between the values of the integrals when integrating over
x1
up to0.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.