module Lacaml_C:sig
..end
Lacaml.C
contains linear algebra routines for
complex numbers (precision: complex32). It is recommended to use this
module by writing
open Lacaml.C
at the top of your file.typeprec =
Bigarray.complex32_elt
typenum_type =
Complex.t
typevec =
(Complex.t, Bigarray.complex32_elt, Bigarray.fortran_layout)
Bigarray.Array1.t
typervec =
(float, Bigarray.float32_elt, Bigarray.fortran_layout) Bigarray.Array1.t
typemat =
(Complex.t, Bigarray.complex32_elt, Bigarray.fortran_layout)
Bigarray.Array2.t
typetrans3 =
[ `C | `N | `T ]
val prec : (Complex.t, Bigarray.complex32_elt) Bigarray.kind
C
. Allows to write precision
independent code.module Vec:sig
..end
module Mat:sig
..end
val pp_num : Format.formatter -> Complex.t -> unit
pp_num ppf el
is equivalent to fprintf ppf "(%G, %Gi)"
el.re el.im
.val pp_vec : (Complex.t, 'a) Lacaml_io.pp_vec
val pp_mat : (Complex.t, 'a) Lacaml_io.pp_mat
val dotu : ?n:int ->
?ofsx:int ->
?incx:int ->
Lacaml_complex32.vec ->
?ofsy:int -> ?incy:int -> Lacaml_complex32.vec -> Lacaml_complex32.num_type
dotu ?n ?ofsx ?incx x ?ofsy ?incy y
see BLAS documentation!n
: default = greater n s.t. ofsx+(n-1)(abs incx) <= dim x
ofsx
: default = 1incx
: default = 1ofsy
: default = 1incy
: default = 1val dotc : ?n:int ->
?ofsx:int ->
?incx:int ->
Lacaml_complex32.vec ->
?ofsy:int -> ?incy:int -> Lacaml_complex32.vec -> Lacaml_complex32.num_type
dotc ?n ?ofsx ?incx x ?ofsy ?incy y
see BLAS documentation!n
: default = greater n s.t. ofsx+(n-1)(abs incx) <= dim x
ofsx
: default = 1incx
: default = 1ofsy
: default = 1incy
: default = 1val lansy_min_lwork : int -> Lacaml_common.norm4 -> int
lansy_min_lwork m norm
lansy
-function.val lansy : ?n:int ->
?up:bool ->
?norm:Lacaml_common.norm4 ->
?work:Lacaml_complex32.rvec ->
?ar:int -> ?ac:int -> Lacaml_complex32.mat -> float
lansy ?n ?up ?norm ?work ?ar ?ac a
see LAPACK documentation!n
: default = number of columns of matrix a
up
: default = true (reference upper triangular part of a
)norm
: default = `Owork
: default = allocated work space for norm `Ival gecon_min_lwork : int -> int
gecon_min_lwork n
gecon
-function.val gecon_min_lrwork : int -> int
gecon_min_lrwork n
gecon
-function.val gecon : ?n:int ->
?norm:Lacaml_common.norm2 ->
?anorm:float ->
?work:Lacaml_complex32.vec ->
?rwork:Lacaml_complex32.rvec ->
?ar:int -> ?ac:int -> Lacaml_complex32.mat -> float
gecon ?n ?norm ?anorm ?work ?rwork ?ar ?ac a
a
n
: default = available number of columns of matrix a
norm
: default = 1-normanorm
: default = norm of the matrix a
as returned by lange
work
: default = automatically allocated workspacerwork
: default = automatically allocated workspacear
: default = 1ac
: default = 1val sycon_min_lwork : int -> int
sycon_min_lwork n
sycon
-function.val sycon : ?n:int ->
?up:bool ->
?ipiv:Lacaml_common.int32_vec ->
?anorm:float ->
?work:Lacaml_complex32.vec ->
?ar:int -> ?ac:int -> Lacaml_complex32.mat -> float
sycon ?n ?up ?ipiv ?anorm ?work ?ar ?ac a
a
n
: default = available number of columns of matrix a
up
: default = upper triangle of the factorization of a
is storedipiv
: default = vec of length n
anorm
: default = 1-norm of the matrix a
as returned by lange
work
: default = automatically allocated workspaceval pocon_min_lwork : int -> int
pocon_min_lwork n
pocon
-function.val pocon_min_lrwork : int -> int
pocon_min_lrwork n
pocon
-function.val pocon : ?n:int ->
?up:bool ->
?anorm:float ->
?work:Lacaml_complex32.vec ->
?rwork:Lacaml_complex32.rvec ->
?ar:int -> ?ac:int -> Lacaml_complex32.mat -> float
pocon ?n ?up ?anorm ?work ?rwork ?ar ?ac a
a
n
: default = available number of columns of matrix a
up
: default = upper triangle of Cholesky factorization
of a
is storedanorm
: default = 1-norm of the matrix a
as returned by lange
work
: default = automatically allocated workspacerwork
: default = automatically allocated workspaceval gees : ?n:int ->
?jobvs:Lacaml_common.schur_vectors ->
?sort:Lacaml_common.eigen_value_sort ->
?w:Lacaml_complex32.vec ->
?vsr:int ->
?vsc:int ->
?vs:Lacaml_complex32.mat ->
?work:Lacaml_complex32.vec ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat -> int * Lacaml_complex32.vec * Lacaml_complex32.mat
gees ?n ?jobvs ?sort ?w ?vsr ?vsc ?vs ?work ?ar ?ac a
See gees
-function for details about arguments.val gesvd_min_lwork : m:int -> n:int -> int
gesvd_min_lwork ~m ~n
gesvd
-function for matrices with m
rows and n
columns.val gesvd_lrwork : m:int -> n:int -> int
gesvd_lrwork m n
gesvd
-function.val gesvd_opt_lwork : ?m:int ->
?n:int ->
?jobu:Lacaml_common.svd_job ->
?jobvt:Lacaml_common.svd_job ->
?s:Lacaml_complex32.rvec ->
?ur:int ->
?uc:int ->
?u:Lacaml_complex32.mat ->
?vtr:int ->
?vtc:int ->
?vt:Lacaml_complex32.mat -> ?ar:int -> ?ac:int -> Lacaml_complex32.mat -> int
val gesvd : ?m:int ->
?n:int ->
?jobu:Lacaml_common.svd_job ->
?jobvt:Lacaml_common.svd_job ->
?s:Lacaml_complex32.rvec ->
?ur:int ->
?uc:int ->
?u:Lacaml_complex32.mat ->
?vtr:int ->
?vtc:int ->
?vt:Lacaml_complex32.mat ->
?work:Lacaml_complex32.vec ->
?rwork:Lacaml_complex32.rvec ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
Lacaml_complex32.rvec * Lacaml_complex32.mat * Lacaml_complex32.mat
val geev_min_lwork : int -> int
geev_min_lwork n
geev
-function.val geev_min_lrwork : int -> int
geev_min_lrwork n
geev
-function.val geev_opt_lwork : ?n:int ->
?vlr:int ->
?vlc:int ->
?vl:Lacaml_complex32.mat option ->
?vrr:int ->
?vrc:int ->
?vr:Lacaml_complex32.mat option ->
?ofsw:int ->
?w:Lacaml_complex32.vec -> ?ar:int -> ?ac:int -> Lacaml_complex32.mat -> int
geev ?work ?rwork ?n ?vlr ?vlc ?vl ?vrr ?vrc ?vr ?ofsw w ?ar ?ac a
See geev
-function for details about arguments.val geev : ?n:int ->
?work:Lacaml_complex32.vec ->
?rwork:Lacaml_complex32.vec ->
?vlr:int ->
?vlc:int ->
?vl:Lacaml_complex32.mat option ->
?vrr:int ->
?vrc:int ->
?vr:Lacaml_complex32.mat option ->
?ofsw:int ->
?w:Lacaml_complex32.vec ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
Lacaml_complex32.mat * Lacaml_complex32.vec * Lacaml_complex32.mat
geev ?work ?rwork ?n
?vlr ?vlc ?vl
?vrr ?vrc ?vr
?ofsw w
?ar ?ac a
Failure
if the function fails to converge(lv, w, rv)
, where lv
and rv
correspond to the left and
right eigenvectors respectively, w
to the eigenvalues. lv
(rv
)
is the empty matrix if vl
(vr
) is set to None
.n
: default = available number of columns of matrix a
work
: default = automatically allocated workspacerwork
: default = automatically allocated workspacevl
: default = Automatically allocated left eigenvectors.
Pass None
if you do not want to compute them,
Some lv
if you want to provide the storage.
You can set vlr
, vlc
in the last case.
(See LAPACK GEEV docs for details about storage of complex eigenvectors)vr
: default = Automatically allocated right eigenvectors.
Pass None
if you do not want to compute them,
Some rv
if you want to provide the storage.
You can set vrr
, vrc
in the last case.w
: default = automatically allocate eigenvaluesval swap : ?n:int ->
?ofsx:int ->
?incx:int ->
x:Lacaml_complex32.vec ->
?ofsy:int -> ?incy:int -> Lacaml_complex32.vec -> unit
swap ?n ?ofsx ?incx ~x ?ofsy ?incy y
see BLAS documentation!n
: default = greater n s.t. ofsx+(n-1)(abs incx) <= dim x
ofsx
: default = 1incx
: default = 1ofsy
: default = 1incy
: default = 1val scal : ?n:int ->
Lacaml_complex32.num_type ->
?ofsx:int -> ?incx:int -> Lacaml_complex32.vec -> unit
scal ?n alpha ?ofsx ?incx x
see BLAS documentation!n
: default = greater n s.t. ofsx+(n-1)(abs incx) <= dim x
ofsx
: default = 1incx
: default = 1val copy : ?n:int ->
?ofsy:int ->
?incy:int ->
?y:Lacaml_complex32.vec ->
?ofsx:int -> ?incx:int -> Lacaml_complex32.vec -> Lacaml_complex32.vec
copy ?n ?ofsy ?incy ?y ?ofsx ?incx x
see BLAS documentation!y
, which is overwritten.n
: default = greater n s.t. ofsx+(n-1)(abs incx) <= dim x
ofsy
: default = 1incy
: default = 1y
: default = new vector with ofsy+(n-1)(abs incy)
rowsofsx
: default = 1incx
: default = 1val nrm2 : ?n:int -> ?ofsx:int -> ?incx:int -> Lacaml_complex32.vec -> float
nrm2 ?n ?ofsx ?incx x
see BLAS documentation!n
: default = greater n s.t. ofsx+(n-1)(abs incx) <= dim x
ofsx
: default = 1incx
: default = 1val axpy : ?alpha:Lacaml_complex32.num_type ->
?n:int ->
?ofsx:int ->
?incx:int ->
Lacaml_complex32.vec ->
?ofsy:int -> ?incy:int -> Lacaml_complex32.vec -> unit
axpy ?alpha ?n ?ofsx ?incx x ?ofsy ?incy y
see BLAS documentation!alpha
: default = { re = 1.; im = 0. }
n
: default = greater n s.t. ofsx+(n-1)(abs incx) <= dim x
ofsx
: default = 1incx
: default = 1ofsy
: default = 1incy
: default = 1val iamax : ?n:int -> ?ofsx:int -> ?incx:int -> Lacaml_complex32.vec -> int
iamax ?n ?ofsx ?incx x
see BLAS documentation!n
: default = greater n s.t. ofsx+(n-1)(abs incx) <= dim x
ofsx
: default = 1incx
: default = 1val amax : ?n:int ->
?ofsx:int -> ?incx:int -> Lacaml_complex32.vec -> Lacaml_complex32.num_type
amax ?n ?ofsx ?incx x
x
.n
: default = greater n s.t. ofsx+(n-1)(abs incx) <= dim x
ofsx
: default = 1incx
: default = 1val gemv : ?m:int ->
?n:int ->
?beta:Lacaml_complex32.num_type ->
?ofsy:int ->
?incy:int ->
?y:Lacaml_complex32.vec ->
?trans:Lacaml_complex32.trans3 ->
?alpha:Lacaml_complex32.num_type ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?ofsx:int -> ?incx:int -> Lacaml_complex32.vec -> Lacaml_complex32.vec
gemv ?m ?n ?beta ?ofsy ?incy ?y ?trans ?alpha ?ar ?ac a ?ofsx ?incx x
see BLAS documentation! BEWARE that the 1988 BLAS-2 specification
mandates that this function has no effect when n=0
while the
mathematically expected behabior is y ← beta * y
.y
, which is overwritten.m
: default = number of available rows in matrix a
n
: default = available columns in matrix a
beta
: default = { re = 0.; im = 0. }
ofsy
: default = 1incy
: default = 1y
: default = vector with minimal required length (see BLAS)trans
: default = `Nalpha
: default = { re = 1.; im = 0. }
ar
: default = 1ac
: default = 1ofsx
: default = 1incx
: default = 1val gbmv : ?m:int ->
?n:int ->
?beta:Lacaml_complex32.num_type ->
?ofsy:int ->
?incy:int ->
?y:Lacaml_complex32.vec ->
?trans:Lacaml_complex32.trans3 ->
?alpha:Lacaml_complex32.num_type ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
int ->
int -> ?ofsx:int -> ?incx:int -> Lacaml_complex32.vec -> Lacaml_complex32.vec
gbmv
?m ?n ?beta ?ofsy ?incy ?y ?trans ?alpha ?ar ?ac a kl ku ?ofsx ?incx x
see BLAS documentation!y
, which is overwritten.m
: default = same as n
(i.e., a
is a square matrix)n
: default = available number of columns in matrix a
beta
: default = { re = 0.; im = 0. }
ofsy
: default = 1incy
: default = 1y
: default = vector with minimal required length (see BLAS)trans
: default = `Nalpha
: default = { re = 1.; im = 0. }
ar
: default = 1ac
: default = 1ofsx
: default = 1incx
: default = 1val symv : ?n:int ->
?beta:Lacaml_complex32.num_type ->
?ofsy:int ->
?incy:int ->
?y:Lacaml_complex32.vec ->
?up:bool ->
?alpha:Lacaml_complex32.num_type ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?ofsx:int -> ?incx:int -> Lacaml_complex32.vec -> Lacaml_complex32.vec
symv ?n ?beta ?ofsy ?incy ?y ?up ?alpha ?ar ?ac a ?ofsx ?incx x
see BLAS documentation!y
, which is overwritten.n
: default = dimension of symmetric matrix a
beta
: default = { re = 0.; im = 0. }
ofsy
: default = 1incy
: default = 1y
: default = vector with minimal required length (see BLAS)up
: default = true (upper triangular portion of a
is accessed)alpha
: default = { re = 1.; im = 0. }
ar
: default = 1ac
: default = 1ofsx
: default = 1incx
: default = 1val trmv : ?n:int ->
?trans:Lacaml_complex32.trans3 ->
?diag:Lacaml_common.diag ->
?up:bool ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?ofsx:int -> ?incx:int -> Lacaml_complex32.vec -> unit
trmv ?n ?trans ?diag ?up ?ar ?ac a ?ofsx ?incx x
see BLAS documentation!n
: default = dimension of triangular matrix a
trans
: default = `Ndiag
: default = false (not a unit triangular matrix)up
: default = true (upper triangular portion of a
is accessed)ar
: default = 1ac
: default = 1ofsx
: default = 1incx
: default = 1val trsv : ?n:int ->
?trans:Lacaml_complex32.trans3 ->
?diag:Lacaml_common.diag ->
?up:bool ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?ofsx:int -> ?incx:int -> Lacaml_complex32.vec -> unit
trsv ?n ?trans ?diag ?up ?ar ?ac a ?ofsx ?incx x
see BLAS documentation!n
: default = dimension of triangular matrix a
trans
: default = `Ndiag
: default = false (not a unit triangular matrix)up
: default = true (upper triangular portion of a
is accessed)ar
: default = 1ac
: default = 1ofsx
: default = 1incx
: default = 1val tpmv : ?n:int ->
?trans:Lacaml_complex32.trans3 ->
?diag:Lacaml_common.diag ->
?up:bool ->
?ofsap:int ->
Lacaml_complex32.vec ->
?ofsx:int -> ?incx:int -> Lacaml_complex32.vec -> unit
tpmv ?n ?trans ?diag ?up ?ofsap ap ?ofsx ?incx x
see BLAS documentation!n
: default = dimension of packed triangular matrix ap
trans
: default = `Ndiag
: default = false (not a unit triangular matrix)up
: default = true (upper triangular portion of ap
is accessed)ofsap
: default = 1ofsx
: default = 1incx
: default = 1val tpsv : ?n:int ->
?trans:Lacaml_complex32.trans3 ->
?diag:Lacaml_common.diag ->
?up:bool ->
?ofsap:int ->
Lacaml_complex32.vec ->
?ofsx:int -> ?incx:int -> Lacaml_complex32.vec -> unit
tpsv ?n ?trans ?diag ?up ?ofsap ap ?ofsx ?incx x
see BLAS documentation!n
: default = dimension of packed triangular matrix ap
trans
: default = `Ndiag
: default = false (not a unit triangular matrix)up
: default = true (upper triangular portion of ap
is accessed)ofsap
: default = 1ofsx
: default = 1incx
: default = 1val gemm : ?m:int ->
?n:int ->
?k:int ->
?beta:Lacaml_complex32.num_type ->
?cr:int ->
?cc:int ->
?c:Lacaml_complex32.mat ->
?transa:Lacaml_complex32.trans3 ->
?alpha:Lacaml_complex32.num_type ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?transb:Lacaml_complex32.trans3 ->
?br:int -> ?bc:int -> Lacaml_complex32.mat -> Lacaml_complex32.mat
gemm ?m ?n ?k ?beta ?cr ?cc ?c ?transa ?alpha ?ar ?ac a ?transb ?br ?bc b
see BLAS documentation!c
, which is overwritten.m
: default = number of rows of a
(or tr a
) and c
n
: default = number of columns of b
(or tr b
) and c
k
: default = number of columns of a
(or tr a
) and
number of rows of b
(or tr b
)beta
: default = { re = 0.; im = 0. }
cr
: default = 1cc
: default = 1c
: default = matrix with minimal required dimensiontransa
: default = `Nalpha
: default = { re = 1.; im = 0. }
ar
: default = 1ac
: default = 1transb
: default = `Nbr
: default = 1bc
: default = 1val symm : ?m:int ->
?n:int ->
?side:Lacaml_common.side ->
?up:bool ->
?beta:Lacaml_complex32.num_type ->
?cr:int ->
?cc:int ->
?c:Lacaml_complex32.mat ->
?alpha:Lacaml_complex32.num_type ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?br:int -> ?bc:int -> Lacaml_complex32.mat -> Lacaml_complex32.mat
symm ?m ?n ?side ?up ?beta ?cr ?cc ?c ?alpha ?ar ?ac a ?br ?bc b
see BLAS documentation!c
, which is overwritten.m
: default = number of rows of c
n
: default = number of columns of c
side
: default = `L (left - multiplication is a
b
)up
: default = true (upper triangular portion of a
is accessed)beta
: default = { re = 0.; im = 0. }
cr
: default = 1cc
: default = 1c
: default = matrix with minimal required dimensionalpha
: default = { re = 1.; im = 0. }
ar
: default = 1ac
: default = 1br
: default = 1bc
: default = 1val trmm : ?m:int ->
?n:int ->
?side:Lacaml_common.side ->
?up:bool ->
?transa:Lacaml_complex32.trans3 ->
?diag:Lacaml_common.diag ->
?alpha:Lacaml_complex32.num_type ->
?ar:int ->
?ac:int ->
a:Lacaml_complex32.mat -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
trmm ?m ?n ?side ?up ?transa ?diag ?alpha ?ar ?ac ~a ?br ?bc b
see BLAS documentation!m
: default = number of rows of b
n
: default = number of columns of b
side
: default = `L (left - multiplication is a
b
)up
: default = true (upper triangular portion of a
is accessed)transa
: default = `Ndiag
: default = `N (non-unit)alpha
: default = { re = 1.; im = 0. }
ar
: default = 1ac
: default = 1br
: default = 1bc
: default = 1val trsm : ?m:int ->
?n:int ->
?side:Lacaml_common.side ->
?up:bool ->
?transa:Lacaml_complex32.trans3 ->
?diag:Lacaml_common.diag ->
?alpha:Lacaml_complex32.num_type ->
?ar:int ->
?ac:int ->
a:Lacaml_complex32.mat -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
trsm ?m ?n ?side ?up ?transa ?diag ?alpha ?ar ?ac ~a ?br ?bc b
see BLAS documentation!b
, which is overwritten.m
: default = number of rows of b
n
: default = number of columns of b
side
: default = `L (left - multiplication is a
b
)up
: default = true (upper triangular portion of a
is accessed)transa
: default = `Ndiag
: default = `N (non-unit)alpha
: default = { re = 1.; im = 0. }
ar
: default = 1ac
: default = 1br
: default = 1bc
: default = 1val syrk : ?n:int ->
?k:int ->
?up:bool ->
?beta:Lacaml_complex32.num_type ->
?cr:int ->
?cc:int ->
?c:Lacaml_complex32.mat ->
?trans:Lacaml_common.trans2 ->
?alpha:Lacaml_complex32.num_type ->
?ar:int -> ?ac:int -> Lacaml_complex32.mat -> Lacaml_complex32.mat
syrk ?n ?k ?up ?beta ?cr ?cc ?c ?trans ?alpha ?ar ?ac a
see BLAS documentation!c
, which is overwritten.n
: default = number of rows of a
(or a
'), c
k
: default = number of columns of a
(or a
')up
: default = true (upper triangular portion of c
is accessed)beta
: default = { re = 0.; im = 0. }
cr
: default = 1cc
: default = 1c
: default = matrix with minimal required dimensiontrans
: default = `Nalpha
: default = { re = 1.; im = 0. }
ar
: default = 1ac
: default = 1val syr2k : ?n:int ->
?k:int ->
?up:bool ->
?beta:Lacaml_complex32.num_type ->
?cr:int ->
?cc:int ->
?c:Lacaml_complex32.mat ->
?trans:Lacaml_common.trans2 ->
?alpha:Lacaml_complex32.num_type ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?br:int -> ?bc:int -> Lacaml_complex32.mat -> Lacaml_complex32.mat
syr2k ?n ?k ?up ?beta ?cr ?cc ?c ?trans ?alpha ?ar ?ac a ?br ?bc b
see BLAS documentation!c
, which is overwritten.n
: default = number of rows of a
(or a
'), c
k
: default = number of columns of a
(or a
')up
: default = true (upper triangular portion of c
is accessed)beta
: default = { re = 0.; im = 0. }
cr
: default = 1cc
: default = 1c
: default = matrix with minimal required dimensiontrans
: default = `Nalpha
: default = { re = 1.; im = 0. }
ar
: default = 1ac
: default = 1br
: default = 1bc
: default = 1val lacpy : ?uplo:[ `L | `U ] ->
?m:int ->
?n:int ->
?br:int ->
?bc:int ->
?b:Lacaml_complex32.mat ->
?ar:int -> ?ac:int -> Lacaml_complex32.mat -> Lacaml_complex32.mat
lacpy ?uplo ?m ?n ?br ?bc ?b ?ar ?ac a
copy a (triangular)
(sub-)matrix a
(to an optional (sub-)matrix b
).uplo
: default = whole matrixval lassq : ?n:int ->
?scale:float ->
?sumsq:float ->
?ofsx:int -> ?incx:int -> Lacaml_complex32.vec -> float * float
lassq ?n ?ofsx ?incx ?scale ?sumsq
(scl, ssq)
, where
scl
is a scaling factor and ssq
the sum of squares of vector
x
starting at ofs
and using increment incx
and initial
scale
and sumsq
. The following equality holds:
scl**2. *. ssq = x.{1}**2. +. ... +. x.{n}**2. +. scale**2. *. sumsq
.
See LAPACK-documentation for details!n
: default = greater n s.t. ofsx+(n-1)(abs incx) <= dim x
scale
: default = 0.sumsq
: default = 1.ofsx
: default = 1incx
: default = 1val larnv : ?idist:[ `Normal | `Uniform0 | `Uniform1 ] ->
?iseed:Lacaml_common.int32_vec ->
?n:int ->
?ofsx:int -> ?x:Lacaml_complex32.vec -> unit -> Lacaml_complex32.vec
larnv ?idist ?iseed ?n ?ofsx ?x ()
idist
, random seed iseed
, vector offset
ofsx
and optional vector x
.idist
: default = `Normal
iseed
: default = integer vector of size 4 with all ones.n
: default = dim x - ofsx + 1
if x
is provided, 1
otherwise.ofsx
: default = 1
x
: default = vector of length ofsx - 1 + n
if n
is provided.val lange_min_lwork : int -> Lacaml_common.norm4 -> int
lange_min_lwork m norm
lange
-function.val lange : ?m:int ->
?n:int ->
?norm:Lacaml_common.norm4 ->
?work:Lacaml_complex32.rvec ->
?ar:int -> ?ac:int -> Lacaml_complex32.mat -> float
lange ?m ?n ?norm ?work ?ar ?ac a
norm = `O
), or the Frobenius norm (norm = `F
), or the infinity
norm (norm = `I
), or the element of largest absolute value
(norm = `M
) of a real matrix a
.m
: default = number of rows of matrix a
n
: default = number of columns of matrix a
norm
: default = `O
work
: default = allocated work space for norm `I
ar
: default = 1ac
: default = 1val lauum : ?up:bool -> ?n:int -> ?ar:int -> ?ac:int -> Lacaml_complex32.mat -> unit
lauum ?up ?n ?ar ?ac a
computes the product U * U**T or L**T * L,
where the triangular factor U or L is stored in the upper or lower
triangular part of the array a
. The upper or lower part of a
is overwritten.up
: default = true
n
: default = minimum of available number of rows/columns in matrix a
ar
: default = 1ac
: default = 1val getrf : ?m:int ->
?n:int ->
?ipiv:Lacaml_common.int32_vec ->
?ar:int -> ?ac:int -> Lacaml_complex32.mat -> Lacaml_common.int32_vec
getrf ?m ?n ?ipiv ?ar ?ac a
computes an LU factorization of a
general m
-by-n
matrix a
using partial pivoting with row
interchanges. See LAPACK documentation.Failure
if the matrix is singular.ipiv
, the pivot indices.m
: default = number of rows in matrix a
n
: default = number of columns in matrix a
ipiv
: = vec of length min(m, n)
ar
: default = 1ac
: default = 1val getrs : ?n:int ->
?ipiv:Lacaml_common.int32_vec ->
?trans:Lacaml_complex32.trans3 ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
getrs ?n ?ipiv ?trans ?ar ?ac a ?nrhs ?br ?bc b
solves a system
of linear equations a
* X = b
or a
' * X = b
with a general
n
-by-n
matrix a
using the LU factorization computed by
Lacaml_C.getrf
.
Note that matrix a
will be passed to Lacaml_C.getrf
if ipiv
was not
provided.Failure
if the matrix is singular.n
: default = number of columns in matrix a
ipiv
: default = result from getrf
applied to a
trans
: default = `Nar
: default = 1ac
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1val getri_min_lwork : int -> int
getri_min_lwork n
Lacaml_C.getri
-function if the matrix has n
columns.val getri_opt_lwork : ?n:int -> ?ar:int -> ?ac:int -> Lacaml_complex32.mat -> int
getri_opt_lwork ?n ?ar ?ac a
Lacaml_C.getri
-function.n
: default = number of columns of matrix a
ar
: default = 1ac
: default = 1val getri : ?n:int ->
?ipiv:Lacaml_common.int32_vec ->
?work:Lacaml_complex32.vec ->
?ar:int -> ?ac:int -> Lacaml_complex32.mat -> unit
getri ?n ?ipiv ?work ?ar ?ac a
computes the inverse of a matrix
using the LU factorization computed by Lacaml_C.getrf
. Note that matrix
a
will be passed to Lacaml_C.getrf
if ipiv
was not provided.Failure
if the matrix is singular.n
: default = number of columns in matrix a
ipiv
: default = vec of length m
from getriwork
: default = vec of optimum lengthar
: default = 1ac
: default = 1val sytrf_min_lwork : unit -> int
sytrf_min_lwork ()
Lacaml_C.sytrf
-function.val sytrf_opt_lwork : ?n:int -> ?up:bool -> ?ar:int -> ?ac:int -> Lacaml_complex32.mat -> int
sytrf_opt_lwork ?n ?up ?ar ?ac a
Lacaml_C.sytrf
-function.n
: default = number of columns of matrix a
up
: default = true (store upper triangle in a
)ar
: default = 1ac
: default = 1val sytrf : ?n:int ->
?up:bool ->
?ipiv:Lacaml_common.int32_vec ->
?work:Lacaml_complex32.vec ->
?ar:int -> ?ac:int -> Lacaml_complex32.mat -> Lacaml_common.int32_vec
sytrf ?n ?up ?ipiv ?work ?ar ?ac a
computes the factorization of
the real symmetric matrix a
using the Bunch-Kaufman diagonal
pivoting method.Failure
if D in a
= U*D*U' or L*D*L' is singular.n
: default = number of columns in matrix a
up
: default = true (store upper triangle in a
)ipiv
: = vec of length nwork
: default = vec of optimum lengthar
: default = 1ac
: default = 1val sytrs : ?n:int ->
?up:bool ->
?ipiv:Lacaml_common.int32_vec ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
sytrs ?n ?up ?ipiv ?ar ?ac a ?nrhs ?br ?bc b
solves a system of
linear equations a
*X = b
with a real symmetric matrix a
using the factorization a
= U*D*U**T or a
= L*D*L**T computed
by Lacaml_C.sytrf
. Note that matrix a
will be passed to Lacaml_C.sytrf
if
ipiv
was not provided.Failure
if the matrix is singular.n
: default = number of columns in matrix a
up
: default = true (store upper triangle in a
)ipiv
: default = vec of length n
ar
: default = 1ac
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1val sytri_min_lwork : int -> int
sytri_min_lwork n
Lacaml_C.sytri
-function if the matrix has n
columns.val sytri : ?n:int ->
?up:bool ->
?ipiv:Lacaml_common.int32_vec ->
?work:Lacaml_complex32.vec ->
?ar:int -> ?ac:int -> Lacaml_complex32.mat -> unit
sytri ?n ?up ?ipiv ?work ?ar ?ac a
computes the inverse of the
real symmetric indefinite matrix a
using the factorization a
=
U*D*U**T or a
= L*D*L**T computed by Lacaml_C.sytrf
. Note that matrix
a
will be passed to Lacaml_C.sytrf
if ipiv
was not provided.Failure
if the matrix is singular.n
: default = number of columns in matrix a
up
: default = true (store upper triangle in a
)work
: default = vec of optimum lengthar
: default = 1ac
: default = 1val potrf : ?n:int ->
?up:bool ->
?ar:int ->
?ac:int -> ?jitter:Lacaml_complex32.num_type -> Lacaml_complex32.mat -> unit
potrf ?n ?up ?ar ?ac ?jitter a
factorizes symmetric positive
definite matrix a
(or the designated submatrix) using Cholesky
factorization.
Due to rounding errors ill-conditioned matrices may actually appear
as if they were not positive definite, thus leading to an exception.
One remedy for this problem is to add a small jitter
to the
diagonal of the matrix, which will usually allow Cholesky to complete
successfully (though at a small bias). For extremely ill-conditioned
matrices it is recommended to use (symmetric) eigenvalue decomposition
instead of this function for a numerically more stable factorization.
Raises Failure
if the matrix is singular.
n
: default = number of columns in matrix a
up
: default = true (store upper triangle in a
)ar
: default = 1ac
: default = 1jitter
: default = nothingval potrs : ?n:int ->
?up:bool ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?nrhs:int ->
?br:int ->
?bc:int ->
?factorize:bool ->
?jitter:Lacaml_complex32.num_type -> Lacaml_complex32.mat -> unit
potrs ?n ?up ?ar ?ac a ?nrhs ?br ?bc ?factorize ?jitter b
solves
a system of linear equations a
*X = b
, where a
is symmetric
positive definite matrix, using the Cholesky factorization a
=
U**T*U or a
= L*L**T computed by Lacaml_C.potrf
.Failure
if the matrix is singular.n
: default = number of columns in matrix a
up
: default = truear
: default = 1ac
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1factorize
: default = true (calls Lacaml_C.potrf
implicitly)jitter
: default = nothingval potri : ?n:int ->
?up:bool ->
?ar:int ->
?ac:int ->
?factorize:bool ->
?jitter:Lacaml_complex32.num_type -> Lacaml_complex32.mat -> unit
potri ?n ?up ?ar ?ac ?factorize ?jitter a
computes the inverse
of the real symmetric positive definite matrix a
using the
Cholesky factorization a
= U**T*U or a
= L*L**T computed by
Lacaml_C.potrf
.Failure
if the matrix is singular.n
: default = number of columns in matrix a
up
: default = true (upper triangle stored in a
)ar
: default = 1ac
: default = 1factorize
: default = true (calls Lacaml_C.potrf
implicitly)jitter
: default = nothingval trtrs : ?n:int ->
?up:bool ->
?trans:Lacaml_complex32.trans3 ->
?diag:Lacaml_common.diag ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
trtrs ?n ?up ?trans ?diag ?ar ?ac a ?nrhs ?br ?bc b
solves a
triangular system of the form a
* X = b
or a
**T * X = n
,
where a
is a triangular matrix of order n
, and b
is an
n
-by-nrhs
matrix.Failure
if the matrix a
is singular.n
: default = number of columns in matrix a
up
: default = truetrans
: default = `Ndiag
: default = `Nar
: default = 1ac
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1val tbtrs : ?n:int ->
?kd:int ->
?up:bool ->
?trans:Lacaml_complex32.trans3 ->
?diag:Lacaml_common.diag ->
?abr:int ->
?abc:int ->
Lacaml_complex32.mat ->
?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
tbtrs ?n ?kd ?up ?trans ?diag ?abr ?abc ab ?nrhs ?br ?bc b
solves a triangular system of the form a
* X = b
or a
**T * X = b
,
where a
is a triangular band matrix of order n
, and b
is
an n
-by-nrhs
matrix.Failure
if the matrix a
is singular.n
: default = number of columns in matrix ab
kd
: default = number of rows in matrix ab
- 1up
: default = truetrans
: default = `Ndiag
: default = `Nabr
: default = 1abc
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1val trtri : ?n:int ->
?up:bool ->
?diag:Lacaml_common.diag ->
?ar:int -> ?ac:int -> Lacaml_complex32.mat -> unit
trtri ?n ?up ?diag ?ar ?ac a
computes the inverse of a real
upper or lower triangular matrix a
.Failure
if the matrix a
is singular.n
: default = number of columns in matrix a
up
: default = true (upper triangle stored in a
)diag
: default = `Nar
: default = 1ac
: default = 1val geqrf_opt_lwork : ?m:int -> ?n:int -> ?ar:int -> ?ac:int -> Lacaml_complex32.mat -> int
geqrf_opt_lwork ?m ?n ?ar ?ac a
Lacaml_C.geqrf
-function given matrix
a
and optionally its logical dimensions m
and n
.m
: default = number of rows in matrix a
n
: default = number of columns in matrix a
ar
: default = 1ac
: default = 1val geqrf_min_lwork : n:int -> int
geqrf_min_lwork ~n
Lacaml_C.geqrf
-function if the matrix has n
columns.val geqrf : ?m:int ->
?n:int ->
?work:Lacaml_complex32.vec ->
?tau:Lacaml_complex32.vec ->
?ar:int -> ?ac:int -> Lacaml_complex32.mat -> Lacaml_complex32.vec
geqrf ?m ?n ?work ?tau ?ar ?ac a
computes a QR factorization of
a real m
-by-n
matrix a
. See LAPACK documentation.tau
, the scalar factors of the elementary reflectors.m
: default = number of rows in matrix a
n
: default = number of columns in matrix a
work
: default = vec of optimum lengthtau
: default = vec of required lengthar
: default = 1ac
: default = 1val gesv : ?n:int ->
?ipiv:Lacaml_common.int32_vec ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
gesv ?n ?ipiv ?ar ?ac a ?nrhs ?br ?bc b
computes the solution to
a real system of linear equations a
* X = b
, where a
is an
n
-by-n
matrix and X and b
are n
-by-nrhs
matrices. The
LU decomposition with partial pivoting and row interchanges is
used to factor a
as a
= P * L * U, where P is a permutation
matrix, L is unit lower triangular, and U is upper triangular.
The factored form of a
is then used to solve the system of
equations a
* X = b
. On exit, b
contains the solution matrix X.Failure
if the matrix a
is singular.n
: default = available number of columns in matrix a
ipiv
: default = vec of length n
ar
: default = 1ac
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1val gbsv : ?n:int ->
?ipiv:Lacaml_common.int32_vec ->
?abr:int ->
?abc:int ->
Lacaml_complex32.mat ->
int -> int -> ?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
gbsv ?n ?ipiv ?abr ?abc ab kl ku ?nrhs ?br ?bc b
computes the
solution to a real system of linear equations a
* X = b
, where
a
is a band matrix of order n
with kl
subdiagonals and ku
superdiagonals, and X and b
are n
-by-nrhs
matrices. The LU
decomposition with partial pivoting and row interchanges is used
to factor a
as a
= L * U, where L is a product of permutation and
unit lower triangular matrices with kl
subdiagonals, and U is
upper triangular with kl+ku
superdiagonals. The factored form of
a
is then used to solve the system of equations a
* X = b
.Failure
if the matrix a
is singular.n
: default = available number of columns in matrix ab
ipiv
: default = vec of length n
abr
: default = 1abc
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1val gtsv : ?n:int ->
?ofsdl:int ->
Lacaml_complex32.vec ->
?ofsd:int ->
Lacaml_complex32.vec ->
?ofsdu:int ->
Lacaml_complex32.vec ->
?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
gtsv ?n ?ofsdl dl ?ofsd d ?ofsdu du ?nrhs ?br ?bc b
solves the
equation a
* X = b
where a
is an n
-by-n
tridiagonal
matrix, by Gaussian elimination with partial pivoting. Note that
the equation A
'*X = b
may be solved by interchanging the order
of the arguments du
and dl
.Failure
if the matrix is singular.n
: default = available length of vector d
ofsdl
: default = 1ofsd
: default = 1ofsdu
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1val posv : ?n:int ->
?up:bool ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
posv ?n ?up ?ar ?ac a ?nrhs ?br ?bc b
computes the solution to a
real system of linear equations a
* X = b
, where a
is an
n
-by-n
symmetric positive definite matrix and X and b
are
n
-by-nrhs
matrices. The Cholesky decomposition is used to
factor a
as
a
= U**T * U, if up = true
, or
a
= L * L**T, if up = false
,
where U is an upper triangular matrix and L is a lower triangular
matrix. The factored form of a
is then used to solve the system
of equations a
* X = b
.Failure
if the matrix is singular.n
: default = available number of columns in matrix a
up
: default = true i.e., upper triangle of a
is storedar
: default = 1ac
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1val ppsv : ?n:int ->
?up:bool ->
?ofsap:int ->
Lacaml_complex32.vec ->
?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
ppsv ?n ?up ?ofsap ap ?nrhs ?br ?bc b
computes the solution to
the real system of linear equations a
* X = b
, where a
is an
n
-by-n
symmetric positive definite matrix stored in packed
format and X and b
are n
-by-nrhs
matrices. The Cholesky
decomposition is used to factor a
as
a
= U**T * U, if up = true
, or
a
= L * L**T, if up = false
,
where U is an upper triangular matrix and L is a lower triangular
matrix. The factored form of a
is then used to solve the system
of equations a
* X = b
.Failure
if the matrix is singular.n
: default = the greater n s.t. n(n+1)/2 <= Vec.dim ap
up
: default = true i.e., upper triangle of ap
is storedofsap
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1val pbsv : ?n:int ->
?up:bool ->
?kd:int ->
?abr:int ->
?abc:int ->
Lacaml_complex32.mat ->
?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
pbsv ?n ?up ?kd ?abr ?abc ab ?nrhs ?br ?bc b
computes the
solution to a real system of linear equations a
* X = b
, where
a
is an n
-by-n
symmetric positive definite band matrix and X
and b
are n
-by-nrhs
matrices. The Cholesky decomposition is
used to factor a
as
a
= U**T * U, if up = true
, or
a
= L * L**T, if up = false
,
where U is an upper triangular band matrix, and L is a lower
triangular band matrix, with the same number of superdiagonals or
subdiagonals as a
. The factored form of a
is then used to
solve the system of equations a
* X = b
.Failure
if the matrix is singular.n
: default = available number of columns in matrix ab
up
: default = true i.e., upper triangle of ab
is storedkd
: default = available number of rows in matrix ab
- 1abr
: default = 1abc
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1val ptsv : ?n:int ->
?ofsd:int ->
Lacaml_complex32.vec ->
?ofse:int ->
Lacaml_complex32.vec ->
?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
ptsv ?n ?ofsd d ?ofse e ?nrhs ?br ?bc b
computes the solution to
the real system of linear equations a
*X = b
, where a
is an
n
-by-n
symmetric positive definite tridiagonal matrix, and X
and b
are n
-by-nrhs
matrices. A is factored as a
=
L*D*L**T, and the factored form of a
is then used to solve the
system of equations.Failure
if the matrix is singular.n
: default = available length of vector d
ofsd
: default = 1ofse
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1val sysv_opt_lwork : ?n:int ->
?up:bool ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> int
sysv_opt_lwork ?n ?up ?ar ?ac a ?nrhs ?br ?bc b
sysv
-function given matrix
a
, optionally its logical dimension n
and given right hand side
matrix b
with an optional number nrhs
of vectors.n
: default = available number of columns in matrix a
up
: default = true i.e., upper triangle of a
is storedar
: default = 1ac
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1val sysv : ?n:int ->
?up:bool ->
?ipiv:Lacaml_common.int32_vec ->
?work:Lacaml_complex32.vec ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
sysv ?n ?up ?ipiv ?work ?ar ?ac a ?nrhs ?br ?bc b
computes the
solution to a real system of linear equations a
* X = b
, where
a
is an N-by-N symmetric matrix and X and b
are n
-by-nrhs
matrices. The diagonal pivoting method is used to factor a
as
a
= U * D * U**T, if up = true
, or
a
= L * D * L**T, if up = false
,
where U (or L) is a product of permutation and unit upper (lower)
triangular matrices, and D is symmetric and block diagonal with
1-by-1 and 2-by-2 diagonal blocks. The factored form of a
is
then used to solve the system of equations a
* X = b
.Failure
if the matrix is singular.n
: default = available number of columns in matrix a
up
: default = true i.e., upper triangle of a
is storedipiv
: default = vec of length n
work
: default = vec of optimum length (-> sysv_opt_lwork
)ar
: default = 1ac
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1val spsv : ?n:int ->
?up:bool ->
?ipiv:Lacaml_common.int32_vec ->
?ofsap:int ->
Lacaml_complex32.vec ->
?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
spsv ?n ?up ?ipiv ?ofsap ap ?nrhs ?br ?bc b
computes the
solution to the real system of linear equations a
* X = b
,
where a
is an n
-by-n
symmetric matrix stored in packed
format and X and b
are n
-by-nrhs
matrices. The diagonal
pivoting method is used to factor a
as
a
= U * D * U**T, if up = true
, or
a
= L * D * L**T, if up = false
,
where U (or L) is a product of permutation and unit upper (lower)
triangular matrices, D is symmetric and block diagonal with 1-by-1
and 2-by-2 diagonal blocks. The factored form of a
is then used
to solve the system of equations a
* X = b
.Failure
if the matrix is singular.n
: default = the greater n s.t. n(n+1)/2 <= Vec.dim ap
up
: default = true i.e., upper triangle of ap
is storedipiv
: default = vec of length n
ofsap
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1val gels_min_lwork : m:int -> n:int -> nrhs:int -> int
gels_min_lwork ~m ~n ~nrhs
gels
-function if the logical dimensions
of the matrix are m
rows and n
columns and if there are nrhs
right hand side vectors.val gels_opt_lwork : ?m:int ->
?n:int ->
?trans:Lacaml_common.trans2 ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> int
gels_opt_lwork ?m ?n ?trans ?ar ?ac a ?nrhs ?br ?bc b
gels
-function given
matrix a
, optionally its logical dimensions m
and n
and given
right hand side matrix b
with an optional number nrhs
of vectors.m
: default = available number of rows in matrix a
n
: default = available number of columns in matrix a
trans
: default = `Nar
: default = 1ac
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1val gels : ?m:int ->
?n:int ->
?work:Lacaml_complex32.vec ->
?trans:Lacaml_common.trans2 ->
?ar:int ->
?ac:int ->
Lacaml_complex32.mat ->
?nrhs:int -> ?br:int -> ?bc:int -> Lacaml_complex32.mat -> unit
gels ?m ?n ?work ?trans ?ar ?ac a ?nrhs ?br ?bc b
see
LAPACK documentation!m
: default = available number of rows in matrix a
n
: default = available number of columns of matrix a
work
: default = vec of optimum length (-> Lacaml_C.gels_opt_lwork
)trans
: default = `Nar
: default = 1ac
: default = 1nrhs
: default = available number of columns in matrix b
br
: default = 1bc
: default = 1