def knapsack(self, f, d, b, quantities):
     """
     max z: f1X1 + ... + frXr
            d1X1 + ... + frXr <= b
            X1 .. Xr >=0, integer
            
     @param f, list of parameters to be maximized
     @param d, list of objective parameters
     @param b, int boundary of the objective
     @return (x, z)
              x list of values
              z, the maximized value
     """
     problem = pulp.LpProblem("Knapsakc", pulp.LpMaximize)
     
     nrCols = len(f)       
     
     x = []
     for r in range(nrCols):
         # Create variables Xi, int, >=0
         x.append(pulp.LpVariable("x%d"%r , 0, quantities[r], pulp.LpInteger))
     
     problem += sum( d[r] * x[r] for r in range(nrCols)) <= b
     problem += sum( f[r] * x[r] for r in range(nrCols))
     
     #status = problem.solve(pulp.GLPK(msg = 0))
     #problem.writeLP("/tmp/knapsack.lp")
     status = problem.solve()
     if self.LOG:
         print problem
        
     return ([pulp.value(a) for a in x], pulp.value(problem.objective))

def LPboost_solver():
    d = pp.LpVariable.dicts("d", range(3), 0, 1)
    prob = pp.LpProblem("LPboost", pp.LpMinimize)
    prob += pp.lpSum(d) == 1
    print pp.LpStatus[prob.solve()]
    for i in range(3):
        print pp.value(d[i])

def get_optimal_allocations(domain, criterion, batch_size, method):
    candidates = get_candidates(domain, criterion, method, limit=1000)
    n = len(candidates)
    ilp_vars = [pulp.LpVariable('x%s' % x, 0, batch_size, pulp.LpInteger) \
                  for x in range(n)]
    prices = [candidates[i]['price'] for i in range(n)]
    confs = [candidates[i]['conf'] for i in range(n)]
    dists = [candidates[i]['dist'] for i in range(n)]
    histo = get_availability_histo(domain)
    ilp = None
    if method == 'min_price':
        ilp = pulp.LpProblem('minprice', pulp.LpMinimize)
        ilp += pulp.lpDot(prices, ilp_vars)
        ilp += pulp.lpDot(confs, ilp_vars) >= decimal.Decimal(batch_size * criterion)
    elif method == 'max_conf':
        ilp = pulp.LpProblem('maxconf', pulp.LpMaximize)
        ilp += pulp.lpDot(confs, ilp_vars)
        ilp += pulp.lpDot(prices, ilp_vars) <= decimal.Decimal(criterion)
    ilp += pulp.lpSum(ilp_vars) == batch_size
    for level in histo:
        ilp += pulp.lpDot([dists[i][level - 1] for i in range(n)], \
                          ilp_vars) <= histo[level]
    status = ilp.solve(pulp.GLPK(msg=0))
    allocs = []
    if pulp.LpStatus[status] == 'Optimal':
        for i in range(n):
            x = ilp_vars[i]
            if pulp.value(x) > 0:
                allocs.append((pulp.value(x), candidates[i]))
    return allocs

def vector_solver():
    d = pp.LpVariable.dicts("x", range(2), 0, 1)
    prob = pp.LpProblem("myV", pp.LpMinimize)
    prob += d[0] + d[1] == 1
    print pp.LpStatus[prob.solve()]
    print pp.value(d[0])
    print pp.value(d[1])

def solveProblemWithRestr(weigthLimit):
	prob = pulp.LpProblem("myProblem", pulp.LpMaximize)
	
	# Only one of each piece
	prob += sum([a.var for a in arms]) <= 1.0
	prob += sum([a.var for a in chests]) <= 1.0
	prob += sum([a.var for a in heads]) <= 1.0
	prob += sum([a.var for a in legs]) <= 1.0

	# Maximum Weigth
	restr = sum(getRestr(armors))
	prob += restr <= weigthLimit

	obj = sum(getObj(armors))
	prob += obj


	status = prob.solve()
	
	print("| Name | Armor Rating |  Weigth |")
	print("|:-----|------:|------:|")
	
	for a in armors:
		if pulp.value(a.var) > 0:
			print("|", a.name, "|", a.ar, "|", a.w, "|")

	print("| Total | ", pulp.value(obj), "|", pulp.value(restr), "|")

def print_model(model):
    print "STATUS: ", LpStatus[model.status]
    print "TOTAL COST: ", value(model.total_cost)
    print "FUEL COSTS: ", value(model.fuel_costs)
    print "STOP COSTS: ", value(model.stop_costs)
    print "CONTRACT COSTS: ", value(model.contract_costs)

    for yard in model.d.YARDS:
        print "Yard ", yard, value(model.v_contract[yard])

def basic_pulp_solver():
    x = pp.LpVariable("x", 0, 3)
    y = pp.LpVariable("y", 0, 1)
    prob = pp.LpProblem("myP", pp.LpMinimize)
    prob += x + y <= 2
    prob += -4*x + y
    status = prob.solve()
    print pp.LpStatus[status]
    print pp.value(x), pp.value(y)

