10 1 1 391 1581
10 1 1 391 1581
10 1 1 391 1581
Hieu D. Nguyen
Copyright 2011
Table of Contents
Preface
Chapter 1. The Joy of Number Patterns
Exercises
Preface
Mathematics is not a careful march down a well-cleared highway, but a journey into a strange wilderness, where the
explorers often get lost. Rigour should be a signal to the historian that the maps have been made, and the real
explorers have gone elsewhere.
W. S. Anglin, "Mathematics and History", Mathematical Intelligencer, v. 4, no. 4.
Welcome to the jungle! This book is an attempt to excite the reader about the beauty of mathematics through exploration and
experimentation with number patterns arising from integer sequences. It seeks to demonstrate how these fascinating patterns can
be discovered with the help technology, in particular computer algebra systems such as Stephen Wolfram's Mathematica and
Internet resources such as Neil Sloan's Online Encyclopedia of Integer Sequences (OEIS). Those who find joy in playing with
integer sequences and hunting for formulas experimentally by computer will hopefully find this book worthwhile to read.
More concretely, this book seeks to expose the reader to a wide range of number patterns that can be extracted from experimental
data and to intuitively convince (or deceive at times) the reader of their validity before subjecting them to rigorous proof. This
requires very good bookkeeping and presentation of numerical data. Thus, the author has put considerable effort into grooming
numerical data so that patterns become clear to reader and formulas can be easily conjectured, or if no pattern is evident, then use
tools at our disposable to systematically determine if such a formula exists.
Mathematics by Experiment
Mathematics is not a deductive science -- that's a cliche. When you try to prove a theorem, you don't just list the
hypotheses, and then start to reason. What you do is trial and error, experimentation, guesswork.
Paul R. Halmos, I Want to be a Mathematician, Washington: MAA Spectrum, 1985.
According to the author it seems that little emphasis is placed on exploration and discovery in mathematics education today. No
systematic approach is formally taught to our students, who are typically asked to solve problems that have been so well paved or
narrowly constructed that there is little room for creativity and experimentation. If students are merely taught to prove theorems,
then who will be the explorers of new conjectures?
Mathematical experimentation has largely been overshadowed by mathematical proof. Certainly, developing students’ ability to
understand and formulate mathematical arguments is a worthy goal in any mathematics course; unfortunately, many students
easily get frustrated when they are first exposed to mathematical proofs and asked to construct them independently on their own.
The road to mathematical “nirvana” can be long arduous for these students. However, by instilling in them a deep appreciation
and ownership of the material through exploration and experimentation, the author believes that students will be much more
motivated to tackle mathematical proofs. This issue not only affects college students majoring in mathematics, but more broadly
and more importantly, it also affects liberal arts students and high school students alike.
Experimental Mathematics
Fortunately, there is a new movement to regain the flavor of mathematical exploration and experimentation by computer as
championed by Jonathan Borwein and David Bailey in their series of books on experimental mathematics [References]. Of
course, the intuitive or investigative approach to learning and discovering mathematics is not new. Some middle-school text-
books already emphasize such an approach, for example the Connected Mathematics curriculum, and are divided into units
consisting of “Investigations”. However, the availability of powerful computers, advanced software such as computer algebra
systems, and online mathematical databases to perform numerical experiments and search for deep patterns has given birth to the
digital world of experimental mathematics.
i
Fortunately, there is a new movement to regain the flavor of mathematical exploration and experimentation by computer as
Mathematics by Experiment
championed by Jonathan Borwein and David Bailey in their series of books on experimental mathematics [References]. Of
course, the intuitive or investigative approach to learning and discovering mathematics is not new. Some middle-school text-
books already emphasize such an approach, for example the Connected Mathematics curriculum, and are divided into units
consisting of “Investigations”. However, the availability of powerful computers, advanced software such as computer algebra
systems, and online mathematical databases to perform numerical experiments and search for deep patterns has given birth to the
digital world of experimental mathematics.
It should be emphasized that experimental mathematics is not the same as statistical data analysis or data mining although they do
share some similarities and overlap in certain topics. To this end, Borwein and Bailey give an accurate working definition of
experimental mathematics in [BB]:
Experimental mathematics is the methodology of doing mathematics that includes the use of computations for:
1. Gaining insight and intuition.
2. Discovering new patterns and relationships.
3. Using graphical displays to suggest underlying mathematical principles.
4. Testing and especially falsifying conjectures.
5. Exploring a possible result to see if it is worth formal proof.
6. Suggesting approaches for formal proof.
7. Replacing lengthy hand derivations with computer-based derivations.
8. Confirming analytically derived results.
Thus, experimental mathematics refers to an approach of studying mathematics where a field can be effectively studied using
advanced computing technology such as computer algebra systems. This approach developed and grown so quickly as the
present. Number crunching is what computers do best; therefore, any field that requires significant computation or symbolic
manipulation has or will be touched by experimental mathematics. However, it has not always been this way, as observed by
Borwein and Bailey in [BB]:
One of the greatest ironies of the information technology revolution is that while the computer was
conceived and born in the field of pure mathematics, through the genius of the giants such as John von
Neumann and Alan Turing, until recently this marvelous technology had only a minor impact within the
field that gave it birth.
identify it as such and provide formulas, either as a recurrence Fn+1 = Fn + Fn-1 or explicitly as Fn = (Binet’s
formula), but also to automatically generate other related identities, such as Únk=0 Fk = Fn+2 - 1, Únm=1 Fk2 =
2n 5
Fn Fn+1 , and
Ln = Fn + 2 Fn-1 , where 8Ln < is the Lucas sequence, a cousin of the Fibonacci sequence that satisfies the same recurrence but
begins with the values L0 = 2 and L1 = 1. What is needed then is a study of certain relations common to many integer sequences
and program them for detection.
Thus, the last chapter in this book represents an attempt towards this next step, namely to employ a CAS such as Mathematica to
automatically generate formulas and identities from a given set of integer sequences and to employ databases such as OEIS to
recognize them. We believe the perfect technological storm has arrived to achieve through a convergence of powerful comput-
ers, large online databases, and sophisticated computer algebra systems. These advances will allow us to leapfrog from merely
identifying a single integer sequence to data mining an entire database of integer sequences to reveal all the number patterns and
mathematical identities that lie inside it.
Acknowledgements
In writing this book the author was inspired by Robert M. Young’s classic text, Excusions in Calculus: An Interplay of the
Continuous and the Discrete (MAA, 1992), to foster the same spirit of mathematical exploration by which one is guided by
George Polya’s craft of discovery and the insights of the great masters. Thus, this book seeks by analogy to provide a 21st-
century adaptation of Young’s work by describing the interplay between the analog (human mind) and the digital (computer) in
the exploration of integer sequences and their number patterns.
Jonathan Borwein and David Bailey’s series of books on experimental mathematics (see[BBC], [BB], [BBG], [DB]) greatly
influenced the author’s experimental approach to writing this book.
Neil Sloan’s Online Encyclopedia of Integer Sequences (OEIS) has been an invaluable resource. Its overwhelming collection of
integer sequences has allowed the author to cherry-pick some of the most fascinating integer sequences to share with the reader.
Many of the examples discussed in this book were taken from a variety of expository math journals, such as those published by
the Mathematical Association of America (MAA) and The Fibonaaci Quarterly. Other classic textbooks that have been helpful
and should be mentioned include Concrete Mathematics: A Foundation for Computer Science (Ronald Graham, Donald Knuth,
and Oren Patashnik), The Art of Computer Programming (Donald Knuth), and Unsolved Problems in Number Theory (Richard
Guy).
Lastly, I would like to thank my young daughter, Ai Lan, for her help in typing the title of this book while sitting on her father’s
lap. There are few simple joys in life more memorable than this.
Mathematica Notes
This book was completely typeset using Mathematica 8.0, including all the computations and data that appear in it. Why choose
Mathematica instead of other computer algebra systems (CAS) that are available, such as Maple, Matlab, and Sage? In the
author’s opinion all CAS are essentially comparable in terms of features and their library of mathematical functions. Each has its
advantages and disadvantages, a debate that the author will not touch on, but refer the reader to Internet discussion boards. The
author has chosen Mathematica primarily because it is the CAS of choice at his institution and one that he is most comfortable
iii presented in this textbook can be easily ported over to other
with, having used it for the last fifteen years. Most of the examples
CAS. The point of this book is to demonstrate that modern computer algebra systems such as Mathematica are powerful experi-
mental tools for investigating number patterns of integer sequences.
Mathematics by Experiment
This book was completely typeset using Mathematica 8.0, including all the computations and data that appear in it. Why choose
Mathematica instead of other computer algebra systems (CAS) that are available, such as Maple, Matlab, and Sage? In the
author’s opinion all CAS are essentially comparable in terms of features and their library of mathematical functions. Each has its
advantages and disadvantages, a debate that the author will not touch on, but refer the reader to Internet discussion boards. The
author has chosen Mathematica primarily because it is the CAS of choice at his institution and one that he is most comfortable
with, having used it for the last fifteen years. Most of the examples presented in this textbook can be easily ported over to other
CAS. The point of this book is to demonstrate that modern computer algebra systems such as Mathematica are powerful experi-
mental tools for investigating number patterns of integer sequences.
There are two types of Mathematica commands that the reader should be aware of:
Built-In Commands: These commands are built into Mathematica. Some commands require that the add-on Mathematica
package Combinatorica be loaded. This is automatically done when one loads the Mathematica package MathematicsbyExperi-
ment.m created solely for this book, which can also be downloaded from the book’s website at:
http://www.rowan.edu/colleges/las/departments/math/facultystaff/nguyen/experimentalmath/index.html
It is recommended that you save this file to the same location as the Mathematica notebook containing the online edition of this
book. Then choose ‘Yes’ to evaluate all initialization cells when you first evaluate an input cell in this notebook.
In-House Commands: These commands were programmed by the author. Code for these commands (written as Mathemat-
ica modules) are contained in the same Mathematica package, MathematicsbyExperiment.m. Follow the same instructions as
above.
iv
1
The Joy of Number Patterns
Euler's work reminds us of what we should strive for in the teaching of our students and ourselves--to observe, to
experiment, to search for patterns, and above all, to rejoice not just in the final rigor, but in that wondrous process of
discovery that must always precede it.
Robert M. Young, Excursions in Calculus, 1992.
Blake
Notes on Reynold's Discourses, c. 1808.
Mathematics plays an important role in the sciences, a notion that has been firmly established since the Scientific Revolution.
Much has been written about the impact of mathematics on the sciences. However, the reverse impact has received far less
attention. One aspect that the author believes where modern science can help mathematics is through the former’s most basic
creed, the Scientific Method. Ask a random sample of well educated people to explain the Scientific Method as it was taught to
them school and some might remember it correctly as a systematic approach to performing scientific research. Then ask the same
group if they recall such a method for mathematics. Most likely there will be little or no response. This is because no such
method is taught let alone emphasized throughout their childhood instruction of mathematics. Why is it that scientific experimen-
tation is systemically taught to children in schools at an early age while mathematical experimentation is largely ignored? The
instinct for mathematical exploration remains undeveloped and fails to see light in many children's minds even as they reach
adulthood.
2
Chapter 1
Scientific Method
Make
Observation
Formulate
Try Again
Hypothesis
Conduct
Experiment
Accept Reject
Hypothesis Hypothesis
Mathematical Method
Make
Observation
Formulate
Try Again
Conjecture
Prove
Conjecture
Proof Counter-Example
Found Found
The emphasis on rigorous proof is what sets mathematics apart from the sciences, the gold standard by which every theorem is
measured by. While mathematical proof has received much attention in undergraduate education, mathematical experimentation
has been largely ignored in all mathematics courses. Math instructors spend little time on the discovery aspects of the material
that they teach; most of their time is spent on proof demonstrations.
3 While the author strongly believes in the importance of
proof, he also believes in the power of discovery to motivate students.
Mathematics by Experiment
The emphasis on rigorous proof is what sets mathematics apart from the sciences, the gold standard by which every theorem is
measured by. While mathematical proof has received much attention in undergraduate education, mathematical experimentation
has been largely ignored in all mathematics courses. Math instructors spend little time on the discovery aspects of the material
that they teach; most of their time is spent on proof demonstrations. While the author strongly believes in the importance of
proof, he also believes in the power of discovery to motivate students.
The real danger is not that computers will begin to think like men, but that men will begin to think like computers.
Harris, Sydney J. (In H. Eves, Return to Mathematical Circles, Boston: Prindle, Weber and Schmidt, 1988)
4
Chapter 1
A master of inductive research in mathematics, he made important discoveries (on infinite series, in the
Theory of Numbers, and in other branches of mathematics) by induction, that is, by observation, daring
guess, and shrewd verification. In this respect, however, Euler is not unique; other mathematicians,
great and small, used induction extensively in their work.
Yet Euler seems to me almost unique in one aspect: he takes pains to present the relevant
inductive evidence carefully, in detail, in good order. He presents it convincingly but honestly, as a
genuine scientist should do. He presentation is "the candid exposition" of the ideas that led him to those
discoveries" and has a distinctive charm. Naturally enough, as any other author, he tries to impress his
readers, but, as a really good author, he tries to impress his readers only by such things as have
genuinely impressed himself.
In other words Euler had talent for good bookkeeping. His discovery of the recurrence for the sum of divisors functions and its
connection to the Pentagonal Number Theorem is a classic example of first-rate bookkeeping.
In addition to good bookkeeping, the ability (or patience) to compute mathematical values to great accuracy is another trademark
of great mathematicians. Their computations were amazingly done by hand with mere pencil and paper. Of course these tools
have been replaced by modern calculators and computers, which has allowed even students to compute values quickly and with
extremely high precision.
Each chapter in this book ends with a mathematical story of a great bookkeeper, including the four mentioned above. The author
hopes that the reader will enjoy these stories, revel in their achievements, and be inspired as the author was to explore the
mathematical jungle of integer sequences.
5
Mathematics by Experiment
6
Chapter 1
http://www-history.mcs.st-and.ac.uk/Mathematicians/Euler.html
Leonard Euler is one of the greatest mathematicians of all time and certainly the most prolific. His works span over 50 years and
are accessible online at the Euler Archive (http://www.math.dartmouth.edu/~euler/). Among his many achievements, Euler’s
application of the Pentagonal number theorem to the divisor function stands out as a mathematical gem, which he recounts in his
paper Discovery of a most extraordinary law of the numbers concerning the sum of their divisors, written in French and pub-
lished in 1751 (see [Eu]). We retell the story based on Young’s treatment in [Yo], which includes excerpts of Polya’s translation
of Euler’s paper into English (see [Po]).
First studied by Euler, the sum of divisors (or sigma) function is one of the most interesting mathematical objects in number
theory. Define ΣHnL to be sum of the divisors of n, i.e., ΣHnL = Úd n d. For example, ΣH12L = 1 + 2 + 3 + 4 + 6 + 12 = 28. Here is
a table listing values for Σ HnL for n = 1, ..., 100:
Sum of Divisors ΣHnL
n ΣHnL n ΣHnL n ΣHnL n ΣHnL n ΣHnL
1 1 11 12 21 32 31 32 41 42
2 3 12 28 22 36 32 63 42 96
3 4 13 14 23 24 33 48 43 44
4 7 14 24 24 60 34 54 44 84
5 6 15 24 25 31 35 48 45 78
6 12 16 31 26 42 36 91 46 72
7 8 17 18 27 40 37 38 47 48
8 15 18 39 28 56 38 60 48 124
9 13 19 20 29 30 39 56 49 57
10 18 20 42 30 72 40 90 50 93
The tables above show no obvious pattern for the values of ΣHnL, even for n prime, and leads Euler to give the following bleak
diagnosis:
If we examine a little the sequence of these numbers, we are almost driven to despair. We cannot hope to
discover the least order. The irregularity of the primes is so deeply involved in it that we must think it
impossible to disentangle any law governing this sequence, unless we know the sequence of the primes
itself. It could appear even that the sequence before us is still more mysterious than the sequence of the
primes.
However, after such dire words, Euler, like a mathemagician, reveals that he has cracked the code to the sigma function:
Nevertheless, I observed that this sequence is subject to a completely definite law and could even be
regarded as a recurring sequence. This mathematical expression means that each term can be computed
from the foregoing terms, according to an invariable rule.
To imitate Euler’s discovery, we begin by exploring recursive patterns between consecutive elements ΣHnL, ΣHn - 1L, and
ΣHn - 2L:
7
Mathematics by Experiment
Unfortunately, this pattern fails to hold for ΣH5L through ΣH10L. When a pattern fails, it’s useful to determine the way in which it
failed in hopes of adapting the pattern as Euler did. In this case we find that
ΣH5L = 6 = 7 + 4 - 5 = ΣH4L + ΣH3L - 5
Does this new pattern hold for ΣH6L by subtracting ΣH5L + ΣH4L by 6? Unfortunately, no. However, what we find is that
ΣH6L = 12 = 6 + 7 - 1 = ΣH5L + ΣH4L - ΣH1L
where we have interpreted the subtraction of 1 as subtraction by ΣH1L.
What about ΣH7L? In this case we observe that
ΣH7L = 8 = 12 + 6 - 10 = ΣH6L + ΣH5L - 10
How can we interpret the subtraction by 10? It doesn't seem to follow either of the previous patterns of subtraction by 7 or
subtraction by ΣH2L = 3? Aha! What about subtraction by both terms:
ΣH7L = 8 = 12 + 6 - 3 - 7 = ΣH6L + ΣH5L - ΣH2L - 7
Examining further values, we discover that
ΣH8L = 15 = 8 + 12 - 4 - 1 = ΣH7L + ΣH6L - ΣH3L - ΣH1L
Further experimentation shows that the sigma function obeys a remarkable recurrence, which can be inferred from the table
below:
8
Chapter 1
ΣH1L = 1
ΣH2L = ΣH1L + 2
The following is Euler’s remarkable explanation on how he recognized the divisor sequence:
I confess that I did not hit on this discovery by mere chance, but another proposition opened the path to this
beautiful property...In considering the partition of numbers, I examined, a long time ago, the expression
H1 - xL I1 - x2 M I1 - x3 M I1 - x4 M I1 - x5 M I1 - x6 M I1 - x7 M I1 - x8 M × × ×
in which the product is assumed to be infinite. In order to see what kind of series will result, I multiplied
actually a great number of factors and found
1 - x - x2 + x5 + x7 - x12 - x15 + x22 + x26 - x35 - x40 × × ×
The exponents of x are the same which enter into the above formulas; also the signs + and - arise twice in
succession.
This result known as the Pentagonal Number Theorem. Thus, not only did Euler discover a wonderful pattern involving the sum
of divisors function, but made a profound connection with a generating function for the generalized pentagonal numbers defined
pHnL:
9
Mathematics by Experiment
This result known as the Pentagonal Number Theorem. Thus, not only did Euler discover a wonderful pattern involving the sum
of divisors function, but made a profound connection with a generating function for the generalized pentagonal numbers defined
pHnL:
Theorem (Sum of Divisors): Let n be a positive integer and p(K) be the largest generalized pentagonal number less than n.
Then Σ(n) satisfies the recurrence
H-1LdHN-1L2D n, if n = pHNL
where dHnL = : .
0, otherwise
What a remarkable discovery!
10
Chapter 1
11
Mathematics by Experiment
2
Exploring Patterns of Integer Sequences
If mathematics describes an objective world just like physics, there is no reason why inductive methods should not be
applied in mathematics just the same as in physics.
Kurt Godel, Some Basic Theorems on the Foundations, 1951
The Fibonacci sequence 8Fn < = 80, 1, 1, 2, 3, 5, 8, 13, ...< is one of the most recognized integer sequences in the world (recorded
as entry A000045 in the Online Encyclopedia of Integer Sequences (OEIS)). It is named after Fibonacci (also known as
Leonardo of Pisa), who first wrote on such numbers in his most famous work, the Liber Abaci (Book of Calculations), published
in 1202 AD. However, there is evidence that indicates Fibonacci numbers were studied as early as 700 AD by Indian mathemati-
cians (see [Fi]). It is the only sequence to sport its own publication, the 48-year old journal The Fibonacci Quarterly
(http://www.fq.math.ca/), which publishes mathematical results having connections to the Fibonacci sequence.
The Fibonacci sequence appears in the Liber Abaci [Fi] as a solution to a counting problem involving a population of rabbits:
How Many Pairs of Rabbits Are Created by One Pair in One Year.
A certain man had one pair of rabbits in a certain enclosed place, and one wishes to know how many
are created from the pair in one year when it is the nature of them in a single month to bear another
pair, and in the second month those born to bear also. Because the abovewritten pair in the first month
bore, you will double it; there will be two pairs in one month. One of these, namely the first, bears in the
second month, and thus there are in the second month 3 pairs; of these in one month two are pregnant,
and in the third month 2 pairs of rabbits are born, and thus there are 5 pairs in the month; in this month
3 pairs are pregnant, and in the fourth month there are 8 pairs, of which 5 pairs bear another 5 pairs;
these are added to the 8 pairs making 13 pairs in the fifth month; those 5 pairs that are born in this
month do not mate in this month, but another 8 pairs are pregnant, and thus there are in the sixth month
21 pairs; to these are added 13 pairs that are born in the seventh month; there will be 34 pairs in this
month;...there will be 377 pairs, and this many pairs are produced from the abovewritten pair in the
mentioned place at the end of one year.
You can see in the margin [see below] how we operated, namely that we added the first number to the
second, namely the 1 to the 2, and the second to the third, and the third to the fourth, and the fourth to
the fifth, and thus one after another until we added the tenth to the eleventh, namely the 144 to the 233,
and we had the abovewritten sum of rabbits, namely 377, and thus you can in order to find it for an
unending number of months.
beginning 1
first 2
second 3
third 5
fourth 8
fifth 13
sixth 21 12
seventh 34
eighth 55
beginning 1
Chapter 2
first 2
second 3
third 5
fourth 8
fifth 13
sixth 21
seventh 34
eighth 55
ninth 89
tenth 144
eleventh 233
end 377
Thus, Fibonacci recognized the simple recurrence pattern, Fn+1 = Fn + Fn-1 , that governs these numbers. And yet from simple
recurrence, thousands of number patterns, formulas and identities involving the Fibonacci sequence have been discovered,
including the well-known identity discovered by Cassini,
Fn-1 Fn+1 - Fn2 = H-1Ln (2.1)
â
¥ 1 7- 5
= (2.2)
n=0 F2 n 2
Even to this day, the Fibonacci sequence continues to fascinate mathematicians and be a fruitful subject of investigation.
So how doe one join the hunt to experimentally discover patterns of integer sequences such as the Fibonacci sequence? As a
start, given a sequence 8an <, one should generate enough data, say be calculating a sufficient number of terms of the sequence, to
allow a number pattern to emerge. If necessary, massage the data and systematically present it in such a way that either an
explicit formula, recurrence, or algorithm becomes evident. Additional patterns can then be obtained by considering subse-
quences and transformations of 8an <. It is the myriad of such transformations and the creative techniques for exttracting patterns
that will make our exploration interesting.
Thus, this book discusses methods and tools used to transform integer sequences and analyze them in order to obtain interesting
formulas, identities and connections to other sequences.The following algorithm describes the approach that we shall take to
explore integer sequences:
13
Mathematics by Experiment
Number Pattern
Search Algorithm
Generate
Sequence
Generate Apply
Subsequence Transformation
Analyze
Sequence
Pattern No Pattern
Found Found
n Fn
0 0
1 1
2 1
3 2
4 3
5 5
6 8
7 13
8 21
9 34
2. Generate Subsequence and/or Apply Transformation: We apply the partial sums transformation to 8Fn < to obtain a new
sequence 8Sn < defined by
Sn = âFk
n
(2.3)
k=0
14
Chapter 2
nMax = 9;
dataFibonaccinumbersandpartialsums =
Table@8n, Fibonacci@nD, Sum@Fibonacci@kD, 8k, 0, n<D<, 8n, 0, nMax<D;
ColumnDataDisplayBdataFibonaccinumbersandpartialsums,
k=0
n Fn Sn =Únk=0 Fn
0 0 0
1 1 1
2 1 2
3 2 4
4 3 7
5 5 12
6 8 20
7 13 33
8 21 54
9 34 88
Unfortunately, the two columns 8Fn < and 8Sn < do not appear to be correlated if we make a direct comparison between the values
in each row. However, if we realign these two rows by shifting 8Fn < up by two rows, i.e. replacing 8Fn < with 8Fn+2 <, then a
pattern emerges. Can you describe it?
nMax = 9;
dataFibonaccinumbersandpartialsums =
Table@8n, Fibonacci@n + 2D, Sum@Fibonacci@kD, 8k, 0, n<D<, 8n, 0, nMax<D;
ColumnDataDisplayBdataFibonaccinumbersandpartialsums,
k=0
n Fn+2 Sn =Únk=0 Fn
0 1 0
1 2 1
2 3 2
3 5 4
4 8 7
5 13 12
6 21 20
7 34 33
8 55 54
9 89 88
4. Pattern Found or No Pattern Found: In this case, a pattern was found: 8Fn+2 < and 8Sn < differ by 1. Thus, we’ve discovered
the classic identity
âFk = Fn+2 - 1
n
(2.4)
k=0
Of course, we’ve only experimentally verified that this identity is true and only for a finite number of rows. It remains to prove
(2.4) deductively, say by mathematical induction. Those acquainted with techniques of proof should have no difficulty in proving
this identity. For those new or still transitioning to mathematical proofs (see Appendix A for an introduction to techniques of
proof), below is a proof by mathematical induction. 15
Mathematics by Experiment
Of course, we’ve only experimentally verified that this identity is true and only for a finite number of rows. It remains to prove
(2.4) deductively, say by mathematical induction. Those acquainted with techniques of proof should have no difficulty in proving
this identity. For those new or still transitioning to mathematical proofs (see Appendix A for an introduction to techniques of
proof), below is a proof by mathematical induction.
Proof of (2.4): The base case n = 0 is clearly true:
âFk = F0 = 0
0
k=0
F0+2 - 1 = 1 - 1 = 0
As for the inductive step, assume that the n-th case holds. Then merely add Fn+1 to both sides of (2.4) and apply the Fibonacci
recurrence Fn+3 = Fn+2 + Fn+1 to demonstrate that the Hn + 1L-th case is true:
k=0 k=0
This completes the proof and guarantees that (2.4) holds for all n.
A word of caution is in order. The example presented above involved a pattern that was relatively easy to detect. However, not
all patterns will be this easy. Some will force us to jump through many hoops and hurdles before revealing themselves. Of
course, there will be integer sequences where no patterns exist at all (for now). Fortunately, there are tools at our disposable to
increase our chances of success.
2.1.1.1 Mathematica
Mathematica is a powerful computer algebra system capable of symbolic and high-precision computation. Its large library of
built-in functions, including many that are useful for performing numerical experiments, makes it an indispensable tool for the
experimental mathematician.
Below we give a brief description of some Mathematica commands, based on version 8.0, that will be useful for generating
sequence data and detecting their number patterns. Information about a Mathematica command can be displayed by inserting a
question mark (?) before the command and evaluating it as input. Complete documentation on all Mathematica commands is
available through the program’s help menu or online at:
http://reference.wolfram.com/mathematica/guide/Mathematica.html
Those with access to Mathematica will benefit from the Mathematica version of this book, which allows the user to evaluate all
of the Mathematica programs that appear in it. Please download and run the Mathematica package MathematicsbyExperiment.m
to load those commands not built into Mathematica. This package is available at http://www.rowan.edu/math/facultystaff/nguyen/-
experimentalmath/index.html.
16
Chapter 2
Those with access to Mathematica will benefit from the Mathematica version of this book, which allows the user to evaluate all
of the Mathematica programs that appear in it. Please download and run the Mathematica package MathematicsbyExperiment.m
to load those commands not built into Mathematica. This package is available at http://www.rowan.edu/math/facultystaff/nguyen/-
experimentalmath/index.html.
FindInstance
? FindInstance
FindInstance@expr, varsD finds an instance of vars that makes the statement expr be True.
FindInstance@expr, vars, domD finds an instance over the
domain dom. Common choices of dom are Complexes, Reals, Integers and Booleans.
FindInstance@expr, vars, dom, nD finds n instances.
NOTE: FindInstance is not the same command as Solve; the former attempts to find particular solutions to a set of equations
over a specified domain while the latter attempts to find the most general solution.
Any Diophantine equation of the form x2 - n × y2 = 1, where n is not a perfect square, is called Pell’s equation. These equations
were originally studied by ancient Greek and Indian mathematicians who recognized their usefulness in approximating square
roots. In this example we shall investigate positive integer solutions to the most elementary case, x2 - 2 y2 = 1.
NOTE: If n = m2 is a perfect square, then x2 - n × y2 = 1 reduces to finding Pythagorean triples of the form H1, my, xL.
a) Using Mathematica's FindInstance command, we find that there are five solutions over the domain 0 £ x £ 1000 and
0 £ y £ 1000.
Pellsolutions = Sort@ FindInstance@
x ^ 2 - 2 y ^ 2 == 1 && 0 < x < 1000 && 0 < y < 1000, 8x, y<, Integers, 10DD
88x ® 3, y ® 2<, 8x ® 17, y ® 12<, 8x ® 99, y ® 70<, 8x ® 577, y ® 408<<
Let’s display these solutions in column form to help facilitate our recognition of patterns:
Solutions to Pell's Equation x2 -2y2 =1
x y
3 2
17 12
99 70
577 408
Do you recognize a pattern for the solutions Hx, yL? The shrewd pattern hunter might recognize a recurrence formula from just
these four solutions, from which an explicit formula can be obtain by applying the theory of recursive sequences. For those new
to this game, there is fortunately a Mathematica command that can easily find this explicit formula.
FindSequenceFunction
? FindSequenceFunction
NOTE: To increase Mathematica’s chances of finding a formula for a given sequence using FindSequenceFunction, it may be
necessary to input up to ten terms depending on the complexity of the formula (if it exists). Beware however that it may give
incorrect results; some examples of this will be discussed later in this section. Thus, it is good practice to always confirms
formulas generated from this command.
17
Mathematics by Experiment
NOTE: To increase Mathematica’s chances of finding a formula for a given sequence using FindSequenceFunction, it may be
necessary to input up to ten terms depending on the complexity of the formula (if it exists). Beware however that it may give
incorrect results; some examples of this will be discussed later in this section. Thus, it is good practice to always confirms
formulas generated from this command.
Example 2.2
a) Consider the finite sequence 81, 4, 9, 16, 25, 36, 49, 64, 81<. Applying the FindSequenceFunction command to it yields the
formula for the perfect squares, n2 (A000290), as expected.
FindSequenceFunction@81, 4, 9, 16, 25, 36, 49, 64, 81<, nD
n2
b) However, Mathematica does not always produce a formula that one expects from a given pattern. Consider the sequence
consisting of all natural numbers NOT divisible by 3: 8an < = 81, 2, 4, 5, 7, 8, 10, 11, 13, 14, 16, ...< (A001651). Applying the
FindSequenceFunction command to the first ten terms of this sequence yields the formula:
FindSequenceFunction@81, 2, 4, 5, 7, 8, 10, 11, 13, 14, 16<, nD
H- 3 - H- 1Ln + 6 nL
1
4
This formula reveals little of the fact that it excludes multiples of 3, but is in fact correct. Can you prove that it generates an ?
NOTE: Another formula for an is given by
an = n + dHn + 1L 2t - 1
Here, dxt refers to the floor function (also called the greatest integer function), defined to be the greatest integer less than or
equal to x. Can you prove that this formula and the one given by Mathematica are equivalent? HINT: Equate the two formulas
and simplify to obtain
H- 3 - H- 1Ln + 6 nL == n + Floor@Hn + 1L 2D - 1F
1
SimplifyB
4
Unfortunately, it appears that Mathematica did not evaluate our FindSequenceFunction commands; however, this only indicates
that it was not able to find formulas for x and y, most likely because we did not provide Mathematica with enough terms. Thus,
we’ll need to enlarge our solution set, say double the number of solutions to eight. As the output below shows, this requires
enlarging the domain that needs to be search by several orders of magnitude. Fortunately, this poses no difficulty for Mathemat-
ica given the computer powerful of today’s desktop computers.
18
Chapter 2
Feeding these eight solutions into the FindSequenceFunction now yields the following formulas:
dataPellsolutionsx =
Table@Pellsolutionsdouble@@k, 1, 2DD, 8k, 1, Length@PellsolutionsdoubleD<D
dataPellsolutionsy = Table@Pellsolutionsdouble@@k, 2, 2DD,
8k, 1, Length@PellsolutionsdoubleD<D
83, 17, 99, 577, 3363, 19 601, 114 243, 665 857<
Clear@x0, y0D;
x0 = FindSequenceFunction@dataPellsolutionsx, nD
y0 = FindSequenceFunction@dataPellsolutionsy, nD
II3 - 2 2 M + I3 + 2 2M M
1 n n
I4 + 3 2 M II3 - 2 2 M - I3 + 2 2M M
n n
4 I3 + 2 2M
-
This teaches us a lesson: don’t give up at the first try with FindSequenceFunction -- try longer sequences.
It remains to verify that these formulas indeed satisfy the Pell equation x2 - 2 y2 = 1:
Simplify@x0 ^ 2 - 2 y0 ^ 2 1D
True
NOTE: We have not demonstrated that these formulas generate ALL positive integer solutions. For a proof of this fact and a
complete treatment of Pell equations, see [Ba]. The solutions for x and y appear as entries A001541 and A001542 in OEIS,
respectively.
FURTHER EXPLORATION: Find formulas for positive integer solutions to the negative Pell equation x2 - 2 y2 = -1 and
confirm that your formulas are valid.
Consider the following “Coconut” problem that first appeared in The Saturday Evening Post (October 9, 1926), written by Ben
Ames Williams:
Five men and a monkey were shipwrecked on a desert island, and they spent the first day gathering
coconuts for food. Piled them all up together and then went to sleep for the night.
But when they were all asleep one man woke up, and he thought there might be a row about dividing
the coconuts in the morning, so he decided to take his share. So he divided the coconuts into five piles.
He had once coconut left over, and he gave that to the monkey, and he hid his pile and put the rest all
back together.
By and by the next man woke up and did the same thing. And he had one left over, and he gave it to
the monkey. And all five of the men did the same thing, one after the other; each one taking a fifth of
the coconuts in the pile when he woke up, and each one having one left over for the monkey. And in
the morning they divided what coconuts were left, and they came out in five equal shares. Of course,
each one must have known there were coconuts missing; but each one was guilty as the others, so they
didn’t say anything. How many coconuts were there in the beginning?
19
Mathematics by Experiment
By and by the next man woke up and did the same thing. And he had one left over, and he gave it to
the monkey. And all five of the men did the same thing, one after the other; each one taking a fifth of
the coconuts in the pile when he woke up, and each one having one left over for the monkey. And in
the morning they divided what coconuts were left, and they came out in five equal shares. Of course,
each one must have known there were coconuts missing; but each one was guilty as the others, so they
didn’t say anything. How many coconuts were there in the beginning?
An interesting discussion of the solution by Martin Gardner can be found in [Ga], p.3. Here we shall take a more experimental
approach to determine the number of coconuts in the original pile.
Let c be the total number of coconuts and sk be the number of coconuts that the k-th sailor receives from his division of the
cocunts during the night, and r the number of coconuts that each sailor receives after the final division in the morning. The
problem can then be described by the following system of Diophantine equations:
c = 5 s1 + 1
4 s1 = 5 s2 + 1
4 s2 = 5 s3 + 1
(2.5)
4 s3 = 5 s4 + 1
4 s4 = 5 s5 + 1
4 s5 = 5 r + 1
Eliminating the intermediate variables leads to a single Diophantine equation:
Clear@c, sD;
Eliminate@8c 5 s@1D + 1, 4 s@1D == 5 s@2D + 1, 4 s@2D == 5 s@3D + 1, 4 s@3D == 5 s@4D + 1,
4 s@4D == 5 s@5D + 1, 4 s@5D == 5 r + 1<, 8s@1D, s@2D, s@3D, s@4D, s@5D<D
11 529 + 15 625 r 1024 c
Thus, c = 15 621, which of course is an unrealistically large number of coconuts that the sailors would have discovered. To
obtain other solutions, we argue as follows: since c is divided six times into 5 piles, then 56 can be added to any solution to obtain
the next highest solution. It follows that 15621 must be the smallest positive integer solution.
What if we now generalize the problem to n sailors? What is the smallest position integer solution for c? Here is a table listing
the values of r and c for n = 1 up to n = 5:
Solutions to the Coconut Problem
n r c
1 3 21
2 15 121
3 63 621
4 255 3121
5 1023 15 621
- 4 + 51+n
The sequence 81, 21, 121, 621, 3121, 15 621, ...< is entry A164785 in OEIS.
FURTHER EXPLORATION: Suppose m coconuts remain and are given to the monkey (instead of one coconut) each time the
coconuts are divided into five piles by each of the five sailors. How many coconuts were in the original pile? What if there were
n sailors?
20
Chapter 2
a) WARNING: The FindSequenceFunction command does not alway produce the desired answer since it is sensitive to the
number of terms used. Consider the finite sequence 82, 3, 5, 6, 7, 8, 10, 11, 12, 13, 14, 15, 17< consisting of positive integers
that are NOT squares (A000037). Applying the FindSequenceFunction to it yields
FindSequenceFunction@82, 3, 5, 6, 7, 8, 10, 11, 12, 13, 14, 15, 17<, nD
In the output above notice that the entry 14, which appears in the original sequence, is missing and that the entry 16 should not be
included. Here's a formula, n + n+ n , that correctly generates the original finite sequence:
At this point the reader might object by arguing that the FindSequenceFunction command would have produced the formula
n+ n+ n if it had been given more terms, say 20 terms instead of 13? Unfortunately, even then Mathematica fails to
give this formula. Instead, Mathematica gives a formula involving the DifferenceRoot function:
Table@Floor@n + Sqrt@Sqrt@nD + nDD, 8n, 1, 20<D
82, 3, 5, 6, 7, 8, 10, 11, 12, 13, 14, 15, 17, 18, 19, 20, 21, 22, 23, 24<
Clear@formulaD;
formula@n_D =
FindSequenceFunction@Table@Floor@n + Sqrt@Sqrt@nD + nDD, 8n, 1, 20<D, nD
.=, 9- 1 - y
.@n
.D + y
.@1 + n
.D 0, y
.@1D 2,
. . . . . . .
DifferenceRootAFunctionA9y ., n
.@2D 3, y
.@3D 5, y
.@4D 6, y
.@5D 7, y.@6D 8, y .@7D 10, y
.@8D 11,
. . . . . . .
y
.@9D 12, y
.@10D 13, y
.@11D 14, y .@12D 15, y .@13D 17=EE@nD
. . . . .
y
Let’s confirm if this formula is correct. Here’s a list of the first 20 terms generated from it:
Table@formula@nD, 8n, 1, 20<D
82, 3, 5, 6, 7, 8, 10, 11, 12, 13, 14, 15, 17, 18, 19, 20, 21, 22, 23, 24<
From the output above it does seems that Mathematica‘s formula gives the same terms as those generated from n + n+ n ,
but alas, if we compare the two formulas out to 30 terms, then we find that they disagree starting with the 21th term (the value at
this term should be 26, not 25):
Table@8Floor@n + Sqrt@Sqrt@nD + nDD, formula@nD<, 8n, 1, 30<D
882, 2<, 83, 3<, 85, 5<, 86, 6<, 87, 7<, 88, 8<, 810, 10<, 811, 11<, 812, 12<,
813, 13<, 814, 14<, 815, 15<, 817, 17<, 818, 18<, 819, 19<, 820, 20<,
821, 21<, 822, 22<, 823, 23<, 824, 24<, 826, 25<, 827, 26<, 828, 27<,
829, 28<, 830, 29<, 831, 30<, 832, 31<, 833, 32<, 834, 33<, 835, 34<<
21
Mathematics by Experiment
H- 1 + nL2
On the other hand, Mathematica fails to do this for the Fibonacci sequence. For example, inputting the finite sequence
81, 1, 2, 3, 5, 8, 13, 21<, which begins with the term F1 = 1, yields the Fibonacci sequence as expected:
FindSequenceFunction@81, 1, 2, 3, 5, 8, 13, 21<, nD
Fibonacci@nD
However, prepending the element F0 = 0 doesn’t produce the desired shift in the formula, Fn-1 , that we would expect:
FindSequenceFunction@80, 1, 1, 2, 3, 5, 8, 13, 21<, nD
H- Fibonacci@nD + LucasL@nDL
1
2
This is not necessary a bad thing. In this case, equating the two formulas leads to an identity between the Fibonacci sequence and
its cousin, the Lucas sequence Ln (A000032), defined by the same recurrence, Ln+1 = Ln + Ln-1 , but with different initial values,
namely L0 = 2 and L1 = 1:
Fn-1 = HLn - Fn L 2 (2.6)
or equivalently,
Ln = Fn + 2 Fn-1 (2.7)
The next Mathematica command will be helpful in simplifying complicated expressions.
Simplify
? Simplify
Simplify@exprD performs a sequence of algebraic and other transformations on expr, and returns the simplest form it finds.
Simplify@expr, assumD does simplification using assumptions.
Example 2.5
a) Here’s an example where the Simplify command comes to the rescue to simplify a trigonometric identity that Mathematica
does not recognize:
Sum@Cos@n * Pi 100D, 8n, 0, 100<D
I- 1 - 5 M+ I1 - 5 M+ I- 1 + 5 M+ I1 + 5M
1 1 1 1
4 4 4 4
Simplify@%D
0
b) However, Mathematica isn't always saavy when it comes to reducing formulas, even with the Simplify command. For
example, the following double sum reduces to a constant value, namely 1, which we easily discover by computing some test
values. However, Mathematica does not seem to recognize this or at least cannot prove it to be true based on its methods.
22
Chapter 2
b) However, Mathematica isn't always saavy when it comes to reducing formulas, even with the Simplify command. For
example, the following double sum reduces to a constant value, namely 1, which we easily discover by computing some test
values. However, Mathematica does not seem to recognize this or at least cannot prove it to be true based on its methods.
TimeConstrained@Simplify@Sum@H- 1L ^ j Hi ! j !L, 8i, 0, n<, 8j, 0, n - i<DD, 30D
$Aborted
Mathematica is not a perfect CAS as demonstrated with the FindSequenceFunction command In fact, many of its errors have
been documented, most due to differences in assumptions between Mathematica and the user. Here are some examples to keep
in mind:
Example 2.6
The function f HxL = x0 is well defined only for x > 0 and in that case equals 1. Thus, be careful with improper evaluation such as
f H0L:
f@x_D = x ^ 0;
f@0D
1
Indeterminate
Example 2.7
163 Π
Consider Ramanujan’s constant e , which is an almost integer (a value that is very close to an integer).
R = Exp@Pi Sqrt@163DD
N@R, 32D
163 Π
ã
2.6253741264076874399999999999925 ´ 1017
However, if compute its decimal part, R - dRt, using Mathematica’s approximation command, N, but without specifying the n-
digit precision, then we obtain the wrong answer (recall that the decimal part of every real value is bounded between 0 and 1):
N@R - Floor@RDD
- 480.
This serves as a warning that the reader should be careful with approximating extremely large values when Mathematica.
Example 2.8
23
Mathematics by Experiment
Example 2.8
Solve@x ^ 2 + 1 x x + 1 x, xD
88x ® 0<, 8x ® 1<<
Example 2.9
In the evaluation below, Mathematica claims that the two sums are equal, where the first defines the harmonic number
1 1 1
Hn = 1 + 2
+ 3
+ ... + n , i.e., sum of reciprocals of the first n natural numbers (A001008 and A002805), and the second involves
a harmonic sum of the number of divisors of the natural numbers up to n.
Sum@1 k, 8k, 1, n<D
Sum@1 Hk * Length@Divisors@kDDL, 8k, 1, n<D
HarmonicNumber@nD
HarmonicNumber@nD
:1, >
3 11 25 137
, , ,
2 6 12 60
:1, >
5 17 3 8
, , ,
4 12 2 5
The Online Encyclopedia of Integer Sequences (OEIS) is an online database created by Neil Sloan in 1965 that as of January
2011 contains a collection of almost 200,000 integer sequences and allows the user to search whether a given integer sequence
appears in it (see [O]). Each sequence entry in OEIS (or the Encyclopedia as we’ll sometimes refer to it) contains a host of
information about the sequence, including a list of initial terms, formulas, known results, and references. Here’s an example
using an in-house command called OEIS, i.e. a command created by the author and not built into Mathematica (see Mathematics-
byExperiment.m package file for Mathematica code) to query the OEIS website.
OEIS
? OEIS
a) Suppose we searched the Encyclopedia to see if it recognizes the finite sequence 81, 1, 2, 3, 5, 8, 13, 21, 34, 55<.
OEIS@81, 1, 2, 3, 5, 8, 13, 21, 34, 55<D
24
Chapter 2
99A000045,
Fibonacci numbers: FHnL = FHn- 1 L + FHn- 2 L, FH0L = 0, FH 1 L = 1 , FH
2 L = 1 , ... HFormerly M0692 N0256L , +180 2136 =,
80, 1, 1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, 1597,
2584, 4181, 6765, 10 946, 17 711, 28 657, 46 368, 75 025,
121 393, 196 418, 317 811, 514 229, 832 040, 1 346 269, 2 178 309,
3 524 578, 5 702 887, 9 227 465, 14 930 352, 24 157 817, 39 088 169<=
Not only does OEIS recognize this sequence as being part of the Fibonacci sequence (A000045), but also part of 3185 other
sequences. This shows that the Fibonacci numbers, one of the most recognized sequences in the world, arise in many other
patterns. To learn more about the Fibonacci sequence, just click on the link for A000045, which will open up your web browser
to the OEIS web page describing them. To see the other 3185 sequences, click on the button ‘Go to OEIS complete search
results’, which will open up your web browser and show the full search results directly from OEIS.
b) Here's an example of a finite sequence where one would think that because of the extremely large last element, OEIS would
find either an exact match or no match at all:
OEIS@81, 1, 3, 21, 987, 2 178 309, 10 610 209 857 723<, 2D
OEIS Query: 81, 1, 3, 21, 987, 2 178 309, 10 610 209 857 723<
The On-Line Encyclopedia of Integer
Sequences, published electronically at http:oeis.org, 2010
Suprisingly, there are seven match results with the first two shown above. This shows that there integer sequences that OEIS
cannot match exactly because they are either embedded in other integer sequences or have gone through simplifications during
calculation; but even then partial matches can still provide clues to help us find their formulas. This is demonstrated in the next
example.
25
Mathematics by Experiment
Suppose we wish to find an explicit formula for the recursive rational sequence defined by
u@1D = 1;
u@n_D := H3 u@n - 1D + 1L H5 u@n - 1D + 3L
Below is a list of the first twenty values:
datau = Table@u@nD, 8n, 1, 15<D
:1,
1 5 13 17 89 233 305 1597
, , , , , , , ,
2 11 29 38 199 521 682 3571
>
4181 5473 28 657 75 025 98 209 514 229
, , , , ,
9349 12 238 64 079 167 761 219 602 1 149 851
Observe that applying the FindSequenceFunction command to this sequence yields a complicated formula:
1 5 13 17 89 233
FindSequenceFunctionB:1, , , , , , ,
11 29 38 199 521 2
>, nF
305 1597 4181 5473 28 657 75 025 98 209 514 229
, , , , , , ,
682 3571 9349 12 238 64 079 167 761 219 602 1 149 851
- J- 370 248 451 2-1+n + 165 580 141 ´ 2-1+n 5 - 969 323 029 I7 - 3 5M
-1+n
+
Let’s divide the problem by analyzing the numerators, {1,1,5,13,17,89,...}, and denominators, {1,2,11,29,38,199,...}, separately.
Again, we find that the FindSequenceFunction command fails to give an elementary pattern for the numerators, expressing them
in terms of Mathematica’s DifferenceRoot function.
datanumerator = Numerator@datauD
81, 1, 5, 13, 17, 89, 233, 305, 1597, 4181, 5473, 28 657, 75 025, 98 209, 514 229<
FindSequenceFunction@datanumerator, nD
n=, 9y.@n
.D - 18 y
.@3 + .
nD + y
.@6 + .
nD 0,
. . . . . . . .
DifferenceRootAFunctionA9y ., .
.@1D 1, y
.@2D 1, y
.@3D 5, .@4D 13, y .@5D 17, y.@6D 89=EE@nD
. . . . . .
y y
On the other hand, feeding the same sequence into OEIS yields three possible matches involving sequences that are related to the
Fibonacci sequence:
OEIS@datanumerator, 10D
26
Chapter 2
OEIS Query: 81, 1, 5, 13, 17, 89, 233, 305, 1597, 4181, 5473, 28 657, 75 025, 98 209, 514 229<
The On-Line Encyclopedia of Integer
Sequences, published electronically at http:oeis.org, 2010
With this as a clue, the keen bookkeeper will notice that most of the numerators are in fact Fibonacci numbers, namely
81, 1, 5, 13, 89, 233, 1597, 4181, ...<, and in fact those that are not Fibonacci numbers, 817, 305, 5473, ...<, occur at every third
entry starting with 17.
Using the following in-house command, which identifies the positions of elements within a given integer sequence, we verify that
the numerators consist of Fibonnacci numbers of the form F2 n-1 :
MatchPosition
? MatchPosition
MatchPosition@subset,sequence,n,first,lastD displays the positions of the elements specified by subset within the list generated
by the function sequence@nD from n=first to n=last. Requires loading of MathematicsbyExperimentPackage.m package.
27
Mathematics by Experiment
As for the corresponding list of denominators, we follow the trail above and eliminate every third entry starting with the denomina-
tor 38 from consideration. Then matching these denominators with the Lucas numbers yields the same pattern:
datadenominator = Denominator@datauD
81, 2, 11, 29, 38, 199, 521, 682, 3571, 9349, 12 238, 64 079, 167 761, 219 602, 1 149 851<
The formula for un in terms of Fibonacci and Lucas numbers is now clear:
F2 n-1
un =
L2 n-1
where Fn and Ln are Fibonacci and Lucas numbers, respectively. The following table confirms this:
Table@Fibonacci@2 n - 1D LucasL@2 n - 1D, 8n, 1, 10<D
:1, >
1 5 13 17 89 233 305 1597 4181
, , , , , , , ,
2 11 29 38 199 521 682 3571 9349
The Inverse Symbolic Calculator (ISC), created by Simon Plouffe in 1995, is useful for searching real numbers that represent the
exact value of a given approximation (see [I]). Here is an example using the in-house command, ISC, to query the ISC website
(see Mathematics by Experiment package file for Mathematica code) .
28
Chapter 2
The Inverse Symbolic Calculator (ISC), created by Simon Plouffe in 1995, is useful for searching real numbers that represent the
exact value of a given approximation (see [I]). Here is an example using the in-house command, ISC, to query the ISC website
(see Mathematics by Experiment package file for Mathematica code) .
ISC
? ISC
For example, feeding the approximation 3.14159 for Π into ISC yields many matches (289). The first five matches are shown
below:
ISC@"3.14159", 5D
Inverse Symbolic Calculator, published electronically
at http:oldweb.cecm.sfu.caprojectsISCISCmain.html
We can, of course, reduce the number of matches by feeding ISC a better approximation:
ISC@"3.141592653589793", 5D
Inverse Symbolic Calculator, published electronically
at http:oldweb.cecm.sfu.caprojectsISCISCmain.html
Surprisingly, we still found 43 matches, which goes to show that there are many interesting values that approximate Π.
The Pell equation, x2 - 2 y2 = 1 (discussed in Example 1.1), and its negative, x2 - 2 y2 = -1, share a common pattern. Let’s
combine their first four solutions:
29
Mathematics by Experiment
To obtain an even better approximation of the limiting ratio, we can extend the solution set by using the FindSequenceFunction
to find formulas for xn and yn :
x@n_D = FindSequenceFunction@81, 3, 7, 17, 41, 99, 239, 577<, nD
y@x_D = FindSequenceFunction@81, 2, 5, 12, 29, 70, 169, 408<, nD
II1 - 2 M + I1 + 2M M
1 n n
2
Fibonacci@n, 2D
The output above refers to the Fibonacci polynomial Fn H2L. Here is a list of values for rn accurate to 12 decimal places:
30
Chapter 2
It is now clear that these ratios converge to limit 2 . Those that did not recognize this value can feed 1.414213562373 into ISC:
ISC@"1.414213562373", 5D
Inverse Symbolic Calculator, published electronically
at http:oldweb.cecm.sfu.caprojectsISCISCmain.html
This confirms 2 as the limit. Of course we can obtain the same answer by also using Mathematica’s Limit command since
exact formulas for xn and yn are known:
Limit@x@nD y@nD, n ® InfinityD
Suppose we wanted to find the limit of the sequence un in Example 1.10 as n ® ¥. Here is a table of the first 30 values of un :
31
Mathematics by Experiment
n un n un n un
1 1.000000000 11 0.4472361809 21 0.4472135970
2 0.3333333333 12 0.4472049689 22 0.4472135949
3 0.5000000000 13 0.4472168906 23 0.4472135957
4 0.4285714286 14 0.4472123369 24 0.4472135954
5 0.4545454545 15 0.4472140762 25 0.4472135955
6 0.4444444444 16 0.4472134119 26 0.4472135955
7 0.4482758621 17 0.4472136656 27 0.4472135955
8 0.4468085106 18 0.4472135687 28 0.4472135955
9 0.4473684211 19 0.4472136057 29 0.4472135955
10 0.4471544715 20 0.4472135916 30 0.4472135955
This sequence appears to stabilize out to 10 decimal places to the value 0.4472135955. However, feeding this value into ISC
fails to give an answer.
ISC@"0.4472135955", 5D
Inverse Symbolic Calculator, published electronically
at http:oldweb.cecm.sfu.caprojectsISCISCmain.html
This is because Mathematica had rounded off the last significant digit (10-th decimal place). A more precise calculation shows
that this digit should be 4 not 5.
Table@N@Fibonacci@nD LucasL@nD, 15D, 8n, 1, 40<D
81.00000000000000, 0.333333333333333, 0.500000000000000, 0.428571428571429,
0.454545454545455, 0.444444444444444, 0.448275862068966, 0.446808510638298,
0.447368421052632, 0.447154471544715, 0.447236180904523, 0.447204968944099,
0.447216890595010, 0.447212336892052, 0.447214076246334, 0.447213411871319,
0.447213665639877, 0.447213568708896, 0.447213605733234, 0.447213591591195,
0.447213596992973, 0.447213594929677, 0.447213595717786, 0.447213595416755,
0.447213595531739, 0.447213595487819, 0.447213595504595, 0.447213595498187,
0.447213595500634, 0.447213595499700, 0.447213595500057, 0.447213595499920,
0.447213595499972, 0.447213595499952, 0.447213595499960, 0.447213595499957,
0.447213595499958, 0.447213595499958, 0.447213595499958, 0.447213595499958<
1
Feeding the correct value 0.4472135954 into ISC now yields the exact limit, .
5
ISC@"0.4472135954", 5D
32
Chapter 2
NOTE: Of course we can obtain the exact answer by also using Mathematica's Limit command since exact formulas for Fn and
Ln are known:
Limit@Fibonacci@nD LucasL@nD, n ® InfinityD
1
You can use all the quantitative data you can get, but you still have to distrust it and use your own intelligence and
judgment.
-- Alvin Toffler (www.quotationspage.com)
In this modern age and with tools such as computer algebra systems, it is quite easy to generate large data sets for analysis. In
this section we discuss frequently used techniques (tricks of the trade) for manipulating and transforming integer sequences to
obtain number patterns. The examples presented demonstrate that although the computer can be programmed to detect most of
these patterns, it still pays to do good bookkeeping and have a keen eye.
I never guess. It is a capital mistake to theorize before one has data. Insensibly one begins to twist facts to suit theories,
instead of theories to suit facts.
-- Sir Arthur Conan Doyle (1859 - 1930), The Sign of Four, A Scandal in Bohemia
Good bookkeeping begins with arranging data in an orderly fashion that allows patterns to reveal themselves either through
symmetry or periodicity.
Partitioning a sequence into a collection of subsets or a two-dimensional data set is an effective method to reveal patterns that
would be difficult to recognize otherwise.
33
Mathematics by Experiment
Define a sequence 8an <, called the diatomic sequence (A002487), whose generation proceeds as follows:
1. Begin with the finite sequence 80, 1<.
2. Insert between every two consecutive terms of the sequence a term equal to the sum of those two entries
3. Iterate step 2.
The sequence in the limit is the diatomic sequence. The table below describes the first five iterations.
Generation of Diatomic Sequence
80, 1<
Iteration Sequence
80, 1, 1<
0
80, 1, 1, 2, 1<
1
80, 1, 1, 2, 1, 3, 2, 3, 1<
2
80, 1, 1, 2, 1, 3, 2, 3, 1, 4, 3, 5, 2, 5, 3, 4, 1<
3
4
5 80, 1, 1, 2, 1, 3, 2, 3, 1, 4, 3, 5, 2, 5,
3, 4, 1, 5, 4, 7, 3, 8, 5, 7, 2, 7, 5, 8, 3, 7, 4, 5, 1<
NOTE: Observe that this construction has memory in the sense that each iteration preserves the sequence from the previous
iteration.
The sequence 8an < can be programmed in Mathematica as follows:
Clear@a, sternD;
nMax = 11;
a@0D = 0;
a@1D = 1;
Do@
stern =
Riffle@Table@a@kD, 8k, 0, 2 ^ i - 1<D, Table@a@k - 1D + a@kD, 8k, 1, 2 ^ Hi - 1L<DD;
Do@a@2 ^ Hi - 1L + kD = stern@@2 ^ Hi - 1L + k + 1DD, 8k, 1, 2 ^ Hi - 1L<D,
8i, 1, nMax<D;
datastern = Table@8n, a@nD<, 8n, 0, 39<D;
ColumnDataDisplay@datastern, 10,
8"n", "aHnL"<, "Generation of Diatomic Sequence", LeftD
Generation of Diatomic Sequence
n aHnL n aHnL n aHnL n aHnL
0 0 10 3 20 3 30 4
1 1 11 5 21 8 31 5
2 1 12 2 22 5 32 1
3 2 13 5 23 7 33 6
4 1 14 3 24 2 34 5
5 3 15 4 25 7 35 9
6 2 16 1 26 5 36 4
7 3 17 5 27 8 37 11
8 1 18 4 28 3 38 7
9 4 19 7 29 7 39 10
a) Observe that the value 1 seems to repeat itself and occur at positions equal to a power of 2. In other words, we have a2n = 1.
Let’s use these positions as markers then to create an array (ignore the first element a0 ) whose rows begin with 1 and whose
lengths increase by powers of two:
34
Chapter 2
1
1 2
1 3 2 3
1 4 3 5 2 5 3 4
1 5 4 7 3 8 5 7 2 7 5 8 3 7 4 5
1 6 5 9 4 11 7 10 3 11 8 13 5 12 7 9 2 9 7 12 5 13 8 11 3 10
â ak = 3n .
2n+1 -1
(2.8)
k=2n
b) Next, notice that the gaps along each column is constant, i.e. the difference between any two terms in each column is the same.
Gaps Along Each Column
Column Gap
1 0
2 1
3 1
4 2
5 1
6 3
7 2
8 3
9 1
10 4
Do you recognize this sequence of gaps? You should since it is just the diatomic sequence itself. Thus, we’ve found the identity
a2n+1 +n - a2n +n = an (2.9)
c) Also, observe that values in certain columns above repeat themselves. For example, column 2 is repeated as colunns 3, 5, 9,
17, etc. (positions that are one more than a power of 2) as shown in red below:
Grid@Table@If@n 2 ^ p + 1 ÈÈ n == 2 ^ p + 2 ÈÈ n 2 ^ p + 4 ÈÈ n 2 ^ p + 8 ÈÈ n 2 ^ p + 16,
Style@a@nD, RedD, a@nDD, 8p, 0, 5<, 8n, 2 ^ p, 2 ^ Hp + 1L - 1<D,
Frame ® All, Background ® 8None, 88None, LightGray<<<D
1
1 2
1 3 2 3
1 4 3 5 2 5 3 4
1 5 4 7 3 8 5 7 2 7 5 8 3 7 4 5
1 6 5 9 4 11 7 10 3 11 8 13 5 12 7 9 2 9 7 12 5 13 8 11 3 10
35
Mathematics by Experiment
NOTE: It seems that each entry is the sum of two entries above similar to Pascal's triangle; however, this pattern fails for those
entries that are repeated from the previous column. Therefore, if we consider instead a diatomic array where we pad 1's at the
end of each row to make it palindromic and arrange the entries into columns as shown below for the first five rows:
diatomicarray@nMax_D := Table@ReplacePart@
PadRight@8<, 2 ^ HnMax - 1L + 1, ""D, H1 + 2 ^ HnMax - nL * ð ® a@2 ^ Hn - 1L + ðDL &
HRange@0, 2 ^ Hn - 1LDLD, 8n, 1, nMax<D Grid
diatomicarray@5D
1 1
1 2 1
1 3 2 3 1
1 4 3 5 2 5 3 4 1
1 5 4 7 3 8 5 7 2 7 5 8 3 7 4 5 1
Thus we see that each entry is either the sum of the two entries above it or if there is an entry directly above it, then it takes on the
same value as that entry. Recursively, the next row in the array can be constructed by copying the previous row and inserting the
sum of every two consecutive entries in the previous row as a new entry.
A sequence 8an < is recursive if it satisfies a recurrence, a relationship where each element an is expressed in terms of any of the
previous elements an-1 , an-2 , ..., a1 , a0 . A linear homogeneous recurrence with constant coefficients is one of the form
an = c1 × an-1 + c2 × an-2 + ... + cn-1 × a1 + cn × a0
where 8c1 , c2 , ..., cn < are constant real values. This recurrence will be called a 8c1 , c2 , ..., cn <-recurrence. For example, recall
that the Fibonacci sequence 80, 1, 1, 2, 3, 5, 8, 13, ...< has a 81, 1<-recurrence: Fn+1 = Fn + Fn-1 .
Generating Recursive Sequences
In many cases generating integer sequences by recursion is more efficient than using an explicit formula. There are many
methods to generate recursive sequences in Mathematica. We demonstrate four different methods to illustrate their strengths and
weaknesses.
Example: Let’s generate the sequence an defined by a 85, -1<-recurrence:
an = 5 an-1 - an-2
METHOD 1 - We define the sequence an as a delayed function a[n] using the SetDelayed assignment specified by the symbol
“:=” (colon-equal sign):
Clear@aD;
a@0D = 2;
a@1D = 5;
a@n_D := 5 a@n - 1D - a@n - 2D
Timing@Table@a@nD, 8n, 0, 25<DD
83.291, 82, 5, 23, 110, 527, 2525, 12 098, 57 965, 277 727, 1 330 670, 6 375 623, 30 547 445,
146 361 602, 701 260 565, 3 359 941 223, 16 098 445 550, 77 132 286 527, 369 562 987 085,
1 770 682 648 898, 8 483 850 257 405, 40 648 568 638 127, 194 758 992 933 230,
933 146 396 028 023, 4 470 972 987 206 885, 21 421 718 540 006 402, 102 637 619 712 825 125<<
The Timing command shows that it took Mathematica about 3 seconds to generate the first 25 terms, which is quite slow. This is
because previous terms, instead of being stored to memory, are always recalculated (due to the SetDelayed assignment) in
generating the next term of the sequence.
36
Chapter 2
METHOD 2 - To force Mathematica to store previous terms, we adapt our definition of a[n] as follows:
ClearAll@aD;
a@0D = 2;
a@1D = 5;
a@n_D := a@nD = 5 a@n - 1D - a@n - 2D
Timing@Table@a@nD, 8n, 0, 25<DD
80., 82, 5, 23, 110, 527, 2525, 12 098, 57 965, 277 727, 1 330 670, 6 375 623, 30 547 445,
146 361 602, 701 260 565, 3 359 941 223, 16 098 445 550, 77 132 286 527, 369 562 987 085,
1 770 682 648 898, 8 483 850 257 405, 40 648 568 638 127, 194 758 992 933 230,
933 146 396 028 023, 4 470 972 987 206 885, 21 421 718 540 006 402, 102 637 619 712 825 125<<
Observe that this methods drastically reduces the run-time to almost zero. In fact, calculating the first 1000 terms requires less
than a second (output is suppressed):
ClearAll@aD;
a@0D = 2;
a@1D = 5;
a@n_D := a@nD = 5 a@n - 1D - a@n - 2D
Timing@Table@a@nD, 8n, 0, 1000<D;D
80.032, Null<
METHOD 3 - An equivalent approach to Method 2 is to program a loop to generate an inside a Mathematica module (subroutine):
ClearAll@aD;
sequence@nMax_D := Module@8a, n<,
a@0D = 2;
a@1D = 5;
Do@
a@nD = 5 a@n - 1D - a@n - 2D,
8n, 2, nMax<
D;
Table@a@nD, 8n, 0, nMax<D
D
Timing@sequence@1000D;D
80.032, Null<
METHOD 4
One can also use Mathematica’s LinearRecurrence command, which is useful for generating linear recursive sequences.
LinearRecurrence
? LinearRecurrence
37
Mathematics by Experiment
Besides being simple to use, timing shows that LinearRecurrence is twice as fast as Methods 2 and 3 in generating the first
1000 terms.
NOTE: We can of course use the FindSequenceFunction command to find an explicit formula for an based on say the first ten
terms:
Simplify@FindSequenceFunction@LinearRecurrence@85, - 1<, 82, 5<, 10D, nDD
- I- 5 + 21 M I5 + 21 M I5 - 21 M I5 + 21 M
1 1 n 1 n
+
2 2 2
However, we find that this formula is much slower in generating the first 1000 terms:
ClearAll@aD
- J- 5 + 21 N J5 + 21 N + J5 - 21 N J5 + 21 N ;
1 1 n 1 n
a@n_D =
85.163, Null<
This is because considerable computation is required to simplify the binomial expressions for an to an integer value.
Finding Linear Recurrences
If a recurrence is not known for a finite sequence, then the following command attempts to find such a recurrence.
FindLinearRecurrence
? FindLinearRecurrence
FindLinearRecurrence@listD finds if possible the minimal linear recurrence that generates list.
FindLinearRecurrence@list, dD finds if possible the linear recurrence of maximum order d that generates list.
Suppose we wish to find a linear recurrrence for the set of solutions rHnL to the Coconut Problem discussed in Example :
FindLinearRecurrence@821, 121, 621, 3121, 15 621<D
86, - 5<
Thus, rHnL satisfies the recurrence rHnL = 6 rHn - 1L - 5 rHn - 2L, which we now use to generate the next five solutions:
LinearRecurrence@86, 5<, 821, 121<, 10D
821, 121, 831, 5591, 37 701, 254 161, 1 713 471, 11 551 631, 77 877 141, 525 021 001<
The integers, commonly expressed as a sequence which runs in both directions, 8 ..., -3, -2, -1, 0, 1, 2, 3, ...<, satisfy a simple
recurrence: n = Hn - 1L + 1. In this example we will show that rearranging them can produce some rather interesting recurrences
(see [MP]).
a) Let’s rearranging the integers into a sequence that begins with 0 and then alternatees between n and -n:
80, 1, -1, 2, -2, 3, -3, 4, -4, ...<. What recurrence does this sequence satisfy?
FindLinearRecurrence@80, 1, - 1, 2, - 2, 3, - 3, 4, - 4, 5, - 5<D
8- 1, 1, 1<
38
Chapter 2
b) How about the recurring sequence in which every integer appears exactly twice? Let’s try
80, 0, 1, 1, -1, -1, 2, 2, -2, -2, 3, 3, -3, -3, 4, 4, -4, -4, ...<:
FindLinearRecurrence@
80, 0, 1, 1, - 1, - 1, 2, 2, - 2, - 2, 3, 3, - 3, - 3, 4, 4, - 4, - 4, 5, 5, - 5, - 5<D
81, - 2, 2, - 1, 1<
c) Find a formula for the recurrence of a corresponding sequence in which every integer occurs exactly n-times?
NOTE: As its name implies, FindLinearRecurrence is only useful in finding constant-coefficient linear recurrences (of homoge-
neous type). Thus, it is not able to find more general linear recurrences satisfied by sequences such the factorials aHnL = n !,
which has the recurrence
aH1L = 1; aHnL = n × aHn - 1L
or nonlinear recurrences such as aHnL = aHn - 1L aHn - 2L with initial values aH0L = 1 and aH1L = 2.. Try it for yourself to see how
FindLinearRecurrence evaluates these examples.
A classic problem asks for a house number in Louvain, Belgium where the sum of the house numbers before it equals the sum of
the house numbers after it (see [Pa]). It is assumed that house numbers take on integers from 1 to k. For example, if k = 8, then
the desired house number is 6 since
1+2+3+4+5=7+8
However, not all integer values of k will yield an integer solution for h (the reader can easily check this for k = 5). In fact, as we
shall see later, such solutions are sparse.
More generally, given a positive integer k, we seek to find another positive integer h such that
1 + 2 + ... + Hh - 1L = Hh + 1L + Hh + 2L + ... + k (2.11)
Simplifying the equation above in Mathematica yields
Simplify@Sum@i, 8i, 1, h - 1<D Sum@i, 8i, h + 1, k<DD
2 h2 k + k2
Therefore, it suffices to solve for h and check which values of k will yield integers solutions for h. Here is a table listing the first
six solutions where k is less than 10,000:
Solutions to the House Number Problem
h k
1 1
6 8
35 49
204 288
1189 1681
6930 9800
simplify equation (2.11) above. Without these formulas, the run time required to generate the table of solutions above would be
much longer. However, even with these formulas, finding the next five solutions, which grow relatively quickly, will require an
exceedingly long time.
A more efficient approach to finding additional solutions is to find a recursion satisfied by the first six values for h:
39
Mathematics by Experiment
This recursion allows us to now quickly generate the first ten solutions:
LinearRecurrence@86, - 1<, 81, 6<, 10D
81, 6, 35, 204, 1189, 6930, 40 391, 235 416, 1 372 105, 7 997 214<
NOTES:
1. The equation 2 h2 = k 2 + k can be transformed into a Pell equation by completing the square in k and making an appropriate
change of variables:
2 h2 = k 2 + k = Hk + 1 2L2 - 1 4 Þ 8 h2 = H2 k + 1L2 - 1
Setting x = 2 k + 1 and y = 2 h leads to the Pell equation
x2 - 2 y2 = 1
which we discussed in Example (2.1). The solutions for Hx, yL correspond precisely to solutions for Hh, kL as seen in the follow-
ing table:
ClearAll@x, yD; Pellsolutionsplus = FindInstance@
Hx ^ 2 - 2 y ^ 2 == 1L && 0 < x < 1000 && 0 < y < 1000, 8x, y<, Integers, 10D
dataPellsolutionsplus = Table@8Sort@PellsolutionsplusD@@k, 1, 2DD,
Sort@PellsolutionsplusD@@k, 2, 2DD<, 8k, 1, Length@PellsolutionsplusD<D
ColumnDataDisplayATable@8n, dataPellsolutionsplus@@n, 1DD,
dataPellsolutionsplus@@n, 2DD, dataPellsolutionsplus@@n, 2DD 2,
HdataPellsolutionsplus@@n, 1DD - 1L 2<, 8n, 1, Length@dataPellsolutionsplusD<D,
10, 8"n", "xn ", "yn ", "hn =yn 2", "kn =Hxn -1L2"<, "Pell Equations x2 -2y2 =1"E
88x ® 3, y ® 2<, 8x ® 17, y ® 12<, 8x ® 99, y ® 70<, 8x ® 577, y ® 408<<
2. The values for h also correspond to the products xn yn , where Hxn yn L is a solution to the Pell equations x2 - 2 y2 = ± 1:
88x ® 1, y ® 1<, 8x ® 7, y ® 5<, 8x ® 41, y ® 29<, 8x ® 239, y ® 169<,
8x ® 3, y ® 2<, 8x ® 17, y ® 12<, 8x ® 99, y ® 70<, 8x ® 577, y ® 408<<
40
Chapter 2
Can you explain why? In fact, the sequence 81, 6, 35, 204, 1189, 6930, ...< are square-triangular numbers (A001109), which we
discuss further in the next chapter. Note also that the solutions for xn and yn share the same recurrence:
FindLinearRecurrence@Table@dataPellsolutionsplusminus@@n, 1DD,
8n, 1, Length@dataPellsolutionsplusminusD<DD
FindLinearRecurrence@Table@dataPellsolutionsplusminus@@n, 2DD,
8n, 1, Length@dataPellsolutionsplusminusD<DD
82, 1<
82, 1<
FURTHER EXPLORATION: Find a formula involving xn , yn , and y2 n . Then prove your formula.
Sequences sometimes reveal themselves through patterns involving their subsequences. Given a sequence 8an <, a subsequence
8bn < is a sequence whose elements come from 8an < and are listed in the same order. More formally, 8bn < is called a subsequence
of 8an < if bn = a f HnL where f HnL is an increasing integer function, i.e., f HnL < f Hn + 1L. Here are some common subsequences that
we shall consider in this book:
S1. Even Elements: 8a2 n < = 8a0 , a2 , a4 , a6 , ...<
S2. Odd Elements: 8a2 n+1 < = 8a1 , a3 , a5 , a7 , ...<
S3. Elements at Square Positions: 8an2 < = 8a1 , a4 , a9 , a16 , ...<
S4. Elements at Powers of 2: 8a2n < = 8a1 , a2 , a4 , a8 , ...<
S5. Elements at Prime Positions: 8aPrimeHnL < = 8a2 , a3 , a5 , a7 , ...<, where PrimeHnL refers to the n-th prime, starting with
PrimeH1L = 2.
Let’s continue our search for patterns with the diatomic sequence 8an < discussed in Example 2.13 and compare the sequence itself
with some of its subsequences, for example, those at even positions and at odd positions.
41
Mathematics by Experiment
Patterns for these subsequences now become clear: a2 n = an and a2 n+1 = an + an+1 .
Thus, the diatomic sequence can be alternatively defined through the recurrence
:
a2 n = an
(2.12)
a2 n+1 = an + an+1
with a0 = 0, a1 = 1. We generate the first 10 terms using this recurrence to confirm that it agrees with the original definition.
a@0D = 0;
a@1D = 1;
a@n_D := If@Mod@n, 2D 0, a@n 2D, a@Hn - 1L 2D + a@Hn + 1L 2DD
Table@a@nD, 8n, 0, 9<D
80, 1, 1, 2, 1, 3, 2, 3, 1, 4<
Special Transformations
a) At the beginning of this chapter we had found a pattern for the partial sums of the Fibonacci sequence 8Fn } based on the
following table of values:
42
Chapter 2
Shifting the column for 8Fn < up two rows revealed the following identity:
âFk = Fn+2 - 1
n
(2.13)
k=0
For this simple example, we can avoid the above analysis altogether by asking Mathematica to evaluate the partial sums transfor-
mation symbolically, which yields the same identity:
Sum@Fibonacci@kD, 8k, 0, n<D
- 1 + Fibonacci@2 + nD
It is clear that these sums are given by the Fibonacci numbers at even positions (A001906):
NOTE: If one uses the index 2 k - 1 instead of 2 k + 1 to represent the odd Fibonacci numbers, then Mathematica doesn’t
recognize their partial sums as an identity, but instead gives an explicit formula.
43
Mathematics by Experiment
4n I1 + 5M J- 1 + I 2 I1 + 5 MM N
-2 n 1 4n
5
This follows from Mathematica’s knowledge of Binet's formula for the Fibonacci numbers:
J1 + 5 N - J1 - 5N
n n
Fn = (2.15)
n
2 5
To obtain identity (2.14), it suffices to verify that Mathematica‘s explicit formula is equal to F2 n , which follows from replacing n
by 2 n in Binet’s formula.
c) What do you expect the formula to be for the sum of the even Fibonacci numbers (make a guess before performing your
experiment)? Let's verify your conjecture:
Partial Sums of Odd Fibonacci Numbers
n Fn Únk=0 F2 k
0 0 0 = 0
1 1 1 = 0+1
2 1 4 = 0+1+3
3 2 12 = 0+1+3+8
4 3 33 = 0+1+3+8+21
5 5 88 = 0+1+3+8+21+55
6 8 232 = 0+1+3+8+21+55+144
7 13 609 = 0+1+3+8+21+55+144+377
8 21 1596 = 0+1+3+8+21+55+144+377+987
9 34 4180 = 0+1+3+8+21+55+144+377+987+2584
This time a slight twist emerges in the pattern: the sum of the first n even Fibonacci numbers is one less than the Fibonacci
number F2 n+1 :
â F2 k = F2 n+1 - 1
n
(2.16)
m=1
d) How about the partial sum of the squares of the Fibonacci numbers?
44
Chapter 2
k=0 k=0
81, 2, 3,
2 1 2
81, 3, 5,
3 2 6 6<
81, 2, 4,
4 3 15 15<
81, 2, 4,
5 5 40 5, 8, 10, 20, 40<
81, 3, 7,
6 8 104 8, 13, 26, 52, 104<
81, 2, 3,
7 13 273 13, 21, 39, 91, 273<
81, 2, 5,
8 21 714 6, 7, 14, 17, 21, 34, 42, 51, 102, 119, 238, 357, 714<
9 34 1870 10, 11, 17, 22, 34, 55, 85, 110, 170, 187, 374, 935, 1870<
Yes, in this case each sum is equal to the product of two consecutive Fibonacci numbers, also known as the golden rectangle
numbers (A001654):
Sum@Fibonacci@kD ^ 2, 8k, 1, n<D
Fibonacci@nD Fibonacci@1 + nD
â Fk2 = Fn Fn+1
n
(2.17)
m=1
NOTE: Observe in the table above that the two Fibonacci divisors Fn and Fn+1 seem to always be consecutive in the list of
divisors for Únm=1 Fk2 . In other words, there are no divisors of Únm=1 Fk2 that lie between Fn and Fn+1 . Can you explain why?
e) Let’s apply the determinant transformation to the Fibonacci sequence:
45
Mathematics by Experiment
The pattern in this case is quite simple and yields Cassini’s identity:
2
Fn Fn+2 - Fn+1 = H-1Ln+1 (2.18)
NOTE: Cassini’s identity is the basis for the following visual deception known as the Missing Square Puzzle: Take a square
whose sides are 8 units in length. Divide the square into four pieces according to the illustration below and rearrange them to
form a 5 x 13 rectangle:
5 8
3 5
However, the area the square is 82 = 64 square units, whereas the area of the rectangle is 5 13 = 65 square units. What hap-
pened to the missing square? The answer lies in the fact that the four pieces do not form an exact rectangle; their inclines have
different slopes. There is a gap in the shape of a parallelogram in the middle of the rectangle. This becomes more visible if we
zoom in:
FURTHER EXPLORATION: Explore on your own to find similar patterns for the determinant transformation of even Fibonacci
numbers. Repeat for the odd Fibonacci numbers.
46
Chapter 2
John Wallis, in his investigation of areas under polynomial functions, discovered some remarkable simple patterns regarding
ratios of power sums, i.e., sums of consecutive integers of the form 1 p + 2 p + ... + n p . Define
1 p + 2 p + ... + n p Únk=1 k p
RHn, pL = = (2.19)
n p + n p + ... + n p n p+1
where the denominator is the sum of n copies of n p . Wallis was able to detect patterns for RHn, pL, which we demonstrate for
RHn, 2L. Towards this end, let’s examine the values of RHn, 2L for 1 £ n £ 10 as shown in the table below.
n RHn,2L n RHn,2L
02 +12 1 02 +12 +22 +32 +42 +52 +62 13
1 2 6 62 +62 +62 +62 +62 +62 +62
36
12 +12
02 +12 +22 5 02 +12 +22 +32 +42 +52 +62 +72 5
2 12 7 72 +72 +72 +72 +72 +72 +72 +72
14
22 +22 +22
02 +12 +22 +32 7 02 +12 +22 +32 +42 +52 +62 +72 +82 17
3 18 8 82 +82 +82 +82 +82 +82 +82 +82 +82
48
32 +32 +32 +32
02 +12 +22 +32 +42 3 02 +12 +22 +32 +42 +52 +62 +72 +82 +92 19
4 8 9 92 +92 +92 +92 +92 +92 +92 +92 +92 +92
54
42 +42 +42 +42 +42
02 +12 +22 +32 +42 +52 11 02 +12 +22 +32 +42 +52 +62 +72 +82 +92 +102 7
5 52 +52 +52 +52 +52 +52
30 10
102 +102 +102 +102 +102 +102 +102 +102 +102 +102 +102 20
Visual inspection suggests that the numerators in the fractions above consist of the odd integers and the denominators consists of
multiples of 6, although the pattern is not consistent due to cancellation of factors. Mathematica confirms that this is indeed the
pattern:
FindSequenceFunction@81 2, 5 12, 7 18, 3 8, 11 30, 13 36, 5 14<, nD
1+2n
6n
Thus,
1+2n
RHn, 2L = (2.20)
6n
NOTE : For those with a background in calculus, observe that the limiting value of RHn, 2L as n ® ¥ equals 1/3, which gives the
exact area under the parabola f HxL = x2 along the interval @0, 1D. This is because RHn, 2L represents an approximation of this area
by n rectangles, which can be demonstrated as follows: partition the interval @0, 1D into n subintervals having a uniform width of
1 n. Each subinterval @Hk - 1L n, k nD then defines the base of a rectangle Rk whose height is specified by the value of f HxL at
the right endpoint of the subinterval. The area of each rectangle is given by
areaHRk L =
1 k k2
×f = (2.21)
n n n3
and the total of these areas, called a Riemann sum, is precisely RHn, 2L:
The method of finite differences is very powerful tool for determining formulas of sequences generated by a polynomial. In a
collection of his 50 most interesting columns (as judged by reader response), The Colossal Book of Mathematics, Martin Gardner
writes (p. 15):
47
Mathematics by Experiment
The method of finite differences is very powerful tool for determining formulas of sequences generated by a polynomial. In a
collection of his 50 most interesting columns (as judged by reader response), The Colossal Book of Mathematics, Martin Gardner
writes (p. 15):
Recreational problems involving permutations and combinations often contain low-order formulas
that can be correctly guessed by the method of finite differences and later (one hopes ) proved.
Given a sequence 8a1 , a2 , a3 , ...<, we define its first differences (or gap) to be
D an = an+1 - an (2.23)
and more generally its d-th differences (or differences to level d) to be
Dd an = Dd-1 an+1 - Dd-1 an (2.24)
for p ³ 1 where D0 an = an .
If a formula is known for the m-th differences, then this yields a recurrence for the original sequence by backwards computation:
Dd an = Dd-1 an+1 - Dd-1 an
= Dd-2 an+2 - 2 Dd-2 an+1 + Dd-2 an
= Dd-3 an+3 - 3 Dd-3 an+2 + 3 Dd-3 an+1 - Dd-3 an
... (2.25)
= âH-1Li
d
d
an+d-i
i
i=0
On the other hand, if all d-differences Dd ak are known for some fixed k, then an explicit formula for the original sequence is
given by
an+k = â K O Dm ak
n
n
m (2.26)
m=0
? Differences
Example 2.21
Consider the sequence 8an < consisting of sums of squares as defined by an = 12 + 22 + ... + n2 :
48
Chapter 2
Sums of Squares
n an
1 1 12
2 5 12 + 22
3 14 12 + 22 + 32
4 30 12 + 22 + 32 + 42
5 55 12 + 22 + 32 + 42 + 52
6 91 12 + 22 + 32 + 42 + 52 + 62
7 140 12 + 22 + 32 + 42 + 52 + 62 + 72
8 204 12 + 22 + 32 + 42 + 52 + 62 + 72 + 82
9 285 12 + 22 + 32 + 42 + 52 + 62 + 72 + 82 + 92
10 385 12 + 22 + 32 + 42 + 52 + 62 + 72 + 82 + 92 + 102
The following array lists the successive differences of 8an < up to level 4:
d Dd an
84, 9,
85, 7,
1 16, 25, 36, 49, 64, 81, 100<
82, 2,
2 9, 11, 13, 15, 17, 19<
80, 0,
3 2, 2, 2, 2, 2<
4 0, 0, 0, 0<
H- 1 + nL n + H- 2 + nL H- 1 + nL n
3 1
n+
2 3
Simplifying the expression above now yields the classic formula for sums of squares:
Factor@Simplify@%DD
n H1 + nL H1 + 2 nL
1
6
NOTE: Mathematica employs similar techniques to obtain the same answer when we ask it to evaluate an :
a@nD
n H1 + nL H1 + 2 nL
1
6
49
Mathematics by Experiment
Suppose we revise our definition of d-differences to ignore sign and consider only the absolute value of the differences, which
we refer to as d¤-differences:
Then computing d¤-differences for the sequence of primes reveals an interesting pattern known as Gilbreath’s conjecture: each
sequence of d¤-differences, except for the first (d = 1), begins with 1. This is indicated by the red entries in the array below.
Clear@sD
s@0, n_D := Prime@nD
s@m_, n_D := Abs@s@m - 1, n + 1D - s@m - 1, nDD
TableForm@
Table@If@m > 0 && n 1, Style@s@m, nD, RedD, s@m, nDD, 8m, 0, 11<, 8n, 1, 12 - m<DD
2 3 5 7 11 13 17 19 23 29 31 37
1 2 2 4 2 4 2 4 6 2 6
1 0 2 2 2 2 2 2 4 4
1 2 0 0 0 0 0 2 0
1 2 0 0 0 0 2 2
1 2 0 0 0 2 0
1 2 0 0 2 2
1 2 0 2 0
1 2 2 2
1 0 0
1 0
1
In 1993, this pattern was verified by Odlyzko [Od] for all primes up to 1013 ; however, no proof is known of Gilbreath’s conjec-
ture (see [Gi]).
bHnL = â
n
n
aHnL (2.27)
k
k=0
n
Recall that are binomial coefficients defined earlier.
k
The inversion of a sequence 8an < (also referred to as the binomial transform) is defined as
cHnL = âH-1Lk
n
n
aHnL (2.28)
k
k=0
Inversion refers to the fact that this transformation is equal to its inverse transform so that one can recover 8an < using the same
formula:
aHnL = âH-1Lk
n
n
cHnL (2.29)
k
k=0
Example 2.23
50
Chapter 2
Clear@bD;
b@n_D = Sum@Binomial@n, kD * Fibonacci@kD, 8k, 0, n<D
- I 2 I3 - 5 MM + I 2 I3 + 5 MM
1 n 1 n
5
Although Mathematica does not recognize it, a comparison of the first several values bn reveals the pattern to be nothing more
than the even Fibonacci numbers, F2 n .
Thus,
â
n
n
F = F2 n (2.30)
k k
k=0
Do you think this pattern is special to the Fibonacci sequence? What about other sequences having the same 81, 1< recurrence,
but with different initial values, such as the Lucas sequence 8Ln <?
Ln+1 = Ln + Ln-1 ; L0 = 2, L1 = 1 (2.31)
Here are the first ten Lucas numbers, generated using the Mathematica function LucasL.
Table@LucasL@nD, 8n, 0, 9<D
82, 1, 3, 4, 7, 11, 18, 29, 47, 76<
Applying the binomial transform shows that it yields a similar formula as the Fibonacci numbers:
Clear@cD;
c@n_D = Sum@Binomial@n, kD * LucasL@kD, 8k, 0, n<D
I3 - 5M I3 + 5M
1 n 1 n
+
2 2
Again, comparing the first several values of cn shows that indeed they are the Lucas numbers at even positions (A005248):
51
Mathematics by Experiment
â
n
n
L = L2 n (2.32)
k k
k=0
FURTHER EXPLORATION: Explore whether this identity holds for all 81, 1<-recurrences.
b) Now consider the inversion of the Fibonacci sequence:
Clear@bD;
b@n_D = Sum@H- 1L ^ k * Binomial@n, kD * Fibonacci@kD, 8k, 0, n<D
I 2 I1 - 5 MM - I 2 I1 + 5 MM
1 n 1 n
5
Do you recognize this formula? You should since it is nothing more than (the negative of) Binet’s formula for the Fibonacci
numbers as seen from its initial values given below.
Table@Simplify@b@nDD, 8n, 0, 10<D
80, - 1, - 1, - 2, - 3, - 5, - 8, - 13, - 21, - 34, - 55<
Thus,
âH-1Lk
n
n
F = -Fn (2.33)
k k
k=0
NOTE: Observe that Mathematica fails to recognize bn as the negative Fibonacci numbers, even if we apply the FindSequence-
Function command to its first ten values:
Table@Simplify@b@nDD, 8n, 0, 10<D
FindSequenceFunction@Table@Simplify@b@nDD, 8n, 0, 10<D, nD
80, - 1, - 1, - 2, - 3, - 5, - 8, - 13, - 21, - 34, - 55<
HFibonacci@nD - LucasL@nDL
1
2
Let’s proceed to the inversion of the Lucas numbers:
52
Chapter 2
Clear@cD;
c@n_D = Sum@H- 1L ^ k * Binomial@n, kD * LucasL@kD, 8k, 0, n<D
I1 - 5M I1 + 5M
1 n 1 n
+
2 2
Here, we find that the Lucas sequence is preserved under inversion:
Table@Simplify@c@nDD, 8n, 0, 10<D
82, 1, 3, 4, 7, 11, 18, 29, 47, 76, 123<
Thus,
âH-1Lk
n
n
L = Ln (2.34)
k k
k=0
NOTE: For those with a background in linear algebra, inversion is an involution, i.e., a linear transformation I whose square is
the identity transformation, I 2 = Id. It follows that 1 and -1 are its only eigenvalues. The example above demonstrates then that
the Lucas and Fibonacci sequences are eigenvectors corresponding to these eigenvalues, respectively. Can you find other
sequences that are eigenvectors of the inversion transform?
Convolution
Given two sequences 8an < and 8bn <, we define their convolution to be the sequence
cn = âak bn-k
n
(2.35)
k=0
Convolution arises natural in many applications such as signal processing and statistics. Later in this chapter, we’ll describe how
convolution formulas can be easily derived for many different sequences using generating functions.
Let’s investigate how sequences like the Fibonacci sequence behave under self-convolution. Mathematica gives the following
unappealing convolution formula:
Clear@cD;
c@n_D = Sum@Fibonacci@kD * Fibonacci@n - kD, 8k, 0, n<D
1
5 I5 + 5M
Instead, we’ll call on Mathematica to find a recurrence for 8cn < (A001629):
Table@Simplify@c@nDD, 8n, 1, 10<D
80, 1, 2, 5, 10, 20, 38, 71, 130, 235<
Let’s compare this recurrence with that of the convolved Lucas numbers.
53
Mathematics by Experiment
Clear@cD;
c@n_D = Sum@LucasL@kD * LucasL@n - kD, 8k, 0, n<D
I- 2 I1 + 5 MM
1 -n
5+ 5
I4 I9 + 5 M + I- 2 I3 + 5 MM I11 + 3 5 M + I5 + 5 M I4n + I- 2 I3 + 5 MM M nM
n n n
FindLinearRecurrence@84, 13, 22, 45, 82, 152, 274, 491, 870, 1531<D
82, 1, - 2, - 1<
Thus, we see that both recurrences are the same. As a last experiment, let’s convolve the Fibonacci sequence with the Lucas
sequence.
Clear@cD;
c@n_D = Sum@Fibonacci@kD * LucasL@n - kD, 8k, 0, n<D
I- 1 - 5M I- 3 - 5M H1 + nL
1 1 -n 1 n
-1 +
5 2 2
FindLinearRecurrence@dataConvolveFibonacciLucasD
82, 1, - 2, - 1<
Again, we obtain the same recurrence for 8cn < = 80, 2, 3, 8, 18, 47, ...< (A099920). We also observe the following: the first four
initial values of 8cn < are Fibonacci numbers. This leads us to conjecture that perhaps Fibonacci numbers appear as factors.
Indeed, by generating a table of values for the ratio cn Fn (except for n = 0), we discover the following pattern:
n cn Fn
1 2
2 3
3 4
4 5
5 6
6 7
7 8
8 9
9 10
10 11
âFk Ln-k = Hn + 1L Fn
n
(2.36)
k=0
FURTHER EXPLORATION:
a) Investigate whether convolutions of any two 81, 1<-sequences must satisfy a 82, 1, -2, -1<-recurrence. Can you prove this?
54
Chapter 2
b) Investigate recurrences for the convolution of 8m, 1<-sequences by experimenting with various positive integer values for m.
Can you find a general formula for the recurrence?
In 1751, Leonard Euler proposed to Christian Goldbach the problem of finding the number of ways that a n-sided regular polygon
can be divided into triangles using its diagonals, which we denote by En . Here are plots of the solutions for the first few values of
n:
n Polygon Divisions
n Number of solutions
3 1
4 2
5 5
6 14
7 42
8 132
Thus, En is given by the Catalan number Cn-2 (A000108), which is known to have formula
H2 nL !
n ! Hn + 1L !
Cn = (2.37)
Thus,
55
Mathematics by Experiment
H2 Hn - 2LL !
Hn - 2L ! Hn - 1L !
En = (2.38)
The Catalan numbers arise in many other applications in combinatorics (see [St] for 66 different ways of defining the Catalan
numbers) and have many interesting properties. For example, suppose we convolve the sequence of Catalan numbers with itself
by defining
cn = âCk Cn-k
n
(2.39)
k=0
n CHnL cHnL
0 1 1
1 1 2
2 2 5
3 5 14
4 14 42
5 42 132
6 132 429
7 429 1430
8 1430 4862
9 4862 16 796
This shows that convolution shifts the Catalan numbers (as well as the values En ) and yields the following recurrence formula:
Generating Functions
An extremely useful analytic approach to studying convolution of sequences is via generating functions. A function GHxL is
called a generating function for a sequence 8an < if its elements describe the power series coefficients of GHxL, i.e.
GHxL = âan xn
¥
(2.41)
n=0
The connection with convolution arises when the generating functions for two sequences 8an < and 8bn < are multiplied together;
their product yields a function whose power series coefficients give precisely the convolution of 8an < and 8bn <:
There are many techniques for finding the generating function for a given sequence. If the sequence has a recurrence, then a
formula for the generating function can typically be found by manipulating its power series. We demonstrate this with some
examples.
Let Fn denote the Fibonacci sequence as usual and define its generating function by
GHxL = âFn xn
¥
n=0
We then translate the recurrence for the Fibonacci numbers, Fn = Fn-1 + Fn-2 , into an identity involving GHxL as follows:
56
Chapter 2
GHxL = âFn xn = 0 + x + âFn xn = x + âHFn-1 + Fn-2 L xn = x + x âFn-1 xn-1 + x2 âFn-2 xn-2 = x + x GHxL + x2 GHxL
¥ ¥ ¥ ¥ ¥
? GeneratingFunction
A related command, FindGeneratingFunction, which combines the two functions GeneratingFunction and FindSequenceFunc-
tion, is useful when a formula for 8an < is not immediately available.
FindGeneratingFunction
? FindGeneratingFunction
FindGeneratingFunction@8a1 , a2 , …<, xD attempts to find a simple generating function in x whose nth series coefficient is an .
FindGeneratingFunction@88n1 , a1 <, 8n2 , a2 <, …<, xD
attempts to find a simple generating function whose ni th series coefficient is ai .
Generating functions are useful for counting partitions of integers. To demonstrate this, consider the following product of two
power series where we multiply term by term and carefully keep track of exponents:
âcn xn = âxi âx2 j = Ix0 + x1 + x2 + x3 + ...M Ix0 + x2×1 + x2×2 + x2×3 + ...M
¥ ¥ ¥
= x0+0 + x1+0 + Ix2+0 + x0+2×1 M + Ix3+0 + x1+2×1 M + Ix4+0 + x2+2×1 + x0+2×2 M + ...
(2.44)
= x0 + x1 + 2 x2 + 2 x3 + 3 x4 + ...
Thus, we see that the coefficient cn counts the number of ways that the integer n can be partitioned as a sum of i copies of 1 and j
copies of 2. For example, if n = 4, then the expansion above shows that there are three such partitions:
4=1+1+1+1
4=1+1+2
4=2+2
57
Mathematics by Experiment
Here's a table listing the number of partitions for n ranging from 1 to 10:
n Number of Partitions
1 1
2 2
3 2
4 3
5 3
6 4
7 4
8 5
9 5
10 6
The pattern is clear: the values of 8cn < consist of positive integers listed in order, but repeated (except for the first value)
(A008619). To prove this pattern, we compute the generating functions for the two given power series Ú¥ i=0 x and Új=0 x
i ¥ 2j
by
evaluating them directly in Mathematica: (the GeneratingFunction command is not necessary in this case):
G1 = Sum@x ^ i, 8i, 0, Infinity<D
1
1-x
G2 = Sum@x ^ H2 iL, 8i, 0, Infinity<D
1
1 - x2
Multiplying these two functions together gives us the generating function for Ú¥ n
n=0 cn x . To obtain a formula for cn , it suffices to
use Mathematica’s SeriesCoefficient command
SeriesCoefficient
? SeriesCoefficient
SeriesCoefficient@series, nD finds the coefficient of the nth -order term in a power series in the form generated by Series.
SeriesCoefficient@ f , 8x, x0 , n<D finds the coefficient of Hx - x0 Ln in the expansion of f about the point x = x0 .
SeriesCoefficientA f , 8x, x0 , n x <, 9y, y0 , n y =, …E finds a coefficient in a multivariate series.
H3 + H- 1Ln + 2 nL n ³ 0
1
4
0 True
Thus,
H3 + H-1Ln + 2 nL
1
cn = (2.45)
4
This confirms the pattern that we had observed in the table above. Of course, we could have obtained the same answer by
computing the coefficient sequence 8cn < as the convolution of the coefficient sequences 8an < and 8bn < corresponding to the series
n=0 an x = Úi=0 x and Ún=0 bn x = Új=0 x , respectively. Since an = 1 and bn = H1 + H-1L L 2, we find their convolution gives
Ú¥ n ¥ i ¥ n ¥ 2j n
58
Chapter 2
H3 + H- 1Ln + 2 nL
1
4
NOTE: It makes sense to assume in our definition of a generating function that GHxL = Ú¥ n
n=0 an x as a power series should
converge on some interval of non-zero length; however, this is not always necessary. Even if Ú¥ n
n=0 an x is a divergent series, we
can still algebraically manipulate Ún=0 an x as a formal power series and hope to get useful results (see [GKP], p. 346 for an
¥ n
Let a0 be an integer and 8a1 , a2 , ..., an < a finite set of positive integers. The expression
pn 1
= a0 +
qn 1
a1 + 1
(2.46)
a2 + 1
...+
an
is called a finite continued fraction having 8a1 , a2 , ..., an < as quotients. It is also denoted by @a0 ; a1 , ..., an D. If the quotients
8a1 , a2 , ...< is infinite, then the corresponding continued fraction @a0 ; a1 , a2 , ...D is said to be infinite, viewed as a limit of
@a0 ; a1 , a2 , ..., an D (called its convergents), i.e.,
@a0 ; a1 , a2 , ...D = lim @a0 ; a1 , a2 , ..., an D (2.47)
n®¥
Thus, the convergents pn qn = @a0 ; a1 , a2 , ..., an D provide rational approximations of x = @a0 ; a1 , a2 , ...D.
Generating Continued Fractions
The continued fraction representation of a real value x can be generated from its integer part dxt and fractional part x - dxt and
repeating this process as follows:
a0 = dxt; f1 = x - dxt
; f2 = a1 - da1 t
1
a1 =
f1
; f3 = a2 - da2 t
1
a2 =
f2
...
; f4 = a3 - da3 t
1
an =
fn
...
If fn = 0 for some n, then this process is terminated and the resulting continued fraction for x = @a0 ; a1 , a2 , ..., an D is finite.
Let x = Π. Then the first three quotients of the continued fraction for Π are given by
a0 = dΠt = 3; f1 = Π - dΠt = 0.14159 ...
59
Mathematics by Experiment
Convergents@Pi, 10D
:3, >
22 333 355 103 993 104 348 208 341 312 689 833 719 1 146 408
, , , , , , , ,
7 106 113 33 102 33 215 66 317 99 532 265 381 364 913
Unfortunately, these quotients do not follow a recognizable pattern. However, there are generalized continued fractions of Π that
do follow regular patterns (see [La]).
NOTE: Observe that the convergents above provide better and better approximations of Π, as they should since the sequence of
convergents must converge to Π:
dataconvergentspi =
Table@8n, Row@8Convergents@Pi, 10D@@nDD, "=", N@Convergents@Pi, 10D@@nDD, 9D<D<,
8n, 1, 10<D;
ColumnDataDisplay@dataconvergentspi, 10, 8"n", "pn qn "<, "Convergents of Π"D
Convergents of Π
n pn qn
1 3=3.00000000
22
2 7
=3.14285714
333
3 106
=3.14150943
355
4 113
=3.14159292
103 993
5 33 102
=3.14159265
104 348
6 33 215
=3.14159265
208 341
7 66 317
=3.14159265
312 689
8 99 532
=3.14159265
833 719
9 265 381
=3.14159265
1 146 408
10 364 913
=3.14159265
60
Chapter 2
Convergents@Sqrt@2D, 10D
:1, >
3 7 17 41 99 239 577 1393 3363
, , , , , , , ,
2 5 12 29 70 169 408 985 2378
Do you recognize the numerators 81, 3, 7, 17, 41, 99, 239, ...< and denominators 81, 2, 5, 12, 29, 70, 169, ...< in the convergents
above? They are precisely solutions to the Pell equations x2 - 2 y2 = ± 1 discussed in Example (2.1).
FURTHER EXPLORATION:
1. Explain why the quotients above are all the same, i.e., an = 2 for all integers n ³ 1.
2. Repeat this example by computing the quotients and convergents of x = 3 . Find patterns for them and determining whether
the numerators and denominators of the convergents satisfy Pell equations of any kind.
3. What about x = n where n is a positive integer?
An integer d is called a divisor (or factor) of n if d divides n. For example, 4 is a divisor of 12 since 4 divides 12 (i.e. 12 4 = 3).
The Mathematica command Divisors will generate a list of all the divisors of a given integer. For example, all the divisors of 20
are
Divisors@20D
81, 2, 4, 5, 10, 20<
A positive integer p > 1 is called prime if the only divisors are 1 and itself; otherwise, it is said to be composite. For example, 13
is prime since the only divisors are 1 and 13; however, 15 is composite since 15 = 3 × 5. One can use the Mathematica command
PrimeQ to determine whether an integer is prime.
PrimeQ@13D
PrimeQ@15D
True
False
n Prime Factors of 2n -1
1 1 11
2 3 31
3 7 71
4 15 31 × 51
5 31 311
6 63 32 × 71
7 127 1271
8 255 31 × 51 × 171
9 511 71 × 731
10 1023 31 × 111 × 311
61
Mathematics by Experiment
For n = 2, 4, 6, and 10, we observe that n + 1 (namely 3, 5, 7, and 11) appears as a prime factor of 2n - 1. Let’s confirm this for
the first twenty prime values of n + 1:
n Prime Factors of 2n -1
1 1 11
2 3 31
4 15 31 × 51
6 63 32 × 71
10 1023 31 × 111 × 311
12 4095 32 × 51 × 71 × 131
16 65 535 31 × 51 × 171 × 2571
18 262 143 33 × 71 × 191 × 731
22 4 194 303 31 × 231 × 891 × 6831
28 268 435 455 31 × 51 × 291 × 431 × 1131 × 1271
2340 -1 = 31 × 52 × 111 × 311 × 411 × 1371 × 9531 × 10211 × 44211 × 26 3171 × 43 6911 × 131 0711 × 550 8011 ×
23 650 0611 × 7 226 904 352 843 746 8411 × 9 520 972 806 333 758 4311 × 26 831 423 036 065 352 6111
Lagrange’s Four-Square Theorem states that any natural number n can be written as a sum of four squares of integers:
n21 + n22 + n23 + n24 = n (2.48)
But how many ways are there of doing this for each n, i.e. how many solutions are of the form 8n1 , n2 , n3 , n4 < where we distin-
guish sign and order? Here are several solutions for n = 2:
12 + 12 + 02 + 02 = 2
02 + 02 + 12 + 12 = 2
H-1L2 + H-1L2 + 02 + 02 = 2
We denote by rd HnL the number of ways of writing n as a sum of d squares. Then r4 H2L = 24. Here's a Mathematica command
called SumsOfSquaresRepresentations for generating all such representations (see http://mathworld.wolfram.com/-
SumofSquaresFunction.html):
62
Chapter 2
SumOfSquaresRepresentations@4, 2D
88- 1, - 1, 0, 0<, 8- 1, 0, - 1, 0<, 8- 1, 0, 0, - 1<, 8- 1, 0, 0, 1<,
8- 1, 0, 1, 0<, 8- 1, 1, 0, 0<, 80, - 1, - 1, 0<, 80, - 1, 0, - 1<, 80, - 1, 0, 1<,
80, - 1, 1, 0<, 80, 0, - 1, - 1<, 80, 0, - 1, 1<, 80, 0, 1, - 1<, 80, 0, 1, 1<,
80, 1, - 1, 0<, 80, 1, 0, - 1<, 80, 1, 0, 1<, 80, 1, 1, 0<, 81, - 1, 0, 0<,
81, 0, - 1, 0<, 81, 0, 0, - 1<, 81, 0, 0, 1<, 81, 0, 1, 0<, 81, 1, 0, 0<<
? SquaresR
SquaresR@d, nD gives the number of ways rd HnL to represent the integer n as a sum of d squares.
SquaresR@4, 2D
24
n r4 HnL
1 8
2 24
3 32
4 24
5 48
6 96
7 64
8 24
9 104
10 144
Let’s now compare these values with the divisors of n, and in particular with the sum of the divisors of n:
81, 3<
2 24 3
81, 2, 4<
3 32 4
81, 5<
4 24 7
81, 2, 3, 6<
5 48 6
81, 7<
6 96 12
81, 2, 4, 8<
7 64 8
81, 3, 9<
8 24 15
81, 2, 5, 10<
9 104 13
10 144 18
Do you see a connection r4 HnL and the sum of the divisors of n? It seems that the former is 8 times the latter, although this
relationship fails for n = 4 and n = 8. However, if we ignore the divisors 4 and 8, then it is true. Let’s confirm this for higher
values of n:
63
Mathematics by Experiment
n r4 HnL
81, 11<
Divisors of n Sum of Divisors of n
81, 2, 3, 4, 6, 12<
11 96 12
81, 13<
12 96 28
81, 2, 7, 14<
13 112 14
81, 3, 5, 15<
14 192 24
81, 2, 4, 8, 16<
15 192 24
81, 17<
16 24 31
81, 2, 3, 6, 9, 18<
17 144 18
81, 19<
18 312 39
Thus, we see that the PHnL depends only on the sum of those divisors that are NOT divisible by 4.
Theorem: The number of representations of n ³ 1 as the sum of four squares is given by
r4 HnL = 8 â d
(2.49)
d n,4Id
? DivisorSigma
For example, the divisors of n = 6 are 81, 2, 3, 6<. Thus, ΣH6L = 12 and so 6 is a perfect number.
DivisorSigma@1, 6D
12
Some higher perfect numbers are 28, 496, and 8128. It is unknown whether there are infinite many perfect numbers.
64
Chapter 2
No pattern emerges for ΣHnL until we restrict our attention to certain subsequences. The following table, which lists ΣHnL for n
prime, clearly shows that ΣHnL = n + 1 in this case:
n ΣHnL
2 3
3 4
5 6
7 8
11 12
13 14
17 18
19 20
23 24
29 30
This of course is due to the fact that any prime integer n has only two divisors, 1 and n (itself).
The next table lists the values of ΣHnL for those values of n equal to a power of 2:
n ΣHnL
2 3
4 7
8 15
16 31
32 63
It appears that ΣH2n L = 2n+1 - 1 and follows from the fact that the divisors of 2n are all the non-negative powers of 2 less than or
equal to 2n . Thus, ΣH2n L = 20 + 21 + ... + 2n = 2n+1 - 1.
Euler's totient (phi) function:
Two positive integers m and n are said to be relatively prime if they have no commond divisor other than 1. For example, 5 and
11 are relatively prime. However, 6 and 15 are NOT relatively prime, having 3 as a common divisor.
The totient (phi) function ΦHnL is defined to be the number of positive integers less than or equal to n that are relatively prime to
n. The corresponding command in Mathematica is EulerPhi:
EulerPhi
65
Mathematics by Experiment
? EulerPhi
For example, ΦH6L = 2 since there are two positive integers less than or equal to 6 that are relatively prime to 6, namely 1 and 5:
EulerPhi@6D
2
a) The following table lists the values of ΦHnL for the first thirty positive integers:
Euler's Φ function
n ΦHnL n ΦHnL n ΦHnL
1 1 11 10 21 12
2 1 12 4 22 10
3 2 13 12 23 22
4 2 14 6 24 8
5 4 15 8 25 20
6 2 16 8 26 12
7 6 17 16 27 18
8 4 18 6 28 12
9 6 19 18 29 28
10 4 20 8 30 8
If we again focus only on primes values of n, then the following pattern emerges:
ΦHnL for n prime
n ΦHnL
2 1
3 2
5 4
7 6
11 10
13 12
17 16
19 18
23 22
29 28
Thus, ΣHnL = n - 1 whenever n is prime. This is because every positive integer less than n is relatively prime to n (there are n - 1
such integers).
b) It's clear from the results that we’ve obtained so far that ΦHnL + ΣHnL = 2 n whenever n is prime.
66
Chapter 2
n ΣHnL+ΦHnL
2 4
3 6
5 10
7 14
11 22
13 26
17 34
19 38
23 46
29 58
Let’s generalize this problem and consider those values for n where ΦHnL + ΣHnL = k n for some fixed positive integer k > 2. For
example, here are some solutions for k = 3:
SigmaPhiSum@k_, number_D := Module@8i = 0, n = 1, list = 8<<,
While@i < number,
If@DivisorSigma@1, nD + EulerPhi@nD k * n, AppendTo@list, nD; i ++D; n ++D; list
D
SigmaPhiSum@3, 6D
8312, 560, 588, 1400, 85 632, 147 492<
Here are some interesting open questions involving solutions to ΦHnL + ΣHnL = k n (see [Gu]):
1. Are there infinitely many solutions for each k?
2. Is there an odd solution?
3. The solutions above are of the form 4 m? Are there solutions of the form 4 m + 2?
c) What about the product ΣHnL ΦHnL? Do you recognize a pattern?
n ΣHnL×ΦHnL
1 1
2 3
3 8
4 14
5 24
6 24
7 48
8 60
9 78
10 72
Again, if we focus on the primes, then we observe that Σ HnL Φ HnL is always one less than a square:
67
Mathematics by Experiment
n Prime
n ΣHnL×ΦHnL
2 3
3 8
5 24
7 48
11 120
13 168
17 288
19 360
23 528
29 840
But note that there are composite integers n where ΣHnL ΦHnL is also one less than a square.
SigmaPhiProduct@number_D := Module@8i = 0, n = 1, list = 8<<,
While@i < number, If@Sqrt@DivisorSigma@1, nD * EulerPhi@nD + 1D
IntegerPart@Sqrt@DivisorSigma@1, nD * EulerPhi@nDD + 1D,
AppendTo@list, nD; i ++D; n ++D; list
D
SigmaPhiProduct@10D
82, 3, 5, 6, 7, 11, 13, 17, 19, 22<
Here are some open problems involving Σ HnL Φ HnL (see [Gu]):
1. Are there infinitely values for n where Σ HnL Φ HnL is a perfect square?
2. Are there infinitely many composite values for n where Σ HnL Φ HnL is one less than a square?
2. Characterize those values for n where Σ HnL Φ HnL is divisible by n.
n 2n -1 2n -1 prime n 2n -1 2n -1 prime
1 1 False 11 2047 False
2 3 True 12 4095 False
3 7 True 13 8191 True
4 15 False 14 16 383 False
5 31 True 15 32 767 False
6 63 False 16 65 535 False
7 127 True 17 131 071 True
8 255 False 18 262 143 False
9 511 False 19 524 287 True
10 1023 False 20 1 048 575 False
68
Chapter 2
If we examine the values for n in which is 2n - 1 is prime (called Mersenne primes), we find that they are also prime, namely
n = 2, 3, 5, 7, 13, 17, 19. Thus, we’ve discovered the result that
Theorem: If 2n - 1 is prime, then n is prime.
NOTE:
1. Observe that the converse is false since n = 11 is prime, but 211 - 1 = 2047 = 23 × 89 is NOT prime.
2. Euclid proved that if 2n - 1 is prime, then 2n-1 H2n - 1L is an even perfect number. Euler showed that the converse is also true.
However, it is not known whether there exists odd perfect numbers. See Wikipedia entry [Wi-Pe].
Every integer n, when divided by a positive integer p, leaves a remainder r. In other words, n can be written in the form
n = q× p + r (2.50)
with 0 £ r £ n - 1. Using this fact, we define the congruence of n modulo p to be equal to r (called the residue) and write this in
shorthand as r º n Hmod pL. In the special case where p = 2, then r = 0 or 1 and determines the parity of n, i.e. n is even if r = 0
and n is odd if r = 1.
In this example we investigate which positive integers n can be expressed as a sum of two square integers, i.e. n21 + n22 = n, where
signs and order are distinguished. Then recall that r2 HnL denotes the number of solutions for 8n1 , n2 <. Then r2 H5L = 8. Here are
the 8 different representations obtained using the SumOfSquaresRepresentations command discussed in Example 2.31.
SumOfSquaresRepresentations@2, 5D
88- 2, - 1<, 8- 2, 1<, 8- 1, - 2<, 8- 1, 2<, 81, - 2<, 81, 2<, 82, - 1<, 82, 1<<
Also recall that we can compute r2 H5L using the SquaresR command.
SquaresR@2, 5D
8
n r2 HnL n r2 HnL
1 4 11 0
2 4 12 0
3 0 13 8
4 4 14 0
5 8 15 0
6 0 16 4
7 0 17 8
8 4 18 4
9 4 19 0
10 8 20 8
No pattern is evident from the table above. Let’s restrict our attention then to prime integers:
69
Mathematics by Experiment
n r2 HnL n r2 HnL
2 4 31 0
3 0 37 8
5 8 41 8
7 0 43 0
11 0 47 0
13 8 53 8
17 8 59 0
19 0 61 8
23 0 67 0
29 8 71 0
We find that starting with n = 3, the values of r2 HnL are nonzero for n = 5, 13, 17, 29, 37, 41, 53, 61. Is their a pattern to these
values? To discover the answer, let’s compute the residues of n (modulo 4).
Thus, it is now clear that odd prime p is expressible as a sum of two squares if and only if p º 1 Hmod 4L.
Let’s consider gaps (successive differences) along each column. Define D rd HnL = Hn + 1L Hmod dL - n Hmod dL. Then we see that
the values of D rd HnL are either 1 or 1 - d (except for the first column):
Table 2.3: Values of D rd HnL
70
Chapter 2
d=1 d=2 d=3 d=4 d=5 d=6 d=7 d=8 d=9 d=10
n=1 0 -1 1 1 1 1 1 1 1 1
n=2 0 1 -2 1 1 1 1 1 1 1
n=3 0 -1 1 -3 1 1 1 1 1 1
n=4 0 1 1 1 -4 1 1 1 1 1
n=5 0 -1 -2 1 1 -5 1 1 1 1
n=6 0 1 1 1 1 1 -6 1 1 1
n=7 0 -1 1 -3 1 1 1 -7 1 1
n=8 0 1 -2 1 1 1 1 1 -8 1
n=9 0 -1 1 1 -4 1 1 1 1 -9
n=10 0 1 1 1 1 1 1 1 1 1
A little analysis reveals that the choice of value depends on whether or not d divides Hn + 1L. Thus, we have the following
formula:
1 - d, if d Hn + 1L
D rd HnL = : (2.51)
1, otherwise
b) Define the sum of remainders function ΡHnL = Únd=1 Hn mod dL which sums the first n elements along each row in Table 2.2.
Here’s a listing of the first fifty values of ΡHnL (A004125):
Table 2.4: Sums of Remainders Function ΡHnL
Ρ@n_D := Sum@Mod@n, dD, 8d, 1, n<D
Observe that certain values repeat:ΡH1L = ΡH2L = 0, ΡH3L = ΡH4L = 1, ΡH7L = ΡH8L = 8, ΡH15L = ΡH16L = 36, ΡH31L = ΡH32L = 167.
These positions correspond to powers of 2. Does this pattern hold for all powers of 2? Let's verify for larger values:
Table 2.5: Values of ΡH2n L
71
Mathematics by Experiment
Are there any other patterns besides D ΡH2n - 1L = 0? Suppose we investigate those values for n where D ΡHnL = 1. However, an
experimental search of the first 10,000 values shows that there are only three cases where D ΡHnL = 1, namely n = 2, 9, 135.
DΡ@n_D := Ρ@n + 1D - Ρ@nD;
Do@
If@Ρ@n + 1D - Ρ@nD 1, Print@nDD, 8n, 1, 10 000<D
2
9
135
Thus, no pattern seems to exist. On the other hand, let’s consider when D ΡHnL = -1. Here, we find that for the first 10,000
values, there are four cases:
Do@
If@Ρ@n + 1D - Ρ@nD - 1, Print@nDD, 8n, 1, 10 000<D
72
Chapter 2
5
27
495
8127
A pattern now emerges: adding 1 to these values for n give the perfect numbers 6, 28, 496, and 8128. This leads us to the
following result:
Theorem: n is perfect if and only if D Ρ Hn - 1L = -1.
See [Sp] for a proof of this theorem.
NOTE: Observe the false pattern in Table 2.6: if n = 10 p, then ΡHnL = 10 q + 9 for some integer q. For example, D ΡH50L = 29.
Alas, this fails for n = 80 since D ΡH80L = 40.
Example 2.37
n 2n -1 mod n n 2n -1 mod n
1 0 11 1
2 1 12 3
3 1 13 1
4 3 14 3
5 1 15 7
6 3 16 15
7 1 17 1
8 7 18 9
9 7 19 1
10 3 20 15
It appears that 2n - 1 º 1 Hmod nL for n = 2, 3, 5, 7, 11, 13, 17, 19. These are prime numbers. We check this for the first twenty
primes:
n 2n -1 mod n n 2n -1 mod n
2 1 31 1
3 1 37 1
5 1 41 1
7 1 43 1
11 1 47 1
13 1 53 1
17 1 59 1
19 1 61 1
23 1 67 1
29 1 71 1
73
Mathematics by Experiment
Do@
If@Mod@2 ^ n - 1, nD 1 && ! PrimeQ@nD, Print@nDD, 8n, 1, 1000<D
341
561
645
Observe that these counterexamples correspond to the same counterexamples mentioned at the end of Example 2.30. This is
because Theorem 2.6 is equivalent to Theorem 2.2.
We shall refer to the sum S p HnL = 1n + 2n + ... + n p as a power sum. Here is table of power sums for 1 £ n £ 4 and 1 £ p £ 4:
Table 2.8: Values of 1n + 2n + ... + n p
p=1 p=2 p=3 p=4
n=1 11 1 12 1 13 1 14 1
n=2 11 + 21 3 12 + 22 5 13 + 23 9 14 + 24 17
n=3 11 + 21 + 31 6 12 + 22 + 32 14 13 + 23 + 33 36 14 + 24 + 34 98
n=4 11 + 21 + 31 + 41 10 12 + 22 + 32 + 42 30 13 + 23 + 33 + 43 100 14 + 24 + 34 + 44 354
74
Chapter 2
The following pattern emerges from the table above for odd primes p:
Theorem: S p H4 n + 2L = 2 n + 1
This is because of Fermat's Little Theorem:
Fermat’s Little Theorem: Let p be a prime. Then p I a implies a p-1 º 1 Hmod pL
p I a Þ a p-1 º 1 Hmod pL (2.52)
See if you can prove Theorem using Fermat’s Little Theorem.
b) Let's now restrict m to being primes:
Table 2.11:
p=1 p=2 p=3 p=4 p=5 p=6 p=7 p=8 p=9 p=10 p=11 p=12 p=13 p=14 p=15
n=2 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1
n=3 0 2 0 2 0 2 0 2 0 2 0 2 0 2 0
n=5 0 0 0 4 0 0 0 4 0 0 0 4 0 0 0
n=7 0 0 0 0 0 6 0 0 0 0 0 6 0 0 0
n=11 0 0 0 0 0 0 0 0 0 10 0 0 0 0 0
n=13 0 0 0 0 0 0 0 0 0 0 0 12 0 0 0
n=17 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
n=19 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
n=23 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
n=29 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
a) Let’s investigate the residues (remainders) of the Fibonacci numbers modulo 2. Below is a table listing the residues of the first
twenty Fibonacci numbers:
Table 2.12:
n Fn mod 2 n Fn mod 2
0 0 10 1
1 1 11 1
2 1 12 0
3 0 13 1
4 1 14 1
5 1 15 0
6 0 16 1
7 1 17 1
8 1 18 0
9 0 19 1
By examining these residues, it’s clear that they repeat every third entry and can be described using the formula
75
Mathematics by Experiment
0 if n º 0 mod 3
Fn mod 2 = : 1 if n º 1 mod 3 (2.54)
1 if n º 2 mod 3
ModAn2 , 3E
To prove this, we consider the three possible remainders for n when divided by 3: 0, 1, 2. This corresponds to n = 3 q,
n = 3 q + 1, and n = 3 q + 2, respectively. Then substituting each of these forms for n into n2 yields
Expand@n ^ 2 . n ® 83 q, 3 q + 1, 3 q + 2<D
99 q2 , 1 + 6 q + 9 q2 , 4 + 12 q + 9 q2 =
To obtain the congruence of each of these expressions, we apply modular arithmetic to obtain the desired answer. For example,
if n = 3 q + 1, then
1 + n + n2 º 1 + 6 q + 9 q2 º 1 mod 3 (2.55)
since both 6 q and 9 q2 are divisible by 3. This proves Fn mod 2 º 1. Another option is to have Mathematica perform the
modular arithmetic:
Simplify@Mod@1 + 6 q + 9 q ^ 2, 3D, Element@q, IntegersDD
1
The reader is encouraged to verify the other two cases by carrying out the same calculations.
b) Are the residues of the Fibonacci numbers periodic for other moduli? To answer this, we make a table of the residues
Fn Hmod mL, with n ranging from 1 to 20 and m ranging from 1 to 10:
Table 2.13:
Fn mod m Hn=1,...,15L
n=1 n=2 n=3 n=4 n=5 n=6 n=7 n=8 n=9 n=10 n=11 n=12 n=13 n=14 n=15
m=1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
m=2 1 1 0 1 1 0 1 1 0 1 1 0 1 1 0
m=3 1 1 2 0 2 2 1 0 1 1 2 0 2 2 1
m=4 1 1 2 3 1 0 1 1 2 3 1 0 1 1 2
m=5 1 1 2 3 0 3 3 1 4 0 4 4 3 2 0
m=6 1 1 2 3 5 2 1 3 4 1 5 0 5 5 4
m=7 1 1 2 3 5 1 6 0 6 6 5 4 2 6 1
m=8 1 1 2 3 5 0 5 5 2 7 1 0 1 1 2
m=9 1 1 2 3 5 8 4 3 7 1 8 0 8 8 7
76
Chapter 2
Fn mod m Hn=16,...,30L
n= n= n= n= n= n= n= n= n= n= n= n= n= n= n=
16 17 18 19 20 21 22 23 24 25 26 27 28 29 30
m=1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0
m=2 1 1 0 1 1 0 1 1 0 1 1 0 1 1 0
m=3 0 1 1 2 0 2 2 1 0 1 1 2 0 2 2
m=4 3 1 0 1 1 2 3 1 0 1 1 2 3 1 0
m=5 2 2 4 1 0 1 1 2 3 0 3 3 1 4 0
m=6 3 1 4 5 3 2 5 1 0 1 1 2 3 5 2
m=7 0 1 1 2 3 5 1 6 0 6 6 5 4 2 6
m=8 3 5 0 5 5 2 7 1 0 1 1 2 3 5 0
m=9 6 4 1 5 6 2 8 1 0 1 1 2 3 5 8
Looking at the first several rows of this array seem to indicate that the residues are periodic for every modulus m. Indeed, this is
true (can you prove why?). Here is a table listing the periods for m ranging from 1 to 20:
Table 2.14:
No formula in terms of m is known for these periods. However, there are other interesting patterns that can be extracted. Denote
the period of Fn Hmod mL by PHmL. Then observe that for m > 2, PHmL is always even. Another pattern involves the sequence of
residues 81, 0, 1<, which seems to appear with relative high frequency, especially in the second array above. Let’s isolate this
sequence in our table:
Table 2.15:
Fn mod m Hn=1,...,15L - Residue Pattern 81,0,1<
n=1 n=2 n=3 n=4 n=5 n=6 n=7 n=8 n=9 n=10 n=11 n=12 n=13 n=14 n=15
m=1 - - - - - - - - - - - - - - -
m=2 1 1 0 1 1 0 1 1 0 1 1 0 1 1 0
m=3 1 - - - - - 1 0 1 - - - - - 1
m=4 1 - - - 1 0 1 - - - 1 0 1 - -
m=5 1 - - - - - - - - - - - - - -
m=6 1 - - - - - - - - - - - - - -
m=7 1 - - - - - - - - - - - - - 1
m=8 1 - - - - - - - - - 1 0 1 - -
m=9 1 - - - - - - - - - - - - - -
77
Mathematics by Experiment
A discernible pattern now emerges. The locations of 81, 0, 1< seem to be regularly spaced. Here’s a table listing the positions of
the residue 0 at these locations in comparison to the periods PHmL.
Table 2.16:
820<
4 6
824<
5 20
816<
6 24
812, 24<
7 16
8 12
The pattern is now clear. The locations of 81, 0, 1< occur at positions n that are multiples of PHmL. Thus, we’ve discovered the
following result:
Theorem: Let i = 0 or 1. Then Fn±i = i Hmod mL if and only if P HmL n.
Can you prove this theorem?
2.2.5.1 Permutations
A permutation of a set of n elements is an ordering of its elements. For example, 81, 2, 3< and 82, 1, 3< are two different permuta-
tions of 81, 2, 3<. The Mathematica command Permutations will generate all permutations of a given set of elements:
? Permutations
The set 81, 2, ..., n< has n ! permutations because of the Multiplication Principle.
2.2.5.2 Inversions
78
Chapter 2
2.2.5.2 Inversions
An inversion of a permutation Σ is a pair of elements of Σ that is ‘out of order’. For example, the permutation 82, 3, 1< has two
inversions: 82, 1< and 83, 1<. Denote by iHΣL to be the total number of inversions of Σ. The Mathematica command Inversions
will calculate iHΣL:
? Inversions
In this example we analyze how the total number of inversions are distributed among all permutations on n elements. For
example, the following table lists iHΣL for all permutations of the set 81, 2, 3< corresponding to n = 3:
Total Number of Inversions
Σ
81,
IHΣL
81,
2, 3< 0
82,
3, 2< 1
82,
1, 3< 1
83,
3, 1< 2
83,
1, 2< 2
2, 1< 3
Thus, we see that the total number of inversions are distributed as follows:
Distribution of iHpL for permutations on 3 elements
k Number of Permutations with iHΣL=k
0 1
1 2
2 2
3 1
More generally, we define IHn, kL to be the number of permutations Σ on n elements having iHΣL = k. Below is a table of values
for IHn, kL:
IHn,kL
IHn,kL k=0 k=1 k=2 k=3 k=4 k=5 k=6 k=7 k=8 k=9 k=10
n=1 1 0 0 0 0 0 0 0 0 0 0
n=2 1 1 0 0 0 0 0 0 0 0 0
n=3 1 2 2 1 0 0 0 0 0 0 0
n=4 1 3 5 6 5 3 1 0 0 0 0
n=5 1 4 9 15 20 22 20 15 9 4 1
Observe that each row in the table above is symmtric. Also, we find that the following recurrence holds:
IHn, kL = IHn - 1, kL + IHn - 1, k - 1L (2.56)
FURTHER EXPLORATION: Can you find other patterns in the table for IHn, kL?
2.2.5.3 Derangements
79
Mathematics by Experiment
2.2.5.3 Derangements
Derangements are permutations with no fixed points, i.e., no elements remains in place. For example, there are two derange-
ments of the set 81, 2, 3<: 82, 3, 1< and 83, 1, 2<. Here is the Mathematica command for generating derangements:
? Derangements
The number of derangements of the set 81, ..., n< is called the subfactorial of n and denoted by ni (or ! n). Thus, 3i = 2. Here is
the corresponding Mathematica command:
? Subfactorial
Subfactorial@nD gives the number of permutations of n objects that leave no object fixed.
Subfactorial@4D
9
Subfactorial
n n¡
1 0
2 1
3 2
4 9
5 44
6 265
7 1854
8 14 833
9 133 496
Do you see a connection between any two consecutive terms? Let’s look at their ratios, ni Hn - 1Li .
80
Chapter 2
nMax = 9;
datasubfactorialsratios2elements = Table@
8n, Subfactorial@nD, If@n < 3, "-", Row@8Subfactorial@nD Subfactorial@n - 1D,
"»", N@Subfactorial@nD Subfactorial@n - 1D, 5D<DD<, 8n, 1, nMax<D;
ColumnDataDisplayAdatasubfactorialsratios2elements, 10,
9"n", "n¡ ", "n¡ Hn-1L¡ "=, "Subfactorial "E
Subfactorial
n n¡ n¡ Hn-1L¡
1 0 -
2 1 -
3 2 2»2.0000
9
4 9 2
»4.5000
44
5 44 9
»4.8889
265
6 265 44
»6.0227
1854
7 1854 265
»6.9962
14 833
8 14 833 1854
»8.0005
133 496
9 133 496 14 833
»8.9999
Observe that the ratios above are almost integers. If we round them off, then a simple linear progression emerges as a pattern.
Denote by @ni Hn - 1Li D to be the round of ni Hn - 1Li .
nMax = 9;
datasubfactorialsratios2elementsround = Table@8n, Subfactorial@nD,
If@n < 3, "-", Round@Subfactorial@nD Subfactorial@n - 1DDD<, 8n, 1, nMax<D;
ColumnDataDisplayAdatasubfactorialsratios2elementsround, 10,
9"n", "n¡ ", "@n¡ Hn-1L¡ D"=, "Subfactorial "E
Subfactorial
n n¡ @n¡ Hn-1L¡ D
1 0 -
2 1 -
3 2 2
4 9 4
5 44 5
6 265 6
7 1854 7
8 14 833 8
9 133 496 9
Since it appears that @ni Hn - 1Li D = n, we should consider scaling subfactorials by n, namely n × Hn - 1Li .
81
Mathematics by Experiment
nMax = 9;
datasubfactorialsscaled =
Table@8n, Subfactorial@nD, If@n < 3, "-", n * Subfactorial@n - 1DD,
If@n < 3, "-", Subfactorial@nD - n * Subfactorial@n - 1DD<, 8n, 1, nMax<D;
ColumnDataDisplayAdatasubfactorialsscaled, 10,
9"n", "n¡ ", "nHn-1L¡ ", "n¡ -nHn-1L¡ "=, "Subfactorial "E
Subfactorial
n n¡ nHn-1L¡ n¡ -nHn-1L¡
1 0 - -
2 1 - -
3 2 3 -1
4 9 8 1
5 44 45 -1
6 265 264 1
7 1854 1855 -1
8 14 833 14 832 1
9 133 496 133 497 -1
A permutation is called an even (or odd) if it can be expressed as an even number of of even (or odd) number of transpositions,
i.e., exchanges of two elements, respectively. An even (or odd) derangement is one that is also an even (or odd) permutation,
respectively. Denote by de HnL and do HnL the number of even and odd derangements, respectively. Here’s a table comparing de HnL
(A003221)and do HnL (A145221):
Table 2.17:
82
Chapter 2
In this example we demonstrate how to approximate 2 using a two-dimensional sequence (see [LP]). To start, consider the
sequence
8an < = :J 2 N > = :1,
n
2 , 2, 2 2 , 4, 4 2 , ...>
consisting of powers of 2 . Since 2 is unknown, let’s replace it with a crude approximation, say 1, to obtain a new sequence
8bn < = 81, 1, 2, 2, 4, 4, ...<
whose formula is given by
FindSequenceFunction@81, 2, 2, 4, 4, 8, 8<, nD
2- 2 + 2 I1 - H- 1Ln + 2 + H- 1Ln 2M
3 n
a) Let’s find a formula for bn, k . Unfortunately, the FindSequenceFunction command cannot be used directly on two-dimen-
sional sequences. Instead, we’ll find formulas for the entries in each row and paste them together to get our two-dimensional
formula. Define fk HnL =bn, k . The following table gives fomulas for fk HnL for k ranging from 1 to 10:
Table 2.19:
83
Mathematics by Experiment
Observe that each formula involves a linear combination of H-1Ln and 2 . Thus, we’ll need to find formulas for their corre-
sponding coefficients (as functions of k):
81, 3, 7, 17, 41, 99, 239, 577, 1393, 3363< and 81, 2, 5, 12, 29, 70, 169, 408, 985, 2378<
Entering these coefficients into FindSequenceFunction yields
FindSequenceFunction@83, 7, 17, 41, 99, 239, 577, 1393, 3363<, kD
- I1 - 2 M + 3 I1 + 2 M +2 2 I1 + 2M
k k k
2 I1 + 2M
2 I1 - 2 M + 4 I1 + 2 M +3 2 I1 + 2M
k k k
4 I1 + 2M
Substituting these formulas for the coefficients gives us the following formula for bn, k :
b@n_, k_D =
H-3+nL
KKK- J1 - 2 N + 3 J1 + 2 N +2 2 J1 + 2 N O J2 J1 + 2 NNO -
1 k k k
SimplifyB2 2
KK 2 J1 - 2 N + 4 J1 + 2 N +3 2 J1 + 2 N O J4 J1 + 2 NNO
k k k
2 +
KK 2 J1 - 2 N + 4 J1 + 2 N +3 2 J1 + 2 N O J4 J1 + 2 NNO H- 1Ln 2 OF
k k k
H-3+nL
I- 1 + 2 M JH- 1Ln I1 - 2 M + I1 + 2 M I3 + 2 2 MN
1
k k
22
b) What do you notice about the ratios bk1 bk0 as k ® ¥? Here’s a table listing the first twenty ratios:
Table 2.20:
84
Chapter 2
2
This is easy to derive from the formula for bk1 bk0 :
Simplify@b@1, kD b@0, kDD
2 J- I1 - 2 M + I1 + 2 M I3 + 2 2 MN
k k
I1 - 2 M + I1 + 2 M I3 + 2 2M
k k
O + J3 + 2 2N
k
1- 2
2 -K
bk1 bk0
1+ 2
= (2.61)
K O + J3 + 2 2N
k
1- 2
1+ 2
2 J3 + 2 2N
bk1 bk0 ®
J3 + 2 2N
= 2 (2.62)
Thus, the ratios bk1 bk0 provide a rational approximation of 2 . What about the ratios bkn+1 bkn as k ® ¥ for arbitrary n? See if
you can determine a pattern for their limits.
Simplify@b@2, kD b@1, kDD
2 I1 - 2 M + I1 + 2 M I4 + 3 2M
k k
- I1 - 2 M + I1 + 2 M I3 + 2 2M
k k
NOTE: Observe that each row of bn, k satisfies the recurrence b2 n, k = 2 bn, k . A recursive formula for bn, k in terms of the ele-
ments in the first row is given by
bkn = â
k
k 0
b (2.63)
j n+ j
j=0
85
Mathematics by Experiment
FURTHER EXPLORATION: Explore for powers of 3 instead of 2 as initial values for b0n .
86
Chapter 2
87
Mathematics by Experiment
http://www-history.mcs.st-and.ac.uk/Mathematicians/Wallis.html
John Wallis, Savilian Professor of Geometry and founding member of the Royal Society, first honed his pattern detection abilities
from his experience as a cryptologist during the English Civil War (1642-1651) in which he decoded Royalist messages for the
victorious Parliamentarians. In his most important work, Arithmetica Infinitorum, published in 1656, Wallis demonstrated how
he was able to derive his famous infinite product formula for Π by interpolating number patterns obtained by quadrature (area
enclosed by a curve). His journey begins with the area of a quarter unit circle, Π/4, bounded by the graph of y = 1 - x2 on the
interval @0, 1D.
y 1 - x2
1
-1 0 1
In the language of calculus, this shaded area can be represented by a definite integral and denoted by
à I1 - x M
1
12 Π
2
âx= (2.64)
0 4
Unfortunately, the definite integral above allows Wallis little room for algebraic manipulation. Instead, he considers the definite
integral Ù0 I1 - x12 M â x (where the exponents are replaced by their reciprocals) and more generally those of the form
1 2
Ù0 I1 - x1 p M â x for arbitrary integers p and q. Now, Wallis’ contemporaries, namely Robertval, Fermat and Pascal, had already
1 q
à x âx=
1 1
pq + 1
pq
(2.66)
0
As a result, he is now able to calculate Ù0 I1 - x1 p M â x using binomial expansion and the formula above. For example, since
1 q
I1 - x12 M = 1 - 2 x12 + x
2
(2.67)
it follows that the area under I1 - x12 M can be calculated as the sum of the areas under the three curves 1, 2 x12 , and x, which
2
it follows that the area under I1 - x12 M can be calculated as the sum of the areas under the three curves 1, 2 x12 , and x, which
2
à I1 - x M â x = à I1 - 2 x + xM â x = à 1 â x - à 2 x â x + à x â x = 1 - 2 ×
1
2
1 1 1 1 1 1 1
12 + 1
12 12 12
+ = (2.68)
0 0 0 0 0 1+1 6
With this technique in mind, Wallis then makes a tabulation of the values of Ù0 I1 - x1 p M â x using other positive integers for p
1 q
and q. This generates a symmetrical table given below consisting of values that Wallis recognized to be reciprocals of figurate
numbers, a generalization of the triangular numbers to higher-dimensions (discussed in detail in Chapter 3).
Ù0 H1-xp Lq âx
1
Ù0 I1 - x1 p M â x
gH p, qL = (2.69)
1 q
then
Hq + 1L Hq + 2L × × × Hq + pL
gH p, qL = f pq+1 = (2.70)
1×2× × × p
are figurate numbers whenever p and q are positive integers. The values of gH p, qL together form the Figurate Triangle, now
commonly referred to as Pascal's triangle (see Chapter 3), with the exception that the entries equal to 1 are missing:
gHp,qL
q=1 q=2 q=3 q=4 q=5
p=1 2 3 4 5 6
p=2 3 6 10 15 21
p=3 4 10 20 35 56
p=4 5 15 35 70 126
p=5 6 21 56 126 252
Next, to obtain a formula for w = gH1 2, 1 2L, which corresponds to the reciprocal of the original definite integral
Ù0 I1 - x2 M â x, Wallis interpolates the values of gHp, qL by assuming that the following formula for figurate numbers continues
1 12
gH p, qL = :
1 if p = 0
Hq+1L Hq+2L×××Hq+ pL (2.71)
1×2××× p
if p ¹ 0
Wallis it seems is not bothered by the fact that gH0, qL is not formally defined according to definition 1.6; the scent of the figurate
numbers must have been too strong for him not to follow. As a result of formula 1.8, Wallis is able to produce the following table
containing column values for half-integer values of q:
89
Mathematics by Experiment
gHp,qL
1 3 5 7 9
q=0 q= 2 q=1 q= 2 q=2 q= 2 q=3 q= 2 q=4 q= 2 q=5
1 3 5 7 9 11
p=1 2
1 2
2 2
3 2
4 2
5 2
6
3 15 35 63 99 143
p=2 8
1 8
3 8
6 8
10 8
15 8
21
5 35 105 231 429 715
p=3 16
1 16
4 16
10 16
20 16
35 16
56
To complete the table for half-integer values of p, he takes advantage of the following recurrence satisfied by figurate numbers
and again assumes that it continues to hold for all non-negative values of p and q, including p = -1 2:
q+ p+1
gH p + 1, qL = gH p, qL × (2.72)
p+1
This, coupled with the assumption that the values must be symmetric, i.e. gH p, qL = gHq, pL, allows him to fill in all the missing
values by expressing them in terms of w = gH1 2, 1 2L:
gHp,qL
1 1 3 5 7
q=- 2 q=0 q= 2 q=1 q= 2 q=2 q= 2 q=3 q= 2 q=4
1 w 1 w 3 4w 5 8w 35
p=- 2 ¥ 1 2 2 3 8 15 16 35 128
p=0 1 1 1 1 1 1 1 1 1 1
1 w 3 4w 15 8w 35 64 w 315
p= 2 2
1 w 2 3 8 5 16 35 128
1 3 5 7 9
p=1 2
1 2
2 2
3 2
4 2
5
H128 wL 21
3 w 4w 5 8w 35 64 w 105 1155
p= 2 3
1 3 2 3 8 15 16 128
3 15 35 63 99
p=2 8
1 8
3 8
6 8
10 8
15
H128 wL 15 H512 wL 35
5 4w 8w 7 64 w 63 231 3003
p= 2 15
1 5 2 15 8 16 128
5 35 105 231 429
p=3 16
1 16
4 16
10 16
20 16
35
Again we note that gH-1 2, pL is undefined, but Wallis does not seem to care, as by now he has committed himself to following
the trail marked by the pattern of figurate numbers.
With his table complete, Wallis focuses in on the row of values for p = 1 2 and factors them according to the pattern:
gH12,qL
1 1 3 5 7
q=- 2 q=0 q= 2 q=1 q= 2 q=2 q= 2 q=3 q= 2 q=4
H3 × 5 × 7L H4 × 6 × 8L H3 × 5 × 7 × 9L
1 3 4 3×5 4×6
H2 × 4 × 6L H3 × 5 × 7Lw H2 × 4 × 6 × 8L
p=12 2
w 1 w 2 3
w 2×4 3×5
w
Lastly, Wallis assumes that the ratios of consecutive terms in each row decrease montonically, i.e.,
90
Chapter 2
p+1 p+2
gJ 2
, qN gJ 2
, qN
p
> (2.73)
p+1
gI 2 , qM gJ , qN
2
3 × 3 × 5 × 5 × × × H2 n + 1L H2 n + 1L 2n+3 3 × 3 × 5 × 5 × × × H2 n + 1L H2 n + 1L 2n+2
2 × 4 × 4 × 6 × × × H2 nL H2 n + 2L 2 × 4 × 4 × 6 × × × H2 nL H2 n + 2L
<w< (2.76)
2n+2 2n+1
2 n+3 2 n+2
Since the factors 2 n+2
and 2 n+1
both converge to 1 in the limit as n ® ¥, we thus obtain Wallis’ formula by setting
w = 4 Π:
4 3×3×5×5×7×7× × ×
= (2.77)
Π 2×4×4×6×6×8× × ×
A brilliant mathematical tour-de-force!
91
Mathematics by Experiment
92
Chapter 2
Exercises
Readers can hunt on their own for patterns by solving the following exercises.
1. Products of Consecutive Integers
a. Find a formula for the product of four consecutive integers beginning with n, i.e., nHn + 1L Hn + 2L Hn + 3L in terms of
perfect squares. For example, here are the products for n ranging from 1 to 5:
s = Table@n * Hn + 1L Hn + 2L Hn + 3L, 8n, 1, 5<D
824, 120, 360, 840, 1680<
b. Prove your formula algebraically.
2. Sums of Squares of Fibonacci Numbers (see [Ho])
a. Find a pattern for the sum of squares of two consecutive Fibonacci numbers:
Sums of Squares
n FHnL2 +FHn+1L2
0 1
1 2
2 5
3 13
4 34
5 89
6 233
7 610
8 1597
9 4181
NOTE: Mathematica's FindSequenceFunction gives a direct formula for the sequence above, but not a very interesting one.
Find a more interesting formula involving the Fibonacci numbers.
b. Can you find patterns involving sums of squares of non-consecutive Fibonacci numbers, e.g., the even Fibonacci numbers?
c. What about sums of square of three consecutive Fibonacci numbers?
93
Mathematics by Experiment
881<<
n Partions of n
5 885<, 84, 1<, 83, 2<, 83, 1, 1<, 82, 2, 1<, 82, 1, 1, 1<, 81, 1, 1, 1, 1<<
4
Now define bHnL to be the number of partitions of n that consist only of 2's or 3's.
Number of Partitions
Consisting of 2's or 3's
n bHnL n bHnL n bHnL
1 0 11 2 21 4
2 1 12 3 22 4
3 1 13 2 23 4
4 1 14 3 24 5
5 1 15 3 25 4
6 2 16 3 26 5
7 1 17 3 27 5
8 2 18 4 28 5
9 2 19 3 29 5
10 2 20 4 30 6
94
Chapter 2
a. For which integer values of n does this equation yield rational solutions for x? Hint: Use the Solve command to obtain the
solution formula for x and then apply the FindInstance command to find particular solutions for n.
b. Find formulas for the two rational solutions for x. Then prove your formulas.
6. Divisibility of the Perrin sequence by primes (see [PS]
Define a sequence an by a0 = 3, a1 = 0, a2 = 2, and an+3 = an+1 + an . Find a divisibility pattern involving an .
Answer: an is divisible by prime integers n.
a@0D = 3; a@1D = 0; a@2D = 2;
a@n_D := a@n - 2D + a@n - 3D
NOTE: This result generalizes to sequences a0 = k, a1 = a2 =. .. ak-2 = 0, ak-1 = k - 1, and an+k = an+1 + an .
7. Lacunary Recurrences (see [Yo-P])
a. Consider the Fibonacci sequence F0 = 0, F1 = 1, Fn+1 = Fn + Fn-1 .
i. Find a recurrence for the subsequence 8F2 m <. Prove your answer.
ii. Find a recurrence for the subsequence 8Fa m < where a is a positive integer:
iii. Find a recurrence for the subsequence 8Fa m+b < where a and b are positive integers.
b. Next consider the sequence x0 = 0, x1 = 1, xn+1 = 2 xn + xn-1 .
i. Find a recurrence for the subsequence 8xa m < where a is a positive integer:
ii. Find a recurrence for the subsequence 8xa m+b < where a is a positive integer:
i=0
i=0 j=0
c. Denote the roots of the quartic equation x3 - x2 - x - 1 = 0 by a, b, c, d. Conjecture a recurrence for the trinomial sum
95
Mathematics by Experiment
and experimentally verify your conjecture. NOTE: You may find that Mathematica will have difficulty evaluating uHnL for large
values of n using the trinomial sum formula above. This shows that recurrences can be more effective in generating sequences
than direct formulas.
10. Fibonacci to Lucas (see [ES3])
Find a recurrence for the ratios of Fibonacci to Lucas numbers:
Table@Fibonacci@nD LucasL@nD, 8n, 1, 10<D
:1, >
1 1 3 5 4 13 21 17 55
, , , , , , , ,
3 2 7 11 9 29 47 38 123
11. Non-Totients (see [Pu])
A non-totient is an integer n for which there is no solution to the equation jHxL = n, where jHxL is Euler's totient function. Find a
number pattern among the set of non-totients and prove that it is true.
12. Continued Fractions
Find a formula for the convergents of the continued fraction (1,1,1,...):
Convergents@PadRight@8<, 10, 1DD
:1, 2, >
3 5 8 13 21 34 55 89
, , , , , , ,
2 3 5 8 13 21 34 55
Give a proof of your formula.
13. Subsets with no adjacent elements
Let sn denote the number of subsets of 81, 2, ..., n< such that no two elements are adjacent. For example, if n = 4, then the
subsets with no adjacent elements are: {{1,3},{2,4}, {1,4}}
For n = 5, we have: {{1,3},{1,4},{1,5},{2,4},{2,5}}
Recursion: Fn = Fn-1 + Fn-2 (non-adjacent subsets of {1,2,3,4,6} = non-adjacent subsets of {1,2,3,4,5} + non-adjacent subsets of
{1,2,3,4} (with 6 added as an element since any non-adjacent subset containing 6 will be a non-adjacent subset of {1,2,3,4} with
6 deleted))
14. Counting Triangles in a Square (see [BoKo])
Let T HnL denote the number of lattice triangles lying inside the region @0, nD ´@0, nD of the Cartesian plane whose sides lie on
lines of slope 0, ¥, 1, or -1. Determine a closed formula for T HnL.
15. Diophantine Triplets (see [DeBr])
A Diophantine triplet is a set of three positive integers Ha, b, cL such that a < b < c and a b + 1, b c + 1, and a c + 1 are all perfect
squares. Find patterns describing Diophantine triplets.
96
Chapter 3
3
Classical Number Patterns
In this chapter we describe some classical numerical experiments involving integer sequences to reveal the myriad of number
patterns that can arise. We make two disclaimers about these experiments. First they focus more on analysis of data and formula-
tion of mathematical conjectures as opposed to rigorous proofs of the results obtained, although proofs of certain conjectures are
given if they are short and elementary; otherwise, more complicated proofs are given in the appendix or referenced elsewhere in
the mathematics literature. Second we provide litte historical background behind the numerical experiments discussed,, some of
which have a rich history and date as far back as the Greeks, and instead provide references where appropriate.
97
Mathematics by Experiment
Triangular Numbers
n tn =Únk=1 k
1 1 = 1
2 3 = 1+2
3 6 = 1+2+3
4 10 = 1+2+3+4
5 15 = 1+2+3+4+5
6 21 = 1+2+3+4+5+6
7 28 = 1+2+3+4+5+6+7
8 36 = 1+2+3+4+5+6+7+8
9 45 = 1+2+3+4+5+6+7+8+9
10 55 = 1+2+3+4+5+6+7+8+9+10
STEP 1
The first goal in understanding any sequence of values defined by summation is to find an efficient formula for calculating it
(besides by brute force). For example, can you compute t100 without summing all 100 terms? Since no obvious multiplicative
formula appears in the table above, we shall demonstrate four different methods for extracting such a formula for tn :
Method 1:
This classic technique is a variation of the “summing in pairs” trick commonly told as an anecdote involving Carl Friedrich
Gauss, who at the age of 10, was able to quickly sum the first 100 positive integers to the surprise of his schoolmaster (see [Hay]).
Suppose we sum two copies of tn by adding corresponding terms in pairs, but reverse the order of the terms in the second copy:
Triangular Numbers
n tn +tn =Únk=1 k+Ú1k=n k n tn +tn =Únk=1 k+Ú1k=n k
1 = 1 21 = 1+2+3+4+5+6
1 = 1 21 = 6+5+4+3+2+1
1 6
------------- -------------
2 = 2 42 = 7+7+7+7+7+7
3 = 1+2 28 = 1+2+3+4+5+6+7
3 = 2+1 28 = 7+6+5+4+3+2+1
2 7
------------- -------------
6 = 3+3 56 = 8+8+8+8+8+8+8
6 = 1+2+3 36 = 1+2+3+4+5+6+7+8
6 = 3+2+1 36 = 8+7+6+5+4+3+2+1
3 8
------------- -------------
12 = 4+4+4 72 = 9+9+9+9+9+9+9+9
10 = 1+2+3+4 45 = 1+2+3+4+5+6+7+8+9
10 = 4+3+2+1 45 = 9+8+7+6+5+4+3+2+1
4 9
------------- -------------
20 = 5+5+5+5 90 = 10+10+10+10+10+10+10+10+10
15 = 1+2+3+4+5 55 = 1+2+3+4+5+6+7+8+9+10
15 = 5+4+3+2+1 55 = 10+9+8+7+6+5+4+3+2+1
5 10
------------- -------------
30 = 6+6+6+6+6 110 = 11+11+11+11+11+11+11+11+11+11
Aha! It is now clear from the table above that 2 tn equals the sum of n copies of n + 1, or equivalently, nHn + 1L. Assuming this,
we conclude that
nHn + 1L
tn = (3.2)
2
Method 2:
98
Chapter 3
n tn Divisors of tn
1 1 81<
2 4 81, 3<
3 10 81, 2, 3, 6<
4 20 81, 2, 5, 10<
5 35 81, 3, 5, 15<
6 56 81, 3, 7, 21<
7 84 81, 2, 4, 7, 14, 28<
8 120 81, 2, 3, 4, 6, 9, 12, 18, 36<
9 165 81, 3, 5, 9, 15, 45<
10 220 81, 5, 11, 55<
It appears that for n is a divisor of tn for n odd. This suggests that we should examine the values of tn n:
n tn n
1 1
3
2 2
3 2
5
4 2
5 3
7
6 2
7 4
9
8 2
9 5
11
10 2
99
Mathematics by Experiment
Triangular Numbers
n tn =Únk=1 k on en
1 1 = 1 1 = 1 0 = 8<
2 3 = 1+2 1 = 1 2 = 2
3 6 = 1+2+3 4 = 1+3 2 = 2
4 10 = 1+2+3+4 4 = 1+3 6 = 2+4
5 15 = 1+2+3+4+5 9 = 1+3+5 6 = 2+4
6 21 = 1+2+3+4+5+6 9 = 1+3+5 12 = 2+4+6
7 28 = 1+2+3+4+5+6+7 16 = 1+3+5+7 12 = 2+4+6
8 36 = 1+2+3+4+5+6+7+8 16 = 1+3+5+7 20 = 2+4+6+8
9 45 = 1+2+3+4+5+6+7+8+9 25 = 1+3+5+7+9 20 = 2+4+6+8
10 55 = 1+2+3+4+5+6+7+8+9+10 25 = 1+3+5+7+9 30 = 2+4+6+8+10
Clear patterns now emerge. Depending on whether n is odd or even, which we re-index as 2 n - 1 and 2 n, respectively, we find
that
o2 n-1 = o2 n = n2
As for en , observe that we can re-express it as a sum of odd integers by subtracting 1 from each even integer. For example,
e5 = 2 + 4 = H1 + 1L + H3 + 1L = H1 + 3L + 2 = o4 + 2
e6 = 2 + 4 + 6 = H1 + 1L + H3 + 1L + H5 + 1L = H1 + 3 + 5L + 3 = o6 + 3
More generally, we have
e2 n = e2 n+1 = o2 n + n = n2 + n = nHn + 1L
It follows that
d Dd tn
0 80, 1, 3, 6, 10, 15, 21, 28, 36, 45, 55<
1 81, 2, 3, 4, 5, 6, 7, 8, 9, 10<
2 81, 1, 1, 1, 1, 1, 1, 1, 1<
3 80, 0, 0, 0, 0, 0, 0, 0<
tn = â K O Dm a0 =
n
n n 0 n 1 n 2 nHn - 1L n2 + n nHn + 1L
D a0 + D a0 + D a0 = 1 × 0 + n × 1 + ×1 = = (3.4)
m 0 1 2 2 2 2
m=0
100
Chapter 3
n H1 + nL
1
2
2. A visual “proof without words” of (3.2) is given below ([]):
tn tn nHn+1L
n tn tn +tn+1
1 1 4
2 3 9
3 6 16
4 10 25
5 15 36
6 21 49
7 28 64
8 36 81
9 45 100
10 55 121
H1 + nL2
Can you sketch a similar “proof without words” to demonstrate this formula?
STEP 3
What about a formula for Tn = t1 + t2 + ... + tn , the sum of the first n triangular numbers?
101
Mathematics by Experiment
n tn t1 +...+tn+1
1 1 1
2 3 4
3 6 10
4 10 20
5 15 35
6 21 56
7 28 84
8 36 120
9 45 165
10 55 220
The pattern is not so simple here. Let’s try the same trick as before by examining the values of Tn n:
n tn Tn n
1 1 1
2 3 2
10
3 6 3
4 10 5
5 15 7
28
6 21 3
7 28 12
8 36 15
55
9 45 3
10 55 22
n tn 3Tn n
1 1 3
2 3 6
3 6 10
4 10 15
5 15 21
6 21 28
7 28 36
8 36 45
9 45 55
10 55 66
This shows that values for 3 Tn n are the same as the triangular numbers, except shifted by one position. Thus,
3 Tn n = tn+1 = Hn + 1L Hn + 2L 2, or equivalently,
nHn + 1L Hn + 2L
Tn = (3.6)
6
Again, we can obtain the same answer using Mathematica:
102
Chapter 3
n H1 + nL H2 + nL
1
6
NOTE: The numbers Tn are also referred to as tetrahedral (or pyramidal) numbers. This can be seen geometrically by stacking
the triangular numbers t1 , t2 ,...,tn on top of each other as layers to form the tetrahedron corresponding to Tn :
Tetrahedral Numbers
n=1 n=2 n=3 n=4
FURTHER EXPLORATION:
1. Can you detect any patterns for the weighted sum Únk=1 k tn+1 = t1 + 2 t2 + ... + n tn ? Use Mathematica to verify your conjec-
tures.
2. Leonard Euler observed that if n is a triangular number, then 9 n + 1, 25 n + 3, 49 n + 6, 81 n + 10, and so forth, are also
triangular numbers involving odd square multiples of n. Are there other triangular numbers that involve odd square multiples of
n?
103
Mathematics by Experiment
n rn
1 1
2 36
3 1225
4 41 616
5 1 413 721
A further search for larger square-triangular numbers becomes too exhaustive; thus, it is desirable to find a more efficient
formula for generating them.
Unfortunately, FindSequenceFunction fails to determine an explicit formula for square-triangular numbers:
Simplify@FindSequenceFunction@datasquaretriangularnumbers@@All, 2DD, nDD
FindSequenceFunction@81, 36, 1225, 41 616, 1 413 721<, nD
However, since these values are integer squares, perhaps we should employ their radicals. Here’s a table listing the square roots
of the six square-triangular numbers given above:
n rn
1 1
2 6
3 35
4 204
5 1189
Perhaps more terms are needed. Fortunately, we now have a recurrence formula, which we can use to quickly generate additional
terms:
datasquaretriangularnumbersmoreterms = LinearRecurrence@86, - 1<, 81, 6<, 10D
81, 6, 35, 204, 1189, 6930, 40 391, 235 416, 1 372 105, 7 997 214<
Feeding this longer list of terms into FindSequenceFunction yields the following formula for tn :
104
Chapter 3
Simplify@FindSequenceFunction@datasquaretriangularnumbersmoreterms, nDD
- I3 - 2 2 M + I3 + 2 2M
n n
4 2
Squaring this formula yields a corresponding formula for tn , originally discovered by Leonard Euler in 1788:
J3 + 2 2 N - J3 - 2 2N
n n 2
rn = (3.10)
4 2
tHt+1L
NOTE: A more novel approach to studying square-triangular numbers is to describe them algebraically, i.e., rn = s2 = 2
for
some positive integers s and t. Multiplying both sides of this equation by 8 and completing the square on the left-hand side yields
H2 t + 1L2 - 1 = 8 s2 . This is equivalent to Pell’s equation x2 - 2 y2 = 1 with x = 2 t + 1 and y = 2 s, which we discussed in the
previous chapter.
In fact, the values for rn correspond to the product x × y, where Hx, yL represent solutions to the more general Pell equation
x2 - 2 y2 = ± 1.
105
Mathematics by Experiment
n pn =Únk=1 H3k-2L
1 1 = 1
2 5 = 1+4
3 12 = 1+4+7
4 22 = 1+4+7+10
5 35 = 1+4+7+10+13
6 51 = 1+4+7+10+13+16
7 70 = 1+4+7+10+13+16+19
8 92 = 1+4+7+10+13+16+19+22
9 117 = 1+4+7+10+13+16+19+22+25
10 145 = 1+4+7+10+13+16+19+22+25+28
STEP 1
Conjecture a formula for the pentagonal numbers pn in terms of n. HINT: Try either Gauss’s technique of summing in pairs or
else examine their divisors. Does your formula match Mathematica's formula?
STEP 2
What connection do pentagonal numbers have with triangular numbers. HINT: Partition the pentagon that corresponds to each
pentagonal number into an appropriate number of triangles. Prove your conjecture algebraically.
“
Step 1
Based on the pattern, can you guess what the first five values for some higher-order polygonal numbers, say hexagonal, heptago-
nal, and octagonal?
Polygonal Numbers Pdn
Pdn n=1 n=2 n=3 n=4 n=5
d=3 1 3 6 10 15
d=4 1 4 9 16 25
d=5 1 5 12 22 35
d=6 1 6 15 28 45
d=7 1 7 18 34 55
d=8 1 8 21 40 65
106
Chapter 3
nHn + 1L
fn2 = f11 + f21 + ... + fn1 = 1 + 2 + ... + n = tn = (3.14)
2
nHn + 1L Hn + 2L
fn3 = f12 + f22 + ... + fn2 = t1 + t2 + ... + tn = Tn = (3.15)
6
...
fnd = f1d-1 + f2d-1 ... + fnd-1 (3.16)
Based on known formulas for fn1 , fn2 , and fn3 , can you conjecture a formula for fnd ?
The following is a tabulation of fnd :
Figurate Numbers fdn
fdn n=1 n=2 n=3 n=4 n=5
d=0 1 1 1 1 1
d=1 1 2 3 4 5
d=2 1 3 6 10 15
d=3 1 4 10 20 35
d=4 1 5 15 35 70
It is clear from this table that the following recursive identity holds:
fnd = fnd-1 + fn-1
d
(3.17)
This recurrence is the basis for many other fascinating patterns in the table, which in the context of figurate numbers is referred
to as the Figurate Triangle. However, it is more well known as Pascal’s triangle, the next topic in this chapter.
107
Mathematics by Experiment
Hx + yLn = y + ... + K O yn
n n n n-1 n n-2 2 n
x + x y+ x (3.18)
0 1 2 n
= K O = 1.
n n
It is clear that
0 n
5 5!
For example, = = 10. Can you prove formula (3.20)? A proof is given in Appendix A.1.
2 2!×3!
The most easily recognized patttern involving Pascal’s triangle is the following: every interior entry is the sum of the two entries
above it (e.g. the entry 3 in the third row is the sum of entries 1 and 2 in the second row). This fundamental relationship is
expressed mathematically by the identity (called Pascal’s identity)
n n-1 n-1
= + (3.21)
k k-1 k
Observe that if we merely input this identity into Mathematica, then Mathematica is not able to verify it:
108
Chapter 3
However, if we use the FullSimplify command and assume n and k to be non-negative integers, then Mathematica is now able to
confirm it:
FullSimplify@Binomial@n, kD Binomial@n - 1, k - 1D + Binomial@n - 1, kD,
Element@n, IntegersD && n >= 0 && Element@k, IntegersD && k >= 0D
True
NOTE: The equilateral configuration of Pascal’s triangle is quite arbitrary and is not the only one that is useful. A more natural
configuration which aligns coefficients corresponding to the monomials xk yn-k (see 1.1-1.5) is that of a right triangle (column
justified), a form used originally by Michael Stifel and other when it made it appearance in Western mathematical texts in the
1500's:
Binomial Coefficents H L
n
k
H L k=0 k=1 k=2 k=3 k=4 k=5 k=6 k=7 k=8 k=9
n
k
n=0 1
n=1 1 1
n=2 1 2 1
n=3 1 3 3 1
n=4 1 4 6 4 1
n=5 1 5 10 10 5 1
n=6 1 6 15 20 15 6 1
n=7 1 7 21 35 35 21 7 1
n=8 1 8 28 56 70 56 28 8 1
n=9 1 9 36 84 126 126 84 36 9 1
n
We note that Pascal's triangle can also be displayed in array form by including those values of for 1 £ n < k , a form that is
k
referred to as Pascal's matrix (or square).
Binomial Coefficents J N
n
k
k=0 k=1 k=2 k=3 k=4 k=5 k=6 k=7 k=8 k=9
n=0 1 0 0 0 0 0 0 0 0 0
n=1 1 1 0 0 0 0 0 0 0 0
n=2 1 2 1 0 0 0 0 0 0 0
n=3 1 3 3 1 0 0 0 0 0 0
n=4 1 4 6 4 1 0 0 0 0 0
n=5 1 5 10 10 5 1 0 0 0 0
n=6 1 6 15 20 15 6 1 0 0 0
n=7 1 7 21 35 35 21 7 1 0 0
n=8 1 8 28 56 70 56 28 8 1 0
n=9 1 9 36 84 126 126 84 36 9 1
109
Mathematics by Experiment
Step 1
Let's consider the sum of the terms in each row of Pascal's triangle. For example, the sum of the entries in row n = 3 equals 8
(entires in this row are shown in red below):
Binomial Coefficents J N
n
k
k=0 k=1 k=2 k=3 k=4 k=5
n=0 1
n=1 1 1
n=2 1 2 1
n=3 1 3 3 1
n=4 1 4 6 4 1
n=5 1 5 10 10 5 1
The following table lists the sums of the first 6 rows (n = 0 corresponds to the first row):
Únk=0 J N
n
n
k
0 1 = 1
1 2 = 1+1
2 4 = 1+2+1
3 8 = 1+3+3+1
4 16 = 1+4+6+4+1
5 32 = 1+5+10+10+5+1
The table reveals that each row sums to a power of 2. Thus we have discovered the another well known fundamental identity
involving Pascal's triangle:
â
n
n
= 2n (3.22)
k
k=0
â
n
n
= 1 + 1 + ... + 1 = n (3.23)
0
k=0
For example, the partial sum of the first four entries in the first colum (corresponding to n = 3) equals 4. Here is a table of partial
sums for n ranging from 1 to 5:
110
Chapter 3
Únk=0 J N
n
n
k
0 1 = 1
1 2 = 1+1
2 3 = 1+1+1
3 4 = 1+1+1+1
4 5 = 1+1+1+1+1
5 6 = 1+1+1+1+1+1
The interesting observation is that these partial sums are recorded precisely in the second column of Pascal's triangle.
Binomial Coefficents J N
n
k
k=0 k=1 k=2 k=3 k=4 k=5
n=0 1
n=1 1 1
n=2 1 2 1
n=3 1 3 3 1
n=4 1 4 6 4 1
n=5 1 5 10 10 5 1
To test if this is a coincidence, we compute the partial sums of the second column, given by
â = ân = 1 + 2 + 3 + ... + n
n n
n
1 (3.24)
k=1 k=1
to see if a similar pattern holds (the entries corresponding case n = 4 is displayed in red):
Únk=0 J N
n
n
k
1 1 = 1
2 3 = 1+2
3 6 = 1+2+3
4 10 = 1+2+3+4
5 15 = 1+2+3+4+5
Indeed, we find that the partial sums 81, 3, 6, 10, 15, ...< are again given precisely by the third column in Pascal's triangle, which
the reader should recognize as the triangular numbers tn discussed in the previous section.
111
Mathematics by Experiment
Binomial Coefficents J N
n
k
k=0 k=1 k=2 k=3 k=4 k=5
n=0 1
n=1 1 1
n=2 1 2 1
n=3 1 3 3 1
n=4 1 4 6 4 1
n=5 1 5 10 10 5 1
Verify on your own that this relationship continues to hold for any column so that the sum of the first n entries in the m-th column
of Pascal's triangle is given by the Hn + 1L-th entry in the Hm + 1L-th column:
Thus, we’ve discovered classic binomial formula:
â
n
k n+1
= (3.25)
m m+1
k=m
HPascal's TriangleL
n Sum of n-th Diagonal
0 1 = 1
1 1 = 1+0
2 2 = 1+1+0
3 3 = 1+2+0+0
4 5 = 1+3+1+0+0
5 8 = 1+4+3+0+0+0
It appears that the diagonals sum to the Fibonacci numbers F1 = 1, F2 = 1, F3 = 2, F4 = 3, F5 = 5, F6 = 8. Thus we have
discovered the identity
112
Chapter 3
âK O = Fn+1
n
n-m
m (3.26)
m=0
Can you prove this? HINT: Use Pascal’s Identity to prove that the left hand side of (3.26) satisfies the same recurrence as the
Fibonacci numbers.
-n
What patterns do you observe? Can you find a formula for where n > 0 and k > 0?
k
113
Mathematics by Experiment
nMax = 11;
databinomialmodulon =
Table@If@n 0, "-", Mod@Binomial@n, kD, nDD, 8n, 0, nMax<, 8k, 0, nMax<D;
ColumnB:"H L Hmod nL",
n
k
Grid@Prepend@HPrepend@databinomialmodulon@@ðDD, "n=" ~~ ToString@ð - 1DDL &
Range@1, Length@databinomialmodulonDD,
Prepend@H"k=" ~~ ToString@ðDL & Range@0, nMaxD, " "DD, Frame ® All,
Alignment ® Center, Background ® 881 ® LightBrown<, 81 ® LightBrown<<D>, CenterF
H L Hmod nL
n
k
k=0 k=1 k=2 k=3 k=4 k=5 k=6 k=7 k=8 k=9 k=10 k=11
n=0 - - - - - - - - - - - -
n=1 0 0 0 0 0 0 0 0 0 0 0 0
n=2 1 0 1 0 0 0 0 0 0 0 0 0
n=3 1 0 0 1 0 0 0 0 0 0 0 0
n=4 1 0 2 0 1 0 0 0 0 0 0 0
n=5 1 0 0 0 0 1 0 0 0 0 0 0
n=6 1 0 3 2 3 0 1 0 0 0 0 0
n=7 1 0 0 0 0 0 0 1 0 0 0 0
n=8 1 0 4 0 6 0 4 0 1 0 0 0
n=9 1 0 0 3 0 0 3 0 0 1 0 0
n=10 1 0 5 0 0 2 0 0 5 0 1 0
n=11 1 0 0 0 0 0 0 0 0 0 0 1
Do you see any patterns involving the rows? Columns? HINT: Consider the rows and columns at prime positions.
FURTHER EXPLORATION: Binomials to Binomials (see [Os])
Recall that binomial coefficients describe the expansion of H1 + xLn . Let’s fix x = 2 and consider the expansion of the binomial
K1 + 2 O = an + bn
n
2 (3.28)
K1 + 2 O =1
0
K1 + 2 O =1+
1
2
K1 + 2 O =1+2 2 +K 2 O =3+2
2 2
2
114
Chapter 3
a) Can you find both recursive and explicit formulas for an and bn ?
c) Find recursive formulas for an and bn (in terms of c and d) for the binomial expansion Jc + d N = an + bn
n
2.
a2 3 + 2 a
115
Mathematics by Experiment
We rule out the negative solution a = -1; thus, the only solution is a = 3, corresponding to the triple 83, 4, 5<.
Step 2
Let’s relax the restriction that all three sides must be consecutive integers by considering triples 8a, b, c< where say, the longer
leg and hypotenuse, are consecutive integers, i.e., c = b + 1; one example, besides 83, 4, 5<, is the triple 85, 12, 13<. Now, not all
Pythagorean triples have this property, for example, 88, 15, 17<. On the other hand, are there more such triples? Are there
infinitely many such triples? If so, is there a formula for generating them?
To answer these questions, again it suffices to substitute c = b + 1 into the equation a2 + b2 = c2 and solve for b in terms of a:
Clear@a, b, eqD;
eq = Simplify@a ^ 2 + b ^ 2 == Hb + 1L ^ 2D
Solve@eq, bD
a2 1 + 2 b
::b ® I- 1 + a2 M>>
1
2
We now argue as follows: if b is required to a positive integer, then a2 - 1 must be even (divisible by 2). It follows that a2 must
be an odd integer greater than 1 and hence a must be odd integer greater than 1. This allows us to index the solutions as
Clear@a, b , cD;
a@n_D := 2 n + 1;
b@n_D := Ha@nD ^ 2 - 1L 2;
c@n_D := b@nD + 1;
8a@nD, Simplify@b@nDD, Simplify@c@nD, Assumptions ® Element@n, IntegersD && n > 0D<
91 + 2 n, 2 n H1 + nL, 1 + 2 n + 2 n2 =
where n ranges over the positive integers. Here is a table listing the first ten solutions:
Pythagorean Triples 8a,b,c< with c=b+1
n an =2n+1 bn =2nHn+1L cn =2n2 +2n+1
1 3 4 5
2 5 12 13
3 7 24 25
4 9 40 41
5 11 60 61
6 13 84 85
7 15 112 113
8 17 144 145
9 19 180 181
10 21 220 221
FURTHER EXPLORATION: Can you find formulas to describe Pythagorean triples where the hypotenuse and longer leg
differ by 2?
Step 3
A more challenging problem is to find Pythagorean triples where both legs are consecutive integers, i.e., triples 8a, b, c< with
b = a + 1. Substituting this restriction into a2 + b2 = c2 yields
116
Chapter 3
Simplify@a ^ 2 + b ^ 2 c ^ 2 . 8b ® a + 1<D
a2 + H1 + aL2 c2
This equation is not as useful as in previous cases. All we can conclude here is that c must be odd. This follows from the fact if
a2 is even (or odd), then Ha + 1L2 is odd (or even, respectively); thus, c2 = a2 + b2 must be odd, or equivalently, c must be odd.
Here is a table of the first seven solutions obtained through a brute force search:
Pythagorean Triples 8a,b,c< with b=a+1
n an bn =an +1 cn = a2 + Ha + 1L2
1 3 4 5
2 20 21 29
3 119 120 169
4 696 697 985
5 4059 4060 5741
6 23 660 23 661 33 461
7 137 903 137 904 195 025
These are enough solutions for Mathematica to find a formula and recurrence for an :
Simplify@FindSequenceFunction@dataconsecutivelegs@@All, 1DD, nDD
I- 6 - 4 2 - I3 - 2 2 M I1 + 2 M + I3 + 2 2 M I7 + 5 2 MM
1 n n
12 + 8 2
FindLinearRecurrence@dataconsecutivelegs@@All, 1DDD
87, - 7, 1<
As for cn , we find
Simplify@FindSequenceFunction@dataconsecutivelegs@@All, 3DD, nDD
II3 - 2 2 M I2 + 2 M + I3 + 2 2 M I10 + 7 2 MM
1 n n
12 + 8 2
FindLinearRecurrence@dataconsecutivelegs@@All, 3DDD
86, - 1<
NOTE: Mathematica doesn’t seem to recognize that an satisfies an even more simple, although non-homogeneous, recurrence:
an+1 = 6 an - an-1 + 2.
RecurrenceTable@8a@n + 1D 6 a@nD - a@n - 1D + 2, a@1D 3, a@2D 20<, a, 8n, 1, 7<D
83, 20, 119, 696, 4059, 23 660, 137 903<
117
Mathematics by Experiment
a=2n
b = n2 - 1
c = n2 + 1
Step 1
Let’s complete the analysis for all higher values of H. Substituting c = b + H into a2 + b2 = c2 yields
Clear@a, b, cD;
sol = Solve@Simplify@a ^ 2 + b ^ 2 c ^ 2 . 8c ® b + H<D, bD
::b ® >>
a2 - H2
2H
We now use the formula above for b to tabulate values of Pythagorean triples for H ranging from 1 to 10:
An analysis of the differences between consecutive values of an in the first few tables seems to suggest the following:
an+1 - an = :
H if H even
2 H if H odd
However, this pattern fails for H = 8 and H = 9. A check of the higher values of H (up to H = 50) shows that this pattern also
fails for the following values: 16, 18, 24, 25, 27, 32, 36, 40, 45, 48, and 50 as shown in the table below:
118
Chapter 3
Formula for an
H an H an
1+2n 11 11 H1 + 2 nL
2 H1 + nL
1
2 12 12 H1 + nL
3 3 H1 + 2 nL 13 13 H1 + 2 nL
4 4 H1 + nL 14 14 H1 + nL
5 5 H1 + 2 nL 15 15 H1 + 2 nL
6 6 H1 + nL 16 8 H2 + nL
7 7 H1 + 2 nL 17 17 H1 + 2 nL
8 4 H2 + nL 18 6 H3 + nL
9 3 H3 + 2 nL 19 19 H1 + 2 nL
10 10 H1 + nL 20 20 H1 + nL
H an H an H an
21 21 H1 + 2 nL 31 31 H1 + 2 nL 41 41 H1 + 2 nL
22 22 H1 + nL 32 8 H4 + nL 42 42 H1 + nL
23 23 H1 + 2 nL 33 33 H1 + 2 nL 43 43 H1 + 2 nL
24 12 H2 + nL 34 34 H1 + nL 44 44 H1 + nL
25 5 H5 + 2 nL 35 35 H1 + 2 nL 45 15 H3 + 2 nL
26 26 H1 + nL 36 12 H3 + nL 46 46 H1 + nL
27 9 H3 + 2 nL 37 37 H1 + 2 nL 47 47 H1 + 2 nL
28 28 H1 + nL 38 38 H1 + nL 48 24 H2 + nL
29 29 H1 + 2 nL 39 39 H1 + 2 nL 49 7 H7 + 2 nL
30 30 H1 + nL 40 20 H2 + nL 50 10 H5 + nL
Observe that this list includes the odd perfect squares 9, 25, and 49. Let’s focus then on just these types of values:
Formula for an
H an H an
1+2n 21 21 H21 + 2 nL
3 H3 + 2 nL
1
3 23 23 H23 + 2 nL
5 5 H5 + 2 nL 25 25 H25 + 2 nL
7 7 H7 + 2 nL 27 27 H27 + 2 nL
9 9 H9 + 2 nL 29 29 H29 + 2 nL
11 11 H11 + 2 nL 31 31 H31 + 2 nL
13 13 H13 + 2 nL 33 33 H33 + 2 nL
15 15 H15 + 2 nL 35 35 H35 + 2 nL
17 17 H17 + 2 nL 37 37 H37 + 2 nL
19 19 H19 + 2 nL 39 39 H39 + 2 nL
an = H K H + 2 nO = H2 h + 1L HH2 h + 1L + 2 nL
Let's now check to see if this pattern holds for the even perfect squares:
119
Mathematics by Experiment
Formula for an
H an H an
2 4 H1 + nL 22 44 H11 + nL
4 8 H2 + nL 24 48 H12 + nL
6 12 H3 + nL 26 52 H13 + nL
8 16 H4 + nL 28 56 H14 + nL
10 20 H5 + nL 30 60 H15 + nL
12 24 H6 + nL 32 64 H16 + nL
14 28 H7 + nL 34 68 H17 + nL
16 32 H8 + nL 36 72 H18 + nL
18 36 H9 + nL 38 76 H19 + nL
20 40 H10 + nL 40 80 H20 + nL
If we distribute a factor of 2 through each formula for an , then indeed the same pattern holds in the case H = H2 hL2 .
What about those other exceptional values of H that are not perfect squares such as 8, 18, 24, 27, 32, etc. It appears that these
values are multiples of perfect squares. As a start, let’s focus on double perfect squares, i.e., H = 2 h2 :
Formula for an
h= H 2 an h= H 2 an
1 2 H1 + nL 11 22 H11 + nL
2 4 H2 + nL 12 24 H12 + nL
3 6 H3 + nL 13 26 H13 + nL
4 8 H4 + nL 14 28 H14 + nL
5 10 H5 + nL 15 30 H15 + nL
6 12 H6 + nL 16 32 H16 + nL
7 14 H7 + nL 17 34 H17 + nL
8 16 H8 + nL 18 36 H18 + nL
9 18 H9 + nL 19 38 H19 + nL
10 20 H10 + nL 20 40 H20 + nL
an = 2H K H 2 + 2 nO = 2 hHh + 2 nL
120
Chapter 3
Formula for an
h= H 3 an h= H 3 an
1 3 H1 + 2 nL 11 33 H11 + 2 nL
2 12 H1 + nL 12 72 H6 + nL
3 9 H3 + 2 nL 13 39 H13 + 2 nL
4 24 H2 + nL 14 84 H7 + nL
5 15 H5 + 2 nL 15 45 H15 + 2 nL
6 36 H3 + nL 16 96 H8 + nL
7 21 H7 + 2 nL 17 51 H17 + 2 nL
8 48 H4 + nL 18 108 H9 + nL
9 27 H9 + 2 nL 19 57 H19 + 2 nL
10 60 H5 + nL 20 120 H10 + nL
If we manipulate each formula so that the factor Hh + 2 nL appears, then we find that a similar pattern holds:
an = 3H K H 2 + 2 nO = 2 hHh + 2 nL
The reader is encouraged to verify this pattern for higher multiples of perfect squres. Thus, we conclude with the following
theorem:
Theorem: The following formulas generate Pythagorean triples 8a, b, c< of height H=c-b:
a) If H is NOT a multiple of a perfect square, then
an = :
2 hHn + 1L if H = 2 h even
H2 h + 1L H2 n + 1L if H = 2 h + 1 odd
=:
a2n - H 2 h nHn + 2L if H = 2 h even
2 H2 h + 1L nHn + 1L if H = 2 h + 1 odd
bn = (3.29)
2H
cn = bn + H = :
hIn2 + 2 n + 2M if H = 2 h even
H2 h + 1L I2 n + 2 n + 1M if H = 2 h + 1 odd
2
Step 2
Let’s investigate whether Pythagorean triples satisfy any recurrences. It suffices to consider the different cases as in Step 1
depending on whether the height H equals a multiple of a perfect square. Let’s begin with H an even integer, but not a multiple
of a perfect square:
121
Mathematics by Experiment
By examining the sums an + bn , it appears that the values for bn in the tables above satisfy the following recurrence:
bn+1 = an + bn + H 2
Does a similar recurrence hold for H an odd integer, but not a multiple of a perfect square? Let’s analyze the tables below:
nMax = 500;
Hmax = 5;
dataPythagoreantriplesheightsumab =
Map@Take@DeleteCases@Table@8ð, a, sol@@1, 1, 2DD . H ® ð,
Simplify@sol@@1, 1, 2DD + HD . H ® ð<, 8a, 1, nMax<D,
x__ ; H! HIntegerQ@x@@3DDD && x@@3DD > 0LLD, 5D &, 83, 5, 7, 11, 13<D;
Row@Table@ColumnDataDisplay@Table@ReplacePart@
dataPythagoreantriplesheightsumab@@HDD@@nDD, 1 -> nD,
8n, 1, Length@dataPythagoreantriplesheightsumab@@HDDD<D,
10, 8"n", "an ", "bn ", "cn "<, "
H=c-b=" <> ToString@dataPythagoreantriplesheightsumab@@HDD@@1, 1DDD <> "",
88Left, Left, Left, Left<, Automatic<D,
8H, 1, Length@dataPythagoreantriplesheightsumabD<D, " "D
We now move on to the case where H is a multiple of a perfect square. Let’s begin with H = h2 :
122
Chapter 3
It appears that for H = 1, we have bn+1 = 2 an + bn + 2, but for H = 4, we have bn+1 = an + bn + 2. No recurrence seems to exist
for H = 9, 16, 25. However, in these cases, if we restrict to triples that are spaced H rows apart, then an+k - an is a multiple
of H and the pattern holds. For example, if H = 9, then every third triple satisfies the recurrence bn+3 = 2 an + bn + 18, for
example,
8a1 , b1 , c1 <, 8a4 , b4 , c4 <, 8a7 , b7 , c7 <, ...
On the other hand, if H = 16, then every fourth triple satisfies the recurrence bn+4 = an + bn + 8, for example,
123
Mathematics by Experiment
3.4 Permutations
Recall that a permutation of a set is a just an ordering of the set. The following experiments reveal the wide range of properties
of permutations.
c@n_, k_D :=
Count@n - HammingDistance@ð, Range@1, nDD & Permutations@Range@1, nDD, kD
nMax = 5;
Table@c@n, kD, 8n, 0, nMax<, 8k, 0, n<D Grid
1
0 1
1 0 1
2 3 0 1
9 8 6 0 1
44 45 20 10 0 1
We immediately observe that cHn, n - 1L = 0 and cHn, nL = 1. The latter identity is clear: there is only one permutation, namely
the identity permutation, in which every student catches his or her own cap. The former identity can be explained as follows: if
Hn - 1L students catch their own caps, then the remaining student must also catch his or her own cap since it is the only cap
remaining; thus, cHn, n - 1L = 0.
124
Chapter 3
We immediately observe that cHn, n - 1L = 0 and cHn, nL = 1. The latter identity is clear: there is only one permutation, namely
the identity permutation, in which every student catches his or her own cap. The former identity can be explained as follows: if
Hn - 1L students catch their own caps, then the remaining student must also catch his or her own cap since it is the only cap
remaining; thus, cHn, n - 1L = 0.
Observe that since the number of permutations equals n ! and grows rapidly, a brute force computation of cHn, kL is impractical
here for large n, as demonstrated by the amount of time (in seconds) required to compute cHn, 0L for n = 1, 2, ... 10:
temp = Table@8n, Timing@c@n, 0DD<, 8n, 1, 10<D Grid
80., 0<
80., 1<
1
80., 2<
2
80., 9<
3
80., 44<
4
80.015, 265<
5
80.11, 1854<
6
80.343, 14 833<
7
Thus it is useful to find a more efficient formula for cHn, kL. Let’s start by trying to find a pattern for the values cHn, 0L in the first
column:
datahn0 = Table@c@n, 0D, 8n, 1, 8<D
80, 1, 2, 9, 44, 265, 1854, 14 833<
FindSequenceFunction@datahn0, nD
Subfactorial@nD
This shows that cHn, 0L represents the number of permutations where no student retrieves his or her own cap, i.e., the number of
derangements (discussed in Chapter 2). Thus,
cHn, 0L = ni (3.32)
where ni is the subfactorial function. This leads us to suspect that cHn, kL should involve the subfactorial Hn - kLi since in this
case Hn - kL students will NOT have retrieved his or her own cap. Thus, we consider values of cHn, kL divided by Hn - kLi :
Quiet@Table@c@n, kD Subfactorial@n - kD, 8n, 0, 5<, 8k, 0, n<D GridD
1
Indeterminate 1
1 Indeterminate 1
1 3 Indeterminate 1
1 4 6 Indeterminate 1
1 5 10 10 Indeterminate 1
We quickly recognize this table as Pascal’s triangle consisting of binomial coefficients (the Indeterminate values resulted from
division of cHn, n - 1L by Hn - Hn - 1LLi = 1i = 0). Thus, the formula for cHn, kL in terms of the factorial function becomes clear:
n n
cHn, kL = cHn - k, 0L = cHn - k, 0L (3.33)
k k
NOTE: We can of course have reasoned further by viewing cHn, kL as the number of ways of choosing k students who end up
retrieving their own cap multiplied with the number of ways in which the remaining Hn - kL students who do NOT retrieve their
own cap, i.e., cHn - k, 0L. Thus,
Hn - kLi
n
cHn, kL = (3.34)
k
125
Mathematics by Experiment
Thus, the probability that k students will catch their own caps equals cHn, kL n !.
cHn, kL n Hn - kLi n ! Hn - kLi Hn - kLi
k ! Hn - kL ! n ! k ! Hn - kL !
= = = (3.35)
n! k n!
To answer the original problem, the probability that exactly 5 out of 20 students will catch their own caps equals
cH20, 5L H15Li 48 066 515 734
= = » 0.003,
20 ! 5 ! × 15 ! 120 × 1 307 674 368 000
which is quite remote.
NOTE: Much more likely is when no student catches his or her own cap. The probability of this occuring out of 20 students
equals
cH20, 0L H20Li 895 014 631 192 902 121
= = » 0.368
20 ! 20 ! 2 432 902 008 176 640 000
or more than a third chance.
FURTHER EXPLORATION: Can you determine the limiting value of cHn, kL n ! as n ® ¥? If necessary, use the in-house
Mathematica command ISC or go directly to the Inverse Symbolic Calculator (ISC) website.
3.4.2 Runs
ACP Vol. 3, p.34
An ascending run of a permutation is an increasing contiguous subsequence of the permutation that cannot be extended at either
end. For example, the permutation 82, 4, 1, 3, 5< contains the run 81, 3, 5<. To find all runs, we use the Mathematica comm-
mand Runs.
Runs@82, 4, 1, 3, 5<D
882, 4<, 81, 3, 5<<
Thus, 82, 4, 1, 3, 5< has two runs, 82, 4< and 81, 3, 5<. Note that runs specify a partition for the permutation.
Denote by [ _ to be the number of permutations of length n that have k runs each. For example, here is a table listing the runs
n
k
for each permutation of length 3:
Runs for Permutations of Length 3
Σ Runs of Σ
81, 881, 2, 3<<
81, 2< 881, 3<, 82<<
2, 3<
126
Chapter 3
datanumberofruns = Table@
8k, Count@Runs@ðD & Permutations@81, 2, 3<D, t_ ; Length@tD kD<, 8k, 1, 3<D;
ColumnDataDisplayBdatanumberofruns, 10, :"k", "X \">, "Distribution of Runs"F
3
k
Distribution of Runs
k X \
3
k
1 1
2 4
3 1
X \ k=1 k=2
n
k=3 k=4 k=5 k=6 k=7
k
n=1 1 0 0 0 0 0 0
n=2 1 1 0 0 0 0 0
n=3 1 4 1 0 0 0 0
n=4 1 11 11 1 0 0 0
n=5 1 26 66 26 1 0 0
n=6 1 57 302 302 57 1 0
n=7 1 120 1191 2416 1191 120 1
Using Mathematica, we find that the second column satisfies the formula
FindSequenceFunction@80, 1, 4, 11, 26, 57, 120<, nD
- 1 + 2n - n
127
replaced ascending runs by descending runs. However, it is also possible to consider ascending and descending runs that alter-
nate in a given permutation as follows. Let Σ = 8a1 , a2 , ..., an < be a permutation of 81, 2, ..., n<. We define the first run to begin
at a1 and prescribe it to be ascending or descending based on whether a1 < a2 or a1 > a2 , respectively. Suppose the first run ends
Mathematics by Experiment
at ai before changing its climb, i.e., changes from ascending to descending or vice versa. Then ai becomes the start of the second
run and its ascension or descension is opposite that of the first run. By continuing this process down the last element of Σ, we
obtain a sequence of alternating runs. For example, the permutation 82, 5, 3, 1, 4< has three alternating runs: 82, 5<, 85, 3, 1<, and
81, 4<. The in-house Mathematica command AlternatingRuns will generate alternating runs of a permutation.
AlternatingRuns@82, 5, 3, 1, 4<D
882, 5<, 85, 3, 1<, 81, 4<<
Let’s now generate a table of alternating runs for all permutations of length 3.
Alternating Runs for Permutations of Length 3
Σ Alternating Runs of Σ Number of Alternating Runs
81, 881, 2, 3<<
81, 881, 3<, 83, 2<<
2, 3< 1
To investigate further, denote by [[ __ to be the number of permutations of length n which have exactly k alternating runs. For
n
k
example, we see from the table above that [[ __ = 2 and [[ __ = 4. Here is an array of values for [[ __.
3 3 n
1 2 k
Can you find any patterns in the array above for [[ __?
n
k
3.4.4 Involutions
ACP Vol. 3, p. 65
A permutation Σ of length n can also be thought of as a mapping of the set 81, 2, ..., n<. For example, if Σ = 82, 3, 1<, then we
can view it as the function
ΣH1L = 2
ΣH2L = 3
ΣH3L = 1
More generally, if Σ = 8k1 , k2 , ..., kn <, then ΣHiL = ki , i.e., the integer i is mapped to ki . Every permutation has an inverse,
denoted by Σ-1 , which maps ki back to i. For example, if Σ = 82, 3, 1<, then Σ-1 = 83, 1, 2< since we require
Σ-1 H1L = 3
Σ-1 H2L = 1
Σ-1 H3L = 2
Thus, Σ-1 satisfies the properties Σ-1 HΣHiLL = i and ΣIΣ-1 HiLM = i.
An involution is a permutation which is equal to its inverse, namely Σ = Σ-1 . For example, Σ = 83, 2, 1< is an involution. We
confirm this using the Mathematica command InversePermutation.
128
Chapter 3
InversePermutation@83, 2, 1<D
83, 2, 1<
Let’s investigate the numer of involutions of a given length. Towards this end, denote by IHnL to be the number of involutions of
length n. Here is a table of values of IHnL:
Number of Involutions of Length n
n IHnL
1 1
2 2
3 4
4 10
5 26
6 76
7 232
8 764
3.5 Partitions
Let’s consider the partition congruences, i.e., the congruence of pHnL mod q, where q is a positive integer. The table below gives
congruences for the first thirty values of pHnL mod 2.
Partition Congruences
n pHnL pHnL mod 2 n pHnL pHnL mod 2 n pHnL pHnL mod 2
1 1 1 11 56 0 21 792 0
2 2 0 12 77 1 22 1002 0
3 3 1 13 101 1 23 1255 1
4 5 1 14 135 1 24 1575 1
5 7 1 15 176 0 25 1958 0
6 11 1 16 231 1 26 2436 0
7 15 1 17 297 1 27 3010 0
8 22 0 18 385 1 28 3718 0
9 30 0 19 490 0 29 4565 1
10 42 0 20 627 1 30 5604 0
129
Mathematics by Experiment
Unfortunately, no pattern is evident. A similar dead-end arises if we consider congruences of pHnL mod q for q = 3 and q = 4.
However, when q = 5, we find an interesting pattern:
Partition Congruences
n pHnL pHnL mod 2 n pHnL pHnL mod 2 n pHnL pHnL mod 2
1 1 1 11 56 1 21 792 2
2 2 2 12 77 2 22 1002 2
3 3 3 13 101 1 23 1255 0
4 5 0 14 135 0 24 1575 0
5 7 2 15 176 1 25 1958 3
6 11 1 16 231 1 26 2436 1
7 15 0 17 297 2 27 3010 0
8 22 2 18 385 0 28 3718 3
9 30 0 19 490 0 29 4565 0
10 42 2 20 627 2 30 5604 4
Observe that the table above includes those integers n ending in 4 or 9, i.e., integers of the form 5 n + 4. Thus, we’ve discovered
the first of Ramanujan’s congruences for the partition function.
pH5 n + 4L º 0 mod 5 (3.36)
FURTHER INVESTIGATION: There are two other Ramanujan congruences. See if you can discover them by testing different
values of the modulus q.
3.5.1.2 Number of Smallest Parts
http://www.math.psu.edu/vstein/alg/antheory/preprint/andrews/17.pdf
The total number of smallest parts appearing in all the partitions of a positive integer n is defined to be sptHnL. For example, the
partitions of n = 4 are:
Replace@ð, ð@@Length@ðDDD ® Style@ð@@Length@ðDDD, UnderlinedD, 2D &
IntegerPartitions@4D
884<, 83, 1<, 82, 2<, 82, 1, 1<, 81, 1, 1, 1<<
The smallest parts in each partition are underlined in the above output. The table below lists the number of smallest parts for
each partition:
130
Chapter 3
84<
Partition Number of Smallest Parts
83, 1<
1
82, 2<
1
82, 1, 1<
2
81, 1, 1, 1<
2
4
Thus, we have sptH4L = 10. Let’s define spt as the command spt in Mathematica:
Clear@sptD
spt@n_D := Total@Count@ð, ð@@Length@ðDDDD & IntegerPartitions@nDD
The following table lists the first 30 values of spt.
We follow the trail blazed in the previous part by considering congruences of sptHnL mod 5.
We find that the same congruence holds for sptHnL as it does for the partition function, namely
sptH5 n + 4L º 0 mod 5 (3.37)
FURTHER EXPLORATION: Find two other congruence relations for sptHnL.
3.5.1.3 Palindromic Compositions
http://www.fq.math.ca/Scanned/41-3/heubach.pdf
A composition of n is an ordered sequence of positive integers whose sum is n. Thus, a composition is an ordered partition where
the order of the terms is taken into account. For example, there are five partitions of n = 4.
131
Mathematics by Experiment
IntegerPartitions@4D
884<, 83, 1<, 82, 2<, 82, 1, 1<, 81, 1, 1, 1<<
A palindromic composition is a composition that reads the same forwards and backwards (palindrome). For example, n = 4 has
three palindromic compositions:
4
1+2+1
1+1+1+1
Let’s investigate the number of palindromic compositions of an arbitrary positive integer n.
n Bn mod 2 n Bn mod 2
1 1 11 0
2 0 12 1
3 1 13 1
4 1 14 0
5 0 15 1
6 1 16 1
7 1 17 0
8 0 18 1
9 1 19 1
10 1 20 0
The pattern is clear: the residues above are periodic with period 3, i.e.,
BHn + 3L º BHnL mod 2 (3.38)
Let’s continue this trail by considering congruences of Bn mod 3.
132
Chapter 3
Here, the residues seem to have period 13 over a range of 2 periods. We should experimentally verify this over a wider range.
Towards this end, we use our in-house Mathematica command SequencePeriod to experimentally determine whether this period
holds for the first 100 terms:
? SequencePeriod
The reader can confirm this over an even wider range, say the first 1000 terms. Therefore, we’ve discovered the congruence
relation
BHn + 13L º BHnL mod 3 (3.39)
What about other values for q? Is the sequence of congruences Bn mod 3 always periodic for every q? Define Πq to be the period
of Bn mod 3 q. Below is a table of values for Πq for q ranging from 1 to 6:
q Πq
1 1
2 3
3 13
4 12
5 781
6 39
Observe that the periods in the table above are relatively small except for q = 5, which jumps to a period 781.
dataperiodBellnumbers = 81, 3, 13, 12, 781, 39, 137 257, 24,
39, 2343, 28 531 167 061, 156, 25 239 592 216 021, 411 771, 10 153, 48,
51 702 516 367 896 047 761, 117, 109 912 203 092 239 643 840 221, 9372, 1 784 341,
85 593 501 183, 949 112 181 811 268 728 834 319 677 753, 312, 3905, 117<
81, 3, 13, 12, 781, 39, 137 257, 24, 39, 2343, 28 531 167 061, 156, 25 239 592 216 021,
411 771, 10 153, 48, 51 702 516 367 896 047 761, 117, 109 912 203 092 239 643 840 221,
9372, 1 784 341, 85 593 501 183, 949 112 181 811 268 728 834 319 677 753, 312, 3905, 117<
133
Mathematics by Experiment
FURTHER EXPLORATION:
1. Make a table of periods for Bn mod 3 for q ranging from 1 to 20. NOTE: Some of these periods will be extremely large.
2. Do you see a pattern for those values of q where the period jumps to a relatively large value (in comparison to its immediate
neighbors)?
3. Can you find a formula for the relatively large periods mentioned in part 2.
3.6 Hyper-Polyhedra
In this section we investigate patterns involving the number of vertices, edges, and faces of higher-dimensional polyhedra.
Euler’s polyhedron formula describes how V , E, and F are related for all convex polyhedra:
V -E+F=2 (3.40)
3.6.2 Hypertetrahedron
A hypertetrahedron (also called a simplex or pentatope) is generalization of a tetrahedron to four dimensions (4-D). More
generally, an n-tetrahedron is a generalization of a tetrahedron to n dimensions. An n-dimensional unit hypertrahedron is defined
to be the object obtained by inserting a vertex along the n-th dimension and forming edges of length 1 between it and all vertices
of an Hn - 1L-dimensional unit hypertetrahedron. Essentially an n-dimensional hypertetrahedron can be represented by the
complete graph Kn on Hn + 1L vertices where any two vertices are connected by an edge.
Step 1
Let's determine the number of vertices and edges for a hypertetrahedron:
134
Chapter 3
Hypertetrahedra
Dimension ð Vertices ð Edges
0 1 0
1 2 1
2 3 3
3 4 6
4 5 10
The pattern for the number of vertices is clear. Do you recognize the values for number of edges?
3.6.3 Hypercube
A hypercube is generalization of a cube to four dimensions (4-D). More generally, an n-cube is a generalization of a cube to n
dimensions. An n-dimensional unit hypercube is defined to be the object obtained by taking two copies of an Hn - 1L-dimensional
unit hypercube, parallel to each other along the n-th dimension, and forming edges of length 1 between corresponding vertices of
the two copies.
Step 2
The following table lists the number of vertices and edges of an n-dimensional hypercube for n = 1, 2, 3, 4:
Hypercube
Dimension ð Vertices ð Edges
0 1 0
1 2 1
2 4 4
3 8 12
4 16 32
5 32 80
Do you recognize the pattern for the number of vertices and edges? What about faces of a n-cube? See if you can find a formula
for the number of faces.
http://www-history.mcs.st-and.ac.uk/Mathematicians/Newton.html
135
Mathematics by Experiment
Binomial Theorem:
Newton's Discovery of the General Binomial Theorem
D. T. Whiteside
The Mathematical Gazette
Vol. 45, No. 353 (Oct., 1961), pp. 175-180
(article consists of 6 pages)
Published by: The Mathematical Association
Stable URL: http://www.jstor.org/stable/3612767
Newton, having studied Wallis’ Arithmetica Infinitorum, knew that Wallis had interpolated the definite integral Ù0 I1 - x12 M â x
1 2
that led to his infinite product representation of Π. By following the same trail, Newton was able to generalize the Binomial
He originally succeeded in obtaining infinite series expressions for Ù0 Ia2 - x2 M
x 12
Theorem to non-integer exponents. â x,
Ù0 Ia2 + x2 M â x, and Ù0 a2 Hb + xL-1 â x, representing the areas the circle y2 = a2 - x2 and hyperbolas yH b + xL = a2 , respectively.
x 12 x
As he explains it in [], Newton in each case reduced the problem to interpolating the integrand, say a2 Hb + xL-1 , from the
sequence of expansions a2 Hb + xLn for non-negative integers n, whose coefficients were already known to be given by Pascal’s
triangle:
Binomial Expansion of a2 Hb+xLn
n a2 Hb+xLn
0 a2
1 a2 b + a2 x
2 a2 b2 + 2 a2 b x + a2 x2
3 a2 b3 + 3 a2 b2 x + 3 a2 b x2 + a2 x3
4 a2 b4 + 4 a2 b3 x + 6 a2 b2 x2 + 4 a2 b x3 + a2 x4
5 a2 b5 + 5 a2 b4 x + 10 a2 b3 x2 + 10 a2 b2 x3 + 5 a2 b x4 + a2 x5
n
Now, by creating a table of these coefficients, called binomial coefficients (see Chapter 2) and denoted by to refer to the
k
coefficient of a2 bn-k xk in the expansion of a2 Hb + xLn , Newton was able to extend the recursive pattern for binomial coefficients
to negative integer exponents:
Binomial Coefficients H L
n
k
k=0 k=1 k=2 k=3 k=4 k=5
n=-5 1 - 5 15 - 35 70 - 126
n=-4 1 - 4 10 - 20 35 - 56
n=-3 1 -3 6 - 10 15 - 21
n=-2 1 -2 3 -4 5 -6
n=-1 1 -1 1 -1 1 -1
n=0 1 0 0 0 0 0
n=1 1 1 0 0 0 0
n=2 1 2 1 0 0 0
n=3 1 3 3 1 0 0
n=4 1 4 6 4 1 0
n=5 1 5 10 10 5 1
As a result, this led him to the series expansion for a2 Hb + xL-1 corresponding to n = -1:
a2 Hb + xL-1 =
a2 a2 a2 a2
- x+ x2 - x3 + ... (3.41)
b b2 b3 b4
136
Chapter 3
Of course, the more difficult problem was to interpolate the binomial coefficients of a2 Hb + xLn for fractional exponents n. To
achieve this, Newton made the spectacular observation that the values in each column followed a linear progression of the form
Coefficients for xk in expansion of a2 Hb+xLn
k=0 k=1 k=2 k=3 k=4 k=5
n=-5 a b - 5 c d - 5 e + 15 f g - 5 h + 15 i - 35 j k - 5 l + 15 m - p - 5 q + 15 r -
35 n + 70 o 35 s + 70 t - 126 u
n=-4 a b - 4 c d - 4 e + 10 f g - 4 h + 10 i - 20 j k - 4 l + 10 m - p - 4 q + 10 r -
20 n + 35 o 20 s + 35 t - 56 u
n=-3 a b-3c d-3e+6f g - 3 h + 6 i - 10 j k-3l+6m- p-3q+6r-
10 n + 15 o 10 s + 15 t - 21 u
n=-2 a b-2c d-2e+3f g-2h+3i-4j k-2l+ p-2q+3r-
3m-4n+5o 4s+5t-6u
n=-1 a b-c d-e+f g-h+i-j k-l+m-n+o p-q+r-s+t-u
n=0 a b d g k p
n=1 a b+c d+e g+h k+l p+q
n=2 a b+2c d+2e+f g+2h+i k+2l+m p+2q+r
n=3 a b+3c d+3e+3f g+3h+3i+j k+3l+3m+n p+3q+3r+s
n=4 a b+4c d+4e+6f g+4h+6i+4j k+4l+6m+4n+o p+4q+6r+4s+t
n=5 a b+5c d + 5 e + 10 f g + 5 h + 10 i + 10 j k + 5 l + 10 m + p + 5 q + 10 r +
10 n + 5 o 10 s + 5 t + u
where a = 1, b = 0, c = 1, d = 0, e = 0, f = 1, etc.
Newton then boldly assumed that a similar linear progression would continue to hold when values for half-integer exponents, i.e.,
n = m 2, were inserted into the Table? above (indicated by * since they are unknown for the moment):
Coefficients for xk in expansion of a2 Hb+xLn
k=0 k=1 k=2 k=3 k=4
n=-4 a b-8c d - 8 e + 36 f g - 8 h + 36 i - 120 j k - 8 l + 36 m - 120 n + 330 o
n=-72 a b-7c d - 7 e + 28 f g - 7 h + 28 i - 84 j k - 7 l + 28 m - 84 n + 210 o
n=-3 a b-6c d - 6 e + 21 f g - 6 h + 21 i - 56 j k - 6 l + 21 m - 56 n + 126 o
n=-52 a b-5c d - 5 e + 15 f g - 5 h + 15 i - 35 j k - 5 l + 15 m - 35 n + 70 o
n=-2 a b-4c d - 4 e + 10 f g - 4 h + 10 i - 20 j k - 4 l + 10 m - 20 n + 35 o
n=-32 a b-3c d-3e+6f g - 3 h + 6 i - 10 j k - 3 l + 6 m - 10 n + 15 o
n=-1 a b-2c d-2e+3f g-2h+3i-4j k-2l+3m-4n+5o
n=-12 a b-c d-e+f g-h+i-j k-l+m-n+o
n=0 a b d g k
n=12 a b+c d+e g+h k+l
n=1 a b+2c d+2e+f g+2h+i k+2l+m
n=32 a b+3c d+3e+3f g+3h+3i+j k+3l+3m+n
n=2 a b+4c d+4e+6f g+4h+6i+4j k+4l+6m+4n+o
n=52 a b+5c d + 5 e + 10 f g + 5 h + 10 i + 10 j k + 5 l + 10 m + 10 n + 5 o
n=3 a b+6c d + 6 e + 15 f g + 6 h + 15 i + 20 j k + 6 l + 15 m + 20 n + 15 o
n=72 a b+7c d + 7 e + 21 f g + 7 h + 21 i + 35 j k + 7 l + 21 m + 35 n + 35 o
n=4 a b+8c d + 8 e + 28 f g + 8 h + 28 i + 56 j k + 8 l + 28 m + 56 n + 70 o
Since these values are known for those rows where n is an integer exponent, it suffices to solve the infinite family of equations
(assuming they are all consistent) for the variables a, b, c, d, .... For example, the equations derived from the third column are
... (3.42)
d - 4 e + 10 f = 3
d -2e+3 f =1
(3.43)
137
Mathematics by Experiment
d =0
d +2e+ f =0
d +4e+6 f =1
...
It follows that d = 0, e = -1 8, and f = 1 4. By applying this to every column, Newton was able to fill in his table for half-
integer exponents:
Coefficients for xk in expansion of a2 Hb+xLn
k=0 k=1 k=2 k=3 k=4
n=-4 1 - 4 10 - 20 35
7 63 231 3003
n=-72 1 -2 8
- 16 128
n=-3 1 -3 6 - 10 15
5 35 105 1155
n=-52 1 - 2 8
- 16 128
n=-2 1 -2 3 -4 5
3 15 35 315
n=-32 1 - 2 8
- 16 128
n=-1 1 -1 1 -1 1
1 3 5 35
n=-12 1 - 2 8
- 16 128
n=0 1 0 0 0 0
1 1 1 5
n=12 1 2
- 8 16
- 128
n=1 1 1 0 0 0
3 3 1 3
n=32 1 2 8
- 16 128
n=2 1 2 1 0 0
5 15 5 5
n=52 1 2 8 16
- 128
n=3 1 3 3 1 0
7 35 35 35
n=72 1 2 8 16 128
n=4 1 4 6 4 1
In particular,
Ia2 - x2 M
12 2 3 4
12 x2 1 x2 1 x2 1 x2 5 x2
=a 1- = aB1 + - + - + ...F (3.44)
a2 2 a2 8 a2 16 a2 128 a2
Lastly, it remains to find a general formula for generating binomial coefficients for all fractional exponents without having to
interpolate each row individually. This is where Newton must have been influenced by Wallis’ use of infinite products for he
wrote down the formula
pq 1 ´ p ´Hp - qL ´H p - 2 qL ´Hp - 3 qL ´ × × ×
= (3.45)
k 1 ´k ´2 k ´3 k ´4 k ´ × × ×
leading to the modern version of the Binomial Theorem, which holds for all real exponents n:
Theorem: If n is real-valued and |x|<1, then
nHn - 1L Hn - 2L nHn - 1L Hn - 2L Hn - 3L
H1 + xLn = 1 +
n nHn - 1L
x+ x2 + x3 + x4 + ... (3.46)
1 1×2 1×2×3 1×2×3×4
(3.43)
138
Chapter 3
139
Appendix
i
Bibliography
[Ba] Edward J. Barbeau, Pell's Equation, Problem Books in Mathematics, Springer-Verlag, 2003.
[BB] Jonathan Borwein and David Bailey, Mathematics by Experiment: Plausible Reasoning for the 21st Century, A K Peters,
Wellesley, MA, 2008.
[BBC] David Bailey, Jonathan M. Borwein, Neil J. Calkin, Roland Girgensohn, D. Russell Luke, and Victor H. Moll, Experimen-
tal Mathematics in Action, A K Peters, Wellesley, MA, 2007.
[BoCu] John P. Bonomo and Carolyn K. Cuff, How Do You Stack Up?, The College Mathematics Journal, Vol. 35, No. 5 (Nov.,
2004), pp. 351-361: http://www.jstor.org/stable/4146844.
[BBG] Jonathan Borwein, David Bailey and Roland Girgensohn, Experimentation in Mathematics: Computational Paths to
Discovery, A K Peters, Natick, MA, 2004.
[BBN] Arthur T. Benjamin, Curtis T. Bennett, Florence Newberger, Recounting the Odds of an Even Derangement, Mathematics
Magazine, Vol. 78, No. 5, Dec., 2005, pp. 387-390.
[BoKo] Douglas Bowman and Nasha Komanda, E3450, The American Mathematical Monthly, Vol. 100, No. 3 (Mar., 1993), p.
300: http://www.jstor.org/stable/2324473.
[Ch] Hongwei Chen, Excursions in Classical Analysis, Mathematical Association of America, 2010.
[DB] Keith Devlin and Jonathan Borwein, The Computer as Crucible: An Introduction to Experimental Mathematics, A K
Peters, Wellesley, MA, 2008.
[DeBr] M. N. Deshpande and Ezra Brown, Diophantine Triplets and the Pell Sequence, The Fibonacci Quarterly, Vol. 39, No. 3,
2001: http://www.fq.math.ca/Scanned/39-3/deshpande.pdf
[El] Niles Eldridge, The Pattern of Evolution, W. H. Freeman, 1998.
[Eu] Leonard Euler, Decouverte d’une loi tout extraordinaire des nombres par rapport a la somme de leurs diviseurs (Discovery
of a most extraordinary law of the numbers concerning the sum of their divisors), Bibliotheque impartiale 3, 1751, pp. 10-31
(Opera Omnia: Series 1, Volume 2, pp. 241 - 253), Eneström index: E175. Available on the Euler Archive: http://www.math.-
dartmouth.edu/~euler/.
[ES] Russ Euler and Jawad Sadek, Elementary Problems and Solutions (B951), The Fibonacci Quarterly, Vol. 41, No. 5, 2003.
http://www.fq.math.ca/Scanned/41-5/elementary41-5.pdf
[ES2] Russ Euler and Jawad Sadek, Elementary Problems and Solutions (B959), The Fibonacci Quarterly, Vol. 41, No. 2, 2003:
http://www.fq.math.ca/Scanned/41-2/elementary41-2.pdf
[ES3] Russ Euler and Jawad Sadek, Elementary Problems and Solutions (B959), The Fibonacci Quarterly, Vol. 42, No. 2, 2004:
http://www.fq.math.ca/Papers1/42-2/May2004elementary.pdf
[EwSc] J. Ewing and G. Schober, The area of the Mandelbrot set, Numerische Mathematik, Vol. 61 (1992), 59 - 72.
[Fi] Fibonacci, Liber Abaci. A Translation into Modern English of Leonardo Pisano’s Book of Calculuation, L. E. Sigler
(translator), Springer-Verlag, NY 2002.
[Ga] Martin Gardner, The Book of Colossal Mathematics, W. W. Norton, New York, 2001.
[Go] L. J. Goldstein, A History of the Prime Number Theorem, The American Mathematical Monthly, Vol. 80, No. 6 (Jun. - Jul.,
1973), pp. 599-615: http://www.jstor.org/stable/2319162
[Gr] Andrew Granville, Zaphod Beeblebrox’s Brian and the Fifty-ninth Row of Pascal’s Triangle, The American Mathematical
Monthly
Vol. 99, No. 4 (Apr., 1992), pp. 318-331: http://www.jstor.org/stable/2324898.
i
Mathematics by Experiment
[Gr] Andrew Granville, Zaphod Beeblebrox’s Brian and the Fifty-ninth Row of Pascal’s Triangle, The American Mathematical
Monthly
Vol. 99, No. 4 (Apr., 1992), pp. 318-331: http://www.jstor.org/stable/2324898.
[Gu] Richard Guy, Divisors and Desires, The American Mathematical Monthly, Vol. 104, No. 4 (Apr., 1997), pp. 359-360:
http://www.jstor.org/stable/2974586
[GKP] Ronald L. Graham, Donald E. Knuth, and Oren Patashnik, Concrete Mathematics, Addison-Wesley, 1994.
[Ha] Paul Halmos, I Want to be a Mathematician: An Automathography, Mathematical Association of America, 1985.
[Hay] Brian Hayes, Gauss’s Day of Reckoning, American Scientist, Volume 94, Number 3 (May-June 2006), p. 200.
[Ho] F. T. Howard, The Sum of the Squares of Two Generalized Fibonacci Numbers, The Fibonacci Quarterly, Vol. 41, No. 1,
2003: http://www.fq.math.ca/Scanned/41-1/howard.pdf.
[I] Inverse Symbolic Calculator: Wikipedia - the free encyclopedia: http://en.wikipedia.org/wiki/Inverse_symbolic_calculator
[Ki] Clark Kimberling, Fractal Sequences, http://faculty.evansville.edu/ck6/integer/fractals.html.
[Kn] Donald E. Knuth, The Art of Computer Programming: Volume 3 - Sorting and Searching, Addison-Wesley, 1973.
[La] L. J. Lange, An Elegant Continued Fraction for Π, The American Mathematical Monthly, Vol. 106, No. 5 (May, 1999), pp.
456-458: http://www.jstor.org/stable/2589152
[LP] Thomas M. Liggett and Peter Petersen, ,The Law of Large Numbers and $\sqrt 2$, The American Mathematical Monthly,
Vol. 102, No. 1 (Jan., 1995), pp. 31-35: http://www.jstor.org/stable/2974854
[LLHC] John F. Loase, David Lansing, Cassie Hryczaniuk and Jamie Cahoon, A Variant of the Partition Function, The College
Mathematics Journal, Vol. 36, No. 4 (Sep., 2005), pp. 320-321: http://www.jstor.org/stable/30044876
[LoHe] Calvin T. Long and Boyd Henry, Linear Relations Between Powers of Terms in Arithmetic Progression, The College
Mathematics Journal, Vol. 32, No. 2 (Mar., 2001), pp. 135-137: http://www.jstor.org/stable/2687121
[MT] Colin L. Mallows and John W. Tukey, An overview of techniques of data analysis, in Some Recent Advances in Statistics,
de Oliviera and Epstein (editors), Academic Press, London, 1982.
[Me] R. S. Melham, A Three-Variable Identity Involving Cubes of Fibonacci Numbers, The Fibonacci Quarterly, Vol. 41, No. 3,
2003: http://www.fq.math.ca/Scanned/41-3/melham.pdf.
[MP] G. Myerson and A. J. van der Poorten, Some Problems Concerning Recurrence Sequences, The American Mathematical
Monthly, Vol. 102, No. 8 (Oct., 1995), pp. 698-705: http://www.jstor.org/stable/2974639
[MW] Darryl McCullough and Elizabeth Wade, Recursive Enumeration of Pythagorean Triples, The College Mathematics
Journal Vol. 34, No. 2 (Mar., 2003), pp. 107-111: http://www.jstor.org/stable/3595782
[O] On-Line Encyclopedia of Integer Sequences: Wikipedia - the free encyclopedia: http://en.wikipedia.org/wiki/On-Line_Ency-
clopedia _of _Integer _Sequences
[Od] Odlyzko, A. M. (1993), Iterated absolute values of differences of consecutive primes, Mathematics of Computation 61:
373–380, http://www.dtc.umn.edu/~odlyzko/doc/arch/gilbreath.conj.ps
[Os] Thomas J. Osler, Binomials to Binomials, The College Mathematics Journal, Vol. 31, No. 3 (May, 2000), pp. 211-212:
http://www.jstor.org/stable/2687492
[Pa] Poo-Sung Park, Ramanujan’s Continued Fraction for a Puzzle, The College Mathematics Journal, Vol. 36, No. 5 (Nov.,
2005), pp. 363-365: http://www.jstor.org/stable/30044886
[PS] Problems and Solutions 653, The College Mathematics Journal, Vol. 31, No. 3 (May, 2000), pp. 219-226: http://www.js-
tor.org/stable/2687495.
[Pu] Jan-Christoph Puchta, On the Distribution of Totients, The Fibonacci Quarterly, Vol. 40, No. 1, 2002: http://www.fq.math.-
ca/Scanned/40-1/puchta2.pdf
[Wi-Fi] Fibonacci number: Wikipedia - the free encyclopedia: http://en.wikipedia.org/wiki/Fibonacci_number
[Wi-Gi] Gilbreath’s Conjecture: Wikipedia - the free encyclopedia: http://en.wikipedia.org/wiki/Gilbreath%27 s_conjecture
[Wi-Pe] Perfect Number: Wikipedia - the free encyclopedia: http://en.wikipedia.org/wiki/Perfect_number
[Po] George Polya, Mathematics and Plausible Reasoning, Volume 1: Induction and Analogy in Mathematics, Princeton Univer-
sity Press, 1954.
[ShSh] Harris S. Shultz and Ray C. Shiflett, M&m Sequences, The College Mathematics Journal, Vol. 36, No. 3 (May, 2005),
pp. 191-198: http://www.jstor.org/stable/30044851.
ii
Preface
[Sp] Michael Z. Spivey, The Humble Sum of Remainders Function, Mathematics Magazine, Vol. 78, No. 4, Oct., 2005, pp.
300-305: http://www.jstor.org/stable/30044176
[St1] Richard P. Stanley, Enumerative Combinatorics: Volume 1,
[St2] Richard P. Stanley, Enumerative Combinatorics: Volume 2,
[WW] Peter W. Wade and William R. Wade, Recursions That Produce Pythagorean Triples, The College Mathematics Journal,
Vol. 31, No. 2 (Mar., 2000), pp. 98-101: http://www.jstor.org/stable/2687578
[Yo-P] Paul T. Young, On Lacunary Recurrences, The Fibonacci Quarterly, Volume 41 (2003), No. 1, pp.41-47: http://www.fq.-
math.ca/Scanned/41-1/young.pdf
[Yo-T] Robert M. Young, Excursions in Calculus: An Interplay of the Continuous and the Discrete, The Mathematical Associa-
tion of America, 1992.
iii