(Discrete Mathematics Techniques III)
1st ed. Apr 2nd, 2010
Abstract
This is the last in the 3-part series of articles on finding for oneself the solution to the sum of integer power problem, and in the process discovering the Bernoulli numbers. In Part 3 (this paper), we find a direct closed-form solution, i.e. one that does not require iteration, for the general case of the finite-summation-of-integer-powers problem . Having established in Part 2 that the closed-form solution is a polynomial, the summation is here rewritten as the sum of the
independent monomials
(
), where the
are unknown coefficients. Using the recurrence relation
, we obtain a linear combination of the monomials, which reduces to an easily solvable
-by-
triangular linear system in the unknown coefficients
of the closed-form polynomial solution. Maxima and Octave/Matlab codes for directly computing the closed-form solutions are included in the Appendices.
A lovely paper by Bearden (March 1996, American Mathematical Monthly), which was shared with me by a reader, tells the mathematical story nicely, with much of the history filled in.
NOTE: This paper contains many formulas. As such, you may prefer reading the paper as a typeset PDF.
Finding a Closed-Form solution for
using a Direct Method
Motivation
Our goal is to obtain, without using iteration, a closed form solution for the general case of the finite-summation-of-integer-powers problem:
Recall that in Part 2 of this paper , we used a method of recurrence relations to obtain the following results about :
-
may be expressed as a linear recurrence relation that uses the
lower order formulas
and
as follows:
.
Thein (2) are the coefficients from the
-st polynomial solution
, while the terms
are defined in terms of binomial coefficients1 as:
- The closed form solution of
is a polynomial in
of degree
with rational coefficients and a constant term equal to zero.
While the expression for in (2) is indeed accurate, it requires repeated iteration to obtain the closed form expression for any particular
. Such an approach without the assistance of a computer algebra system such as Maxima2 would be prohibitively time-consuming and prone to error.
We are hence left with the following question:
Can we find a closed-form solution for the coefficients of the general solution polynomial that can be obtained directly, i.e. that does not require iteration?
This paper uses linear algebra and matrices to achieve precisely this goal.
The Direct Approach
Polynomial with Undetermined Coefficients
Since we have established in (2) that is a polynomial of order
with no constant term, we may write down the closed form solution of
in the form:
where , are coefficients that have yet to be determined. As a first step towards determining these coefficients, observe that every finite summation can be written as a first order recurrence in
by peeling off the last summand:
Substituting (3) into (4) gives:
Summation Manipulations
What follows is a sequence of summation manipulations3 and simplifications of (5) aimed at obtaining an expression from which a closed-form solution for the coefficicents becomes transparent:
- Expand the binomial powers using the binomial formula on both the left-hand side (LHS) and right-hand side (RHS) of (5):
- Modify the LHS and RHS of (6) as follows, exploiting the fact that
:
- LHS: Extend upper range of inner sum to match upper range of outer sum.
- RHS: Peel off 0-index term and add a vacuous
-index term to last sum.
- Result:
-
- LHS: Interchange order of summations.
- RHS: Combine terms.
- Result:
-
- LHS: Peel off 0-index term from the (now outside) summation on
.
- RHS: Change dummy index variable.
- Result:
- LHS: Peel off 0-index term from the (now outside) summation on
- Bring everything over to the left-hand side to obtain a homogeneous linear combination of the monomials
:
Linear Independence Reveals Equations for Unknown Coefficients
Since the monomials are a basis for the vector space
of polynomials of degree less than or equal to
, any linear combination of the
will equal zero if and only if each coefficient equals zero.
Therefore, setting the coefficients from (7) equal to zero gives us a system of simultaneous linear equations in the
unknown coefficients
, one equation for each of the coefficients of
:
Observe that (*) reduces to , since
for all
, and
when
. Thus we are left with a square system of
equations (
) in
unknowns
, (
).
By observing that for
, we can combine (8) and (9) into the single set of linear equations:
Taking into account the inequality , and noting once again that
when
, (10) simplifies to a triangular linear system:
and hence all equations in the system are linearly independent.
From Linear System to Matrix Equation
The triangular linear system (11) may be readily solved by combining the summations into a single matrix equation and using any one of a number of matrix solver packages. But this requires that the system (11) be re-indexed4 so that the start value of index
matches that of index
. Following the usual matrix convention, referred to as 1-indexing, we choose the start indices to be
.
Maxima code for solving the 0-indexed5 triangular linear system (13) is given in Appendix B. Octave/Matlab code for solving the 1-indexed triangular linear system (12) is given in Appendix C.
The matrix equation equivalent to (12) is: , where
Since is an upper triangular matrix with all non-zero upper triangular entries, we know we can readily solve this for any given values of
and
using a matrix solver system such as Maxima/Mathematica, Octave/Matlab, or Maple, among others. Source code for solutions in Maxima and Octave/Matlab are given in Appendix B and Appendix C, respectively.
Computing the General, Closed-Form Solution
We know from (2) that the general, closed-form solution to (1) is a -degree polynomial (3) in
with rational coefficients
given by solving the triangular linear system (12), or equivalently (14).
The highest few coefficients can be readily computed by hand. For all
, we have:
-
-
-
-
-
-
-
.
- In particular, we claim (without proof) that
are rational coefficients independent of
, and
for all even positive
. Hence we have:
Substituting the above coefficients into the general closed-form polynomial solution gives us:
A listing of the solution formulas for the first ten values is given in Appendix A.
Conclusion
Reviewing the path we have taken in this 3-part paper: Part 1 motivated the problem and illustrated a recurrence relation method for obtaining the solution for small . Part 2 generalized this method and found a
-th order recurrence relation for the solution. Part 2 also illustrated, by an induction argument, that the closed form solutions for all
are polynomials of degree
with rational coefficients and no constant term.
This paper used the closed-form polynomial expression motivated in Part 2 to obtain a direct solution. By writing out the polynomial form with undetermined coefficients, the recurrence relation (4) was manipulated into a linear combination of polynomials (7). Using a linear independence argument, (7) was reduced to a triangular linear system (12) and associated matrix equation (14).The general closed-form solution is (16).
Solving the finite-summation-of-integer-powers problem , for arbitrary positive integers
and
has provided a natural setting to use a variety of techniques from discrete mathematics and linear algebra, and poses additional interesting questions (characterization of the denominators in the coefficients of (16) (the
in (15)) and divisibility properties of
for various
). In particular, we have obtained a direct method for solving the finite-summation-of-integer-powers problem that invokes a matrix solution and does not require iteration.
Appendix A: Listing of Solutions to (1)
A listing of the solution formulas for the first ten values is as follows:
Appendix B: Maxima Source Code for Solving (14)
For particular , the coefficients
for the closed form solutions
are immediate. Maxima code, is given below, for determining this solution using the 0-indexed triangular linear system (13):
solution(p):= block([a, eq], /* give subroutine variables local scope */ v : makelist(a[i], i, 0, p), /* create list of unknowns (0-indexed) */ /* create list of equations (0-indexed) */ eq : makelist(sum(binom(j+1,i)*a[j],j,i,p) = binom(p,i), i, 0, p), linsolve(eq, v) )$
A more elaborate function that computes the series form and factored forms of the solution is given below. This was used to generate the first ten solution listings.
SpN_mat(p):= block([a, eq], /* give subroutine variables local scope */ v : makelist(a[i], i, 0, p), /* create list of unknowns (0-indexed) */ /* create list of equations (0-indexed) */ eq : makelist(sum(binom(j+1,i)*a[j],j,i,p) = binom(p,i), i, 0, p), /* find coefficients of solution polynomial by solving linear system */ sol : linsolve(eq, v), /* create polynomial: inner product of {N^i} with coefficients */ pol : makelist(N^(i+1),i,0,p), /* monomials {N^i} */ /* closed form formula */ cff : rhs(sol.pol), /* inner product of {N^i} with coefficients */ cff /* return closed form formula in series */ )$
Download the complete Maxima source code listing here.
Appendix C: Octave/Matlab Source Code for Solving (14)
For particular , the coefficients
for the closed form solutions
are immediate. Octave/Matlab code, is given below, for determining this solution using the 1-indexed triangular linear system (12):
% Functions listing coefficients from a_p+1 to a_1 as fractions function c = sumkp_matrix(p) M=zeros(p+1); % initialize matrix M for i=1:p+1 for j=i:p+1 % set upper triangular elements M(i,j)=nchoosek(j,i-1); % Equation (12) end; end; b=zeros(p+1,1); % initialize column matrix b for i=1:p+1 b(i)=nchoosek(p,i-1); % Equation (12) end c=inv(M)*b; % Solve to obtain coefficients c=fliplr(c'); % List from highest index to lowest sml=abs(c)<1e-8; % Numerical correction: set tiny coeffients to exactly zero c(find(sml==1))=0; disp("Solution Coefficients (High to Low order):") c=rats(c); % return coefficients as fractions end;
Download the complete Octave source code listing here.
This article is available for download as a PDF here.
The Maxima and Octave complete source code listings implementing the direct matrix solution are available for download here.
References
- [Ebr10a]
-
Assad Ebrahim.
Finite summation of integer powers xp, part 1.
January 2010. - [EO10b]
-
Assad Ebrahim.
Finite summation of integer powers xp, part 2.
February 2010. - [GKP]
-
Ron Graham, Donald Knuth, and Oren Patashnik.
Concrete Mathematics: A Foundation for Computer Science.
Addison Wesley. - [Knu]
-
Donald Knuth.
The Art of Computer Programming (3 Volumes).
Footnotes
- The notation
used above denotes binomial coefficients, often expressed verbally as "n choose k". Other representations of these coefficients include
and
. ↩
- Maxima code for automatically computing solutions to (1) using solution formula (2) is given in Part 2 of this paper . ↩
- Such manipulations are taught in material such as (GKP) and (Knu). ↩
- Index manipulation is covered in (GKP) and (Knu). ↩
- 0-indexing is more convenient for various computational software packages and programming langauges, including Maxima. The 0-indexed system is:
[…] If you’re not sure what LaTeX / TeX is, start with Introduction to LaTeX. If you already have a TeX platform set-up, skip to Modularity in TeX. To see an equation heavy paper on an elementary topic (accessible to undergraduates), check out Finite Summation of Integer Powers (Part 3). […]
[…] from this site are: Finite Summations of Integers (1-Mathematical, 2-Maxima, 3-Octave), Maxima for Symbolic Computation, Mathematics of Duelling (R […]