convert math.blas.matrices to use fortran calls

db4
Joe Groff 2009-02-09 16:59:26 -06:00
parent 35c54a91ac
commit d160b80dac
2 changed files with 54 additions and 58 deletions

View File

@ -14,34 +14,34 @@ ARTICLE: "math.blas-types" "BLAS interface types"
"BLAS vectors come in single- and double-precision, real and complex flavors:" "BLAS vectors come in single- and double-precision, real and complex flavors:"
{ $subsection float-blas-vector } { $subsection float-blas-vector }
{ $subsection double-blas-vector } { $subsection double-blas-vector }
{ $subsection float-complex-blas-vector } { $subsection complex-float-blas-vector }
{ $subsection double-complex-blas-vector } { $subsection complex-double-blas-vector }
"These vector types all follow the " { $link sequence } " protocol. In addition, there are corresponding types for matrix data:" "These vector types all follow the " { $link sequence } " protocol. In addition, there are corresponding types for matrix data:"
{ $subsection float-blas-matrix } { $subsection float-blas-matrix }
{ $subsection double-blas-matrix } { $subsection double-blas-matrix }
{ $subsection float-complex-blas-matrix } { $subsection complex-float-blas-matrix }
{ $subsection double-complex-blas-matrix } { $subsection complex-double-blas-matrix }
"There are BOA constructors for all vector and matrix types, which provide the most flexibility in specifying memory layout:" "There are BOA constructors for all vector and matrix types, which provide the most flexibility in specifying memory layout:"
{ $subsection <float-blas-vector> } { $subsection <float-blas-vector> }
{ $subsection <double-blas-vector> } { $subsection <double-blas-vector> }
{ $subsection <float-complex-blas-vector> } { $subsection <complex-float-blas-vector> }
{ $subsection <double-complex-blas-vector> } { $subsection <complex-double-blas-vector> }
{ $subsection <float-blas-matrix> } { $subsection <float-blas-matrix> }
{ $subsection <double-blas-matrix> } { $subsection <double-blas-matrix> }
{ $subsection <float-complex-blas-matrix> } { $subsection <complex-float-blas-matrix> }
{ $subsection <double-complex-blas-matrix> } { $subsection <complex-double-blas-matrix> }
"For the simple case of creating a dense, zero-filled vector or matrix, simple empty object constructors are provided:" "For the simple case of creating a dense, zero-filled vector or matrix, simple empty object constructors are provided:"
{ $subsection <empty-vector> } { $subsection <empty-vector> }
{ $subsection <empty-matrix> } { $subsection <empty-matrix> }
"BLAS vectors and matrices can also be constructed from other Factor sequences:" "BLAS vectors and matrices can also be constructed from other Factor sequences:"
{ $subsection >float-blas-vector } { $subsection >float-blas-vector }
{ $subsection >double-blas-vector } { $subsection >double-blas-vector }
{ $subsection >float-complex-blas-vector } { $subsection >complex-float-blas-vector }
{ $subsection >double-complex-blas-vector } { $subsection >complex-double-blas-vector }
{ $subsection >float-blas-matrix } { $subsection >float-blas-matrix }
{ $subsection >double-blas-matrix } { $subsection >double-blas-matrix }
{ $subsection >float-complex-blas-matrix } { $subsection >complex-float-blas-matrix }
{ $subsection >double-complex-blas-matrix } ; { $subsection >complex-double-blas-matrix } ;
ARTICLE: "math.blas.matrices" "BLAS interface matrix operations" ARTICLE: "math.blas.matrices" "BLAS interface matrix operations"
"Transposing and slicing matrices:" "Transposing and slicing matrices:"
@ -87,8 +87,8 @@ HELP: blas-matrix-base
{ $list { $list
{ { $link float-blas-matrix } } { { $link float-blas-matrix } }
{ { $link double-blas-matrix } } { { $link double-blas-matrix } }
{ { $link float-complex-blas-matrix } } { { $link complex-float-blas-matrix } }
{ { $link double-complex-blas-matrix } } { { $link complex-double-blas-matrix } }
} }
"All of these subclasses share the same tuple layout:" "All of these subclasses share the same tuple layout:"
{ $list { $list
@ -104,14 +104,14 @@ HELP: float-blas-matrix
{ $class-description "A matrix of single-precision floating-point values. For details on the tuple layout, see " { $link blas-matrix-base } "." } ; { $class-description "A matrix of single-precision floating-point values. For details on the tuple layout, see " { $link blas-matrix-base } "." } ;
HELP: double-blas-matrix HELP: double-blas-matrix
{ $class-description "A matrix of double-precision floating-point values. For details on the tuple layout, see " { $link blas-matrix-base } "." } ; { $class-description "A matrix of double-precision floating-point values. For details on the tuple layout, see " { $link blas-matrix-base } "." } ;
HELP: float-complex-blas-matrix HELP: complex-float-blas-matrix
{ $class-description "A matrix of single-precision floating-point complex values. Complex values are stored in memory as two consecutive float values, real part then imaginary part. For details on the tuple layout, see " { $link blas-matrix-base } "." } ; { $class-description "A matrix of single-precision floating-point complex values. Complex values are stored in memory as two consecutive float values, real part then imaginary part. For details on the tuple layout, see " { $link blas-matrix-base } "." } ;
HELP: double-complex-blas-matrix HELP: complex-double-blas-matrix
{ $class-description "A matrix of double-precision floating-point complex values. Complex values are stored in memory as two consecutive float values, real part then imaginary part. For details on the tuple layout, see " { $link blas-matrix-base } "." } ; { $class-description "A matrix of double-precision floating-point complex values. Complex values are stored in memory as two consecutive float values, real part then imaginary part. For details on the tuple layout, see " { $link blas-matrix-base } "." } ;
{ {
float-blas-matrix double-blas-matrix float-complex-blas-matrix double-complex-blas-matrix float-blas-matrix double-blas-matrix complex-float-blas-matrix complex-double-blas-matrix
float-blas-vector double-blas-vector float-complex-blas-vector double-complex-blas-vector float-blas-vector double-blas-vector complex-float-blas-vector complex-double-blas-vector
} related-words } related-words
HELP: Mwidth HELP: Mwidth
@ -272,7 +272,7 @@ HELP: cmatrix{
{ 0.0 0.0 -1.0 3.0 } { 0.0 0.0 -1.0 3.0 }
{ 0.0 0.0 0.0 C{ 0.0 -1.0 } } { 0.0 0.0 0.0 C{ 0.0 -1.0 } }
} "> } } "> }
{ $description "Construct a literal " { $link float-complex-blas-matrix } ". Note that although BLAS matrices are stored in column-major order, the literal is specified in row-major order." } ; { $description "Construct a literal " { $link complex-float-blas-matrix } ". Note that although BLAS matrices are stored in column-major order, the literal is specified in row-major order." } ;
HELP: zmatrix{ HELP: zmatrix{
{ $syntax <" zmatrix{ { $syntax <" zmatrix{
@ -281,7 +281,7 @@ HELP: zmatrix{
{ 0.0 0.0 -1.0 3.0 } { 0.0 0.0 -1.0 3.0 }
{ 0.0 0.0 0.0 C{ 0.0 -1.0 } } { 0.0 0.0 0.0 C{ 0.0 -1.0 } }
} "> } } "> }
{ $description "Construct a literal " { $link double-complex-blas-matrix } ". Note that although BLAS matrices are stored in column-major order, the literal is specified in row-major order." } ; { $description "Construct a literal " { $link complex-double-blas-matrix } ". Note that although BLAS matrices are stored in column-major order, the literal is specified in row-major order." } ;
{ {
POSTPONE: smatrix{ POSTPONE: dmatrix{ POSTPONE: smatrix{ POSTPONE: dmatrix{

View File

@ -1,11 +1,13 @@
USING: accessors alien alien.c-types arrays byte-arrays combinators USING: accessors alien alien.c-types arrays byte-arrays combinators
combinators.short-circuit fry kernel locals macros combinators.short-circuit fry kernel locals macros
math math.blas.cblas math.blas.vectors math.blas.vectors.private math math.blas.ffi math.blas.vectors math.blas.vectors.private
math.complex math.functions math.order functors words math.complex math.functions math.order functors words
sequences sequences.merged sequences.private shuffle sequences sequences.merged sequences.private shuffle
specialized-arrays.direct.float specialized-arrays.direct.double specialized-arrays.direct.float specialized-arrays.direct.double
specialized-arrays.float specialized-arrays.double specialized-arrays.float specialized-arrays.double
parser prettyprint.backend prettyprint.custom ; specialized-arrays.direct.complex-float specialized-arrays.direct.complex-double
specialized-arrays.complex-float specialized-arrays.complex-double
parser prettyprint.backend prettyprint.custom ascii ;
IN: math.blas.matrices IN: math.blas.matrices
TUPLE: blas-matrix-base underlying ld rows cols transpose ; TUPLE: blas-matrix-base underlying ld rows cols transpose ;
@ -25,7 +27,7 @@ GENERIC: n*M.M+n*M! ( alpha A B beta C -- C=alpha*A.B+beta*C )
<PRIVATE <PRIVATE
: (blas-transpose) ( matrix -- integer ) : (blas-transpose) ( matrix -- integer )
transpose>> [ CblasTrans ] [ CblasNoTrans ] if ; transpose>> [ "T" ] [ "N" ] if ;
GENERIC: (blas-matrix-like) ( data ld rows cols transpose exemplar -- matrix ) GENERIC: (blas-matrix-like) ( data ld rows cols transpose exemplar -- matrix )
@ -38,19 +40,18 @@ GENERIC: (blas-matrix-like) ( data ld rows cols transpose exemplar -- matrix )
unless ; unless ;
:: (prepare-gemv) :: (prepare-gemv)
( alpha A x beta y >c-arg -- order A-trans m n alpha A-data A-ld x-data x-inc beta y-data y-inc ( alpha A x beta y -- A-trans m n alpha A-data A-ld x-data x-inc beta y-data y-inc
y ) y )
A x y (validate-gemv) A x y (validate-gemv)
CblasColMajor
A (blas-transpose) A (blas-transpose)
A rows>> A rows>>
A cols>> A cols>>
alpha >c-arg call alpha
A underlying>> A underlying>>
A ld>> A ld>>
x underlying>> x underlying>>
x inc>> x inc>>
beta >c-arg call beta
y underlying>> y underlying>>
y inc>> y inc>>
y ; inline y ; inline
@ -64,13 +65,12 @@ GENERIC: (blas-matrix-like) ( data ld rows cols transpose exemplar -- matrix )
unless ; unless ;
:: (prepare-ger) :: (prepare-ger)
( alpha x y A >c-arg -- order m n alpha x-data x-inc y-data y-inc A-data A-ld ( alpha x y A -- m n alpha x-data x-inc y-data y-inc A-data A-ld
A ) A )
x y A (validate-ger) x y A (validate-ger)
CblasColMajor
A rows>> A rows>>
A cols>> A cols>>
alpha >c-arg call alpha
x underlying>> x underlying>>
x inc>> x inc>>
y underlying>> y underlying>>
@ -89,21 +89,20 @@ GENERIC: (blas-matrix-like) ( data ld rows cols transpose exemplar -- matrix )
unless ; unless ;
:: (prepare-gemm) :: (prepare-gemm)
( alpha A B beta C >c-arg -- order A-trans B-trans m n k alpha A-data A-ld B-data B-ld beta C-data C-ld ( alpha A B beta C -- A-trans B-trans m n k alpha A-data A-ld B-data B-ld beta C-data C-ld
C ) C )
A B C (validate-gemm) A B C (validate-gemm)
CblasColMajor
A (blas-transpose) A (blas-transpose)
B (blas-transpose) B (blas-transpose)
C rows>> C rows>>
C cols>> C cols>>
A Mwidth A Mwidth
alpha >c-arg call alpha
A underlying>> A underlying>>
A ld>> A ld>>
B underlying>> B underlying>>
B ld>> B ld>>
beta >c-arg call beta
C underlying>> C underlying>>
C ld>> C ld>>
C f >>transpose ; inline C f >>transpose ; inline
@ -250,16 +249,18 @@ FUNCTOR: (define-blas-matrix) ( TYPE T U C -- )
VECTOR IS ${TYPE}-blas-vector VECTOR IS ${TYPE}-blas-vector
<VECTOR> IS <${TYPE}-blas-vector> <VECTOR> IS <${TYPE}-blas-vector>
>ARRAY IS >${TYPE}-array >ARRAY IS >${TYPE}-array
TYPE>ARG IS ${TYPE}>arg XGEMV IS ${T}GEMV
XGEMV IS cblas_${T}gemv XGEMM IS ${T}GEMM
XGEMM IS cblas_${T}gemm XGERU IS ${T}GER${U}
XGERU IS cblas_${T}ger${U} XGERC IS ${T}GER${C}
XGERC IS cblas_${T}ger${C}
MATRIX DEFINES-CLASS ${TYPE}-blas-matrix MATRIX DEFINES-CLASS ${TYPE}-blas-matrix
<MATRIX> DEFINES <${TYPE}-blas-matrix> <MATRIX> DEFINES <${TYPE}-blas-matrix>
>MATRIX DEFINES >${TYPE}-blas-matrix >MATRIX DEFINES >${TYPE}-blas-matrix
XMATRIX{ DEFINES ${T}matrix{
t [ T >lower ]
XMATRIX{ DEFINES ${t}matrix{
WHERE WHERE
@ -277,21 +278,16 @@ M: MATRIX (blas-vector-like)
drop <VECTOR> ; drop <VECTOR> ;
: >MATRIX ( arrays -- matrix ) : >MATRIX ( arrays -- matrix )
[ >ARRAY underlying>> ] (>matrix) [ >ARRAY underlying>> ] (>matrix) <MATRIX> ;
<MATRIX> ;
M: VECTOR n*M.V+n*V! M: VECTOR n*M.V+n*V!
[ TYPE>ARG ] (prepare-gemv) (prepare-gemv) [ XGEMV ] dip ;
[ XGEMV ] dip ;
M: MATRIX n*M.M+n*M! M: MATRIX n*M.M+n*M!
[ TYPE>ARG ] (prepare-gemm) (prepare-gemm) [ XGEMM ] dip ;
[ XGEMM ] dip ;
M: MATRIX n*V(*)V+M! M: MATRIX n*V(*)V+M!
[ TYPE>ARG ] (prepare-ger) (prepare-ger) [ XGERU ] dip ;
[ XGERU ] dip ;
M: MATRIX n*V(*)Vconj+M! M: MATRIX n*V(*)Vconj+M!
[ TYPE>ARG ] (prepare-ger) (prepare-ger) [ XGERC ] dip ;
[ XGERC ] dip ;
: XMATRIX{ \ } [ >MATRIX ] parse-literal ; parsing : XMATRIX{ \ } [ >MATRIX ] parse-literal ; parsing
@ -304,12 +300,12 @@ M: MATRIX pprint-delims
: define-real-blas-matrix ( TYPE T -- ) : define-real-blas-matrix ( TYPE T -- )
"" "" (define-blas-matrix) ; "" "" (define-blas-matrix) ;
: define-complex-blas-matrix ( TYPE T -- ) : define-complex-blas-matrix ( TYPE T -- )
"u" "c" (define-blas-matrix) ; "U" "C" (define-blas-matrix) ;
"float" "s" define-real-blas-matrix "float" "S" define-real-blas-matrix
"double" "d" define-real-blas-matrix "double" "D" define-real-blas-matrix
"float-complex" "c" define-complex-blas-matrix "complex-float" "C" define-complex-blas-matrix
"double-complex" "z" define-complex-blas-matrix "complex-double" "Z" define-complex-blas-matrix
>> >>