Python pulp.LpProblem() Examples

def create_model(self):
        def distances(assignment):
            return l2_distance(self.data[assignment[0]], self.centroids[assignment[1]])

        clusters = list(range(self.k))
        assignments = [(i, j)for i in range(self.n) for j in range(self.k)]

        # outflow variables for data nodes
        self.y = pulp.LpVariable.dicts('data-to-cluster assignments',
                                  assignments,
                                  lowBound=0,
                                  upBound=1,
                                  cat=pulp.LpInteger)

        # outflow variables for cluster nodes
        self.b = pulp.LpVariable.dicts('cluster outflows',
                                  clusters,
                                  lowBound=0,
                                  upBound=self.n-self.min_size,
                                  cat=pulp.LpContinuous)

        # create the model
        self.model = pulp.LpProblem("Model for assignment subproblem", pulp.LpMinimize)

        # objective function
        self.model += pulp.lpSum(distances(assignment) * self.y[assignment] for assignment in assignments)

        # flow balance constraints for data nodes
        for i in range(self.n):
            self.model += pulp.lpSum(self.y[(i, j)] for j in range(self.k)) == 1

        # flow balance constraints for cluster nodes
        for j in range(self.k):
            self.model += pulp.lpSum(self.y[(i, j)] for i in range(self.n)) - self.min_size == self.b[j]

        # flow balance constraint for the sink node
        self.model += pulp.lpSum(self.b[j] for j in range(self.k)) == self.n - (self.k * self.min_size) 

def find_constrained_input(self, problem: pulp.LpProblem,
                             metric: PulpLinearMetric,
                             x: np.ndarray,
                             name_prefix = 'x_0_0'):

    d_var = LpVariable(metric.dist_var_name,
                       lowBound = metric.draw_lower_bound (),
                       upBound = metric.upper_bound)
    problem += d_var

    in_vars = self.input_layer_encoder.pulp_in_vars ()
    assert (in_vars.shape == x.shape)
    metric.pulp_constrain (problem, d_var, in_vars, x, name_prefix)

    tp1 ('LP solving: {} constraints'.format(len(problem.constraints)))
    problem.solve (self.solver)
    tp1 ('Solved!')

    if LpStatus[problem.status] != 'Optimal':
      return None

    res = np.zeros(in_vars.shape)
    for idx, var in np.ndenumerate (in_vars):
      res[idx] = pulp.value (var)

    return pulp.value(problem.objective), res


# --- 

def define_lp(fares, product_names):
    """Set up LP.

    Parameters
    ----------
    fares: 2D np array
            contains fares for products, size n_classes*n_products
    products_names: list
            product names (typically relation/class combinations)

    Returns
    -------
    tuple of the form (LP problem, decision variables)
    """
    prob = pulp.LpProblem('network_RM', pulp.LpMaximize)

    # decision variables: available seats per product
    x = pulp.LpVariable.dicts('x', product_names, lowBound=0)

    # objective function
    revenue = pulp.lpSum([x[it]*fares.ravel()[i] for i, it in
                          enumerate(product_names)])
    prob += revenue

    return prob, x 

def update_objective(self, input_data, dmu_code, input_variables,
                         output_variables, lp_model):
        ''' Updates coefficients of the objective function of a given model.

            Args:
                input_data (InputData): object that stores input data.
                dmu_code (str): DMU code.
                input_variables (dict of str to pulp.LpVariable): dictionary
                    that maps variable name to pulp variable corresponding
                    to input categories.
                output_variables (dict of str to pulp.LpVariable): dictionary
                    that maps variable name to pulp variable corresponding
                    to output categories.
                lp_model (pulp.LpProblem): linear programming model.
        '''
        for category in input_data.output_categories:
            lp_model.objective[output_variables[category]] = input_data.coefficients[
                dmu_code, category] 

