Apéry was a French mathematician
https://www.wikiwand.com/en/Roger...
For a pretty robust list of mathematical constants, see this Wikipe...
The Golden ratio is an important number in fields like art, archite...
Because of the importance of this problem, the Clay Institute a mil...
The Basel problem is asking to calculate the value of the Riemann Z...
A rational number is a number than can be written as a quotient of ...
The Gompertz constant (or Euler-Gompertz constant) is
$\delta=\fra...
The Catalan constant $G$ is defined as $G=\frac{1}{1^2}-\frac{1}{3^...
The Euler-Mascheroni constant $\gamma$ is
$\gamma = \int_1 ^\inft...
The irrationality measure, also called 'Liouville number' is given ...
Continued fractions have a rich mathematical history, with ancient ...
Polynomial continued fractions have a remarkable property of encomp...
The (standard) Bessel functions, denoted $J_\alpha, Y_\alpha$, are ...
The gamma function $\Gamma(x)$ generalizes the factorial operation ...
The hypergeometric function $_2 F_1$ generalizes many functions, su...
The Ramanujan Machine project, aside from being a tribute to Ramanu...
"Experimental mathematics is that branch of mathematics that concer...
Traditionally, mathematics is approached from a 'top-down' perspect...
In the original Ramanujan Machine paper, we aspired to have a futur...
This was very unexpected to the team who worked on the Ramanujan Ma...
The plot presented here takes inspiration from a real-world analog,...
BOINC is a platform that allows volunteers' computers to contribute...
The PSLQ algorithm was selected as one of the top 10 algorithms of ...
The Lerch zeta function is a generalization of the Riemann Zeta fun...
A conservative matrix field can also be viewed as a method for acce...
The process performed here is analogous to defining a potential on ...
At first glance, this condition may seem unrelated to fundamental m...
The analogy between conservative matrix fields and conservative vec...
This observation is truly remarkable. We are essentially working wi...
The [Liouville–Roth](https://en.wikipedia.org/wiki/Liouville_number...
Factorial reduction can play a crucial part in achieving a substant...
Even if we were to discover a continued fraction within the ZigZagZ...
These are good problems to have. If anything - this is a mark of th...
Since practically all interesting constants and closed-form functio...
Numerical validation challenges are in principle "just a technical ...
This information-theoretical approach in a sense attempts to 'quant...
The quality of PSLQ approximation is related (among other things) t...
If true, this would mean that each constant has a basic property - ...
It also provides a natural tool to measure the complexity and irrat...
Personally, the vision of man and AI, working together - each contr...
Algorithm-assisted discovery of an intrinsic order among
mathematical constants
Rotem Elimelech, Ofir David, Carlos De la Cruz Mengual, Rotem Kalisch, Wolfgang
Berndt, Michael Shalyt, Mark Silberstein, Yaron Hadad, and Ido Kaminer
Technion - Israel Institute of Technology, Haifa 3200003, Israel
August 24, 2023
Abstract
In recent decades, a growing number of discoveries in fields of mathematics have been assisted by computer
algorithms, primarily for exploring large parameter spaces that humans would take too long to investigate
[1–7]. As computers and algorithms become more powerful, an intriguing possibility arises—the interplay
between human intuition and computer algorithms can lead to discoveries of novel mathematical concepts
that would otherwise remain elusive [8–13]. To realize this perspective, we have developed a massively
parallel computer algorithm that discovers an unprecedented number of continued fraction formulas for
fundamental mathematical constants. The sheer number of formulas discovered by the algorithm unveils a
novel mathematical structure that we call the conservative matrix field. Such matrix fields (1) unify thousands
of existing formulas, (2) generate infinitely many new formulas, and most importantly, (3) lead to unexpected
relations between different mathematical constants, including multiple integer values of the Riemann zeta
function. Conservative matrix fields also enable new mathematical proofs of irrationality. In particular,
we can use them to generalize the celebrated proof by Apéry for the irrationality of ζ(3) [14, 15]. Utilizing
thousands of personal computers worldwide [16], our computer-supported research strategy demonstrates
the power of experimental mathematics, highlighting the prospects of large-scale computational approaches
to tackle longstanding open problems and discover unexpected connections across diverse fields of science.
1 Mathematical constants
Mathematical constants emerge naturally in various fields of science and branches of mathematics, such
as geometry, combinatorics, number theory, and probability theory [17]. Well-known examples of constants
are the ratio between the circumference and the diameter of a circle π, the base of the natural logarithm e,
the golden ratio φ, and the values of the infinite sum for integers n ≥ 2:
ζ(n)
:
=
1
1
n
+
1
2
n
+
1
3
n
+ ···
The latter are the integer values of the Riemann zeta function ζ, a cornerstone of mathematical exploration
in the field of number theory, due to the remarkable connection between its complex zeros and the asymp-
totic distribution of prime numbers. The precise distribution of these zeros of ζ constitutes the Riemann
hypothesis, arguably the most important unsolved question in pure mathematics to this day.
The occurrence of a constant in diverse mathematical contexts reflects its significance and often leads to
the discovery of fundamental connections. A remarkable example is given by the so-called Basel problem,
posed in 1650, yet only solved in the mid-18th century by Euler [18] (see also [19]). His solution established
the relation ζ(2) = π
2
/6 and, in general, expressed all positive even values ζ(2n) as rational multiples of
even powers of π. Euler’s discovery also hinted for the first time at the profound, yet a priori not obvious,
relevance of π in questions related to the distribution of prime numbers. Notably, no analogous relation has
been found for the odd values ζ(2n + 1) to date, three centuries after Euler’s work.
1
arXiv:2308.11829v1 [cs.AI] 22 Aug 2023
Despite their tremendous impact, discoveries of new relations between mathematical constants and other
mathematical structures are sporadic events that stem from strong intuitive leaps or great strokes of creativ-
ity. In light of this, a captivating possibility is that a deeper underlying concept exists, encompassing all or a
significant group of constants and providing a framework for classifying and ordering them. Here we propose
such a mathematical structure—called conservative matrix field—and show that it not only recreates known
relations between constants, but also reveals new ones. Through these relations, the matrix fields help to
order constants within a cohesive framework that may shed new light on their interrelationships, intrinsic
properties, and complexity.
The complexity of mathematical constants
Central to understanding the intrinsic nature of a mathematical constant is the question of its irrationality,
which can be understood as a measure of its complexity. A rational number can be regarded as possessing
the most simple representation, being expressible as a quotient of two integers. In turn, an irrational number
exhibits a higher level of complexity, for it cannot be represented by such a quotient.
Even though the question of whether a number is rational or irrational may appear innocent at first
glance, for numerous fundamental constants, determining their rationality can be a highly non-trivial task.
The situation for the odd values of the Riemann zeta function serves to illustrate this claim. Indeed, the
proof of the irrationality of ζ(3) by Apéry in 1978 [14, 15] stands as the only definitive accomplishment in
this regard, while the matter remains unsettled for all of the remaining odd values ζ(2n + 1); see [20–22] for
partial results in this direction. The odd values of the zeta function, along with the Catalan, Gompertz, and
Euler–Mascheroni constants, represent a large collection of prominent mathematical constants for which the
question of irrationality is still open.
The distinction rational vs. irrational offers a first kind of hierarchy that orders numbers according to their
complexity. An attempt to refine this binary hierarchy arises from the concept of irrationality measure [23], a
numerical value that quantifies how fast a constant may be approximated by an infinite sequence of distinct
rational numbers. In this sense, any number that admits a fast enough rational approximation must be
irrational. Rational numbers are characterized by an irrationality measure of zero, while irrational numbers
are characterized by an irrationality measure of at least one. Although there are numbers with arbitrary
irrationality measures greater than one, it is noteworthy that most real numbers possess an irrationality
measure of one [24]. Consequently, there is more to be desired from a further-refined dichotomy between
rational and irrational numbers.
Our endeavor to establish a finer hierarchy revolves around quantifying the complexity of the actual
formulas that may represent each number, and clustering them accordingly. We now explain our notion
of formulas and elaborate on why this specific class of formulas is highly suitable for computer-supported
research.
Algorithm-driven discovery of formulas for mathematical constants
Apéry’s celebrated proof of the irrationality of ζ(3) relies on producing an efficient approximation of this
number via the continued fraction formula
6
ζ(3)
= 5 −
1
6
117 −
2
6
535 −
3
6
.
.
.
−
n
6
17(n
3
+ (n + 1)
3
) − 12(2n + 1) −
.
.
.
(1)
Truncating such an infinite continued fraction at a finite number of steps yields a sequence of rational
numbers approximating the target constant. In this way, the exploration of continued fractions offers an
approach for estimating irrationality measures.
Within the class of continued fractions, formulas like Eq. 1 constitute the distinguished subclass of
polynomial continued fractions, for which the partial numerators and denominators are defined by polyno-
mials with integer coefficients as a function of the depth n. The study of polynomial continued fractions is
2
particularly attractive due to their broad mathematical scope and computational simplicity [25, §VI] (see
also [26–29]). Remarkably, these continued fractions encompass a broad range of mathematical functions
including trigonometric, Bessel, gamma, and hypergeometric functions; they represent a large collection
of constants [30]; they generalize a wide class of infinite sums [31]; and they represent linear recurrence
formulas with polynomial coefficients [25]. The computational simplicity of polynomial continued fractions
stems from two key aspects. Firstly, the space of polynomial continued fractions is countable, which enables
its systematic exploration and analysis. Secondly, the rational approximation by a polynomial continued
fraction can be evaluated very efficiently [32, §5.2]. Those two facts make this space a fitting candidate for
computer-supported research.
The Ramanujan Machine project [6], named to honor Srinivasa Ramanujan’s unique contributions to
mathematics [33], proposed to automate the process of discovery of formulas via an algorithmic approach.
The project developed algorithms to find formulas for constants numerically, yielding many new ones such as
notable formulas for Catalan’s constant. The most successful algorithm employed there relied on a brute-force
search over the space of polynomial continued fractions.
To broaden the exploration beyond the space that can be covered by brute-force searches, one is confronted
with the necessity to develop a more sophisticated exploration strategy. It was through the analysis of the
database of formulas obtained in [6] that we were able to identify an intriguing numerical property in the
space of formulas. This property, called factorial reduction (defined in Section 2), guided the development
of a more efficient algorithm for formula discovery.
Our novel algorithm is thus designed to search for formulas based on factorial reduction, discovering
hundreds of polynomial continued fraction formulas for well-known mathematical constants. This output
significantly exceeds what has been manually discovered in the past (e.g., [28]). The algorithm is heuristic
in nature, as it prunes the search space based on the conjecture that continued fractions representing known
constants tend to have factorial reduction. At any rate, the actual formulas that the algorithm discovers can
be and are independently verified. A search utilizing the factorial reduction property makes the algorithm
conducive to distributed computing, enabling us to implement a massively parallel search, with which we
covered a significantly larger space of formulas than what was possible before.
Our algorithm-driven approach exemplifies the vision of experimental mathematics [11,34,35]—leveraging
algorithms at scale to generate new insights in mathematics. Through the generation of large amounts of
data, algorithms offer valuable hints about patterns and behaviors of mathematical structures that would be
otherwise imperceptible within a reasonable timeframe. These hints guide mathematicians in formulating
stronger conjectures and even provide a blueprint for rigorous proofs and generalizations. This process
can be regarded as a cycle: the results of algorithmic experiments create new conjectures that themselves
inspire new algorithmic experiments, leading to further discoveries and a deeper understanding of the original
problem. Our work executes two cycles of an algorithmic-driven scientific method (Fig. 1), leading to the
discovery of an underlying model—a novel mathematical structure.
We see a parallel between the relationship of algorithms to mathematics and that of experiments to the-
oretical physics. Computational frameworks can serve as virtual laboratories for mathematicians, enabling
systematic exploration of mathematical hypotheses on a large scale. The greater computational power avail-
able, the more precise and comprehensive the experiments, providing potential breakthrough opportunities.
This collaborative approach between algorithms and experts facilitates the exploration and testing of hy-
potheses, leading to an accelerated understanding of unsolved problems and augments the serendipitous
insights of brilliant mathematicians.
3
Figure 1: Algorithmic-driven scientific method in experimental mathematics and its realization in our work. Our work
encompasses two cycles of the scientific method, leading to the discovery of a novel mathematical concept.
Discovery of a novel mathematical structure: conservative matrix field
The insight that led us to find the concept of conservative matrix fields resulted exclusively from having
the large dataset of polynomial continued fraction formulas discovered by our novel algorithm [36]. When
examining these formulas, certain recurring symmetries and patterns arise. As an example, note the resem-
blance between the following three formulas for π:
(2)
The first one was discovered by Brouncker and proven by Euler (see [37]), the second one by the Ramanujan
Machine [6], and the third one by Pickett et al. [38]. Such striking similarities were found in many more
examples in our dataset of formulas, sparking the question of whether there exists a governing structure
within the space of formulas for each given constant.
The massive number of formulas produced by our algorithm allowed us to identify several infinite para-
metric families of formulas representing constants such as e, ln(2), π and ζ(3) (see Section 3). A study of
the relationships between the formulas revealed a mathematical structure, coined conservative matrix field
(Section 3), which unifies the infinite families and even derives new formulas with faster rates of convergence.
This structure was identified by observing a kind of “conservation law”, reminiscent of the defining property
of conservative vector fields in physics.
A remarkable application of conservative matrix fields is identifying promising paths for proving the
irrationality of constants, as we exemplify for ζ(3) [39]. This approach now directly applies to more cases
(Appendix I).
4
We put forward the conjecture that these matrix fields have the potential to unify all polynomial continued
fraction formulas associated with a specific constant. Interestingly, matrix fields also unveil connections be-
tween different constants, for instance, between π
2
and Catalan’s constant, π
3
and ζ(3), or e and Gompertz’s
constant. Based on these findings, we hypothesize that many constants, including values of the Riemann
zeta function, can be interconnected through the framework of conservative matrix fields, supporting the
idea that they give rise to a new hierarchy of mathematical constants.
2 The Distributed Factorial Reduction algorithm
Our algorithm aims to discover formulas that equate linear fractional transformations of mathematical
constants and polynomial continued fractions, both with integer coefficients. Due to the nature of this prob-
lem, the search space is intrinsically infinite. Even if we limit our search to a small list of constants, the search
space encompasses all their combinations with arbitrary coefficients, and the number of combinations grows
exponentially with the maximal coefficients allowed. A greater difficulty lies in the fact that an attempt to
explore the search space in parallel by distributing it across different workers leads to redundant computa-
tions of the same expressions. Specifically, division of the search space between workers over the space of
constants requires repeated computing of the same continued fraction by different workers, each equating
the result to different constants or different transformations of constants. Moreover, sharing the already
computed results among workers to avoid such redundancy leads to a prohibitively high communication
overhead between them.
Our algorithm overcomes these challenges by identifying continued fractions that possess a novel property—
which we call factorial reduction. This identification enables us to pinpoint which continued fractions relate
to (well-known) mathematical constants without having to equate them to the (infinite) space of possible
constants and transformations of constants. Our algorithm thus avoids the need to scan the space of con-
stants, which substantially reduces its complexity, and most importantly, makes it amenable to large-scale
distributed computing. This way, we scan the search space across non-communicating workers without
causing redundant calculations and avoid the communication bottleneck.
We conjecture that the property of factorial reduction is a signature of recurrence sequences that converge
to mathematical constants. Consequently, our algorithm not only searches for conjectured formulas but is
itself based on a conjecture. Importantly, every formula generated by the algorithm is independently verified,
thereby liberating it from the necessity of justifying the factorial reduction algorithm, making it a stand-alone
conjecture awaiting formal proof.
This conjecture-based algorithm proved to be extremely effective, leading to the discovery of hundreds
of new formulas for mathematical constants, many of which would have been hard or impossible to discover
by other means. Identifying this conjectured property was enabled by the abundance of formulas discovered
by the older algorithms of the Ramanujan Machine project [36]. In retrospect, the substantial number of
independently verified formulas exhibiting factorial reduction reinforces the belief that factorial reduction is
a distinctive feature of formulas converging to constants of interest.
Factorial reduction: an observation from the ‘laboratory’ of experimental math-
ematics
At the core of our algorithm lies a novel observation about the greatest common divisors (GCDs) g
n
of the
numerators p
n
and denominators q
n
of the convergents of polynomial continued fractions. The convergents
are the rational numbers p
n
/q
n
found by truncating an infinite continued fraction. More precisely, if a
n
and
b
n
are the partial denominators and numerators of a continued fraction respectively (as in the rightmost
term of Eq. 2), then
p
n
q
n
= a
0
+
b
1
a
1
+
b
2
a
2
+
.
.
.
+
b
n
a
n
. (3)
These p
n
and q
n
are defined by the recursion u
n
= a
n
u
n−1
+ b
n
u
n−2
with initial conditions p
−1
= 1, p
0
= a
0
and q
−1
= 0, q
0
= 1. The growth rate of the resulting convergent numerator p
n
and denominator q
n
is a
power of a factorial, i.e. p
n
, q
n
∼ (n!)
d
for some positive integer d. See [6,25,40,41] for additional background
on continued fractions.
5
Our algorithm relies on the observation that polynomial continued fractions that converge to mathemat-
ical constants exhibit a peculiar behavior: Denoting g
n
= GCD(p
n
, q
n
) as the greatest common divisor of
the convergents, the reduced numerator and denominator grow exponentially at most, i.e.
p
n
g
n
,
q
n
g
n
∼ s
n
, (4)
instead of the much faster factorial growth rate of p
n
and q
n
before reduction. We name this observed
property factorial reduction. See our complementary work [42] for additional details about this concept.
Factorial reduction is found to be extremely rare: a wide search hints that the space of continued fractions
with factorial reduction is a “lower-dimensional” subspace inside the space of all polynomial continued frac-
tions. However, surprisingly, we observed factorial reduction in every polynomial continued fraction formula
found by the Ramanujan Machine using previous algorithms (see Appendix A) for a wide range of constants,
including π, ζ(3), and Catalan’s constant. In fact, we examined formulas converging to these constants in
literature spanning over hundreds of years (e.g., [18, 25–28,31, 33,41,43, 44]). Interestingly, all the continued
fractions and infinite sums [45] we checked displayed the property of factorial reduction.
The algorithm presented in this work relies on this unexpected, albeit simple property. Concretely, the
algorithm numerically tests the growth rate of reduced convergents of polynomial continued fractions, saving
the ones that have exponential rates (rather than factorial rates). This approach alleviates the need to
identify a limit constant in advance. Thus, the property of factorial reduction enables efficient identification
of continued fractions that converge to mathematical constants. Fig. 2 displays the test for several polynomial
continued fractions, showcasing the efficacy of our algorithm.
Figure 2: Observation of the factorial reduction (FR) property - its connection to different constants and its sparsity
in the wide space of continued fractions. (a) Calculation of the growth rate s =
n
p
|p
n
/g
n
| of the reduced convergents
for different continued fractions, detailed in the table. We plot a negative s value to denote a negative limit to p
n
/q
n
.
Even very similar continued fractions can have vastly different growth rates. We observe exponential growth whenever
the continued fraction converges to a mathematical constant. This approach identifies formulas for a wide range of
mathematical constants. (b) Comparison of the growth rate of continued fractions of the same form shows that the only
one having factorial reduction is precisely the one found by Apéry [14, 15] to converge to ζ(3). The factorial reduction
property is extremely rare, shown by its sensitivity to slight changes in the parameters of the continued fraction. (c) Table
of the polynomial continued fractions calculated in panel (a), denoting which ones have or do not have factorial reduction,
and which ones converge to a known constant or to a constant with no known closed form. Here, G is Catalan’s constant,
and
ˆ
ζ(5, 1) = ζ(5) − ζ(4) + ζ(3) − ζ(2) + 1 is the continued fraction from Eq. 8 with s = 5, R = 1.
6
The concept of factorial reduction not only provides a useful identification tool for formulas of constants
but is also a mathematical curiosity on its own [42]. This concept is relevant to Apéry-like irrationality
proofs [46], and can help prove the irrationality of other constants (see also the notion of “integer-ating
factor” in [47]).
Distributing the search for formulas
The identification of factorial reduction made possible the formulation of our Distributed Factorial Re-
duction (DFR) algorithm, whose key steps are (Fig. 3):
• We create a scheme on-site that defines the search space of polynomial continued fractions.
• This scheme is sent to a server, where it is parsed into smaller chunks, designed to compute in reasonable
time on a standard home computer. The distribution of the computing chunks is done using the
Berkeley Open Infrastructure for Network Computing (BOINC) [16].
• Each chunk is computed by off-site worker computers donated by volunteers, checking every polynomial
continued fraction in the chunk for the factorial reduction property. Each worker returns whether
factorial reduction was identified (a very rare event).
• Every identified case of factorial reduction is verified automatically on-site.
• We then attempt to match each verified case to specific mathematical constants using the PSLQ
algorithm [48,49]. We save all the verified cases, with and without matching constants. Each continued
fraction for which matching constants were found is further validated by calculating it to a greater depth
and comparing more digits, before storing the match as a new conjectured formula.
Figure 3: Implementation of the Distributed Factorial Reduction (DFR) algorithm. On-site scheme creation: Our
database collects formulas from the literature and from past results of algorithms. This collection helps us define schemes
of parameterized polynomial continued fractions. The range of parameters of the search space is decided based on our
computational power. In this work, our search spaces are families of polynomial continued fractions. Off-site factorial
reduction testing: Our BOINC resource manager server distributes the search space between the workers, who evaluate the
polynomial continued fractions and test each one for factorial reduction. This step is the most computationally intensive
part of the algorithm. The computing power that enables this effort is donated by volunteers from the BOINC community,
and, at the time of writing, involves over 5000 personal computers. The rare events of positive identification are sent back
for further analysis. On-site result verification: each identified continued fraction is first independently verified for factorial
reduction. We then attempt to match each verified continued fraction to constants using PSLQ [48, 49]. Each PSLQ
match is then validated to greater precision. We also collect the results for which a PSLQ match was not found, to enable
future tests with additional constants, refined applications of PSLQ, and potentially superior algorithms that may emerge.
7
Although we have been able to identify constants for a significant majority of factorial reduction cases,
we have encountered certain instances where a corresponding constant could not be identified. This is
not surprising, as PSLQ is only able to test against the finite list of constants provided to it as an input.
Nevertheless, we also store as formula-candidates all the results for which PSLQ failed, i.e. found to have
factorial reduction but without a connection to concrete constants (examples in Table 4 in Appendix C).
The most intensive part of the search is the identification of factorial reduction, due to the large space of
candidate polynomial continued fractions. For example, consider a search over polynomials a
n
, b
n
of degrees
5 and 10 respectively. This specific subspace of polynomial continued fractions is of interest because it was
found to contain new formulas for ζ(5), a prominent constant for which the irrationality question remains
open. For an intuitive estimate of the size of this search space, consider that even only testing polynomials
with coefficients between -10 to 10 already contains 3 × 10
22
(74 bit) candidate continued fractions, an
enormous effort for existing computational capabilities. We reduce the search space by selecting specific
(non-canonical) representations of the polynomials, as described in Appendix C. Even this limited space
requires immense computational resources.
Distributed computing is used in this part of the algorithm, to enable a search over such large spaces.
We executed the distributed algorithm first on a Technion computing cluster (Zeus), and then distributed
it on individual computers with the help of the BOINC community [50]. The thousands of new formulas
discovered in this project are in large part thanks to the contribution of the BOINC community. The choice
of distributed computing via BOINC is especially valuable as it allows everyone, novices or experts, the
opportunity to contribute to the discovery of new mathematics.
Through the support of the BOINC community, we are able to run our algorithm on extremely large search
spaces, making it the largest computational experiment on formulas of mathematical constants. Specifically,
our Distributed Factorial Reduction algorithm has been running online since October 2021. At the time of
writing, the algorithm is distributed over >5000 computers, and has been supported by >1000 volunteers
throughout its operation.
Example formulas discovered by our algorithm
Our distributed algorithm led to the discovery of an unprecedented number of continued fraction formulas,
creating a dataset that serves as the basis for the discoveries presented in this work. This section presents
selected examples of formulas from this dataset, all of which were unknown before this work, and most
remain so far unproven. A review of the results can be found on our website [36], and the full list is being
prepared in an online library format [51], importable for future projects in experimental mathematics.
We first present selected formulas for values of the Riemann zeta function. Our Distributed Factorial
Reduction algorithm found results such as:
72
72ζ(2) − 115
= 21 −
1
25 −
16
.
.
.
−
n
4
n
2
+ (n + 1)
2
+ 20 −
.
.
.
. (5)
This formula and others can be generalized to infinite families that connect to the Lerch zeta function Φ
(Appendix C). This function has explicit expressions for some of its values that combine more than one
constant. For example, the following formula connects ζ(3) and π
3
:
9
Φ(1, 3,
7
3
)
=
9
13ζ(3) −
1755
64
+
2π
3
3
√
3
= 65 −
81 · 1
6
249 −
81 · 2
6
.
.
.
−
81 · n
6
9(n
3
+ (n + 1)
3
) + 56(2n + 1) −
.
.
.
. (6)
The complexity of this formula and the examples below emphasize the prospects of the Distributed Factorial
Reduction algorithm, since existing algorithms could not have discovered them.
8
Our algorithm discovered formulas that mix different zeta values, e.g.,
1
ζ(3) − ζ(2) + 1
= 2 −
2
13 −
96
.
.
.
−
n
5
(n + 1)
n
3
+ (n + 1)
3
+ (n + 1)
2
−
.
.
.
, (7)
1
ζ(7) − 4ζ(3) + 4
= 5 −
1
333 −
16384
.
.
.
−
n
14
n
7
+ (n + 1)
7
+ 8(n
5
+ (n + 1)
5
) − 8(n
3
+ (n + 1)
3
) + 4(2n + 1) −
.
.
.
.
The search also discovered infinite families of formulas that connect an arbitrary number of integer zeta
values, with arguments up to half the degree of the polynomial b
n
in the continued fraction. We present one
of these families, with polynomial degree s and a root shifted by R, denoting its limit by 1/
ˆ
ζ(s, R):
1
ζ(s) + α
s−1
ζ(s − 1) + ... + α
1
= 1 + R −
1
2s−1
(1 + R)
1 + 2
s
+ 2
s−1
R −
2
2s−1
(2 + R)
.
.
.
−
n
2s−1
(n + R)
n
s
+ (n + 1)
s
+ R(n + 1)
s−1
−
.
.
.
, (8)
where α
1
, ..., α
s−1
∈ Q can be derived for every R ∈ Q (Appendix C). This general family of continued
fractions is merely a subset of a larger infinite family that constructs a parametric set of a
n
polynomials for
each b
n
with rational roots. We explore this family and provide a proof in Appendix C.
Our algorithmic-driven research led to the discovery of formulas for other constants, such as the Catalan
constant G
1
2G
= 1 −
2
7 −
32
.
.
.
−
2n
4
3n
2
+ 3n + 1 −
.
.
.
.
The automated search led to many formulas for algebraic numbers, including ones with degrees greater than
2, such as
3
√
4 + 1
3
√
4 − 1
= 5 −
8
15 −
35
.
.
.
−
(3n)
2
− 1
5 + 10n −
.
.
.
.
This continued fraction was also recently discovered in [29], and is found here to be part of a parametric
family of formulas for arbitrary-degree roots
c
√
c + 1 + 1
c
√
c + 1 − 1
= (2 + c) −
c
2
− 1
3(2 + c) −
(2c)
2
− 1
.
.
.
−
(cn)
2
− 1
(1 + 2n)(2 + c) −
.
.
.
.
Appendix H further shows that this parametric family is itself a particular case of a wider family of continued
fractions whose limit is a quotient of hypergeometric functions
2
F
1
[25, §49].
9
3 Conservative matrix fields
After implementing our Distributed Factorial Reduction algorithm, we seek to extract insights from
the numerous resulting formulas. One such approach involves identifying clusters of formulas that exhibit
common patterns and relate to distinct mathematical constants.
Remarkably, such clusters emerge naturally (see for example Fig. 4), and exploring their properties
leads to novel insights. Specifically, we present a mathematical structure that generalizes infinitely many
polynomial continued fractions. This structure possesses intriguing mathematical properties and enables
deriving additional polynomial continued fractions with their own unique characteristics. Notably, this
structure leads to the polynomial continued fraction that Apéry employed to demonstrate the irrationality
of ζ(3) and further enables the generalization of his approach to other constants.
Infinite families of formulas found by the algorithm
The Distributed Factorial Reduction algorithm produces a multitude of formulas, which we group into
infinite parametric families. One illustrative example of an infinite family of formulas found by our algorithm
is related to the constant ζ(3) (Fig. 4). This family of formulas is parameterized by α. For any rational
value of α, the limit was computationally observed to be a value of the polygamma function ψ
(2)
. Each
integer value of α provides a formula for ζ(3) = −ψ
(2)
(1)/2 up to a linear fractional transformation.
Applying the corresponding inverse transformation to each of the formulas for integer values of α yields
an infinite family of sequences that all converge to ζ(3). We construct a new sequence of rational numbers
that converge to ζ(3) by selecting the n-th element from the n-th sequence (illustrated in Fig. 4). The
new sequence exhibits a faster convergence rate than any of the continued fractions that constructed it,
providing a more efficient representation of ζ(3). This sequence is exactly the one used by Apéry to prove
the irrationality of ζ(3).
Figure 4: Generating an efficient converging sequence from an infinite family of continued fractions. The table presents a
parametric family of continued fractions. The rows correspond to integer values of α, which provide formulas converging to
linear fractional transformations of ζ(3). Inverting this transformation for every element in each sequence creates sequences
that all approach ζ(3). Then, we construct a new sequence by sampling the "diagonal" trajectory along the 2-dimensional
grid of sequences, i.e., we select each consequent element from the consequent sequence. This constructed sequence
creates an efficient approximation of ζ(3) that proves its irrationality.
The implication of this observation is that it is possible to create efficient formulas by building on infinite
families of continued fractions. The next section formalizes this procedure, revealing a novel mathematical
structure that generalizes our findings and even connects different mathematical constants. Further examples
of such infinite families are presented in Appendix C.
10
Unifying formulas using the conservative matrix field
We propose a mathematical structure that unites families of continued fractions representing a specific
constant. To explain the origin and implications of this concept, we pass to a matrix representation of
continued fractions and use the formulas of ζ(3) (Fig. 4) as the leading example.
Being a composition of linear fractional transformations, any continued fraction can be represented as a
multiplication of 2 × 2 matrices. The recurrence formula Eq. 3 satisfied by the convergent numerators p
n
and denominators q
n
can be translated into:
p
n−1
p
n
q
n−1
q
n
| {z }
V (n+1)
=
1 a
0
0 1
| {z }
V (1)
·
0 b
1
1 a
1
| {z }
M
X
(1)
· ··· ·
0 b
n
1 a
n
| {z }
M
X
(n)
. (9)
The matrix V (n) is understood as the “state matrix” of the system at the integer point n = 1, 2, 3, . . .
The matrix M
X
(n) will represent the “work matrix” required for a one-step displacement from n to n + 1,
changing states from V (n) to V (n + 1). Eq. 9 then describes the displacement from the integer point 1 to
n + 1. At the limit n → ∞, both ratios p
n−1
/q
n−1
and p
n
/q
n
of the column vectors of the matrix V (n + 1)
tend to the limit of the continued fraction.
Consider now the infinite family of formulas, indexed by the integer parameter α, listed in the table of
Fig. 4. Each formula converges to a fractional linear transformation of ζ(3). Adding the dependency on the
index α to Eq. 9 for each individual formula gives rise to the equation:
V (n + 1, α) = V (1, α) · M
X
(1, α) · ··· · M
X
(n, α), (10)
where the matrix M
X
(x, y) is defined as
M
X
(x, y) =
0 −x
6
1 x
3
+ (x + 1)
3
+ 2y(y − 1)(2x + 1)
. (11)
Eq. 10 can be interpreted graphically as consequent displacements rightward along horizontal trajectories.
We can arrange the horizontal trajectories corresponding to α = 1, 2, 3, . . . as parallel lines ordered from
bottom to top. This arrangement completes a 2-dimensional grid. For convenience, we relabel the parameters
(n, α) by (x, y).
We now wish to examine the work necessary for a vertical unit displacement in the grid, from V (x, y)
to V (x, y + 1). We denote the work matrix for such a displacement by M
Y
(x, y). For this notation to be
well-defined, the horizontal and vertical displacements must commute in the sense of satisfying
M
X
(x, y) · M
Y
(x + 1, y) = M
Y
(x, y) · M
X
(x, y + 1) (12)
for every pair of integers x, y ≥ 1. We illustrate this requirement in Fig. 5(b). In general, given an arbitrary
M
X
(x, y), there is a priori no reason why the connecting matrices M
Y
(x, y) should have a polynomial
dependency on (x, y). Surprisingly, in our example, M
Y
(x, y) takes the following form:
M
Y
(x, y) =
−(x − y)(x
2
− xy + y
2
) −x
6
1 (x + y)(x
2
+ xy + y
2
)
. (13)
The existence of polynomial matrices M
Y
(x, y) satisfying this equation in this example is quite remarkable.
A conservative matrix field is defined as a choice of 2×2 matrices M
X
(x, y) and M
Y
(x, y) with integer
polynomials in x and y, such that the condition in Eq. 12 holds for a regime of x, y in which the determinants
of M
X
, M
Y
are non-zero. We call Eq. 12 the conservative property (or the cocycle equation). The conservative
property guarantees that the work for a displacement along any trajectory in the 2-dimensional grid is given
by a matrix that only depends on its initial and final coordinates. This path-invariance feature is analogous to
the inherent property of conservative vector fields, with path integration replaced by matrix multiplications.
Fig. 5 displays this striking analogy.
11
Figure 5: Comparison of conservative vector fields and conservative matrix fields. (a) Conservative vector fields exhibit
path-independence, meaning that integrals taken over different paths yield the same result. (b) Analogously, conservative
matrix fields also offer path-independence, but for multiplications along the paths. The resulting matrices only depend
on the initial and final positions. Interestingly, conservative matrix fields are found to exhibit additional features: paths
going to infinity represent continued fraction formulas. (c) Each presented path provides a different formula of the same
mathematical constant. The example continued fractions here are derived from the conservative matrix field of e in Appendix
D. Such continued fractions can be derived for any constant that has a conservative matrix field. We hypothesize that all
continued fractions and infinite sums converging to a certain constant can be derived from a single conservative matrix
field of that constant.
Drawing on this analogy, the state V (x, y) of the conservative matrix field can be understood as a “matrix
potential”, expressing the work matrix of any trajectory from a fixed origin (e.g., (1, 1)) to the point (x, y).
Concretely, V (x, y) can be calculated as a successive multiplication of matrices M
X
and M
Y
along any
trajectory connecting the origin to the point (x, y). For example, two possible trajectories are
V (x, y) = M
X
(1, 1) · M
X
(2, 1) · ··· · M
X
(x − 1, 1) · M
Y
(x, 1) · M
Y
(x, 2) · ··· · M
Y
(x, y − 1)
V (x, y) = M
Y
(1, 1) · M
Y
(1, 2) · ··· · M
Y
(1, y − 1) · M
X
(1, y) · M
X
(2, y) · ··· ·M
X
(x − 1, y)
i.e., V (x, y) can be computed either by first performing a rightward displacement of x −1 units and then an
upward displacement of y − 1 units, or vice versa.
We define the limit along a trajectory starting at the origin as the ratio of the entries in the right column
of the potential matrix V as we proceed on the trajectory toward infinity. In the special case of a horizontal
trajectory in the example matrix field of Eqs. 11,13, this limit corresponds to the continued fraction in
the first row of Fig. 4. It is a remarkable observation that the limit along any trajectory of our example
conservative matrix field converges to 1/ζ(3). This observation is proven in the complementary article [39].
Building on this remarkable property, we can derive an abundance of new formulas of ζ(3), each from a
different trajectory direction in the conservative matrix field. These formulas serve as the first outcomes
from exploring the implication of the concept of conservative matrix fields.
12
It is even more remarkable that the property of a trajectory-invariant limit prevails across multiple
conservative matrix fields, spawning infinite families of formulas for each mathematical constant, each tied
to a distinct trajectory.
We produced a rich collection of conservative matrix fields (summarized in Appendix D) that were all
originally discovered via the Distributed Factorial Reduction algorithm, originating from infinite families of
formulas analogous to that in Fig. 4. These matrix fields in turn led us to discover an analytic approach
to constructing additional conservative matrix fields purely from the conservation property of Eq. 12 (sum-
marized in Appendix H). This analytical construction shows that conservative matrix fields generalize the
notion of rational Wilf–Zeilberger pairs [34]. Interestingly, our constructed matrix fields for constants such
as π and ζ(2) also demonstrated the property of having the same limit along all infinite trajectories, just
as the matrix fields discovered by the Distributed Factorial Reduction algorithm. The occurrence of this
peculiar phenomenon suggests the existence of a fundamental principle underlying a distinguished class of
conservative matrix fields with well-defined limits. The axiomatization of this class could broaden the ap-
plications of conservative matrix fields, from refining methods for numerical approximations and aiding in
symbolic computations, to providing new paths for irrationality proofs.
The definition of a conservative matrix field requires two matrices M
X
(x, y) and M
Y
(x, y) that commute
in the sense of Eq. 12. This definition can be naturally generalized to encompass a set of d pairwise
commuting matrices
M
X
1
(x
1
, . . . , x
d
), . . . , M
X
d
(x
1
, . . . , x
d
).
We call such a structure a d-dimensional conservative matrix field (each matrix is still 2 × 2). Appendix E
presents an example of a 4-dimensional conservative matrix field that converges to ζ(2) along every infinite
trajectory.
Connections between constants via conservative matrix fields
So far we have seen that each conservative matrix field generates infinite different formulas of a single
constant by following different straight trajectories. Each such formula can be directly translated to a
continued fraction form [52], and interestingly, the property of factorial reduction appears to be preserved
in all of them (once a single direction had it). These observations suggest an appealing hypothesis — that
all continued fraction formulas of a certain constant can be derived from a single conservative matrix field
(of two or more dimensions).
Further computational experiments reveal that conservative matrix fields can establish new relationships
between various mathematical constants. Specifically, by initializing the trajectory at a point other than
(1, 1)—for example, starting at (1/2, 1) (with the same M
X
and M
Y
)—produces a sequence that converges to
a different limit than the one starting at (1, 1), but arises from the same conservative matrix field. If this shift
is integer, then the resulting limit will be related to the same constant by a linear fractional transformation.
Conversely, when the shift is a non-integer rational number, in most cases the trajectory produces a formula
that converges to a different limit that is not related to the original constant. For example, π and ln 2 both
emerge from the same conservative matrix field via a shift of 1/2 of the trajectory initialization in either
axis. A similar connection exists for other constants such as ζ(2) and Catalan’s constant, or ζ(3) and π
3
.
These connections transfer formulas that were developed for one constant to the other constant, and may
help translate other properties between the constants. e.g. the 4-dimensional conservative matrix field of
ζ(2) in Appendix E directly gives rise to a similar matrix field for G. Constants that share their matrix field
up to such rational shifts can be placed at the same level of the hierarchy — the matrix field guarantees that
they are derived from recurrence formulas with the same structure, with polynomials of the same degree.
Different trajectories in a conservative matrix field can express the same constant by formulas with vastly
different attributes. Altering the trajectory may enable generating formulas with desired attributes like more
efficient convergence, crucial for proving the irrationality of various constants. In the section below, we follow
this observation and analyze straight trajectories with different slopes in the 2-dimensional grid. We show
how to provide quantitative formulas that yield different lower bounds on the irrationality measure of the
target constant.
13
Proofs of irrationality
Ever since Apéry’s breakthrough proof of the irrationality of ζ(3), there have been various attempts to
generalize his argument to provide irrationality proofs for other constants. In this section, we will show that
the conservative matrix field of ζ(3) can be used to provide a systematic approach to prove its irrationality.
The same approach extends to different conservative matrix fields, opening a pathway for irrationality proofs
for other constants.
All these proofs are based on the premise that any continued fraction that converges to a constant
produces a sequence of rational approximations, whose fast enough convergence implies the irrationality
of the constant. This argument can be quantified by the Liouville–Roth irrationality measure [23]. The
irrationality measure of a real number L is the largest possible number δ for which there exists a sequence
of distinct rational numbers {p
n
/q
n
} that converges to L and satisfies the inequality
p
n
q
n
− L
<
1
q
1+δ
n
(14)
for sufficiently large n. Rational numbers have irrationality measure δ = 0 and irrational numbers have
irrationality measure δ ≥ 1 by a classical result of Dirichlet. These two facts combined give rise to the next
irrationality criterion: if an infinite sequence of rational numbers converging to L has a positive irrationality
measure, then L is irrational.
It is possible to define the irrationality measure δ of a rational sequence as the largest value that satisfies
Eq. 14 for a sub-sequence. This δ bounds the irrationality measure of the constant to which the sequence
converges, and can thus be used to prove the irrationality of the constant. For example, Apéry’s continued
fraction in Eq. 1 has an irrationality measure δ ≈ 0.08, which is how Apéry proved the irrationality of
ζ(3) [14, 15].
Conservative matrix fields enable deriving a closed-form formula for δ that generalizes the work by Apéry
to arbitrary polynomial continued fractions, providing a parametric family of rational approximations for
each constant. The formula followed from our observation of a strong form of factorial reduction that
appears in matrix fields (after balancing their degrees [53]): Let g(x, y) = GCD(V (x, y)) be the greatest
common divisor calculated over the four entries of the matrix V . Numerical experiments show that most
of the conservative matrix fields that we found [54] have their reduced quotients V
ij
(x, y)/g(x, y) all grow
exponentially rather than factorially, for each straight trajectory t(n) = (x(n), y(n)) that goes to infinity in
the 2-dimensional grid. We hypothesize that this remarkable strong factorial reduction property is intrinsic
to all conservative matrix fields in which the polynomials in the four quotients have the same degree. The
base of the exponential growth s is trajectory dependent, and can be derived from any of the matrix entries
ij, as s
n
≈ V
ij
(t(n))/g(t(n)).
Our closed-form δ formula thus provides a different value for each trajectory:
1 + δ =
ln |e
max
/e
min
|
ln s
, (15)
where |e
min
| ≤ |e
max
| are the eigenvalues of the matrix that represents a step along the trajectory t(n) =
(x(n), y(n)) [55] in the limit of x, y going to infinity.
The irrationality of ζ(3)
To find the trajectory that offers the best δ in a conservative matrix field, we extract a numerical value
δ = δ(x, y) at every point in the 2-dimensional grid. δ(x, y) can be extracted from the following equation,
which arises from imposing equality in (14):
˜
V
12
(x, y)
˜
V
22
(x, y)
− L
=
1
˜
V
22
(x, y)
1+δ
. (16)
Here,
˜
V (x, y) = V (x, y)/g(x, y) is the reduced potential for every location (x, y) in the matrix field. Doing
so produces a function δ : N
2
→ R on the 2-dimensional grid. The limit of this function at infinity along
each direction provides the value of the irrationality measure of the rational approximation corresponding
14
to that direction. We thus search for the direction at which this limit is maximal. The resulting δ(x, y) for
a few examples are presented in Fig. 6. Every location where δ > 0 is colored in red. Hence, assembling
a sequence along a trajectory contained in the red regime of the figure yields a sequence that demonstrates
the irrationality of the constant.
Figure 6: Using conservative matrix fields to extract irrationality measures. (a) The extracted irrationality measures from
the conservative matrix field of ζ(3). Every infinite trajectory that remains in the red area yields a sequence with δ > 0,
which can prove the irrationality of ζ(3). The optimal δ is along the x = y trajectory, which is equivalent to Apery’s
continued fraction [14, 15, 39]. (b) The conservative matrix field of ζ(3) shifted by 1/3 along the x axis. The result is a
conservative matrix field that converges to a different constant: one that combines ζ(3) and π
3
; it is not known whether
this constant is rational or irrational. The optimal δ is no longer along the x = y trajectory. We present two methods to
detect the optimal δ trajectory through optimization algorithms, depicted in Appendix F and marked here by colored dotted
trajectories. Although this conservative matrix field does not contain a trajectory with δ > 0, the same methods provide
δ > 0 for other constants (Appendix I), and may yet provide a positive δ for this constant in a higher dimensional extension
of this matrix field. (c) The extracted δ for a conservative matrix field corresponding to
3
√
4. The optimal trajectory is
x = y, offering an irrationality-proving δ > 0.
To find an exact expression for the trajectory x = y, we define M
traj
(n) = M
X
(n, n) ·M
Y
(n + 1, n). The
resulting matrix (after balancing the degree [53]) is
M
traj
(n) = n
3
·
−(n + 1)
3
−6n
3
− 9n
2
− 5n − 1
6(n + 1)
3
35n
3
+ 54n
2
+ 30n + 6
, (17)
for which the eigenvalues are (for n → ∞) e
±
= 17 ± 12
√
2. A direct numerical test shows that log s ≈ 6.5,
and thus δ ≈ 0.08, matching the irrationality measure of the continued fraction found by Apéry [14, 15],
which was used in the first proof of the irrationality of ζ(3). Our complementary work [39] shows how to
use this conservative matrix field to prove the irrationality of ζ(3). In Appendix F, we show the equivalence
between this matrix representation and Apéry’s continued fraction.
Consequently, the conservative matrix field can be seen as an underlying algebraic structure that lies
behind Apéry’s continued fraction and constitutes his irrationality proof. This finding further motivates
searching for similar conservative matrix fields for other values of the Riemann zeta function and other
mathematical constants. As an example, the same method, applied on a conservative matrix field for ζ(2)
(Appendix D) also provides a way to prove the irrationality of this constant, providing δ ≈ 0.09. Generalizing
beyond ζ(2) and ζ(3), it remains to be discovered whether similar matrix fields exist for other ζ(n), and
whether they could help prove the irrationality of these constants. Appendix G shows a similar method
applied for higher zeta values, yielding sequences with non-trivial irrationality measures, though not positive
values so far.
15
Results for ζ(5) and the quest to prove its irrationality
It would be of great interest to create a conservative matrix field for ζ(5), for which irrationality is a long-
standing open question [21, 22]. Although we have not yet found such a matrix field, this section describes
several formulas for ζ(5) (Table 1) that may serve as the foundation for constructing a conservative matrix
field. Part of the formulas in Table 1 combine several values of the Riemann zeta function, hinting that they
belong to different conservative matrix fields (or to different rational shifts of the same matrix field). Three
of the discovered formulas in the table did not fit any combination of zetas (in numerical tests using PSLQ)
and are thus presented with their numerical values.
a
n
b
n
Formula
n
5
+ (n + 1)
5
+ 6(n
3
+ (n + 1)
3
) −n
10
2
2ζ(5)+6ζ(3)−9
n
5
+ (n + 1)
5
+ 6(n
3
+ (n + 1)
3
) − 4(2n + 1) −n
10
2
2ζ(5)−2ζ(3)−1
n
5
+ (n + 1)
5
+ 16(n
3
+ (n + 1)
3
) − 4(2n + 1) −n
10
64
64ζ(5)+176ζ(3)−273
8(n
5
+ (n + 1)
5
) − 15(n
3
+ (n + 1)
3
) + 9(2n + 1) −64n
10
1.20426...
8(n
5
+ (n + 1)
5
) − 12(n
3
+ (n + 1)
3
) + 7(2n + 1) −64n
10
2.45174...
8(n
5
+ (n + 1)
5
) + 20(n
3
+ (n + 1)
3
) − 5(2n + 1) −64n
10
22.8410...
Table 1: Sporadic results suspected as related to ζ(5), found via the Distributed Factorial Reduction algorithm for
polynomial continued fractions of degree 5. Part of the formulas were linked to a combination of ζ(5) and ζ(3).
The Distributed Factorial Reduction algorithm discovered infinite families of formulas encompassing
different combinations of zeta function values, such as Eq. 8. Although we did not identify a matrix
field that generalizes these continued fractions, we propose a method that uses them to extract meaningful
irrationality measures. We demonstrate this method on ζ(5). We construct new sequences that converge to
ζ(5) by assembling linear combinations of continued fractions. To find the coefficients of the combinations
c
i
, we use the closed-form expression of the limit of the continued fraction from Eq. 8, denoted by
ˆ
ζ(s, R).
For each set of parameters R
i
, we can determine rational coefficients c
i
such that:
ζ(5) = c
2
·
ˆ
ζ(s = 2, R
2
) + c
3
·
ˆ
ζ(s = 3, R
3
) + c
4
·
ˆ
ζ(s = 4, R
4
) + c
5
·
ˆ
ζ(s = 5, R
5
) (18)
To construct a sequence of rational numbers that converge to ζ(5), we substitute each
ˆ
ζ(s, R) by the sequence
of convergents of its continued fraction form (Eq. 8). Varying the R
i
parameters, we can generate infinitely
many sequences for ζ(5), organized in a structure reminiscent of Fig. 4. The analysis of the convergence
rates for different R
i
values yields nontrivial irrationality measures for ζ(5) (Appendix G). This method may
help improve the irrationality measure beyond the recently achieved record value [22].
4 Discussion
The incorporation of the distributed algorithm caused a dramatic increase in the number of formula
candidates obtained by the Ramanujan Machine, which, in turn, posed new kinds of algorithmic challenges.
Specifically, we are in need of an automated method to identify the relevant constants that connect to
each formula candidate, as well as a method to generalize a set of formulas to parametric families and to
conservative matrix fields. These challenges are discussed in the next subsections.
Algorithms for finding connections between continued fractions and mathemat-
ical constants
The Distributed Factorial Reduction algorithm identifies promising candidate continued fractions but
does not identify the mathematical constants to which they converge. To find a closed-form formula for
the continued fraction, we try to match a selection of candidate constants to the numerical value of the
discovered continued fraction using the integer-relation algorithm PSLQ [48,49]. Appendix A compares this
approach to other algorithms that we had previously developed [6] for this purpose.
16
The PSLQ algorithm accepts a vector of real numbers z
i
, and yields a vector of integers c
i
such that
P
c
i
z
i
= 0 up to a preset error. In our simplest use of PSLQ, for a continued fraction that evaluates to a
numerical value v and a candidate constant η, the input vector is (1, η, −v, −vη). If the algorithm finds a
set of integers c
0
, c
1
, c
2
, c
3
, then the formula
c
0
+c
1
η
c
2
+c
3
η
= v holds. Most of the results found in this work were
based on this application of PSLQ, resulting in conjectured formulas such as Eq. 5.
Importantly, the input vector of PSLQ can be extended to include several constants and powers thereof,
expanding the scope of our findings to more complex formulas. Successful examples of operating PSLQ with
such extended vectors led to the most intriguing results of this work, including the conjectural formula of
Eqs. 7 and 8 (more information in Appendix C). Our use of PSLQ can be taken further to include other
functions of constants (e.g. square roots, exponents) and even to connect several continued fractions. We
hypothesize that eventually, a refined version of PSLQ will find a closed-form formula for all the continued
fractions that were found in this work.
Result verification
Conjectures found by data-driven algorithms in experimental mathematics are tested via their numerical
approximations and are hence bound to the possibility of being false positives (until they are proven [56,57]).
In order to reduce the occurrence of false detection, we calculate the continued fraction to a greater depth [58]
and compare it with additional digits of the constants. Most formulas reported in this work have been tested
to at least 100 digits of precision. However, the slow convergence rate of certain formulas (e.g., the first few
members of the ζ(3) infinite family in Fig. 4) makes the task of guaranteeing that precision in these cases
computationally infeasible.
To acquire additional confidence in the slowly converging formulas, we are again inspired by experimental
physics, where one can mitigate errors by identifying patterns that support an underlying model. In such
cases, multiple observations—each possibly with low precision—can together lead to the discovery of a
model with high confidence, if a shared pattern supports that model. The same process translates to
our investigation in experimental mathematics. We compensate for the low precision in slowly converging
continued fractions by identifying shared patterns, like parametric families that generalize several formulas
(e.g., Fig. 4 and Appendix C). In such a case, even though every formula is verified only up to a limited
precision, the appearance of a family of formulas provides more confidence in the validity of each of its
members. For example, the first two members of the ζ(3) infinite family from Fig. 4 acquire less than 50
correct digits when calculated to a depth of 10
5
. Still, the faster convergence of other formulas in that family
adds to our confidence in the slower converging members.
In physics, for additional verification, candidate models can be used to generate new predictions for
targeted experiments in previously untested parameters. The successful validation of these predictions
serves to enhance our confidence in the physics model. Translating this approach to our case in experimental
mathematics, we tested the validity of every candidate parametric family of formulas by substituting new
(rational) parameters and evaluating the formula numerically. Breaking the analogy with physics, it is not
straightforward to assess the numerical error for particular continued fractions. We estimate the error based
on the rate of convergence of the continued fraction (Appendix B).
Lastly, we consider the innate beauty of mathematical representations as a tool for verification. PSLQ is
more successful when finding a formula that holds to a large precision using only small integers and only a few
constants. When the formula suggested by PSLQ is more compact, we usually consider it more “beautiful”,
perceiving the formula as more likely to be correct. A quantitative rule of thumb for success in PSLQ is
searching for results in which the total number of digits in all the integers found by PSLQ is significantly
smaller than the number of input digits (resulting from the numerical evaluation of a continued fraction).
From the point of view of the compactness of the representation, smaller integers imply that most of the
information provided by the continued fraction is held by the constant rather than by the integers. We also
gain increased confidence in a new candidate result of PSLQ when the resulting integers remain unchanged
after running PSLQ again with more digits as input.
17
Outlook and open questions
The new algorithmic approach shown here is the most successful ever in finding conjectured formulas for
mathematical constants. The sheer number of new formulas led to new challenges and motivated clustering
of formulas of the same constant. In certain cases, we identified such connections, which in turn unveiled
a novel mathematical structure—the conservative matrix field. This structure has created infinite families
of continued fractions, enabled irrationality proofs, and revealed a hidden hierarchy connecting different
mathematical constants.
Our algorithm-assisted research has inspired numerous intriguing, open questions. Many of these ques-
tions concern the properties of conservative matrix fields. The foremost question is whether indeed a unique
conservative matrix field exists for each constant, a premise that our findings appear to support so far. If
corroborated, this concept could have intriguing implications, suggesting that a unified matrix field may
encapsulate all the formulas corresponding to each specific constant. We must also establish the maximum
dimension (number of matrices) within the matrix field of each mathematical constant, as this attribute
could become a fundamental characteristic of each constant.
Capturing additional constants may be possible by generalizing the concept of a conservative matrix field
to larger matrices beyond 2 ×2, or to other analogous mathematical structures. It is also important to find
whether 2 × 2 conservative matrix fields exist at all for an arbitrary polynomial degree (Appendix H shows
matrix fields up to degree 3).
The hierarchy uncovered using conservative matrix fields reveals a classification of constants based on the
complexity of the formulas producing them. This ordered classification offers a valuable avenue for exploring
shared properties between constants of the same order. Consequently, conservative matrix fields emerge as
a potent and promising tool for identifying invariants among mathematical constants, paving the way for
novel insights into their fundamental nature.
At present, algorithmic-assisted research strategies are still at their infancy. Our work has detailed the
process of large-scale experimental exploration and discovery, unveiling many numerically-validated formulas,
which we subsequently generalized into a unified mathematical structure. Traditionally, such generalizations
are attributed to moments of inspiration—when a mathematician discerns a pattern or creates a novel
structure. We foresee this process of generalization becoming algorithmically-assisted in the future. This
can be achieved through techniques such as feature extraction and pattern recognition applied to the wealth
of data generated by large-scale mathematical experiments. Such algorithms could generate high-quality
hypotheses, automatically proposing generalizations for researchers to review, validate, and ultimately prove.
We could consider the heuristic of factorial reduction as an intriguing case study of automating the
act of generalization and hypothesis. This heuristic, which led to most of the results presented here, was
discovered by (manually) noticing an emergent pattern. In future works, identifying statistical anomalies or
surprising order within experimental results (like the unexpected prevalence of factorial reduction in formulas
of constants) could be done algorithmically as well—potentially identifying effects and structures that may
otherwise be overlooked.
Looking ahead, algorithmic approaches in experimental mathematics are set to provide ever-stronger
methods to study long-standing open questions. In the coming years, more sophisticated algorithms will be
tailored to generate conjectures of growing complexity. Such conjectures, formulated automatically, could
provide new clues and accelerate progress across various fields of mathematics, tackling profound questions
such as the structure and properties of fundamental constants.
Acknowledgements
This research is supported by the generosity of Eric and Wendy Schmidt by recommendation of the
Schmidt Futures Polymaths program. We are grateful to the volunteers in the BOINC community whose
contribution made this discovery possible.
18
References
[1] Siemion Fajtlowicz. On conjectures of Graffiti. Annals of Discrete Mathematics, 38:113–118, 1988.
[2] Kenneth I. Appel and Wolfgang Haken. Every planar map is four colorable, volume 98. Am. Math.
Soc., 1989.
[3] Doron Zeilberger. A holonomic systems approach to special functions identities. Journal of computa-
tional and applied mathematics, 32:321–368, 1990.
[4] Doron Zeilberger. A fast algorithm for proving terminating hypergeometric identities. Discrete mathe-
matics, 80:207–211, 1990.
[5] Herbert S Wilf and Doron Zeilberger. An algorithmic proof theory for hypergeometric (ordinary and
“q”) multisum/integral identities. Inventiones mathematicae, 108:575–633, 1992.
[6] Gal Raayoni, Shahar Gottlieb, Yahel Manor, George Pisha, Yoav Harris, Uri Mendlovic, Doron Haviv,
Yaron Hadad, and Ido Kaminer. Generating conjectures on fundamental constants with the Ramanujan
Machine. Nature, 590:67–73, 2021.
[7] Alhussein Fawzi, Matej Balog, Aja Huang, Thomas Hubert, Bernardino Romera-Paredes, Moham-
madamin Barekatain, Alexander Novikov, Francisco J R Ruiz, Julian Schrittwieser, Grzegorz Swirszcz,
et al. Discovering faster matrix multiplication algorithms with reinforcement learning. Nature, 610:47–
53, 2022.
[8] Hao Wang. Toward mechanical mathematics. IBM J.Res.Dev., 4:2–22, 1960.
[9] Randall Davis and Douglas B Lenat. Knowledge-Based Systems in Artificial Intelligence: 2 Case Studies.
McGraw-Hill, Inc., 1982.
[10] Douglas B Lenat and John Seely Brown. Why AM and EURISKO appear to work. Artificial intelligence,
23:269–294, 1984.
[11] Stephen Wolfram et al. A new kind of science, volume 5. Wolfram media Champaign, 2002.
[12] Bruno Buchberger, Adrian Crˇaciun, Tudor Jebelean, Laura Kovács, Temur Kutsia, Koji Nakagawa,
Florina Piroi, Nikolaj Popov, Judit Robu, Markus Rosenkranz, et al. Theorema: Towards computer-
aided mathematical theory exploration. Journal of applied logic, 4:470–504, 2006.
[13] Alex Davies, Petar Veličković, Lars Buesing, Sam Blackwell, Daniel Zheng, Nenad Tomašev, Richard
Tanburn, Peter Battaglia, Charles Blundell, András Juhász, et al. Advancing mathematics by guiding
human intuition with AI. Nature, 600:70–74, 2021.
[14] Roger Apéry. Irrationalité de ζ(2) et ζ(3). Astérisque, 61:1, 1979.
[15] Alf van der Poorten. A proof that Euler missed. Math. Intelligencer, 1:195–203, 1979.
[16] David P Anderson. Boinc: A system for public-resource computing and storage. Fifth IEEE/ACM
international workshop on grid computing, pages 4–10, 2004.
[17] Steven R. Finch. Mathematical constants. Cambridge University Press, 2003.
[18] Leonhard Euler. De summis serierum reciprocarum. Commentarii Academiae Scientiarum, 7:123–134,
1740.
[19] Raymond Ayoub. Euler and the Zeta Function. Am. Math. Mon., 81:1067–1086, 1974.
[20] Keith Ball and Tanguy Rivoal. Irrationalité d’une infinité de valeurs de la fonction zêta aux entiers
impairs. Inventiones mathematicae, 146:193–207, 2001.
[21] Wadim Zudilin. One of the numbers ζ(5), ζ(7), ζ(9), ζ(11) is irrational. Uspekhi Mat. Nauk, 56:149–150,
2001.
19
[22] Francis Brown and Wadim Zudilin. On cellular rational approximations to ζ(5). arXiv:2210.03391,
2022.
[23] Godfrey Harold Hardy, Edward Maitland Wright, et al. An introduction to the theory of numbers.
Oxford university press, 1979.
[24] Aleksander Ya. Khinchin. Continued fractions. The University of Chicago Press, 1964.
[25] Oskar Perron. Die Lehre von den Kettenbrüchen. Band II. Analytisch–funktionentheoretische Ketten-
brüche. 3. verbesserte und erweiterte Auflage. Teubner, 1957.
[26] Douglas Bowman and James McLaughlin. Polynomial continued fractions. Acta Arith., 103:329–342,
2002.
[27] James McLaughlin and Nancy J. Wyshinski. Real numbers with polynomial continued fraction expan-
sions. Acta Arith., 116:63–79, 2005.
[28] Annie Cuyt, Vigdis Petersen, Brigitte Verdonk, Haakon Waadeland, and William B. Jones. Handbook
of continued fractions for special functions. Springer Science & Business Media, 2008.
[29] Dzmitry Badziahin. On effective irrationality exponents of cubic irrationals. arXiv:2301.02391, 2023.
[30] By default, we shall use the term “formula” to refer to the polynomial continued fraction formulas, such
as (1). Our approach covers other types of formulas as well.
[31] Leonhard Euler. Introductio in analysin infinitorum, volume 1,2. Apud Marcum-Michaelem Bousquet
& Socios, 1748.
[32] William H Press, Saul A Teukolsky, William T Vettering, and Brian P Flannery. Numerical recipes in
C++ : the art of scientific computing. Cambridge University Press, 2002.
[33] Bruce C Berndt. Ramanujan’s notebooks: Part III. Springer Science & Business Media, 2012.
[34] Marko Petkovsek, Herbert S. Wilf, and Doron Zeilberger. A=B. AK Peters Ltd, 1996.
[35] David Bailey, Jonathan Borwein, Neil Calkin, Russell Luke, Roland Girgensohn, and Victor Moll.
Experimental mathematics in action. CRC press, 2007.
[36] Results found by the Ramanujan Machine. http://www.ramanujanmachine.com/results/.
[37] Jesus Guillera. History of the formulas and algorithms for pi. arXiv:0807.0872, 2009.
[38] Thomas J Pickett and Ann Coleman. Another continued fraction for π. Am. Math. Mon., 115:930–933,
2008.
[39] Ofir David. The conservative matrix field. arXiv:2303.09318, 2023.
[40] Hubert Stanley Wall. Analytic theory of continued fractions. Courier Dover Publications, 2018.
[41] Oskar Perron. Die Lehre von den Kettenbrüchen. Teubner, 1913.
[42] Nadav Ben David, Guy Nimri, Uri Mendlovic, Yahel Manor, and Ido Kaminer. On the Connection
Between Irrationality Measures and Polynomial Continued Fractions. arXiv:2111.04468, 2021.
[43] David Naccache and Ofer Yifrach-Stav. On Catalan Constant Continued Fractions. arXiv:2210.15669,
2022.
[44] David Naccache and Ofer Yifrach-Stav. The balkans continued fraction. arXiv:2308.06291, 2023.
[45] We can test an infinite sum for the factorial reduction property by converting it to a continued fraction
using Euler’s continued fraction formula, or by a direct examination of the growth of the reduced
numerators and denominators of the partial sums.
20
[46] Doron Zeilberger and Wadim Zudilin. Automatic discovery of irrationality proofs and irrationality
measures. Int. J. Number Theory, 17:815–825, 2021.
[47] Robert Dougherty-Bliss, Christoph Koutschan, and Doron Zeilberger. Tweaking the Beukers integrals
in search of more miraculous irrationality proofs a la Apéry. The Ramanujan Journal, 58:973–994, 2022.
[48] Helaman RP Ferguson, Daivd H Bailey, and Paul Kutler. A polynomial time, numerically stable integer
relation algorithm. Technical report, 1998.
[49] Helaman Ferguson, David Bailey, and Steve Arno. Analysis of PSLQ, an integer relation finding algo-
rithm. Mathematics of Computation, 68:351–369, 1999.
[50] The Ramanujan Machine BOINC community. https://www.ramanujanmachine.com/run-the-
ramanujan-machine/.
[51] Library of Integer RElations and Constants. https://github.com/RamanujanMachine/LIReC.
[52] Ofir Razon, Yoav Harris, Shahar Gottlieb, Dan Carmon, Ofir David, and Ido Kaminer. Automated
Search for Conjectures on Mathematical Constants using Analysis of Integer Sequences. International
Conference on Machine Learning, 202:28809–28842, 2023.
[53] We transform the matrices to reduce their maximal degree. This transform is motivated by the compu-
tational effort of matrix multiplication, which is dominated by the largest element in the matrix. For
example, the matrix in equation 11 contains a
n
∝ n
3
and b
n
∝ n
6
. By splitting the roots of b
n
to two
groups, and multiplying one group to the preceding matrix, we decrease the largest matrix entry to be
O(n
3
) instead of O(n
6
). This process offers a significant computational advantage, while only changing
the produced expressions by a linear fractional transform.
[54] Except for the conservative matrix field of e, where the convergence is also faster and thus a nontrivial
δ is still extracted.
[55] Let M
t
(n) be a matrix that represents the direction of the trajectory. For example, a diagonal going in
45
◦
is represented by M
t
(n) = M
Y
(n, n) · M
X
(n, n + 1), and a diagonal going in 30
◦
relative to the x
axis is represented by M
t
(n) = M
X
(2n, n) · M
X
(2n + 1, n) · M
Y
(2n + 1, n + 1).
[56] Robert Dougherty-Bliss and Doron Zeilberger. Automatic conjecturing and proving of exact values of
some infinite families of infinite continued fractions. The Ramanujan Journal, pages 1–17, 2020.
[57] Eric Brier, David Naccache, and Ofer Yifrach-Stav. A Note on the Ramanujan Machine.
arXiv:2211.01058, 2022.
[58] The number of terms is decided as the calculation takes place, depending on the rate of convergence,
which determines the (increaseing) depths at which we sample the decimal approximations.
21
A Past algorithms in the Ramanujan Machine project
This section describes the algorithms for discovery of formulas for mathematical constants that were
developed as part of the Ramanujan Machine project prior to this work.
• Meet-in-the-Middle [6] - Suggests conjectured formulas in the form of polynomial continued frac-
tions based on brute-force matching of decimal approximations. The algorithm computes polynomial
continued fractions from a user-defined scheme to a predetermined depth, creating a decimal approxi-
mation of the continued fraction. The decimal approximation is then matched against a pre-calculated
list of decimal approximations of rational functions of fundamental constants to identify potential for-
mulas. To ensure efficient execution, the algorithm incorporates bloom filters, caching of the sequences
a
n
and b
n
, and other optimizations. The following section provides a more comprehensive description
of this algorithm.
• Descent and Repel [6] - Proposes conjectured formulas in the form of polynomial continued fractions
(but applicable to arbitrary forms) by solving an optimization problem via an altered version of gradient
descent combined with rounding to the nearest integer parameters. The gradient descent is executed
with several points in parallel, incorporating a Coulomb-like repulsion between the points to enable a
better cover of the search space.
• Enumerated Signed-Continued-Fraction Massey Approve (ESMA) [52] - Finds conjectured
formulas in the form of signed continued fractions (where b
n
is a sequence made from ±1) by identifying
patterns in integer sequences. For every sequence of ± signs, each mathematical constant has a unique
signed continued fraction expansion in the form of an integer sequence. ESMA searches for patterns
in such integer sequences using a modification of the Berlekamp–Massey algorithm.
All of these algorithms are limited by the requirement to specify the target constants. i.e., one must guess
what constants are expected to be the exact convergence limit of the continued fractions. The choice of the
said constants is often based on intuition, or extracted from within a prescribed finite list. The Distributed
Factorial Reduction algorithm bypasses this limitation, as it searches for continued fractions without relying
on a convergence to any specific constant. The following section provides a brief review of the Meet-in-the-
Middle algorithm, exemplifying the above-mentioned challenge of choosing the target constant. The section
then concludes by explaining how factorial reduction serves to bypass this challenge.
The Meet-in-the-Middle algorithm
This algorithm receives a predetermined constant η for every execution and aims to find expressions of
the form:
a + b · η
c + d · η
= a
0
+
b
1
a
1
+
b
2
a
2
+
b
3
.
.
.
+
b
n
a
n
+
.
.
.
(19)
with a, b, c, d ∈ Z. More general versions of the algorithm enumerate over more complex expressions in the
left-hand-side. Either way, the algorithm has to enumerate over the spaces of values for each side of the
equation separately and to compare them.
First, the algorithm calculates many left-hand-side expressions. These expressions are stored in a hash
table, and the key assigned to every entry is the first 10 digits of the decimal approximation. The keys are
stored in a bloom filter, which allows for quick queries to determine whether an element is in the filter.
Next, the algorithm iterates over all polynomial continued fractions with numerators and denominators
with fixed polynomial degrees, and with integer coefficients in a pre-specified interval. Each continued fraction
is calculated to depth 30, and the first 10 digits of the resulting decimal approximation are computed. This
value is then searched for in the bloom filter. If the value is found, the algorithm indicates a candidate result
that is later tested further. Variants of the algorithm operate with a different number of stored digits and a
different depth of calculation.
22
Finally, to eliminate false positives and identify the remaining formulas, the algorithm calculates the
candidate polynomial continued fractions to higher precision by evaluating them to depth 10,000. The
algorithm then compares the first 100 digits to the exact formula.
Limitations of Meet-in-the-Middle that were alleviated in the Distributed Fac-
torial Reduction algorithm
While being very successful, Meet-in-the-Middle had several drawbacks that limited our ability to dis-
tribute it, and restricted the continued fractions that it can identify.
Firstly, the algorithm requires the storage of a bloom filter in memory and a hash table on disk. This
not only restricts the number of potential expressions on the left-hand-side of Eq. 19 that we can check, but
also complicates the distribution process, requiring the download of additional data for every off-site worker.
This limitation is not present in the Distributed Factorial Reduction algorithm, since it identifies continued
fractions without requiring prescribing the constant in the left-hand-side.
Secondly, the algorithm assumes that the continued fraction converges quickly enough to stabilize its
first few digits after a small number of iterations. Consequently, continued fractions that converge slowly
will not be detected under this assumption. For example, the second formula in the infinite family of ζ(3)
presented in Fig. 4 only gets 11 correct digits after calculating the first 100k terms of the continued fraction.
Therefore, such a formula cannot be detected by this algorithm. In contrast, the Distributed Factorial
Reduction algorithm does not rely on decimal approximations to detect potential conjectures, instead using
the greatest common divisor of p
n
and q
n
. Thus, the convergence rate of a continued fraction does not
prevent its detection. i.e., this algorithm can identify slow and fast continued fractions equally well.
B Verification of results
This section discusses the challenges of using finite numerical approximations when finding conjectures,
and the mechanisms we use to address these challenges. As mentioned in Sec. 4, these approximations
are prone to introduce errors to the process and consequently lead to false positive results. These false
positives are eliminated by the mechanisms detailed below, allowing us to estimate the magnitude of errors
and translate that error into a confidence measure in each result.
Approximation error in the evaluation of a continued fraction
Our algorithms use decimal representations for comparing the limit of a polynomial continued fraction
and the fundamental constant expression L. Hence, when discussing the precision of a continued fraction
PCF
n
at a finite depth n, we refer to the number of digits in the fraction’s approximation that are identical
to its limit. A precision of K means an approximation error ϵ
n
= |PCF
n
− L| < 10
−K
.
We observes in numerical experiments that ϵ
n
monotonically decreases with increasing depth for continued
fractions of the type used in this work. This observation suggests that for every required precision K, there
exists an N such that for every n > N the first K digits of the continued fraction calculated to depth n are
identical to the first K digits of the limit L:
⌊PCF
n
· 10
K
⌋ = ⌊L · 10
K
⌋
Note that for some edge cases, this equation will not be satisfied. For example, if the limit is an integer and
the sequence approaches it from below, we can get
P
n
i=1
9 · 10
−i
→ 1, while at no point the digits match.
However, since comparing digits is efficient computationally, and almost all limits of interest are non-integer,
this method satisfies our needs.
To find how many digits we can trust as already accurate, our algorithm attempts to find the largest K for
a finite predetermined n. To achieve this, we calculate the decimal approximation of the continued fraction
for two different depths, n
1
< n
2
. Since the error monotonically decreases, the first digit where PCF
n
1
and
PCF
n
2
do not match identifies a bound from above on the precision K offered by PCF
n
1
. By a careful choice
of n
1
and n
2
, we can decrease the upper bound calculated for K and get a better approximation for it.
Our algorithm finds this upper bound of K and halts the calculation when the upper bound exceeds the
desired precision K by a confidence margin, or when we reach the maximal allowed depth. In the latter case,
we return the best attainable K at this depth, allowing us to collect some useful data even for continued
fractions with slow convergence.
23
The choice of n
1
and n
2
is done by a form of exponential backoff algorithm. Most of our continued
fractions converge at least at a polynomial rate, with the slowest convergence rate thus being ϵ
n
∼ 1/n.
Such a convergence rate implies that some continued fractions will present the same incorrect digits over a
wide range of terms. This phenomenon makes it more challenging to distinguish digits that match the limit
from those that do not, and as a result makes evaluating K harder. To obtain a good approximation of
K, we want a significant difference between PCF
n
1
and PCF
n
2
, since this provides higher confidence in the
accuracy of the digits of PCF
n
1
that are still equal to PCF
n
2
, and therefore a tighter bound on K.
We sample the continued fraction at depths n
1
, n
2
, ..., choosing n
l+1
−n
l
∝ l, which matches the sampling
rate to the convergence rate in the slowest case of ϵ
n
∼ 1/n. This sampling rate makes the calculated decimal
approximation satisfies the requirement for substantial difference between subsequent PCF
n
l
, resulting in a
tight bound of K. As a result of this process, we attain a decimal approximation of the continued fraction’s
limit, along with the number of digits we refer to as “correct”.
Identifying over-fitted formulas discovered by PSLQ
To match a linear fractional transformation (i.e., a Mobïus transform) of a constant to a decimal
approximation, we utilize PSLQ to find integer coefficients ⃗a = (a
0
, a
1
, a
2
, a
3
) that satisfy the equation
a
0
+a
1
η
a
2
+a
3
η
= PCF
n
for a given constant η up to a given precision. Using the method described in the previous
section, we employ finite approximations of the continued fraction to evaluate the first K correct digits of
the limit. Given the approximation error of the continued fraction (∼ 10
−K
), we require PSLQ to find ⃗a
that satisfies the formula to the same accuracy.
The structure of linear fractional transformations is highly expressive. By using large enough coefficients
⃗a, the algorithm can always find a formula that matches any constant to any continued fraction up to any
finite precision. We should thus distinguish between true results and false-positives. A true result implies
that the formula will remain accurate when extended to more digits. A false-positive result arises from
over-fitting. For example, the formula
314157+e
100000
= π is correct up to 10
−5
, but is obviously false.
To distinguish over-fitted results, we consider that a decimal approximation using K digits can represent
10
K
different numbers. This amount sets a limit on the size of the PSLQ coefficients ⃗a. If the total number
of digits in the ⃗a found by PSLQ is also K, it indicates that the digits alone can already express 10
K
different numbers. Thus, the match is most-likely a coincidence, indicating a false-positive. Conversely, if
the collective number of digits in ⃗a is significantly smaller then K, it means that the result is more likely
true, as the probability for an accidental match is very small. Specifically, to identify false results, we use
the difference between the number of digits given by the decimal approximation K and the total number of
digits l in the vector ⃗a of PSLQ. A larger value of K − l implies high confidence in the formula (with the
chance of an accidental match scaling as ∼ 10
l−k
), while K ≈ l or K < l lead to disqualifying the conjecture.
This condition can be explained intuitively as expressing the same information of K digits using less
data. Unlike the decimal approximation, which conveys information solely through its digits, the formula
uses the mathematical constant for extra representation capability (added to the digits of ⃗a that define the
transformation). The mathematical constant thus enables a more efficient representation of the information.
In other words, true results can be understood as a means of data compression of the numerically computed
digits using the mathematical constant to provide a compact form. Consider for example the special case
where the limit of the continued fraction is exactly the constant itself. Then, the linear fractional trans-
formation does not need to provide any digits (⃗a = (0, 1, 1, 0)), hence we say that the representation is
maximally compressed - all the information is in the constant. The smaller the coefficients of the linear
fractional transformation, the closer it is to maximal compression.
C Selected families of continued fractions found by the Distributed
Factorial Reduction algorithm
In this section, we show a collection of infinite families of continued fractions discovered by the Distributed
Factorial Reduction algorithm. Some of the families presented in this section are used in Appendix D to
create conservative matrix fields for the corresponding constants.
24
ZigZagZeta - formulas mixing several values of the Riemann zeta function
The following continued fractions (Table 2) with a
n
of degree 5 and b
n
of degree 10 were found to have
the factorial reduction property. These results led us to discover the ZigZagZeta family, converging to a mix
of integer values of the Riemann zeta function, up to half the degree of b
n
.
a
n
b
n
Formula
n
5
+ (n + 1)
4
(n + 2) −n
9
(n + 1)
1
ζ(5)−ζ(4)+ζ(3)−ζ(2)+1
n
5
+ (n + 1)
4
(n + 3) −n
9
(n + 2)
32
32ζ(5)−48ζ(4)+56ζ(3)−60ζ(2)+61
n
5
+ (n + 1)
4
(n + 4) −n
9
(n + 3)
7776
7776ζ(5)−14256ζ(4)+18360ζ(3)−20700ζ(2)+21317
n
5
+ (n + 1)
4
(n + 1 + α) −n
9
(n + α) f(ζ(5), ζ(4), ζ(3), ζ(2))
n
4
(n + 1) + (n + 1)
5
−n
9
(n + 1)
1
ζ(4)
n
4
(n + 2) + (n + 1)
5
−n
9
(n + 2)
2
ζ(4)+ζ(3)
n
4
(n + 3) + (n + 1)
5
−n
9
(n + 3)
6
2ζ(4)+3ζ(3)+ζ(2)
n
4
(n + α) + (n + 1)
5
−n
9
(n + α) g(ζ(4), ζ(3), ζ(2), α)
n
4
(n + 1) + (n + 1)
4
(n + 2) −n
4
(n + 1) × n
4
(n + 1)
1
ζ(4)−ζ(3)+ζ(2)−1
n
4
(n + 1) + (n + 1)
4
(n + 3) −n
4
(n + 1) × n
4
(n + 2)
16
16ζ(4)−24ζ(3)+28ζ(2)−29
n
4
(n + 2) + (n + 1)
4
(n + 4) −n
4
(n + 2) × n
4
(n + 1)
2
ζ(4)
5
Q
i=1
(n + r
i
) +
5
Q
i=1
(n + 1 + s
i
) −
5
Q
i=1
(n + r
i
)
5
Q
i=1
(n + s
i
) t(ζ(4), ζ(3), ζ(2), ⃗r,⃗s)
2(n
5
+ (n + 1)
4
(n + 1 + 1/2)) −4n
9
(n + 1/2)
3
6
F
5
[
1,1,1,1,1,1
2,2,2,2,5/2;1
]
3(n
5
+ (n + 1)
4
(n + 1 + 1/3)) −9n
9
(n + 1/3)
4
6
F
5
[
1,1,1,1,1,1
2,2,2,2,7/3;1
]
4(n
5
+ (n + 1)
4
(n + 1 + 1/4)) −16n
9
(n + 1/4)
5
6
F
5
[
1,1,1,1,1,1
2,2,2,2,9/4;1
]
C(n
5
+ (n + 1)
4
(n + 1 + 1/C) −C
2
n
9
(n + 1/C)
C+1
6
F
5
[
1,1,1,1,1,1
2,2,2,2,(2C+1)/C;1
]
Table 2: Results for degree 5 polynomial continued fractions and their generalizations to ZigZagZeta families. These
parametric families were manually generalized based on examples discovered by the BOINC community, executing the
Distributed Factorial Reduction algorithm. The function f that appears in the general form is defined recursively (see
below). These results can be converted to infinite sums [39] using Euler’s continued fraction formula and the decomposition
b
n
= h
1
(n) · h
2
(n), a
n
= h
1
(n) + h
2
(n + 1). This method also provides a closed form for the functions g and t. Similar
families exist for other degrees of polynomial continued fractions.
For the case of a
n
= n
s
+ (n + 1)
s−1
(n + 1 + R) and b
n
= −n
2s−1
(n + R) with values s ≥ 2 and R ∈ Q
(Eq. 8), we denote the limit of the continued fraction by 1/
ˆ
ζ(s, R). It was found by Wolfgang Berndt that
for R ∈ N,
ˆ
ζ(s, R) can be calculated by the recursion
ˆ
ζ(s, R) =
ˆ
ζ(s, R − 1) −
1
R
ˆ
ζ(s − 1, R),
using initial conditions
ˆ
ζ(2, R) =
P
∞
j=R+1
1/j
2
and
ˆ
ζ(s, 0) = ζ(s).
25
Infinite family of continued fractions for ζ(2)
The following is an infinite family of continued fractions that converge to a linear fractional transformation
of ζ(2), and serves as a basis for the conservative matrix field of ζ(2) shown in Appendix D.
α a
n
b
n
Formula
1 n
2
+ (n + 1)
2
−n
4
1
ζ(2)
2 n
2
+ (n + 1)
2
+ 2 −n
4
1
2−ζ(2)
3 n
2
+ (n + 1)
2
+ 6 −n
4
2
3−2ζ(2)
4 n
2
+ (n + 1)
2
+ 12 −n
4
18
31−18ζ(2)
...
α n
2
+ (n + 1)
2
+ α(α − 1) −n
4
1
2Φ(−1,2,α)
Table 3: An infinite family of ζ(2) formulas. These formulas for ζ(2) can be expressed through the Lerch zeta function
Φ(z, s, α) =
P
∞
n=0
z
n
(n+α)
s
for z = −1, s = 2, and integer values of α. The same formula also holds for non-integer values
of α. Using non-integer rational values of α yields formulas that involve more constants other than ζ(2).
We observed that the limit of the continued fractions can be expressed in terms of the Lerch zeta
function Φ(−1, 2, α). Substituting non-integer rational values in α gives rise to constants different from ζ(2).
Specifically, this parametric family generalizes several formulas for Catalan’s constant found by the first
algorithm of the Ramanujan Machine project [6,36] and by more recent algorithms [43] that followed on the
Ramanujan Machine work.
The family presented here is remarkably similar to the one of ζ(3) (see Fig. 4), as both can be expressed
through a
n
= n
d
+ (n + 1)
d
+ c(n
d−2
+ (n + 1)
d−2
) and b
n
= −n
2d
for some rational value of c. This
pattern did not give rise to other infinite families of continued fractions for degrees d higher than 3, but it
did produce some sporadic results that we describe in the next subsection.
Results for different values of the Riemann zeta function
In an attempt to generalize the formulas found for ζ(2) (Table 3) and ζ(3) (Fig. 4) to higher ζ values, we
executed the Distributed Factorial Reduction algorithm over a specially constructed sub-space of polynomial
continued fractions. Specifically, we searched over the following parametric family with integer parameters
B, c
d
, c
d−2
, c
d−4
... (down to c
1
or c
0
depending on the parity of d):
a
n
= c
d
σ(n, d) + c
d−2
σ(n, d − 2) + c
d−4
σ(n, d − 4) + . . . b
n
= −Bn
2d
, (20)
defining σ(n, d) = n
d
+ (n + 1)
d
. The resulting formulas do not belong to the infinite families described
in other sections and are instead sporadic formulas (Table 4). It is currently unknown whether any of the
sporadic formulas can be generalized to form a new infinite family.
Infinite family of continued fractions for values of the Polylogarithm function
Attempting to discover new formulas that converge to expressions involving ζ(4), we executed the Dis-
tributed Factorial Reduction algorithm scanning continued fractions of degree 4. While the search did not
yield many formulas explicitly including ζ(4), numerous formulas presenting constants from the same level of
the hierarchy emerged. Among these results, a subset of continued fractions yielded values of the fourth order
Polylogarithm function Li
4
(x). Subsequently, these findings were rigorously proven using Euler’s continued
fraction formula and generalized to any order of the Polylogarithm function. This generalized continued
fraction family is defined through the polynomials:
a
n
= n
d
+ c(n + 1)
d
b
n
= −cn
2d
.
The resulting continued fractions converge to Li
d
(1/c).
26
a
n
b
n
Formula
n
4
+ (n + 1)
4
+ 2(n
2
+ (n + 1)
2
) −n
8
1
ζ(4)+2ζ(2)−8
n
5
+ (n + 1)
5
+ 6(n
3
+ (n + 1)
3
) −n
10
2
2ζ(5)+6ζ(3)−9
n
5
+ (n + 1)
5
+ 6(n
3
+ (n + 1)
3
) − 4(2n + 1) −n
10
2
2ζ(5)−2ζ(3)−1
n
5
+ (n + 1)
5
+ 16(n
3
+ (n + 1)
3
) − 4(2n + 1) −n
10
64
64ζ(5)+176ζ(3)−273
8(n
5
+ (n + 1)
5
) + 20(n
3
+ (n + 1)
3
) − 5(2n + 1) −64n
10
22.8410...
8(n
5
+ (n + 1)
5
) − 15(n
3
+ (n + 1)
3
) + 9(2n + 1) −64n
10
1.20426...
8(n
5
+ (n + 1)
5
) − 12(n
3
+ (n + 1)
3
) + 7(2n + 1) −64n
10
2.45174...
n
7
+ (n + 1)
7
+ 8(n
5
+ (n + 1)
5
) − 8(n
3
+ (n + 1)
3
) + 4(2n + 1) −n
14
1
ζ(7)−4ζ(3)+4
Table 4: Sporadic results for higher values (≥ 4) of the Riemann zeta function. All the continued fractions in this table
exhibit factorial reduction. Lines 5 to 7 are only formula candidates as we did not find a closed-form expression for them,
despite their recursion formulas being strikingly similar to the others. We would expect them to yield a formula involving
constants from the same level of the hierarchy as ζ(5).
Inspired by this parametric family, we found that the ZigZagZeta family can be similarly expanded.
Indeed, adding a parameter c to the general formula from Table 2 creates the following continued fractions
that have the factorial reduction property:
a
n
=
d
Y
i=1
(n + r
i
) + c
d
Y
i=1
(n + 1 + s
i
)
b
n
= −c
d
Y
i=1
(n + r
i
) ×
d
Y
i=1
(n + s
i
).
D Discovered conservative matrix fields
This section presents examples of conservative matrix fields for several important mathematical constants.
The conservative matrix field of e
The Ramanujan Machine [6] generated hundreds of continued fraction formulas for e. Generalizing
these formulas, we constructed parametric families of continued fractions, which led us to find the following
conservative matrix field:
M
X
=
0 x + 1
1 −(x + y + 1)
, M
Y
=
−1 x + 1
1 −(x + y + 2)
This matrix field can be generated through the degree 1 analytical construction of matrix fields presented
in Appendix H by setting ⃗c = (0, −1, −1, 0). By applying rational indents, this conservative matrix field can
also generate sequences that converge to linear fractional transformations of the Euler–Gompertz constant
and to other values of the
1
F
1
hyper-geometric function.
The conservative matrix field of π
Through an identical process to the e matrix field, we created the following matrix field for π:
M
X
=
0 −(2x + 1)x
1 y + 3x + 2
, M
Y
=
y −x −(2x + 1)x
1 2x + 2y + 1
This matrix field was later identified as a member of the degree 1 analytical matrix field presented in
Appendix H by setting ⃗c = (1, 2, 0, 1).
27
The conservative matrix field of ζ(2)
The infinite family of continued fractions shown in Table 3 enabled us to construct a conservative matrix
field that converges to linear fractional transformations of ζ(2). The matrices that define the matrix field
are:
M
X
=
0 −x
2
(x + 1)
2
x
2
+ (x + 1)
2
+ y(y − 1)
M
Y
=
−x
2
+ xy − y
2
/2 −x
2
x
2
x
2
+ xy + y
2
/2
Section 3 discusses the effect of applying a rational shift to the x or y variables. Applying such a shift
in this case leads to a conservative matrix field that converges to expressions involving Catalan’s constant
G and other constants. Naccache et al. [43,44] have recently generated an infinite family of formulas for the
Catalan constant by a new approach following the Ramanujan Machine. The conservative matrix field of
ζ(2) (with its extra dimensions in the next section) seems to generate these infinite families and generalize
them to higher-dimensional parametric spaces.
E Multidimensional conservative matrix field of degree 2
The definition of a conservative matrix field can be modified to allow for grids of dimension d > 2,
specified by d matrices that satisfy the pairwise conservative property:
M
U
i
(u
0
, ..., u
i
, .., u
j
, ..., u
d
)M
U
j
(u
0
, ..., u
i
+ 1, .., u
j
, ..., u
d
) =
M
U
j
(u
0
, ..., u
i
, .., u
j
, ..., u
d
)M
U
i
(u
0
, ..., u
i
, .., u
j
+ 1, ..., u
d
)
(21)
Multidimensional matrix fields support more trajectories, and with them a possibility for faster converging
sequences and larger irrationality measures.
We formulated a set of four matrices that satisfy the generalized conservative property, using second
degree polynomials as matrix entries. Each of these four matrices can be constructed through a single
matrix (with an arbitrary rational C):
M =
C(xy + xz + xw + yz + yw + zw) −1
C
2
xyzw 0
+ Cx
2
· I,
which enables composing the 4D matrix field:
M
X
(x, y, z, w) = M(x, y, z, w)
M
Y
(x, y, z, w) = M(y, z, w, x)
M
Z
(x, y, z, w) = M(z, w, x, y)
M
W
(x, y, z, w) = M(w, x, y, z).
(22)
Below, we examine the case of C = 1.
Surprisingly, this 4D matrix field can generate continued fractions from the ZigZagZeta family (Appendix
C). To achieve this intriguing result, a specific transformation must be applied to the matrices, forcing them
into a form that represents a continued fraction (as presented in Eq. 9). This transformation involves
taking a co-boundary on the matrices, which can be achieved by using a matrix denoted as U(x, y, z, w) and
applying it to the matrices in each direction according to the following rule: M
s
→ U
−1
M
s
U(s → s + 1),
where s represents the direction index, and U is defined as:
U =
1 xy + xz + xw + yz + yw + zw
0 2xyzw
.
This transformation retrieves the ZigZagZeta family of formulas along the X direction:
M
x
=
0 −(x + 1)(x + y)(x + z)(x + w)
1 (x + 1)
2
+ (x + 1)(x + y + z + w) + (yz + zw + wy)
.
Such a matrix is the direct representation of the polynomial continued fraction with a
n
= (x + 1)
2
+ (x +
1)(x + y + z + w) + (yz + zw + wy) and b
n
= −(x + 1)(x + y)(x + z)(x + w), which exactly matches the
ZigZagZeta family shown in Appendix C.
28
F Methods and algorithms for the conservative matrix fields
The conservative matrix fields are novel mathematical structures with various applications, such as the
generation of infinitely many recurrence formulas that correspond to fundamental constants. This new
mathematical structure calls for developing new algorithmic tools tailored to its unique properties. This
section presents ways to extract optimal trajectories in matrix fields and explains how to transform each
trajectory to a continued fraction.
Algorithms to identify optimal trajectories for maximal irrationality measures
Different trajectories on the conservative matrix field offer different values of δ. In some cases, as depicted
in Fig. 6, the trajectory x = y exhibits the largest δ. But this phenomenon is not consistent across
other matrix fields, leading to a new challenge in identifying the optimal trajectory on a matrix field. The
challenge becomes more significant when dealing with higher-dimension matrix fields (Appendix E), or when
attempting to search a large number of matrix fields for trajectories that yield positive δ.
To search for the optimal trajectory, the first algorithm that we apply is a version of gradient decent. We
begin at the origin and calculate δ(x, y) in neighboring locations using Eq. 16. We then choose the step that
maximizes δ(x, y), and repeat the process. As an example, Fig. 6(b) presents the trajectory found by this
algorithm in black. This algorithm can often miss the optimal trajectory when the function δ(x, y) is not
smooth enough. Various generalizations of gradient-descent algorithms can be directly applied to address
this issue and yield better results.
The second algorithm that we apply uses linear least squares (LLS) to find an optimal trajectory. For
any trajectory initiated at coordinate (1, 1), the endpoint after n steps would be a coordinate (x, y) where
x+y −2 = n. We analyze the δ at each coordinate satisfying this condition to determine the best achievable δ
at the nth step. We then repeat the process for increasing n values and use LLS to fit a linear path of optimal
δ. For the case depicted in Fig. 6(b), the trajectory obtained via this method provides a larger irrationality
measure relative to gradient descent methods, reaching δ = −0.4741 for this constant that combines ζ(3)
and π
3
and is not known to be rational or irrational.
Equivalence between continued fractions and conservative matrix fields
The matrix attained by choosing the x = y trajectory is
M
n
(n) = M
x
(n, n)M
y
(n + 1, n).
This matrix does not generally have a continued fraction form (like M
x
in Eq. 9). However, it is always
possible to transform this matrix to take the form of a continued fraction. In fact, there are certain transfor-
mations that convert an arbitrary sequence of matrix multiplications to a continued fraction form [39] [52].
i.e., it is always possible to extract effective polynomials a
n
, b
n
for any sequence of 2 × 2 matrix multipli-
cations. For example, calculating the conservative matrix field of π along the trajectory x = y creates the
following formula for π:
10
π −4
= −12 +
− 1
−238 +
− 16
−968 +
− 81
.
.
.
+
− n
2
(2n + 1)
2
(4n − 3)(4n + 5)
−2(4n + 3)(6n
2
+ 9n + 2) +
.
.
.
(23)
The conversion from the polynomial continued fraction to a matrix multiplication changes the initial
conditions by a linear fractional transformation with integer coefficients, which does not affect the properties
of the sequence. This behavior is part of a more general property of conservative matrix fields. We can always
change the initial condition of a trajectory (the position or initial matrix) and the resulting constant only
changes by a linear fractional transformation with integer coefficients. We refer to such cases as converging
to the “same” constant.
An implementation of this algorithm can be found in our git repository https://github.com/RamanujanMachine/
ResearchTools
29
G Using families of polynomial continued fractions to derive non-
trivial irrationality measures for ζ(5)
Figure 7: Example delta approximations for
ζ(5), acquired by taking different trajectories
on a 2D space of ZigZagZeta continued frac-
tion formulas.
Table 2 presented several ZigZagZeta families that converge
to expressions involving multiple constants. In this section, we
combine these formulas to generate sequences converging to a sin-
gle constant, by linear combinations that eliminate unwanted con-
stants.
To exemplify this concept, consider a special case of Eq. 8,
where we substitute R = 1. Two members of this family are:
PCF[s = 4, R = 1] =
1
ˆ
ζ(4, 1)
=
1
ζ(4) − ζ(3) + ζ(2) − 1
,
PCF[s = 5, R = 1] =
1
ˆ
ζ(5, 1)
=
1
ζ(5) − ζ(4) + ζ(3) − ζ(2) + 1
.
Summing the inverse of each of these polynomial continued frac-
tions yields a sequence that converges to ζ(5):
ζ(5) =
1
PCF[s = 5, R = 1]
+
1
PCF[s = 4, R = 1]
.
By taking the convergents of the continued fractions, this formula
produces a rational sequence converging to ζ(5).
We can directly generalize this approach by taking linear combinations of formulas with additional values
of R, producing a multi-dimensional infinite family of formulas of any value of the Riemann zeta function.
For example, we take the following formulas with the same R for all powers s up to 5:
1
PCF[s = 2, R]
=
ˆ
ζ(2, R) = ζ(2) + α
2
1
PCF[s = 3, R]
=
ˆ
ζ(3, R) = ζ(3) + β
3
ζ(2) + α
3
1
PCF[s = 4, R]
=
ˆ
ζ(4, R) = ζ(4) + γ
4
ζ(3) + β
4
ζ(2) + α
4
1
PCF[s = 5, R]
=
ˆ
ζ(5, R) = ζ(5) + δ
5
ζ(4) + γ
5
ζ(3) + β
5
ζ(2) + α
5
.
(24)
These linear equations form a matrix that can be inverted to find the coefficients c
i
, providing an explicit
expression for ζ(5) = c
2
·
ˆ
ζ(2, R) + c
3
·
ˆ
ζ(3, R) + c
4
·
ˆ
ζ(4, R) + c
5
·
ˆ
ζ(5, R).
Substituting the corresponding continued fraction into each
ˆ
ζ yields an infinite family of rational sequences
that all converge to ζ(5). We place these sequences as lines of a 2D grid and search for optimal sequences
along different trajectories on it, assessing the irrationality measure provided by each trajectory. Following
this procedure for the specific 2D grid presented here already provides a non-trivial irrationality measure
for ζ(5) of δ ≈ −0.717. Applying this procedure on the much richer multi-dimensional infinite family of
formulas that we found could help optimize δ, potentially increasing it to above 0, and thus providing a path
for proving the irrationality of ζ(5).
30
H Analytical construction of conservative matrix fields
In other sections of this paper (Section 3, Appendix D) we have derived matrices M
X
and M
Y
that
define a conservative matrix field by generalizing an infinite family of continued fractions. In this section,
we present another approach, which involves finding a general expression for the elements of M
X
and M
Y
that satisfy the conservative property. A detailed review of this process is presented in [39]. The outcome of
this analysis is the realization that an infinite family of matrix fields can be constructed using two functions,
f(x, y) and
¯
f(x, y), that satisfy the following two conditions:
• Linear condition: f(x, y) − f(x + 1, y − 1) =
¯
f(x + 1, y) −
¯
f(x, y − 1)
• Quadratic condition: (f
¯
f)(x, y) + (f
¯
f)(0, 0) = (f
¯
f)(x, 0) + (f
¯
f)(0, y)
Note that the linear condition is a modified version of the defining property of Wilf–Zeilberger pairs [34].
Given such f(x, y) and
¯
f(x, y), we define:
a(x, y) = f(x, y) −
¯
f(x + 1, y),
b(x) = (f
¯
f)(x, 0) − (f
¯
f)(0, 0).
We can then obtain the matrices M
X
and M
Y
as follows:
M
X
(x, y) =
0 b(x)
1 a(x, y)
, M
Y
(x, y) =
¯
f(x, y) b(x)
1 f (x, y)
(25)
By setting f (x, y) and
¯
f(x, y) to be polynomials of the same degree, we were able to find closed-form
solutions for the linear and quadratic conditions for f,
¯
f of degrees 1 and 2. The solution takes the form of a
parametric family of conservative matrix fields. For f,
¯
f of degree 1, each conservative matrix field is defined
by four parameters ⃗c = (c
0
, c
1
, c
2
, c
3
), and the solutions are of the form:
f(x, y) = c
0
+ c
1
(x + y)
¯
f(x, y) = c
2
+ c
3
(x − y)
(26)
The horizontal trajectories in this matrix field correspond to polynomial continued fractions with partial
numerator b
n
and denominator a
n
of degrees at most 2 and 1, respectively. For such continued fractions,
the limit can be computed in terms of the hypergeometric functions
1
F
1
and
2
F
1
in [25, §49].
In particular, the degree 1 matrix field can produce exponential, logarithms, and algebraic numbers
as their limit, among other values. We have identified families of ⃗c parameters that create matrix fields
converging to predictable limits:
⃗c = (0, k, l, 0) → e
l
k
⃗c = (0, k, 0, m) → ln
k + m
k
⃗c = (0, k, l, m) →
k + m
k
l
m
These formulas have been empirically verified (each tested numerically along the x = y diagonal). The
obtained constants coincide with the limits along horizontal directions of each of the corresponding matrix
fields, which can be computed in terms of the hypergeometric functions
1
F
1
and
2
F
1
(see [25, §48, 49]).
Interestingly, the parametrization of Eq. 26 encapsulates various properties of conservative matrix fields,
such as the shifts to the initial location of different trajectories as described in Sec. 3. The parameters c
0
and c
2
can be thought of as shifts along the x + y and x − y respectively, by
c
0
c
1
and
c
2
c
3
(assuming c
1
̸= 0 or
c
3
̸= 0). Since we take the limit along the x = y axis, a shift by an integer along it is equivalent to changing
the starting point of the trajectory, which only introduces a linear fractional transform on the limit, but
does not change the fundamental constant underlying the limit.
31
One consequence of the above is that the following matrix fields converge to the same constant (up to a
linear fractional transform), for any integer N:
(c
0
, c
1
, c
2
, c
3
) ⇐⇒ (c
0
+ N, c
1
, c
1
, c
2
, c
3
)
An example of this effect can be seen in the entry ⃗c = (2, 2, 0, 1) in Table 5, which converges to a linear
fractional transform of ln
3
2
, the same constant as in the equation above for k = 2, m = 1, i.e., ⃗c = (0, 2, 0, 1).
Another example is the table entry ⃗c = (3, 3, 2, 4), which converges to a linear fractional transform of
√
21,
the same constant as in the equation above for k = 3, l = 2, m = 4, i.e., ⃗c = (0, 3, 2, 4).
The degree 1 matrix field captures additional cases that we discovered using algorithmic approaches
(Appendix D). For example, setting c
1
= (1, 2, 0, 1) and transforming the matrix to its balanced form ( [53])
results in the conservative matrix field of π.
For f,
¯
f of degree 2, each conservative matrix field is defined by four parameters (c
0
, c
1
, c
2
, c
3
), and the
solutions are of the form:
f(x, y) = (2c
1
+ c
2
)(c
1
+ c
2
) − c
3
c
0
− c
3
((2c
1
+ c
2
)(x + y) + (c
1
+ c
2
)(2x + y)) + c
2
3
(2x
2
+ 2xy + y
2
),
¯
f(x, y) = c
3
(c
0
+ c
2
x + c
1
y −c
3
(2x
2
− 2xy + y
2
)).
One special case is found by setting ⃗c = (0, 0, 0, 1), yielding the conservative matrix field for ζ(2) shown in
Appendix D.
As the degree of f and
¯
f increases, the linear and quadratic conditions evolve into a larger number of
more complicated conditions. We leave this challenge to future works, and present below an instructive
special case. The degenerate case in which f (0, 0) =
¯
f(0, 0) = 0 includes exactly three families of f and
¯
f
polynomials of degree 3. Each of these families is defined by two parameters (c
0
, c
1
):
• f
1
(x, y) = −((c
0
+ c
1
(x + y))(c
0
(x + 2y) + c
1
(x
2
+ xy + y
2
)))
¯
f
1
(x, y) = (c
0
+ c
1
(−x + y))(c
0
(x −2y) −
c
1
(x
2
− xy + y
2
))
• f
2
(x, y) = (c
0
+ c
1
(−x + y))(c
0
(x − 2y) − c
1
(x
2
− xy + y
2
))
¯
f
2
(x, y) = (c
0
+ c
1
(−x + y))(c
0
(x − 2y) −
c
1
(x
2
− xy + y
2
))
• f
3
(x, y) = (x + y)(c
2
0
−c
0
c
1
(x −y) −2c
2
1
(x
2
+ xy + y
2
))
¯
f
3
(x, y) = (c
0
+ c
1
(x −y))(3c
0
(x −y) + 2c
1
(x
2
−
xy + y
2
))
Choosing f
1
,
¯
f
1
and c = (0, 1) results in the conservative matrix field for ζ(3) (Fig. 4, Eqs. 11,13).
In Appendix I, we evaluate irrationality measures of the approximations produced by diagonal trajectories
inside these families of conservative matrix fields, identifying ones that can be used to prove the irrationality
of certain constants.
32
I Using conservative matrix fields to derive irrationality measures
of different constants
The family of conservative matrix fields presented in Appendix H generates infinitely many options to
extract rational sequences with accelerated convergence that can each help proving the irrationality of a
certain constant. Different members of the family converge to different constants, providing a systematic
methodology of creating potentially irrationality-proving sequences to different constants.
By calculating sequences generated from the x = y trajectory of each conservative matrix field, we
constructed dozens of irrationality-proving sequences. Table 5 shows a set of such sequences and their
respective irrationality measures. We note that some of the constants presented here are quite exotic and
hard to fit with PSLQ. Thus, we also used inverse-calculators as Wolfram’s closed-forms detection to connect
these sequences to their matching limits.
Constant CMF degree
CMF parameters
(c
0
, c
1
, c
2
, c
3
)
Numerically
optimized δ
4
1/3
1 (0, 1, 1, 3) 0.0938
14
1/3
1 (0, 4, 4, 3) 0.376
2
1/3
1 (3, 3, 2, 3) 0.321
5
1/3
1 (0, 5, 1, 3) 0.405
e 1 (0, 2, 2, 0) 1.0
√
e 1 (0, 4, 2, 0) 1.0
e
2
1 (0, 1, 2, 0) 1.0
ln(2) 1 (0, 1, 0, 1) 0.288
√
2 1 (0, 4, 2, 4) 1.0
√
3 1 (0, 1, 1, 2) 1.0
√
5 1 (0, 5, 2, 4) 1.0
√
6 1 (0, 4, 1, 2) 1.0
1
4
−5
√
21 +
q
2(21 + 5
√
21
1 (0, 3, 1, 4) 0.0137
9
4
√
3 · 7
3
− 3
1 (0, 3, 3, 4) 0.0098
1
3
−5 − 2
3
√
14 +
3
√
14
2
1 (0, 4, 1, 3) 0.3621
5 − 12 log 3 + log 4096
log 27/8 − 1
1 (2, 2, 0, 1) 0.4041
1
25
21
√
21 − 19
1 (3, 3, 2, 4) 0.9830
ζ(2) 2 (0, 0, 0, 1) 0.0988
ζ(3) 3 (0, 1) using f
1
,
¯
f
1
0.0833
Table 5: Irrationality measures for various constants, using the x = y trajectory in different conservative matrix fields.
33
Because of the importance of this problem, the Clay Institute a million dollar to anyone who would shed light on this question
https://www.claymath.org/millennium/riemann-hypothesis/
For a pretty robust list of mathematical constants, see this Wikipedia link
https://www.wikiwand.com/en/List_of_mathematical_constants
Apéry was a French mathematician
https://www.wikiwand.com/en/Roger_Ap%C3%A9ry
who proved that the Riemann zeta function $\zeta(s)$ evaluated at $s=3$ is irrational.
Apéry's proof is available here
https://www.wikiwand.com/en/Ap%C3%A9ry%27s_theorem
and was simplified over the years. After he proved his theorem, the constant $\zeta(3)$ was named after him - Apéry's constant, and it has a bunch of cool properties https://www.wikiwand.com/en/Ap%C3%A9ry%27s_constant.
The Basel problem is asking to calculate the value of the Riemann Zeta function values at $s=2$, i.e. calculating the numerical value of $\zeta(2)$. It was solved by Euler when he was 28 at the time. He became well known due to this result.
https://www.wikiwand.com/en/Basel_problem
A rational number is a number than can be written as a quotient of two integers, i.e. in the form $\frac{p}{q}$ for $p$ and $q\neq 0$ integers.
The Catalan constant $G$ is defined as $G=\frac{1}{1^2}-\frac{1}{3^2}+\frac{1}{5^2}-\frac{1}{7^2}+\dots$. It is not known whether it irrational or transcendental. The Catalan constant $G$ has been called "arguably the most basic constant whose irrationality and transcendence remain unproven." https://www.wikiwand.com/en/Catalan%27s_constant
The irrationality measure, also called 'Liouville number' is given by bounding a real number $x$ with this inequality
$0<|x-\frac{p}{q}|<\frac{1}{q^n}$
Read more here
https://www.wikiwand.com/en/Liouville_number
The Euler-Mascheroni constant $\gamma$ is
$\gamma = \int_1 ^\infty \left(-\frac{1}{x}+\frac{1}{\lfloor x \rfloor}\right)dx$
and it arises in domino tilings, spanning tress, grid graphs, the asymptotic of the prime numbers, mass distribution of spiral galaxies, and more.
https://www.wikiwand.com/en/Catalan%27s_constant
You gotta love Euler =)
The Gompertz constant (or Euler-Gompertz constant) is
$\delta=\frac{1}{2-\frac{1}{4-\frac{4}{6-\frac{9}{8-\dots}}}}$
It can be represented more simply as an integral $\delta=\int_0 ^\infty \ln(1+x)e^{-x}dx$
https://www.wikiwand.com/en/Gompertz_constant
"Experimental mathematics is that branch of mathematics that concerns itself ultimately with the codification and transmission of insights within the mathematical community through the use of experimental (in either the Galilean, Baconian, Aristotelian or Kantian sense) exploration of conjectures and more informal beliefs and a careful analysis of the data acquired in this pursuit."
(definition taken from https://www.wikiwand.com/en/Experimental_mathematics)
Polynomial continued fractions have a remarkable property of encompassing a wide variety of mathematical functions beyond just rational functions. They provide representations for trigonometric functions, Bessel functions used in physics and engineering, gamma functions with applications in complex analysis and statistics, and hypergeometric functions that arise in differential equations and special functions theory.
Continued fractions have a rich mathematical history, with ancient roots dating back to ancient Greek mathematicians like Euclid. The Indian mathematician Srinivasa Ramanujan, after whom the Ramanujan Machine project was named, made notable contributions to the field of continued fractions, particularly in the area of mock theta functions, which have connections to modular forms and number theory.
The Ramanujan Machine project, aside from being a tribute to Ramanujan's legacy, reflects the growing role of computational methods in mathematical research. The use of algorithms to discover new formulas for constants showcases the synergy between human intuition and computational power.
This was very unexpected to the team who worked on the Ramanujan Machine.
BOINC is a platform that allows volunteers' computers to contribute to large-scale computations when the computers are idle. This concept has been used for projects like the search for extraterrestrial intelligence (SETI) and protein folding simulations.
https://www.seti.org/
The Lerch zeta function is a generalization of the Riemann Zeta function, see the definition here
https://www.wikiwand.com/en/Lerch_zeta_function
The analogy between conservative matrix fields and conservative vector fields is the reason for the term "conservative matrix field". It'd be intriguing to see conservative matrix fields also appearing in the context of theoretical physics.
In the original Ramanujan Machine paper, we aspired to have a future generation of the algorithm that can run in a distributed manner over many computers without redundancy. This was not possible until the discovery of factorial reduction in this work.
The Golden ratio is an important number in fields like art, architecture, and even biology. In math two quantities are in the golden ratio if their ratio is the same as the ratio of their sum to the larger of the two quantities.
https://www.wikiwand.com/en/Golden_ratio
The Golden ratio is ubiquitous in nature
https://www.wikiwand.com/en/Golden_ratio#/Applications_and_observations
The PSLQ algorithm was selected as one of the top 10 algorithms of the century in 2000, and can be used to solve integer equations. Read more here.
https://www.wikiwand.com/en/Integer_relation_algorithm
These are good problems to have. If anything - this is a mark of the broad success of the Factorial Reduction heuristic.
The process performed here is analogous to defining a potential on a conservative vector field.
This information-theoretical approach in a sense attempts to 'quantify beauty' - which is pretty cool in its own right.
A conservative matrix field can also be viewed as a method for accelerating sequences. In the figure provided, the convergence rate increases as the parameter α grows. Ideally, to obtain a sequence with the optimal convergence rate, we would choose the entry corresponding to '$\alpha=\infty$', but this isn't practically feasible. To leverage the increasing convergence rate, we select terms from different sequences. At each step, we sample from a sequence that converges faster. With careful selection, we can construct a new sequence that converges to the same number as the original set of sequences, but does so much faster. In general, faster sequences tend to provide a larger irrationality measure.
It also provides a natural tool to measure the complexity and irrationality of numbers by the complexity of their matrix field.
Personally, the vision of man and AI, working together - each contributing their own ideas and expertise, makes me optimistic about the future.
The current algorithm-assisted research is one of the precursers of human-AI collaboration, and we intend to push the envelop further - applying AI to higher levels on conjecture generation and hypothesization.
Traditionally, mathematics is approached from a 'top-down' perspective. This means starting with a formal definition and progressing through multiple steps and manipulations to gain novel insights. In contrast, experimental science adopts a 'bottom-up' approach. Scientists gather observations through laboratory experiments, from which they detect underlying patterns. These patterns are then leveraged to predict the outcomes of future experiments.
Since our "mathematical experiments" are being conducted by computers, the capacity for collecting observations on a very large scale has opened new avenues for discovering intricate and complex patterns.
Since practically all interesting constants and closed-form functions have continued fractions which exhibit Factorial Reduction, we hypothesize that it's a fundamental property of numbers - meaning that any continued fraction that exhibits Factorial Reduction is interesting/important in some deep way, even if we didn't recognize all of them yet.
Additional research is required.
The quality of PSLQ approximation is related (among other things) to rational approximations of real numbers.
If true, this would mean that each constant has a basic property - previously unknown to mathematicians: its own conservative matrix field.
Assuming this discovery is as broad as we believe, one might compare it to the discovery of the spin of elementary particles: a new fundamental property.
Numerical validation challenges are in principle "just a technical issue". We could just run the calculations longer and get more precision - right?
A - yes, but we have many results exhibiting this issue, so the cost quickly becomes prohibitive.
B - using patterns and the discovered internal structure of continued fraction families is an approach that 'talks in the native language of the numbers'. One might say that we actually want to verify the whole family, and using single continued fraction validations just as a method for achieving the deeper goal.
It should be noted that mathematical discoveries and developments often fall under the bottom-up definition, contrary to tradition. An example of this is complex numbers, which were first used to express the roots of general cubic equations. https://www.wikiwand.com/en/Complex_number#History
The hypergeometric function $_2 F_1$ generalizes many functions, such as the natural logarithm and the geometric sum. It is a part of a larger family of generalized hypergeometric functions $_p F_q$.
https://www.wikiwand.com/en/Hypergeometric_function
https://www.wikiwand.com/en/Generalized_hypergeometric_function
The (standard) Bessel functions, denoted $J_\alpha, Y_\alpha$, are the solutions to the differential equation $x^2 y''+x y'+(x^2-\alpha^2)y=0$. This equation often appears in physics, such as when describing an electromagnetic wave in a cylinder.
https://www.wikiwand.com/en/Bessel_function
Note that this paper defines irrationality measure as $n-1$ instead of $n$. This should explain any "off-by-one" discrepancies between this paper and other sources.
The gamma function $\Gamma(x)$ generalizes the factorial operation $x!$. It is used in physics and in number theory, particularly when describing the Riemann Zeta function.
https://www.wikiwand.com/en/Gamma_function
At first glance, this condition may seem unrelated to fundamental mathematical constants. However, all of the matrices (with polynomial values) we've found that satisfy this condition, are linked to a fundamental constant.
This observation suggests a potential connection between fundamental constants and various symmetries, offering valuable insights into how these constants are intricately tied to different mathematical structures.
The [Liouville–Roth](https://en.wikipedia.org/wiki/Liouville_number) theorem explores the limits of approximating a number using a rational sequence.
For any given denominator $q$, there will always be a numerator $p$ such that:
$$ \frac{p-1}{q} < L < \frac{p}{q}$$
A pair of sequence $p_n$ and $q_n$ that satisfies this equation will achieve $\delta=0$. If we have a _faster_ sequence that approaches $L$, we can get a larger $\delta$.
In practice, finding a formula for $p_n$ and $q_n$ that gives $\delta>-1$ can be an exceptionally challenging task. For certain constants, even marginal enhancements in the best possible irrationality measure hold significant meaning.
This observation is truly remarkable. We are essentially working with the same formula, making only minor adjustments to the initial conditions. The fact that such a minor change yields two sequences that converge to apparently distinct constants is very surprising.
Factorial reduction can play a crucial part in achieving a substantial irrationality measure $\delta$.
When we reduce the common divisor between $p_n$ and $q_n$, the approximation error (left-hand side of Equation 14) does not change. However, we can opt for the _reduced_ denominator $q_n$ in the right-hand side, and thus allow for a larger $\delta$.
In fact, continued fractions like Apery's and other continued fractions that prove irrationality present strong factorial reduction.
Even if we were to discover a continued fraction within the ZigZagZeta family (as shown in Equation 8) that contains $\zeta(5)$ and yields a positive value for $\delta$, this would not proof the irrationality of $\zeta(5)$. This is because values of the continued fractions in the ZigZagZeta family contain other constants we know to be irrational.
To prove the irrationality of $\zeta(5)$, our aim is to find sequences that converge solely to a rational transformation of $\zeta(5)$. Our algorithm achieves this by using other continued fractions from the same family to eliminating unwanted constants and yield an infinite set of sequences that converge exclusively to $\zeta(5).
The plot presented here takes inspiration from a real-world analog, the [bubble chamber](https://en.wikipedia.org/wiki/Bubble_chamber). In a bubble chamber, the trajectory of a particle passing through is influenced by its charge and mass, allowing scientists to identify and distinguish different particles based on their trajectories.
Our 'Mathematical bubble chamber' performs a similar task within the realm of mathematics. In this context, the trajectory of continued fractions passing through our 'Mathematical bubble chamber' is determined by the common divider of their partial numerators and denominators. By analyzing these trajectories, we can effectively identify and differentiate continued fractions that undergo factorial reduction from those that do not