Matrix  $Rev: 3071 $ at $LastChangedDate: 2015-03-26 15:35:47 +0100 (Thu, 26 Mar 2015) $
cs.h
Go to the documentation of this file.
1 #ifndef _CS_H
2 #define _CS_H
3 #include <stdlib.h>
4 #include <limits.h>
5 #include <math.h>
6 // needed for FILE:
7 #include <stdio.h>
8 #include <stddef.h>
9 // For use with R package 'Matrix'
10 # include <R_ext/Print.h>
11 # include <R_ext/Random.h>
12 # define printf Rprintf
13 #ifdef MATLAB_MEX_FILE
14 #include "mex.h"
15 #endif
16 #define CS_VER 3 /* CSparse Version */
17 #define CS_SUBVER 1
18 #define CS_SUBSUB 2
19 #define CS_DATE "April 16, 2013" /* CSparse release date */
20 #define CS_COPYRIGHT "Copyright (c) Timothy A. Davis, 2006-2013"
21 
22 #ifdef MATLAB_MEX_FILE
23 #undef csi
24 #define csi mwSignedIndex
25 #endif
26 // Matrix pkg:
27 #define csi int
28 #ifndef csi
29 #define csi int
30 #endif
31 
32 /* --- primary CSparse routines and data structures ------------------------- */
33 typedef struct cs_sparse /* matrix in compressed-column or triplet form */
34 {
35  csi nzmax ; /* maximum number of entries */
36  csi m ; /* number of rows */
37  csi n ; /* number of columns */
38  csi *p ; /* column pointers (size n+1) or col indices (size nzmax) */
39  csi *i ; /* row indices, size nzmax */
40  double *x ; /* numerical values, size nzmax */
41  csi nz ; /* # of entries in triplet matrix, -1 for compressed-col */
42 } cs ;
43 
44 cs *cs_add (const cs *A, const cs *B, double alpha, double beta) ;
45 csi cs_cholsol (csi order, const cs *A, double *b) ;
46 cs *cs_compress (const cs *T) ;
47 csi cs_dupl (cs *A) ;
48 csi cs_entry (cs *T, csi i, csi j, double x) ;
49 csi cs_gaxpy (const cs *A, const double *x, double *y) ;
50 cs *cs_load (FILE *f) ;
51 csi cs_lusol (csi order, const cs *A, double *b, double tol) ;
52 cs *cs_multiply (const cs *A, const cs *B) ;
53 double cs_norm (const cs *A) ;
54 csi cs_print (const cs *A, csi brief) ;
55 csi cs_qrsol (csi order, const cs *A, double *b) ;
56 cs *cs_transpose (const cs *A, csi values) ;
57 /* utilities */
58 void *cs_calloc (csi n, size_t size) ;
59 void *cs_free (void *p) ;
60 void *cs_realloc (void *p, csi n, size_t size, csi *ok) ;
61 cs *cs_spalloc (csi m, csi n, csi nzmax, csi values, csi triplet) ;
62 cs *cs_spfree (cs *A) ;
63 csi cs_sprealloc (cs *A, csi nzmax) ;
64 void *cs_malloc (csi n, size_t size) ;
65 
66 /* --- secondary CSparse routines and data structures ----------------------- */
67 typedef struct cs_symbolic /* symbolic Cholesky, LU, or QR analysis */
68 {
69  csi *pinv ; /* inverse row perm. for QR, fill red. perm for Chol */
70  csi *q ; /* fill-reducing column permutation for LU and QR */
71  csi *parent ; /* elimination tree for Cholesky and QR */
72  csi *cp ; /* column pointers for Cholesky, row counts for QR */
73  csi *leftmost ; /* leftmost[i] = min(find(A(i,:))), for QR */
74  csi m2 ; /* # of rows for QR, after adding fictitious rows */
75  double lnz ; /* # entries in L for LU or Cholesky; in V for QR */
76  double unz ; /* # entries in U for LU; in R for QR */
77 } css ;
78 
79 typedef struct cs_numeric /* numeric Cholesky, LU, or QR factorization */
80 {
81  cs *L ; /* L for LU and Cholesky, V for QR */
82  cs *U ; /* U for LU, R for QR, not used for Cholesky */
83  csi *pinv ; /* partial pivoting for LU */
84  double *B ; /* beta [0..n-1] for QR */
85 } csn ;
86 
87 typedef struct cs_dmperm_results /* cs_dmperm or cs_scc output */
88 {
89  csi *p ; /* size m, row permutation */
90  csi *q ; /* size n, column permutation */
91  csi *r ; /* size nb+1, block k is rows r[k] to r[k+1]-1 in A(p,q) */
92  csi *s ; /* size nb+1, block k is cols s[k] to s[k+1]-1 in A(p,q) */
93  csi nb ; /* # of blocks in fine dmperm decomposition */
94  csi rr [5] ; /* coarse row decomposition */
95  csi cc [5] ; /* coarse column decomposition */
96 } csd ;
97 
98 csi *cs_amd (csi order, const cs *A) ;
99 csn *cs_chol (const cs *A, const css *S) ;
100 csd *cs_dmperm (const cs *A, csi seed) ;
101 csi cs_droptol (cs *A, double tol) ;
102 csi cs_dropzeros (cs *A) ;
103 csi cs_happly (const cs *V, csi i, double beta, double *x) ;
104 csi cs_ipvec (const csi *p, const double *b, double *x, csi n) ;
105 csi cs_lsolve (const cs *L, double *x) ;
106 csi cs_ltsolve (const cs *L, double *x) ;
107 csn *cs_lu (const cs *A, const css *S, double tol) ;
108 cs *cs_permute (const cs *A, const csi *pinv, const csi *q, csi values) ;
109 csi *cs_pinv (const csi *p, csi n) ;
110 csi cs_pvec (const csi *p, const double *b, double *x, csi n) ;
111 csn *cs_qr (const cs *A, const css *S) ;
112 css *cs_schol (csi order, const cs *A) ;
113 css *cs_sqr (csi order, const cs *A, csi qr) ;
114 cs *cs_symperm (const cs *A, const csi *pinv, csi values) ;
115 csi cs_updown (cs *L, csi sigma, const cs *C, const csi *parent) ;
116 csi cs_usolve (const cs *U, double *x) ;
117 csi cs_utsolve (const cs *U, double *x) ;
118 /* utilities */
119 css *cs_sfree (css *S) ;
120 csn *cs_nfree (csn *N) ;
121 csd *cs_dfree (csd *D) ;
122 
123 /* --- tertiary CSparse routines -------------------------------------------- */
124 csi *cs_counts (const cs *A, const csi *parent, const csi *post, csi ata) ;
125 double cs_cumsum (csi *p, csi *c, csi n) ;
126 csi cs_dfs (csi j, cs *G, csi top, csi *xi, csi *pstack, const csi *pinv) ;
127 csi cs_ereach (const cs *A, csi k, const csi *parent, csi *s, csi *w) ;
128 csi *cs_etree (const cs *A, csi ata) ;
129 csi cs_fkeep (cs *A, csi (*fkeep) (csi, csi, double, void *), void *other) ;
130 double cs_house (double *x, double *beta, csi n) ;
131 csi cs_leaf (csi i, csi j, const csi *first, csi *maxfirst, csi *prevleaf,
132  csi *ancestor, csi *jleaf) ;
133 csi *cs_maxtrans (const cs *A, csi seed) ;
134 csi *cs_post (const csi *parent, csi n) ;
135 csi *cs_randperm (csi n, csi seed) ;
136 csi cs_reach (cs *G, const cs *B, csi k, csi *xi, const csi *pinv) ;
137 csi cs_scatter (const cs *A, csi j, double beta, csi *w, double *x, csi mark,
138  cs *C, csi nz) ;
139 csd *cs_scc (cs *A) ;
140 csi cs_spsolve (cs *G, const cs *B, csi k, csi *xi, double *x,
141  const csi *pinv, csi lo) ;
142 csi cs_tdfs (csi j, csi k, csi *head, const csi *next, csi *post,
143  csi *stack) ;
144 /* utilities */
145 csd *cs_dalloc (csi m, csi n) ;
146 csd *cs_ddone (csd *D, cs *C, void *w, csi ok) ;
147 cs *cs_done (cs *C, void *w, void *x, csi ok) ;
148 csi *cs_idone (csi *p, cs *C, void *w, csi ok) ;
149 csn *cs_ndone (csn *N, cs *C, void *w, void *x, csi ok) ;
150 
151 #define CS_MAX(a,b) (((a) > (b)) ? (a) : (b))
152 #define CS_MIN(a,b) (((a) < (b)) ? (a) : (b))
153 #define CS_FLIP(i) (-(i)-2)
154 #define CS_UNFLIP(i) (((i) < 0) ? CS_FLIP(i) : (i))
155 #define CS_MARKED(w,j) (w [j] < 0)
156 #define CS_MARK(w,j) { w [j] = CS_FLIP (w [j]) ; }
157 #define CS_CSC(A) (A && (A->nz == -1))
158 #define CS_TRIPLET(A) (A && (A->nz >= 0))
159 #endif
csi cs_happly(const cs *V, csi i, double beta, double *x)
Definition: cs.c:902
void * cs_calloc(csi n, size_t size)
Definition: cs.c:1142
css * cs_sqr(csi order, const cs *A, csi qr)
Definition: cs.c:1740
csi cs_ereach(const cs *A, csi k, const csi *parent, csi *s, csi *w)
Definition: cs.c:811
cs * cs_load(FILE *f)
Definition: cs.c:971
csi cs_pvec(const csi *p, const double *b, double *x, csi n)
Definition: cs.c:1397
csi cs_entry(cs *T, csi i, csi j, double x)
Definition: cs.c:799
cs * cs_transpose(const cs *A, csi values)
Definition: cs.c:1830
csn * cs_chol(const cs *A, const css *S)
Definition: cs.c:393
csi * cs_pinv(const csi *p, csi n)
csi * pinv
Definition: cs.h:83
cs * cs_spalloc(csi m, csi n, csi nzmax, csi values, csi triplet)
Definition: cs.c:1907
csi cs_print(const cs *A, csi brief)
Definition: cs.c:1359
csi n
Definition: cs.h:37
csn * cs_ndone(csn *N, cs *C, void *w, void *x, csi ok)
Definition: cs.c:2009
struct cs_dmperm_results csd
csi cs_ipvec(const csi *p, const double *b, double *x, csi n)
Definition: cs.c:942
csi * cp
Definition: cs.h:72
csi * cs_counts(const cs *A, const csi *parent, const csi *post, csi ata)
Definition: cs.c:510
void * cs_realloc(void *p, csi n, size_t size, csi *ok)
Definition: cs.c:1155
double * B
Definition: cs.h:84
csn * cs_nfree(csn *N)
Definition: cs.c:1946
cs * cs_spfree(cs *A)
Definition: cs.c:1936
csi * cs_maxtrans(const cs *A, csi seed)
Definition: cs.c:1204
csi m2
Definition: cs.h:74
csd * cs_dalloc(csi m, csi n)
Definition: cs.c:1969
csi m
Definition: cs.h:36
csi * q
Definition: cs.h:90
csn * cs_qr(const cs *A, const css *S)
Definition: cs.c:1405
double cs_cumsum(csi *p, csi *c, csi n)
Definition: cs.c:556
double cs_house(double *x, double *beta, csi n)
Definition: cs.c:921
csd * cs_ddone(csd *D, cs *C, void *w, csi ok)
Definition: cs.c:2018
csi * cs_etree(const cs *A, csi ata)
Definition: cs.c:833
csi * p
Definition: cs.h:38
cs * U
Definition: cs.h:82
csi cs_gaxpy(const cs *A, const double *x, double *y)
Definition: cs.c:886
void * cs_free(void *p)
Definition: cs.c:1148
csi cs_droptol(cs *A, double tol)
Definition: cs.c:753
csi * i
Definition: cs.h:39
struct cs_numeric csn
cs * cs_done(cs *C, void *w, void *x, csi ok)
Definition: cs.c:1993
double lnz
Definition: cs.h:75
Definition: cs.h:79
cs * cs_symperm(const cs *A, const csi *pinv, csi values)
Definition: cs.c:1769
csi * parent
Definition: cs.h:71
csi cs_qrsol(csi order, const cs *A, double *b)
Definition: cs.c:1477
csi * leftmost
Definition: cs.h:73
double unz
Definition: cs.h:76
csd * cs_dfree(csd *D)
Definition: cs.c:1982
csi rr[5]
Definition: cs.h:94
struct cs_sparse cs
struct cs_symbolic css
csi cc[5]
Definition: cs.h:95
csi cs_lusol(csi order, const cs *A, double *b, double tol)
Definition: cs.c:1104
csi * cs_randperm(csi n, csi seed)
Definition: cs.c:1531
csi * pinv
Definition: cs.h:69
cs * cs_add(const cs *A, const cs *B, double alpha, double beta)
Definition: cs.c:3
double cs_norm(const cs *A)
Definition: cs.c:1288
csi cs_cholsol(csi order, const cs *A, double *b)
Definition: cs.c:450
csi cs_leaf(csi i, csi j, const csi *first, csi *maxfirst, csi *prevleaf, csi *ancestor, csi *jleaf)
Definition: cs.c:950
csi nzmax
Definition: cs.h:35
csi cs_dropzeros(cs *A)
Definition: cs.c:761
csi * p
Definition: cs.h:89
Definition: cs.h:33
csi cs_lsolve(const cs *L, double *x)
Definition: cs.c:985
csi * q
Definition: cs.h:70
css * cs_schol(csi order, const cs *A)
Definition: cs.c:1631
csi * cs_post(const csi *parent, csi n)
Definition: cs.c:1336
cs * cs_multiply(const cs *A, const cs *B)
Definition: cs.c:1254
csi cs_reach(cs *G, const cs *B, csi k, csi *xi, const csi *pinv)
Definition: cs.c:1553
void * cs_malloc(csi n, size_t size)
Definition: cs.c:1136
cs * cs_compress(const cs *T)
Definition: cs.c:475
csd * cs_scc(cs *A)
Definition: cs.c:1591
#define csi
Definition: cs.h:27
csi * cs_amd(csi order, const cs *A)
Definition: cs.c:45
csi cs_dupl(cs *A)
Definition: cs.c:766
csi nz
Definition: cs.h:41
csi cs_fkeep(cs *A, csi(*fkeep)(csi, csi, double, void *), void *other)
Definition: cs.c:862
csi cs_usolve(const cs *U, double *x)
Definition: cs.c:1890
cholmod_common c
Definition: chm_common.c:15
cs * cs_permute(const cs *A, const csi *pinv, const csi *q, csi values)
Definition: cs.c:1302
csi cs_spsolve(cs *G, const cs *B, csi k, csi *xi, double *x, const csi *pinv, csi lo)
Definition: cs.c:1656
double * x
Definition: cs.h:40
css * cs_sfree(css *S)
Definition: cs.c:1957
csi cs_updown(cs *L, csi sigma, const cs *C, const csi *parent)
Definition: cs.c:1854
csi * r
Definition: cs.h:91
csn * cs_lu(const cs *A, const css *S, double tol)
Definition: cs.c:1019
csi * s
Definition: cs.h:92
csi * cs_idone(csi *p, cs *C, void *w, csi ok)
Definition: cs.c:2001
csi cs_ltsolve(const cs *L, double *x)
Definition: cs.c:1002
csd * cs_dmperm(const cs *A, csi seed)
Definition: cs.c:672
csi cs_sprealloc(cs *A, csi nzmax)
Definition: cs.c:1922
cs * L
Definition: cs.h:81
Definition: cs.h:67
csi cs_dfs(csi j, cs *G, csi top, csi *xi, csi *pstack, const csi *pinv)
Definition: cs.c:572
csi cs_utsolve(const cs *U, double *x)
Definition: cs.c:2025
csi cs_scatter(const cs *A, csi j, double beta, csi *w, double *x, csi mark, cs *C, csi nz)
Definition: cs.c:1570
csi cs_tdfs(csi j, csi k, csi *head, const csi *next, csi *post, csi *stack)
Definition: cs.c:1807