[Arjen Markus] (28 december 2009) I have been contemplating this for several years now, but I finally made a 
start with a wrapper extension to get the functionality of LAPACK into Tcl.

I started off with two simple BLAS routines (BLAS is the more basic library that underlies LAPACK - '''b'''asic '''l'''inear 
'''a'''lgebra '''s'''ystem, if I am not mistaken) and the Fortran wrapper you can find in my ftcl library (see [http://ftcl.sf.net]).
Immediately I ran into design issues:

Take the routine ''dapxy'', which computes a*X + Y (a a scalar, X and Y vectors) and returns the result in Y. The interface is such
that X and Y could also be two-dimensional arrays (matrices) so that you could use this routine to update the rows and
columns of a matrix Y using a vector X ''o'' a matrix X.

Questions that arise:

   * Is that flexibility something we want in the Tcl command as well?
   * Do we want to update a vector or should the Tcl command simpy return the result?
   * Do we want to expose the full interface? That is the most flexible, but it is also distinctly non-Tcl

Right now, my answers are:

   * We merely want to pass two vectors X and Y, not columns or rows of a matrix
   * We do not want to update any arguments - so the Tcl command will return the resulting vector and it is up to the user to put it into the original Y argument, if he/she wants to.
   * The interface should be as simple as:

   set result [daxpy $a $x $y]

(The Fortran interface is:

   call dapxy( n, a, dx, incx, dy, incy )

Because of incx and incy, you can pass two-dimensional arrays such that only a part is updated.)

Another question is: should we wrap all routines (some 1500!) or just the more commonly used ones? (If so, which ones?)

----
[AM] (31 december 2009) A discussion in the Tcl chatroom brought back an old idea: use byte arrays (binary strings) to store the vector or matrix data for a Fortran or C routine directly. Here is a small Tcl program that illustrates the basic storage mechanism (for vectors only at the moment):

======
# vecmat.tcl --
#     Utilities to manipulate Tcl lists that hold
#     binary data that can directly act as Fortran or C arrays
#
#     TODO: error checking
#

# Vecmat --
#     Namespace to hold these utilities
#
namespace eval ::Vecmat {
    namespace export vector matrix
}

# vector --
#     Utility to manipulate vectors
#
# Arguments:
#     method     Method in question
#     args       Any arguments for the method
#
# Returns:
#     Depends on the method
#
proc ::Vecmat::vector {method args} {

    switch -- $method {
        "create" {
             if { [llength $args] != 1 } {
                 return -code error "Invalid number of arguments - should be: vector create number"
             }
             return [::Vecmat::Vcreate {*}$args]
        }
        "set" {
             if { [llength $args] != 3 } {
                 return -code error "Invalid number of arguments - should be: vector set name index number"
             }
             return [::Vecmat::Vset {*}$args]
        }
        "tolist" {
             if { [llength $args] != 1 } {
                 return -code error "Invalid number of arguments - should be: vector tolist vector-value"
             }
             return [::Vecmat::Vtolist {*}$args]
        }
        "fromlist" {
             if { [llength $args] != 1 } {
                 return -code error "Invalid number of arguments - should be: vector fromlist list"
             }
             return [::Vecmat::Vfromlist {*}$args]
        }
        default {
             return -code error "Unknown method: $method"
        }
    }
}

# Vcreate --
#     Create a vector of doubles of a given length
#
# Arguments:
#     number     Number of elements
#
# Returns:
#     Vector with values initialised to 0.0
#
proc ::Vecmat::Vcreate {number} {

    set data [binary format d* [lrepeat $number 0.0]]
    return [list VECTOR $number $data]
}

# Vset --
#     Set a single element in a vector
#
# Arguments:
#     name       Name of the vector
#     index      Index of the element
#     value      The new value
#
# Returns:
#     Vector with new value at position index
#
# Note:
#     This updates the vector variable whose name is $name
#
proc ::Vecmat::Vset {name index value} {
    upvar 2 $name name_

    set data [lindex $name_ 2]
    set data [binary format a*@[expr {8*$index}]d $data $value]
    lset name_ 2 $data
    return $name_
}

# Vtolist --
#     Convert a vector into a regular Tcl list
#
# Arguments:
#     vector     value of the vector
#
# Returns:
#     List containing the values of the vector
#
proc ::Vecmat::Vtolist {vector} {

    binary scan [lindex $vector 2] d* data
    return $data
}

# Vfromlist --
#     Create a vector from a regular Tcl list of doubles
#
# Arguments:
#     list       list of values to be converted
#
# Returns:
#     Vector containing the values of the list
#
proc ::Vecmat::Vfromlist {list} {

    set data [binary format d* $list]
    return [list VECTOR [llength $list] $data]
}

# main --
#     Test the vector and matrix utilities
#
set v [::Vecmat::vector create 10]

foreach index {0 4} {
    ::Vecmat::vector set v $index [expr {cos($index)}]
}

puts [::Vecmat::vector tolist $v]

set w [::Vecmat::vector fromlist {1.0 .2 3.0}]
puts [::Vecmat::vector tolist $w] 
======

<<categories>>Category Mathematics|Category Numerical Analysis