# tensors.tcl --
#     Use the Tensor package to solve a simple diffusion problem
#
package require Tensor

# initialise --
#     Create the tensors
#
# Arguments:
#     value           Value at cell x=0, t=0
#     noCells         Number of cells
#
# Result:
#     None
#
# Side effects:
#     Creates the tensors solution, work and deriv
#
proc initialise {value noCells} {

    ::tensor::create solution -type double -size [expr {$noCells+2}]
    ::tensor::create deriv    -type double -size [expr {$noCells+2}]
    ::tensor::create work     -type double -size [expr {$noCells+2}]

    set zero [expr {$noCells/2}]
    solution = scalar 0.0
    solution section $zero = scalar $value

    work  = scalar 0.0
    deriv = scalar 0.0
}

# nextTime --
#     Compute the solution for the next time
#
# Arguments:
#     diff            Diffusion coefficient (m2/s)
#     dx              Size of a grid cell
#     dt              Time step (s)
#
# Result:
#     None
#
# Side effects:
#     An updated solution
#
proc nextTime {diff dx dt} {

    set maxCell   [expr {[solution dimensions] - 1}]
    set maxCellM1 [expr {$maxCell - 1}]
    set maxCellM2 [expr {$maxCell - 2}]

    ::tensor::expr "work(i=1:$maxCellM1) = -2.0 * solution(i=1:$maxCellM1)
                                                  + solution(i=0:$maxCellM2)
                                                  + solution(i=2:$maxCell)"
    deriv = scalar 0.0
    ::tensor::expr "deriv(i) = $diff * work(i) / $dx / $dx"
    ::tensor::expr "solution(i) = solution(i) + $dt * deriv(i)"
}

# main --
#     Run it
#
initialise 1.0 10

set dx   0.1
set dt   0.01
set diff 0.1

for { set i 0 } { $i < 10 } { incr i } {
   nextTime $diff $dx $dt
}

#solution writechannel stdout text ;# Did not print anything useful ...

solution put into a
parray a