1 /* ========================================================================== */
2 /* === Core/cholmod_factor ================================================== */
3 /* ========================================================================== */
4
5 /* -----------------------------------------------------------------------------
6 * CHOLMOD/Core Module. Copyright (C) 2005-2013,
7 * Univ. of Florida. Author: Timothy A. Davis
8 * -------------------------------------------------------------------------- */
9
10 /* Core utility routines for the cholmod_factor object:
11 *
12 * The data structure for an LL' or LDL' factorization is too complex to
13 * describe in one sentence. This object can hold the symbolic analysis alone,
14 * or in combination with a "simplicial" (similar to a sparse matrix) or
15 * "supernodal" form of the numerical factorization. Only the routine to free
16 * a factor is primary, since a factor object is created by the factorization
17 * routine (cholmod_factorize). It must be freed with cholmod_free_factor.
18 *
19 * Primary routine:
20 * ----------------
21 * cholmod_free_factor free a factor
22 *
23 * Secondary routines:
24 * -------------------
25 * cholmod_allocate_factor allocate a symbolic factor (LL' or LDL')
26 * cholmod_reallocate_factor change the # entries in a factor
27 * cholmod_change_factor change the type of factor (e.g., LDL' to LL')
28 * cholmod_pack_factor pack the columns of a factor
29 * cholmod_reallocate_column resize a single column of a factor
30 * cholmod_factor_to_sparse create a sparse matrix copy of a factor
31 * cholmod_copy_factor create a copy of a factor
32 *
33 * Note that there is no cholmod_sparse_to_factor routine to create a factor
34 * as a copy of a sparse matrix. It could be done, after a fashion, but a
35 * lower triangular sparse matrix would not necessarily have a chordal graph,
36 * which would break the many CHOLMOD routines that rely on this property.
37 *
38 * The cholmod_factor_to_sparse routine is provided so that matrix operations
39 * in the MatrixOps module may be applied to L. Those operations operate on
40 * cholmod_sparse objects, and they are not guaranteed to maintain the chordal
41 * property of L. Such a modified L cannot be safely converted back to a
42 * cholmod_factor object.
43 */
44
45 #include "cholmod_internal.h"
46 #include "cholmod_core.h"
47
48
49 /* ========================================================================== */
50 /* === cholmod_allocate_factor ============================================== */
51 /* ========================================================================== */
52
53 /* Allocate a simplicial symbolic factor, with L->Perm and L->ColCount allocated
54 * and initialized to "empty" values (Perm [k] = k, and ColCount[k] = 1).
55 * The integer and numerical parts of L are not allocated. L->xtype is returned
56 * as CHOLMOD_PATTERN and L->is_super are returned as FALSE. L->is_ll is also
57 * returned FALSE, but this may be modified when the matrix is factorized.
58 *
59 * This is sufficient (but far from ideal) for input to cholmod_factorize,
60 * since the simplicial LL' or LDL' factorization (cholmod_rowfac) can
61 * reallocate the columns of L as needed. The primary purpose of this routine
62 * is to allocate space for a symbolic factorization, for the "expert" user to
63 * do his or her own symbolic analysis. The typical user should use
64 * cholmod_analyze instead of this routine.
65 *
66 * workspace: none
67 */
68
CHOLMOD(allocate_factor)69 cholmod_factor *CHOLMOD(allocate_factor)
70 (
71 /* ---- input ---- */
72 size_t n, /* L is n-by-n */
73 /* --------------- */
74 cholmod_common *Common
75 )
76 {
77 Int j ;
78 Int *Perm, *ColCount ;
79 cholmod_factor *L ;
80 int ok = TRUE ;
81
82 RETURN_IF_NULL_COMMON (FALSE) ;
83 Common->status = CHOLMOD_OK ;
84
85 /* ensure the dimension does not cause integer overflow */
86 (void) CHOLMOD(add_size_t) (n, 2, &ok) ;
87 if (!ok || n > Int_max)
88 {
89 ERROR (CHOLMOD_TOO_LARGE, "problem too large") ;
90 return (NULL) ;
91 }
92
93 L = CHOLMOD(malloc) (sizeof (cholmod_factor), 1, Common) ;
94 if (Common->status < CHOLMOD_OK)
95 {
96 return (NULL) ; /* out of memory */
97 }
98 L->n = n ;
99 L->is_ll = FALSE ;
100 L->is_super = FALSE ;
101 L->is_monotonic = TRUE ;
102 L->itype = ITYPE ;
103 L->xtype = CHOLMOD_PATTERN ;
104 L->dtype = DTYPE ;
105
106 /* allocate the purely symbolic part of L */
107 L->ordering = CHOLMOD_NATURAL ;
108 L->Perm = CHOLMOD(malloc) (n, sizeof (Int), Common) ;
109 L->IPerm = NULL ; /* only created by cholmod_solve2 when needed */
110 L->ColCount = CHOLMOD(malloc) (n, sizeof (Int), Common) ;
111
112 /* simplicial part of L is empty */
113 L->nzmax = 0 ;
114 L->p = NULL ;
115 L->i = NULL ;
116 L->x = NULL ;
117 L->z = NULL ;
118 L->nz = NULL ;
119 L->next = NULL ;
120 L->prev = NULL ;
121
122 /* supernodal part of L is also empty */
123 L->nsuper = 0 ;
124 L->ssize = 0 ;
125 L->xsize = 0 ;
126 L->maxesize = 0 ;
127 L->maxcsize = 0 ;
128 L->super = NULL ;
129 L->pi = NULL ;
130 L->px = NULL ;
131 L->s = NULL ;
132 L->useGPU = 0;
133
134 /* L has not been factorized */
135 L->minor = n ;
136
137 if (Common->status < CHOLMOD_OK)
138 {
139 CHOLMOD(free_factor) (&L, Common) ;
140 return (NULL) ; /* out of memory */
141 }
142
143 /* initialize Perm and ColCount */
144 Perm = L->Perm ;
145 for (j = 0 ; j < ((Int) n) ; j++)
146 {
147 Perm [j] = j ;
148 }
149 ColCount = L->ColCount ;
150 for (j = 0 ; j < ((Int) n) ; j++)
151 {
152 ColCount [j] = 1 ;
153 }
154
155 return (L) ;
156 }
157
158
159 /* ========================================================================== */
160 /* === cholmod_free_factor ================================================== */
161 /* ========================================================================== */
162
163 /* Free a factor object.
164 *
165 * workspace: none
166 */
167
CHOLMOD(free_factor)168 int CHOLMOD(free_factor)
169 (
170 /* ---- in/out --- */
171 cholmod_factor **LHandle, /* factor to free, NULL on output */
172 /* --------------- */
173 cholmod_common *Common
174 )
175 {
176 Int n, lnz, xs, ss, s ;
177 cholmod_factor *L ;
178
179 RETURN_IF_NULL_COMMON (FALSE) ;
180
181 if (LHandle == NULL)
182 {
183 /* nothing to do */
184 return (TRUE) ;
185 }
186 L = *LHandle ;
187 if (L == NULL)
188 {
189 /* nothing to do */
190 return (TRUE) ;
191 }
192
193 n = L->n ;
194 lnz = L->nzmax ;
195 s = L->nsuper + 1 ;
196 xs = (L->is_super) ? ((Int) (L->xsize)) : (lnz) ;
197 ss = L->ssize ;
198
199 /* symbolic part of L */
200 CHOLMOD(free) (n, sizeof (Int), L->Perm, Common) ;
201 CHOLMOD(free) (n, sizeof (Int), L->IPerm, Common) ;
202 CHOLMOD(free) (n, sizeof (Int), L->ColCount, Common) ;
203
204 /* simplicial form of L */
205 CHOLMOD(free) (n+1, sizeof (Int), L->p, Common) ;
206 CHOLMOD(free) (lnz, sizeof (Int), L->i, Common) ;
207 CHOLMOD(free) (n, sizeof (Int), L->nz, Common) ;
208 CHOLMOD(free) (n+2, sizeof (Int), L->next, Common) ;
209 CHOLMOD(free) (n+2, sizeof (Int), L->prev, Common) ;
210
211 /* supernodal form of L */
212 CHOLMOD(free) (s, sizeof (Int), L->pi, Common) ;
213 CHOLMOD(free) (s, sizeof (Int), L->px, Common) ;
214 CHOLMOD(free) (s, sizeof (Int), L->super, Common) ;
215 CHOLMOD(free) (ss, sizeof (Int), L->s, Common) ;
216
217 /* numerical values for both simplicial and supernodal L */
218 if (L->xtype == CHOLMOD_REAL)
219 {
220 CHOLMOD(free) (xs, sizeof (double), L->x, Common) ;
221 }
222 else if (L->xtype == CHOLMOD_COMPLEX)
223 {
224 CHOLMOD(free) (xs, 2*sizeof (double), L->x, Common) ;
225 }
226 else if (L->xtype == CHOLMOD_ZOMPLEX)
227 {
228 CHOLMOD(free) (xs, sizeof (double), L->x, Common) ;
229 CHOLMOD(free) (xs, sizeof (double), L->z, Common) ;
230 }
231
232 *LHandle = CHOLMOD(free) (1, sizeof (cholmod_factor), (*LHandle), Common) ;
233 return (TRUE) ;
234 }
235
236
237 /* ========================================================================== */
238 /* === cholmod_reallocate_factor ============================================ */
239 /* ========================================================================== */
240
241 /* Change the size of L->i and L->x, or allocate them if their current size
242 * is zero. L must be simplicial.
243 *
244 * workspace: none
245 */
246
CHOLMOD(reallocate_factor)247 int CHOLMOD(reallocate_factor)
248 (
249 /* ---- input ---- */
250 size_t nznew, /* new # of entries in L */
251 /* ---- in/out --- */
252 cholmod_factor *L, /* factor to modify */
253 /* --------------- */
254 cholmod_common *Common
255 )
256 {
257 /* ---------------------------------------------------------------------- */
258 /* get inputs */
259 /* ---------------------------------------------------------------------- */
260
261 RETURN_IF_NULL_COMMON (FALSE) ;
262 RETURN_IF_NULL (L, FALSE) ;
263 RETURN_IF_XTYPE_INVALID (L, CHOLMOD_REAL, CHOLMOD_ZOMPLEX, FALSE) ;
264 PRINT1 (("realloc factor: xtype %d\n", L->xtype)) ;
265 if (L->is_super)
266 {
267 /* L must be simplicial, and not symbolic */
268 ERROR (CHOLMOD_INVALID, "L invalid") ;
269 return (FALSE) ;
270 }
271 Common->status = CHOLMOD_OK ;
272 PRINT1 (("realloc factor %g to %g\n", (double) L->nzmax, (double) nznew)) ;
273
274 /* ---------------------------------------------------------------------- */
275 /* resize (or allocate) the L->i and L->x components of the factor */
276 /* ---------------------------------------------------------------------- */
277
278 CHOLMOD(realloc_multiple) (nznew, 1, L->xtype, &(L->i), NULL,
279 &(L->x), &(L->z), &(L->nzmax), Common) ;
280 return (Common->status == CHOLMOD_OK) ;
281 }
282
283
284 /* ========================================================================== */
285 /* === cholmod_reallocate_column =========================================== */
286 /* ========================================================================== */
287
288 /* Column j needs more space, reallocate it at the end of L->i and L->x.
289 * If the reallocation fails, the factor is converted to a simplicial
290 * symbolic factor (no pattern, just L->Perm and L->ColCount).
291 *
292 * workspace: none
293 */
294
CHOLMOD(reallocate_column)295 int CHOLMOD(reallocate_column)
296 (
297 /* ---- input ---- */
298 size_t j, /* the column to reallocate */
299 size_t need, /* required size of column j */
300 /* ---- in/out --- */
301 cholmod_factor *L, /* factor to modify */
302 /* --------------- */
303 cholmod_common *Common
304 )
305 {
306 double xneed ;
307 double *Lx, *Lz ;
308 Int *Lp, *Lprev, *Lnext, *Li, *Lnz ;
309 Int n, pold, pnew, len, k, tail ;
310
311 /* ---------------------------------------------------------------------- */
312 /* get inputs */
313 /* ---------------------------------------------------------------------- */
314
315 RETURN_IF_NULL_COMMON (FALSE) ;
316 RETURN_IF_NULL (L, FALSE) ;
317 RETURN_IF_XTYPE_INVALID (L, CHOLMOD_REAL, CHOLMOD_ZOMPLEX, FALSE) ;
318 if (L->is_super)
319 {
320 ERROR (CHOLMOD_INVALID, "L must be simplicial") ;
321 return (FALSE) ;
322 }
323 n = L->n ;
324 if (j >= L->n || need == 0)
325 {
326 ERROR (CHOLMOD_INVALID, "j invalid") ;
327 return (FALSE) ; /* j out of range */
328 }
329 Common->status = CHOLMOD_OK ;
330
331 DEBUG (CHOLMOD(dump_factor) (L, "start colrealloc", Common)) ;
332
333 /* ---------------------------------------------------------------------- */
334 /* increase the size of L if needed */
335 /* ---------------------------------------------------------------------- */
336
337 /* head = n+1 ; */
338 tail = n ;
339 Lp = L->p ;
340 Lnz = L->nz ;
341 Lprev = L->prev ;
342 Lnext = L->next ;
343
344 ASSERT (Lnz != NULL) ;
345 ASSERT (Lnext != NULL && Lprev != NULL) ;
346 PRINT1 (("col %g need %g\n", (double) j, (double) need)) ;
347
348 /* column j cannot have more than n-j entries if all entries are present */
349 need = MIN (need, n-j) ;
350
351 /* compute need in double to avoid integer overflow */
352 if (Common->grow1 >= 1.0)
353 {
354 xneed = (double) need ;
355 xneed = Common->grow1 * xneed + Common->grow2 ;
356 xneed = MIN (xneed, n-j) ;
357 need = (Int) xneed ;
358 }
359 PRINT1 (("really new need %g current %g\n", (double) need,
360 (double) (Lp [Lnext [j]] - Lp [j]))) ;
361 ASSERT (need >= 1 && need <= n-j) ;
362
363 if (Lp [Lnext [j]] - Lp [j] >= (Int) need)
364 {
365 /* no need to reallocate the column, it's already big enough */
366 PRINT1 (("colrealloc: quick return %g %g\n",
367 (double) (Lp [Lnext [j]] - Lp [j]), (double) need)) ;
368 return (TRUE) ;
369
370 }
371
372 if (Lp [tail] + need > L->nzmax)
373 {
374 /* use double to avoid integer overflow */
375 xneed = (double) need ;
376 if (Common->grow0 < 1.2) /* fl. pt. compare, false if NaN */
377 {
378 /* if grow0 is less than 1.2 or NaN, don't use it */
379 xneed = 1.2 * (((double) L->nzmax) + xneed + 1) ;
380 }
381 else
382 {
383 xneed = Common->grow0 * (((double) L->nzmax) + xneed + 1) ;
384 }
385 if (xneed > Size_max ||
386 !CHOLMOD(reallocate_factor) ((Int) xneed, L, Common))
387 {
388 /* out of memory, convert to simplicial symbolic */
389 CHOLMOD(change_factor) (CHOLMOD_PATTERN, L->is_ll, FALSE, TRUE,
390 TRUE, L, Common) ;
391 ERROR (CHOLMOD_OUT_OF_MEMORY, "out of memory; L now symbolic") ;
392 return (FALSE) ; /* out of memory */
393 }
394 PRINT1 (("\n=== GROW L from %g to %g\n",
395 (double) L->nzmax, (double) xneed)) ;
396 /* pack all columns so that each column has at most grow2 free space */
397 CHOLMOD(pack_factor) (L, Common) ;
398 ASSERT (Common->status == CHOLMOD_OK) ;
399 Common->nrealloc_factor++ ;
400 }
401
402 /* ---------------------------------------------------------------------- */
403 /* reallocate the column */
404 /* ---------------------------------------------------------------------- */
405
406 Common->nrealloc_col++ ;
407
408 Li = L->i ;
409 Lx = L->x ;
410 Lz = L->z ;
411
412 /* remove j from its current position in the list */
413 Lnext [Lprev [j]] = Lnext [j] ;
414 Lprev [Lnext [j]] = Lprev [j] ;
415
416 /* place j at the end of the list */
417 Lnext [Lprev [tail]] = j ;
418 Lprev [j] = Lprev [tail] ;
419 Lnext [j] = n ;
420 Lprev [tail] = j ;
421
422 /* L is no longer monotonic; columns are out-of-order */
423 L->is_monotonic = FALSE ;
424
425 /* allocate space for column j */
426 pold = Lp [j] ;
427 pnew = Lp [tail] ;
428 Lp [j] = pnew ;
429 Lp [tail] += need ;
430
431 /* copy column j to the new space */
432 len = Lnz [j] ;
433 for (k = 0 ; k < len ; k++)
434 {
435 Li [pnew + k] = Li [pold + k] ;
436 }
437
438 if (L->xtype == CHOLMOD_REAL)
439 {
440 for (k = 0 ; k < len ; k++)
441 {
442 Lx [pnew + k] = Lx [pold + k] ;
443 }
444 }
445 else if (L->xtype == CHOLMOD_COMPLEX)
446 {
447 for (k = 0 ; k < len ; k++)
448 {
449 Lx [2*(pnew + k) ] = Lx [2*(pold + k) ] ;
450 Lx [2*(pnew + k)+1] = Lx [2*(pold + k)+1] ;
451 }
452 }
453 else if (L->xtype == CHOLMOD_ZOMPLEX)
454 {
455 for (k = 0 ; k < len ; k++)
456 {
457 Lx [pnew + k] = Lx [pold + k] ;
458 Lz [pnew + k] = Lz [pold + k] ;
459 }
460 }
461
462 DEBUG (CHOLMOD(dump_factor) (L, "colrealloc done", Common)) ;
463
464 /* successful reallocation of column j of L */
465 return (TRUE) ;
466 }
467
468
469 /* ========================================================================== */
470 /* === cholmod_pack_factor ================================================== */
471 /* ========================================================================== */
472
473 /* Pack the columns of a simplicial LDL' or LL' factor. This can be followed
474 * by a call to cholmod_reallocate_factor to reduce the size of L to the exact
475 * size required by the factor, if desired. Alternatively, you can leave the
476 * size of L->i and L->x the same, to allow space for future updates/rowadds.
477 *
478 * Each column is reduced in size so that it has at most Common->grow2 free
479 * space at the end of the column.
480 *
481 * Does nothing and returns silently if given any other type of factor.
482 *
483 * Does NOT force the columns of L to be monotonic. It thus differs from
484 * cholmod_change_factor (xtype, -, FALSE, TRUE, TRUE, L, Common), which
485 * packs the columns and ensures that they appear in monotonic order.
486 */
487
CHOLMOD(pack_factor)488 int CHOLMOD(pack_factor)
489 (
490 /* ---- in/out --- */
491 cholmod_factor *L, /* factor to modify */
492 /* --------------- */
493 cholmod_common *Common
494 )
495 {
496 double *Lx, *Lz ;
497 Int *Lp, *Li, *Lnz, *Lnext ;
498 Int pnew, j, k, pold, len, n, head, tail, grow2 ;
499
500 /* ---------------------------------------------------------------------- */
501 /* get inputs */
502 /* ---------------------------------------------------------------------- */
503
504 RETURN_IF_NULL_COMMON (FALSE) ;
505 RETURN_IF_NULL (L, FALSE) ;
506 RETURN_IF_XTYPE_INVALID (L, CHOLMOD_PATTERN, CHOLMOD_ZOMPLEX, FALSE) ;
507 Common->status = CHOLMOD_OK ;
508 DEBUG (CHOLMOD(dump_factor) (L, "start pack", Common)) ;
509 PRINT1 (("PACK factor %d\n", L->is_super)) ;
510
511 if (L->xtype == CHOLMOD_PATTERN || L->is_super)
512 {
513 /* nothing to do unless L is simplicial numeric */
514 return (TRUE) ;
515 }
516
517 /* ---------------------------------------------------------------------- */
518 /* pack */
519 /* ---------------------------------------------------------------------- */
520
521 grow2 = Common->grow2 ;
522 PRINT1 (("\nPACK grow2 "ID"\n", grow2)) ;
523
524 pnew = 0 ;
525 n = L->n ;
526 Lp = L->p ;
527 Li = L->i ;
528 Lx = L->x ;
529 Lz = L->z ;
530 Lnz = L->nz ;
531 Lnext = L->next ;
532
533 head = n+1 ;
534 tail = n ;
535
536 for (j = Lnext [head] ; j != tail ; j = Lnext [j])
537 {
538 /* pack column j */
539 pold = Lp [j] ;
540 len = Lnz [j] ;
541 ASSERT (len > 0) ;
542 PRINT2 (("col "ID" pnew "ID" pold "ID"\n", j, pnew, pold)) ;
543 if (pnew < pold)
544 {
545 PRINT2 ((" pack this column\n")) ;
546
547 for (k = 0 ; k < len ; k++)
548 {
549 Li [pnew + k] = Li [pold + k] ;
550 }
551
552 if (L->xtype == CHOLMOD_REAL)
553 {
554 for (k = 0 ; k < len ; k++)
555 {
556 Lx [pnew + k] = Lx [pold + k] ;
557 }
558 }
559 else if (L->xtype == CHOLMOD_COMPLEX)
560 {
561 for (k = 0 ; k < len ; k++)
562 {
563 Lx [2*(pnew + k) ] = Lx [2*(pold + k) ] ;
564 Lx [2*(pnew + k)+1] = Lx [2*(pold + k)+1] ;
565 }
566 }
567 else if (L->xtype == CHOLMOD_ZOMPLEX)
568 {
569 for (k = 0 ; k < len ; k++)
570 {
571 Lx [pnew + k] = Lx [pold + k] ;
572 Lz [pnew + k] = Lz [pold + k] ;
573 }
574 }
575
576 Lp [j] = pnew ;
577 }
578 len = MIN (len + grow2, n - j) ;
579 pnew = MIN (Lp [j] + len, Lp [Lnext [j]]) ;
580 }
581 PRINT2 (("final pnew = "ID"\n", pnew)) ;
582 return (TRUE) ;
583 }
584
585
586 /* ========================================================================== */
587 /* === cholmod_factor_to_sparse ============================================= */
588 /* ========================================================================== */
589
590 /* Constructs a column-oriented sparse matrix containing the pattern and values
591 * of a simplicial or supernodal numerical factor, and then converts the factor
592 * into a simplicial symbolic factor. If L is already packed, monotonic,
593 * and simplicial (which is the case when cholmod_factorize uses the simplicial
594 * Cholesky factorization algorithm) then this routine requires only O(1)
595 * memory and takes O(1) time.
596 *
597 * Only operates on numeric factors (real, complex, or zomplex). Does not
598 * change the numeric L->xtype (the resulting sparse matrix has the same xtype
599 * as L). If this routine fails, L is left unmodified.
600 */
601
CHOLMOD(factor_to_sparse)602 cholmod_sparse *CHOLMOD(factor_to_sparse)
603 (
604 /* ---- in/out --- */
605 cholmod_factor *L, /* factor to copy, converted to symbolic on output */
606 /* --------------- */
607 cholmod_common *Common
608 )
609 {
610 cholmod_sparse *Lsparse ;
611
612 /* ---------------------------------------------------------------------- */
613 /* get inputs */
614 /* ---------------------------------------------------------------------- */
615
616 RETURN_IF_NULL_COMMON (NULL) ;
617 RETURN_IF_NULL (L, NULL) ;
618 RETURN_IF_XTYPE_INVALID (L, CHOLMOD_REAL, CHOLMOD_ZOMPLEX, NULL) ;
619 Common->status = CHOLMOD_OK ;
620 DEBUG (CHOLMOD(dump_factor) (L, "start convert to matrix", Common)) ;
621
622 /* ---------------------------------------------------------------------- */
623 /* convert to packed, monotonic, simplicial, numeric */
624 /* ---------------------------------------------------------------------- */
625
626 /* leave as LL or LDL' */
627 if (!CHOLMOD(change_factor) (L->xtype, L->is_ll, FALSE, TRUE, TRUE, L,
628 Common))
629 {
630 ERROR (CHOLMOD_INVALID, "cannot convert L") ;
631 return (NULL) ;
632 }
633
634 /* ---------------------------------------------------------------------- */
635 /* create Lsparse */
636 /* ---------------------------------------------------------------------- */
637
638 /* allocate the header for Lsparse, the sparse matrix version of L */
639 Lsparse = CHOLMOD(malloc) (sizeof (cholmod_sparse), 1, Common) ;
640 if (Common->status < CHOLMOD_OK)
641 {
642 return (NULL) ; /* out of memory */
643 }
644
645 /* transfer the contents from L to Lsparse */
646 Lsparse->nrow = L->n ;
647 Lsparse->ncol = L->n ;
648 Lsparse->p = L->p ;
649 Lsparse->i = L->i ;
650 Lsparse->x = L->x ;
651 Lsparse->z = L->z ;
652 Lsparse->nz = NULL ;
653 Lsparse->stype = 0 ;
654 Lsparse->itype = L->itype ;
655 Lsparse->xtype = L->xtype ;
656 Lsparse->dtype = L->dtype ;
657 Lsparse->sorted = TRUE ;
658 Lsparse->packed = TRUE ;
659 Lsparse->nzmax = L->nzmax ;
660 ASSERT (CHOLMOD(dump_sparse) (Lsparse, "Lsparse", Common) >= 0) ;
661
662 /* ---------------------------------------------------------------------- */
663 /* convert L to symbolic, but do not free contents transfered to Lsparse */
664 /* ---------------------------------------------------------------------- */
665
666 L->p = NULL ;
667 L->i = NULL ;
668 L->x = NULL ;
669 L->z = NULL ;
670 L->xtype = CHOLMOD_PATTERN ;
671 CHOLMOD(change_factor) (CHOLMOD_PATTERN, FALSE, FALSE, TRUE, TRUE, L,
672 Common) ;
673
674 return (Lsparse) ;
675 }
676
677
678 /* ========================================================================== */
679 /* === cholmod_copy_factor ================================================== */
680 /* ========================================================================== */
681
682 /* Create an exact copy of a factor, with one exception:
683 *
684 * Entries in unused space are not copied (they might not be initialized,
685 * and copying them would cause program checkers such as purify and
686 * valgrind to complain).
687 *
688 * Note that a supernodal L cannot be zomplex.
689 */
690
CHOLMOD(copy_factor)691 cholmod_factor *CHOLMOD(copy_factor)
692 (
693 /* ---- input ---- */
694 cholmod_factor *L, /* factor to copy */
695 /* --------------- */
696 cholmod_common *Common
697 )
698 {
699 cholmod_factor *L2 ;
700 double *Lx, *L2x, *Lz, *L2z ;
701 Int *Perm, *ColCount, *Lp, *Li, *Lnz, *Lnext, *Lprev, *Lsuper, *Lpi, *Lpx,
702 *Ls, *Perm2, *ColCount2, *L2p, *L2i, *L2nz, *L2next, *L2prev, *L2super,
703 *L2pi, *L2px, *L2s ;
704 Int n, j, p, pend, s, xsize, ssize, nsuper ;
705
706 /* ---------------------------------------------------------------------- */
707 /* get inputs */
708 /* ---------------------------------------------------------------------- */
709
710 RETURN_IF_NULL_COMMON (NULL) ;
711 RETURN_IF_NULL (L, NULL) ;
712 RETURN_IF_XTYPE_INVALID (L, CHOLMOD_PATTERN, CHOLMOD_ZOMPLEX, NULL) ;
713 Common->status = CHOLMOD_OK ;
714 DEBUG (CHOLMOD(dump_factor) (L, "start copy", Common)) ;
715
716 n = L->n ;
717
718 /* ---------------------------------------------------------------------- */
719 /* allocate a simplicial symbolic factor */
720 /* ---------------------------------------------------------------------- */
721
722 /* allocates L2->Perm and L2->ColCount */
723 L2 = CHOLMOD(allocate_factor) (n, Common) ;
724 if (Common->status < CHOLMOD_OK)
725 {
726 return (NULL) ; /* out of memory */
727 }
728 ASSERT (L2->xtype == CHOLMOD_PATTERN && !(L2->is_super)) ;
729
730 Perm = L->Perm ;
731 ColCount = L->ColCount ;
732 Perm2 = L2->Perm ;
733 ColCount2 = L2->ColCount ;
734 L2->ordering = L->ordering ;
735
736 for (j = 0 ; j < n ; j++)
737 {
738 Perm2 [j] = Perm [j] ;
739 }
740 for (j = 0 ; j < n ; j++)
741 {
742 ColCount2 [j] = ColCount [j] ;
743 }
744 L2->is_ll = L->is_ll ;
745
746 /* ---------------------------------------------------------------------- */
747 /* copy the rest of the factor */
748 /* ---------------------------------------------------------------------- */
749
750 if (L->xtype != CHOLMOD_PATTERN && !(L->super))
751 {
752
753 /* ------------------------------------------------------------------ */
754 /* allocate a simplicial numeric factor */
755 /* ------------------------------------------------------------------ */
756
757 /* allocate L2->p, L2->nz, L2->prev, L2->next, L2->i, and L2->x.
758 * packed = -1 so that cholmod_change_factor allocates space of
759 * size L2->nzmax */
760 L2->nzmax = L->nzmax ;
761 if (!CHOLMOD(change_factor) (L->xtype, L->is_ll, FALSE, -1, TRUE,
762 L2, Common))
763 {
764 CHOLMOD(free_factor) (&L2, Common) ;
765 return (NULL) ; /* out of memory */
766 }
767 ASSERT (MAX (1, L->nzmax) == L2->nzmax) ;
768
769 /* ------------------------------------------------------------------ */
770 /* copy the contents of a simplicial numeric factor */
771 /* ------------------------------------------------------------------ */
772
773 Lp = L->p ;
774 Li = L->i ;
775 Lx = L->x ;
776 Lz = L->z ;
777 Lnz = L->nz ;
778 Lnext = L->next ;
779 Lprev = L->prev ;
780
781 L2p = L2->p ;
782 L2i = L2->i ;
783 L2x = L2->x ;
784 L2z = L2->z ;
785 L2nz = L2->nz ;
786 L2next = L2->next ;
787 L2prev = L2->prev ;
788 L2->xtype = L->xtype ;
789 L2->dtype = L->dtype ;
790
791 for (j = 0 ; j <= n ; j++)
792 {
793 L2p [j] = Lp [j] ;
794 }
795
796 for (j = 0 ; j < n+2 ; j++)
797 {
798 L2prev [j] = Lprev [j] ;
799 }
800
801 for (j = 0 ; j < n+2 ; j++)
802 {
803 L2next [j] = Lnext [j] ;
804 }
805
806 for (j = 0 ; j < n ; j++)
807 {
808 L2nz [j] = Lnz [j] ;
809 }
810
811 for (j = 0 ; j < n ; j++)
812 {
813 p = Lp [j] ;
814 pend = p + Lnz [j] ;
815 for ( ; p < pend ; p++)
816 {
817 L2i [p] = Li [p] ;
818 }
819 p = Lp [j] ;
820
821 if (L->xtype == CHOLMOD_REAL)
822 {
823 for ( ; p < pend ; p++)
824 {
825 L2x [p] = Lx [p] ;
826 }
827 }
828 else if (L->xtype == CHOLMOD_COMPLEX)
829 {
830 for ( ; p < pend ; p++)
831 {
832 L2x [2*p ] = Lx [2*p ] ;
833 L2x [2*p+1] = Lx [2*p+1] ;
834 }
835 }
836 else if (L->xtype == CHOLMOD_ZOMPLEX)
837 {
838 for ( ; p < pend ; p++)
839 {
840 L2x [p] = Lx [p] ;
841 L2z [p] = Lz [p] ;
842 }
843 }
844
845 }
846
847 }
848 else if (L->is_super)
849 {
850
851 /* ------------------------------------------------------------------ */
852 /* copy a supernodal factor */
853 /* ------------------------------------------------------------------ */
854
855 xsize = L->xsize ;
856 ssize = L->ssize ;
857 nsuper = L->nsuper ;
858
859 L2->xsize = xsize ;
860 L2->ssize = ssize ;
861 L2->nsuper = nsuper ;
862
863 /* allocate L2->super, L2->pi, L2->px, and L2->s. Allocate L2->x if
864 * L is numeric */
865 if (!CHOLMOD(change_factor) (L->xtype, TRUE, TRUE, TRUE, TRUE, L2,
866 Common))
867 {
868 CHOLMOD(free_factor) (&L2, Common) ;
869 return (NULL) ; /* out of memory */
870 }
871
872 ASSERT (L2->s != NULL) ;
873
874 /* ------------------------------------------------------------------ */
875 /* copy the contents of a supernodal factor */
876 /* ------------------------------------------------------------------ */
877
878 Lsuper = L->super ;
879 Lpi = L->pi ;
880 Lpx = L->px ;
881 Ls = L->s ;
882 Lx = L->x ;
883
884 L2super = L2->super ;
885 L2pi = L2->pi ;
886 L2px = L2->px ;
887 L2s = L2->s ;
888 L2x = L2->x ;
889
890 L2->maxcsize = L->maxcsize ;
891 L2->maxesize = L->maxesize ;
892
893 for (s = 0 ; s <= nsuper ; s++)
894 {
895 L2super [s] = Lsuper [s] ;
896 }
897 for (s = 0 ; s <= nsuper ; s++)
898 {
899 L2pi [s] = Lpi [s] ;
900 }
901 for (s = 0 ; s <= nsuper ; s++)
902 {
903 L2px [s] = Lpx [s] ;
904 }
905
906 L2s [0] = 0 ;
907 for (p = 0 ; p < ssize ; p++)
908 {
909 L2s [p] = Ls [p] ;
910 }
911
912 if (L->xtype == CHOLMOD_REAL)
913 {
914 for (p = 0 ; p < xsize ; p++)
915 {
916 L2x [p] = Lx [p] ;
917 }
918 }
919 else if (L->xtype == CHOLMOD_COMPLEX)
920 {
921 for (p = 0 ; p < 2*xsize ; p++)
922 {
923 L2x [p] = Lx [p] ;
924 }
925 }
926 }
927
928 L2->minor = L->minor ;
929 L2->is_monotonic = L->is_monotonic ;
930
931 DEBUG (CHOLMOD(dump_factor) (L2, "L2 got copied", Common)) ;
932 ASSERT (L2->xtype == L->xtype && L2->is_super == L->is_super) ;
933 return (L2) ;
934 }
935