def update_equality_constraint(self, input_data, dmu_code, input_variables,
                                   output_variables, lp_model):
        ''' Updates coefficients of the equality constraint of a given model.

            Args:
                input_data (InputData): object that stores input data.
                dmu_code (str): DMU code.
                input_variables (dict of str to pulp.LpVariable): dictionary
                    that maps variable name to pulp variable corresponding
                    to input categories.
                output_variables (dict of str to pulp.LpVariable): dictionary
                    that maps variable name to pulp variable corresponding
                    to output categories.
                lp_model (pulp.LpProblem): linear programming model.
        '''
        for category, var in input_variables.items():
            lp_model.constraints['equality_constraint'][var] = input_data.coefficients[
                dmu_code, category] 

def update_objective(self, input_data, dmu_code, input_variables,
                         output_variables, lp_model):
        ''' Updates coefficients of the objective function of a given model.

            Args:
                input_data (InputData): object that stores input data.
                dmu_code (str): DMU code.
                input_variables (dict of str to pulp.LpVariable): dictionary
                    that maps variable name to pulp variable corresponding
                    to input categories.
                output_variables (dict of str to pulp.LpVariable): dictionary
                    that maps variable name to pulp variable corresponding
                    to output categories.
                lp_model (pulp.LpProblem): linear programming model.
        '''
        for category in input_data.input_categories:
            lp_model.objective[input_variables[category]] = input_data.coefficients[
                dmu_code, category] 

def update_equality_constraint(self, input_data, dmu_code, input_variables,
                                   output_variables, lp_model):
        ''' Updates coefficients of the equality constraint of a given model.

            Args:
                input_data (InputData): object that stores input data.
                dmu_code (str): DMU code.
                input_variables (dict of str to pulp.LpVariable): dictionary
                    that maps variable name to pulp variable corresponding
                    to input categories.
                output_variables (dict of str to pulp.LpVariable): dictionary
                    that maps variable name to pulp variable corresponding
                    to output categories.
                lp_model (pulp.LpProblem): linear programming model.
        '''
        for category, var in output_variables.items():
            lp_model.constraints['equality_constraint'][var] = input_data.coefficients[
                dmu_code, category] 

def make_into_lp_problem(good_for, N, add_noise=False):
    """
    Helper function for solve_with_lp_and_reduce()

    N --- number of isoform sequences
    good_for --- dict of <isoform_index> --> list of matched paths index
    """
    prob = LpProblem("The Whiskas Problem",LpMinimize)

    # each good_for is (isoform_index, [list of matched paths index])
    # ex: (0, [1,2,4])
    # ex: (3, [2,5])
    used_paths = []
    for t_i, p_i_s in good_for:
        used_paths += p_i_s
    used_paths = list(set(used_paths))

    variables = [LpVariable(str(i),0,1,LpInteger) for i in used_paths]
    #variables = [LpVariable(str(i),0,1,LpInteger) for i in xrange(N)]

    # objective is to minimize sum_{Xi}
    prob += sum(v for v in variables)

    already_seen = set()
    # constraints are for each isoform, expressed as c_i * x_i >= 1
    # where c_i = 1 if x_i is matched for the isoform
    # ex: (0, [1,2,4]) becomes t_0 = x_1 + x_2 + x_4 >= 1
    for t_i, p_i_s in good_for:
        #c_i_s = [1 if i in p_i_s else 0 for i in xrange(N)]
        #prob += sum(variables[i]*(1 if i in p_i_s else 0) for i in xrange(N)) >= 1
        p_i_s.sort()
        pattern = ",".join(map(str,p_i_s))
        #print >> sys.stderr, t_i, p_i_s, pattern
        if pattern not in already_seen:
            if add_noise:
                prob += sum(variables[i]*(1+random.random() if p in p_i_s else 0) for i,p in enumerate(used_paths)) >= 1
            else:
                prob += sum(variables[i]*(1 if p in p_i_s else 0) for i,p in enumerate(used_paths)) >= 1
        already_seen.add(pattern)
    prob.writeLP('cogent.lp')
    return prob 

def label_prop(C, nt, Dct, lp="linear"):
    
