Arjen Markus (22 november 2002) It must be the time of year, but I felt very mathematical this week. So I sat down and wrote a little script that does a mathematical job: it allows you to define polynomial functions and (truncated) power series. Obviously, this is just a quick version and a full-blown package should include such things as:
Nevertheless, if you have a need for other functions than the ones offered by expr, this may give you some idea. The proof of the pudding is usually in the eating - so the example includes the zeroth order Bessel function of the first kind, J0. It is my favourite :)
# polynoms.tcl -- # # Package for handling polynoms and power series # (sample Workbench module) # # Notes: # This package is a quick hack to get started only # # Version information: # version 0.1: initial implementation, november 2002
package provide Polynoms 0.1
namespace eval ::polynoms {
namespace export polynom powerseries
# powerseries -- # Construct a procedure that implements a power series # and return its name # # Arguments: # name Name of the new procedure # nterms Number of terms # calcterm Procedure to calculate a term # # Result: # Name of procedure that will calculate the power series # proc powerseries {name nterms calcterm} {
set terms {}
for { set n 0 } { $n <= $nterms } { incr n } {
lappend terms [$calcterm $n]
}
polynom $name $terms}
# polynom -- # Construct a procedure that implements a polynom # and return its name # # Arguments: # name Name of the new procedure # terms Coefficients (from degree 0 to degree n) # # Result: # Name of procedure that will calculate the polygon # proc polynom {name terms} {
set first [lindex $terms end]
set rterms {}
for { set n [expr {[llength $terms]-2}] } { $n >= 0 } { incr n -1 } {
lappend rterms [lindex $terms $n]
}
interp alias {} $name {} [namespace current]::EvalPolynom \
$first $rterms
return $name}
# EvalPolynom -- # Calculate the value of a polynom # # Arguments: # rterms Coefficients from degree n to degree 0 # x Value of the argument # # Result: # Value of the polynom at x # proc EvalPolynom {first rterms x} {
set value $first
foreach term $rterms {
set value [expr {$value*$x+$term}]
}
return $value}
} ;# End of namespace
# # Run the program # namespace import ::polynoms::*
puts "Simple parabola" polynom parabola {0.0 0.0 1.0} for { set i 0 } { $i < 10 } { incr i } {
puts [parabola $i]
}
puts "Exponential function"
proc exp_terms {n} {
set fac 1
for {set i 1} {$i <= $n} {incr i} {
set fac [expr {$fac*$i}]
}
expr {1.0/$fac}} powerseries exp 10 exp_terms for { set i 0 } { $i < 10 } { incr i } {
puts "[exp $i] [expr {exp($i)}]"}
puts "Bessel function J0"
proc J0_terms {n} {
if { $n%2 == 1 } {
return 0.0
} else {
set fac 1
set sign [expr {$n%4==2?-1:1}]
for {set i 1} {$i <= $n/2} {incr i} {
set fac [expr {$fac*$i}]
set sign [expr {$sign/4.0}]
}
expr {$sign/($fac*$fac)}
}} powerseries J0 10 J0_terms for { set x 0.0 } { $x < 2.0 } { set x expr {$x+0.2} } {
puts "$x [J0 $x]"
}