Search code examples
linear-programmingscip

Best way to obtain the LP relaxation of a MIP using SCIP


I've a C++ app which can use indistinctively CPLEX or COIN, and now I want to add SCIP as a choice. As a part of the functionality, it's required to run the LP relaxation of a MIP problem. In CPLEX or COIN that's simple, I just create the problem and then solve the LP using:

  • CPLEX: int const rval = CPXXlpopt(m_env, m_lp);
  • COIN: m_modelIP.solver()->getModelPtr().initialSolve(m_lpOptions);

And the solution to that LP is in the same place the solution to the regular MIP would be, there's no "special" nor "separate" solution object I need to look at, just the regular one.

However, I haven't been able to figure out the equivalent using SCIP. So far, I've seen these two possible solutions to achieve that:

    int nVars = SCIPgetNVars(m_scip);
    std::vector<SCIP_VARTYPE> varTypes(nVars, SCIP_VARTYPE_CONTINUOUS);
    SCIP_VAR **pVars = SCIPgetVars(m_scip);
    
    /* Move all binary and integral variables to continuous */
    for (size_t i = 0; i < static_cast<size_t>(nVars); i++)
    {
        SCIP_VARTYPE const varType = SCIPvarGetType(pVars[i]);
        if (varType == SCIP_VARTYPE_CONTINUOUS)
            continue;

        varTypes[i] = varType;
    
        SCIP_Bool infeasible;
        if (SCIPchgVarType(m_scip, pVars[i], SCIP_VARTYPE_CONTINUOUS, &infeasible) != SCIP_OKAY)
            throw std::runtime_error("Error changing integrality type to SCIP variable " + std::to_string(i) + " to continuos to solve LP relaxation");
    }

    /* Store parameters current value */
    SCIP_Longint oldNodes;
    if (SCIPgetLongintParam(m_scip, "limits/nodes", &oldNodes) != SCIP_OKAY)
        throw std::runtime_error("Couldn't get limits/nodes parameter in SCIP");

    int oldMaxRounds;
    if (SCIPgetIntParam(m_scip, "presolving/maxrounds", &oldMaxRounds) != SCIP_OKAY)
        throw std::runtime_error("Couldn't get presolving/maxrounds parameter in SCIP");
        
    /* Change parameters to create the "LP-like" scenario */
    if (SCIPsetLongintParam(m_scip, "limits/nodes", 1) != SCIP_OKAY)
        throw std::runtime_error("Couldn't set limits/nodes parameter in SCIP");

    if (SCIPsetIntParam(m_scip, "presolving/maxrounds", 0) != SCIP_OKAY)
        throw std::runtime_error("Couldn't set presolving/maxrounds parameter in SCIP");

    /* Solve the supposedly now LP problem */
    if (SCIPsolve(m_scip) != SCIP_OKAY)
        throw std::runtime_error("Couldn't solve LP problem in SCIP");

    /* Restore the original state of the parameters */
    if (SCIPsetLongintParam(m_scip, "limits/nodes", oldNodes) != SCIP_OKAY)
        throw std::runtime_error("Couldn't restore limits/nodes parameter in SCIP");

    if (SCIPsetIntParam(m_scip, "presolving/maxrounds", oldMaxRounds) != SCIP_OKAY)
        throw std::runtime_error("Couldn't restore presolving/maxrounds parameter in SCIP");

    /* Restore the original state of the variables */
    for (size_t i = 0; i < static_cast<size_t>(nVars); i++)
    {
        SCIP_VARTYPE const varType = varTypes[i];
        if (varType == SCIP_VARTYPE_CONTINUOUS)
            continue;

        SCIP_Bool infeasible;
        if (SCIPchgVarType(m_scip, pVars[i], varType, &infeasible) != SCIP_OKAY)
            throw std::runtime_error("Error restoring integrality type to SCIP variable " + std::to_string(i) + " from continuos after solving LP relaxation");
    }
  • Method 2: Use the LP interface of SCIP.
    /* Presolve the original MIP problem. Without this step, the following SCIPgetLPI() call throws a SIGSEGV */
    if (SCIPpresolve(m_scip) != SCIP_OKAY)
        throw std::runtime_error("Couldn't presolve SCIP");

    /* Access the LP interface of SCIP */
    SCIP_LPI *pLpi = nullptr;
    if (SCIPgetLPI(m_scip, &pLpi) != SCIP_OKAY)
        throw std::runtime_error("Couldn't obtain access to LP interface of SCIP");

    /* And solve it */
    if (SCIPlpiSolvePrimal(pLpi) != SCIP_OKAY)
        throw std::runtime_error("Couldn't solve LP interface of SCIP");

Although method 1 seems to work, I have the (uncorroborated) feeling that method 2 is more correct. However, I'm not sure that the method used to solve the problem (SCIPlpiSolvePrimal()) is indeed the correct one, or if the solution provided by that method is located in such a place that I can simply retrieve it using SCIPgetBestSol() (My library doesn't have separated methods to obtain the IP or the LP solution, since neither CPLEX nor COIN really need it).

SO my main question is: Is my approach #1 correct? Or is it too slow/not the most appropiate one, and instead #2 is the one that should be used instead? If so, is the solution provided by SCIPlpiSolvePrimal() accessible from the regular MIP interface? Or if not, is there a way to move it into that?

Or is there even a third way to tackle this case? I've read here that maybe using SCIP*Dive() is a better approach, but I haven't found any samples regarding that.

Thanks a lot for your help.


Solution

  • If it is possible for you to execute your code inside a SCIP plugin, then that would be the most efficient/elegant way to do it. (You could, e.g. add an event handler that is executed after the first LP solve and just gets the LP solution). This is also where you could use the diving mode to change some aspects of the problem, just solve the (changed) LP and then resume solving the original problem. You cannot do without using plugins, because you need to be in the solving stage to execute diving mdoe. Method 2 that you mentioned has the same issue, the LPI will only have the correct information if you are already in solving mode. Of course you could create a new LPI struct and basically copy the whole problem, but that is probably too inefficient if you are solving large problems.

    If you need your code structure to be the same as it is (call LP solve, do something with the information, continue solving) you are stuck with option 1 for now (the changing of variable types will not be the performance problem, but you will essentially restart the solve after changing the variable types).

    Unfortunately, SCIP currently does not have the functionality to just solve the relaxation and then stop after, with the option of continuing the whole MIP solve.

    I hope this helps, let me know if you something with my answer is not clear.