#Inputs:
#  C      :    Number of share classes between src and tar
#  nt     :    Number of target domain samples
#  Dct    :    All d_ct in matrix form, nt * C
#  lp     :    Type of linear programming: linear (default) | binary
#Outputs:
#  Mcj    :    all M_ct in matrix form, m * C
    
    Dct = abs(Dct)
    model = pulp.LpProblem("Cost minimising problem", pulp.LpMinimize)
    Mcj = pulp.LpVariable.dicts("Probability",
                                ((i, j) for i in range(C) for j in range(nt)),
                                lowBound=0,
                                upBound=1,
                                cat='Continuous')
    
    # Objective Function
    model += (
    pulp.lpSum([Dct[j, i]*Mcj[(i, j)] for i in range(C) for j in range(nt)])
    )
    
    # Constraints
    for j in range(nt):
        model += pulp.lpSum([Mcj[(i, j)] for i in range(C)]) == 1
    for i in range(C):
        model += pulp.lpSum([Mcj[(i, j)] for j in range(nt)]) >= 1
    
    # Solve our problem
    model.solve()
    pulp.LpStatus[model.status]
    Output = [[Mcj[i, j].varValue for i in range(C)] for j in range(nt)]
    
    return np.array(Output) 

def __min_one_norm(B, initial_seed, seed):
    weight_initial = 1 / float(len(initial_seed))
    weight_later_added = weight_initial / float(0.5)
    difference = len(seed) - len(initial_seed)
    [r, c] = B.shape
    prob = pulp.LpProblem("Minimum one norm", pulp.LpMinimize)
    indices_y = range(0, r)
    y = pulp.LpVariable.dicts("y_s", indices_y, 0)
    indices_x = range(0, c)
    x = pulp.LpVariable.dicts("x_s", indices_x)

    f = dict(zip(indices_y, [1.0] * r))

    prob += pulp.lpSum(f[i] * y[i] for i in indices_y)  # objective function

    prob += pulp.lpSum(y[s] for s in initial_seed) >= 1

    prob += pulp.lpSum(y[r] for r in seed) >= 1 + weight_later_added * difference

    for j in range(r):
        temp = dict(zip(indices_x, list(B[j, :])))
        prob += pulp.lpSum(y[j] + (temp[k] * x[k] for k in indices_x)) == 0

    prob.solve()

    result = []
    for var in indices_y:
        result.append(y[var].value())

    return result 

def solve_ilp(self, N):
        # build the A matrix: a_ij is 1 if j-th gram appears in the i-th sentence

        A = np.zeros((len(self.sentences_idx), len(self.ref_ngrams_idx)))
        for i in self.sentences_idx:
            sent = self.sentences[i].untokenized_form
            sngrams = list(extract_ngrams2([sent], self.stemmer, self.LANGUAGE, N))
            for j in self.ref_ngrams_idx:
                if self.ref_ngrams[j] in sngrams:
                    A[i][j] = 1

        # Define ILP variable, x_i is 1 if sentence i is selected, z_j is 1 if gram j appears in the created summary
        x = pulp.LpVariable.dicts('sentences', self.sentences_idx, lowBound=0, upBound=1, cat=pulp.LpInteger)
        z = pulp.LpVariable.dicts('grams', self.ref_ngrams_idx, lowBound=0, upBound=1, cat=pulp.LpInteger)

        # Define ILP problem, maximum coverage of grams from the reference summaries
        prob = pulp.LpProblem("ExtractiveUpperBound", pulp.LpMaximize)
        prob += pulp.lpSum(z[j] for j in self.ref_ngrams_idx)

        # Define ILP constraints, length constraint and consistency constraint (impose that z_j is 1 if j
        # appears in the created summary)
        prob += pulp.lpSum(x[i] * self.sentences[i].length for i in self.sentences_idx) <= self.sum_length

        for j in self.ref_ngrams_idx:
            prob += pulp.lpSum(A[i][j] * x[i] for i in self.sentences_idx) >= z[j]

        # Solve ILP problem and post-processing to get the summary
        try:
            print('Solving using CPLEX')
            prob.solve(pulp.CPLEX(msg=0))
        except:
            print('Fall back to GLPK')
            prob.solve(pulp.GLPK(msg=0))
                

        summary_idx = []
        for idx in self.sentences_idx:
            if x[idx].value() == 1.0:
                summary_idx.append(idx)

        return summary_idx 

