Python 和 CSharp 中 Google OrTools 的不同结果

问题描述 投票:0回答:1

我需要您帮助解决有关 Google OrTools LinearSolver 的问题。我想在 python 和 CSharp 中计算一些数据。 状态如下:

我将两个 double[,] 和一个 double[] 数组赋予一个函数。然后我将数据转换为所需的格式并将其发送到求解器中。我在两种语言的函数中获得的数据是相同的,我逐行校对了。

有关数据的一些信息:

- fsi: A matrix with 1 on the main diagonal, otherwise 0.
- type(fsi): <class 'numpy.ndarray'> 
- fsi.shape: (300,300)
- type(zsl): <class 'numpy.ndarray'> 
- zsl.shape: (300, 322)
- type(ol):  <class 'numpy.ndarray'> 
- len(ol):   (322)

这是 python func 的代码:

def master_problem(fsi, zsl, ol, tlim=None, relaxation=True, enable_output=False):
    """
    Solve the master problem, either in its full version (MP)
    or in its relaxed version (RMP). Returns the following:
    - Objective value: minimization of sum(alpha[l] * h[l]), with h heights and l layer
    - Alpha values: alpha[l] represents layer selection
    - [RMP] Duals: one dual for each item
    """
    #logger.info("RMP defining variables and constraints")
    
    # Solver
    if relaxation:
        slv = pywraplp.Solver("RMP", pywraplp.Solver.GLOP_LINEAR_PROGRAMMING)
    else:
        slv = pywraplp.Solver("MP", pywraplp.Solver.BOP_INTEGER_PROGRAMMING)

    # Enable verbose output from solver
    if enable_output:
        slv.EnableOutput()

    # Utility
    infinity = slv.infinity()
    n_layers = 300
    n_items = fsi.shape[-1]

    # Variables
    if relaxation:
        al = [slv.NumVar(0, infinity, f"alpha_{l}") for l in range(n_layers)]
    else:
        al = [slv.BoolVar(f"alpha_{l}") for l in range(n_layers)]

    # Constraints
    constraints = []
    coefficients = np.matmul(fsi.T, zsl)

    # Select each item at least once
    # sum(al[l] * zsl[s, l] * fsi[s, i])
    for i in range(n_items):
        c = slv.Constraint(1, infinity, f"c_{i}")
        for l in range(n_layers):
            if coefficients[i, l] > 0:
                c.SetCoefficient(al[l], float(coefficients[i, l]))
        constraints += [c]

    # Objective
    obj = slv.Objective()
    for l, h in enumerate(ol):
        #logger.warning(str(ol) + " | " + str(h))
        obj.SetCoefficient(al[l], float(h))
    obj.SetMinimization()

    # Set a time limit in milliseconds
    if tlim is not None:
        slv.SetTimeLimit(1000 * tlim)

    # Solve
    #logger.debug(f"RMP variables: {slv.NumVariables()}")
    #logger.debug(f"RMP constraints: {slv.NumConstraints()}")
    status = slv.Solve()
    #logger.debug(f"RMP iterations: {slv.iterations()}")

    # Extract results
    duals, alphas = None, None
    objective = float("inf")
    logger.debug(f"RMP alphas: {alphas}")
    if status in (slv.OPTIMAL, slv.FEASIBLE):
        logger.info(f"RMP solved")

        # Extract alpha values
        alphas = [al[l].solution_value() for l in range(n_layers)]
        logger.debug(f"RMP alphas: {alphas}")
        if not all(alphas[l] in (0, 1) for l in range(n_layers)):
            logger.debug("RMP solution not feasible (at least one alpha value is not binary)")

        # Extract objective value
        objective = slv.Objective().Value()
        #logger.debug(f"RMP objective: {objective}")

        # Extract duals
        if relaxation:
            duals = np.array([c.DualValue() for c in constraints])
            #logger.debug(f"RMP duals: {duals}")
    else:
        logger.warning("RMP unfeasible")
    
    return objective, alphas, duals

以下是 C-Sharp 版本:

