Math 236  •  12 April 2001

Factorial function: an introduction to Mathematica programming paradigms

Introduction

Recall the mathematical definition of the factorial function:  For a positive integer n, the value n! is the product 1 · 2 · 
        3 · ⋯ · n of the first n positive integers, and 0!  =  1.  In particular, 1!  =  1, 2!  =  2, 3!  =  6, and 4!  =  24.

The purpose of this notebook is to show various ways of expressing the factorial function in Mathematica.  These ways can serve as introductions to some of Mathematica's programming paradigms.

Mathematica's Factorial function

The factorial function is built into Mathematica and, in its abbreviated input form, uses the very same notation as in mathematics itself:

4!

The form n! is just an abbreviation for the full form Factorial[n].  For example:

Factorial[4]

FullForm of n!

To see that that n! is in general (for nonegative integers n) an abbreviation for Factorial[n].   Look at the FullForm of n!—the way Mathematica internally represents the abbreviated expression n!:

FullForm[n!]

FullForm[4!]

Now why wasn't the result here Factorial[4] rather than the computed value 24?  The reason is that Mathematica evaluated 4! first, before even evaluting FullForm of its argument 4!.  You can force Mathematica to delay evaluating the argument 4! by Holding the entire expression:

Hold[FullForm[4!]]

Exercise.  (a) What is the FullForm of:   2 * n?   of 2 n?   of 2n?  of n2? 
(b) What is the FullForm of 2 + n?   of 2 + 5?  
(c) What is the FullForm of {1, 2, 3}?
(d) What is the FullForm of 1/2?

FullForm[{1, 2, 3}]

FullForm[1/2]

Factorial is Listable

One of the attributes of Factorial is that it is listable …

Attributes[Factorial]

… which means that Factorial automatically distributes over an entire list of nonnegative integers given as its argument.  For example:

Factorial[{1,  3,  5}]

To get a list of consecutive integers for the argument to Factorial, use the Range function:

Range[5]
Factorial[ Range[5] ]

Range[2,  6]
Factorial[ Range[2,  6] ]

Exercise.  (a) Produce as simply as you can the list of factorials of the integers from 0 to 50, inclusive.  
(c) Produce as simply as you can the list of factorials of just the even integers from 0 to 50, inclusive.
(c) Produce as simply as you can the list of factorials of just the odd integers from 1 to 49, inclusive.

To see both the inputs to Factorial and the corresponding values together, use the Table function:

Table[{n,  n!},  {n,  1,  5}]

To display the result of such a Table expression as an actual two-dimensional table, use the TableForm function in the same sort of way that you earlier used MatrixForm to display a matrix:

Table[{n,  n!},  {n,  1,  5}]  //  TableForm

Exercise.  After consulting the Help Browser, cause TableForm to put column headers "n" and "n!" over the table (and no row headings}.

Exercise.  Make a single table that shows the positive integers n from 1 to 25 along with the corresponding values of (2^n) and n!.

Exercise.  For what values of n is (2^n)  =  n!?  For what values of n is (2^n)  <  n!?  For what values of n is (2^n)  >  n!?

Factorial[n] as a product

The definition n!  =  
        1 · 2 · 
          3 · ⋯ · n (and 0!  =  1) can be expressed by using the product notation ∏ that is analogous to the familiar sum notation ∑:
	n!  =  (∏_(i = 1)^n) i and 0!  =  1.