def for_layer(self, cl) -> Tuple[pulp.LpProblem, PulpVarMap]:
    index = cl.layer_index + (0 if activation_is_relu (cl.layer) else 1)
    return self.base_constraints[index].copy () 

def _solve_balancing_ilp_pulp(A):
    import pulp
    x = [pulp.LpVariable('x%d' % i, lowBound=1, cat='Integer') for i in range(A.shape[1])]
    prob = pulp.LpProblem("chempy balancing problem", pulp.LpMinimize)
    prob += reduce(add, x)
    for expr in [pulp.lpSum([x[i]*e for i, e in enumerate(row)]) for row in A.tolist()]:
        prob += expr == 0
    prob.solve()
    return [pulp.value(_) for _ in x] 

def add_capacity_constraints(prob, x, A, product_names, capacities,
                             leg_names=None):
    """Add capacity contraints as upper bound on the segments/legs.

  Parameters
  ----------
    prob: pulp.LpProblem
        the LP problem
    x: dict
        contains decision variables (each of which has type pulp.pulp.LpVariable)
    A: 2D np array
        incidence matrix, size n_relations*n_legs
    products_names: list
        list of product names (e.g. relation/class combinations)
    capacities: Iterable
        list of capacity on each leg/segment
    leg_names: list
        list of leg names

    Returns
    -------
    list of tuples of the form (constraint, constraint_name)
    where `constraint` is of type pulp.LpConstraint
    and `constraint_name` is of type str
    """
    n_trips, n_legs = A.shape
    n_classes = int(len(product_names) / n_trips)
    A_ = np.tile(A, (n_classes, 1))

    capacity_constraints = []
    for leg, leg_name in enumerate(leg_names):
        leg_load = pulp.lpSum([x[it]*A_[i, leg] for i, it
                               in enumerate(product_names)])
        capacity_constraint = (leg_load <= capacities[leg],
                               "cap_{}".format(leg_name))
        prob += capacity_constraint
        capacity_constraints.append(capacity_constraint)

    return capacity_constraints 

def solve(g):
    el = g.get_edge_list()
    nl = g.get_node_list()
    p = LpProblem('min_cost', LpMinimize)
    capacity = {}
    cost = {}
    demand = {}
    x = {}
    for e in el:
        capacity[e] = g.get_edge_attr(e[0], e[1], 'capacity')
        cost[e] = g.get_edge_attr(e[0], e[1], 'cost')
    for i in nl:
        demand[i] = g.get_node_attr(i, 'demand')
    for e in el:
        x[e] = LpVariable("x"+str(e), 0, capacity[e])
    # add obj
    objective = lpSum (cost[e]*x[e] for e in el)
    p += objective
    # add constraints
    for i in nl:
        out_neig = g.get_out_neighbors(i)
        in_neig = g.get_in_neighbors(i)
        p += lpSum(x[(i,j)] for j in out_neig) -\
             lpSum(x[(j,i)] for j in in_neig)==demand[i]
    p.solve()
    return x, value(objective) 

