• Home
  • Features
  • Pricing
  • Docs
  • Announcements
  • Sign In

dance858 / PSLP / 22554923004

01 Mar 2026 11:02PM UTC coverage: 89.536% (+0.5%) from 89.013%
22554923004

Pull #38

github

web-flow
Merge 73b50cb64 into e13878911
Pull Request #38: [WIP] Speedup

1272 of 1418 branches covered (89.7%)

Branch coverage included in aggregate %.

267 of 287 new or added lines in 12 files covered. (93.03%)

2 existing lines in 1 file now uncovered.

3956 of 4421 relevant lines covered (89.48%)

6275.8 hits per line

Source File
Press 'n' to go to next uncovered line, 'b' for previous

90.34
/src/explorers/Parallel_cols.c
1
/*
2
 * Copyright 2025 Daniel Cederberg
3
 *
4
 * This file is part of the PSLP project (LP Presolver).
5
 *
6
 * Licensed under the Apache License, Version 2.0 (the "License");
7
 * you may not use this file except in compliance with the License.
8
 * You may obtain a copy of the License at
9
 *
10
 *     http://www.apache.org/licenses/LICENSE-2.0
11
 *
12
 * Unless required by applicable law or agreed to in writing, software
13
 * distributed under the License is distributed on an "AS IS" BASIS,
14
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
15
 * See the License for the specific language governing permissions and
16
 * limitations under the License.
17
 */
18

19
#include "Parallel_cols.h"
20
#include "Activity.h"
21
#include "Bounds.h"
22
#include "Constraints.h"
23
#include "CoreTransformations.h"
24
#include "Debugger.h"
25
#include "Matrix.h"
26
#include "Numerics.h"
27
#include "Parallel_rows.h"
28
#include "Problem.h"
29
#include "RowColViews.h"
30
#include "SimpleReductions.h"
31
#include "State.h"
32
#include "Tags.h"
33
#include "Workspace.h"
34

35
static inline PresolveStatus process_single_bin(const Problem *prob, const int *bin,
18✔
36
                                                int bin_size, int *rows_to_recompute)
