1 /* ========================================================================== */
2 /* === Core/cholmod_dense =================================================== */
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_dense object:
11 *
12 * The solve routines and some of the MatrixOps and Modify routines use dense
13 * matrices as inputs. These are held in column-major order. With a leading
14 * dimension of d, the entry in row i and column j is held in x [i+j*d].
15 *
16 * Primary routines:
17 * -----------------
18 * cholmod_allocate_dense allocate a dense matrix
19 * cholmod_free_dense free a dense matrix
20 *
21 * Secondary routines:
22 * -------------------
23 * cholmod_zeros allocate a dense matrix of all zeros
24 * cholmod_ones allocate a dense matrix of all ones
25 * cholmod_eye allocate a dense identity matrix
26 * cholmod_sparse_to_dense create a dense matrix copy of a sparse matrix
27 * cholmod_dense_to_sparse create a sparse matrix copy of a dense matrix
28 * cholmod_copy_dense create a copy of a dense matrix
29 * cholmod_copy_dense2 copy a dense matrix (pre-allocated)
30 *
31 * All routines in this file can handle the real, complex, and zomplex cases.
32 * Pattern-only dense matrices are not supported. cholmod_sparse_to_dense can
33 * take a pattern-only input sparse matrix, however, and cholmod_dense_to_sparse
34 * can generate a pattern-only output sparse matrix.
35 */
36
37 #include "cholmod_internal.h"
38 #include "cholmod_core.h"
39
40 /* ========================================================================== */
41 /* === TEMPLATE ============================================================= */
42 /* ========================================================================== */
43
44 #define PATTERN
45 #include "t_cholmod_dense.c"
46 #define REAL
47 #include "t_cholmod_dense.c"
48 #define COMPLEX
49 #include "t_cholmod_dense.c"
50 #define ZOMPLEX
51 #include "t_cholmod_dense.c"
52
53
54 /* ========================================================================== */
55 /* === cholmod_allocate_dense =============================================== */
56 /* ========================================================================== */
57
58 /* Allocate a dense matrix with leading dimension d. The space is not
59 * initialized.
60 */
61
CHOLMOD(allocate_dense)62 cholmod_dense *CHOLMOD(allocate_dense)
63 (
64 /* ---- input ---- */
65 size_t nrow, /* # of rows of matrix */
66 size_t ncol, /* # of columns of matrix */
67 size_t d, /* leading dimension */
68 int xtype, /* CHOLMOD_REAL, _COMPLEX, or _ZOMPLEX */
69 /* --------------- */
70 cholmod_common *Common
71 )
72 {
73 cholmod_dense *X ;
74 size_t nzmax, nzmax0 ;
75 int ok = TRUE ;
76
77 /* ---------------------------------------------------------------------- */
78 /* get inputs */
79 /* ---------------------------------------------------------------------- */
80
81 RETURN_IF_NULL_COMMON (NULL) ;
82 if (d < nrow)
83 {
84 ERROR (CHOLMOD_INVALID, "leading dimension invalid") ;
85 return (NULL) ;
86 }
87 if (xtype < CHOLMOD_REAL || xtype > CHOLMOD_ZOMPLEX)
88 {
89 ERROR (CHOLMOD_INVALID, "xtype invalid") ;
90 return (NULL) ;
91 }
92
93 /* ensure the dimensions do not cause integer overflow */
94 (void) CHOLMOD(add_size_t) (ncol, 2, &ok) ;
95
96 /* nzmax = MAX (1, d*ncol) ; */
97 nzmax = CHOLMOD(mult_size_t) (d, ncol, &ok) ;
98 nzmax = MAX (1, nzmax) ;
99
100 if (!ok || nrow > Int_max || ncol > Int_max || nzmax > Int_max)
101 {
102 ERROR (CHOLMOD_TOO_LARGE, "problem too large") ;
103 return (NULL) ;
104 }
105 Common->status = CHOLMOD_OK ;
106
107 /* ---------------------------------------------------------------------- */
108 /* allocate header */
109 /* ---------------------------------------------------------------------- */
110
111 X = CHOLMOD(malloc) (sizeof (cholmod_dense), 1, Common) ;
112 if (Common->status < CHOLMOD_OK)
113 {
114 return (NULL) ; /* out of memory */
115 }
116
117 PRINT1 (("cholmod_allocate_dense %d-by-%d nzmax %d xtype %d\n",
118 nrow, ncol, nzmax, xtype)) ;
119
120 X->nrow = nrow ;
121 X->ncol = ncol ;
122 X->nzmax = nzmax ;
123 X->xtype = xtype ;
124 X->dtype = DTYPE ;
125 X->x = NULL ;
126 X->z = NULL ;
127 X->d = d ;
128
129 /* ---------------------------------------------------------------------- */
130 /* allocate the matrix itself */
131 /* ---------------------------------------------------------------------- */
132
133 nzmax0 = 0 ;
134 CHOLMOD(realloc_multiple) (nzmax, 0, xtype, NULL, NULL, &(X->x), &(X->z),
135 &nzmax0, Common) ;
136
137 if (Common->status < CHOLMOD_OK)
138 {
139 CHOLMOD(free_dense) (&X, Common) ;
140 return (NULL) ; /* out of memory */
141 }
142
143 return (X) ;
144 }
145
146
147 /* ========================================================================== */
148 /* === cholmod_zeros ======================================================== */
149 /* ========================================================================== */
150
151 /* Allocate a dense matrix and set it to zero */
152
CHOLMOD(zeros)153 cholmod_dense *CHOLMOD(zeros)
154 (
155 /* ---- input ---- */
156 size_t nrow, /* # of rows of matrix */
157 size_t ncol, /* # of columns of matrix */
158 int xtype, /* CHOLMOD_REAL, _COMPLEX, or _ZOMPLEX */
159 /* --------------- */
160 cholmod_common *Common
161 )
162 {
163 cholmod_dense *X ;
164 double *Xx, *Xz ;
165 Int i, nz ;
166
167 /* ---------------------------------------------------------------------- */
168 /* allocate a dense matrix and set it to zero */
169 /* ---------------------------------------------------------------------- */
170
171 RETURN_IF_NULL_COMMON (NULL) ;
172 X = CHOLMOD(allocate_dense) (nrow, ncol, nrow, xtype, Common) ;
173 if (Common->status < CHOLMOD_OK)
174 {
175 return (NULL) ; /* NULL Common, out of memory, or inputs invalid */
176 }
177
178 Xx = X->x ;
179 Xz = X->z ;
180 nz = MAX (1, X->nzmax) ;
181
182 switch (xtype)
183 {
184 case CHOLMOD_REAL:
185 for (i = 0 ; i < nz ; i++)
186 {
187 Xx [i] = 0 ;
188 }
189 break ;
190
191 case CHOLMOD_COMPLEX:
192 for (i = 0 ; i < 2*nz ; i++)
193 {
194 Xx [i] = 0 ;
195 }
196 break ;
197
198 case CHOLMOD_ZOMPLEX:
199 for (i = 0 ; i < nz ; i++)
200 {
201 Xx [i] = 0 ;
202 }
203 for (i = 0 ; i < nz ; i++)
204 {
205 Xz [i] = 0 ;
206 }
207 break ;
208 }
209
210 return (X) ;
211 }
212
213
214 /* ========================================================================== */
215 /* === cholmod_ones ========================================================= */
216 /* ========================================================================== */
217
218 /* Allocate a dense matrix and set it to zero */
219
CHOLMOD(ones)220 cholmod_dense *CHOLMOD(ones)
221 (
222 /* ---- input ---- */
223 size_t nrow, /* # of rows of matrix */
224 size_t ncol, /* # of columns of matrix */
225 int xtype, /* CHOLMOD_REAL, _COMPLEX, or _ZOMPLEX */
226 /* --------------- */
227 cholmod_common *Common
228 )
229 {
230 cholmod_dense *X ;
231 double *Xx, *Xz ;
232 Int i, nz ;
233
234 /* ---------------------------------------------------------------------- */
235 /* allocate a dense matrix and set it to all ones */
236 /* ---------------------------------------------------------------------- */
237
238 RETURN_IF_NULL_COMMON (NULL) ;
239 X = CHOLMOD(allocate_dense) (nrow, ncol, nrow, xtype, Common) ;
240 if (Common->status < CHOLMOD_OK)
241 {
242 return (NULL) ; /* NULL Common, out of memory, or inputs invalid */
243 }
244
245 Xx = X->x ;
246 Xz = X->z ;
247 nz = MAX (1, X->nzmax) ;
248
249 switch (xtype)
250 {
251 case CHOLMOD_REAL:
252 for (i = 0 ; i < nz ; i++)
253 {
254 Xx [i] = 1 ;
255 }
256 break ;
257
258 case CHOLMOD_COMPLEX:
259 for (i = 0 ; i < nz ; i++)
260 {
261 Xx [2*i ] = 1 ;
262 Xx [2*i+1] = 0 ;
263 }
264 break ;
265
266 case CHOLMOD_ZOMPLEX:
267 for (i = 0 ; i < nz ; i++)
268 {
269 Xx [i] = 1 ;
270 }
271 for (i = 0 ; i < nz ; i++)
272 {
273 Xz [i] = 0 ;
274 }
275 break ;
276 }
277
278 return (X) ;
279 }
280
281
282 /* ========================================================================== */
283 /* === cholmod_eye ========================================================== */
284 /* ========================================================================== */
285
286 /* Allocate a dense matrix and set it to the identity matrix */
287
CHOLMOD(eye)288 cholmod_dense *CHOLMOD(eye)
289 (
290 /* ---- input ---- */
291 size_t nrow, /* # of rows of matrix */
292 size_t ncol, /* # of columns of matrix */
293 int xtype, /* CHOLMOD_REAL, _COMPLEX, or _ZOMPLEX */
294 /* --------------- */
295 cholmod_common *Common
296 )
297 {
298 cholmod_dense *X ;
299 double *Xx, *Xz ;
300 Int i, n, nz ;
301
302 /* ---------------------------------------------------------------------- */
303 /* allocate a dense matrix and set it to the identity matrix */
304 /* ---------------------------------------------------------------------- */
305
306 RETURN_IF_NULL_COMMON (NULL) ;
307 X = CHOLMOD(zeros) (nrow, ncol, xtype, Common) ;
308 if (Common->status < CHOLMOD_OK)
309 {
310 return (NULL) ; /* NULL Common, out of memory, or inputs invalid */
311 }
312
313 nz = MAX (1, nrow*ncol) ;
314 Xx = X->x ;
315 Xz = X->z ;
316
317 n = MIN (nrow, ncol) ;
318
319 switch (xtype)
320 {
321 case CHOLMOD_REAL:
322 case CHOLMOD_ZOMPLEX:
323 for (i = 0 ; i < n ; i++)
324 {
325 Xx [i + i*nrow] = 1 ;
326 }
327 break ;
328
329 case CHOLMOD_COMPLEX:
330 for (i = 0 ; i < n ; i++)
331 {
332 Xx [2 * (i + i*nrow)] = 1 ;
333 }
334 break ;
335 }
336
337 return (X) ;
338 }
339
340 /* ========================================================================== */
341 /* === cholmod_free_dense =================================================== */
342 /* ========================================================================== */
343
344 /* free a dense matrix
345 *
346 * workspace: none
347 */
348
CHOLMOD(free_dense)349 int CHOLMOD(free_dense)
350 (
351 /* ---- in/out --- */
352 cholmod_dense **XHandle, /* dense matrix to deallocate, NULL on output */
353 /* --------------- */
354 cholmod_common *Common
355 )
356 {
357 cholmod_dense *X ;
358
359 RETURN_IF_NULL_COMMON (FALSE) ;
360
361 if (XHandle == NULL)
362 {
363 /* nothing to do */
364 return (TRUE) ;
365 }
366 X = *XHandle ;
367 if (X == NULL)
368 {
369 /* nothing to do */
370 return (TRUE) ;
371 }
372
373 switch (X->xtype)
374 {
375 case CHOLMOD_REAL:
376 X->x = CHOLMOD(free) (X->nzmax, sizeof (double), X->x, Common) ;
377 break ;
378
379 case CHOLMOD_COMPLEX:
380 X->x = CHOLMOD(free) (X->nzmax, 2*sizeof (double), X->x, Common) ;
381 break ;
382
383 case CHOLMOD_ZOMPLEX:
384 X->x = CHOLMOD(free) (X->nzmax, sizeof (double), X->x, Common) ;
385 X->z = CHOLMOD(free) (X->nzmax, sizeof (double), X->z, Common) ;
386 break ;
387 }
388
389 *XHandle = CHOLMOD(free) (1, sizeof (cholmod_dense), (*XHandle), Common) ;
390 return (TRUE) ;
391 }
392
393 /* ========================================================================== */
394 /* === cholmod_ensure_dense ================================================= */
395 /* ========================================================================== */
396
397 /* Ensure that the input matrix has a certain size and type. If not, free
398 * the existing matrix and reallocate one of the right size and type.
399 * Returns a pointer to the cholmod_dense matrix, possibly reallocated.
400 * Also modifies the input matrix handle, XHandle, if necessary.
401 */
402
CHOLMOD(ensure_dense)403 cholmod_dense *CHOLMOD(ensure_dense)
404 (
405 /* ---- input/output ---- */
406 cholmod_dense **XHandle, /* matrix handle to check */
407 /* ---- input ---- */
408 size_t nrow, /* # of rows of matrix */
409 size_t ncol, /* # of columns of matrix */
410 size_t d, /* leading dimension */
411 int xtype, /* CHOLMOD_REAL, _COMPLEX, or _ZOMPLEX */
412 /* --------------- */
413 cholmod_common *Common
414 )
415 {
416 cholmod_dense *X ;
417
418 RETURN_IF_NULL_COMMON (NULL) ;
419 if (XHandle == NULL)
420 {
421 ERROR (CHOLMOD_INVALID, "matrix invalid") ;
422 return (NULL) ;
423 }
424
425 X = *XHandle ;
426 if (X == NULL || X->nrow != nrow || X->ncol != ncol
427 || X->d != d || X->xtype != xtype)
428 {
429 /* Matrix X is not allocated, or has the wrong size. Free it and
430 * reallocate it in the right size and shape. If an error occurs
431 * (out of memory or inputs nrow, etc invalid), then the error is
432 * set in cholmod_allocate_dense and X is returned as NULL. */
433 CHOLMOD(free_dense) (XHandle, Common) ;
434 X = CHOLMOD(allocate_dense) (nrow, ncol, d, xtype, Common) ;
435 *XHandle = X ;
436 }
437 return (X) ;
438 }
439
440
441 /* ========================================================================== */
442 /* === cholmod_sparse_to_dense ============================================== */
443 /* ========================================================================== */
444
445 /* Convert a sparse matrix to a dense matrix.
446 * The output dense matrix has the same xtype as the input sparse matrix,
447 * except that a pattern-only sparse matrix A is converted into a real dense
448 * matrix X, with 1's and 0's. All xtypes are supported.
449 */
450
CHOLMOD(sparse_to_dense)451 cholmod_dense *CHOLMOD(sparse_to_dense)
452 (
453 /* ---- input ---- */
454 cholmod_sparse *A, /* matrix to copy */
455 /* --------------- */
456 cholmod_common *Common
457 )
458 {
459 cholmod_dense *X = NULL ;
460
461 /* ---------------------------------------------------------------------- */
462 /* check inputs */
463 /* ---------------------------------------------------------------------- */
464
465 RETURN_IF_NULL_COMMON (NULL) ;
466 RETURN_IF_NULL (A, NULL) ;
467 RETURN_IF_XTYPE_INVALID (A, CHOLMOD_PATTERN, CHOLMOD_ZOMPLEX, NULL) ;
468 if (A->stype && A->nrow != A->ncol)
469 {
470 ERROR (CHOLMOD_INVALID, "matrix invalid") ;
471 return (NULL) ;
472 }
473 Common->status = CHOLMOD_OK ;
474 ASSERT (CHOLMOD(dump_sparse) (A, "A", Common) >= 0) ;
475
476 /* ---------------------------------------------------------------------- */
477 /* convert the matrix, using template routine */
478 /* ---------------------------------------------------------------------- */
479
480 switch (A->xtype)
481 {
482 case CHOLMOD_PATTERN:
483 X = p_cholmod_sparse_to_dense (A, Common) ;
484 break ;
485
486 case CHOLMOD_REAL:
487 X = r_cholmod_sparse_to_dense (A, Common) ;
488 break ;
489
490 case CHOLMOD_COMPLEX:
491 X = c_cholmod_sparse_to_dense (A, Common) ;
492 break ;
493
494 case CHOLMOD_ZOMPLEX:
495 X = z_cholmod_sparse_to_dense (A, Common) ;
496 break ;
497 }
498 return (X) ;
499 }
500
501
502 /* ========================================================================== */
503 /* === cholmod_dense_to_sparse ============================================== */
504 /* ========================================================================== */
505
506 /* Convert a dense matrix to a sparse matrix, similar to the MATLAB statements:
507 *
508 * C = sparse (X) values = TRUE
509 * C = spones (sparse (X)) values = FALSE
510 *
511 * except that X must be double (it can be of many different types in MATLAB)
512 *
513 * The resulting sparse matrix C has the same numeric xtype as the input dense
514 * matrix X, unless "values" is FALSE (in which case C is real, where C(i,j)=1
515 * if (i,j) is an entry in X.
516 */
517
CHOLMOD(dense_to_sparse)518 cholmod_sparse *CHOLMOD(dense_to_sparse)
519 (
520 /* ---- input ---- */
521 cholmod_dense *X, /* matrix to copy */
522 int values, /* TRUE if values to be copied, FALSE otherwise */
523 /* --------------- */
524 cholmod_common *Common
525 )
526 {
527 cholmod_sparse *C = NULL ;
528
529 DEBUG (CHOLMOD(dump_dense) (X, "X", Common)) ;
530
531 /* ---------------------------------------------------------------------- */
532 /* check inputs */
533 /* ---------------------------------------------------------------------- */
534
535 RETURN_IF_NULL_COMMON (NULL) ;
536 RETURN_IF_NULL (X, NULL) ;
537 RETURN_IF_XTYPE_INVALID (X, CHOLMOD_REAL, CHOLMOD_ZOMPLEX, NULL) ;
538 if (X->d < X->nrow)
539 {
540 ERROR (CHOLMOD_INVALID, "matrix invalid") ;
541 return (NULL) ;
542 }
543 Common->status = CHOLMOD_OK ;
544
545 /* ---------------------------------------------------------------------- */
546 /* convert the matrix, using template routine */
547 /* ---------------------------------------------------------------------- */
548
549 switch (X->xtype)
550 {
551 case CHOLMOD_REAL:
552 C = r_cholmod_dense_to_sparse (X, values, Common) ;
553 break ;
554
555 case CHOLMOD_COMPLEX:
556 C = c_cholmod_dense_to_sparse (X, values, Common) ;
557 break ;
558
559 case CHOLMOD_ZOMPLEX:
560 C = z_cholmod_dense_to_sparse (X, values, Common) ;
561 break ;
562 }
563 return (C) ;
564 }
565
566
567 /* ========================================================================== */
568 /* === cholmod_copy_dense2 ================================================== */
569 /* ========================================================================== */
570
571 /* Y = X, where X and Y are both already allocated. The leading dimensions of
572 * X and Y may differ, but both must be >= the # of rows in X and Y.
573 * Entries in rows nrow to d-1 are not copied from X, since the space might not
574 * be initialized. Y->nzmax is unchanged. X->nzmax is typically
575 * (X->d)*(X->ncol), but a user might modify that condition outside of any
576 * CHOLMOD routine.
577 *
578 * The two dense matrices X and Y must have the same numeric xtype.
579 */
580
CHOLMOD(copy_dense2)581 int CHOLMOD(copy_dense2)
582 (
583 /* ---- input ---- */
584 cholmod_dense *X, /* matrix to copy */
585 /* ---- output --- */
586 cholmod_dense *Y, /* copy of matrix X */
587 /* --------------- */
588 cholmod_common *Common
589 )
590 {
591 /* ---------------------------------------------------------------------- */
592 /* check inputs */
593 /* ---------------------------------------------------------------------- */
594
595 RETURN_IF_NULL_COMMON (FALSE) ;
596 RETURN_IF_NULL (X, FALSE) ;
597 RETURN_IF_NULL (Y, FALSE) ;
598 RETURN_IF_XTYPE_INVALID (X, CHOLMOD_REAL, CHOLMOD_ZOMPLEX, FALSE) ;
599 RETURN_IF_XTYPE_INVALID (Y, CHOLMOD_REAL, CHOLMOD_ZOMPLEX, FALSE) ;
600 if (X->nrow != Y->nrow || X->ncol != Y->ncol || X->xtype != Y->xtype)
601 {
602 ERROR (CHOLMOD_INVALID, "X and Y must have same dimensions and xtype") ;
603 return (FALSE) ;
604 }
605 if (X->d < X->nrow || Y->d < Y->nrow
606 || (X->d * X->ncol) > X->nzmax || (Y->d * Y->ncol) > Y->nzmax)
607 {
608 ERROR (CHOLMOD_INVALID, "X and/or Y invalid") ;
609 return (FALSE) ;
610 }
611 Common->status = CHOLMOD_OK ;
612
613 /* ---------------------------------------------------------------------- */
614 /* copy the matrix, using template routine */
615 /* ---------------------------------------------------------------------- */
616
617 switch (X->xtype)
618 {
619 case CHOLMOD_REAL:
620 r_cholmod_copy_dense2 (X, Y) ;
621 break ;
622
623 case CHOLMOD_COMPLEX:
624 c_cholmod_copy_dense2 (X, Y) ;
625 break ;
626
627 case CHOLMOD_ZOMPLEX:
628 z_cholmod_copy_dense2 (X, Y) ;
629 break ;
630 }
631 return (TRUE) ;
632 }
633
634
635 /* ========================================================================== */
636 /* === cholmod_copy_dense =================================================== */
637 /* ========================================================================== */
638
639 /* Y = X, copy a dense matrix */
640
CHOLMOD(copy_dense)641 cholmod_dense *CHOLMOD(copy_dense)
642 (
643 /* ---- input ---- */
644 cholmod_dense *X, /* matrix to copy */
645 /* --------------- */
646 cholmod_common *Common
647 )
648 {
649 cholmod_dense *Y ;
650
651 /* ---------------------------------------------------------------------- */
652 /* check inputs */
653 /* ---------------------------------------------------------------------- */
654
655 RETURN_IF_NULL_COMMON (NULL) ;
656 RETURN_IF_NULL (X, NULL) ;
657 RETURN_IF_XTYPE_INVALID (X, CHOLMOD_REAL, CHOLMOD_ZOMPLEX, NULL) ;
658 Common->status = CHOLMOD_OK ;
659
660 /* ---------------------------------------------------------------------- */
661 /* allocate result */
662 /* ---------------------------------------------------------------------- */
663
664 Y = CHOLMOD(allocate_dense) (X->nrow, X->ncol, X->d, X->xtype, Common) ;
665 if (Common->status < CHOLMOD_OK)
666 {
667 return (NULL) ; /* out of memory or X invalid */
668 }
669
670 /* ---------------------------------------------------------------------- */
671 /* Y = X */
672 /* ---------------------------------------------------------------------- */
673
674 /* This cannot fail (X and Y are allocated, and have the same nrow, ncol
675 * d, and xtype) */
676 CHOLMOD(copy_dense2) (X, Y, Common) ;
677
678 /* ---------------------------------------------------------------------- */
679 /* return result */
680 /* ---------------------------------------------------------------------- */
681
682 return (Y) ;
683 }
684