diff --git a/CMakeLists.txt b/CMakeLists.txt index 41984b28..d6a7da31 100644 --- a/CMakeLists.txt +++ b/CMakeLists.txt @@ -3,7 +3,7 @@ project(PSLP LANGUAGES C CXX) # The library is written in C, but we also have so set(PSLP_VERSION_MAJOR 0) set(PSLP_VERSION_MINOR 0) -set(PSLP_VERSION_PATCH 10) +set(PSLP_VERSION_PATCH 11) set(PSLP_VERSION "${PSLP_VERSION_MAJOR}.${PSLP_VERSION_MINOR}.${PSLP_VERSION_PATCH}") add_compile_definitions(PSLP_VERSION="${PSLP_VERSION}") diff --git a/include/PSLP/PSLP_API.h b/include/PSLP/PSLP_API.h index 40f08362..b065965d 100644 --- a/include/PSLP/PSLP_API.h +++ b/include/PSLP/PSLP_API.h @@ -124,6 +124,14 @@ extern "C" presolver contains the presolved problem. */ PresolveStatus run_presolver(Presolver *presolver); + /* Free the presolved problem data owned by 'presolver' after the caller has + copied all arrays in 'reduced_prob'. This does not free 'presolver'; the + caller must still call 'free_presolver' after postsolve. Normal postsolve + and dual-infeasibility-ray postsolve remain available. Primal- + infeasibility-ray postsolve is no longer available because it requires + the reduced constraint matrix. */ + void free_presolver_reduced_problem(Presolver *presolver); + /* Postsolve the problem given the primal-dual solution (x, y, z) of the reduced problem. The function populates presolver->sol, so if you're looking for the solution to the original problem, you want to look there. diff --git a/src/core/Presolver.c b/src/core/Presolver.c index c5640b65..c0bdc9e7 100644 --- a/src/core/Presolver.c +++ b/src/core/Presolver.c @@ -736,6 +736,12 @@ PresolveStatus run_presolver(Presolver *presolver) stats->nnz_reduced = A->nnz; DEBUG(run_debugger_stats_consistency_check(stats)); populate_presolved_problem(presolver); + + // The transpose is only needed during presolve. Postsolve uses the recorded + // reductions and row/column mappings. + free_matrix(prob->constraints->AT); + prob->constraints->AT = NULL; + clock_gettime(CLOCK_MONOTONIC, &outer_timer.end); stats->time_presolve = GET_ELAPSED_SECONDS(outer_timer); @@ -747,6 +753,50 @@ PresolveStatus run_presolver(Presolver *presolver) return final_status(stats); } +void free_presolver_reduced_problem(Presolver *presolver) +{ + + // 'Ap' is non-null when 'run_presolver' has populated the reduced problem. + // Bailing out here keeps this function safe to call before a successful + // presolve. + if (!presolver || !presolver->reduced_prob || !presolver->reduced_prob->Ap) + { + return; + } + + PresolvedProblem *reduced_prob = presolver->reduced_prob; + Problem *prob = presolver->prob; + if (prob) + { + Constraints *constraints = prob->constraints; + if (constraints) + { + free_matrix(constraints->A); + constraints->A = NULL; + free_matrix(constraints->AT); + constraints->AT = NULL; + PS_FREE(constraints->lhs); + PS_FREE(constraints->rhs); + PS_FREE(constraints->bounds); + PS_FREE(constraints->row_tags); + PS_FREE(constraints->col_tags); + } + if (prob->obj) + { + PS_FREE(prob->obj->c); + } + } + + reduced_prob->Ax = NULL; + reduced_prob->Ai = NULL; + reduced_prob->lhs = NULL; + reduced_prob->rhs = NULL; + reduced_prob->c = NULL; + PS_FREE(reduced_prob->Ap); + PS_FREE(reduced_prob->lbs); + PS_FREE(reduced_prob->ubs); +} + void postsolve(Presolver *presolver, const double *x, const double *y, const double *z) { @@ -798,6 +848,14 @@ void postsolve_primal_infeas_ray(Presolver *presolver, const double *y, size_t len_z = MAX((size_t) 1, stats->n_cols_reduced); double *z; assert(reduced_prob != NULL); + if (!reduced_prob || !reduced_prob->Ax) + { + // The reduced constraint matrix was released by + // 'free_presolver_reduced_problem', so the ray cannot be recovered. + fprintf(stderr, "PSLP warning: postsolve_primal_infeas_ray called after " + "free_presolver_reduced_problem! Don't do this! \n"); + return; + } z = (double *) ps_malloc(len_z, sizeof(double)); assert(z != NULL); diff --git a/tests/test_postsolve.h b/tests/test_postsolve.h index f928363d..97b548f3 100644 --- a/tests/test_postsolve.h +++ b/tests/test_postsolve.h @@ -46,12 +46,26 @@ static char *test_0_postsolve() Constraints *constraints = prob->constraints; Matrix *A = constraints->A; - run_presolver(presolver); + // calling before presolve must be a no-op that keeps the original problem + free_presolver_reduced_problem(presolver); + mu_assert("free before presolve must not touch the problem", + constraints->A != NULL && constraints->AT != NULL); + run_presolver(presolver); + mu_assert("transpose should be released after presolve", + constraints->AT == NULL); Mapping *maps = prob->constraints->state->work->mappings; int *rows_map = maps->rows; int *cols_map = maps->cols; + // should be able to call it multiple times + free_presolver_reduced_problem(presolver); + free_presolver_reduced_problem(presolver); + mu_assert("reduced matrix should be released", constraints->A == NULL); + mu_assert("reduced problem buffers should be released", + presolver->reduced_prob->Ax == NULL && + presolver->reduced_prob->Ap == NULL); + // construct optimal primal solution to reduced problem (computed offline) double x[] = {10., -10., -10., 2.71428571, -10., -6.85714286, -10.}; double y[] = {-2.57142857, 0.14285714}; diff --git a/tests/test_ray_postsolve.h b/tests/test_ray_postsolve.h index 1c9ddaaa..4e665333 100644 --- a/tests/test_ray_postsolve.h +++ b/tests/test_ray_postsolve.h @@ -242,6 +242,7 @@ static char *test_0_dual_ray_postsolve() PresolveStatus status = run_presolver(presolver); mu_assert("dual ray postsolve presolve status", status == UNCHANGED); + free_presolver_reduced_problem(presolver); double x[] = {0.0, 0.0}; double x_orig[] = {0.0, 0.0};