Should I believe this numerical trivariate integration result involving Liouville's constant to fifteen decimal places?

Question

I have two constraints

jj1 = Q1 > 0 && Q2 > 0 && Q3 > 0 && Q1 + 3 Q2 + 2 Q3 < 1 && Q1^2 + 3 Q2 Q1 + (3 Q2 + Q3)^2 < 3 Q2 + 2 Q1 Q3

and

S1 = (256 (Q1 - Q3)^6 (1 - 9 Q2 - 6 Q3 + 
 3 (Q1^2 + 9 Q2^2 + 6 Q2 Q3 + 4 Q3^2 + Q1 (-1 + 3 Q2 + 4 Q3))))/6561 > (8192 Sqrt[2/91])/16209796869

If I issue the command

36 Integrate[Boole[jj1], {Q1, 0, 1}, {Q3, 0, 1}, {Q2, 0, 1}]

I obtain the result

(8 Pi)/(27 Sqrt[3])

I can also obtain this result using the command

36  Integrate[1, {Q1, 0, 1/3}, {Q3, 0, (1 + Q1)/4}, {Q2, 1/6 (1 - Q1 - 2 Q3) - 1/6 Sqrt[1 - 2 Q1 - 3 Q1^2 - 4 Q3 + 12 Q1 Q3], 1/6 (1 - Q1 - 2 Q3) + 1/6 Sqrt[1 - 2 Q1 - 3 Q1^2 - 4 Q3 + 12 Q1 Q3]}]

where the integration limits are obtained from

GenericCylindricalDecomposition[jj1, {Q1, Q3, Q2}]

Now, if I modify the three-dimensional integration problem immediately above by replacing Integrate by NIntegrate and replacing the integrand of 1 by Boole[S1] and use Working Precision->20, I obtain a result of

0.076083637421679183037

which WolframAlpha informs me coincides to fifteen decimal places with \begin{equation} \frac{30 (10 \mathcal{L}_{Li}+1)}{64 \mathcal{L}_{Li}+821} \approx 0.076083637421679151988, \end{equation} where $\mathcal{L}_{Li}$ is Liouville's constant https://mathworld.wolfram.com/LiouvillesConstant.html (I substituted the constant $\approx 0.11000100000000000000000100000000000000000000000000000000000000000000\ 000000000000000000000000000000000000000000000000001$ and the formula does hold.)

The source of my sceptism as to this result (which was totally unexpected/surprising to me) is that I have a third constraint

S2 = (4 Sqrt[(Q1 - Q3)^2] + 
4/3 Sqrt[
 1 - 9 Q2 - 6 Q3 + 
  3 (Q1^2 + 9 Q2^2 + 6 Q2 Q3 + 4 Q3^2 + Q1 (-1 + 3 Q2 + 4 Q3))]) >4/3

Now, if I replace Boole[S1] in the numerical integration by Boole[S2], I get 0.018903515328655050694 which to high precision is $\frac{2}{81} \left(4 \sqrt{3} \pi -21\right) \approx 0.018903515328657140917$. (This $\frac{2}{81} \left(4 \sqrt{3} \pi -21\right)$ is obtainable by exact integration, so I have full confidence in it. I, however, have no comparable exact integration result based on S1.)

Now, what I find immediately puzzling is that I am engaged in a related quasi-Monte Carlo (quasirandom) procedure (discussed in Given measures on sets and on certain Boolean combinations of the sets, can one check their consistency and/or extend them to other combinations? and How can one generate an open-ended sequence of low-discrepancy points in 3D?)

Based on more than twelve billion points, I get an S2-based estimate of 0.018903823--matching the known exact value to six places--while for the S1-based estimate, I currently have 0.07912772, agreeing to only two decimal places with the Liouville's constant-related one of 0.07608363.

Of course, I would be pleased to have an exact computation using Boole[S1]--which would definitively settle the question. Further, subject-matter (quantum-information-theoretic) considerations make an exact computation of Boole[S1||S2] of central interest.

Answer 1

Well, the coincidence between \begin{equation} \frac{30 (10 \mathcal{L}_{Li}+1)}{64 \mathcal{L}_{Li}+821} \approx 0.076083637421679151988472257219968531634399, \end{equation} where $\mathcal{L}_{Li}$ is the Lioville's constant, and 36 NIntegrate[Boole[S1], {Q1, 0, 1/3}, {Q3, 0, (1 + Q1)/4}, {Q2, 1/6 (1 - Q1 - 2 Q3) - 1/6 Sqrt[1 - 2 Q1 - 3 Q1^2 - 4 Q3 + 12 Q1 Q3], 1/6 (1 - Q1 - 2 Q3) + 1/6 Sqrt[1 - 2 Q1 - 3 Q1^2 - 4 Q3 + 12 Q1 Q3]}]

seems to extend to fifteen decimal places, but no more. With WorkingPrecision->32, the ratio of the numerical integration to the Liouville's constant-based formula is given as 1.0000000000000004080857769315461.

So, this just all seems to be a numerical oddity--with nothing particular of a fundamental nature being revealed.

Also, the disparity between the numerical integration result and the current quasirandom estimate (based on 340 million points) of 0.07912772, leads me to doubt the precision of the former.

The NIntegrate commands do yield the message "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." So, such a course should lead to greater agreement with the 0.07912772 quasirandom estimate, it might seem.

So, my basic search for the exact integration results employing Boole[S1] and Boole[S1||S2] ($\approx 0.081502827$) has yet to be completed.