bifurcation

Page started by Lars H, 9 July 2004.

For some time now I've been toying with an idea concerning Tcl and functional programming. Recursion is the fundamental concept in functional programming (in the purer installations it is even the only way to handle iteration), but recursively heavy algorithms are usually not a good idea in Tcl, because the cost of recursion is comparatively large, and certainly larger than in the average functional language. One interpretation of this is that Tcl is useless for functional programming, but since I don't believe in that my conclusion must instead be that straightforward recursion is in Tcl a too primitive control structure to be the single solution to all problems. There has to be other things that can do better; the only problem is figure out what they are!

Tail call optimization is one thing. Even languages with low cost for recursive calls tend to suck if they haven't got this built in; it is what pure functional languages use to achieve decent loop performance. But perhaps it is, for practical purposes, rather a distinct control structure that can be dressed out as a recursion? Consider the following command

proc iterate {vars vals body} {
    while {![
        catch {
            uplevel 1 [list foreach $vars $vals break]
            uplevel 1 [list if 1 then $body]
        } vals
    ]} {} ; # Empty while body
    return $vals
}

Syntax-wise it is like foreach, but it only uses the given list of values on the first iteration. On the second iteration, it uses the result from the first iteration, on the third iteration it uses the result from the second iteration, etc. However if one uses the return command in the body then the iterate command returns with that result (break would have been more natural, but that does not support returning a value).

To see how this can be equivalent to a tail recursion, consider the following tail call optimizable implementation of the factorial function.

proc factorial1 { n { result 1. } } {
    if {$n <= 1} then {return $result}
    factorial1 [expr {$n-1}] [expr { $n * $result }]
}

The corresponding implementation using iterate is

proc factorial2 {n} {
    iterate {n result} [list $n 1.] {
        if {$n <= 1} then {return $result}
        list [expr {$n-1}] [expr { $n * $result }]
    }
}

Is factorial2 more natural than factorial1? That is probably to equal parts a matter of custom and personal taste. For things usually described as functions (e.g. the factorial) the factorial1 style is natural, but for things more commonly described as sequences (or a family of sequences) the factorial2 style has a point.

In terms of speed, and despite what was said above concerning recursion overhead, factorial1 wins by roughly a factor 5. Furthermore if one removes the seemingly pointless "if 1 then" (which is there to force byte-compilation of $body) from the definition of iterate then it will be more like a factor 20! From that one can learn that the overhead of a script level implementation of this kind of control structure is generally much larger than the overhead of recursion, so speed comparisons are probably uninteresting unless one reimplements the control structure commands in C. I'll leave the speed issue at that (unless someone jumps in to show how to make a Critcl implementation of this kind of command). Besides, the fastest Tcl implementation of the factorial is anyway the straightforward iterative one, not any functional programming variant.

iterate is however only a warm-up, the main topic here is bifurcation. A problem with tail call optimization is that there must not be any post-processing of the tail call result, whereas most naive recursive definitions of a function does post-process the result before returning it. The naive definition of a recursive factorial is rather

proc factorial0 {n} {
    if {$n>1} then {
        return [expr {$n * [factorial0 [expr {$n-1}]]}]
    } else {
        return 1
    }
}

The idea of the bifurcation command is to provide something similar to recursion by mimicing the fork/join of parallel processes. The body of iterate is split in two parts: one that spawns "recursive" calls and one that joins the results from these spawned calls back up. The factorial implementation there is

proc factorial3 {x} {
    bifurcation [list $x] n {
        if {$n > 1} then {
            spawn result [expr {$n-1}]
        } else {
            expr 1
        }
    } result {
        expr {$n * $result}
    }
}

but as the name hints, it is really meant to support recursions that branch out to multiple recursive calls at each step.

Syntax-wise, bifurcation is like a foreach that has been outfitted for fork/join branching. The syntax is

bifurcation initial-arguments argument-names fork-body ?branch-results join-body ...?

The first step is simple. The initial-arguments are assigned to the variables listed in argument-names, and then the fork-body is evaluated. If the fork-body doesn't use the companion command spawn then the bifurcation is complete and its result is the result of the fork-body.

However if it uses spawn, whose syntax is

spawn branch-name ?arguments ...

