Matrix r5059
Loading...
Searching...
No Matches
cholmod-api.c File Reference
#include "Mdefines.h"
#include "cholmod-api.h"

Go to the source code of this file.

Macros

#define errorFree(...)
 
#define MAYBE_FREE
 
#define MAYBE_FREE
 
#define MAYBE_FREE
 
#define MAYBE_FREE
 

Functions

cholmod_factor * sexp_as_cholmod_factor (cholmod_factor *L, SEXP from)
 Coerce from sparseCholesky to (cholmod_factor *)
 
cholmod_sparse * sexp_as_cholmod_sparse (cholmod_sparse *A, SEXP from, Rboolean allocUnit, Rboolean sortInPlace)
 Coerce from [CR]sparseMatrix to (cholmod_sparse *)
 
cholmod_triplet * sexp_as_cholmod_triplet (cholmod_triplet *A, SEXP from, Rboolean allocUnit)
 Coerce from TsparseMatrix to (cholmod_triplet *)
 
cholmod_dense * sexp_as_cholmod_dense (cholmod_dense *A, SEXP from)
 Coerce from .geMatrix or vector to (cholmod_dense *)
 
cholmod_dense * numeric_as_cholmod_dense (cholmod_dense *A, double *data, int nrow, int ncol)
 Coerce from (double *) to (cholmod_dense *) with given dimensions.
 
SEXP cholmod_factor_as_sexp (cholmod_factor *L, int doFree)
 Coerce from (cholmod_factor *) to sparseCholesky.
 
SEXP cholmod_sparse_as_sexp (cholmod_sparse *A, int doFree, int ttype, int doLogic, const char *diagString, SEXP dimnames)
 Coerce from (cholmod_sparse *) to CsparseMatrix.
 
SEXP cholmod_triplet_as_sexp (cholmod_triplet *A, int doFree, int ttype, int doLogic, const char *diagString, SEXP dimnames)
 Coerce from (cholmod_triplet *) to TsparseMatrix.
 
SEXP cholmod_dense_as_sexp (cholmod_dense *A, int doFree)
 Coerce from (cholmod_dense *) to [dz]geMatrix.
 
double cholmod_factor_ldetA (cholmod_factor *L)
 Log determinant from Cholesky factorization.
 
cholmod_factor * cholmod_factor_update (cholmod_factor *L, cholmod_sparse *A, double beta)
 Update a Cholesky factorization.
 

Macro Definition Documentation

◆ errorFree

#define errorFree ( ...)
Value:
do { \
Rf_error(__VA_ARGS__); \
} while (0)
#define MAYBE_FREE

Referenced by cholmod_dense_as_sexp(), cholmod_factor_as_sexp(), cholmod_sparse_as_sexp(), and cholmod_triplet_as_sexp().

◆ MAYBE_FREE [1/4]

#define MAYBE_FREE
Value:
do { \
if (doFree != 0) { \
if (doFree < 0) \
R_Free(L); \
else if (L->itype == CHOLMOD_INT) \
cholmod_free_factor(&L, &c); \
else \
cholmod_l_free_factor(&L, &cl); \
} \
} while (0)
cholmod_common c
Definition cholmod-etc.c:5
cholmod_common cl
Definition cholmod-etc.c:6

Referenced by cholmod_dense_as_sexp(), cholmod_factor_as_sexp(), cholmod_sparse_as_sexp(), and cholmod_triplet_as_sexp().

◆ MAYBE_FREE [2/4]

#define MAYBE_FREE
Value:
do { \
if (doFree != 0) { \
if (doFree < 0) \
R_Free(A_); \
else if (A_->itype == CHOLMOD_INT) \
cholmod_free_sparse(&A_, &c); \
else \
cholmod_l_free_sparse(&A_, &cl); \
} \
} while (0)

◆ MAYBE_FREE [3/4]

#define MAYBE_FREE
Value:
do { \
if (doFree != 0) { \
if (doFree < 0) \
R_Free(A); \
else if (A->itype == CHOLMOD_INT) \
cholmod_free_triplet(&A, &c); \
else \
cholmod_l_free_triplet(&A, &cl); \
} \
} while (0)

◆ MAYBE_FREE [4/4]

