math::calculus - Integration and ordinary differential equations
SYNOPSIS
package require Tcl 8
package require math::calculus 0.5
|
This package implements several simple mathematical algorithms:
The package is fully implemented in Tcl. No particular attention has been paid to the accuracy of the calculations. Instead, well-known algorithms have been used in a straightforward manner.
Ordinarily, such an equation would be written as:
d2y dy
a(x)--- + b(x)-- + c(x) y = D(x)
dx2 dx
|
A(x) = a(x)
B(x) = b(x) - a'(x)
C(x) = c(x) - B'(x) = c(x) - b'(x) + a''(x)
|
func(x) = 0
|
Notes:
Several of the above procedures take the names of procedures as arguments. To avoid problems with the visibility of these procedures, the fully-qualified name of these procedures is determined inside the calculus routines. For the user this has only one consequence: the named procedure must be visible in the calling procedure. For instance:
namespace eval ::mySpace {
namespace export calcfunc
proc calcfunc { x } { return $x }
}
#
# Use a fully-qualified name
#
namespace eval ::myCalc {
proc detIntegral { begin end } {
return [integral $begin $end 100 ::mySpace::calcfunc]
}
}
#
# Import the name
#
namespace eval ::myCalc {
namespace import ::mySpace::calcfunc
proc detIntegral { begin end } {
return [integral $begin $end 100 calcfunc]
}
}
|
Enhancements for the second-order boundary value problem:
Integrate x over the interval [0,100] (20 steps):
proc linear_func { x } { return $x }
puts "Integral: [::math::calculus::integral 0 100 20 linear_func]"
|
puts "Integral: [::math::calculus::integralExpr 0 100 20 {$x}]"
|
The differential equation for a dampened oscillator:
x'' + rx' + wx = 0 |
can be split into a system of first-order equations:
x' = y y' = -ry - wx |
Then this system can be solved with code like this:
proc dampened_oscillator { t xvec } {
set x [lindex $xvec 0]
set x1 [lindex $xvec 1]
return [list $x1 [expr {-$x1-$x}]]
}
set xvec { 1.0 0.0 }
set t 0.0
set tstep 0.1
for { set i 0 } { $i < 20 } { incr i } {
set result [::math::calculus::eulerStep $t $tstep $xvec dampened_oscillator]
puts "Result ($t): $result"
set t [expr {$t+$tstep}]
set xvec $result
}
|
Suppose we have the boundary value problem:
Dy'' + ky = 0
x = 0: y = 1
x = L: y = 0
|
This boundary value problem could originate from the diffusion of a decaying substance.
It can be solved with the following fragment:
proc coeffs { x } { return [list $::Diff 0.0 $::decay] }
proc force { x } { return 0.0 }
set Diff 1.0e-2
set decay 0.0001
set length 100.0
set y [::math::calculus::boundaryValueSecondOrder \
coeffs force {0.0 1.0} [list $length 0.0] 100]
|