Python cy.CyClpSimplex类代码示例

    def test(self):
        model = CyLPModel()

        x = model.addVariable('x', 3)

        A = np.matrix([[1,2,3], [1,1,1]])
        b = CyLPArray([5, 3])

        model.addConstraint(A * x == b)
        model.addConstraint(x >= 0)

        model.objective = 1*x[0]  + 1*x[1] + 1.1 * x[2]

        # Solve it a first time
        s = CyClpSimplex(model)
        s.primal()
        sol = s.primalVariableSolution['x']
        self.assertTrue((abs(sol - np.array([1,2,0]) ) <= 10**-6).all())
        # Add a cut
        s.addConstraint(x[0] >= 1.1)
        s.primal()
        sol = s.primalVariableSolution['x']
        self.assertTrue((abs(sol - np.array([1.1, 1.8, 0.1]) ) <= 10**-6).all())

        # Change the objective function
        c = csr_matrixPlus([[1, 10, 1.1]]).T
        s.objective = c.T * x
        s.primal()
        sol = s.primalVariableSolution['x']
        self.assertTrue((abs(sol - np.array([2, 0, 1]) ) <= 10**-6).all())

    def test(self):
        m = CyCoinModel()
        
        m.addVariable(3, np.array(
                        [0, 1, 2], np.int32),
                        np.array([1., 1., 1.], np.double), 0, 10, 5)

        m.addVariable(2, np.array(
                        [1,2], np.int32),
                        np.array([5, 2.], np.double), 0, 10, 2)

        # Add bound for the three constraints (we have two variables)
        m.setConstraintLower(0, 2.3)
        m.setConstraintLower(1, 4.5)
        m.setConstraintLower(0, 1.5)
        
        # Add a 4th constraint
        m.addConstraint(2, 
                            np.array([0, 1], np.int32), 
                            np.array([1., 1.], np.double), 2, 7)
        
        s = CyClpSimplex()
        # Load the problem from the CyCoinModel
        s.loadProblemFromCyCoinModel(m)
        
        s.primal()
         
        self.assertAlmostEqual(s.objectiveValue, 8.7, 7)

    def test_Sparse(self):
        model = CyLPModel()

        x = model.addVariable('x', 3)
        y = model.addVariable('y', 2)

        A = csc_matrixPlus(([1, 2, 1, 1], ([0, 0, 1, 1], [0, 1, 0, 2])),  shape=(2, 3))
        B = csc_matrixPlus(([1, 1], ([0, 1], [0, 2])),  shape=(2, 3))
        D = np.matrix([[1., 2.],[0, 1]])
        a = CyLPArray([5, 2.5])
        b = CyLPArray([4.2, 3])
        x_u= CyLPArray([2., 3.5])

        model.addConstraint(A*x <= a)
        model.addConstraint(2 <= B * x + D * y <= b)
        model.addConstraint(y >= 0)
        model.addConstraint(1.1 <= x[1:3] <= x_u)

        c = CyLPArray([1., -2., 3.])
        model.objective = c * x + 2 * y[0] + 2 * y[1]


        s = CyClpSimplex(model)

        s.primal()
        sol = np.concatenate((s.primalVariableSolution['x'],
                              s.primalVariableSolution['y']))
        self.assertTrue((abs(sol -
                        np.array([0.2, 2, 1.1, 0, 0.9]) ) <= 10**-6).all())

    def QPModel(self, addW=False):
        A = self.A
        c = self.c
        s = CyClpSimplex()
        
        x = s.addVariable('x', self.nCols)
        if addW:
            w = s.addVariable('w', self.nCols)
        
        s += A * x >= 1
        n = self.nCols
        
        if not addW:
            s += 0 <= x <= 1
        else:
            s += x + w == 1
            s += 0 <= w <= 1

##        s += -1 <= x <= 1

        s.objective = c * x
       
        if addW:
            G = sparse.lil_matrix((2*n, 2*n))
            for i in xrange(n/2, n): #xrange(n-1):
                G[i, i] = 1
            G[2*n-1, 2*n-1] = 10**-10
        else:
            G = sparse.lil_matrix((n, n))
            for i in xrange(n/2, n): #xrange(n-1):
                G[i, i] = 1

    
        s.Hessian = G
        return s

def solve(filename, method):

    s = CyClpSimplex()
    s.readMps(filename)

    s.preSolve(feasibilityTolerance=10 ** -8)
    #s.useCustomPrimal(1)

    if method == 'd':
        pivot = DantzigPivot(s)
    elif method == 'l':
        pivot = LIFOPivot(s)
    elif method == 'm':
        pivot = MostFrequentPivot(s)
    elif method == 'p':
        pivot = PositiveEdgePivot(s)
    else:
        print 'Unkown solution method.'
        sys.exit(1)

    s.setPivotMethod(pivot)

    #s.setPerturbation(50)

    start = clock()
    s.primal()
    print 'Problem solved in %g seconds.' % (clock() - start)
    return s.objectiveValue

    def model(self):
        A = self.A
        c = self.c
        s = CyClpSimplex()
        
        x = s.addVariable('x', self.nCols)

        s += A * x >= 1
        s += 0 <= x <= 1
        
        s.objective = c * x

        return s

    def test_multiDim(self):
        from cylp.cy import CyClpSimplex
        from cylp.py.modeling.CyLPModel import CyLPArray
        s = CyClpSimplex()
        x = s.addVariable('x', (5, 3, 6))
        s += 2 * x[2, :, 3].sum() + 3 * x[0, 1, :].sum() >= 5

        s += 0 <= x <= 1
        c = CyLPArray(range(18))

        s.objective = c * x[2, :, :] + c * x[0, :, :]
        s.primal()
        sol = s.primalVariableSolution['x']
        self.assertTrue(abs(sol[0, 1, 0] - 1) <= 10**-6)
        self.assertTrue(abs(sol[2, 0, 3] - 1) <= 10**-6)