def _create_lp(self):
        ''' Creates initial linear program.
        '''
        assert len(self.input_data.DMU_codes) != 0
        self._variables.clear()
        self._constraints.clear()
        # create LP for the first DMU - it is the easiest way to
        # adapt current code
        for elem in self.input_data.DMU_codes:
            dmu_code = elem
            break

        orientation = self._concrete_model.get_orientation()
        obj_type = self._concrete_model.get_objective_type()
        self.lp_model = pulp.LpProblem('Envelopment model: {0}-'
                                       'oriented'.format(orientation),
                                       obj_type)
        ub_obj = self._concrete_model.get_upper_bound_for_objective_variable()
        obj_variable = pulp.LpVariable('Variable in objective function',
                                       self._concrete_model.
                                       get_lower_bound_for_objective_variable(),
                                       ub_obj,
                                       pulp.LpContinuous)
        self._variables['obj_var'] = obj_variable

        self.lp_model += (obj_variable,
                          'Efficiency score or inverse of efficiency score')

        lambda_variables = pulp.LpVariable.dicts('lambda',
                                                 self.input_data.DMU_codes,
                                                 0, None, pulp.LpContinuous)
        self._variables.update(lambda_variables)
        self._add_constraints_for_outputs(lambda_variables, dmu_code,
                                          obj_variable)
        self._add_constraints_for_inputs(lambda_variables, dmu_code,
                                         obj_variable) 

def get_minima(self, phases, bounds):
        """
        Given a set of Phases, get_minima will determine the minimum
        free energy elemental composition as a weighted sum of these
        compounds
        """

        prob = pulp.LpProblem('GibbsEnergyMin', pulp.LpMinimize)
        pvars = pulp.LpVariable.dicts('phase', phases, 0)
        bvars = pulp.LpVariable.dicts('bound', bounds, 0.0, 1.0)
        prob += pulp.lpSum( self.phase_energy(p)*pvars[p] for p in phases ) - \
                pulp.lpSum( self.phase_energy(bound)*bvars[bound] for bound in bounds ), \
                                "Free Energy"
        for elt in self.bound_space:
            prob += sum([ p.unit_comp.get(elt,0)*pvars[p] for p in phases ])\
                        == \
                sum([ b.unit_comp.get(elt, 0)*bvars[b] for b in bounds ]),\
                            'Contraint to the proper range of'+elt
        prob += sum([ bvars[b] for b in bounds ]) == 1, \
                'sum of bounds must be 1'

        if pulp.GUROBI().available():
            prob.solve(pulp.GUROBI(msg=False))
        elif pulp.COIN_CMD().available():
            prob.solve(pulp.COIN_CMD())
        elif pulp.COINMP_DLL().available():
            prob.solve(pulp.COINMP_DLL())
        else:
            prob.solve()

        E = pulp.value(prob.objective)
        xsoln = defaultdict(float,
            [(p, pvars[p].varValue) for p in phases if
                abs(pvars[p].varValue) > 1e-4])
        return xsoln, E 

def _gclp(self, composition={}, mus={}, phases=[]):
        if not qmpy.FOUND_PULP:
            raise Exception('Cannot do GCLP without installing PuLP and an LP',
                    'solver')
        prob = pulp.LpProblem('GibbsEnergyMin', pulp.LpMinimize)
        phase_vars = pulp.LpVariable.dicts('lib', phases, 0.0)
        prob += pulp.lpSum([ (p.energy -
            sum([ p.unit_comp.get(elt,0)*mu
                for elt, mu in mus.items() ])) * phase_vars[p]
            for p in phases]),\
                    "Free Energy"
        for elt, constraint in composition.items():
            prob += pulp.lpSum([
                p.unit_comp.get(elt,0)*phase_vars[p]
                for p in phases ]) == float(constraint),\
                        'Conservation of '+elt
        ##[vh]
        ##print prob
        if pulp.GUROBI().available():
            prob.solve(pulp.GUROBI(msg=False))
        elif pulp.COIN_CMD().available():
            prob.solve(pulp.COIN_CMD())
        else:
            prob.solve()

        phase_comp = dict([ (p, phase_vars[p].varValue)
            for p in phases if phase_vars[p].varValue > 1e-5])
        
        energy = sum( p.energy*amt for p, amt in phase_comp.items() )
        energy -= sum([ a*composition.get(e, 0) for e,a in mus.items()])
        return energy, phase_comp 

