Matrix  $Rev: 3071 $ at $LastChangedDate: 2015-03-26 15:35:47 +0100 (Thu, 26 Mar 2015) $
chm_common.c
Go to the documentation of this file.
1 
3 #include "chm_common.h"
4 // -> Mutils.h
5 
6 Rboolean isValid_Csparse(SEXP x); /* -> Csparse.c */
7 
9  SEXP ans = allocVector(INTSXP, 3);
10  int* version = INTEGER(ans);
11  SuiteSparse_version(version);
12  return ans;
13 }
14 
15 cholmod_common c;
16 cholmod_common cl;
17 
24 
26  SEXP rho = chm_common_env;
27  defineVar(dboundSym, ScalarReal(c.dbound), rho);
28  defineVar(grow0Sym, ScalarReal(c.grow0), rho);
29  defineVar(grow1Sym, ScalarReal(c.grow1), rho);
30  defineVar(grow2Sym, ScalarInteger(c.grow2), rho);
31  defineVar(maxrankSym, ScalarInteger(c.maxrank), rho);
32  defineVar(supernodal_switchSym,
33  ScalarReal(c.supernodal_switch), rho);
34  defineVar(supernodalSym, ScalarInteger(c.supernodal), rho);
35  defineVar(final_asisSym, ScalarLogical(c.final_asis), rho);
36  defineVar(final_superSym, ScalarLogical(c.final_super), rho);
37  defineVar(final_llSym, ScalarLogical(c.final_ll), rho);
38  defineVar(final_packSym, ScalarLogical(c.final_pack), rho);
39  defineVar(final_monotonicSym, ScalarLogical(c.final_monotonic), rho);
40  defineVar(final_resymbolSym, ScalarLogical(c.final_resymbol), rho);
41  defineVar(prefer_zomplexSym, ScalarLogical(c.prefer_zomplex), rho);
42  defineVar(prefer_upperSym, ScalarLogical(c.prefer_upper), rho);
44  ScalarLogical(c.quick_return_if_not_posdef), rho);
45  defineVar(nmethodsSym, ScalarInteger(c.nmethods), rho);
46  defineVar(m0_ordSym, ScalarInteger(c.method[0].ordering), rho);
47  defineVar(postorderSym, ScalarLogical(c.postorder), rho);
48 }
49 
51  SEXP rho = chm_common_env;
52  c.dbound = asReal(findVarInFrame(rho, dboundSym));
53  c.grow0 = asReal(findVarInFrame(rho, grow0Sym));
54  c.grow1 = asReal(findVarInFrame(rho, grow1Sym));
55  c.grow2 = asInteger(findVarInFrame(rho, grow2Sym));
56  c.maxrank = asInteger(findVarInFrame(rho, maxrankSym));
57  c.supernodal_switch = asReal(findVarInFrame(rho, supernodal_switchSym));
58  c.supernodal = asLogical(findVarInFrame(rho, supernodalSym));
59  c.final_asis = asLogical(findVarInFrame(rho, final_asisSym));
60  c.final_super = asLogical(findVarInFrame(rho, final_superSym));
61  c.final_ll = asLogical(findVarInFrame(rho, final_llSym));
62  c.final_pack = asLogical(findVarInFrame(rho, final_packSym));
63  c.final_monotonic = asLogical(findVarInFrame(rho, final_monotonicSym));
64  c.final_resymbol = asLogical(findVarInFrame(rho, final_resymbolSym));
65  c.prefer_zomplex = asLogical(findVarInFrame(rho, prefer_zomplexSym));
66  c.prefer_upper = asLogical(findVarInFrame(rho, prefer_upperSym));
67  c.quick_return_if_not_posdef =
68  asLogical(findVarInFrame(rho, quick_return_if_not_posdefSym));
69  c.nmethods = asInteger(findVarInFrame(rho, nmethodsSym));
70  c.method[0].ordering = asInteger(findVarInFrame(rho, m0_ordSym));
71  c.postorder = asLogical(findVarInFrame(rho, postorderSym));
72 }
73 
74 SEXP CHM_set_common_env(SEXP rho) {
75  if (!isEnvironment(rho))
76  error(_("Argument rho must be an environment"));
77  chm_common_env = rho;
78  dboundSym = install("dbound");
79  grow0Sym = install("grow0");
80  grow1Sym = install("grow1");
81  grow2Sym = install("grow2");
82  maxrankSym = install("maxrank");
83  supernodal_switchSym = install("supernodal_switch");
84  supernodalSym = install("supernodal");
85  final_asisSym = install("final_asis");
86  final_superSym = install("final_super");
87  final_llSym = install("final_ll");
88  final_packSym = install("final_pack");
89  final_monotonicSym = install("final_monotonic");
90  final_resymbolSym = install("final_resymbol");
91  prefer_zomplexSym = install("final_zomplex");
92  prefer_upperSym = install("final_upper");
93  quick_return_if_not_posdefSym = install("quick_return_if_not_posdef");
94  nmethodsSym = install("nmethods");
95  m0_ordSym = install("m0.ord");
96  postorderSym = install("postorder");
98  return R_NilValue;
99 }
100 
114 static int stype(int ctype, SEXP x)
115 {
116  if ((ctype % 3) == 1) return (*uplo_P(x) == 'U') ? 1 : -1;
117  return 0;
118 }
119 
126 static int xtype(int ctype)
127 {
128  switch(ctype / 3) {
129  case 0: /* "d" */
130  case 1: /* "l" */
131  return CHOLMOD_REAL;
132  case 2: /* "n" */
133  return CHOLMOD_PATTERN;
134  case 3: /* "z" */
135  return CHOLMOD_COMPLEX;
136  }
137  return -1;
138 }
139 
140 /* coerce a vector to REAL and copy the result to freshly R_alloc'd memory */
141 static void *RallocedREAL(SEXP x)
142 {
143  SEXP rx = PROTECT(coerceVector(x, REALSXP));
144  int lx = LENGTH(rx);
145  /* We over-allocate the memory chunk so that it is never NULL. */
146  /* The CHOLMOD code checks for a NULL pointer even in the length-0 case. */
147  double *ans = Memcpy((double*)R_alloc(lx + 1, sizeof(double)),
148  REAL(rx), lx);
149  UNPROTECT(1);
150  return (void*)ans;
151 }
152 
153 
154 static void *xpt(int ctype, SEXP x)
155 {
156  switch(ctype / 3) {
157  case 0: /* "d" */
158  return (void *) REAL(GET_SLOT(x, Matrix_xSym));
159  case 1: /* "l" */
160  return RallocedREAL(GET_SLOT(x, Matrix_xSym));
161  case 2: /* "n" */
162  return (void *) NULL;
163  case 3: /* "z" */
164  return (void *) COMPLEX(GET_SLOT(x, Matrix_xSym));
165  }
166  return (void *) NULL; /* -Wall */
167 }
168 
170 {
171  int *Ai = (int*)(A->i), *Ap = (int*)(A->p);
172  int j, p;
173 
174  for (j = 0; j < A->ncol; j++) {
175  int p1 = Ap[j], p2 = Ap[j + 1] - 1;
176  for (p = p1; p < p2; p++)
177  if (Ai[p] >= Ai[p + 1])
178  return FALSE;
179  }
180  return TRUE;
181 }
182 
186 static void chm2Ralloc(CHM_SP dest, CHM_SP src)
187 {
188  int np1, nnz;
189 
190  /* copy all the characteristics of src to dest */
191  memcpy(dest, src, sizeof(cholmod_sparse));
192 
193  /* R_alloc the vector storage for dest and copy the contents from src */
194  np1 = src->ncol + 1;
195  nnz = (int) cholmod_nnz(src, &c);
196  dest->p = (void*) Memcpy((int*)R_alloc(np1, sizeof(int)),
197  (int*)(src->p), np1);
198  dest->i = (void*) Memcpy((int*)R_alloc(nnz, sizeof(int)),
199  (int*)(src->i), nnz);
200  if(src->xtype)
201  dest->x = (void*) Memcpy((double*)R_alloc(nnz, sizeof(double)),
202  (double*)(src->x), nnz);
203 }
204 
208 static void chTr2Ralloc(CHM_TR dest, CHM_TR src)
209 {
210  int nnz;
211 
212  /* copy all the (non-pointer) characteristics of src to dest */
213  memcpy(dest, src, sizeof(cholmod_triplet));
214 
215  /* R_alloc the vector storage for dest and copy the contents from src */
216  nnz = src->nnz;
217  dest->i = (void*) Memcpy((int*)R_alloc(nnz, sizeof(int)),
218  (int*)(src->i), nnz);
219  dest->j = (void*) Memcpy((int*)R_alloc(nnz, sizeof(int)),
220  (int*)(src->j), nnz);
221  if(src->xtype)
222  dest->x = (void*) Memcpy((double*)R_alloc(nnz, sizeof(double)),
223  (double*)(src->x), nnz);
224 }
225 
246  Rboolean check_Udiag, Rboolean sort_in_place)
247 {
248  static const char *valid[] = { MATRIX_VALID_Csparse, ""};
249  int *dims = INTEGER(GET_SLOT(x, Matrix_DimSym)),
250  ctype = R_check_class_etc(x, valid);
251  SEXP islot = GET_SLOT(x, Matrix_iSym);
252 
253  if (ctype < 0) error(_("invalid class of object to as_cholmod_sparse"));
254  if (!isValid_Csparse(x))
255  error(_("invalid object passed to as_cholmod_sparse"));
256  memset(ans, 0, sizeof(cholmod_sparse)); /* zero the struct */
257 
258  ans->itype = CHOLMOD_INT; /* characteristics of the system */
259  ans->dtype = CHOLMOD_DOUBLE;
260  ans->packed = TRUE;
261  /* slots always present */
262  ans->i = INTEGER(islot);
263  ans->p = INTEGER(GET_SLOT(x, Matrix_pSym));
264  /* dimensions and nzmax */
265  ans->nrow = dims[0];
266  ans->ncol = dims[1];
267  /* Allow for over-allocation of the i and x slots. Needed for
268  * sparse X form in lme4. Right now it looks too difficult to
269  * check for the length of the x slot, because of the xpt
270  * utility, but the lengths of x and i should agree. */
271  ans->nzmax = LENGTH(islot);
272  /* values depending on ctype */
273  ans->x = xpt(ctype, x);
274  ans->stype = stype(ctype, x);
275  ans->xtype = xtype(ctype);
276 
277  /* are the columns sorted (increasing row numbers) ?*/
278  ans->sorted = check_sorted_chm(ans);
279  if (!(ans->sorted)) { /* sort columns */
280  if(sort_in_place) {
281  if (!cholmod_sort(ans, &c))
282  error(_("in_place cholmod_sort returned an error code"));
283  ans->sorted = 1;
284  }
285  else {
286  CHM_SP tmp = cholmod_copy_sparse(ans, &c);
287  if (!cholmod_sort(tmp, &c))
288  error(_("cholmod_sort returned an error code"));
289 
290 #ifdef DEBUG_Matrix
291  /* This "triggers" exactly for return values of dtCMatrix_sparse_solve():*/
292  /* Don't want to translate this: want it report */
293  Rprintf("Note: as_cholmod_sparse() needed cholmod_sort()ing\n");
294 #endif
295  chm2Ralloc(ans, tmp);
296  cholmod_free_sparse(&tmp, &c);
297  }
298  }
299 
300  if (check_Udiag && ctype % 3 == 2 // triangular
301  && (*diag_P(x) == 'U')) { /* diagU2N(.) "in place" : */
302  double one[] = {1, 0};
303  CHM_SP eye = cholmod_speye(ans->nrow, ans->ncol, ans->xtype, &c);
304  CHM_SP tmp = cholmod_add(ans, eye, one, one, TRUE, TRUE, &c);
305 
306 #ifdef DEBUG_Matrix_verbose /* happens quite often, e.g. in ../tests/indexing.R : */
307  Rprintf("Note: as_cholmod_sparse(<ctype=%d>) - diagU2N\n", ctype);
308 #endif
309  chm2Ralloc(ans, tmp);
310  cholmod_free_sparse(&tmp, &c);
311  cholmod_free_sparse(&eye, &c);
312  } /* else :
313  * NOTE: if(*diag_P(x) == 'U'), the diagonal is lost (!);
314  * ---- that may be ok, e.g. if we are just converting from/to Tsparse,
315  * but is *not* at all ok, e.g. when used before matrix products */
316 
317  return ans;
318 }
319 
335 SEXP chm_sparse_to_SEXP(CHM_SP a, int dofree, int uploT, int Rkind,
336  const char* diag, SEXP dn)
337 {
338  SEXP ans;
339  char *cls = "";/* -Wall */
340  int *dims, nnz, *ansp, *ansi, *aii = (int*)(a->i), *api = (int*)(a->p),
341  longi = (a->itype) == CHOLMOD_LONG;
342  SuiteSparse_long *ail = (SuiteSparse_long*)(a->i), *apl = (SuiteSparse_long*)(a->p);
343 
344  PROTECT(dn); /* dn is usually UNPROTECTed before the call */
345 
346  /* ensure a is sorted and packed */
347  if (!a->sorted || !a->packed)
348  longi ? cholmod_l_sort(a, &cl) : cholmod_sort(a, &c);
349  /* determine the class of the result */
350 
351 #define DOFREE_MAYBE \
352  if (dofree > 0) \
353  longi ? cholmod_l_free_sparse(&a, &cl) : cholmod_free_sparse(&a, &c); \
354  else if (dofree < 0) Free(a)
355 
356  switch(a->xtype){
357  case CHOLMOD_PATTERN:
358  cls = uploT ? "ntCMatrix": ((a->stype) ? "nsCMatrix" : "ngCMatrix");
359  break;
360  case CHOLMOD_REAL:
361  switch(Rkind) {
362  case 0:
363  cls = uploT ? "dtCMatrix": ((a->stype) ? "dsCMatrix" : "dgCMatrix");
364  break;
365  case 1:
366  cls = uploT ? "ltCMatrix": ((a->stype) ? "lsCMatrix" : "lgCMatrix");
367  break;
368  default:
369  DOFREE_MAYBE;
370  error(_("chm_sparse_to_SEXP(<real>, *): invalid 'Rkind' (real kind code)"));
371  }
372  break;
373  case CHOLMOD_COMPLEX:
374  cls = uploT ? "ztCMatrix": ((a->stype) ? "zsCMatrix" : "zgCMatrix");
375  break;
376  default:
377  DOFREE_MAYBE;
378  error(_("unknown xtype in cholmod_sparse object"));
379  }
380  ans = PROTECT(NEW_OBJECT(MAKE_CLASS(cls)));
381  /* allocate and copy common slots */
382  nnz = longi ? cholmod_l_nnz(a, &cl) : cholmod_nnz(a, &c);
383  dims = INTEGER(ALLOC_SLOT(ans, Matrix_DimSym, INTSXP, 2));
384  dims[0] = a->nrow; dims[1] = a->ncol;
385  ansp = INTEGER(ALLOC_SLOT(ans, Matrix_pSym, INTSXP, a->ncol + 1));
386  ansi = INTEGER(ALLOC_SLOT(ans, Matrix_iSym, INTSXP, nnz));
387  for (int j = 0; j <= a->ncol; j++) ansp[j] = longi ? (int)(apl[j]) : api[j];
388  for (int p = 0; p < nnz; p++) ansi[p] = longi ? (int)(ail[p]) : aii[p];
389  /* copy data slot if present */
390  if (a->xtype == CHOLMOD_REAL) {
391  int i, *m_x;
392  double *a_x = (double *) a->x;
393  switch(Rkind) {
394  case 0:
395  Memcpy(REAL(ALLOC_SLOT(ans, Matrix_xSym, REALSXP, nnz)),
396  a_x, nnz);
397  break;
398  case 1:
399  m_x = LOGICAL(ALLOC_SLOT(ans, Matrix_xSym, LGLSXP, nnz));
400  for (i=0; i < nnz; i++)
401  m_x[i] = ISNAN(a_x[i]) ? NA_LOGICAL : (a_x[i] != 0);
402  break;
403  }
404  }
405  else if (a->xtype == CHOLMOD_COMPLEX) {
406  DOFREE_MAYBE;
407  error(_("complex sparse matrix code not yet written"));
408 /* Memcpy(COMPLEX(ALLOC_SLOT(ans, Matrix_xSym, CPLXSXP, nnz)), */
409 /* (complex *) a->x, nnz); */
410  }
411  if (uploT) { /* slots for triangularMatrix */
412  if (a->stype) error(_("Symmetric and triangular both set"));
413  SET_SLOT(ans, Matrix_uploSym, mkString((uploT > 0) ? "U" : "L"));
414  SET_SLOT(ans, Matrix_diagSym, mkString(diag));
415  }
416  if (a->stype) /* slot for symmetricMatrix */
417  SET_SLOT(ans, Matrix_uploSym,
418  mkString((a->stype > 0) ? "U" : "L"));
419  DOFREE_MAYBE;
420  if (dn != R_NilValue)
421  SET_SLOT(ans, Matrix_DimNamesSym, duplicate(dn));
422 
423  UNPROTECT(2);
424  return ans;
425 }
426 #undef DOFREE_MAYBE
427 
428 
438 Rboolean chm_MOD_xtype(int to_xtype, cholmod_sparse *A, CHM_CM Common) {
439 // *MOD*: shouting, as A is modified in place
440 
441 /* --------------------------------------------------------------------------
442  * cholmod_sparse_xtype: change the xtype of a sparse matrix
443  * --------------------------------------------------------------------------
444  int cholmod_sparse_xtype
445  (
446  // ---- input ----
447  int to_xtype, //
448  // ---- in/out ---
449  cholmod_sparse *A, //
450  // ---------------
451  cholmod_common *Common
452  ) ;
453 
454  int cholmod_l_sparse_xtype (int, cholmod_sparse *, cholmod_common *) ;
455 */
456  if((A->itype) == CHOLMOD_LONG) {
457  return (Rboolean) cholmod_l_sparse_xtype (to_xtype, A, Common);
458  } else {
459  return (Rboolean) cholmod_sparse_xtype (to_xtype, A, Common);
460  }
461 }
462 
463 
480 CHM_TR as_cholmod_triplet(CHM_TR ans, SEXP x, Rboolean check_Udiag)
481 {
482  static const char *valid[] = { MATRIX_VALID_Tsparse, ""};
483  int ctype = R_check_class_etc(x, valid),
484  *dims = INTEGER(GET_SLOT(x, Matrix_DimSym));
485  SEXP islot = GET_SLOT(x, Matrix_iSym);
486  int m = LENGTH(islot);
487  Rboolean do_Udiag = (check_Udiag && ctype % 3 == 2 && (*diag_P(x) == 'U'));
488  if (ctype < 0) error(_("invalid class of object to as_cholmod_triplet"));
489 
490  memset(ans, 0, sizeof(cholmod_triplet)); /* zero the struct */
491 
492  ans->itype = CHOLMOD_INT; /* characteristics of the system */
493  ans->dtype = CHOLMOD_DOUBLE;
494  /* nzmax, dimensions, types and slots : */
495  ans->nnz = ans->nzmax = m;
496  ans->nrow = dims[0];
497  ans->ncol = dims[1];
498  ans->stype = stype(ctype, x);
499  ans->xtype = xtype(ctype);
500  ans->i = (void *) INTEGER(islot);
501  ans->j = (void *) INTEGER(GET_SLOT(x, Matrix_jSym));
502  ans->x = xpt(ctype, x);
503 
504  if(do_Udiag) {
505  /* diagU2N(.) "in place", similarly to Tsparse_diagU2N [./Tsparse.c]
506  (but without new SEXP): */
507  int k = m + dims[0];
508  CHM_TR tmp = cholmod_l_copy_triplet(ans, &c);
509  int *a_i, *a_j;
510 
511  if(!cholmod_reallocate_triplet((size_t) k, tmp, &c))
512  error(_("as_cholmod_triplet(): could not reallocate for internal diagU2N()"
513  ));
514 
515  /* TODO? instead of copy_triplet() & reallocate_triplet()
516  * ---- allocate to correct length + Memcpy() here, as in
517  * Tsparse_diagU2N() & chTr2Ralloc() below */
518  a_i = tmp->i;
519  a_j = tmp->j;
520  /* add (@i, @j)[k+m] = k, @x[k+m] = 1. for k = 0,..,(n-1) */
521  for(k=0; k < dims[0]; k++) {
522  a_i[k+m] = k;
523  a_j[k+m] = k;
524 
525  switch(ctype / 3) {
526  case 0: { /* "d" */
527  double *a_x = tmp->x;
528  a_x[k+m] = 1.;
529  break;
530  }
531  case 1: { /* "l" */
532  int *a_x = tmp->x;
533  a_x[k+m] = 1;
534  break;
535  }
536  case 2: /* "n" */
537  break;
538  case 3: { /* "z" */
539  double *a_x = tmp->x;
540  a_x[2*(k+m) ] = 1.;
541  a_x[2*(k+m)+1] = 0.;
542  break;
543  }
544  }
545  } /* for(k) */
546 
547  chTr2Ralloc(ans, tmp);
548  cholmod_l_free_triplet(&tmp, &c);
549 
550  } /* else :
551  * NOTE: if(*diag_P(x) == 'U'), the diagonal is lost (!);
552  * ---- that may be ok, e.g. if we are just converting from/to Tsparse,
553  * but is *not* at all ok, e.g. when used before matrix products */
554 
555  return ans;
556 }
557 
573 SEXP chm_triplet_to_SEXP(CHM_TR a, int dofree, int uploT, int Rkind,
574  const char* diag, SEXP dn)
575 {
576  SEXP ans;
577  char *cl = ""; /* -Wall */
578  int *dims;
579 
580  PROTECT(dn); /* dn is usually UNPROTECTed before the call */
581  /* determine the class of the result */
582 
583 #define DOFREE_MAYBE \
584  if (dofree > 0) cholmod_free_triplet(&a, &c); \
585  else if (dofree < 0) Free(a)
586 
587  switch(a->xtype) {
588  case CHOLMOD_PATTERN:
589  cl = uploT ? "ntTMatrix" :
590  ((a->stype) ? "nsTMatrix" : "ngTMatrix");
591  break;
592  case CHOLMOD_REAL:
593  switch(Rkind) {
594  case 0:
595  cl = uploT ? "dtTMatrix" :
596  ((a->stype) ? "dsTMatrix" : "dgTMatrix");
597  break;
598  case 1:
599  cl = uploT ? "ltTMatrix" :
600  ((a->stype) ? "lsTMatrix" : "lgTMatrix");
601  break;
602  }
603  break;
604  case CHOLMOD_COMPLEX:
605  cl = uploT ? "ztTMatrix" :
606  ((a->stype) ? "zsTMatrix" : "zgTMatrix");
607  break;
608  default:
609  DOFREE_MAYBE;
610  error(_("unknown xtype in cholmod_triplet object"));
611  }
612  ans = PROTECT(NEW_OBJECT(MAKE_CLASS(cl)));
613  /* allocate and copy common slots */
614  dims = INTEGER(ALLOC_SLOT(ans, Matrix_DimSym, INTSXP, 2));
615  dims[0] = a->nrow; dims[1] = a->ncol;
616  Memcpy(INTEGER(ALLOC_SLOT(ans, Matrix_iSym, INTSXP, a->nnz)),
617  (int *) a->i, a->nnz);
618  Memcpy(INTEGER(ALLOC_SLOT(ans, Matrix_jSym, INTSXP, a->nnz)),
619  (int *) a->j, a->nnz);
620  /* copy data slot if present */
621  if (a->xtype == CHOLMOD_REAL) {
622  int i, *m_x;
623  double *a_x = (double *) a->x;
624  switch(Rkind) {
625  case 0:
626  Memcpy(REAL(ALLOC_SLOT(ans, Matrix_xSym, REALSXP, a->nnz)),
627  a_x, a->nnz);
628  break;
629  case 1:
630  m_x = LOGICAL(ALLOC_SLOT(ans, Matrix_xSym, LGLSXP, a->nnz));
631  for (i=0; i < a->nnz; i++)
632  m_x[i] = ISNAN(a_x[i]) ? NA_LOGICAL : (a_x[i] != 0);
633  break;
634  }
635  }
636  else if (a->xtype == CHOLMOD_COMPLEX) {
637  DOFREE_MAYBE;
638  error(_("complex sparse matrix code not yet written"));
639 /* Memcpy(COMPLEX(ALLOC_SLOT(ans, Matrix_xSym, CPLXSXP, a->nnz)), */
640 /* (complex *) a->x, a->nz); */
641  }
642  if (uploT) { /* slots for triangularMatrix */
643  if (a->stype) error(_("Symmetric and triangular both set"));
644  SET_SLOT(ans, Matrix_uploSym, mkString((uploT > 0) ? "U" : "L"));
645  SET_SLOT(ans, Matrix_diagSym, mkString(diag));
646  }
647  /* set symmetry attributes */
648  if (a->stype)
649  SET_SLOT(ans, Matrix_uploSym,
650  mkString((a->stype > 0) ? "U" : "L"));
651  DOFREE_MAYBE;
652  if (dn != R_NilValue)
653  SET_SLOT(ans, Matrix_DimNamesSym, duplicate(dn));
654  UNPROTECT(2);
655  return ans;
656 }
657 #undef DOFREE_MAYBE
658 
673 {
674 #define _AS_cholmod_dense_1 \
675  static const char *valid[] = { MATRIX_VALID_ge_dense, ""}; \
676  int dims[2], ctype = R_check_class_etc(x, valid), nprot = 0; \
677  \
678  if (ctype < 0) { /* not a classed matrix */ \
679  if (isMatrix(x)) Memcpy(dims, INTEGER(getAttrib(x, R_DimSymbol)), 2); \
680  else {dims[0] = LENGTH(x); dims[1] = 1;} \
681  if (isInteger(x)) { \
682  x = PROTECT(coerceVector(x, REALSXP)); \
683  nprot++; \
684  } \
685  ctype = (isReal(x) ? 0 : \
686  (isLogical(x) ? 2 : /* logical -> default to "l", not "n" */ \
687  (isComplex(x) ? 6 : -1))); \
688  } else Memcpy(dims, INTEGER(GET_SLOT(x, Matrix_DimSym)), 2); \
689  if (ctype < 0) error(_("invalid class of object to as_cholmod_dense")); \
690  memset(ans, 0, sizeof(cholmod_dense)); /* zero the struct */ \
691  \
692  ans->dtype = CHOLMOD_DOUBLE; /* characteristics of the system */ \
693  ans->x = ans->z = (void *) NULL; \
694  /* dimensions and nzmax */ \
695  ans->d = ans->nrow = dims[0]; \
696  ans->ncol = dims[1]; \
697  ans->nzmax = ((size_t)dims[0]) * dims[1]; \
698  /* set the xtype and any elements */ \
699  switch(ctype / 2) { \
700  case 0: /* "d" */ \
701  ans->xtype = CHOLMOD_REAL; \
702  ans->x = (void *) REAL((ctype % 2) ? GET_SLOT(x, Matrix_xSym) : x); \
703  break
704 
706 
707  case 1: /* "l" */
708  ans->xtype = CHOLMOD_REAL;
709  ans->x = RallocedREAL((ctype % 2) ? GET_SLOT(x, Matrix_xSym) : x);
710  break;
711  case 2: /* "n" */
712  ans->xtype = CHOLMOD_PATTERN;
713  ans->x = (void *) LOGICAL((ctype % 2) ? GET_SLOT(x, Matrix_xSym) : x);
714  break;
715 
716 #define _AS_cholmod_dense_2 \
717  case 3: /* "z" */ \
718  ans->xtype = CHOLMOD_COMPLEX; \
719  ans->x = (void *) COMPLEX((ctype % 2) ? GET_SLOT(x, Matrix_xSym) : x); \
720  break; \
721  } \
722  UNPROTECT(nprot); \
723  return ans
724 
726 }
727 
728 /* version of as_cholmod_dense() that produces a cholmod_dense matrix
729  * with REAL 'x' slot -- i.e. treats "nMatrix" as "lMatrix" -- as only difference;
730  * Not just via a flag in as_cholmod_dense() since that has fixed API */
732 {
734 
735  case 1: /* "l" */
736  case 2: /* "n" (no NA in 'x', but *has* 'x' slot => treat as "l" */
737  ans->xtype = CHOLMOD_REAL;
738  ans->x = RallocedREAL((ctype % 2) ? GET_SLOT(x, Matrix_xSym) : x);
739  break;
740 
742 }
743 
744 #undef _AS_cholmod_dense_1
745 #undef _AS_cholmod_dense_2
746 
755 {
756  if (x->xtype != CHOLMOD_REAL)
757  error(_("chm_transpose_dense(ans, x) not yet implemented for %s different from %s"),
758  "x->xtype", "CHOLMOD_REAL");
759  double *xx = x->x, *ansx = ans->x;
760  // Inspired from R's do_transpose() in .../R/src/main/array.c :
761  int i,j, nrow = x->nrow, len = x->nzmax, l_1 = len-1;
762  for (i = 0, j = 0; i < len; i++, j += nrow) {
763  if (j > l_1) j -= l_1;
764  ansx[i] = xx[j];
765  }
766  return;
767 }
768 
769 void R_cholmod_error(int status, const char *file, int line, const char *message)
770 {
771  CHM_restore_common(); /* restore any setting that may have been changed */
772 
773 /* NB: keep in sync with M_R_cholmod_error(), ../inst/include/Matrix_stubs.c */
774 
775  /* From CHOLMOD/Include/cholmod_core.h : ...status values.
776  zero means success, negative means a fatal error, positive is a warning.
777  */
778 #ifndef R_CHOLMOD_ALWAYS_ERROR
779  if(status < 0) {
780 #endif
781  error(_("Cholmod error '%s' at file %s, line %d"), message, file, line);
782 #ifndef R_CHOLMOD_ALWAYS_ERROR
783  }
784  else
785  warning(_("Cholmod warning '%s' at file %s, line %d"),
786  message, file, line);
787 #endif
788 }
789 
790 /* just to get 'int' instead of 'void' as required by CHOLMOD's print_function */
791 static
792 int R_cholmod_printf(const char* fmt, ...)
793 {
794  va_list(ap);
795 
796  va_start(ap, fmt);
797  Rprintf((char *)fmt, ap);
798  va_end(ap);
799  return 0;
800 }
801 
811 {
812  int res;
813  if (!(res = cholmod_start(c)))
814  error(_("Unable to initialize cholmod: error code %d"), res);
815  c->print_function = R_cholmod_printf; /* Rprintf gives warning */
816  /* Since we provide an error handler, it may not be a good idea to allow CHOLMOD printing,
817  * because that's not easily suppressed on the R level :
818  * Hence consider, at least temporarily * c->print_function = NULL; */
819  c->error_handler = R_cholmod_error;
820  return TRUE;
821 }
822 
835 SEXP chm_dense_to_SEXP(CHM_DN a, int dofree, int Rkind, SEXP dn, Rboolean transp)
836 {
837 /* FIXME: should also have args (int uploST, char *diag) */
838  SEXP ans;
839  char *cl = ""; /* -Wall */
840  int *dims, ntot;
841 
842  PROTECT(dn); // <-- no longer protected in caller
843 
844 #define DOFREE_de_MAYBE \
845  if (dofree > 0) cholmod_free_dense(&a, &c); \
846  else if (dofree < 0) Free(a);
847 
848  switch(a->xtype) { /* determine the class of the result */
849 /* CHOLMOD_PATTERN never happens because cholmod_dense can't :
850  * case CHOLMOD_PATTERN:
851  * cl = "ngeMatrix"; break;
852  */
853  case CHOLMOD_REAL:
854  switch(Rkind) { /* -1: special for this function! */
855  case -1: cl = "ngeMatrix"; break;
856  case 0: cl = "dgeMatrix"; break;
857  case 1: cl = "lgeMatrix"; break;
858  default:
860  error(_("unknown 'Rkind'"));
861  }
862  break;
863  case CHOLMOD_COMPLEX:
864  cl = "zgeMatrix"; break;
865  default:
867  error(_("unknown xtype"));
868  }
869 
870  ans = PROTECT(NEW_OBJECT(MAKE_CLASS(cl)));
871  /* allocate and copy common slots */
872  dims = INTEGER(ALLOC_SLOT(ans, Matrix_DimSym, INTSXP, 2));
873  if(transp) {
874  dims[1] = a->nrow; dims[0] = a->ncol;
875  } else {
876  dims[0] = a->nrow; dims[1] = a->ncol;
877  }
878  ntot = ((size_t)dims[0]) * dims[1];
879  if (a->d == a->nrow) { /* copy data slot -- always present in dense(!) */
880  if (a->xtype == CHOLMOD_REAL) {
881  int i, *m_x;
882  double *ansx, *a_x = (double *) a->x;
883  switch(Rkind) {
884  case 0:
885  ansx = REAL(ALLOC_SLOT(ans, Matrix_xSym, REALSXP, ntot));
886  if(transp) {
887  // Inspired from R's do_transpose() in .../R/src/main/array.c :
888  int i,j, nrow = a->nrow, len = ntot, l_1 = len-1;
889  for (i = 0, j = 0; i < len; i++, j += nrow) {
890  if (j > l_1) j -= l_1;
891  ansx[i] = a_x[j];
892  }
893  } else {
894  Memcpy(ansx, a_x, ntot);
895  }
896  break;
897  case -1: /* nge*/
898  case 1: /* lge*/
899  m_x = LOGICAL(ALLOC_SLOT(ans, Matrix_xSym, LGLSXP, ntot));
900  if(transp) {
901  // Inspired from R's do_transpose() in .../R/src/main/array.c :
902  int i,j, nrow = a->nrow, len = ntot, l_1 = len-1;
903  for (i = 0, j = 0; i < len; i++, j += nrow) {
904  if (j > l_1) j -= l_1;
905  m_x[i] = a_x[j];
906  }
907  } else {
908  for (i=0; i < ntot; i++)
909  m_x[i] = ISNAN(a_x[i]) ? NA_LOGICAL : (a_x[i] != 0);
910  }
911  break;
912  }
913  }
914  else if (a->xtype == CHOLMOD_COMPLEX) {
916  error(_("complex sparse matrix code not yet written"));
917 /* Memcpy(COMPLEX(ALLOC_SLOT(ans, Matrix_xSym, CPLXSXP, ntot)), */
918 /* (complex *) a->x, ntot); */
919  }
920  } else {
922  error(_("code for cholmod_dense with holes not yet written"));
923  }
924 
926  if (dn != R_NilValue)
927  SET_SLOT(ans, Matrix_DimNamesSym, duplicate(dn));
928  UNPROTECT(2);
929  return ans;
930 }
931 
942 SEXP chm_dense_to_matrix(CHM_DN a, int dofree, SEXP dn)
943 {
944 #define CHM_DENSE_TYPE \
945  SEXPTYPE typ; \
946  /* determine the class of the result */ \
947  typ = (a->xtype == CHOLMOD_PATTERN) ? LGLSXP : \
948  ((a->xtype == CHOLMOD_REAL) ? REALSXP : \
949  ((a->xtype == CHOLMOD_COMPLEX) ? CPLXSXP : NILSXP)); \
950  if (typ == NILSXP) { \
951  DOFREE_de_MAYBE; \
952  error(_("unknown xtype")); \
953  }
954 
955  PROTECT(dn);
957 
958  SEXP ans = PROTECT(allocMatrix(typ, a->nrow, a->ncol));
959 
960 #define CHM_DENSE_COPY_DATA \
961  if (a->d == a->nrow) { /* copy data slot if present */ \
962  if (a->xtype == CHOLMOD_REAL) \
963  Memcpy(REAL(ans), (double *) a->x, a->nrow * a->ncol); \
964  else if (a->xtype == CHOLMOD_COMPLEX) { \
965  DOFREE_de_MAYBE; \
966  error(_("complex sparse matrix code not yet written")); \
967 /* Memcpy(COMPLEX(ALLOC_SLOT(ans, Matrix_xSym, CPLXSXP, a->nnz)), */ \
968 /* (complex *) a->x, a->nz); */ \
969  } else if (a->xtype == CHOLMOD_PATTERN) { \
970  DOFREE_de_MAYBE; \
971  error(_("don't know if a dense pattern matrix makes sense")); \
972  } \
973  } else { \
974  DOFREE_de_MAYBE; \
975  error(_("code for cholmod_dense with holes not yet written")); \
976  }
977 
979 
981  if (dn != R_NilValue)
982  setAttrib(ans, R_DimNamesSymbol, duplicate(dn));
983  UNPROTECT(2);
984  return ans;
985 }
986 
996 SEXP chm_dense_to_vector(CHM_DN a, int dofree)
997 {
999 
1000  SEXP ans = PROTECT(allocVector(typ, a->nrow * a->ncol));
1003  UNPROTECT(1);
1004  return ans;
1005 }
1006 
1007 CHM_DN numeric_as_chm_dense(CHM_DN ans, double *v, int nr, int nc)
1008 {
1009  ans->d = ans->nrow = nr;
1010  ans->ncol = nc;
1011  ans->nzmax = nr * nc;
1012  ans->x = (void *) v;
1013  ans->xtype = CHOLMOD_REAL;
1014  ans->dtype = CHOLMOD_DOUBLE;
1015  return ans;
1016 }
1017 
1032 {
1033  static const char *valid[] = { MATRIX_VALID_CHMfactor, ""};
1034  int *type = INTEGER(GET_SLOT(x, install("type"))),
1035  ctype = R_check_class_etc(x, valid);
1036  SEXP tmp;
1037 
1038  if (ctype < 0) error(_("invalid class of object to as_cholmod_factor"));
1039  memset(ans, 0, sizeof(cholmod_factor)); /* zero the struct */
1040 
1041  ans->itype = CHOLMOD_INT; /* characteristics of the system */
1042  ans->dtype = CHOLMOD_DOUBLE;
1043  ans->z = (void *) NULL;
1044  ans->xtype = (ctype < 2) ? CHOLMOD_REAL : CHOLMOD_PATTERN;
1045 
1046  ans->ordering = type[0]; /* unravel the type */
1047  ans->is_ll = (type[1] ? 1 : 0);
1048  ans->is_super = (type[2] ? 1 : 0);
1049  ans->is_monotonic = (type[3] ? 1 : 0);
1050  /* check for consistency */
1051  if ((!(ans->is_ll)) && ans->is_super)
1052  error(_("Supernodal LDL' decomposition not available"));
1053  if ((!type[2]) ^ (ctype % 2))
1054  error(_("Supernodal/simplicial class inconsistent with type flags"));
1055  /* slots always present */
1056  tmp = GET_SLOT(x, Matrix_permSym);
1057  ans->minor = ans->n = LENGTH(tmp); ans->Perm = INTEGER(tmp);
1058  ans->ColCount = INTEGER(GET_SLOT(x, install("colcount")));
1059  ans->z = ans->x = (void *) NULL;
1060  if (ctype < 2) {
1061  tmp = GET_SLOT(x, Matrix_xSym);
1062  ans->x = REAL(tmp);
1063  }
1064  if (ans->is_super) { /* supernodal factorization */
1065  ans->xsize = LENGTH(tmp);
1066  ans->maxcsize = type[4]; ans->maxesize = type[5];
1067  ans->i = (int*)NULL;
1068  tmp = GET_SLOT(x, install("super"));
1069  ans->nsuper = LENGTH(tmp) - 1; ans->super = INTEGER(tmp);
1070  /* Move these checks to the CHMfactor_validate function */
1071  if (ans->nsuper < 1)
1072  error(_("Number of supernodes must be positive when is_super is TRUE"));
1073  tmp = GET_SLOT(x, install("pi"));
1074  if (LENGTH(tmp) != ans->nsuper + 1)
1075  error(_("Lengths of super and pi must be equal"));
1076  ans->pi = INTEGER(tmp);
1077  tmp = GET_SLOT(x, install("px"));
1078  if (LENGTH(tmp) != ans->nsuper + 1)
1079  error(_("Lengths of super and px must be equal"));
1080  ans->px = INTEGER(tmp);
1081  tmp = GET_SLOT(x, install("s"));
1082  ans->ssize = LENGTH(tmp); ans->s = INTEGER(tmp);
1083  } else {
1084  ans->nzmax = LENGTH(tmp);
1085  ans->p = INTEGER(GET_SLOT(x, Matrix_pSym));
1086  ans->i = INTEGER(GET_SLOT(x, Matrix_iSym));
1087  ans->nz = INTEGER(GET_SLOT(x, install("nz")));
1088  ans->next = INTEGER(GET_SLOT(x, install("nxt")));
1089  ans->prev = INTEGER(GET_SLOT(x, install("prv")));
1090  }
1091  if (!cholmod_check_factor(ans, &c))
1092  error(_("failure in as_cholmod_factor"));
1093  return ans;
1094 }
1095 
1096 
1106 SEXP chm_factor_to_SEXP(CHM_FR f, int dofree)
1107 {
1108  SEXP ans;
1109  int *dims, *type;
1110  char *class = (char*) NULL; /* -Wall */
1111 
1112 #define DOFREE_MAYBE \
1113  if(dofree) { \
1114  if (dofree > 0) cholmod_free_factor(&f, &c); \
1115  else /* dofree < 0 */ Free(f); \
1116  }
1117 
1118  if(!chm_factor_ok(f)) {
1119  DOFREE_MAYBE;
1120  error(_("CHOLMOD factorization was unsuccessful"));
1121  // error(_("previous CHOLMOD factorization was unsuccessful"));
1122  }
1123 
1124  switch(f->xtype) {
1125  case CHOLMOD_REAL:
1126  class = f->is_super ? "dCHMsuper" : "dCHMsimpl";
1127  break;
1128  case CHOLMOD_PATTERN:
1129  class = f->is_super ? "nCHMsuper" : "nCHMsimpl";
1130  break;
1131  default:
1132  DOFREE_MAYBE;
1133  error(_("f->xtype of %d not recognized"), f->xtype);
1134  }
1135 
1136  ans = PROTECT(NEW_OBJECT(MAKE_CLASS(class)));
1137  dims = INTEGER(ALLOC_SLOT(ans, Matrix_DimSym, INTSXP, 2));
1138  dims[0] = dims[1] = f->n;
1139  /* copy component of known length */
1140  Memcpy(INTEGER(ALLOC_SLOT(ans, Matrix_permSym, INTSXP, f->n)),
1141  (int*)f->Perm, f->n);
1142  Memcpy(INTEGER(ALLOC_SLOT(ans, install("colcount"), INTSXP, f->n)),
1143  (int*)f->ColCount, f->n);
1144  type = INTEGER(ALLOC_SLOT(ans, install("type"), INTSXP, f->is_super ? 6 : 4));
1145  type[0] = f->ordering; type[1] = f->is_ll;
1146  type[2] = f->is_super; type[3] = f->is_monotonic;
1147  if (f->is_super) {
1148  type[4] = f->maxcsize; type[5] = f->maxesize;
1149  Memcpy(INTEGER(ALLOC_SLOT(ans, install("super"), INTSXP, f->nsuper + 1)),
1150  (int*)f->super, f->nsuper+1);
1151  Memcpy(INTEGER(ALLOC_SLOT(ans, install("pi"), INTSXP, f->nsuper + 1)),
1152  (int*)f->pi, f->nsuper + 1);
1153  Memcpy(INTEGER(ALLOC_SLOT(ans, install("px"), INTSXP, f->nsuper + 1)),
1154  (int*)f->px, f->nsuper + 1);
1155  Memcpy(INTEGER(ALLOC_SLOT(ans, install("s"), INTSXP, f->ssize)),
1156  (int*)f->s, f->ssize);
1157  Memcpy(REAL(ALLOC_SLOT(ans, Matrix_xSym, REALSXP, f->xsize)),
1158  (double*)f->x, f->xsize);
1159  } else {
1160  Memcpy(INTEGER(ALLOC_SLOT(ans, Matrix_iSym, INTSXP, f->nzmax)),
1161  (int*)f->i, f->nzmax);
1162  Memcpy(INTEGER(ALLOC_SLOT(ans, Matrix_pSym, INTSXP, f->n + 1)),
1163  (int*)f->p, f->n + 1);
1164  Memcpy(REAL(ALLOC_SLOT(ans, Matrix_xSym, REALSXP, f->nzmax)),
1165  (double*)f->x, f->nzmax);
1166  Memcpy(INTEGER(ALLOC_SLOT(ans, install("nz"), INTSXP, f->n)),
1167  (int*)f->nz, f->n);
1168  Memcpy(INTEGER(ALLOC_SLOT(ans, install("nxt"), INTSXP, f->n + 2)),
1169  (int*)f->next, f->n + 2);
1170  Memcpy(INTEGER(ALLOC_SLOT(ans, install("prv"), INTSXP, f->n + 2)),
1171  (int*)f->prev, f->n + 2);
1172 
1173  }
1174  DOFREE_MAYBE;
1175  UNPROTECT(1);
1176  return ans;
1177 }
1178 #undef DOFREE_MAYBE
1179 
1191 void chm_diagN2U(CHM_SP chx, int uploT, Rboolean do_realloc)
1192 {
1193  int i, n = chx->nrow, nnz = (int)cholmod_nnz(chx, &c),
1194  n_nnz = nnz - n, /* the new nnz : we will have removed n entries */
1195  i_to = 0, i_from = 0;
1196 
1197  if(chx->ncol != n)
1198  error(_("chm_diagN2U(<non-square matrix>): nrow=%d, ncol=%d"),
1199  n, chx->ncol);
1200 
1201  if (!chx->sorted || !chx->packed) cholmod_sort(chx, &c);
1202  /* dimensions and nzmax */
1203 
1204 #define _i(I) ( (int*) chx->i)[I]
1205 #define _x(I) ((double*) chx->x)[I]
1206 #define _p(I) ( (int*) chx->p)[I]
1207 
1208  /* work by copying from i_from to i_to ==> MUST i_to <= i_from */
1209 
1210  if(uploT == 1) { /* "U" : upper triangular */
1211 
1212  for(i = 0; i < n; i++) { /* looking at i-th column */
1213  int j, n_i = _p(i+1) - _p(i); /* = #{entries} in this column */
1214 
1215  /* 1) copy all but the last _above-diagonal_ column-entries: */
1216  for(j = 1; j < n_i; j++, i_to++, i_from++) {
1217  _i(i_to) = _i(i_from);
1218  _x(i_to) = _x(i_from);
1219  }
1220 
1221  /* 2) drop the last column-entry == diagonal entry */
1222  i_from++;
1223  }
1224  }
1225  else if(uploT == -1) { /* "L" : lower triangular */
1226 
1227  for(i = 0; i < n; i++) { /* looking at i-th column */
1228  int j, n_i = _p(i+1) - _p(i); /* = #{entries} in this column */
1229 
1230  /* 1) drop the first column-entry == diagonal entry */
1231  i_from++;
1232 
1233  /* 2) copy the other _below-diagonal_ column-entries: */
1234  for(j = 1; j < n_i; j++, i_to++, i_from++) {
1235  _i(i_to) = _i(i_from);
1236  _x(i_to) = _x(i_from);
1237  }
1238  }
1239  }
1240  else {
1241  error(_("chm_diagN2U(x, uploT = %d): uploT should be +- 1"), uploT);
1242  }
1243 
1244  /* the column pointers are modified the same in both cases :*/
1245  for(i=1; i <= n; i++)
1246  _p(i) -= i;
1247 
1248 #undef _i
1249 #undef _x
1250 #undef _p
1251 
1252  if(do_realloc) /* shorten (i- and x-slots from nnz to n_nnz */
1253  cholmod_reallocate_sparse(n_nnz, chx, &c);
1254  return;
1255 }
1256 
1257 /* Placeholders; TODO: use checks above (search "CHMfactor_validate"): */
1258 
1259 SEXP CHMfactor_validate(SEXP obj) /* placeholder */
1260 {
1261  return ScalarLogical(1);
1262 }
1263 
1264 SEXP CHMsimpl_validate(SEXP obj) /* placeholder */
1265 {
1266  return ScalarLogical(1);
1267 }
1268 
1269 SEXP CHMsuper_validate(SEXP obj) /* placeholder */
1270 {
1271  return ScalarLogical(1);
1272 }
1273 
SEXP Matrix_DimSym
Definition: Syms.h:2
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
#define _i(I)
void CHM_store_common()
Definition: chm_common.c:25
cholmod_sparse * CHM_SP
Definition: chm_common.h:25
SEXP Matrix_xSym
Definition: Syms.h:2
static int stype(int ctype, SEXP x)
stype := "symmetry type".
Definition: chm_common.c:114
static int R_cholmod_printf(const char *fmt,...)
Definition: chm_common.c:792
SEXP Matrix_DimNamesSym
Definition: Syms.h:2
SEXP CHM_set_common_env(SEXP rho)
Definition: chm_common.c:74
static SEXP final_superSym
Definition: chm_common.c:19
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
cholmod_triplet * CHM_TR
Definition: chm_common.h:27
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
static SEXP final_packSym
Definition: chm_common.c:19
cholmod_common cl
Definition: chm_common.c:16
SEXP Matrix_uploSym
Definition: Syms.h:2
#define _p(I)
static SEXP supernodalSym
Definition: chm_common.c:19
void R_cholmod_error(int status, const char *file, int line, const char *message)
Definition: chm_common.c:769
SEXP CHMfactor_validate(SEXP obj)
Definition: chm_common.c:1259
SEXP get_SuiteSparse_version()
Definition: chm_common.c:8
SEXP Matrix_jSym
Definition: Syms.h:2
CHM_DN as_cholmod_x_dense(CHM_DN ans, SEXP x)
Definition: chm_common.c:731
static SEXP nmethodsSym
Definition: chm_common.c:19
static void chm2Ralloc(CHM_SP dest, CHM_SP src)
Copy cholmod_sparse, to an R_alloc()ed version of it.
Definition: chm_common.c:186
int R_cholmod_start(CHM_CM c)
Initialize the CHOLMOD library and replace the print and error functions by R-specific versions...
Definition: chm_common.c:810
static SEXP final_monotonicSym
Definition: chm_common.c:19
static SEXP prefer_upperSym
Definition: chm_common.c:19
#define DOFREE_MAYBE
static SEXP quick_return_if_not_posdefSym
Definition: chm_common.c:19
static int xtype(int ctype)
xtype: the kind of numeric (think "x slot") of Cholmod sparse matrices.
Definition: chm_common.c:126
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 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
static SEXP grow1Sym
Definition: chm_common.c:19
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
static SEXP final_llSym
Definition: chm_common.c:19
static SEXP prefer_zomplexSym
Definition: chm_common.c:19
SEXP Matrix_permSym
Definition: Syms.h:2
SEXP chm_common_env
Definition: chm_common.c:18
#define uplo_P(_x_)
Definition: Mutils.h:174
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
#define _(String)
Definition: Mutils.h:32
Rboolean check_sorted_chm(CHM_SP A)
Definition: chm_common.c:169
cholmod_dense * CHM_DN
Definition: chm_common.h:21
SEXP Matrix_iSym
Definition: Syms.h:2
static void * xpt(int ctype, SEXP x)
Definition: chm_common.c:154
#define MATRIX_VALID_CHMfactor
Definition: Mutils.h:419
cholmod_common * CHM_CM
Definition: chm_common.h:20
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
#define _x(I)
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
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
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
static SEXP final_resymbolSym
Definition: chm_common.c:19
#define _AS_cholmod_dense_1
#define diag_P(_x_)
Definition: Mutils.h:175
static SEXP m0_ordSym
Definition: chm_common.c:19
SEXP Matrix_pSym
Definition: Syms.h:2
static SEXP postorderSym
Definition: chm_common.c:19
static SEXP maxrankSym
Definition: chm_common.c:19
SEXP Matrix_diagSym
Definition: Syms.h:2
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
void CHM_restore_common()
Definition: chm_common.c:50
cholmod_common c
Definition: chm_common.c:15
static SEXP final_asisSym
Definition: chm_common.c:19
static void * RallocedREAL(SEXP x)
Definition: chm_common.c:141
#define MATRIX_VALID_Tsparse
Definition: Mutils.h:390
static void chTr2Ralloc(CHM_TR dest, CHM_TR src)
Copy cholmod_triplet to an R_alloc()ed version of it.
Definition: chm_common.c:208
SEXP CHMsuper_validate(SEXP obj)
Definition: chm_common.c:1269
static SEXP supernodal_switchSym
Definition: chm_common.c:19
#define MATRIX_VALID_Csparse
Definition: Mutils.h:384
SEXP CHMsimpl_validate(SEXP obj)
Definition: chm_common.c:1264
Rboolean isValid_Csparse(SEXP x)
"Cheap" C version of Csparse_validate() - not sorting :
Definition: Csparse.c:11
static SEXP dboundSym
Definition: chm_common.c:19
#define CHM_DENSE_TYPE
#define _AS_cholmod_dense_2
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
#define DOFREE_de_MAYBE
static R_INLINE SEXP ALLOC_SLOT(SEXP obj, SEXP nm, SEXPTYPE type, int length)
Allocate an SEXP of given type and length, assign it as slot nm in the object, and return the SEXP...
Definition: Mutils.h:240
CHM_DN numeric_as_chm_dense(CHM_DN ans, double *v, int nr, int nc)
Definition: chm_common.c:1007
static SEXP grow0Sym
Definition: chm_common.c:19
static SEXP grow2Sym
Definition: chm_common.c:19
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
#define CHM_DENSE_COPY_DATA