I’m going to beat around the bush for a while, but there is a point.
Semigroups
The answer is, the associative property of the binary reduction operation.
That’s pretty abstract, but multiplication is a good example. If x, y and z are some natural numbers (or integers, or rational numbers, or real numbers, or complex numbers, or N×N matrices, or any of a whole bunch more things), then x×y is the same kind of number as both x and y. We started with two numbers, so it’s a binary operation, and got one, so we reduced the count of numbers we had by one, making this a reduction operation. And (x×y)×z is always the same as x×(y×z), which is the associative property.
(If you already know all this, you can skip to the next section.)
A few more things you often see in computer science that work the same way:
- adding any of those kinds of numbers instead of multiplying
- concatenating strings (
"a"+"b"+"c"
is "abc"
whether you start with "ab"+"c"
or "a"+"bc"
)
- Splicing two lists together.
[a]++[b]++[c]
is similarly [a,b,c]
either from back to front or front to back.
cons
on a head and tail, if you think of the head as a singleton list. That’s just concatenating two lists.
- taking the union or the intersection of sets
- Boolean and, Boolean or
- bitwise
&
, |
and ^
- composition of functions: (f∘g)∘h x = f∘(g∘h) x = f(g(h(x)))
- maximum and minimum
- addition modulo p
Some things that don’t:
- subtraction, because 1-(1-2) ≠ (1-1)-2
- x⊕y = tan(x+y), because tan(π/4 + π/4) is undefined
- multiplication over the negative numbers, because -1 × -1 is not a negative number
- division of integers, which has all three problems!
- logical not, because it has only one operand, not two
int print2(int x, int y) { return printf( "%d %d\n", x, y ); }
, as print2( print2(x,y), z );
and print2( x, print2(y,z) );
have different output.
It’s a useful enough concept that we named it. A set with an operation that has these properties is a semigroup. So, the real numbers under multiplication are a semigroup. And your question turns out to be one of the ways this kind of abstraction becomes useful in the real world. Semigroup operations can all be optimized the way you’re asking about.
Try This At Home
As far as I know, this technique was first described in 1974, in Daniel Friedman and David Wise’s paper, “Folding Stylized Recursions into Iterations”, although they assumed a few more properties than it turns out they needed to.
Haskell is a great language to illustrate this in, because it has the Semigroup
typeclass in its standard library. It calls the operation of a generic Semigroup
the operator <>
. Since lists and strings are instances of Semigroup
, their instances define <>
as the concatenation operator ++
, for example. And with the right import, [a] <> [b]
is an alias for [a] ++ [b]
, which is [a,b]
.
But, what about numbers? We just saw that numeric types are semigroups under either addition or multiplication! So which one gets to be <>
for a Double
? Well, either one! Haskell defines the types Product Double
, where (<>) = (*)
(that is the actual definition in Haskell), and also Sum Double
, where (<>) = (+)
.
One wrinkle is that you used the fact that 1 is the multiplicative identity. A semigroup with an identity is called a monoid and is defined in the Haskell package Data.Monoid
, which calls the generic identity element of a typeclass mempty
. Sum
, Product
and list each has an identity element (0, 1 and []
, respectively), so they are instances of Monoid
as well as Semigroup
. (Not to be confused with a monad, so just forget I even brought those up.)
That’s enough information to translate your algorithm into a Haskell function using monoids:
module StylizedRec (pow) where
import Data.Monoid as DM
pow :: Monoid a => a -> Word -> a
{- Applies the monoidal operation of the type of x, whatever that is, by
- itself n times. This is already in Haskell as Data.Monoid.mtimes, but
- let’s write it out as an example.
-}
pow _ 0 = mempty -- Special case: Return the nullary product.
pow x 1 = x -- The base case.
pow x n = x <> (pow x (n-1)) -- The recursive case.
Importantly, note that this is tail recursion modulo semigroup: every case is either a value, a tail-recursive call, or the semigroup product of both. Also, this example happened to use mempty
for one of the cases, but if we hadn’t needed that, we could have done it with the more general typeclass Semigroup
.
Let’s load this program up in GHCI and see how it works:
*StylizedRec> getProduct $ pow 2 4
16
*StylizedRec> getProduct $ pow 7 2
49
Remember how we declared pow
for a generic Monoid
, whose type we called a
? We gave GHCI enough information to deduce that the type a
here is Product Integer
, which is an instance
of Monoid
whose <>
operation is integer multiplication. So pow 2 4
expands recursively to 2<>2<>2<>2
, which is 2*2*2*2
or 16
. So far, so good.
But our function uses only generic monoid operations. Previously, I said that there is another instance of Monoid
called Sum
, whose <>
operation is +
. Can we try that?
*StylizedRec> getSum $ pow 2 4
8
*StylizedRec> getSum $ pow 7 2
14
The same expansion now gives us 2+2+2+2
instead of 2*2*2*2
. Multiplication is to addition as exponentiation is to multiplication!
But I gave one other example of a Haskell monoid: lists, whose operation is concatenation.
*StylizedRec> pow [2] 4
[2,2,2,2]
*StylizedRec> pow [7] 2
[7,7]
Writing [2]
tells the compiler that this is a list, <>
on lists is ++
, so [2]++[2]++[2]++[2]
is [2,2,2,2]
.
Finally, an Algorithm (Two, in Fact)
By simply replacing x
with [x]
, you convert the generic algorithm that uses recursion modulo a semigroup into one that creates a list. Which list? The list of elements the algorithm applies <>
to. Because we used only semigroup operations that lists have too, the resulting list will be isomorphic to the original computation. And since the original operation was associative, we can equally well evaluate the elements from back to front or from front to back.
If your algorithm ever reaches a base case and terminates, the list will be non-empty. Since the terminal case returned something, that will be the final element of the list, so it will have at least one element.
How do you apply a binary reduction operation to every element of a list in order? That’s right, a fold. So you can substitute [x]
for x
, get a list of elements to reduce by <>
, and then either right-fold or left-fold the list:
*StylizedRec> getProduct $ foldr1 (<>) $ pow [Product 2] 4
16
*StylizedRec> import Data.List
*StylizedRec Data.List> getProduct $ foldl1' (<>) $ pow [Product 2] 4
16
The version with foldr1
actually exists in the standard library, as sconcat
for Semigroup
and mconcat
for Monoid
. It does a lazy right fold on the list. That is, it expands [Product 2,Product 2,Product 2,Product 2]
to 2<>(2<>(2<>(2)))
.
This is not efficient in this case because you can’t do anything with the individual terms until you generate all of them. (At one point I had a discussion here about when to use right folds and when to use strict left folds, but it went too far afield.)
The version with foldl1'
is a strictly-evaluated left fold. That is to say, a tail-recursive function with a strict accumulator. This evaluates to (((2)<>2)<>2)<>2
, calculated immediately and not later when it’s needed. (At least, there are no delays within the fold itself: the list being folded is generated here by another function that might contain lazy evaluation.) So, the fold calculates (4<>2)<>2
, then immediately calculates8<>2
, then 16
. This is why we needed the operation to be associative: we just changed the grouping of the parentheses!
The strict left fold is the equivalent of what GCC is doing. The leftmost number in the previous example is the accumulator, in this case a running product. At each step, it’s multiplied by the next number in the list. Another way to express it that is: you iterate over the values to be multiplied, keeping the running product in an accumulator, and on each iteration, you multiply the accumulator by the next value. That is, it’s a while
loop in disguise.
It can sometimes be made just as efficient. The compiler might be able to optimize away the list data structure in memory. In theory, it has enough information at compile time to figure out it should do so here: [x]
is a singleton, so [x]<>xs
is the same as cons x xs
. Each iteration of the function might be able to re-use the same stack frame and update the parameters in place.
Either a right fold or a strict left fold could be more appropriate, in a particular case, so know which one you want. There are also some things only a right fold can do (such as generate interactive output without waiting for all the input, and operate on an infinite list). Here, though, we’re reducing a sequence of operations to a simple value, so a strict left fold is what we want.
So, as you can see, it is possible to automatically optimize tail-recursion modulo any semigroup (one example of which is any of the usual numeric types under multiplication) to either a lazy right fold or a strict left fold, in one line of Haskell.
Generalizing Further
The two arguments of the binary operation don’t have to be the same type, so long as the initial value is the same type as your result. (You can of course always flip the arguments to match the order of the kind of fold you’re doing, left or right.) So you might repeatedly add patches to a file to get an updated file, or starting with an initial value of 1.0, divide by integers to accumulate a floating-point result. Or prepend elements to the empty list to get a list.
Another type of generalization is to apply the folds not to lists but to other Foldable
data structures. Often, an immutable linear linked list is not the data structure you want for a given algorithm. One issue I did not get into above is that it’s a lot more efficient to add elements to the front of a list than to the back, and when the operation is not commutative, applying x
on the left and the right of the operation aren’t the same. So you would need to use another structure, such as a pair of lists or binary tree, to represent an algorithm that could apply x
on the right of <>
as well as to the left.
Also note that the associative property allows you to regroup the operations in other useful ways, such as divide-and-conquer:
times :: Monoid a => a -> Word -> a
times _ 0 = mempty
times x 1 = x
times x n | even n = y <> y
| otherwise = x <> y <> y
where y = times x (n `quot` 2)
Or automatic parallelism, where each thread reduces a subrange to a value that is then combined with the others.
if(n==0) return 0;
(not return 1 like in your question).x^0 = 1
, so that's a bug. Not that it matters for the rest of the question, though; the iterative asm checks for that special case first. But strangely, the iterative implementation introduces a multiply of1 * x
that wasn't present in the source, even if we make afloat
version. gcc.godbolt.org/z/eqwine (and gcc only succeeds with-ffast-math
.) $\endgroup$return 0
has been fixed. The multiplication by 1 is interesting. I'm not sure what to make of it. $\endgroup$float
without-ffast-math
, even though it's the same value being multiplied every time. (Except for the 1.0f` which might be the sticking point?) $\endgroup$