class TestCyCoinIndexedVector(unittest.TestCase):

    def setUp(self):
        self.a = np.array([1, 10.5, -11.3, 100, -50.5, 20], dtype=np.double)
        self.a2 = np.array([1000, 10.5, -11.3, 100, -50.5, 20], dtype=np.double)
        self.s = CyClpSimplex()
    
#    def test_gen(self):
#        w = np.array([3], dtype=np.int32)
#        self.assertEqual(self.s.argWeightedMax(self.a, w, 0.1), 5)
#    
#    def test_empty(self):
#        w = np.array([], dtype=np.int32)
#        self.assertEqual(self.s.argWeightedMax(np.array([]), w, 0.1), 0)
#    
#    def test_first(self):
#        w = np.array([0, 2], dtype=np.int32)
#        self.assertEqual(self.s.argWeightedMax(self.a, w, 99), 3)
#        self.assertEqual(self.s.argWeightedMax(self.a, w, 100), 0)

    def test_argMax4_1(self):
        w_ind = np.array([0, 2, 5], dtype=np.int32)
        self.assertEqual(self.s.argWeightedMax(self.a, 0, 100, w_ind), 5)

    def test_argMax4_2(self):
        w_ind = np.array([0, 2, 5], dtype=np.int32)
        w = np.array([1.5, -10, 4], dtype=np.double)
        self.assertEqual(self.s.argWeightedMax(self.a, 0, w, w_ind), 2)
    
    def test_argMax4_3(self):
        w_ind = np.array([0, 8, 21], dtype=np.int32)
        a_ind = np.array([2, 5, 8, 10, 20, 21] , dtype=np.int32)
        self.assertEqual(self.s.argWeightedMax(self.a, a_ind, 5.1, w_ind), 5)

    def test_argMax4_4(self):
        w_ind = np.array([0, 8, 21], dtype=np.int32)
        w = np.array([100, -10, 4], dtype=np.double)
        a_ind = np.array([2, 5, 8, 10, 20, 21] , dtype=np.int32)
        self.assertEqual(self.s.argWeightedMax(self.a, a_ind, w, w_ind), 2)
    
    def test_argMax_5(self):
        w_ind = np.array([2, 7, 100], dtype=np.int32)
        w = np.array([100, -10, 4], dtype=np.double)
        a_ind = np.array([2, 5, 8, 10, 20, 21] , dtype=np.int32)
        self.assertEqual(self.s.argWeightedMax(self.a2, a_ind, 10, w_ind), 0)

    def test_onlyBounds2(self):
        s = CyClpSimplex()

        x = s.addVariable('x', 3)
        y = s.addVariable('y', 2)

        s += y >= 1
        s += 2 <= x <= 4

        c = CyLPArray([1., -2., 3.])
        s.objective = c * x + 2 * y[0] + 2 * y[1]

        s.primal()

        sol = np.concatenate((s.primalVariableSolution['x'],
                              s.primalVariableSolution['y']))
        self.assertTrue((abs(sol -
                        np.array([2, 4, 2, 1, 1]) ) <= 10**-6).all())

    def test_ArrayIndexing(self):
        from cylp.cy import CyClpSimplex
        from cylp.py.modeling.CyLPModel import CyLPArray
        s = CyClpSimplex()
        x = s.addVariable('x', (5, 3, 6))
        s += 2 * x[2, :, 3].sum() + 3 * x[0, 1, :].sum() >= 5


        s += x[1, 2, [0, 3, 5]] - x[2, 1, np.array([1, 2, 4])] == 1
        s += 0 <= x <= 1
        c = CyLPArray(range(18))

        s.objective = c * x[2, :, :] + c * x[0, :, :]
        s.primal()
        sol = s.primalVariableSolution['x']
        self.assertTrue(abs(sol[1, 2, 0] - 1) <= 10**-6)
        self.assertTrue(abs(sol[1, 2, 3] - 1) <= 10**-6)
        self.assertTrue(abs(sol[1, 2, 5] - 1) <= 10**-6)

    def test_removeVar(self):
        m = self.model
        x = self.x
        A = self.A
        B = self.B
        D = self.D
        b = self.b

        y = m.addVariable('y', 4)
        z = m.addVariable('z', 5)

        m.addConstraint(x >= 0)
        m.addConstraint(y >= -10)
        m.addConstraint(z >= -10)

        m.addConstraint(A * x + D * y + B * z <= b)
        m += x[0] + y[0] + z[0] >= 1.12

        m.objective = x.sum() + y.sum() + z.sum()
        s = CyClpSimplex(m)
        s.primal()

        self.assertTrue('y' in s.primalVariableSolution.keys())

        m.removeVariable('y')
        s = CyClpSimplex(m)
        s.primal()
        self.assertTrue('y' not in s.primalVariableSolution.keys())

    def test_variableBoundSubset(self):
        m = self.model
        x = self.x
        y = m.addVariable('y', 4)
        z = m.addVariable('z', 5)
        A = self.A
        b = self.b

        k = m.addVariable('k', 2)

        s = CyClpSimplex(m)

        s.setColumnLowerSubset(np.array([1, 2], np.int32), np.array([3, 5, 8],
                                        np.int32), np.array([3.2, 3.1, 2.2]))
        self.assertTrue(s.variablesLower[3] != 3.2)
        self.assertTrue(s.variablesLower[5] == 3.1)
        self.assertTrue(s.variablesLower[8] == 2.2)

        s.setColumnUpperSubset(np.array([0, 2], np.int32), np.array([0, 4, 10],
                              np.int32), np.array([3.2, 3.1, 2.2]))
        self.assertTrue(s.variablesUpper[0] == 3.2)
        self.assertTrue(s.variablesUpper[4] != 3.1)
        self.assertTrue(s.variablesUpper[10] == 2.2)