#define MAYBE_FREE
Value:
do { \
if (doFree != 0) { \
if (doFree < 0) \
R_Free(A); \
else \
cholmod_free_dense(&A, &c); \
} \
} while (0)

Function Documentation

◆ cholmod_dense_as_sexp()

SEXP cholmod_dense_as_sexp ( cholmod_dense * A,
int doFree )

Coerce from (cholmod_dense *) to [dz]geMatrix.

Allocates an S4 object of class [dz]geMatrix and copies into the slots from members of a pointed-to cholmod_dense struct. The specific class of the result is determined by struct member xtype.

Parameters
Aa pointer to a cholmod_dense struct.
doFreea flag indicating if and how to free A before returning. (0) don't free, (>0) free with cholmod_free_dense, (<0) free with R_Free.
Returns
A [dz]geMatrix.

Definition at line 851 of file cholmod-api.c.

References _, errorFree, GET_SLOT, Matrix_DimSym, Matrix_xSym, MAYBE_FREE, newObject(), and SET_SLOT.

Referenced by R_init_Matrix().

◆ cholmod_factor_as_sexp()

SEXP cholmod_factor_as_sexp ( cholmod_factor * L,
int doFree )

Coerce from (cholmod_factor *) to sparseCholesky.

Allocates an S4 object inheriting from virtual class sparseCholesky and copies into the slots from members of a pointed-to cholmod_factor struct. The specific class of the result is determined by struct members xtype and is_super.

Parameters
La pointer to a cholmod_factor struct.
doFreea flag indicating if and how to free L before returning. (0) don't free, (>0) free with cholmod_free_factor, (<0) free with R_Free.
Returns
A sparseCholesky.

Definition at line 466 of file cholmod-api.c.

References _, errorFree, GET_SLOT, Matrix_colcountSym, Matrix_DimSym, Matrix_isllSym, Matrix_ismtSym, Matrix_iSym, Matrix_maxcsizeSym, Matrix_maxesizeSym, Matrix_minorSym, Matrix_nextSym, Matrix_nzSym, Matrix_orderingSym, Matrix_permSym, Matrix_piSym, Matrix_prevSym, Matrix_pSym, Matrix_pxSym, Matrix_sSym, Matrix_superSym, Matrix_xSym, MAYBE_FREE, newObject(), and SET_SLOT.

Referenced by R_init_Matrix().

◆ cholmod_factor_ldetA()

double cholmod_factor_ldetA ( cholmod_factor * L)

Log determinant from Cholesky factorization.

Computes log(det(A)) given the Cholesky factorization of A as P1 * A * P1' = L1 * D * L1' = L * L', L = L1 * sqrt(D). The result is computed as sum(log(diag(D))) or 2*sum(log(diag(L))), depending on members is_super and is_ll of the supplied struct. Note that CHOLMOD does not require diag(D) to be positive and that this routine does not check (FIXME).

Parameters
La pointer to a cholmod_factor struct.

Definition at line 911 of file cholmod-api.c.

References _.

Referenced by R_init_Matrix().

◆ cholmod_factor_update()

cholmod_factor * cholmod_factor_update ( cholmod_factor * L,
cholmod_sparse * A,
double beta )

Update a Cholesky factorization.

Updates in-place the Cholesky factorization of a symmetric matrix X+alpha*I with the Cholesky factorization of (1) A+beta*I, where A is a symmetric matrix sharing the nonzero pattern of X, or (2) A*A'+beta*I, where A is a general matrix sharing the nonzero pattern of Y, assuming that X = Y*Y'.

Parameters
La pointer to a cholmod_factor struct, to be modified in-place.
Aa pointer to a cholmod_sparse struct.
betaa multiplier, typically positive, to guarantee strict diagonal dominance.
Returns
L.

Definition at line 958 of file cholmod-api.c.

References _, and c.

Referenced by R_init_Matrix().

◆ cholmod_sparse_as_sexp()

SEXP cholmod_sparse_as_sexp ( cholmod_sparse * A,
int doFree,
int ttype,
int doLogic,
const char * diagString,
SEXP dimnames )

Coerce from (cholmod_sparse *) to CsparseMatrix.

Allocates an S4 object inheriting from virtual class CsparseMatrix and copies into the slots from members of a pointed-to cholmod_sparse struct. The specific class of the result is determined by struct members xtype and stype and by arguments ttype and doLogic.

