File: cholmod_rowadd.c

package info (click to toggle)
suitesparse 1%3A7.11.0%2Bdfsg-2
  • links: PTS, VCS
  • area: main
  • in suites: forky, sid
  • size: 258,172 kB
  • sloc: ansic: 1,153,566; cpp: 48,145; makefile: 4,997; fortran: 2,087; java: 1,826; sh: 1,113; ruby: 725; python: 676; asm: 371; sed: 166; awk: 44
file content (247 lines) | stat: -rw-r--r-- 8,336 bytes parent folder | download | duplicates (5)
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
//------------------------------------------------------------------------------
// CHOLMOD/Modify/cholmod_rowadd: add row/column to an LDL' factorization
//------------------------------------------------------------------------------

// CHOLMOD/Modify Module.  Copyright (C) 2005-2023, Timothy A. Davis,
// and William W. Hager. All Rights Reserved.
// SPDX-License-Identifier: GPL-2.0+

//------------------------------------------------------------------------------

// Adds a row and column to an LDL' factorization, and optionally updates the
// solution to Lx=b.
//
// workspace: Flag (nrow), Head (nrow+1), W (2*nrow), Iwork (2*nrow)
//
// Only real matrices are supported (single or double).  A symbolic L is
// converted into a numeric identity matrix before the row is added.
// The dtypes of all matrices must match, except when L is symbolic (in which
// case it is converted to the dtype of R).

#include "cholmod_internal.h"

#ifndef NGPL
#ifndef NMODIFY

//------------------------------------------------------------------------------
// icomp: for sorting by qsort
//------------------------------------------------------------------------------

static int icomp (Int *i, Int *j)
{
    if (*i < *j)
    {
        return (-1) ;
    }
    else
    {
        return (1) ;
    }
}

//------------------------------------------------------------------------------
// t_cholmod_rowadd_worker
//------------------------------------------------------------------------------

#define DOUBLE
#define REAL
#include "t_cholmod_rowadd_worker.c"

#undef  DOUBLE
#define SINGLE
#define REAL
#include "t_cholmod_rowadd_worker.c"

//------------------------------------------------------------------------------
// cholmod_rowadd
//------------------------------------------------------------------------------

// cholmod_rowadd adds a row to the LDL' factorization.  It computes the kth
// row and kth column of L, and then updates the submatrix L (k+1:n,k+1:n)
// accordingly.  The kth row and column of L should originally be equal to the
// kth row and column of the identity matrix (they are treated as such, if they
// are not).  The kth row/column of L is computed as the factorization of the
// kth row/column of the matrix to factorize, which is provided as a single
// n-by-1 sparse matrix R.  The sparse vector R need not be sorted.

int CHOLMOD(rowadd)
(
    // input:
    size_t k,           // row/column index to add
    cholmod_sparse *R,  // row/column of matrix to factorize (n-by-1)
    // input/output:
    cholmod_factor *L,  // factor to modify
    cholmod_common *Common
)
{
    double bk [2] ;
    bk [0] = 0. ;
    bk [1] = 0. ;
    return (CHOLMOD(rowadd_mark) (k, R, bk, NULL, L, NULL, NULL, Common)) ;
}

//------------------------------------------------------------------------------
// cholmod_rowadd_solve
//------------------------------------------------------------------------------

// Does the same as cholmod_rowadd, and also updates the solution to Lx=b
// See cholmod_updown for a description of how Lx=b is updated.  There is on
// additional parameter:  bk specifies the new kth entry of b.

int CHOLMOD(rowadd_solve)
(
    // input:
    size_t k,           // row/column index to add
    cholmod_sparse *R,  // row/column of matrix to factorize (n-by-1)
    double bk [2],      // kth entry of the right-hand-side b
    // input/output:
    cholmod_factor *L,  // factor to modify
    cholmod_dense *X,   // solution to Lx=b (size n-by-1)
    cholmod_dense *DeltaB,  // change in b, zero on output
    cholmod_common *Common
)
{
    return (CHOLMOD(rowadd_mark) (k, R, bk, NULL, L, X, DeltaB, Common)) ;
}

//------------------------------------------------------------------------------
// cholmod_rowadd_mark
//------------------------------------------------------------------------------

// Does the same as cholmod_rowadd_solve, except only part of L is used in
// the update/downdate of the solution to Lx=b.  This routine is an "expert"
// routine.  It is meant for use in LPDASA only.