def read_instance(module_name = True, file_name = None):

    if module_name:
        lp = CyClpSimplex()

        mip = ilib.import_module(module_name)
            
        A = np.matrix(mip.A)
        #print np.linalg.cond(A)
        b = CyLPArray(mip.b)
        
        #We assume variables have zero lower bounds
        x_l = CyLPArray([0 for _ in range(mip.numVars)])
            
        x = lp.addVariable('x', mip.numVars)
        
        lp += x >= x_l
        try:
            x_u = CyLPArray(getattr(mip, 'x_u'))
            lp += x <= x_u
        except:
            pass
        
        lp += (A * x <= b if mip.sense[1] == '<=' else
               A * x >= b)
        c = CyLPArray(mip.c)
        lp.objective = -c * x if mip.sense[0] == 'Max' else c * x
        
        return lp, x, mip.A, mip.b, mip.sense, mip.integerIndices
    else:
        #TODO Change sense of inequalities so they are all the same
        #     by explicitly checking lp.constraintsUpper and lp.constraintsLower
        #Warning: Reading MP not well tested 
        lp.extractCyLPModel(file_name)
        x = lp.cyLPModel.getVarByName('x')
        sense = ('Min', '>=')
        return lp, x, None, None, sense, integerIndices

    def test_1(self):
        """simplest QP test"""
        s = CyClpSimplex()
        s.readMps(join(currentFilePath, '../input/hs35.qps'))
        #self.assertTrue(abs(cbcModel.objectiveValue - 3089.0) < 10 ** -6)
        
        #print s.Hessian.todense()

        p = WolfePivot(s)
        s.setPivotMethod(p)

        s.primal()
        print s.primalVariableSolution
        print s.objectiveValue

def read_instance(module_name = None, file_name = None):

    if module_name is not None:
        lp = CyClpSimplex()

        mip = ilib.import_module(module_name)
            
        A = np.matrix(mip.A)
        #print np.linalg.cond(A)
        b = CyLPArray(mip.b)
        
        #Warning: At the moment, you must put bound constraints in explicitly for split cuts
        x_l = CyLPArray([0 for _ in range(mip.numVars)])
            
        x = lp.addVariable('x', mip.numVars)
        
        lp += x >= x_l
        try:
            x_u = CyLPArray(getattr(mip, 'x_u'))
            lp += x <= x_u
        except:
            pass
        
        lp += (A * x <= b if mip.sense[1] == '<=' else
               A * x >= b)
        c = CyLPArray(mip.c)
        lp.objective = -c * x if mip.sense[0] == 'Max' else c * x
        return lp, x, mip.A, mip.b, mip.sense[1], mip.integerIndices
    elif file_name is not None:
        lp = CyClpSimplex()
        m = lp.extractCyLPModel(file_name)
        x = m.getVarByName('x')
        integerIndices = [i for (i, j) in enumerate(lp.integerInformation) if j == True]
        infinity = lp.getCoinInfinity()
        sense = None
        for i in range(lp.nRows):
            if lp.constraintsLower[i] > -infinity:
                if sense == '<=':
                    print "Function does not support mixed constraint..."
                    break
                else: 
                    sense = '>='
                    b = lp.constraintsLower
            if lp.constraintsUpper[i] < infinity: 
                if sense == '>=':
                    print "Function does not support mixed constraint..."
                    break
                else: 
                    sense = '<='
                    b = lp.constraintsUpper
        return lp, x, lp.coefMatrix, b, sense, integerIndices
    else:
        print "No file or module name specified..."
        return None, None, None, None, None, None