def main():
    x = LpVariable("x", 0, 3)
    y = LpVariable("y", 0, 1)
    prob = LpProblem("myProblem", LpMinimize)
    prob += x + y <= 2
    prob += -4*x + y
    status = prob.solve(COIN(msg = 0))
    print LpStatus[status]
    print value(x)
    print value(y)

 def get_true_false_classes():
     true_set  = set()
     false_set = set()
     
     for rep, size in classes.items():
         if value(get_var(rep)) == 0:
             false_set.add(rep)
         elif value(get_var(rep)) == size:
             true_set.add(rep)
     return true_set, false_set

    def print_results (self):
        print '<score>%s</score>' % (pulp.value(self.score),)
        print '<total>%s</total>' % fmt_stats(map (lambda x: pulp.value(x), self.total_stats))
        print '<base>%s</base>' % fmt_stats(self.base_stats)

        print '<items>'
        for item_list in self.items.values():
            for item in item_list:
                if pulp.value(self.used[item]) == 1:
                    item.print_results()
        print '</items>'

    def test_expanded_ToyMMB(self):
        self.pkn.compress()
        self.pkn.expand_and_gates()

        model = cno.milp.model.MILPTrain(self.pkn, self.midas)
        model.train()

        eps = 1e-12
        model_error = pulp.value(model.error_objective_expression())
        model_size = pulp.value(model.size_objective_expression())

        assert abs(model_error-10.02) < eps
        assert abs(model_size-9) < eps

 def report_copy_map(self):
     """
     Reports the raw counts seen at each variable. This is normalized by the number of variables in the block.
     """
     copy_map = defaultdict(list)
     for para in self.block_map:
         for i in xrange(len(self.block_map[para])):
             start, stop, var, block = self.block_map[para][i]
             if var is not None:
                 copy_map[para].append([start, stop, pulp.value(var)])
                 prev_var = pulp.value(var)
             else:
                 copy_map[para].append([start, stop, prev_var])
     return copy_map

	def calculate_V(self, s, agent_name):
		""" Calculate V for state s by maximizing it while minimizing opponents actions. Returns the maximium value of V
		"""
		max_v = pulp.LpProblem("Maximize V",  pulp.LpMaximize)

		# Set V as variable to maximize
		v  = pulp.LpVariable("v", 0.0, cat="Continuous")
		max_v += v
		actions = ['West', 'East','North', 'South','Wait']
		# Create policy var for actions
		action_policy_vars = pulp.LpVariable.dicts("A",actions,lowBound =0.0, upBound = 1.0, cat="Continuous")

		# Probabilities sum to 1
		max_v += sum([action_policy_vars[a] for a in actions]) == 1
		for a in actions:
			max_v += action_policy_vars[a] >= 0.000000001 

		# add constraints as summation of actions given an opponent action are bigger than 0
		for o in actions:
			max_v += sum([self.get_qvalue_minimax(s, agent_name, a, o) * action_policy_vars[a] for a in actions]) >= v

		# Solve maximization
		max_v.solve()
		#for i in actions:
		#	if action_policy_vars[i].value() == 1.0:
		#		print i


		return pulp.value(max_v.objective)

