1 // =============================================================================
2 // === spqr_cpack ==============================================================
3 // =============================================================================
4 
5 //  spqr_cpack copies the C matrix from the frontal matrix F and
6 //  stores it in packed form.  F can be overwritten with C (the pack can occur
7 //  in-place), but in that case, the R and H matrices are destroyed.
8 //
9 //  In this example, m = 5, n = 6, npiv = 2, and rank = 2 (the number of good
10 //  pivot columns found; equivalently, the number of rows in the R block).  The
11 //  number of columns in C is cn = n-npiv.  The number of rows is
12 //  cm = MIN (m-rank,cn).  In this example, MIN (m-rank,cn) = 3.
13 //
14 //      . . . . . .
15 //      . . . . . .
16 //      . . c c c c     <- rank = 2
17 //      . . . c c c
18 //      . . . . c c
19 //
20 //  In the next example below, m = 8 instead.  Note that
21 //  cm = MIN (m-rank,cn) = 4.
22 //
23 //      . . . . . .
24 //      . . . . . .
25 //      . . c c c c     <- rank = 2
26 //      . . . c c c
27 //      . . . . c c
28 //      . . . . . c
29 //      . . . . . .
30 //      . . . . . .
31 
32 #include "spqr.hpp"
33 
spqr_cpack(Long m,Long n,Long npiv,Long rank,Entry * F,Entry * C)34 template <typename Entry> Long spqr_cpack     // returns # of rows in C
35 (
36     // input, not modified
37     Long m,                 // # of rows in F
38     Long n,                 // # of columns in F
39     Long npiv,              // number of pivotal columns in F
40     Long rank,              // the C block starts at F (rank,npiv)
41 
42     // input, not modified unless the pack occurs in-place
43     Entry *F,               // m-by-n frontal matrix in column-major order
44 
45     // output, contents not defined on input
46     Entry *C                // packed columns of C, of size cm-by-cn in upper
47                             // trapezoidal form.
48 )
49 {
50     Long i, k, cm, cn ;
51 
52     // -------------------------------------------------------------------------
53     // get inputs
54     // -------------------------------------------------------------------------
55 
56     ASSERT (m >= 0 && n >= 0 && npiv >= 0 && npiv <= n) ;
57 
58     ASSERT (rank >= 0 && rank <= MIN (m,npiv)) ;
59     cn = n - npiv ;                     // number of columns of C
60     cm = MIN (m-rank, cn) ;             // number of rows of C
61     ASSERT (cm <= cn) ;
62     if (cm <= 0 || cn <= 0)
63     {
64         return (0) ;                    // nothing to do
65     }
66 
67     ASSERT (C != NULL && F != NULL) ;
68     ASSERT (C <= F                      // C can be packed in-place, in F
69          || C >= F + m*n) ;             // or C must appear after F
70 
71     F += INDEX (rank,npiv,m) ;          // C starts at F (rank,npiv)
72 
73     // -------------------------------------------------------------------------
74     // pack the upper triangular part of C
75     // -------------------------------------------------------------------------
76 
77     for (k = 0 ; k < cm ; k++)
78     {
79         // pack C (0:k,k)
80         for (i = 0 ; i <= k ; i++)
81         {
82             *(C++) = F [i] ;
83         }
84         F += m ;                        // advance to the next column of F
85     }
86 
87     // -------------------------------------------------------------------------
88     // pack the rectangular part of C
89     // -------------------------------------------------------------------------
90 
91     for ( ; k < cn ; k++)
92     {
93         // pack C (0:cm-1,k)
94         for (i = 0 ; i < cm ; i++)
95         {
96             *(C++) = F [i] ;
97         }
98         F += m ;                        // advance to the next column of F
99     }
100 
101     PR (("Cpack rank %ld cm %ld cn %ld\n", rank, cm, cn)) ;
102     return (cm) ;                       // return # of rows in C
103 }
104 
105 
106 // =============================================================================
107 
108 template Long spqr_cpack <double>     // returns # of rows in C
109 (
110     // input, not modified
111     Long m,                 // # of rows in F
112     Long n,                 // # of columns in F
113     Long npiv,              // number of pivotal columns in F
114     Long rank,              // the C block starts at F (rank,npiv)
115 
116     // input, not modified unless the pack occurs in-place
117     double *F,              // m-by-n frontal matrix in column-major order
118 
119     // output, contents not defined on input
120     double *C               // packed columns of C, of size cm-by-cn in upper
121                             // trapezoidal form.
122 ) ;
123 
124 // =============================================================================
125 
126 template Long spqr_cpack <Complex>    // returns # of rows in C
127 (
128     // input, not modified
129     Long m,                 // # of rows in F
130     Long n,                 // # of columns in F
131     Long npiv,              // number of pivotal columns in F
132     Long rank,              // the C block starts at F (rank,npiv)
133 
134     // input, not modified unless the pack occurs in-place
135     Complex *F,             // m-by-n frontal matrix in column-major order
136 
137     // output, contents not defined on input
138     Complex *C              // packed columns of C, of size cm-by-cn in upper
139                             // trapezoidal form.
140 ) ;
141