3
$\begingroup$

I'm recently studying some topics in analytic number theory and I have encountered results involving the infinite product $$C=\prod_{p}\left(1-\frac{1}{p(p+1)}\right)$$ where $p$ denotes calculating the product over all prime numbers.

Or more intuitively,$$C=\left(1-\frac{1}{2\cdot3}\right)\cdot\left(1-\frac{1}{3\cdot4}\right)\cdot\left(1-\frac{1}{5\cdot6}\right)\cdot\left(1-\frac{1}{7\cdot8}\right)\cdot\left(1-\frac{1}{11\cdot12}\right)\cdots$$ Mathematica can show that $C\approx0.704442.$

However, I need a high precision result, roughly accurate to at least 20 digits, so that I could use the Inverse Symbolic Calculator to make an educated guess of any possible analytic expression for $C$ in forms like $\dfrac{\zeta(2)\zeta(3)}{\zeta(6)}$.

When I use NProduct, I get

NProduct[1 - 1/(Prime[k] (Prime[k] + 1)), {k, 1, Infinity},
         AccuracyGoal -> 20, PrecisionGoal -> 20]

During evaluation of In[14]:= Prime::intpp: Positive integer argument expected in Prime[15.]. >>
During evaluation of In[14]:= Prime::intpp: Positive integer argument expected in Prime[14.]. >>
During evaluation of In[14]:= Prime::intpp: Positive integer argument expected in Prime[13.]. >>
During evaluation of In[14]:= General::stop: Further output of Prime::intpp will be suppressed during this calculation. >>
0.705024

However, when I use N[Product[]], the calculation was so slow that I had to abort the calculation after minutes of calculating.

N[Product[1 - 1/(Prime[k] (Prime[k] + 1)), {k, 1, 100000}], 50]

0.70444223707595873775750824971600704569326374380066

N[Product[1 - 1/(Prime[k] (Prime[k] + 1)), {k, 1, 1000000}], 50]

0.70444220359808164073982578152558851839613830115938

N[Product[1 - 1/(Prime[k] (Prime[k] + 1)), {k, 1, 10000000}], 50]

$Aborted

Is there any way to calculate products like $C$ to high precision in Mathematica? Thanks.

|improve this question
$\endgroup$
  • 1
    $\begingroup$ Page 11 of arxiv.org/pdf/0903.2514v2.pdf gives the expansion, as $Q^{(1)}_1$. It's the "carefree constant", according to Wikipedia. (Thanks to [OEIS] for this.) [OEIS]: oeis.org/A065463 $\endgroup$ – Patrick Stevens Aug 15 '15 at 12:22
  • 1
    $\begingroup$ Did not have any luck with the inverse symbolic calculator though. $\endgroup$ – bobbym Aug 15 '15 at 15:34
  • $\begingroup$ Thanks, @PatrickStevens! @bbgodfrey, $\endgroup$ – Zhenhua Liu Aug 16 '15 at 0:32
  • $\begingroup$ @bbgodfrey, I'll try my best to help, though I'm still a novice in using Mathematica. $\endgroup$ – Zhenhua Liu Aug 16 '15 at 0:41
  • $\begingroup$ @bobbym, Thanks. It seems no analytic expression is available. $\endgroup$ – Zhenhua Liu Aug 16 '15 at 0:41
5
$\begingroup$

Using the formula given in the arXiv preprint Patrick linked to for the "carefree constant" gives:

Exp[NSum[(-1)^k PrimeZetaP[k] (1 - LucasL[k])/k, {k, 2, ∞}, Compiled -> False,
         Method -> "AlternatingSigns", NSumTerms -> 20, WorkingPrecision -> 30]]
   0.704442200999165592738713909247

Note that this agrees with the result in the OEIS up to twenty digits. The only speed-limiting part of this is the calculation of the prime zeta function.

|improve this answer
$\endgroup$
  • $\begingroup$ Thanks a lot! That's great help to me. $\endgroup$ – Zhenhua Liu Aug 16 '15 at 0:42

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.