def opt(C, X):
    orderNum = len(X)
    routeNum = len(C)
    routeIdx = range(routeNum)
    orderIdx = range(orderNum)
    # print routeIdx,orderIdx
    eps = 1.0 / 10 ** 7
    print eps
    var_choice = lp.LpVariable.dicts('route', routeIdx, cat='Binary')
    # var_choice=lp.LpVariable.dicts('route',routeIdx,lowBound=0)#尝试松弛掉01变量
    exceed_labor = lp.LpVariable('Number of routes exceed 1000', 0)
    prob = lp.LpProblem("lastMile", lp.LpMinimize)

    prob += exceed_labor * 100000 + lp.lpSum(var_choice[i] * C[i] for i in routeIdx)

    prob += lp.lpSum(var_choice[i] for i in routeIdx) <= 1000 + exceed_labor + eps
    for i in orderIdx:
        prob += lp.lpSum(var_choice[j] for j in X[i]) >= (1 - eps)

    prob.solve(lp.CPLEX(msg=0))
    print "\n\nstatus:", lp.LpStatus[prob.status]
    if lp.LpStatus[prob.status] != 'Infeasible':
        obj = lp.value(prob.objective)
        print "\n\nobjective:", obj
        sol_list = [var_choice[i].varValue for i in routeIdx]
        print "\n\nroutes:", (sum(sol_list))
        # print "\n\noriginal problem:\n",prob
        return obj, sol_list, lp.LpStatus[prob.status]
    else:
        return None, None, lp.LpStatus[prob.status]

    def export_res(self,queue_b):
        n = self.n
        tn = self.tn

        for i in xrange(n):
            for t in xrange(tn):
                queue_b[i].b[t] = int(math.ceil(pulp.value(self.varb[str(i)][str(t)])))

    def train(self):
        """Initialize and solve the optimization problem.

        The problem is solved in two steps. First, the algorithm searches for
        a network with the best possible fit. Second, the algorithm searches
        for the smallest possible network that returns the best fit found.
        """
        # reset constraint set
        self.model.constraints.clear()

        # add problem constraints
        self.add_problem_constraints()

        # Optimize for error
        error_obj = self.error_objective_expression()
        self.model += error_obj
        sol = self.model.solve()
        if sol is pulp.LpStatusInfeasible:
            raise Exception("Infeasible model.")

        # Optimize for size
        best_fit_value = pulp.value(self.model.objective)
        self.model += error_obj == best_fit_value, 'fit_constraint'
        size_obj = self.size_objective_expression()
        self.model += size_obj
        sol = self.model.solve()
        if sol is pulp.LpStatusInfeasible:
            raise Exception("Infeasible model.")

    def K_dominnace_check_2(self, u_d, v_d, _inequalities):
        """

        :param u_d: a d-dimensional vector(list) like [ 8.53149891  3.36436796]
        :param v_d: tha same list like u_d
        :param _inequalities: list of constraints on d-dimensional Lambda Polytope like
         [[0, 1, 0], [1, -1, 0], [0, 0, 1], [1, 0, -1], [0.0, 1.4770889, -3.1250839]]
        :return: True if u is Kdominance to v regarding given _inequalities otherwise False
        """
        _d = len(u_d)

        prob = LpProblem("Kdominance", LpMinimize)
        lambda_variables = LpVariable.dicts("l", range(_d), 0)

        for inequ in _inequalities:
            prob += lpSum([inequ[j + 1] * lambda_variables[j] for j in range(0, _d)]) + inequ[0] >= 0

        prob += lpSum([lambda_variables[i] * (u_d[i]-v_d[i]) for i in range(_d)])

        #prob.writeLP("show-Ldominance.lp")

        status = prob.solve()
        LpStatus[status]

        result = value(prob.objective)
        if result < 0:
            return False

        return True

    def K_dominance_check(self, _V_best_d, Q_d):
        """
        :param _V_best_d: a list of d-dimension
        :param Q_d: a list of d-dimension
        :return: True if _V_best_d is prefered to Q_d regarding self.Lambda_inequalities and using Kdominance
         other wise it returns False
        """
        _d = len(_V_best_d)

        prob = LpProblem("Ldominance", LpMinimize)
        lambda_variables = LpVariable.dicts("l", range(_d), 0)

        for inequ in self.Lambda_ineqalities:
            prob += lpSum([inequ[j + 1] * lambda_variables[j] for j in range(0, _d)]) + inequ[0] >= 0

        prob += lpSum([lambda_variables[i] * (_V_best_d[i]-Q_d[i]) for i in range(_d)])

        #prob.writeLP("show-Ldominance.lp")

        status = prob.solve()
        LpStatus[status]
        result = value(prob.objective)

        if result < 0:
            return False

        return True

def solve(definition):
    print >> sys.stderr, '-'*20, 'problem'
    print >> sys.stderr, definition

    if debug > 1: print >> sys.stderr, '-'*20, 'parse'
    exs, vars_, target, cat, mm = parse(definition)

    if debug > 1: print >> sys.stderr, '-'*20, 'setup'
    prob, target2, lpvar = setup(exs, vars_, target, cat, mm)

    if debug > 1: print >> sys.stderr, '-'*20, 'solve'
    status = prob.solve(pulp.GLPK(msg=0))

    print >> sys.stderr, '-'*20, 'result'
    print pulp.LpStatus[status]

    if pulp.LpStatus[status] == 'Optimal':
        def sort_var_name(x, y):
            t = cmp(len(x), len(y))
            return t if t != 0 else cmp(x, y)

        res = {}
        for n in sorted(vars_, sort_var_name):
            v = pulp.value(lpvar[n])
            res[n] = v
            print n, ':', v
        print 'target:', eval(target, res), '=', target2

    def pulp_solve(self):
        problem = pulp.LpProblem(pulp.LpMinimize)

        vs = self.get_vars()
        our_vars = dict()
        v = []
        for i in range(len(vs)):
            v.append(pulp.LpVariable('%d' %vs[i][0], vs[i][1], vs[i][2]))
            our_vars[vs[i][0]] = v[i]

        ob = self.get_objective()
        problem.objective = pulp.LpAffineExpression([(our_vars[ob[i][0]], ob[i][1]) for i in range(len(ob))])

        css = self.get_linear_cons()
        for i in range(len(css)):
            ids, coefs, cnst = css[i]
            c = pulp.LpConstraint(
                pulp.LpAffineExpression([(our_vars[ids[j]], coefs[j]) for j in range(len(coefs))], 
                                        constant=cnst)
                , sense=-1)
            problem.constraints[i] = c
            
        problem.solve()
        self.solutions.clear()
        for variable in problem.variables():
            self.solutions[int(variable.name)] = variable.varValue
        self.obj_val = pulp.value(problem.objective)