Parameters
Aa pointer to a cholmod_sparse struct.
doFreea flag indicating if and how to free A before returning. (0) don't free, (>0) free with cholmod_free_sparse, (<0) free with R_Free.
ttypea flag indicating if the result should be a .tCMatrix. (0) not .tCMatrix, (>0) .tCMatrix with uplo="U", (<0) .tCMatrix with uplo="L". If ttype=0, then the result is a .gCMatrix or .sCMatrix depending on stype. (0) .gCMatrix, (>0) .sCMatrix with uplo="U", (<0) .sCMatrix with uplo="L".
doLogica flag indicating if the result should be a l.CMatrix if xtype=CHOLMOD_REAL.
diagStringa null-terminated string or NULL. The diag slot of a .tCMatrix result is "N" if and only if diagString is NULL or diagString[0] is 'N'.
dimnamesan R object specifying the Dimnames slot of the result, unused if not a list of length 2.
Returns
A CsparseMatrix.

Definition at line 625 of file cholmod-api.c.

References _, c, errorFree, GET_SLOT, Matrix_diagSym, Matrix_DimNamesSym, Matrix_DimSym, Matrix_iSym, Matrix_pSym, Matrix_uploSym, Matrix_xSym, MAYBE_FREE, newObject(), SET_SLOT, and TYPEOF.

Referenced by R_init_Matrix().

◆ cholmod_triplet_as_sexp()

SEXP cholmod_triplet_as_sexp ( cholmod_triplet * A,
int doFree,
int ttype,
int doLogic,
const char * diagString,
SEXP dimnames )

Coerce from (cholmod_triplet *) to TsparseMatrix.

Allocates an S4 object inheriting from virtual class TsparseMatrix and copies into the slots from members of a pointed-to cholmod_triplet struct. The specific class of the result is determined by struct members xtype and stype and by arguments ttype and doLogic.

Parameters
Aa pointer to a cholmod_triplet struct.
doFreea flag indicating if and how to free A before returning. (0) don't free, (>0) free with cholmod_free_triplet, (<0) free with R_Free.
ttypea flag indicating if the result should be a .tTMatrix. (0) not .tTMatrix, (>0) .tTMatrix with uplo="U", (<0) .tTMatrix with uplo="L". If ttype=0, then the result is a .gTMatrix or .sTMatrix depending on stype. (0) .gTMatrix, (>0) .sTMatrix with uplo="U", (<0) .sTMatrix with uplo="L".
doLogica flag indicating if the result should be an l.TMatrix if xtype=CHOLMOD_REAL.
diagStringa null-terminated string or NULL. The diag slot of a .tTMatrix result is "N" if and only if diagString is NULL or diagString[0] is 'N'.
dimnamesan R object specifying the Dimnames slot of the result, unused if not a list of length 2.
Returns
A TsparseMatrix.

Definition at line 740 of file cholmod-api.c.

References _, errorFree, GET_SLOT, Matrix_diagSym, Matrix_DimNamesSym, Matrix_DimSym, Matrix_iSym, Matrix_jSym, Matrix_uploSym, Matrix_xSym, MAYBE_FREE, newObject(), SET_SLOT, and TYPEOF.

Referenced by R_init_Matrix().

◆ numeric_as_cholmod_dense()

cholmod_dense * numeric_as_cholmod_dense ( cholmod_dense * A,
double * data,
int nrow,
int ncol )

Coerce from (double *) to (cholmod_dense *) with given dimensions.

An analogue of base::matrix(data, nrow, ncol), where typeof(data)=="double" and length(data)==nrow*ncol.

Parameters
Aa pointer to a cholmod_dense struct, to be modified in-place.
dataa pointer to an nrow*ncol*sizeof(double) block of memory.
nrowthe desired number of rows.
ncolthe desired number of columns.
Returns
A.

Definition at line 436 of file cholmod-api.c.

Referenced by R_init_Matrix().

◆ sexp_as_cholmod_dense()

cholmod_dense * sexp_as_cholmod_dense ( cholmod_dense * A,
SEXP from )

Coerce from .geMatrix or vector to (cholmod_dense *)