then it will be as if the bifurcation command was called recursively, once for each spawn, and these times with the arguments specified there instead of the initial-arguments. The results from these "recursive" (as far as the Tcl interpreter is concerned, it all happens in the same local context) calls are collected and assigned to the respective variables that are the branch-names of the spawn calls. Then the bifurcation command looks for a branch-results argument that lists precisely those branches that were spawned (it is an error if there isn't one, but the order in that the branches are listed is irrelevant) and evaluates the corresponding join-body. The purpose of this join-body is to compute a combined result for the bifurcation invokation, so typically it makes a small calculation and the result of that is the result returned by the bifurcation command. However, a join-body may also call spawn, in which case results from that are computed and subsequently joined back.

As an example, consider the traditional (and brain-dead) recursion for Fibonacci numbers. It can be coded using bifurcation as follows

proc Fibonacci {n} {
    bifurcation [list $n] n {
        # Fork-body
        if {$n>1} then {
            spawn m1 [expr {$n-1}]
            spawn m2 [expr {$n-2}]
        } else {
            expr 1
        }
    } {m1 m2} {
        # Join-body
        expr {$m1 + $m2}
    }
}

One may however take advantage of the fact that it all the bodies are evaluated in the same local context to remember values that have already been computed

proc Fibonacci2 {n} {
    if {![string is integer $n] || $n<0} then {
        error "Invalid input; must be a natural number."
    }
    set F(0) 1
    set F(1) 1
    bifurcation [list [expr {$n+0}]] n {
        if {[info exists F($n)]} then {
            set F($n)
        } else {
            spawn m1 [expr {$n-1}]
            spawn m2 [expr {$n-2}]
        }
    } {m1 m2} {
        set F($n) [expr {$m1 + $m2}]
    }
}

The arguments and branch-result variables are reset to their values for a particular invokation before a body is evaluated, but all other variables are left alone by bifurcation.

Yet another Fibonacci procedure. This one also has a recursion for negative numbers.

proc Fibonacci3 {n} {
    bifurcation [list $n] n {
        if {$n>=2} then {
            spawn m1 [expr {$n-1}]
            spawn m2 [expr {$n-2}]
        } elseif {$n<0} then {
            spawn p1 [expr {$n+1}]
            spawn p2 [expr {$n+2}]
        } else {
            expr 1
        }
    } {m1 m2} {
        expr {$m1 + $m2}
    } {p1 p2} {
        expr {$p2 - $p1}
    }
}

Here is an implementation of these commands:

proc bifurcation {iniargs argnames forkbody args} {
    upvar 1 __bifurcation_spawning spawnA
    if {[llength $args] & 1} then {
        error {Syntax is: bifurcation iniargs argnames forkbody ?joinargs joinbody ...?}
    }
    foreach {vars body} $args {
        set bodyA([lsort $vars]) $body
    }
    set stack [list [list -1 {} $iniargs {}]]
    while {[llength $stack]} {
        array unset spawnA
        uplevel 1 [list ::foreach $argnames [lindex $stack end 2] ::break]
        if {[llength [lindex $stack end]] <= 4} then {
            # Initial call, hence use the forkbody
            set res [uplevel 1 [list ::if 1 then $forkbody]]
        } else {
            # Joining up: Set the joinargs and use the joinbody.
            foreach {var val} [lindex $stack end 4] {
                uplevel 1 [list set $var $val]
            }
            set res [uplevel 1 [list ::if 1 then [lindex $stack end 3]]]
        }
        if {[array size spawnA]} then {
            # The body spawned, so we need to push those spawns onto 
            # the stack and prepare the current stack frame to 
            # recieve result. 
            set case [lsort [array names spawnA]]
            if {![info exists bodyA($case)]} then {
                error "No joinbody for case: $case"
            }
            lset stack end [lreplace [lindex $stack end] 3 end\
              $bodyA($case) [list]]
            set sp [expr {[llength $stack] -1}]
            foreach branch [array names spawnA] {
                lappend stack [list $sp $branch $spawnA($branch) {}]
            }
        } elseif {[lindex $stack end 0] >= 0} then {
            # Return a value to a previous stack frame, and pop this one 
            # off.
            lset stack [lindex $stack end 0] 4 [
              linsert [lindex $stack [lindex $stack end 0] 4] end\
                [lindex $stack end 1] $res
            ]
            set stack [lrange $stack [set stack 0] end-1]
        } else {
            # The last stack frame has been joined up. Return.
            return $res
        }
    }
}

proc spawn {branch args} {
    upvar 1 __bifurcation_spawning spawnA
    set spawnA($branch) $args
}

[To do: Explain data structures. Give more examples.]

First nontrivial example: a procedure for computing the independence number (alpha) of a graph, i.e., the size of the largest independent (no two elements are adjacent) set of vertices. Unlike trivialities like the Fibonacci numbers, this problem is NP-hard in general and thus all known algorithms require exponential execution time. Call syntax is

independence_number vertex-list edge-list

Vertices can be arbitrary strings. Edges are pairs of vertices.

proc independence_number {V E} {
    foreach v $V {set d($v) 0}
    foreach e $E {
        lappend AdjL([lindex $e 0]) [lindex $e 1]
        lappend AdjL([lindex $e 1]) [lindex $e 0]
        incr d([lindex $e 0])
        incr d([lindex $e 1])
    }
    bifurcation [list [array get d]] D {
        set Delta 0
        set 2m 0
        foreach {v deg} $D {
            if {$deg>$Delta} then {
                set Delta $deg
                set vmax $v
            }
            incr 2m $deg
        }
        if {$Delta <= 1} then {
            expr {([llength $D] - ${2m}) / 2}
        } else {
            array unset d
            array set d $D
            unset d($vmax)
            foreach v $AdjL($vmax) {
                if {[info exists d($v)]} then {incr d($v) -1}
            }
            spawn alpha1 [array get d]
            foreach u $AdjL($vmax) {
                if {[info exists d($u)]} then {
                    unset d($u)
                    foreach v $AdjL($u) {
                        if {[info exists d($v)]} then {incr d($v) -1}
                    }
                }
            }
            spawn alpha2 [array get d]
        }
    } {alpha1 alpha2} {
        expr { $alpha2>=$alpha1 ? $alpha2+1 : $alpha1 }
    }
}

Example 1: A 5-cycle.

% independence_number {1 2 3 4 5} {{1 2} {2 3} {3 4} {4 5} {5 1}}
2

Example 2: A 6-cycle.

% independence_number {1 2 3 4 5 6} {{1 2} {2 3} {3 4} {4 5} {5 6} {6 1}}
3

The recursion is based on picking a vertex $vmax and considering two cases: (i) independent sets that don't contain $vmax, (ii) independent sets that do contain $vmax. A largest independent set of type (i) is a largest independent set in the subgraph one gets by deleting $vmax; the alpha1 spawn computes this. If one removes $vmax from a largest independent set of type (ii) then one gets a largest independent set in the subgraph one gets by deleting $vmax and all its neighbours; the alpha2 spawn computes this.

The base for the recursion consists of graphs with maximal degree (Delta) at most 1. In this case the independence number is the number of vertices minus the number of edges.

Second nontrivial example: Rosza Peter's [2 ] function (often misattributed to Ackermann [1 ]). This is the function P of two nonnegative integers which satisfies

    P(0,n) = n+1
  P(m+1,0) = P(m,1)
P(m+1,n+1) = P(m,P(m+1,n))

Does this look harmless? Well, it isn't. The diagonal f(n)=P(n,n) almost surely grows faster anything you can come up with. (Technically it is a recursive (computable) function which is not primitive recursive, so you cannot in terms of primitive recursive functions give a bound on how long time it will take to compute. All functions you have made use of in any practical computation is likely to be primitive recursive.)

What is possible (but still hard) to compute when only standard Tcl integer arithmetic is available is function values with m=3. The straightforward recursive implementation of the above equations is

proc Peter {m n} {
    if {$m>0  && $n>0} then {
        return [Peter [expr {$m-1}] [Peter $m [expr {$n-1}]]]
    } elseif {$m>0} then {
        return [Peter [expr {$m-1}] 1]
    } else {
        return [expr {$n+1}]
    }
}

This is however so heavily recursive that things like [Peter 3 8] almost surely cause stack overflows in Tcl. Defining that procedure makes it very easy to crash the interpreter.

The corresponding bifurcation implementation is far more robust, because its internal stack is much cheaper in terms of memory than the C stack normal Tcl recursion uses.

proc Peter2 {x y} {
    bifurcation [list $x $y] {m n} {
        if {$m>0 && $n>0} then {
            spawn a $m [expr {$n-1}]
        } elseif {$m>0} then {
            spawn b [expr {$m-1}] 1
        } else {
            expr {$n+1}
        }
    } a {
      spawn b [expr {$m-1}] $a
    } b {
      set b
    }
}

NEM How does this compare to the higher-order fold functions (like foldr and foldl in Haskell)? For instance, the factorial function can be defined like this:

 proc fac n { foldr * 1 [.. 1 $n] }

This is functional without resorting to recursion. Of course, we need to define the various helper functions I've used. First, foldr, which is the key:

proc foldr {fun default list} {
    set ret $default
    foreach item [reverse $list] {
        set ret [uplevel 1 $fun [list $ret $item]]
    } 
    return $ret
}

And [reverse] (I could have in-lined this, but the code is clearer without, IMO):

proc reverse {list} {
    set ret [list]
    for {set i 0} {$i < [llength $list]} {incr i} { lappend ret [lindex $list end-$i] }
    return $ret
}

What's still missing? Ah yes, [*] as a function, and the range generator function:

proc * {args} { expr [join [lappend args 1] *] }
proc .. {start end {step 1}} {
    set ret [list]
    for {set i $start} {$i <= $end} {incr i $step} { lappend ret $i }
    return $ret
}

You can similarly define factorial in a variety of different ways using higher order functions or recursion - see [3 ] for an amusing variety. So, what does this bifurcation stuff add?

(NEM realises he probably should have used foldl for shorter code (no reverse)... and being a bit of a lefty!)

MSW Not knowing Haskell, I suppose though that this builds up an intermediate list to fold over in a loop... the recursion on the other hand does just the computation, isn't building up yet another temp to throw away. Ah well...

NEM Yes - it generates a list and then iterates through it (well the Haskell version lazily generates elements of the list as needed). The (non-tail-optimised) recursive version effectively creates the same list on the stack anyway. The point I was making was that most uses of recursive functions are, or can be reduced to, operations on lists (or lists of lists/trees etc), and there are several well-known higher-order functions for dealing with these (map, filter, foldl, foldr...). BTW - the code I gave is Tcl, not Haskell, in case there was some confusion.

MSW Argh, yes. Was just baffled by the [..], but you implemented that, sorry. About recursion mimicking traversal of other data structures via its call-graph, er, that's the point of recursion ? Not having to explicitely handle the data structure... And one can add "...while maintaining constant space usage." if the implementation is tail-recursive/calling.

Now that's elegant. I suppose you can say that you want recursion rather than managing the data structure yourself as long as this management is trivial (like folding and mapping often is), and as long as you don't have the structure handy anyways. Conceptually, there's no difference, except for the explicit / implicit building/traversal of the wanted data structure. (tree/list vs. call graph).

Lars H: To answer NEM's question about how this compares to fold functions, I suppose iterate is near equivalent to those. The main differences are (i) that one doesn't have to name the function being iterated over, which is probably only a stylistic difference, and (ii) that the loop does not have a given length (foldr is to foreach as iterate is to while). bifurcation is another matter though, since it supports "recursions" that branch into several subcases. It's not "constant space usage" in memory as a whole, but the usage of call stack space is constant.

MSW I think you can happily hand an anonymous function (or "block") to fold*/map* in about any programming language having one... Another difference would be that fold usually is a binary operation, working on lists

{a b c d e f}

with e.g. calling '*' like

((((((start * a) * b) * c) * d) * e) * f) (= foldl)

, while with iterate you can construct any arity... ?

Lars H: No, arbitrary arity would rather require bifurcation. iterate, foldl, and foldr all follow a linear path of evaluation, but in order to deal with higher arities it is necessary to traverse a tree. Basically, iterate computes the result of applying a sequence of unary operations to a (typically composite) value. Already generic binary operations require branching, but foldl and foldr restrict themselves to left- or right-associated expressions; the expressions evaluated are always on one of the two forms

((((((start * a) * b) * c) * d) * e) * f)
(a * (b * (c * (d * (e * (f * start))))))

and these can be translated to a sequence of unary operations on start. (In the first case: Multiply by first a, then b, then c, etc. on the right. In the second case: Multiply by first f, then e, then d, etc. on the left.)

As for "happily hand an anonymous function": Yes, I'd expect nothing less. One of the purposes of these experiments is however to see how far one can get without resorting to these traditional functional programming constructions. There could be other things which are more natural in Tcl, and Tcl is one of the few languages that even lets one experiment with inventing new control structures.