Matrix  $Rev: 3071 $ at $LastChangedDate: 2015-03-26 15:35:47 +0100 (Thu, 26 Mar 2015) $
chm_common.h
Go to the documentation of this file.
1 #ifndef CHM_COMMON_H
2 #define CHM_COMMON_H
3 
4 #include "SuiteSparse_config/SuiteSparse_config.h"
5 #include "CHOLMOD/Include/cholmod.h"
6 #include "Mutils.h"
7 // -> R_check_class() et al
8 
9 #ifdef Matrix_with_SPQR
10 # include "SPQR/Include/SuiteSparseQR_C.h"
11 #endif
12 
13 /*
14 typedef struct cholmod_common_struct *CHM_CM ;
15 typedef struct cholmod_dense_struct *CHM_DN ;
16 typedef struct cholmod_factor_struct *CHM_FR ;
17 typedef struct cholmod_sparse_struct *CHM_SP ;
18 typedef struct cholmod_triplet_struct *CHM_TR ;
19 */
20 typedef cholmod_common* CHM_CM;
21 typedef cholmod_dense* CHM_DN;
22 typedef const cholmod_dense* const_CHM_DN;
23 typedef cholmod_factor* CHM_FR;
24 typedef const cholmod_factor* const_CHM_FR;
25 typedef cholmod_sparse* CHM_SP;
26 typedef const cholmod_sparse* const_CHM_SP;
27 typedef cholmod_triplet* CHM_TR;
28 typedef const cholmod_triplet* const_CHM_TR;
29 
30 extern cholmod_common c; /* structure for int CHM routines */
31 extern cholmod_common cl; /* structure for SuiteSparse_long routines */
32 
33 /* NOTE: Versions of these are *EXPORTED* via ../inst/include/Matrix.h
34  * ---- and used e.g., in the lme4 package
35  */
36 CHM_SP as_cholmod_sparse (CHM_SP ans, SEXP x, Rboolean check_Udiag, Rboolean sort_in_place);
37 CHM_TR as_cholmod_triplet(CHM_TR ans, SEXP x, Rboolean check_Udiag);
38 CHM_DN as_cholmod_dense (CHM_DN ans, SEXP x);
39 CHM_DN as_cholmod_x_dense(CHM_DN ans, SEXP x);
40 CHM_DN numeric_as_chm_dense(CHM_DN ans, double *v, int nr, int nc);
41 CHM_FR as_cholmod_factor (CHM_FR ans, SEXP x);
42 
43 #define AS_CHM_DN(x) as_cholmod_dense ((CHM_DN)alloca(sizeof(cholmod_dense)), x )
44 #define AS_CHM_xDN(x) as_cholmod_x_dense ((CHM_DN)alloca(sizeof(cholmod_dense)), x )
45 #define AS_CHM_FR(x) as_cholmod_factor ((CHM_FR)alloca(sizeof(cholmod_factor)), x )
46 #define AS_CHM_SP(x) as_cholmod_sparse ((CHM_SP)alloca(sizeof(cholmod_sparse)), x, TRUE, FALSE)
47 #define AS_CHM_TR(x) as_cholmod_triplet((CHM_TR)alloca(sizeof(cholmod_triplet)),x, TRUE)
48 /* the non-diagU2N-checking versions : */
49 #define AS_CHM_SP__(x) as_cholmod_sparse ((CHM_SP)alloca(sizeof(cholmod_sparse)), x, FALSE, FALSE)
50 #define AS_CHM_TR__(x) as_cholmod_triplet((CHM_TR)alloca(sizeof(cholmod_triplet)), x, FALSE)
51 // optional diagU2N-checking
52 #define AS_CHM_SP2(x,chk) as_cholmod_sparse ((CHM_SP)alloca(sizeof(cholmod_sparse)), x, chk, FALSE)
53 
54 
55 #define N_AS_CHM_DN(x,nr,nc) M_numeric_as_chm_dense((CHM_DN)alloca(sizeof(cholmod_dense)), x , nr, nc )
56 
57 static R_INLINE Rboolean chm_factor_ok(CHM_FR f)
58 {
59  return (Rboolean) (f->minor >= f->n);
60 }
61 
62 Rboolean check_sorted_chm(CHM_SP A);
63 
64 int R_cholmod_start(CHM_CM Common);
65 int R_cholmod_l_start(CHM_CM Common);
66 void R_cholmod_error(int status, const char *file, int line, const char *message);
67 
69 SEXP chm_factor_to_SEXP(CHM_FR f, int dofree);
70 SEXP chm_sparse_to_SEXP(CHM_SP a, int dofree, int uploT, int Rkind,
71  const char *diag, SEXP dn);
72 SEXP chm_triplet_to_SEXP(CHM_TR a, int dofree, int uploT, int Rkind,
73  const char* diag, SEXP dn);
74 SEXP chm_dense_to_SEXP(CHM_DN a, int dofree, int Rkind, SEXP dn, Rboolean transp);
75 /* int uploST, char *diag, SEXP dn); */
76 SEXP chm_dense_to_matrix(CHM_DN a, int dofree, SEXP dn);
77 SEXP chm_dense_to_vector(CHM_DN a, int dofree);
78 
79 Rboolean chm_MOD_xtype(int to_xtype, cholmod_sparse *A, CHM_CM Common);
80 
81 void chm_diagN2U(CHM_SP chx, int uploT, Rboolean do_realloc);
82 void chm_transpose_dense(CHM_DN ans, CHM_DN x);
83 
84 SEXP CHMfactor_validate(SEXP obj);
85 SEXP CHMsimpl_validate(SEXP obj);
86 SEXP CHMsuper_validate(SEXP obj);
87 
88 SEXP CHM_set_common_env(SEXP rho);
89 void CHM_store_common();
90 void CHM_restore_common();
91 #endif
void CHM_restore_common()
Definition: chm_common.c:50
void chm_transpose_dense(CHM_DN ans, CHM_DN x)
Transpose a cholmod_dense matrix ("too trivial to be in CHOLMOD?")
Definition: chm_common.c:754
cholmod_sparse * CHM_SP
Definition: chm_common.h:25
const cholmod_triplet * const_CHM_TR
Definition: chm_common.h:28
SEXP chm_triplet_to_SEXP(CHM_TR a, int dofree, int uploT, int Rkind, const char *diag, SEXP dn)
Copy the contents of a to an appropriate TsparseMatrix object and, optionally, free a or free both a ...
Definition: chm_common.c:573
SEXP chm_sparse_to_SEXP(CHM_SP a, int dofree, int uploT, int Rkind, const char *diag, SEXP dn)
Copy the contents of a to an appropriate CsparseMatrix object and, optionally, free a or free both a ...
Definition: chm_common.c:335
cholmod_triplet * CHM_TR
Definition: chm_common.h:27
CHM_DN as_cholmod_dense(CHM_DN ans, SEXP x)
Populate ans with the pointers from x and modify its scalar elements accordingly. ...
Definition: chm_common.c:672
int R_cholmod_l_start(CHM_CM Common)
SEXP chm_factor_to_SEXP(CHM_FR f, int dofree)
Copy the contents of f to an appropriate CHMfactor object and, optionally, free f or free both f and ...
Definition: chm_common.c:1106
const cholmod_sparse * const_CHM_SP
Definition: chm_common.h:26
const cholmod_factor * const_CHM_FR
Definition: chm_common.h:24
CHM_DN as_cholmod_x_dense(CHM_DN ans, SEXP x)
Definition: chm_common.c:731
void R_cholmod_error(int status, const char *file, int line, const char *message)
Definition: chm_common.c:769
cholmod_common c
Definition: chm_common.c:15
static R_INLINE Rboolean chm_factor_ok(CHM_FR f)
Definition: chm_common.h:57
cholmod_factor * CHM_FR
Definition: chm_common.h:23
SEXP CHMfactor_validate(SEXP obj)
Definition: chm_common.c:1259
SEXP chm_dense_to_vector(CHM_DN a, int dofree)
Copy the contents of a to a numeric R object and, optionally, free a or free both a and its pointer t...
Definition: chm_common.c:996
void chm_diagN2U(CHM_SP chx, int uploT, Rboolean do_realloc)
Drop the (unit) diagonal entries from a cholmod_sparse matrix.
Definition: chm_common.c:1191
cholmod_dense * CHM_DN
Definition: chm_common.h:21
CHM_SP as_cholmod_sparse(CHM_SP ans, SEXP x, Rboolean check_Udiag, Rboolean sort_in_place)
Populate ans with the pointers from x and modify its scalar elements accordingly. ...
Definition: chm_common.c:245
cholmod_common * CHM_CM
Definition: chm_common.h:20
SEXP get_SuiteSparse_version()
Definition: chm_common.c:8
SEXP CHMsuper_validate(SEXP obj)
Definition: chm_common.c:1269
CHM_TR as_cholmod_triplet(CHM_TR ans, SEXP x, Rboolean check_Udiag)
Populate ans with the pointers from x and modify its scalar elements accordingly. ...
Definition: chm_common.c:480
Rboolean chm_MOD_xtype(int to_xtype, cholmod_sparse *A, CHM_CM Common)
Change the "type" of a cholmod_sparse matrix, i.e.
Definition: chm_common.c:438
SEXP CHMsimpl_validate(SEXP obj)
Definition: chm_common.c:1264
SEXP CHM_set_common_env(SEXP rho)
Definition: chm_common.c:74
CHM_FR as_cholmod_factor(CHM_FR ans, SEXP x)
Populate ans with the pointers from x and modify its scalar elements accordingly. ...
Definition: chm_common.c:1031
void CHM_store_common()
Definition: chm_common.c:25
int R_cholmod_start(CHM_CM Common)
Initialize the CHOLMOD library and replace the print and error functions by R-specific versions...
Definition: chm_common.c:810
SEXP chm_dense_to_SEXP(CHM_DN a, int dofree, int Rkind, SEXP dn, Rboolean transp)
Copy the contents of a to an appropriate denseMatrix object and, optionally, free a or free both a an...
Definition: chm_common.c:835
SEXP chm_dense_to_matrix(CHM_DN a, int dofree, SEXP dn)
Copy the contents of a to a matrix object and, optionally, free a or free both a and its pointer to i...
Definition: chm_common.c:942
const cholmod_dense * const_CHM_DN
Definition: chm_common.h:22
CHM_DN numeric_as_chm_dense(CHM_DN ans, double *v, int nr, int nc)
Definition: chm_common.c:1007
Rboolean check_sorted_chm(CHM_SP A)
Definition: chm_common.c:169
cholmod_common cl
Definition: chm_common.c:16