def create_model(self):
        def distances(assignment):
            return l2_distance(self.data[assignment[0]], self.centroids[assignment[1]])

        assignments = [(i, j) for i in range(self.n) for j in range(self.k)]

        # assignment variables
        self.y = pulp.LpVariable.dicts('data-to-cluster assignments',
                                  assignments,
                                  lowBound=0,
                                  upBound=1,
                                  cat=pulp.LpInteger)

        # create the model
        self.model = pulp.LpProblem("Model for assignment subproblem", pulp.LpMinimize)

        # objective function
        self.model += pulp.lpSum([distances(assignment) * self.weights[assignment[0]] * self.y[assignment] for assignment in assignments]), 'Objective Function - sum weighted squared distances to assigned centroid'
        # this is also weighted, otherwise the weighted centroid computation don't make sense.

        # constraints on the total weights of clusters
        for j in range(self.k):
            self.model += pulp.lpSum([self.weights[i] * self.y[(i, j)] for i in range(self.n)]) >= self.min_weight, "minimum weight for cluster {}".format(j)
            self.model += pulp.lpSum([self.weights[i] * self.y[(i, j)] for i in range(self.n)]) <= self.max_weight, "maximum weight for cluster {}".format(j)

        # make sure each point is assigned at least once, and only once
        for i in range(self.n):
            self.model += pulp.lpSum([self.y[(i, j)] for j in range(self.k)]) == 1, "must assign point {}".format(i) 

def create_model(self):
        def distances(assignment):
            return l2_distance(self.data[assignment[0]], self.centroids[assignment[1]])

        clusters = list(range(self.k))
        assignments = [(i, j)for i in range(self.n) for j in range(self.k)]

        # outflow variables for data nodes
        self.y = pulp.LpVariable.dicts('data-to-cluster assignments',
                                  assignments,
                                  lowBound=0,
                                  upBound=1,
                                  cat=pulp.LpInteger)

        # outflow variables for cluster nodes
        self.b = pulp.LpVariable.dicts('cluster outflows',
                                  clusters,
                                  lowBound=0,
                                  upBound=self.n-self.min_size,
                                  cat=pulp.LpContinuous)

        # create the model
        self.model = pulp.LpProblem("Model for assignment subproblem", pulp.LpMinimize)

        # objective function
        self.model += pulp.lpSum([distances(assignment) * self.y[assignment] for assignment in assignments])

        # flow balance constraints for data nodes
        for i in range(self.n):
            self.model += pulp.lpSum(self.y[(i, j)] for j in range(self.k)) == 1

        # flow balance constraints for cluster nodes
        for j in range(self.k):
            self.model += pulp.lpSum(self.y[(i, j)] for i in range(self.n)) - self.min_size == self.b[j]
            
        # capacity constraint on outflow of cluster nodes
        for j in range(self.k):
            self.model += self.b[j] <= self.max_size - self.min_size 

        # flow balance constraint for the sink node
        self.model += pulp.lpSum(self.b[j] for j in range(self.k)) == self.n - (self.k * self.min_size) 

def __init__(self, input_data, input_params):
        self.input_data = input_data
        self.input_params = input_params
        self.model = pulp.LpProblem(name='prod_planning', sense=pulp.LpMinimize)
        self._create_decision_variables()
        self._create_main_constraints()
        self._set_objective_function()

    # ================== Decision variables ================== 

def __init__(self):
        self.prob = LpProblem('Daily Fantasy Sports', LpMaximize) 

