Module Sundials_Matrix.ArrayDense

module ArrayDense: sig .. end

General purpose dense matrix operations on arrays.


type t = Sundials.RealArray2.t 

A dense matrix accessible directly through a Bigarray.

Basic access

val make : int -> int -> float -> t

make m n x returns an m by n array dense matrix with elements set to x.

val create : int -> int -> t

create m n returns an uninitialized m by n array dense matrix.

val size : t -> int * int

m, n = size a returns the numbers of rows m and columns n of a.

val pp : Stdlib.Format.formatter -> t -> unit

Pretty-print an array dense matrix using the Format module.

val get : t -> int -> int -> float

get a i j returns the value at row i and column j of a.

val set : t -> int -> int -> float -> unit

set a i j v sets the value at row i and column j of a to v.

val update : t -> int -> int -> (float -> float) -> unit

update a i j f sets the value at row i and column j of a to f v.

val unwrap : t -> Sundials.RealArray2.data

Direct access to the underlying storage array, which is accessed column first (unlike in Sundials_Matrix.ArrayDense.get).

Operations

val ops : (t, Sundials.RealArray.t)
Sundials_Matrix.matrix_ops

Operations on array-based dense matrices.

val scale_add : float -> t -> t -> unit

scale_add c A B calculates $A = cA + B$.

val scale_addi : float -> t -> unit

scale_addi c A calculates $A = cA + I$.

val matvec : t ->
Sundials.RealArray.t -> Sundials.RealArray.t -> unit

The call matvec a x y computes the matrix-vector product $y = Ax$.

val set_to_zero : t -> unit

Fills the matrix with zeros.

val blit : src:t -> dst:t -> unit

blit ~src ~dst copies the contents of src into dst. Both must have the same size.

val space : t -> int * int

lrw, liw = space a returns the storage requirements of a as lrw realtype words and liw integer words.

Calculations

val add_identity : t -> unit

Increments a square matrix by the identity matrix.

val scale : float -> t -> unit

Multiplies each element by a constant.

val getrf : t -> Sundials.LintArray.t -> unit

getrf a p performs the LU factorization of the square matrix a with partial pivoting according to p. The values in a are overwritten with those of the calculated L and U matrices. The diagonal belongs to U. The diagonal of L is all 1s. Multiplying L by U gives a permutation of a, according to the values of p: p.{k} = j means that rows k and j were swapped (in order, where p.{0} swaps against the original matrix a).

val getrs : t ->
Sundials.LintArray.t -> Sundials.RealArray.t -> unit

getrs a p b finds the solution of ax = b using an LU factorization found by Sundials_Matrix.ArrayDense.getrf. Both p and b must have the same number of rows as a.

val getrs' : t ->
Sundials.LintArray.t -> Sundials.RealArray.t -> int -> unit

Like Sundials_Matrix.ArrayDense.getrs but stores b starting at a given offset.

val potrf : t -> unit

Performs Cholesky factorization of a real symmetric positive matrix.

val potrs : t -> Sundials.RealArray.t -> unit

potrs a b finds the solution of ax = b using the Cholesky factorization found by Sundials_Matrix.ArrayDense.potrf. a must be an n by n matrix and b must be of length n.

val geqrf : t ->
Sundials.RealArray.t -> Sundials.RealArray.t -> unit

geqrf a beta work performs the QR factorization of a. a must be an m by n matrix, where m >= n. The beta vector must have length n. The work vector must have length m.

val ormqr : a:t ->
beta:Sundials.RealArray.t ->
v:Sundials.RealArray.t ->
w:Sundials.RealArray.t -> work:Sundials.RealArray.t -> unit

ormqr q beta v w work computes the product w = qv . Q is an m by n matrix calculated using Sundials_Matrix.ArrayDense.geqrf with m >= n, beta has length n, v has length n, w has length m, and work has length m.

beta : vector passed to Sundials_Matrix.ArrayDense.geqrf
v : vector multiplier
w : result vector
work : temporary vector used in the calculation