This product notation can be entered directly into Mathematica.  The easiest way is by using the BasicInput palette.  This palette of buttons, including (□^□) and (Sqrt[□]) and \(TraditionalForm\`∫□ d□\) should ordinarily be displayed as a separate, untitled window (unless the previous user of Mathematica closed that palette).  If you do not see the palette, from the File menu select Palettes and select BasicInput.  In the following input cell, first select the symbol i.  Then press the \(TraditionalForm\`∏\+\(□ = □\)\%□ □\) button on the BasicInput palette.  The cursor should now be in the leftmost little empty box—a placeholder—below the big product symbol ∏.  Now do this:
• Type i in that box.
• Press ⇥ to move the cursor to the rightmost placeholder box under the ∏ and type 1 there.
• Press  ⇥ to move the cursor to the placeholder box above the ∏ and type, say, 20 there.
You should now have the expression (∏_(i = 1)^20) i in the cell.  Then evaluate the cell.

i                                                    (* 
      you will convert the i to expression (∏_(i = 1)^20) 
          i *)

Exercise.  What is the FullForm of (∏_(i = 1)^n) i?  (Your answer should not be Factorial[n], since that is the FullForm of the product after it is evaluated.  You want the FullForm before the product is evaluated!)

Exercise.  (a) By analogy with the preceding, what expression involving ∑ would yield the sum of the integers 1,  2,  3,  …,  n (for arbitrary, non-specific, n)? 
(b) What is the FullForm of this expression?
(c) Evaluate your expression for an arbitrary, non-specific,  value of n..
(d) Evaluate your expression for a specific value of n, say n  =  10.

⚠  Before proceeding, Clear any value you might have assigned to the symbol n.

Factorials by functional programming

Here's a shortcut that accomplishes the same thing as using the product notation (∏_(i = 1)^n) i:

Apply[Times,  Range[20]]

Check that the result of the preceding input cell is really the value of 20!.

Why does Apply[Times, Range[n]] have the same value as Factorial[n]?  First, what about Times?

Times[2,  7]

FullForm[m n]

Hold[FullForm[2 7]]

Look it up in the HelpBrowser.
Exercise.  Look up Times in the HelpBrowser in order to help understand why the expression in the preceding input cell does actually return result 20!.

Predict the values of each of the following:

Times[{1,  2,  3}]

Times[{1,  2}]

Times[{1}]

Times[{}]

Range[0]

Exercise.  (a) Now look up Apply in the HelpBrowser in order to understand why the expression Apply[Times, Range[20]] has the same value as 20!.
(b) Explain why Apply[Times, {1}] has the same value as 1!.
(c) Explain why Apply[Times, {}] has the same value as 0!.

Let's just make sure again that the names i and n do not have any specific numeric value:

Clear[i,  n]

Recall that an expression of the form fname[var_] := expr is used to define a function of your own.  For example, the expression
	f[x_] := (x^2)
would define the squaring function.
Now define anew the factorial function in terms of Times and Apply:

factorialProduct[n_]  :=  Apply[Times,  Range[n]]

factorialProduct[4]

Exercise.  Evaluate factorialProduct at 2, at 1, and at 0.  Are the values what you expected?  Can you explain why they are correct?

Factorials by recursive programming

In this section you learn how to define the factorial function by a recursive method—a method that, in effect, expresses the factorial function in terms of itself!

The recursive definition is based upon the following two mathematical properties the factorial function:
	\(\(0!\)\  = \ 
            1, \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ \ 
\(n!\)\  = \ 
            n\ \(\((n\  - \ 1)\)!\)\ \ \ \ \ \ \ \ \ \ \(\((n\  = \ 1, \ 2, \ 
                  3, \ …)\)\(.\)\)\)
Apply these properties repeatedly:
	\(TraditionalForm\`\(1!\)\  = \ \(1\ \((\(0!\))\)\  = \ \(1\ ·\ 1\  = \ 
            1\)\), 
\(2!\)\  = \ \(2\ \((\(1!\))\)\  = \ \(2\ ·\ 1\  = \ 
            2\)\), 
\(3!\)\  = \ \(3\ \((\(2!\))\)\  = \ \(3\ ·\ 2\  = \ 6\)\), 
\(\(etc\)\(.\)\)\)

The following definition in Mathematica is a direct realization of the two properties:

fact[0]  =  1;
fact[n_]  :=  n fact[n  -  1]

It really works—try the following:

fact[5]
fact[20]

Another way to view what we have called a recursive definition is a definition by pattern matching.  (Indeed, pattern matching is one of the underlying principles to all of Mathematica.)  The expression fact[0] = 1 is a particular rule.  The expression fact[n_] := n fact[n - 1] is a general rule.  So how does Mathematica know which rule to use for a given input?  It always tries particular rules first.  Only when no particular rule applies does it try a general rule.

The reason that those two simple rules suffice for Mathematica to calculate fact[n] no matter how large n is—provided you wait long enough and your computer doesn't run out of memory—is that when Mathematica evaluates an expression, it repeatedly applies all relevant rule until it cannot further change the value.

You can follow all the steps in applying the two rules by using the Trace function.  For example:

Trace[fact[5]]

Now try the following and note at each step what is known about the function fact.

?fact

fact[1]

?fact

fact[2]

?fact

fact[3]

No matter what argument you supply to fact, Mathematica knows nothing more about fact than when it started.  And that is why this recursive definition, as given above, is terribly inefficient.  Indeed, so inefficient that Mathematica won't even complete the calculation if you uselarge arguments:

fact[100]

fact[250]

fact[300]

To force Mathematica to remember values of fact that it has already calculated—and thereby make the use of fact much more efficient—we use the technique of caching.  Here is the revised definition of fact, using cahcing:

fact[0]  =  1;
fact[n_]  :=  fact[n]  =  n fact[n  -  1]

Try the following and note at each step what is known about the revised function fact.

?fact

fact[1]

?fact

fact[2]

?fact

fact[3]

fact[100]

?fact

fact[300]

fact[500]

Of course there is a downside to caching: a tradeoff of space efficiency for time efficiency.  With the original, uncached version of fact, a large number of calls to fact were requred, one after the other, to find fact[n] for a large n; with the revised, cached version of fact, a large number of particular rules about fact—the values of fact at all numbers up to and including the largest used so far as argument—have to be stored in memory.

Factorials by procedural programming

In this section you learn how to define the factorial function by an iterative procedure.  Of course such an iterative procedure superfluous—because Factorial is a built-in function, and because even if Factorial were not built-in, you could define it by means of the simple "functional programming" expression involving Apply and Times.

First use a Do expression.

?Do

Notice that Do is a full-fledged function, not just a "keyword".  However, it returns a Null result that ordinarily does not display in an Out cell.  To see the final result found by a Do expression, you have to print or otherwise evaluate that final value.

Do[(i^2),  {i,  5}]

Do[Print[(i^2)],  {i,  5}]

Do[temp  =  (i^2),  {i,  5}];  temp

OK, now use Do to calculate Factorial[5]:

runningProduct  =  1;
Do[
        runningProduct  =  i × runningProduct,  {i,  
          5}];
runningProduct

The three separate expressions in the preceding input cell, separated by semicolons, constitute a compound expression.  To make such a compound expression calculate Factorial[n] for arbitrary n and not just the particular n = 5, change 5 to n and encapsulate the entire compound expression in a Module, like this:

factorialDo[n_]  :=  
      Module[{runningProduct,  i},     (* 
          declare local variables *) 
runningProduct  =  1; 
Do[
          runningProduct  =  i × runningProduct,  {i,  
            n}]; 
runningProduct
]

Exercise.  Test the preceding function on various inputs, including 5, 1, and 0.  Explain why it works even when the input is 1 and when the input is 0.

If it makes you more comfortable, you can include in the procedural program factorialDo a separate test for input 0—even though, as you should have observed, such a test is really superfluous.  Here is the (unnecessarily) modified form:

overlyCautiousFactorialDo[
        n_]  :=   Module[{fact,  i},      
If[
          n  ==  0,     
          1,  
  
fact  =  
            1; 
Do[
            fact  =  i × fact,  {i,  
              n}]; 
fact]
]

overlyCautiousFactorialDo[0]

Exercise.  Write a procedural program to define Factorial that use a For expression; write another that uses a While expression.  Begin by looking up For and While in the HelpBrowser.

Type checking of input

Look again at the functional programming version, factorialProduct:

?factorialProduct

What happens if you supply a type of input that was not intended to be used?

factorialProduct[-2]

factorialProduct[1/2]

factorialProduct[0.5]

At the very least, you want the input to be an integer, that is, to be an expression which, when evaluated, has the Head of Integer.

{Head[5],  Head[2  +  3],  Head[-2],  Head[1/2],  
      Head[0.5]}

The way to prevent a function such as factorialProduct from evaluating its input unless that input has a Head of Integer is to append the desired Head to the input variable pattern n_, like this:

factorialProduct[n_Integer]  :=  Apply[Times,  Range[n]]

factorialProduct[5]
factorialProduct[1/2]
factorialProduct[0.5]

Oops!  Why didn't the new definition do what was expected now?  Look at the current rules associated with the name factorialProduct:

?factorialProduct

There are two rules: the more particular one just set down that restricts the input to be an integer; and the more general one from earlier that accepts any input!  Mathematica always tries the more particular rules first and, if it cannot apply those to evaluate the given expression, then it tries the more general rules.  So either you have to use a new name for the intended new function or else Clear the rules associated with the old function.  Let's do the latter:

Clear[factorialProduct]
factorialProduct[n_Integer]  :=  Apply[Times,  Range[n]]

factorialProduct[5]
factorialProduct[1/2]
factorialProduct[0.5]

So far OK.  But what about the following now?

factorialProduct[-2]

Head[-2]

A new way is needed to prevent factorialProduct[n] from being evaluated when n is not a nonnegative integer.  Just examining the Head of the input will not suffice:

Head[-2]

Instead, a new kind of restriction is needed on the input pattern n_.  One way to impose such a restriction is to append an expression of the form ?test, where test is a function that takes only the values True and False.  Here we use the function NonNegative:

Clear[factorialProduct]
factorialProduct[n_?NonNegative]  :=  
      Apply[Times,  Range[n]]

factorialProduct[5]
factorialProduct[-2]

factorialProduct[1/2]

Oops!  The lastest version of factorialProduct checks whether the input is nonnegative, but it no longer checks whether it is integer.  We need to check for both types:

Clear[factorialProduct]
factorialProduct[n_Integer?NonNegative]  :=  
      Apply[Times,  Range[n]]

factorialProduct[5]
factorialProduct[-2]
factorialProduct[1/2]

Exercise.  Include similar type checking via input patterns in the procedural version factorialDo of Factorial.

Exercise.  Include similar type checking via input patterns in the recursive version fact of Factorial.