def get_reaction(self, var, facet=None):
        """
        For a given composition, what is the maximum delta_composition reaction
        on the given facet. If None, returns the whole reaction for the given
        PhaseSpace.

        Examples::

            >>> space = PhaseSpace('Fe2O3-Li2O')
            >>> equilibria = space.hull[0]
            >>> space.get_reaction('Li2O', facet=equilibria)

        """

        if isinstance(var, basestring):
            var = parse_comp(var)

        if facet:
            phases = facet
        else:
            phases = self.stable

        prob = pulp.LpProblem('BalanceReaction', pulp.LpMaximize)
        pvars = pulp.LpVariable.dicts('prod', phases, 0)
        rvars = pulp.LpVariable.dicts('react', phases, 0)
        prob += sum([ p.fraction(var)['var']*pvars[p] for p in phases ])-\
                sum([ p.fraction(var)['var']*rvars[p] for p in phases ]),\
                "Maximize delta comp"
        for celt in self.space:
            prob += sum([ p.fraction(var)[celt]*pvars[p] for p in phases ]) ==\
                    sum([ p.fraction(var)[celt]*rvars[p] for p in phases ]),\
                    'identical %s composition on both sides' % celt
        prob += sum([ rvars[p] for p in phases ]) == 1
        
        if pulp.GUROBI().available():
            prob.solve(pulp.GUROBI(msg=False))
        elif pulp.COIN_CMD().available():
            prob.solve(pulp.COIN_CMD())
        elif pulp.COINMP_DLL().available():
            prob.solve(pulp.COINMP_DLL())
        else:
            prob.solve()

        prods = defaultdict(float,[ (c, pvars[c].varValue) for c in phases
            if pvars[c].varValue > 1e-4 ])
        reacts = defaultdict(float,[ (c, rvars[c].varValue) for c in phases
            if rvars[c].varValue > 1e-4 ])
        n_elt = pulp.value(prob.objective)
        return reacts, prods, n_elt 

def solution(events, slots, objective_function=None, solver=None, **kwargs):
    """Compute a schedule in solution form

    Parameters
    ----------
    events : list or tuple
        of :py:class:`resources.Event` instances
    slots : list or tuple
        of :py:class:`resources.Slot` instances
    solver : pulp.solver
        a pulp solver
    objective_function: callable
        from lp_problem.objective_functions
    kwargs : keyword arguments
        arguments for the objective function

    Returns
    -------
    list
        A list of tuples giving the event and slot index (for the given
        events and slots lists) for all scheduled items.

    Example
    -------
    For a solution where

        * event 0 is scheduled in slot 1
        * event 1 is scheduled in slot 4
        * event 2 is scheduled in slot 5

    the resulting list would be::

        [(0, 1), (1, 4), (2, 5)]
    """
    shape = Shape(len(events), len(slots))
    problem = pulp.LpProblem()
    X = lp.utils.variables(shape)
    beta = pulp.LpVariable("upper_bound")

    for constraint in lp.constraints.all_constraints(
        events, slots, X, beta, 'lpsum'
    ):
        problem += constraint.condition

    if objective_function is not None:
        problem += objective_function(events=events, slots=slots, X=X,
                                      beta=beta,
                                      **kwargs)

    status = problem.solve(solver=solver)
    if status == 1:
        return [item for item, variable in X.items() if variable.value() > 0]
    else:
        raise ValueError('No valid solution found') 

