I can be contacted at mailto:[email protected]

In view of the recent discussion on numerical analysis, (and having been doing such things for a quarter century), I do think it would be useful to have a slightly more friendly syntax. Here is one possible approach, within tcl: a procedure for writing c-like codeproccI don't use tcl for serious number crunching, but I still do a lot of numerics! Something like this would make things a bit easier.

**Advocacy**One of the things I've got from using Tcl is being forced to use other programming paradigms. In particular, the lack of numerically addressed arrays has proved a boon, and pushed me into using regular expressions, lists, associative arrays and other data structures that I was aware of, but tended not to use out of general conservatism. In particular, regexp and regsub provide powerful string manipulation facilities that used to take enormous amounts of work when implemented from scratch.I grew up programming in Fortran, and the only data structure available was the array! So much of my algorithmic thinking has been in terms of manipulating things as a numerically addressed array. If I wanted to build up a list of values, you'd have an array to hold the values, an integer to count the total, and other integers to point into the array when manipulating the data. After moving to c, I tended to keep the same programming style, still using integers to point into an array rather than addressing the thing directly with a pointer.Tcl is FUN. I enjoy tinkering with algorithms and Tcl has sufficient richness as a language to do things in either the traditional way or by some clever lateral thinking.

proc average {theList} { set listLength [llength $theList] set total 0 for {set i 0} {$i < [llength $theList]} {incr i} { set x [lindex $theList $i] set total [expr {$total+$x}] } return [expr $total/$listLength.] } proc tclAverage {theList} { set total [expr [join $theList +]] return [expr $total/[llength $theList].] }

**Tcl Performance**The following program is a good test for pure numerical performance, it calculates the roots for a cubic equation, and is floating point intensive, with transcendentals (trig and inverse trig), powers, square roots and a fair bit of branching.

proc cubic { a1 a2 a3 } { set PI 3.1415926535897932385 set a12 [expr {$a1*$a1}] set Q [expr {($a12-3*$a2)/9.}] set R [expr {(2*$a1*$a12-9*$a1*$a2+27*$a3)/54.}] set Z [expr {$Q*$Q*$Q - $R*$R}] if {$Z == 0.0} { if {$Q==0.0} { set x0 [expr {-$a1/3.}] return [list $x0] } else { if {$R<0} { set x0 [expr {2*sqrt($Q)-$a1/3.}] set x1 [expr {-sqrt($Q)-$a1/3.}] } else { set x0 [expr {sqrt($Q)- $a1/3.}] set x1 [expr {-2*sqrt($Q)-$a1/3.}] } return [list $x0 $x1] } } elseif {$Z<0.0} { set z2 [expr {pow(sqrt(-$Z) + abs($R), 1.0/3.0)}] set z1 [expr {$z2 + $Q/$z2}] if {$R>0} {set z1 [expr {-$z1}]} set x0 [expr {$z1 - $a1/3.}] return [list $x0] } else { set theta [expr {$R/sqrt($Q*$Q*$Q)}] set Q2 [expr {-2.0*sqrt($Q)}] set theta [expr {acos($theta)}] set x0 [expr {$Q2*cos($theta/3.0)-$a1/3.}] set x1 [expr {$Q2*cos(($theta+2*$PI)/3.0)-$a1/3.}] set x2 [expr {$Q2*cos(($theta+4*$PI)/3.0)-$a1/3.}] return [list $x0 $x1 $x2] } }# # And we can test it with the following... #

proc randomCubic {} { set a [expr {rand()*10.0-20.}] set b [expr {rand()*10.0-20.}] set c [expr {rand()*10.0-20.}] cubic $a $b $c }Anyway the proc "cubic" has about forty floating point operations, including several transcendental evaluations. On my 550MHz Celeron based laptop, running

time {randomCubic} 10000gives me 48 microseconds per iteration, so about 20000 iterations per second. This translates into 800,000 floating point operations per second, or 0.8 megaflop.To put this into perspective, the workhorse VAX on which I did much of my Fortran programming 20 odd years ago could knock out 0.1 megaflop if you were the only user on it. It was perfectly adequate for running large computationally intensive finite element programs - at least in its day.Remember, not much more than a generation ago, nearly all engineering design was done using pencil, paper and slide rules. Tcl has perfecly adequate speed for most of the problems the average engineer encounters. You wouldn't want to do real time graphics or model protein folding but for basic engineering modelling, it's ideal.(The cubic procedure above is used in process modelling software to solve for gas compressibility using the Soave-Redlich Kwong and Peng-Robertson models)