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] |