def word_mover_distance_probspec(first_sent_tokens, second_sent_tokens, wvmodel, distancefunc=euclidean, lpFile=None):
    """ Compute the Word Mover's distance (WMD) between the two given lists of tokens, and return the LP problem class.

    Using methods of linear programming, supported by PuLP, calculate the WMD between two lists of words. A word-embedding
    model has to be provided. The problem class is returned, containing all the information about the LP.

    Reference: Matt J. Kusner, Yu Sun, Nicholas I. Kolkin, Kilian Q. Weinberger, "From Word Embeddings to Document Distances," *ICML* (2015).

    :param first_sent_tokens: first list of tokens.
    :param second_sent_tokens: second list of tokens.
    :param wvmodel: word-embedding models.
    :param distancefunc: distance function that takes two numpy ndarray.
    :param lpFile: log file to write out.
    :return: a linear programming problem contains the solution
    :type first_sent_tokens: list
    :type second_sent_tokens: list
    :type wvmodel: gensim.models.keyedvectors.KeyedVectors
    :type distancefunc: function
    :type lpFile: str
    :rtype: pulp.LpProblem
    """
    all_tokens = list(set(first_sent_tokens+second_sent_tokens))
    wordvecs = {token: wvmodel[token] for token in all_tokens}

    first_sent_buckets = tokens_to_fracdict(first_sent_tokens)
    second_sent_buckets = tokens_to_fracdict(second_sent_tokens)

    T = pulp.LpVariable.dicts('T_matrix', list(product(all_tokens, all_tokens)), lowBound=0)

    prob = pulp.LpProblem('WMD', sense=pulp.LpMinimize)
    prob += pulp.lpSum([T[token1, token2]*distancefunc(wordvecs[token1], wordvecs[token2])
                        for token1, token2 in product(all_tokens, all_tokens)])
    for token2 in second_sent_buckets:
        prob += pulp.lpSum([T[token1, token2] for token1 in first_sent_buckets])==second_sent_buckets[token2]
    for token1 in first_sent_buckets:
        prob += pulp.lpSum([T[token1, token2] for token2 in second_sent_buckets])==first_sent_buckets[token1]

    if lpFile!=None:
        prob.writeLP(lpFile)

    prob.solve()

    return prob 

def create_lscp_model(coverage_dict, model_file=None, delineator="$", ):
    """
    Creates a LSCP (Location set covering problem) using the provided coverage and
    parameters. Writes a .lp file which can be solved with Gurobi

    Church, R., & Murray, A. (2009). Coverage Business Site Selection, Location
    Analysis, and GIS (pp. 209-233). Hoboken, New Jersey: Wiley.

    :param coverage_dict: (dictionary) The coverage to use to generate the model
    :param model_file: (string) The model file to output
    :param delineator: (string) The character(s) to use to delineate the layer from the ids
    :return: (Pulp problem) The generated problem to solve
    """
    validate_coverage(coverage_dict, ["coverage"], ["binary"])
    if not isinstance(coverage_dict, dict):
        raise TypeError("coverage_dict is not a dictionary")
    if model_file and not (isinstance(model_file, str)):
        raise TypeError("model_file is not a string")
    if not isinstance(delineator, str):
        raise TypeError("delineator is not a string")
        # create the variables
    demand_vars = {}
    for demand_id in coverage_dict["demand"]:
        demand_vars[demand_id] = pulp.LpVariable("Y{}{}".format(delineator, demand_id), 0, 1, pulp.LpInteger)
    facility_vars = {}
    for facility_type in coverage_dict["facilities"]:
        facility_vars[facility_type] = {}
        for facility_id in coverage_dict["facilities"][facility_type]:
            facility_vars[facility_type][facility_id] = pulp.LpVariable(
                "{}{}{}".format(facility_type, delineator, facility_id), 0, 1, pulp.LpInteger)
    # create the problem
    prob = pulp.LpProblem("LSCP", pulp.LpMinimize)
    # Create objective, minimize number of facilities
    to_sum = []
    for facility_type in coverage_dict["facilities"]:
        for facility_id in coverage_dict["facilities"][facility_type]:
            to_sum.append(facility_vars[facility_type][facility_id])
    prob += pulp.lpSum(to_sum)
    # add coverage constraints
    for demand_id in coverage_dict["demand"]:
        to_sum = []
        for facility_type in coverage_dict["demand"][demand_id]["coverage"]:
            for facility_id in coverage_dict["demand"][demand_id]["coverage"][facility_type]:
                to_sum.append(facility_vars[facility_type][facility_id])
        # Hack to get model to "solve" when infeasible with GLPK.
        # Pulp will automatically add dummy variables when the sum is empty, since these are all the same name,
        # it seems that GLPK doesn't read the lp problem properly and fails
        if not to_sum:
            to_sum = [pulp.LpVariable("__dummy{}{}".format(delineator, demand_id), 0, 0, pulp.LpInteger)]
        prob += pulp.lpSum(to_sum) >= 1, "D{}".format(demand_id)
    if model_file:
        prob.writeLP(model_file)
    return prob