1 /* ========================================================================== */
2 /* === UMFPACK_triplet_to_col =============================================== */
3 /* ========================================================================== */
4 
5 /* -------------------------------------------------------------------------- */
6 /* Copyright (c) 2005-2012 by Timothy A. Davis, http://www.suitesparse.com.   */
7 /* All Rights Reserved.  See ../Doc/License.txt for License.                  */
8 /* -------------------------------------------------------------------------- */
9 
10 /*
11     User callable.  Converts triplet input to column-oriented form.  Duplicate
12     entries may exist (they are summed in the output).  The columns of the
13     column-oriented form are in sorted order.  The input is not modified.
14     Returns 1 if OK, 0 if an error occurred.  See umfpack_triplet_to_col.h for
15     details.
16 
17     If Map is present (a non-NULL pointer to an Int array of size nz), then on
18     output it holds the position of the triplets in the column-form matrix.
19     That is, suppose p = Map [k], and the k-th triplet is i=Ti[k], j=Tj[k], and
20     aij=Tx[k].  Then i=Ai[p], and aij will have been summed into Ax[p].  Also,
21     Ap[j] <= p < Ap[j+1].  The Map array is not computed if it is (Int *) NULL.
22 
23     Dynamic memory usage:
24 
25 	If numerical values are present, then one (two for complex version)
26 	workspace of size (nz+1)*sizeof(double) is allocated via UMF_malloc.
27 	Next, 4 calls to UMF_malloc are made to obtain workspace of size
28 	((nz+1) + (n_row+1) + n_row + MAX (n_row,n_col)) * sizeof(Int).  All of
29 	this workspace (4 to 6 objects) are free'd via UMF_free on return.
30 
31 	For the complex version, additional space is allocated.
32 
33 	An extra array of size nz*sizeof(Int) is allocated if Map is present.
34 */
35 
36 #include "umf_internal.h"
37 #include "umf_malloc.h"
38 #include "umf_free.h"
39 #include "umf_triplet.h"
40 
41 #ifndef NDEBUG
42 PRIVATE Int init_count ;
43 #endif
44 
45 /* ========================================================================== */
46 
UMFPACK_triplet_to_col(Int n_row,Int n_col,Int nz,const Int Ti[],const Int Tj[],const double Tx[],const double Tz[],Int Ap[],Int Ai[],double Ax[],double Az[],Int Map[])47 GLOBAL Int UMFPACK_triplet_to_col
48 (
49     Int n_row,
50     Int n_col,
51     Int nz,
52     const Int Ti [ ],		/* size nz */
53     const Int Tj [ ],		/* size nz */
54     const double Tx [ ],	/* size nz */
55 #ifdef COMPLEX
56     const double Tz [ ],	/* size nz */
57 #endif
58     Int Ap [ ],			/* size n_col + 1 */
59     Int Ai [ ],			/* size nz */
60     double Ax [ ]		/* size nz */
61 #ifdef COMPLEX
62     , double Az [ ]		/* size nz */
63 #endif
64     , Int Map [ ]		/* size nz */
65 )
66 {
67 
68     /* ---------------------------------------------------------------------- */
69     /* local variables */
70     /* ---------------------------------------------------------------------- */
71 
72     Int *RowCount, *Rp, *Rj, *W, nn, do_values, do_map, *Map2, status ;
73     double *Rx ;
74 #ifdef COMPLEX
75     double *Rz ;
76     Int split ;
77 #endif
78 
79 #ifndef NDEBUG
80     UMF_dump_start ( ) ;
81     init_count = UMF_malloc_count ;
82 #endif
83 
84     /* ---------------------------------------------------------------------- */
85     /* check inputs */
86     /* ---------------------------------------------------------------------- */
87 
88     if (!Ai || !Ap || !Ti || !Tj)
89     {
90 	return (UMFPACK_ERROR_argument_missing) ;
91     }
92 
93     if (n_row <= 0 || n_col <= 0)		/* must be > 0 */
94     {
95 	return (UMFPACK_ERROR_n_nonpositive) ;
96     }
97 
98     if (nz < 0)		/* nz must be >= 0 (singular matrices are OK) */
99     {
100 	return (UMFPACK_ERROR_invalid_matrix) ;
101     }
102 
103     nn = MAX (n_row, n_col) ;
104 
105     /* ---------------------------------------------------------------------- */
106     /* allocate workspace */
107     /* ---------------------------------------------------------------------- */
108 
109     Rx = (double *) NULL ;
110 
111     do_values = Ax && Tx ;
112 
113     if (do_values)
114     {
115 #ifdef COMPLEX
116 	Rx = (double *) UMF_malloc (2*nz+2, sizeof (double)) ;
117 	split = SPLIT (Tz) && SPLIT (Az) ;
118 	if (split)
119 	{
120 	    Rz = Rx + nz ;
121 	}
122 	else
123 	{
124 	    Rz = (double *) NULL ;
125 	}
126 #else
127 	Rx = (double *) UMF_malloc (nz+1, sizeof (double)) ;
128 #endif
129 	if (!Rx)
130 	{
131 	    DEBUGm4 (("out of memory: triplet work \n")) ;
132 	    ASSERT (UMF_malloc_count == init_count) ;
133 	    return (UMFPACK_ERROR_out_of_memory) ;
134 	}
135     }
136 
137     do_map = (Map != (Int *) NULL) ;
138     Map2 = (Int *) NULL ;
139     if (do_map)
140     {
141 	DEBUG0 (("Do map:\n")) ;
142 	Map2 = (Int *) UMF_malloc (nz+1, sizeof (Int)) ;
143 	if (!Map2)
144 	{
145 	    DEBUGm4 (("out of memory: triplet map\n")) ;
146 	    (void) UMF_free ((void *) Rx) ;
147 	    ASSERT (UMF_malloc_count == init_count) ;
148 	    return (UMFPACK_ERROR_out_of_memory) ;
149 	}
150     }
151 
152     Rj = (Int *) UMF_malloc (nz+1, sizeof (Int)) ;
153     Rp = (Int *) UMF_malloc (n_row+1, sizeof (Int)) ;
154     RowCount = (Int *) UMF_malloc (n_row, sizeof (Int)) ;
155     W = (Int *) UMF_malloc (nn, sizeof (Int)) ;
156     if (!Rj || !Rp || !RowCount || !W)
157     {
158 	DEBUGm4 (("out of memory: triplet work (int)\n")) ;
159 	(void) UMF_free ((void *) Rx) ;
160 	(void) UMF_free ((void *) Map2) ;
161 	(void) UMF_free ((void *) Rp) ;
162 	(void) UMF_free ((void *) Rj) ;
163 	(void) UMF_free ((void *) RowCount) ;
164 	(void) UMF_free ((void *) W) ;
165 	ASSERT (UMF_malloc_count == init_count) ;
166 	return (UMFPACK_ERROR_out_of_memory) ;
167     }
168 
169     ASSERT (UMF_malloc_count == init_count + 4 +
170 	(Rx != (double *) NULL) + do_map) ;
171 
172     /* ---------------------------------------------------------------------- */
173     /* convert from triplet to column form */
174     /* ---------------------------------------------------------------------- */
175 
176     if (do_map)
177     {
178 	if (do_values)
179 	{
180 	    status = UMF_triplet_map_x (n_row, n_col, nz, Ti, Tj, Ap, Ai, Rp,
181 		Rj, W, RowCount, Tx, Ax, Rx
182 #ifdef COMPLEX
183 		, Tz, Az, Rz
184 #endif
185 		, Map, Map2) ;
186 	}
187 	else
188 	{
189 	    status = UMF_triplet_map_nox (n_row, n_col, nz, Ti, Tj, Ap, Ai, Rp,
190 		Rj, W, RowCount, Map, Map2) ;
191 	}
192     }
193     else
194     {
195 	if (do_values)
196 	{
197 	    status = UMF_triplet_nomap_x (n_row, n_col, nz, Ti, Tj, Ap, Ai, Rp,
198 		Rj, W, RowCount , Tx, Ax, Rx
199 #ifdef COMPLEX
200 		, Tz, Az, Rz
201 #endif
202 		) ;
203 	}
204 	else
205 	{
206 	    status = UMF_triplet_nomap_nox (n_row, n_col, nz, Ti, Tj, Ap, Ai,
207 		Rp, Rj, W, RowCount) ;
208 	}
209     }
210 
211     /* ---------------------------------------------------------------------- */
212     /* free the workspace */
213     /* ---------------------------------------------------------------------- */
214 
215     (void) UMF_free ((void *) Rx) ;
216     (void) UMF_free ((void *) Map2) ;
217     (void) UMF_free ((void *) Rp) ;
218     (void) UMF_free ((void *) Rj) ;
219     (void) UMF_free ((void *) RowCount) ;
220     (void) UMF_free ((void *) W) ;
221     ASSERT (UMF_malloc_count == init_count) ;
222 
223     return (status) ;
224 }
225