37
{
38
    assert(bin_size > 1);
18✔
39
    Constraints *constraints = prob->constraints;
18✔
40
    ColTag *col_tags = constraints->col_tags;
18✔
41
    Bound *bounds = constraints->bounds;
18✔
42
    const Matrix *AT = constraints->AT;
18✔
43
    const RowRange *col_ranges = AT->p;
18✔
44
    const double *c = prob->obj->c;
18✔
45
    iVec *sub_cols_to_delete = constraints->state->sub_cols_to_delete;
18✔
46
    PostsolveInfo *postsolve_info = constraints->state->postsolve_info;
18✔
47

48
    bool recount_ninfs = false;
18✔
49

50
    // column j stays, column k goes
51
    int ii, jj, j, k, len_j;
52
    const int *aj_cols;
53
    double ratio, cj, ck, ak0, lb_j_old, ub_j_old;
54
    const double *aj_vals;
55

56
    for (ii = 0; ii < bin_size - 1; ++ii)
56✔
57
    {
58
        j = bin[ii];
40✔
59

60
        if (HAS_TAG(col_tags[j], C_TAG_INACTIVE))
40✔
61
        {
62
            continue;
15✔
63
        }
64

65
        cj = c[j];
25✔
66
        aj_vals = AT->x + col_ranges[j].start;
25✔
67
        aj_cols = AT->i + col_ranges[j].start;
25✔
68
        len_j = col_ranges[j].end - col_ranges[j].start;
25✔
69
        recount_ninfs = false;
25✔
70

71
        for (jj = ii + 1; jj < bin_size; ++jj)
79✔
72
        {
73
            k = bin[jj];
56✔
74

75
            if (HAS_TAG(col_tags[k], C_TAG_INACTIVE))
56✔
76
            {
77
                continue;
5✔
78
            }
79

80
            ck = c[k];
51✔
81
            ak0 = AT->x[col_ranges[k].start];
51✔
82
            ratio = ak0 / aj_vals[0];
51✔
83
            assert(ratio != 0);
51✔
84

85
            // if column j and column k are parallel, remove column k
86
            // if you run into issues with parallel cols activated, try to
87
            // change this
88
            if (IS_EQUAL_FEAS_TOL(ck * aj_vals[0], cj * ak0))
51✔
89
            {
90
                // mark that we must recount n_infs for rows in which column j
91
                // appears in
92
                recount_ninfs = true;
18✔
93

94
                // ---------------------------------------------------------------------
95
                //                    Update the bounds of column j
96
                // ---------------------------------------------------------------------
97
                lb_j_old = bounds[j].lb;
18✔
98
                ub_j_old = bounds[j].ub;
18✔
99
                if (ratio > 0)
18✔
100
                {
101
                    // update lb
102
                    if (HAS_TAG((col_tags[j] | col_tags[k]), C_TAG_LB_INF))
11✔
103
                    {
104
                        bounds[j].lb = -INF;
2✔
105
                        UPDATE_TAG(col_tags[j], C_TAG_LB_INF);
2✔
106
                    }
107
                    else
108
                    {
109
                        assert(!IS_ABS_INF(bounds[j].lb) &&
9✔
110
                               !IS_ABS_INF(bounds[k].lb));
111
                        bounds[j].lb += bounds[k].lb * ratio;
9✔
112
                    }
113

114
                    // update ub
115
                    if (HAS_TAG((col_tags[j] | col_tags[k]), C_TAG_UB_INF))
11✔
116
                    {
117
                        bounds[j].ub = INF;
2✔
118
                        UPDATE_TAG(col_tags[j], C_TAG_UB_INF);
2✔
119
                    }
120
                    else
121
                    {
122
                        assert(!IS_ABS_INF(bounds[j].ub) &&
9✔
123
                               !IS_ABS_INF(bounds[k].ub));
124
                        bounds[j].ub += bounds[k].ub * ratio;
9✔
125
                    }
126
                }
127
                else
128
                {
129
                    // update lb
130
                    if (HAS_TAG(col_tags[j], C_TAG_LB_INF) ||
7✔
131
                        HAS_TAG(col_tags[k], C_TAG_UB_INF))
6✔
132
                    {
133
                        bounds[j].lb = -INF;
1✔
134
                        UPDATE_TAG(col_tags[j], C_TAG_LB_INF);
1✔
135
                    }
136
                    else
137
                    {
138
                        assert(!IS_ABS_INF(bounds[j].lb) &&
6✔
139
                               !IS_ABS_INF(bounds[k].ub));
140
                        bounds[j].lb += bounds[k].ub * ratio;
6✔
141
                    }
142

143
                    // update ub
144
                    if (HAS_TAG(col_tags[j], C_TAG_UB_INF) ||
7✔
145
                        HAS_TAG(col_tags[k], C_TAG_LB_INF))
7✔
146
                    {
147
                        bounds[j].ub = INF;
2✔
148
                        UPDATE_TAG(col_tags[j], C_TAG_UB_INF);
2✔
149
                    }
150
                    else
151
                    {
152
                        assert(!IS_ABS_INF(bounds[j].ub) &&
5✔
153
                               !IS_ABS_INF(bounds[k].lb));
154
                        bounds[j].ub += bounds[k].lb * ratio;
5✔
155
                    }
156
                }
157

158
                // ---------------------------------------------------------------------
159
                //                  mark col k as substituted
160
                // ---------------------------------------------------------------------
161
                set_col_to_substituted(k, col_tags + k, sub_cols_to_delete);
18✔
162
                save_retrieval_parallel_col(postsolve_info, ub_j_old, lb_j_old,
18✔
163
                                            bounds[k].lb, bounds[k].ub, ratio, j, k,
18✔
164
                                            col_tags[j], col_tags[k]);
18✔
165
            }
166
            else
167
            {
168
                // we could add an implied check here but it doesn't seem to
169
                // make a difference on miplib. If you do it, don't forget that
170
                // the implied free check implicitly corresponds to a bound
171
                // change that must be dual postsolved
172

173
                bool fix_xk_to_lower = false;
33✔
174
                bool fix_xk_to_upper = false;
33✔
175
                bool fix_xj_to_upper = false;
33✔
176
                bool fix_xj_to_lower = false;
33✔
177

178
                assert(!HAS_TAG(col_tags[k], C_TAG_INACTIVE));
33✔
179
                assert(!HAS_TAG(col_tags[j], C_TAG_INACTIVE));
33✔
180

181
                if (ck > ratio * cj)
33✔
182
                {
183
                    if (ratio > 0)
19✔
184
                    {
185
                        fix_xk_to_lower = HAS_TAG(col_tags[j], C_TAG_UB_INF);
11✔
186
                        fix_xj_to_upper = HAS_TAG(col_tags[k], C_TAG_LB_INF);
11✔
187
                    }
188
                    else
189
                    {
190
                        fix_xk_to_lower = HAS_TAG(col_tags[j], C_TAG_LB_INF);
8✔
191
                        fix_xj_to_lower = HAS_TAG(col_tags[k], C_TAG_LB_INF);
8✔
192
                    }
193
                }
194
                else
195
                {
196
                    assert(ck < ratio * cj);
14✔
197
                    if (ratio > 0)
14✔
198
                    {
199
                        fix_xk_to_upper = HAS_TAG(col_tags[j], C_TAG_LB_INF);
10✔
200
                        fix_xj_to_lower = HAS_TAG(col_tags[k], C_TAG_UB_INF);
10✔
201
                    }
202
                    else
203
                    {
204
                        fix_xk_to_upper = HAS_TAG(col_tags[j], C_TAG_UB_INF);
4✔
205
                        fix_xj_to_upper = HAS_TAG(col_tags[k], C_TAG_UB_INF);
4✔
206
                    }
207
                }
208

209
                // check if we can fix xk
210
                if (fix_xk_to_lower)
33✔
211
                {
212
                    if (HAS_TAG(col_tags[k], C_TAG_LB_INF))
3✔
213
                    {
214
                        return UNBNDORINFEAS;
1✔
215
                    }
216

217
                    assert(!HAS_TAG(col_tags[k], C_TAG_INACTIVE));
2✔
218
                    assert(!IS_ABS_INF(bounds[k].lb));
2✔
219
                    fix_col(constraints, k, bounds[k].lb, ck);
2✔
220
                    continue;
2✔
221
                }
222
                else if (fix_xk_to_upper)
30✔
223
                {
224
                    if (HAS_TAG(col_tags[k], C_TAG_UB_INF))
3✔
225
                    {
226
                        return UNBNDORINFEAS;
1✔
227
                    }
228

229
                    assert(!HAS_TAG(col_tags[k], C_TAG_INACTIVE));
2✔
230
                    assert(!IS_ABS_INF(bounds[k].ub));
2✔
231
                    fix_col(constraints, k, bounds[k].ub, ck);
2✔
232
                    continue;
2✔
233
                }
234

235
                // check if we can fix xj
236
                if (fix_xj_to_lower)
27✔
237
                {
238
                    if (HAS_TAG(col_tags[j], C_TAG_LB_INF))
×
239
                    {
240
                        return UNBNDORINFEAS;
×
241
                    }
242

NEW
243
                    assert(!HAS_TAG(col_tags[j], C_TAG_INACTIVE));
×
244
                    assert(!IS_ABS_INF(bounds[j].lb));
×
245
                    fix_col(constraints, j, bounds[j].lb, cj);
×
NEW
246
                    recount_ninfs = false;
×
247
                    //  break from checking if xj is parallel to other cols since
248
                    //  we fix it
UNCOV
249
                    break;
×
250
                }
251
                else if (fix_xj_to_upper)
27✔
252
                {
253
                    if (HAS_TAG(col_tags[j], C_TAG_UB_INF))
×
254
                    {
255
                        return UNBNDORINFEAS;
×
256
                    }
257

NEW
258
                    assert(!HAS_TAG(col_tags[j], C_TAG_INACTIVE));
×
259
                    assert(!IS_ABS_INF(bounds[j].ub));
×
260
                    fix_col(constraints, j, bounds[j].ub, cj);
×
NEW
261
                    recount_ninfs = false;
×
262
                    //  break from checking if xj is parallel to other cols since
263
                    //  we fix it
UNCOV
264
                    break;
×
265
                }
266
            }
267
        }
268

269
        // we might have to recount n_inf_min and n_inf_max for the rows column
270
        // j appears in due to bound change
271
        if (recount_ninfs)
23✔
272
        {
273
            for (jj = 0; jj < len_j; jj++)
52✔
274
            {
275
                rows_to_recompute[aj_cols[jj]] = 1;
38✔
276
            }
277
        }
278
    }
279

280
    return UNCHANGED;
16✔
281
}
282