Sets the members of a pointed-to cholmod_dense struct, using "data" obtained from slots of a .geMatrix. The result should not be freed using cholmod_free_dense, as the resulting members point to memory controlled by R, not by CHOLMOD.

Parameters
Aa pointer to a cholmod_dense struct, to be modified in-place.
froman S4 object inheriting from class .geMatrix or a traditional vector of type "logical", "integer", "double", or "complex" (to be handled as a 1-column matrix if not a matrix).
Returns
A.

Definition at line 356 of file cholmod-api.c.

References ERROR_INVALID_TYPE, GET_SLOT, Matrix_class(), Matrix_DimSym, Matrix_xSym, and TYPEOF.

Referenced by R_init_Matrix().

◆ sexp_as_cholmod_factor()

cholmod_factor * sexp_as_cholmod_factor ( cholmod_factor * L,
SEXP from )

Coerce from sparseCholesky to (cholmod_factor *)

Sets the members of a pointed-to cholmod_factor struct, using "data" obtained from slots of a sparseCholesky. The result should not be freed using cholmod_free_factor, as the resulting members point to memory controlled by R, not by CHOLMOD.

Parameters
La pointer to a cholmod_factor struct, to be modified in-place.
froman S4 object inheriting from virtual class sparseCholesky.
Returns
L.

Definition at line 18 of file cholmod-api.c.

References GET_SLOT, Matrix_class(), Matrix_colcountSym, Matrix_DimSym, Matrix_isllSym, Matrix_ismtSym, Matrix_iSym, Matrix_maxcsizeSym, Matrix_maxesizeSym, Matrix_minorSym, Matrix_nextSym, Matrix_nzSym, Matrix_orderingSym, Matrix_permSym, Matrix_piSym, Matrix_prevSym, Matrix_pSym, Matrix_pxSym, Matrix_sSym, Matrix_superSym, and Matrix_xSym.

Referenced by R_init_Matrix().

◆ sexp_as_cholmod_sparse()

cholmod_sparse * sexp_as_cholmod_sparse ( cholmod_sparse * A,
SEXP from,
Rboolean allocUnit,
Rboolean sortInPlace )

Coerce from [CR]sparseMatrix to (cholmod_sparse *)

Sets the members of a pointed-to cholmod_sparse struct, using "data" obtained from slots of a [CR]sparseMatrix. The result should not be freed using cholmod_free_sparse, as the resulting members point to memory controlled by R, not by CHOLMOD.

Parameters
Aa pointer to a cholmod_sparse struct, to be modified in-place.
froman S4 object inheriting from virtual class [CR]sparseMatrix.
allocUnita boolean indicating if the unit diagonal of formally unit triangular [CR]sparseMatrix should be allocated.
sortInPlacea boolean indicating if unsorted [CR]sparseMatrix should be sorted in place to avoid an allocation.
Returns
A.

Definition at line 127 of file cholmod-api.c.

References _, allocUnit(), c, checkpi(), GET_SLOT, Matrix_class(), Matrix_diagSym, Matrix_DimSym, Matrix_iSym, Matrix_jSym, Matrix_pSym, Matrix_uploSym, Matrix_xSym, TYPEOF, and valid_sparse_compressed.

Referenced by R_init_Matrix().

◆ sexp_as_cholmod_triplet()

cholmod_triplet * sexp_as_cholmod_triplet ( cholmod_triplet * A,
SEXP from,
Rboolean allocUnit )

Coerce from TsparseMatrix to (cholmod_triplet *)

Sets the members of a pointed-to cholmod_triplet struct, using "data" obtained from slots of a TsparseMatrix. The result should not be freed using cholmod_free_sparse, as the resulting members point to memory controlled by R, not by CHOLMOD.

Parameters
Aa pointer to a cholmod_triplet struct, to be modified in-place.
froman S4 object inheriting from virtual class TsparseMatrix.
allocUnita boolean indicating if the unit diagonal of formally unit triangular TsparseMatrix should be allocated.
Returns
A.

Definition at line 247 of file cholmod-api.c.

References allocUnit(), GET_SLOT, Matrix_class(), Matrix_diagSym, Matrix_DimSym, Matrix_iSym, Matrix_pSym, Matrix_uploSym, Matrix_xSym, Matrix_zunit, TYPEOF, and valid_sparse_triplet.

Referenced by R_init_Matrix().