`MonadPlus`in their semantics when tackling the now-trivial pandigital problems!), but I ran into a snag when I was factoring large triangular numbers.

The problem was my naïveté, of course.

My original version of

`factors`

proceeded thusly ...factors ::Integrala ⇒ a → [a]

factors x = filter (λy . x `rem` y ≡ 0) [1..x]

... which is straightforward enough. One can't get any simpler than: "Keep ('

`filter`

') all the numbers from 1 to `x`

that divide `x`

with no `rem`

ainder"!Simple? Yes. Naïve? Yes. For small numbers, this factorization algorithm did just fine, even for the 384

^{th}triangle, 73920, this algorithm found the 112 factors in a blink. The problem appeared when doing a (triangular) search along the Natural number line for an unknown, large, triangle. The factorization algorithm is linear, and combining it with a (nearly) linear search makes the entire activity quadratic. Twenty-six

*hours*after I began my search ("find the first triangular number with more than 500 factors") the algorithm was still going with no solution.

There

*must*be a better way to go about this. And, of course, there are more than a few ways, some of which you can pay for, but none, on the top of the search tree that appeared simple enough to comprehend and to implement in a reasonable time.

Fortunately, there also exists an obvious sub-linear solution: if we are given in the linear progression

`[1..x]`

at the current index, `y`

, that `y`

is a factor of `x`

, then obviously ... `z ``where` z = x / y

... is also a factor. Not only that, but, since we are progressing linearly,

`z`

now becomes the upper bound in the factorization algorithm. The fix-point of the algorithm described gives us a new, much more efficient, factorization ...factors x = factors' 2 (x, [1,x])

wherefactors' y res@(top, ans) | y ≥ top = ans

|otherwise= factors'' (x `divMod` y) y res

factors'' (newtop, 0) y (oldtop, ans)

= factors' (succ y) (newtop, y:oldtop:ans)

factors'' _ y res

= factors' (succ y) res

... Much more efficient? Yes. Self-describing? Not so much. Also, there's now a subtle error introduced by assuming all factors are unique, for ...

`> factors 25 → [5,5,1,25]`

... duplicates are introduced when factoring squares. Testing uniqueness adds its own complexity:

factors x = factors' 2 (x, [1,x])

wherefactors' y res@(top, ans) | y ≥ top = ans

|otherwise= factors'' (x `divMod` y) y res

factors'' (newtop, 0) y (oldtop, ans)

= factors' (succ y) (newtop,addin y oldtop ans)

factors'' _ y res

= factors' (succ y) res

addin x y list|x &equiv y= x:list

|otherwise = x:y:list

Complexity on top of complexity! This goes against my grain (as well as Richard O'Keefe's as you see in the side-bar quote), but is it worth it?

`> sort (factors (triangle 12375)) → [1,2,3,`*...572 others...*,76576500]

Yes, as the above interaction took no noticeable time, whereas before no solution was available after a day's worth of computation. So this optimized factorization algorithm is "good enough" for the task at hand, for I did solve the problem using that algorithm.

But, returning to the fix-point, which more than generalizes recursive algorithms -- it can be used to think of algorithms,

*qua*algorithms, more generally. In this particular case, I have no desire to view every one of the 576 factors of this triangular number. No, I simply wish to ensure that there are more than 500 factors, no matter what those factors are. So, in the large, something like a

`countfactors`

(that simply returns the number of factors, instead the list of all the factors) is more desirable. How are we to go about writing such an algorithm? Most coders, I regret to observe, would punt to the copy-and-paste style of coding, as the desired changed is buried deeply within the algorithm. Not only that, but the type signature of the function varies with our alternatives: in what we've developed so far, the algorithm returns a list, but what we wish to have instead is (only one) number.Fortunately, Dirk Thierbach introduced me to this particular use of the fix-point -- we simply extract the engine of computation into an external function. So now we have two functions that collapse into one line each ...

showfactors ::Integrala ⇒ a → [a]

showfactors x = factors x (λ y z ans .ify ≡ ztheny:anselsey:z:ans) []

countfactors ::Integrala ⇒ a → a

countfactors x = factors x (λ y z ans .ify ≡ zthen1 + anselse2 + ans) 0

... with

`factors`

being (slightly) modified to accept the new functional argument and base case:factors ::Integrala ⇒ a → (a → a → b → b) → b → b

factors x adder ans = factors' 2 (x, adder 1 x ans)

wherefactors' y (tot, ans) | y ≥ tot = ans

|otherwise= factors'' (x `divMod` test) y (tot, ans)

factors'' (tot, 0) y (_, ans)

= factors' (y+1) (tot, adder y tot ans)

factors'' _ _ y ans

= factors' y ans

With this approach, the return type depends on the calling function's use and is threaded throughout the utility function,

`factors`

, by the `adder`

computation engine.One final note: the λ-terms in both calling functions

`showfactors`

and `countfactors`

are also of the same structure, so we can again perform surgery, extracting the engine from the structure:mkAdder ::Integrala ⇒ (a → b) → (b → b → b) → a → a → (b → b)

mkAdder f g y z = g (ify ≡ zthenf yelseg (f y) (f z))

Now, anyone who's been working with monads for more than an brief period will see the

`b`

type as monadic and the functional type `(a → b)`

as a lifting function (`return`) and the composition function,

`g`

, as monadic addition. And, as it turns out, this concept fits very nicely with the new implementation ...showfactors x = factors x (mkAdderreturnmplus) []

countfactors x = factors x (mkadder (const1) (+)) 0

With these new implementations, how much does

`factors`

need to change? Not one bit. Functional purity: why can't all programming languages have this? And with this generalization, I was able to achieve the desired result:`> countfactors (triangle 12375) → 576`

Now it is known that both lists of some type (in this case integers) and also integers under addition are monoids. Given that,

`mkAdder`

can be further simplified. Also, the `if-then-else`is easily replaced by the

`Either`type and proper use of arrows, but these are discussions for another day.

## No comments:

Post a Comment