283
static PresolveStatus process_all_bins(const Problem *prob, const int *parallel_cols,
20✔
284
                                       const iVec *groups, int *rows_to_recompute)
285
{
286
    if (groups->len <= 1)
20✔
287
    {
288
        return UNCHANGED;
4✔
289
    }
290

291
    Constraints *constraints = prob->constraints;
16✔
292
    size_t n_rows = constraints->A->m;
16✔
293
    memset(rows_to_recompute, 0, n_rows * sizeof(int));
16✔
294

295
    size_t n_groups = groups->len - 1;
16✔
296
    PresolveStatus status = UNCHANGED;
16✔
297

298
    for (size_t i = 0; i < n_groups; ++i)
34✔
299
    {
300
        int n_cols_this_group = groups->data[i + 1] - groups->data[i];
18✔
301
        status |= process_single_bin(prob, parallel_cols + groups->data[i],
18✔
302
                                     n_cols_this_group, rows_to_recompute);
303
    }
304

305
    // batch recompute n_infs for all affected rows
306
    const Matrix *A = constraints->A;
16✔
307
    Activity *acts = constraints->state->activities;
16✔
308
    ColTag *col_tags = constraints->col_tags;
16✔
309

310
    for (size_t i = 0; i < n_rows; ++i)
90✔
311
    {
312
        if (rows_to_recompute[i])
74✔
313
        {
314
            recompute_n_infs(acts + i, A->x + A->p[i].start, A->i + A->p[i].start,
23✔
315
                             A->p[i].end - A->p[i].start, col_tags);
23✔
316
        }
317
    }
318

319
    return status;
16✔
320
}
321