public static (double? objective, List<double>? alphas, List<double>? duals) MasterProblem(double[,]fsi, double[,]zsl, double[]ol, double? tlim = null, bool relaxation = true, bool enableOutput = false)
{
            // Prepare data
            var nLayers = layerPool.Count;
            var nItems = fsi.GetLength(0);

            /*Solve the master problem, either in its full version (MP)
            or in its relaxed version (RMP). Returns the following:
            - Objective value: minimization of sum(alpha[l] * h[l]), with h heights and l layer
            - Alpha values: alpha[l] represents layer selection
            - [RMP] Duals: one dual for each item*/

            // Create solver
            Solver solver;
            if (relaxation)
            {
                solver = Solver.CreateSolver("GLOP");
            }
            else
            {
                solver = Solver.CreateSolver("BOP");
            }

            // Enable verbose output from solver
            if (enableOutput)
            {
                solver.EnableOutput();
            }

            // Utility
            // infinity = slv.infinity() <- In CSharp anders gelöst. Wird direkt in Variable verwendet (z.B. double.PositiveInfinity)

            // Variables
            var al = new Variable[nLayers];
            if (relaxation)
            {
                for (int l = 0; l < nLayers; l++)
                {
                    al[l] = solver.MakeNumVar(0.0, double.PositiveInfinity, $"alpha_{l}");
                }
            }
            else
            {
                for (int l = 0; l < nLayers; l++)
                {
                    al[l] = solver.MakeBoolVar($"alpha_{l}");
                }

            }

            // Constraints
            var constraints = new List<Google.OrTools.LinearSolver.Constraint>();
            var fsiMatrix = Matrix<double>.Build.DenseOfArray(fsi);
            var zslMatrix = Matrix<double>.Build.DenseOfArray(zsl);
            var coefficients = fsiMatrix.Transpose() * zslMatrix;
          
            for (int i = 0; i < nItems; i++)
            {
                var c = solver.MakeConstraint(1, double.PositiveInfinity, $"c_{i}");
                for (int l = 0; l < nLayers; l++)
                {
                    if (coefficients[i, l] > 0)
                    {
                        c.SetCoefficient(al[l], coefficients[i, l]);
                    }
                }
                constraints.Add(c);
            }

            // Objective
            var objective = solver.Objective();
            for (int l = 0; l < ol.Count(); l++)
            {
                objective.SetCoefficient(al[l], ol[l]);
            }
            objective.SetMinimization();
            //var objective = solver.Objective();
            //for (int l = 0; l < nLayers; l++)
            //{
            //    objective.SetCoefficient(al[l], ol[l]);
            //}
            //objective.SetMinimization();

            // Set time limit
            if (tlim.HasValue)
            {
                solver.SetTimeLimit((long)(tlim.Value * 1000));
            }

            // Solve
            //Console.WriteLine($"RMP variables: {solver.NumVariables()}");
            //Console.WriteLine($"RMP constraints: {solver.NumConstraints()}");
            var status = solver.Solve();
            //Console.WriteLine($"RMP iterations: {solver.Iterations()}");

            // Extract results
            List<double> alphas = null;
            List<double> duals = null;
            double objectiveValue = double.PositiveInfinity;

            if (status == Solver.ResultStatus.OPTIMAL || status == Solver.ResultStatus.FEASIBLE)
            {
                //Console.WriteLine("RMP solved");

                // Extract alpha values
                alphas = al.Select(a => a.SolutionValue()).ToList();

                // Log alphas
                //Console.WriteLine($"RMP alphas: {string.Join(", ", alphas)}");

                // Check if all alphas are 0 or 1
                if (!alphas.All(alpha => alpha == 0 || alpha == 1))
                {
                    Console.WriteLine("RMP solution not feasible (at least one alpha value is not binary)");
                }

                // Extract objective value
                objectiveValue = objective.Value();

                // Extract duals
                if (relaxation)
                {
                    duals = constraints.Select(c => c.DualValue()).ToList();
                }
            }
            else
            {
                Console.WriteLine("RMP unfeasible");
            }

            return (objectiveValue, alphas, duals);
        }

“objectiveValue”、“alphas”和“duals”的结果在 C-Sharp 中与在 python 中应该不同。我真的不知道问题是什么。有任何想法吗? 如果您需要“zsl”或“ol”值的示例,请告诉我。

python c# matrix solver or-tools
1个回答
0
投票

很可能约束不是完全按照相同的顺序创建的。

MIP 求解器对约束、变量和项的顺序非常敏感。

© www.soinside.com 2019 - 2024. All rights reserved.