package slap

  1. Overview
  2. Docs

Slap.C provides 32-bit complex BLAS and LAPACK functions.

type num_type = Complex.t
type 'a trans3 = 'a Slap_common.trans3

A type of transpose parameters (Slap_common.normal, Slap_common.trans or Slap_common.conjtr).

type (+'n, +'cnt_or_dsc) vec = ('n, num_type, prec, 'cnt_or_dsc) Slap_vec.t

Vectors.

type (+'m, +'n, +'cnt_or_dsc) mat = ('m, 'n, num_type, prec, 'cnt_or_dsc) Slap_mat.t

Matrices.

type (+'n, +'cnt_or_dsc) rvec = ('n, float, rprec, 'cnt_or_dsc) Slap_vec.t

Real vectors. (In Slap.S and Slap.D, rvec is equal to vec.)

val rprec : (float, rprec) Bigarray.kind
module Vec : sig ... end
module Mat : sig ... end
val pp_num : Format.formatter -> num_type -> unit

A pretty-printer for elements in vectors and matrices.

val pp_vec : Format.formatter -> ('n, 'cnt_or_dsc) vec -> unit

A pretty-printer for column vectors.

val pp_mat : Format.formatter -> ('m, 'n, 'cnt_or_dsc) mat -> unit

A pretty-printer for matrices.

BLAS interface

Level 1

val swap : ('n, 'x_cd) vec -> ('n, 'y_cd) vec -> unit

swap x y swaps elements in x and y.

val scal : num_type -> ('n, 'cd) vec -> unit

scal c x multiplies all elements in x by scalar value c, and destructively assigns the result to x.

val copy : ?y:('n, 'y_cd) vec -> ('n, 'x_cd) vec -> ('n, 'y_cd) vec

copy ?y x copies x into y.

  • returns

    vector y, which is overwritten.

val nrm2 : ('n, 'cd) vec -> float

nrm2 x retruns the L2 norm of vector x: ||x||.

val axpy : ?alpha:num_type -> ('n, 'x_cd) vec -> ('n, 'y_cd) vec -> unit

axpy ?alpha x y executes y := alpha * x + y with scalar value alpha, and vectors x and y.

  • parameter alpha

    default = 1.0

val iamax : ('n, 'cd) vec -> int

iamax x returns the index of the maximum value of all elements in x.

val amax : ('n, 'cd) vec -> num_type

amax x finds the maximum value of all elements in x.

Level 2

val gemv : ?beta:num_type -> ?y:('m, 'y_cd) vec -> trans:(('a_m, 'a_n, 'a_cd) mat -> ('m, 'n, 'a_cd) mat) trans3 -> ?alpha:num_type -> ('a_m, 'a_n, 'a_cd) mat -> ('n, 'x_cd) vec -> ('m, 'y_cd) vec

gemv ?beta ?y ~trans ?alpha a x executes y := alpha * OP(a) * x + beta * y.

  • parameter beta

    default = 0.0

  • parameter alpha

    default = 1.0

val gbmv : m:'a_m Slap_size.t -> ?beta:num_type -> ?y:('m, 'y_cd) vec -> trans:(('a_m, 'a_n, 'a_cd) mat -> ('m, 'n, 'a_cd) mat) trans3 -> ?alpha:num_type -> (('a_m, 'a_n, 'kl, 'ku) Slap_size.geband, 'a_n, 'a_cd) mat -> 'kl Slap_size.t -> 'ku Slap_size.t -> ('n, 'x_cd) vec -> ('m, 'y_cd) vec

gbmv ~m ?beta ?y ~trans ?alpha a kl ku x computes y := alpha * OP(a) * x + beta * y where a is a band matrix stored in band storage.

  • returns

    vector y, which is overwritten.

  • parameter beta

    default = 0.0

  • parameter alpha

    default = 1.0

  • parameter kl

    the number of subdiagonals of a

  • parameter ku

    the number of superdiagonals of a

  • since 0.2.0
val symv : ?beta:num_type -> ?y:('n, 'y_cd) vec -> ?up:bool -> ?alpha:num_type -> ('n, 'n, 'a_cd) mat -> ('n, 'x_cd) vec -> ('n, 'y_cd) vec

symv ?beta ?y ?up ?alpha a x executes y := alpha * a * x + beta * y.

  • parameter beta

    default = 0.0

  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
  • parameter alpha

    default = 1.0

val trmv : trans:(('n, 'n, 'a_cd) mat -> ('n, 'n, 'a_cd) mat) trans3 -> ?diag:Slap_common.diag -> ?up:bool -> ('n, 'n, 'a_cd) mat -> ('n, 'x_cd) vec -> unit

trmv ~trans ?diag ?up a x executes x := OP(a) * x.

  • parameter diag

    default = `N

    • If diag = `U, then a is unit triangular;
    • If diag = `N, then a is not unit triangular.
  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
val trsv : trans:(('n, 'n, 'a_cd) mat -> ('n, 'n, 'a_cd) mat) trans3 -> ?diag:Slap_common.diag -> ?up:bool -> ('n, 'n, 'a_cd) mat -> ('n, 'b_cd) vec -> unit

trmv ~trans ?diag ?up a b solves linear system OP(a) * x = b and destructively assigns x to b.

  • parameter diag

    default = `N

    • If diag = `U, then a is unit triangular;
    • If diag = `N, then a is not unit triangular.
  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
val tpmv : trans:(('n, 'n, 'a_cd) mat -> ('n, 'n, 'a_cd) mat) trans3 -> ?diag:Slap_common.diag -> ?up:bool -> ('n Slap_size.packed, Slap_misc.cnt) vec -> ('n, 'x_cd) vec -> unit

tpmv ~trans ?diag ?up a x executes x := OP(a) * x where a is a packed triangular matrix.

  • parameter diag

    default = `N

    • If diag = `U, then a is unit triangular;
    • If diag = `N, then a is not unit triangular.
  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
  • since 0.2.0
val tpsv : trans:(('n, 'n, 'a_cd) mat -> ('n, 'n, 'a_cd) mat) trans3 -> ?diag:Slap_common.diag -> ?up:bool -> ('n Slap_size.packed, Slap_misc.cnt) vec -> ('n, 'x_cd) vec -> unit

tpsv ~trans ?diag ?up a b solves linear system OP(a) * x = b and destructively assigns x to b where a is a packed triangular matrix.

  • parameter diag

    default = `N

    • If diag = `U, then a is unit triangular;
    • If diag = `N, then a is not unit triangular.
  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
  • since 0.2.0

Level 3

val gemm : ?beta:num_type -> ?c:('m, 'n, 'c_cd) mat -> transa:(('a_m, 'a_k, 'a_cd) mat -> ('m, 'k, 'a_cd) mat) trans3 -> ?alpha:num_type -> ('a_m, 'a_k, 'a_cd) mat -> transb:(('b_k, 'b_n, 'b_cd) mat -> ('k, 'n, 'b_cd) mat) trans3 -> ('b_k, 'b_n, 'b_cd) mat -> ('m, 'n, 'c_cd) mat

gemm ?beta ?c ~transa ?alpha a ~transb b executes c := alpha * OP(a) * OP(b) + beta * c.

  • parameter beta

    default = 0.0

  • parameter alpha

    default = 1.0

val symm : side:('k, 'm, 'n) Slap_common.side -> ?up:bool -> ?beta:num_type -> ?c:('m, 'n, 'c_cd) mat -> ?alpha:num_type -> ('k, 'k, 'a_cd) mat -> ('m, 'n, 'b_cd) mat -> ('m, 'n, 'c_cd) mat

symm ~side ?up ?beta ?c ?alpha a b executes

where a is a symmterix matrix, and b and c are general matrices.

  • parameter side

    the side flag to specify direction of multiplication of a and b.

  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
  • parameter beta

    default = 0.0

  • parameter alpha

    default = 1.0

val trmm : side:('k, 'm, 'n) Slap_common.side -> ?up:bool -> transa:(('k, 'k, 'a_cd) mat -> ('k, 'k, 'a_cd) mat) trans3 -> ?diag:Slap_common.diag -> ?alpha:num_type -> a:('k, 'k, 'a_cd) mat -> ('m, 'n, 'b_cd) mat -> unit

trmm ~side ?up ~transa ?diag ?alpha ~a b executes

where a is a triangular matrix, and b is a general matrix.

  • parameter side

    the side flag to specify direction of multiplication of a and b.

  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
  • parameter diag

    default = `N

    • If diag = `U, then a is unit triangular;
    • If diag = `N, then a is not unit triangular.
  • parameter alpha

    default = 1.0

val trsm : side:('k, 'm, 'n) Slap_common.side -> ?up:bool -> transa:(('k, 'k, 'a_cd) mat -> ('k, 'k, 'a_cd) mat) trans3 -> ?diag:Slap_common.diag -> ?alpha:num_type -> a:('k, 'k, 'a_cd) mat -> ('m, 'n, 'b_cd) mat -> unit

trsm ~side ?up ~transa ?diag ?alpha ~a b solves a system of linear equations

where a is a triangular matrix, and b is a general matrix. The solution x is returned by b.

  • parameter side

    the side flag to specify direction of multiplication of a and b.

  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
  • parameter diag

    default = `N

    • If diag = `U, then a is unit triangular;
    • If diag = `N, then a is not unit triangular.
  • parameter alpha

    default = 1.0

val syrk : ?up:bool -> ?beta:num_type -> ?c:('n, 'n, 'c_cd) mat -> trans:(('a_n, 'a_k, 'a_cd) mat -> ('n, 'k, 'a_cd) mat) Slap_common.trans2 -> ?alpha:num_type -> ('a_n, 'a_k, 'a_cd) mat -> ('n, 'n, 'c_cd) mat

syrk ?up ?beta ?c ~trans ?alpha a executes

where a is a general matrix and c is a symmetric matrix.

  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
  • parameter beta

    default = 0.0

  • parameter trans

    the transpose flag for a

  • parameter alpha

    default = 1.0

val syr2k : ?up:bool -> ?beta:num_type -> ?c:('n, 'n, 'c_cd) mat -> trans:(('p, 'q, _) mat -> ('n, 'k, _) mat) Slap_common.trans2 -> ?alpha:num_type -> ('p, 'q, 'a_cd) mat -> ('p, 'q, 'b_cd) mat -> ('n, 'n, 'c_cd) mat

syr2k ?up ?beta ?c ~trans ?alpha a b computes

with symmetric matrix c, and general matrices a and b.

  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
  • parameter beta

    default = 0.0

  • parameter trans

    the transpose flag for a

  • parameter alpha

    default = 1.0

LAPACK interface

Auxiliary routines

val lacpy : ?uplo:[ `L | `U ] -> ?b:('m, 'n, 'b_cd) mat -> ('m, 'n, 'a_cd) mat -> ('m, 'n, 'b_cd) mat

lacpy ?uplo ?b a copies the matrix a into the matrix b.

  • If uplo is omitted, all elements in a is copied.
  • If uplo is `U, the upper trapezoidal part of a is copied.
  • If uplo is `L, the lower trapezoidal part of a is copied.
  • returns

    b, which is overwritten.

  • parameter uplo

    default = all elements in a is copied.

  • parameter b

    default = a fresh matrix.

val lassq : ?scale:float -> ?sumsq:float -> ('n, 'cd) vec -> float * float

lassq ?scale ?sumsq x

  • returns

    (scl, smsq) where scl and smsq satisfy scl^2 * smsq = x1^2 + x2^2 + ... + xn^2 + scale^2 * smsq.

  • parameter scale

    default = 0.0

  • parameter sumsq

    default = 1.0

type larnv_liseed = Slap_size.four
val larnv : ?idist:[ `Normal | `Uniform0 | `Uniform1 ] -> ?iseed:(larnv_liseed, Slap_misc.cnt) Slap_common.int32_vec -> x:('n, Slap_misc.cnt) vec -> unit -> ('n, 'cnt) vec

larnv ?idist ?iseed ~x () generates a random vector with the random distribution specified by idist and random seed iseed.

  • returns

    vector x, which is overwritten.

  • parameter idist

    default = `Normal

  • parameter iseed

    a four-dimensional integer vector with all ones.

type ('m, 'a) lange_min_lwork
val lange_min_lwork : 'm Slap_size.t -> 'a Slap_common.norm4 -> ('m, 'a) lange_min_lwork Slap_size.t

lange_min_lwork m norm computes the minimum length of workspace for lange routine. m is the number of rows in a matrix, and norm is the sort of matrix norms.

val lange : ?norm:'a Slap_common.norm4 -> ?work:('lwork, Slap_misc.cnt) rvec -> ('m, 'n, 'cd) mat -> float

lange ?norm ?work a

  • returns

    the norm of matrix a.

  • parameter work

    default = an optimum-length vector.

val lauum : ?up:bool -> ('n, 'n, 'cd) mat -> unit

lauum ?up a computes

  • U * U^T where U is the upper triangular part of matrix a if up is true.
  • L^T * L where L is the lower triangular part of matrix a if up is false.

The upper or lower triangular part is overwritten.

  • parameter up

    default = true.

Linear equations (computational routines)

val getrf : ?ipiv:(('m, 'n) Slap_size.min, Slap_misc.cnt) Slap_common.int32_vec -> ('m, 'n, 'cd) mat -> (('m, 'n) Slap_size.min, 'cnt) Slap_common.int32_vec

getrf ?ipiv a computes LU factorization of matrix a using partial pivoting with row interchanges: a = P * L * U where P is a permutation matrix, and L and U are lower and upper triangular matrices, respectively. the permutation matrix is returned in ipiv.

  • returns

    vector ipiv, which is overwritten.

  • raises Failure

    if the matrix is singular.

val getrs : ?ipiv:(('n, 'n) Slap_size.min, Slap_misc.cnt) Slap_common.int32_vec -> trans:(('n, 'n, 'a_cd) mat -> ('n, 'n, 'a_cd) mat) trans3 -> ('n, 'n, 'a_cd) mat -> ('n, 'n, 'b_cd) mat -> unit

getrs ?ipiv trans a b solves systems of linear equations OP(a) * x = b where a a 'n-by-'n general matrix, each column of matrix b is the r.h.s. vector, and each column of matrix x is the corresponding solution. The solution x is returned in b.

  • parameter ipiv

    a result of gesv or getrf. It is internally computed by getrf if omitted.

  • raises Failure

    if the matrix is singular.

type 'n getri_min_lwork
val getri_min_lwork : 'n Slap_size.t -> 'n getri_min_lwork Slap_size.t

getri_min_lwork n computes the minimum length of workspace for getri routine. n is the number of columns or rows in a matrix.

val getri_opt_lwork : ('n, 'n, 'cd) mat -> (module Slap_size.SIZE)

getri_opt_lwork a computes the optimal length of workspace for getri routine.

val getri : ?ipiv:(('n, 'n) Slap_size.min, Slap_misc.cnt) Slap_common.int32_vec -> ?work:('lwork, Slap_misc.cnt) vec -> ('n, 'n, 'cd) mat -> unit

getri ?ipiv ?work a computes the inverse of general matrix a by LU-factorization. The inverse matrix is returned in a.

  • parameter ipiv

    a result of gesv or getrf. It is internally computed by getrf if omitted.

  • parameter work

    default = an optimum-length vector.

  • raises Failure

    if the matrix is singular.

type sytrf_min_lwork
val sytrf_min_lwork : unit -> sytrf_min_lwork Slap_size.t

sytrf_min_lwork () computes the minimum length of workspace for sytrf routine.

val sytrf_opt_lwork : ?up:bool -> ('n, 'n, 'cd) mat -> (module Slap_size.SIZE)

sytrf_opt_lwork ?up a computes the optimal length of workspace for sytrf routine.

  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
val sytrf : ?up:bool -> ?ipiv:('n, Slap_misc.cnt) Slap_common.int32_vec -> ?work:('lwork, Slap_misc.cnt) vec -> ('n, 'n, 'cd) mat -> ('n, 'cnt) Slap_common.int32_vec

sytrf ?up ?ipiv ?work a factorizes symmetric matrix a using the Bunch-Kaufman diagonal pivoting method:

  • a = P * U * D * U^T * P^T if up = true;
  • a = P * L * D * L^T * P^T if up = false

where P is a permutation matrix, U and L are upper and lower triangular matrices with unit diagonal, and D is a symmetric block-diagonal matrix. The permutation matrix is returned in ipiv.

  • returns

    vector ipiv, which is overwritten.

  • parameter up

    default = true

  • parameter work

    default = an optimum-length vector.

val sytrs : ?up:bool -> ?ipiv:('n, Slap_misc.cnt) Slap_common.int32_vec -> ('n, 'n, 'a_cd) mat -> ('n, 'nrhs, 'b_cd) mat -> unit

sytrs ?up ?ipiv a b solves systems of linear equations a * x = b where a is a symmetric matrix, each column of matrix b is the r.h.s. vector, and each column of matrix x is the corresponding solution. The solution x is returned in b.

This routine uses the Bunch-Kaufman diagonal pivoting method:

  • a = P * U * D * U^T * P^T if up = true;
  • a = P * L * D * L^T * P^T if up = false

where P is a permutation matrix, U and L are upper and lower triangular matrices with unit diagonal, and D is a symmetric block-diagonal matrix.

  • parameter up

    default = true

  • parameter ipiv

    a result of sytrf. It is internally computed by sytrf if omitted.

type 'n sytri_min_lwork
val sytri_min_lwork : 'n Slap_size.t -> 'n sytri_min_lwork Slap_size.t

sytri_min_lwork () computes the minimum length of workspace for sytri routine.

val sytri : ?up:bool -> ?ipiv:('n, Slap_misc.cnt) Slap_common.int32_vec -> ?work:('lwork, Slap_misc.cnt) vec -> ('n, 'n, 'cd) mat -> unit

sytri ?up ?ipiv ?work a computes the inverse of symmetric matrix a using the Bunch-Kaufman diagonal pivoting method:

  • a = P * U * D * U^T * P^T if up = true;
  • a = P * L * D * L^T * P^T if up = false

where P is a permutation matrix, U and L are upper and lower triangular matrices with unit diagonal, and D is a symmetric block-diagonal matrix.

  • parameter up

    default = true

  • parameter ipiv

    a result of sytrf. It is internally computed by sytrf if omitted.

  • parameter work

    default = an optimum-length vector.

val potrf : ?up:bool -> ?jitter:num_type -> ('n, 'n, 'cd) mat -> unit

potrf ?up ?jitter a computes the Cholesky factorization of symmetrix (Hermitian) positive-definite matrix a:

  • a = U^T * U (real) or a = U^H * U (complex) if up = true;
  • a = L * L^T (real) or a = L * L^H (complex) if up = false

where U and L are upper and lower triangular matrices, respectively. Either of them is returned in the upper or lower triangular part of a, as specified by up.

  • parameter up

    default = true

  • parameter jitter

    default = nothing

val potrs : ?up:bool -> ('n, 'n, 'a_cd) mat -> ?factorize:bool -> ?jitter:num_type -> ('n, 'nrhs, 'b_cd) mat -> unit

potrf ?up a ?jitter b solves systems of linear equations a * x = b using the Cholesky factorization of symmetrix (Hermitian) positive-definite matrix a:

  • a = U^T * U (real) or a = U^H * U (complex) if up = true;
  • a = L * L^T (real) or a = L * L^H (complex) if up = false

where U and L are upper and lower triangular matrices, respectively.

  • parameter up

    default = true

  • parameter factorize

    default = true (potrf is called implicitly)

  • parameter jitter

    default = nothing

val potri : ?up:bool -> ?factorize:bool -> ?jitter:num_type -> ('n, 'n, 'cd) mat -> unit

potrf ?up ?jitter a computes the inverse of symmetrix (Hermitian) positive-definite matrix a using the Cholesky factorization:

  • a = U^T * U (real) or a = U^H * U (complex) if up = true;
  • a = L * L^T (real) or a = L * L^H (complex) if up = false

where U and L are upper and lower triangular matrices, respectively.

  • parameter up

    default = true

  • parameter factorize

    default = true (potrf is called implicitly)

  • parameter jitter

    default = nothing

val trtrs : ?up:bool -> trans:(('n, 'n, 'a_cd) mat -> ('n, 'n, 'a_cd) mat) trans3 -> ?diag:Slap_common.diag -> ('n, 'n, 'a_cd) mat -> ('n, 'nrhs, 'b_cd) mat -> unit

trtrs ?up trans ?diag a b solves systems of linear equations OP(a) * x = b where a is a triangular matrix of order 'n, each column of matrix b is the r.h.s vector, and each column of matrix x is the corresponding solution. The solution x is returned in b.

  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
  • parameter diag

    default = `N

    • If diag = `U, then a is unit triangular;
    • If diag = `N, then a is not unit triangular.
val tbtrs : kd:'kd Slap_size.t -> ?up:bool -> trans:(('n, 'n, 'a_cd) mat -> ('n, 'n, 'a_cd) mat) trans3 -> ?diag:Slap_common.diag -> (('n, 'kd) Slap_size.syband, 'n, 'a_cd) mat -> ('n, 'nrhs, 'b_cd) mat -> unit

tbtrs ~kd ?up ~trans ?diag ab b solves systems of linear equations OP(A) * x = b where A is a triangular band matrix with kd subdiagonals, each column of matrix b is the r.h.s vector, and each column of matrix x is the corresponding solution. Matrix A is stored into ab in band storage format. The solution x is returned in b.

  • parameter kd

    the number of subdiagonals or superdiagonals in A.

  • parameter up

    default = true

    • If up = true, then the upper triangular part of A is used;
    • If up = false, then the lower triangular part of A is used.
  • parameter diag

    default = `N

    • If diag = `U, then A is unit triangular;
    • If diag = `N, then A is not unit triangular.
  • raises Failure

    if the matrix A is singular.

  • since 0.2.0
val trtri : ?up:bool -> ?diag:Slap_common.diag -> ('n, 'n, 'cd) mat -> unit

trtri ?up ?diag a computes the inverse of triangular matrix a. The inverse matrix is returned in a.

  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
  • parameter diag

    default = `N

    • If diag = `U, then a is unit triangular;
    • If diag = `N, then a is not unit triangular.
  • raises Failure

    if the matrix a is singular.

  • since 0.2.0
type 'n geqrf_min_lwork
val geqrf_min_lwork : n:'n Slap_size.t -> 'n geqrf_min_lwork Slap_size.t

geqrf_min_lwork ~n computes the minimum length of workspace for geqrf routine. n is the number of columns in a matrix.

val geqrf_opt_lwork : ('m, 'n, 'cd) mat -> (module Slap_size.SIZE)

geqrf_opt_lwork a computes the optimum length of workspace for geqrf routine.

val geqrf : ?work:('lwork, Slap_misc.cnt) vec -> ?tau:(('m, 'n) Slap_size.min, Slap_misc.cnt) vec -> ('m, 'n, 'cd) mat -> (('m, 'n) Slap_size.min, 'cnt) vec

geqrf ?work ?tau a computes the QR factorization of general matrix a: a = Q * R where Q is an orthogonal (unitary) matrix and R is an upper triangular matrix. R is returned in a. This routine does not generate Q explicitly. It is generated by orgqr.

  • returns

    vector tau, which is overwritten.

  • parameter work

    default = an optimum-length vector.

Linear equations (simple drivers)

val gesv : ?ipiv:('n, Slap_misc.cnt) Slap_common.int32_vec -> ('n, 'n, 'a_cd) mat -> ('n, 'nrhs, 'b_cd) mat -> unit

gesv ?ipiv a b solves a system of linear equations a * x = b where a is a 'n-by-'n general matrix, each column of matrix b is the r.h.s. vector, and each column of matrix x is the corresponding solution.

This routine uses LU factorization: a = P * L * U with permutation matrix P, a lower triangular matrix L and an upper triangular matrix U. By this function, the upper triangular part of a is replaced by U, the lower triangular part by L, and the solution x is returned in b.

  • raises Failure

    if the matrix is singular.

  • since 0.2.0
val gbsv : ?ipiv:('n, Slap_misc.cnt) Slap_common.int32_vec -> (('n, 'n, 'kl, 'ku) Slap_size.luband, 'n, 'a_cd) mat -> 'kl Slap_size.t -> 'ku Slap_size.t -> ('n, 'nrhs, 'b_cd) mat -> unit

gbsv ?ipiv ab kl ku b solves a system of linear equations A * X = B where A is a 'n-by-'n band matrix, each column of matrix B is the r.h.s. vector, and each column of matrix X is the corresponding solution. The matrix A with kl subdiagonals and ku superdiagonals is stored into ab in band storage format for LU factorizaion.

This routine uses LU factorization: A = P * L * U with permutation matrix P, a lower triangular matrix L and an upper triangular matrix U. By this function, the upper triangular part of A is replaced by U, the lower triangular part by L, and the solution X is returned in B.

  • raises Failure

    if the matrix is singular.

  • since 0.2.0
val posv : ?up:bool -> ('n, 'n, 'a_cd) mat -> ('n, 'nrhs, 'b_cd) mat -> unit

posv ?up a b solves systems of linear equations a * x = b where a is a 'n-by-'n symmetric positive-definite matrix, each column of matrix b is the r.h.s vector, and each column of matrix x is the corresponding solution. The solution x is returned in b.

The Cholesky decomposition is used:

  • If up = true, then a = U^T * U (real) or a = U^H * U (complex)
  • If up = false, then a = L^T * L (real) or a = L^H * L (complex)

where U and L are the upper and lower triangular matrices, respectively.

  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
  • raises Failure

    if the matrix is singular.

  • since 0.2.0
val ppsv : ?up:bool -> ('n Slap_size.packed, Slap_misc.cnt) vec -> ('n, 'nrhs, 'b_cd) mat -> unit

ppsv ?up a b solves systems of linear equations a * x = b where a is a 'n-by-'n symmetric positive-definite matrix stored in packed format, each column of matrix b is the r.h.s vector, and each column of matrix x is the corresponding solution. The solution x is returned in b.

The Cholesky decomposition is used:

  • If up = true, then a = U^T * U (real) or a = U^H * U (complex)
  • If up = false, then a = L^T * L (real) or a = L^H * L (complex)

where U and L are the upper and lower triangular matrices, respectively.

  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
  • raises Failure

    if the matrix is singular.

  • since 0.2.0
val pbsv : ?up:bool -> kd:'kd Slap_size.t -> (('n, 'kd) Slap_size.syband, 'n, 'ab_cd) mat -> ('n, 'nrhs, 'b_cd) mat -> unit

pbsv ?up ~kd ab b solves systems of linear equations ab * x = b where ab is a 'n-by-'n symmetric positive-definite band matrix with kd subdiangonals, stored in band storage format, each column of matrix b is the r.h.s vector, and each column of matrix x is the corresponding solution. The solution x is returned in b.

This routine uses the Cholesky decomposition:

  • If up = true, then ab = U^T * U (real) or ab = U^H * U (complex)
  • If up = false, then ab = L^T * L (real) or ab = L^H * L (complex)

where U and L are the upper and lower triangular matrices, respectively.

  • parameter up

    default = true

    • If up = true, then the upper triangular part of ab is used;
    • If up = false, then the lower triangular part of ab is used.
  • parameter kd

    the number of subdiagonals or superdiagonals in ab.

  • raises Failure

    if the matrix is singular.

  • since 0.2.0
val ptsv : ('n, Slap_misc.cnt) vec -> ('n Slap_size.p, Slap_misc.cnt) vec -> ('n, 'nrhs, 'b_cd) mat -> unit

ptsv d e b solves systems of linear equations A * x = b where A is a 'n-by-'n symmetric positive-definite tridiagonal matrix with diagonal elements d and subdiagonal elements e, each column of matrix b is the r.h.s vector, and each column of matrix x is the corresponding solution. The solution x is returned in b.

This routine uses the Cholesky decomposition: A = L^T * L (real) or A = L^H * L (complex) where L is a lower triangular matrix.

  • raises Failure

    if the matrix is singular.

  • since 0.2.0
val sysv_opt_lwork : ?up:bool -> ('n, 'n, 'a_cd) mat -> ('n, 'nrhs, 'b_cd) mat -> (module Slap_size.SIZE)

sysv_opt_lwork ?up a b computes the optimal length of workspace for sysv routine.

val sysv : ?up:bool -> ?ipiv:('n, Slap_misc.cnt) Slap_common.int32_vec -> ?work:('lwork, Slap_misc.cnt) vec -> ('n, 'n, 'a_cd) mat -> ('n, 'nrhs, 'b_cd) mat -> unit

sysv ?up ?ipiv ?work a b solves systems of linear equations a * x = b where a is a 'n-by-'n symmetric matrix, each column of matrix b is the r.h.s. vector, and each column of matrix x is the corresponding solution. The solution x is returned in b.

The diagonal pivoting method is used:

  • If up = true, then a = U * D * U^T
  • If up = false, then a = L * D * L^T where U and L are the upper and lower triangular matrices, respectively.
  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
  • parameter ipiv

    a result of sytrf. It is internally computed by sytrf if omitted.

  • parameter work

    default = an optimum-length vector.

  • raises Failure

    if the matrix is singular.

  • since 0.2.0
val spsv : ?up:bool -> ?ipiv:('n, Slap_misc.cnt) Slap_common.int32_vec -> ('n Slap_size.packed, Slap_misc.cnt) vec -> ('n, 'nrhs, 'b_cd) mat -> unit

spsv ?up a b solves systems of linear equations a * x = b where a is a 'n-by-'n symmetric matrix stored in packed format, each column of matrix b is the r.h.s. vector, and each column of matrix x is the corresponding solution. The solution x is returned in b.

The diagonal pivoting method is used:

  • If up = true, then a = U * D * U^T
  • If up = false, then a = L * D * L^T where U and L are the upper and lower triangular matrices, respectively.
  • parameter up

    default = true

    • If up = true, then the upper triangular part of a is used;
    • If up = false, then the lower triangular part of a is used.
  • parameter ipiv

    a result of sytrf. It is internally computed by sytrf if omitted.

  • raises Failure

    if the matrix is singular.

  • since 0.2.0

Least squares (simple drivers)

type ('m, 'n, 'nrhs) gels_min_lwork
val gels_min_lwork : m:'m Slap_size.t -> n:'n Slap_size.t -> nrhs:'nrhs Slap_size.t -> ('m, 'n, 'nrhs) gels_min_lwork Slap_size.t

gels_min_lwork ~n computes the minimum length of workspace for gels routine.

  • parameter m

    the number of rows in a matrix.

  • parameter n

    the number of columns in a matrix.

  • parameter nrhs

    the number of right hand sides.

val gels_opt_lwork : trans:(('am, 'an, 'a_cd) mat -> ('m, 'n, 'a_cd) mat) Slap_common.trans2 -> ('am, 'an, 'a_cd) mat -> ('m, 'nrhs, 'b_cd) mat -> (module Slap_size.SIZE)

gels_opt_lwork ~trans a b computes the optimum length of workspace for gels routine.

val gels : ?work:('work, Slap_misc.cnt) vec -> trans:(('am, 'an, 'a_cd) mat -> ('m, 'n, 'a_cd) mat) Slap_common.trans2 -> ('am, 'an, 'a_cd) mat -> ('m, 'nrhs, 'b_cd) mat -> unit

gels ?work ~trans a b solves an overdetermined or underdetermined system of linear equations using QR or LU factorization.

  • If trans = Slap_common.normal and 'm >= 'n: find the least square solution to an overdetermined system by minimizing ||b - A * x||^2.
  • If trans = Slap_common.normal and 'm < 'n: find the minimum norm solution to an underdetermined system a * x = b.
  • If trans = Slap_common.trans, and 'm >= 'n: find the minimum norm solution to an underdetermined system a^H * x = b.
  • If trans = Slap_common.trans and 'm < 'n: find the least square solution to an overdetermined system by minimizing ||b - A^H * x||^2.
  • parameter work

    default = an optimum-length vector.

  • parameter trans

    the transpose flag for a.

BLAS interface

Level 1

val dotu : ('n, 'x_cd) vec -> ('n, 'y_cd) vec -> num_type

dotc x y computes x^T y.

  • returns

    an inner product of given two vectors.

  • since 2.0.0
val dotc : ('n, 'x_cd) vec -> ('n, 'y_cd) vec -> num_type

dotc x y computes x^H y.

  • returns

    an inner product of a conjugated vector with another vector.

  • since 2.0.0
OCaml

Innovation. Community. Security.