1 /* ========================================================================== */
2 /* === Modify/t_cholmod_updown ============================================== */
3 /* ========================================================================== */
4
5 /* -----------------------------------------------------------------------------
6 * CHOLMOD/Modify Module. Copyright (C) 2005-2006,
7 * Timothy A. Davis and William W. Hager.
8 * http://www.suitesparse.com
9 * -------------------------------------------------------------------------- */
10
11 /* Updates/downdates the LDL' factorization, by computing a new factorization of
12 *
13 * Lnew * Dnew * Lnew' = Lold * Dold * Lold' +/- C*C'
14 *
15 * This file is not compiled separately. It is included into
16 * cholmod_updown.c. There are no user-callable routines in this file.
17 *
18 * The next include statements, below, create the numerical update/downdate
19 * kernels from t_cholmod_updown_numkr.c. There are 4 compiled versions of this
20 * file, one for each value of WDIM in the set 1, 2, 4, and 8. Each calls
21 * multiple versions of t_cholmod_updown_numkr; the number of versions of each
22 * is equal to WDIM. Each t_cholmod_updown_numkr version is included as a
23 * static function within its t_cholmod_updown.c caller routine. Thus:
24 *
25 * t*_updown.c creates these versions of t_cholmod_updown_numkr.c:
26 * --------- ---------------------------------------------------
27 *
28 * updown_1_r updown_1_1
29 *
30 * updown_2_r updown_2_1 updown_2_2
31 *
32 * updown_4_r updown_4_1 updown_4_2 updown_4_3 updown_4_4
33 *
34 * updown_8_r updown_8_1 updown_8_2 updown_8_3 updown_8_4
35 * updown_8_5 updown_8_6 updown_8_7 updown_8_8
36 *
37 * workspace: Xwork (nrow*wdim)
38 */
39
40 /* ========================================================================== */
41 /* === routines for numeric update/downdate along one path ================== */
42 /* ========================================================================== */
43
44 #undef FORM_NAME
45 #undef NUMERIC
46
47 #define FORM_NAME(k,rank) updown_ ## k ## _ ## rank
48 #define NUMERIC(k,rank) FORM_NAME(k,rank)
49
50 #define RANK 1
51 #include "t_cholmod_updown_numkr.c"
52
53 #if WDIM >= 2
54 #define RANK 2
55 #include "t_cholmod_updown_numkr.c"
56 #endif
57
58 #if WDIM >= 4
59 #define RANK 3
60 #include "t_cholmod_updown_numkr.c"
61 #define RANK 4
62 #include "t_cholmod_updown_numkr.c"
63 #endif
64
65 #if WDIM == 8
66 #define RANK 5
67 #include "t_cholmod_updown_numkr.c"
68 #define RANK 6
69 #include "t_cholmod_updown_numkr.c"
70 #define RANK 7
71 #include "t_cholmod_updown_numkr.c"
72 #define RANK 8
73 #include "t_cholmod_updown_numkr.c"
74 #endif
75
76
77 /* ========================================================================== */
78 /* === numeric update/downdate for all paths ================================ */
79 /* ========================================================================== */
80
NUMERIC(WDIM,r)81 static void NUMERIC (WDIM, r)
82 (
83 int update, /* TRUE for update, FALSE for downdate */
84 cholmod_sparse *C, /* in packed or unpacked, and sorted form */
85 /* no empty columns */
86 Int rank, /* rank of the update/downdate */
87 cholmod_factor *L, /* with unit diagonal (diagonal not stored) */
88 /* temporary workspaces: */
89 double W [ ], /* n-by-WDIM dense matrix, initially zero */
90 Path_type Path [ ],
91 Int npaths,
92 Int mask [ ], /* size n */
93 Int maskmark,
94 cholmod_common *Common
95 )
96 {
97 double Alpha [8] ;
98 double *Cx, *Wpath, *W1, *a ;
99 Int i, j, p, ccol, pend, wfirst, e, path, packed ;
100 Int *Ci, *Cp, *Cnz ;
101
102 /* ---------------------------------------------------------------------- */
103 /* get inputs */
104 /* ---------------------------------------------------------------------- */
105
106 Ci = C->i ;
107 Cx = C->x ;
108 Cp = C->p ;
109 Cnz = C->nz ;
110 packed = C->packed ;
111 ASSERT (IMPLIES (!packed, Cnz != NULL)) ;
112 ASSERT (L->n == C->nrow) ;
113 DEBUG (CHOLMOD(dump_real) ("num_d: in W:", W, WDIM, L->n, FALSE, 1,Common));
114
115 /* ---------------------------------------------------------------------- */
116 /* scatter C into W */
117 /* ---------------------------------------------------------------------- */
118
119 for (path = 0 ; path < rank ; path++)
120 {
121 /* W (:, path) = C (:, Path [path].col) */
122 ccol = Path [path].ccol ;
123 Wpath = W + path ;
124 PRINT1 (("Ordered Columns [path = "ID"] = "ID"\n", path, ccol)) ;
125 p = Cp [ccol] ;
126 pend = (packed) ? (Cp [ccol+1]) : (p + Cnz [ccol]) ;
127 /* column C can be empty */
128 for ( ; p < pend ; p++)
129 {
130 i = Ci [p] ;
131 ASSERT (i >= 0 && i < (Int) (C->nrow)) ;
132 if (mask == NULL || mask [i] < maskmark)
133 {
134 Wpath [WDIM * i] = Cx [p] ;
135 }
136 PRINT1 ((" row "ID" : %g mask "ID"\n", i, Cx [p],
137 (mask) ? mask [i] : 0)) ;
138 }
139 Alpha [path] = 1.0 ;
140 }
141 DEBUG (CHOLMOD(dump_real) ("num_d: W:", W, WDIM, L->n, FALSE, 1,Common)) ;
142
143 /* ---------------------------------------------------------------------- */
144 /* numeric update/downdate of the paths */
145 /* ---------------------------------------------------------------------- */
146
147 /* for each disjoint subpath in Tbar in DFS order do */
148 for (path = rank ; path < npaths ; path++)
149 {
150
151 /* determine which columns of W to use */
152 wfirst = Path [path].wfirst ;
153 e = Path [path].end ;
154 j = Path [path].start ;
155 ASSERT (e >= 0 && e < (Int) (L->n)) ;
156 ASSERT (j >= 0 && j < (Int) (L->n)) ;
157
158 W1 = W + wfirst ; /* pointer to row 0, column wfirst of W */
159 a = Alpha + wfirst ; /* pointer to Alpha [wfirst] */
160
161 PRINT1 (("Numerical update/downdate of path "ID"\n", path)) ;
162 PRINT1 (("start "ID" end "ID" wfirst "ID" rank "ID" ccol "ID"\n", j, e,
163 wfirst, Path [path].rank, Path [path].ccol)) ;
164
165 #if WDIM == 1
166 NUMERIC (WDIM,1) (update, j, e, a, W1, L, Common) ;
167 #else
168
169 switch (Path [path].rank)
170 {
171 case 1:
172 NUMERIC (WDIM,1) (update, j, e, a, W1, L, Common) ;
173 break ;
174
175 #if WDIM >= 2
176 case 2:
177 NUMERIC (WDIM,2) (update, j, e, a, W1, L, Common) ;
178 break ;
179 #endif
180
181 #if WDIM >= 4
182 case 3:
183 NUMERIC (WDIM,3) (update, j, e, a, W1, L, Common) ;
184 break ;
185 case 4:
186 NUMERIC (WDIM,4) (update, j, e, a, W1, L, Common) ;
187 break ;
188 #endif
189
190 #if WDIM == 8
191 case 5:
192 NUMERIC (WDIM,5) (update, j, e, a, W1, L, Common) ;
193 break ;
194 case 6:
195 NUMERIC (WDIM,6) (update, j, e, a, W1, L, Common) ;
196 break ;
197 case 7:
198 NUMERIC (WDIM,7) (update, j, e, a, W1, L, Common) ;
199 break ;
200 case 8:
201 NUMERIC (WDIM,8) (update, j, e, a, W1, L, Common) ;
202 break ;
203 #endif
204
205 }
206 #endif
207
208 }
209 }
210
211 /* prepare for the next inclusion of this file in cholmod_updown.c */
212 #undef WDIM
213