322
PresolveStatus remove_parallel_cols(Problem *prob)
20✔
323
{
324
    assert(prob->constraints->state->ston_rows->len == 0);
20✔
325
    assert(prob->constraints->state->empty_rows->len == 0);
20✔
326
    assert(prob->constraints->state->empty_cols->len == 0);
20✔
327
    DEBUG(verify_problem_up_to_date(prob->constraints));
20✔
328

329
    Constraints *constraints = prob->constraints;
20✔
330
    int *parallel_cols = constraints->state->work->iwork_n_cols;
20✔
331
    int *sparsity_IDs = constraints->state->work->iwork1_max_nrows_ncols;
20✔
332
    int *coeff_hashes = constraints->state->work->iwork2_max_nrows_ncols;
20✔
333
    iVec *group_starts = constraints->state->work->int_vec;
20✔
334

335
    // finding parallel rows of AT is the same as finding parallel cols of A
336
    find_parallel_rows(constraints->AT, constraints->col_tags, group_starts,
20✔
337
                       parallel_cols, sparsity_IDs, coeff_hashes, C_TAG_INACTIVE,
338
                       constraints->state->work->radix_aux);
20✔
339
#ifndef NDEBUG
340
    // there should be no inactive cols in group_starts here
341
    for (int i = 0; i < group_starts->len - 1; ++i)
38✔
342
    {
343
        for (int j = group_starts->data[i]; j < group_starts->data[i + 1]; ++j)
78✔
344
        {
345
            assert(
60✔
346
                !HAS_TAG(constraints->col_tags[parallel_cols[j]], C_TAG_INACTIVE));
347
            assert(constraints->state->col_sizes[parallel_cols[j]] !=
60✔
348
                   SIZE_INACTIVE_COL);
349
        }
350
    }
351
#endif
352

353
    int *rows_to_recompute = sparsity_IDs;
20✔
354
    PresolveStatus status =
355
        process_all_bins(prob, parallel_cols, group_starts, rows_to_recompute);
20✔
356
    assert(constraints->state->empty_rows->len == 0);
20✔
357

358
    delete_fixed_cols_from_problem(prob);
20✔
359
    delete_inactive_cols_from_A_and_AT(prob->constraints);
20✔
360

361
    DEBUG(verify_problem_up_to_date(constraints));
20✔
362

363
    return status;
20✔
364
}
STATUS · Troubleshooting · Open an Issue · Sales · Support · CAREERS · ENTERPRISE · START FREE · SCHEDULE DEMO
ANNOUNCEMENTS · TWITTER · TOS & SLA · Supported CI Services · What's a CI service? · Automated Testing

© 2026 Coveralls, Inc