int CHOLMOD(rowadd_mark)
(
    // input:
    size_t kadd,        // row/column index to add
    cholmod_sparse *R,  // row/column of matrix to factorize (n-by-1)
    double bk [2],      // kth entry of the right hand side, b
    Int *colmark,       // Int array of size 1.  See cholmod_updown.c
    // input/output:
    cholmod_factor *L,  // factor to modify
    cholmod_dense *X,   // solution to Lx=b (size n-by-1)
    cholmod_dense *DeltaB,  // change in b, zero on output
    cholmod_common *Common
)
{

    //--------------------------------------------------------------------------
    // check inputs
    //--------------------------------------------------------------------------

    RETURN_IF_NULL_COMMON (FALSE) ;
    RETURN_IF_NULL (L, FALSE) ;
    RETURN_IF_NULL (R, FALSE) ;
    RETURN_IF_XTYPE_INVALID (L, CHOLMOD_PATTERN, CHOLMOD_REAL, FALSE) ;
    RETURN_IF_XTYPE_INVALID (R, CHOLMOD_REAL, CHOLMOD_REAL, FALSE) ;
    Int k = kadd ;
    if (kadd >= L->n || k < 0)
    {
        ERROR (CHOLMOD_INVALID, "k invalid") ;
        return (FALSE) ;
    }
    if (R->ncol != 1 || R->nrow != L->n)
    {
        ERROR (CHOLMOD_INVALID, "R invalid") ;
        return (FALSE) ;
    }
    if (L->xtype != CHOLMOD_PATTERN && L->dtype != R->dtype)
    {
        ERROR (CHOLMOD_INVALID, "R and L must have the same dtype") ;
        return (FALSE) ;
    }

    if ((X != NULL) && (DeltaB != NULL))
    {
        // also update the solution to Lx=b
        RETURN_IF_XTYPE_INVALID (X, CHOLMOD_REAL, CHOLMOD_REAL, FALSE) ;
        RETURN_IF_XTYPE_INVALID (DeltaB, CHOLMOD_REAL, CHOLMOD_REAL, FALSE) ;
        if (X->nrow != L->n || X->ncol != 1 ||
            DeltaB->nrow != L->n || DeltaB->ncol != 1 ||
            X->dtype != R->dtype || DeltaB->dtype != R->dtype)
        {
            ERROR (CHOLMOD_INVALID, "X and/or DeltaB invalid") ;
            return (FALSE) ;
        }
    }

    Common->status = CHOLMOD_OK ;

    //--------------------------------------------------------------------------
    // allocate workspace
    //--------------------------------------------------------------------------

    // s = 2*n
    int ok = TRUE ;
    size_t s = CHOLMOD(mult_size_t) (L->n, 2, &ok) ;
    if (!ok)
    {
        ERROR (CHOLMOD_TOO_LARGE, "problem too large") ;
        return (FALSE) ;
    }

    CHOLMOD(alloc_work) (L->n, s, s, R->dtype, Common) ;
    if (Common->status < CHOLMOD_OK)
    {
        return (FALSE) ;
    }
    ASSERT (CHOLMOD(dump_work) (TRUE, TRUE, s, R->dtype, Common)) ;

    //--------------------------------------------------------------------------
    // convert to simplicial numeric LDL' factor, if not already
    //--------------------------------------------------------------------------

    if (L->xtype == CHOLMOD_PATTERN || L->is_super || L->is_ll)
    {
        // can only update/downdate a simplicial LDL' factorization
        if (L->xtype == CHOLMOD_PATTERN)
        {
            // L is symbolic; convert it to R->dtype
            L->dtype = R->dtype ;
        }
        CHOLMOD(change_factor) (CHOLMOD_REAL, FALSE, FALSE, FALSE, FALSE, L,
                Common) ;
        if (Common->status < CHOLMOD_OK)
        {
            // out of memory, L is returned unchanged
            return (FALSE) ;
        }
    }

    ASSERT (L->dtype == R->dtype) ;

    //--------------------------------------------------------------------------
    // update L and X
    //--------------------------------------------------------------------------

    float s_bk [2] ;
    s_bk [0] = (float) bk [0] ;
    s_bk [1] = (float) bk [1] ;

    switch (L->dtype & 4)
    {
        case CHOLMOD_SINGLE:
            ok = rs_cholmod_rowadd_worker (k, R, s_bk, colmark, L, X, DeltaB,
                Common) ;
            break ;

        case CHOLMOD_DOUBLE:
            ok = rd_cholmod_rowadd_worker (k, R, bk, colmark, L, X, DeltaB,
                Common) ;
            break ;
    }

    //--------------------------------------------------------------------------
    // return result
    //--------------------------------------------------------------------------

    DEBUG (CHOLMOD(dump_factor) (L, "LDL factorization, L:", Common)) ;
    ASSERT (CHOLMOD(dump_work) (TRUE, TRUE, s, R->dtype, Common)) ;
    return (ok) ;
}

#endif
#endif