Python gurobipy.Model类代码示例

 def __init__(self, name):
     Model.__init__(self, name)
     print str(self)
     self.testVar = 4
     self.createVars()
     self.update()
     self.createObjectiveFunction()
     self.createConstraints()

    def def_PL(self):
        """
			Definie le PL, variables/contraintes/objectifs
		"""
        model = Model("MR")
        self.set_vars(model)
        self.set_constraints(model)
        self.set_objectif(model)
        model.setParam("OutputFlag", False)
        return model

def create_problem(cobra_model, objective_sense="maximize"):
    lp = Model("cobra")
    lp.Params.OutputFlag = 0

    if objective_sense == "maximize":
        objective_sign = -1.0
    elif objective_sense == "minimize":
        objective_sign = 1.0
    else:
        raise ValueError("objective_sense must be 'maximize' or 'minimize'")

    # create metabolites/constraints
    metabolite_constraints = {}
    for metabolite in cobra_model.metabolites:
        metabolite_constraints[metabolite] = lp.addConstr(
            0.0, sense_dict[metabolite._constraint_sense], metabolite._bound, metabolite.id
        )
    lp.update()

    # create reactions/variables along with S matrix
    for j, reaction in enumerate(cobra_model.reactions):
        constraints = [metabolite_constraints[i] for i in reaction._metabolites]
        stoichiometry = reaction._metabolites.values()
        lp.addVar(
            lb=float(reaction.lower_bound),
            ub=float(reaction.upper_bound),
            obj=objective_sign * reaction.objective_coefficient,
            name=reaction.id,
            vtype=variable_kind_dict[reaction.variable_kind],
            column=Column(stoichiometry, constraints),
        )
    lp.update()
    return lp

 def __loadModel(self):
             
     # Create optimization model
     self.__model = Model('PathBackup')
 
     # Auxiliary variables for SOCP reformulation
     U = {}
     R = {}
     
     # Create variables
     for i,j in self.__links:
         self.__BackupCapacity[i,j] = self.__model.addVar(lb=0, obj=1, name='Backup_Capacity[%s,%s]' % (i, j))
     self.__model.update()
      
     for p in self.__paths:
         #LP Relaxation
         #self.__bPath[self.__paths.index(p)] = self.__model.addVar(lb=0,obj=1,name='Backup_Path[%s]' % (self.__paths.index(p)))
         self.__bPath[self.__paths.index(p)] = self.__model.addVar(vtype=GRB.BINARY,obj=1,name='Backup_Path[%s]' % (self.__paths.index(p)))
     self.__model.update()
     
     for i,j in self.__links:
         U[i,j] = self.__model.addVar(obj=1,name='U[%s,%s]' % (i, j))
     self.__model.update()
     
     for i,j in self.__links:
         for s,d in self.__links:
             R[i,j,s,d] = self.__model.addVar(obj=1,name='R[%s,%s,%s,%s]' % (i,j,s,d))
     self.__model.update()
     
     self.__model.modelSense = GRB.MINIMIZE
     self.__model.setObjective(quicksum(self.__BackupCapacity[i,j] for i,j in self.__links))
     self.__model.update()
     
        
     #------------------------------------------------------------------------#
     #                    Constraints definition                              #
     #                                                                        #
     #                                                                        #
     #------------------------------------------------------------------------#
      
     # Link capacity constraints
     for i,j in self.__links:
         self.__model.addConstr(self.__BackupCapacity[i,j] >= quicksum(self.__mean[s,d]*quicksum(self.__bPath[self.__paths.index(p)] for p in self.__Pij[i,j,s,d]) for s,d in self.__links) + U[i,j]*self.__invstd,'[CONST]Link_Cap[%s][%s]' % (i, j))
     self.__model.update()
         
     # SCOP Reformulation Constraints
     for i,j in self.__links:
         self.__model.addConstr(quicksum(R[i,j,s,d]*R[i,j,s,d] for s,d in self.__links) <= U[i,j]*U[i,j],'[CONST]SCOP1[%s][%s]' % (i, j))
     self.__model.update()
         
     # SCOP Reformulation Constraints    
     for i,j in self.__links:
         for s,d in self.__links:
             self.__model.addConstr(self.__std[s,d]*quicksum(self.__bPath[self.__paths.index(p)] for p in self.__Pij[i,j,s,d]) == R[i,j,s,d],'[CONST]SCOP2[%s][%s][%s][%s]' % (i,j,s,d))
     self.__model.update()
     
     # Unique path 
     for s,d in self.__links:
         self.__model.addConstr(quicksum(self.__bPath[self.__paths.index(p)] for p in self.__Psd[s,d]) == 1,'UniquePath[%s,%s]' % (s, d))
     self.__model.update()

def create_problem(cobra_model, quadratic_component=None, **kwargs):
    """Solver-specific method for constructing a solver problem from
    a cobra.Model.  This can be tuned for performance using kwargs


    """
    lp = Model("")

    the_parameters = parameter_defaults
    if kwargs:
        the_parameters = parameter_defaults.copy()
        the_parameters.update(kwargs)

    # Set verbosity first to quiet infos on parameter changes
    if "verbose" in the_parameters:
        set_parameter(lp, "verbose", the_parameters["verbose"])
    for k, v in iteritems(the_parameters):
        set_parameter(lp, k, v)


    # Create variables
    #TODO:  Speed this up
    variable_list = [lp.addVar(_float(x.lower_bound),
                               _float(x.upper_bound),
                               float(x.objective_coefficient),
                               variable_kind_dict[x.variable_kind],
                               str(i))
                     for i, x in enumerate(cobra_model.reactions)]
    reaction_to_variable = dict(zip(cobra_model.reactions,
                                    variable_list))
    # Integrate new variables
    lp.update()

    #Constraints are based on mass balance
    #Construct the lin expression lists and then add
    #TODO: Speed this up as it takes about .18 seconds
    #HERE
    for i, the_metabolite in enumerate(cobra_model.metabolites):
        constraint_coefficients = []
        constraint_variables = []
        for the_reaction in the_metabolite._reaction:
            constraint_coefficients.append(_float(the_reaction._metabolites[the_metabolite]))
            constraint_variables.append(reaction_to_variable[the_reaction])
        #Add the metabolite to the problem
        lp.addConstr(LinExpr(constraint_coefficients, constraint_variables),
                     sense_dict[the_metabolite._constraint_sense.upper()],
                     the_metabolite._bound,
                     str(i))

    # Set objective to quadratic program
    if quadratic_component is not None:
        set_quadratic_objective(lp, quadratic_component)

    lp.update()
    return(lp)

    def __init_prob(self):  # Modify parameters
        logging.info("Generating initial marking problem")
        self.prob = Model("SC1")

        # Gurobi parameters:
        if not self.verbose:
            self.prob.params.OutputFlag = 0
            try:
                os.remove("gurobi.log")
            except OSError:
                pass
        self.prob.params.Threads = 4

    def _init_LP(self):
        if self._lp_inited:
            return

        logging.debug('Init LP')
        self.lp = LPModel('estep')
        self.lp.setAttr("modelSense", 1)    # minimzation

        self.alpha = {}
        beta2 = {}
        beta3 = {}
        # instantiate vars
        logging.debug('Init LP - create vars')
        for a in self.author_graph:
            self.alpha[a] = {}
            for p in self.parts:
                self.alpha[a][p] = self.lp.addVar(lb=0.0)
        for a, b in self.author_graph.edges():
            beta2[(a, b)] = self.lp.addVar()
            beta3[(a, b)] = {}
            for p in self.parts:
                beta3[(a, b)][p] = self.lp.addVar(lb=0.0)
        # integrate added variables into the model
        self.lp.update()
        # add constraints once during this init
        # alphas are indicator vars
        logging.debug('Init LP - indiv constraints')
        ones_arr = [1.0] * len(self.parts)
        for a in self.author_graph:
            self.lp.addConstr(LinExpr(ones_arr, self.alpha[a].values()), GRB.EQUAL, 1.0)
        # beta2 is the sum of beta3s
        logging.debug('Init LP - pair constraints')
        pt_five_array = [0.5] * len(self.parts)
        for a, b in self.author_graph.edges():
            self.lp.addConstr(LinExpr(pt_five_array, beta3[(a, b)].values()), GRB.EQUAL, beta2[(a, b)])
            for p in self.parts:
                self.lp.addConstr(LinExpr([1.0, -1.0], [self.alpha[a][p], self.alpha[b][p]]), GRB.LESS_EQUAL, beta3[(a, b)][p])
                self.lp.addConstr(LinExpr([-1.0, 1.0], [self.alpha[a][p], self.alpha[b][p]]), GRB.LESS_EQUAL, beta3[(a, b)][p])
        self.lp.update()
        # calculate pairwise potentials part of the objective
        # the optimization is to minimize negated log-likelihood = maximize the log-likelihood
        logging.debug('Obj func - pair potentials')
        s = np.log(1 - self.TAU) - np.log(self.TAU)
        lpcoeffs, lpvars = [], []
        for a, b in self.author_graph.edges():
            lpcoeffs.append(-self.author_graph[a][b]['weight'] * s)
            lpvars.append(beta2[(a, b)])
        self.objF_pair = LinExpr(list(lpcoeffs), list(lpvars))

        self._lp_inited = True
        logging.debug('Init LP Done')

 def __loadModel(self,imp_samp, backup_link,link_capacity):
              
     # Create optimization model
     self.__model = Model('Backup')
  
     for i,j in self.__links:
         for k in range(self.__N):
             self.__z[k,i,j] = self.__model.addVar(lb=0,name='z[%s][%s][%s]' % (k,i,j))
     self.__model.update()
     
     for i,j in self.__links: 
         self.__z0[i,j] = self.__model.addVar(lb=-GRB.INFINITY,name='z0[%s][%s]' %(i,j))
     self.__model.update()
      
     for i,j in self.__links: 
         self.__q[i,j] = self.__model.addVar(lb=-GRB.INFINITY,name='q[%s][%s]' %(i,j))
     self.__model.update()
     
     self.__model.modelSense = GRB.MINIMIZE
     
     self.__model.setObjective(quicksum(self.__q[i,j] for i,j in self.__links))
     
     self.__model.update()
      
         
     #------------------------------------------------------------------------#
     #                    Constraints definition                              #
     #                                                                        #
     #                                                                        #
     #------------------------------------------------------------------------#
       
     # Buffer probability I
     for i,j in self.__links:
         self.__model.addConstr(self.__z0[i,j] + 1/(self.__N*self.__epsilon)*quicksum(self.__z[k,i,j]*imp_samp[k] for (k) in range(self.__N)) <= self.__q[i,j],'[CONST]Buffer_Prob_I[%s][%s]'%(i,j))
     self.__model.update()
      
     # Link capacity constraints
     for i,j in self.__links:
         for k in range(self.__N):
             self.__model.addConstr((quicksum(backup_link[i,j,s,d]*self.__capacity[k,s,d] for s,d in self.__links) - link_capacity[i,j] - self.__z0[i,j]) <= self.__z[k,i,j],'[CONST]Buffer_Prob_II[%s][%s][%s]' % (k,i,j))
     self.__model.update()
     
     # Link capacity constraints
     for i,j in self.__links:
         for k in range(self.__N):
             self.__model.addConstr(self.__z[k,i,j] >= 0,'[CONST]Buffer_Prob_III[%s][%s][%s]' % (k,i,j))
     self.__model.update()

    def __init__(self, length):

        """Initialize model, note and chord variables, and constraints for notes
        and chord constraints."""

        self.model = Model("harmony")
        self.length = length

        config = load_source("config",
                join(abspath(dirname(__file__)), "config.py"))

        config.length = self.length

        config.model = self.model

        self.notes = Notes(config)
        config.notes = self.notes

        self.chords = Chords(config)

    def __init__(self, json_data):
        self.m = Model('mip1.log')
        self.solution_matrix = None
        
        # Load from JSON
        self.all_data = json_data
        self.student_ids = [int(k) for k in json_data['students'].keys()]
        self.course_ids = [int(k) for k in json_data['courses'].keys()]
        self.semester_ids = [int(k) for k in json_data['semesters'].keys()]
        self.dependencies = [(d['first_course'], d['second_course']) for d in json_data['course_dependencies']]

        self.student_demand = defaultdict(lambda: False)
        for sd in json_data['student_demand']:
            self.student_demand[sd['student_id'], sd['course_id'], sd['semester_id']] = True

        self.instructor_availability = defaultdict(lambda: False)
        for ia in json_data['instructor_pool']:
            self.instructor_availability[ia['instructor_id'], ia['course_id'], ia['semester_id']] = True

        self.instructor_ids = [s['instructor_id'] for s in json_data['instructor_pool']]

def build_model(data, n_cliques = 0, verbose = True):
    
    # Load Data Format
    n = data['n']
    r = data['r']
    p = data['p']
    s = data['s']
    c = data['c']
    h = data['h']
    w = data['w']
    location = data['location']
    conflicts = data['conflicts']
    locking_times = data['locking_times']
    T = data['T']
    
    model = Model("ExaminationScheduling")
    
    
    if verbose:
        print("Building variables...")
    
    # x[i,k,l] = 1 if exam i is at time l in room k
    x = {}
    for k in range(r):
        for l in range(p):
            if T[k][l] == 1:
                for i in range(n):
                    if location[k] in w[i]:
                        x[i,k,l] = model.addVar(vtype=GRB.BINARY, name="x_%s_%s_%s" % (i,k,l))
    
    # y[i,l] = 1 if exam i is at time l
    y = {}
    for i in range(n):
        for l in range(p):
            y[i, l] = model.addVar(vtype=GRB.BINARY, name="y_%s_%s" % (i,l))
    

    # integrate new variables
    model.update() 

    start = timeit.default_timer()

    # not very readable but same constraints as in GurbiLinear_v_10: speeded up model building by 2 for small problems (~400 exams) and more for huger problem ~1500 exams
    if verbose:
        print("Building constraints...")    
    
    obj = LinExpr()
    sumconflicts = {}
    maxrooms = {}
    for i in range(n):
        sumconflicts[i] = sum(conflicts[i])
        if s[i] <= 50:
            maxrooms[i] = 1
        elif s[i] <= 100:
            maxrooms[i] = 2
        elif s[i] <= 400:
            maxrooms[i] = 7
        elif s[i] <= 700:
            maxrooms[i] = 9
        else:
            maxrooms[i] = 12
        c2 = LinExpr()
        c4 = LinExpr()
        for l in range(p):
            c1 = LinExpr()
            c1 = LinExpr()
            c3 = LinExpr()
            for k in range(r):
                if T[k][l] == 1 and location[k] in w[i]:
                    c1.addTerms(1, x[i, k, l])
                    c4.addTerms(c[k],x[i,k,l])
            obj += c1
            model.addConstr(c1 <= maxrooms[i]* y[i,l], "c1a")
            model.addConstr(c1 >= y[i,l], "C1b")

            for j in conflicts[i]:
                c3.addTerms(1,y[j,l])
            model.addConstr(c3 <= (1 - y[i,l])*sumconflicts[i], "c3")

            c2.addTerms(1,y[i,l])
        model.addConstr( c2 == 1 , "c2")
        model.addConstr(c4 >= s[i], "c4")

    sumrooms = {}
    for l in range(p):
        sumrooms[l] = 0
        cover_inequalities = LinExpr()
        for k in range(r):   
            if T[k][l] == 1:
                sumrooms[l] += 1
                c5 = LinExpr()
                for i in range(n):
                    if location[k] in w[i]:
                        c5.addTerms(1,x[i,k,l])
                model.addConstr( c5 <= 1, "c5")  
                cover_inequalities += c5
        model.addConstr(cover_inequalities <= sumrooms[l], "cover_inequalities")


    model.setObjective( obj, GRB.MINIMIZE)
#.........这里部分代码省略.........

    def _initialize_model(self):

        problem = Model()
        problem.setParam('OutputFlag', False)

        # edges from source to reviewers, capacity controls maximum reviewer load
        self._source_vars = problem.addVars(self.numrev, vtype=GRB.CONTINUOUS, lb=0.0,
                                            ub=self.ability, name='reviewers')

        # edges from papers to sink, capacity controls a number of reviewers per paper
        self._sink_vars = problem.addVars(self.numpapers, vtype=GRB.CONTINUOUS, lb=0.0,
                                          ub=self.demand, name='papers')

        # edges between reviewers and papers. Initially capacities are set to 0 (no edge is added in the network)
        self._mix_vars = problem.addVars(self.numrev, self.numpapers, vtype=GRB.CONTINUOUS,
                                         lb=0.0, ub=0.0, name='assignment')
        problem.update()

        # flow balance equations for reviewers' nodes
        self._balance_reviewers = problem.addConstrs((self._source_vars[i] == self._mix_vars.sum(i, '*')
                                                      for i in range(self.numrev)))

        # flow balance equations for papers' nodes
        self._balance_papers = problem.addConstrs((self._sink_vars[i] == self._mix_vars.sum('*', i)
                                                   for i in range(self.numpapers)))
        problem.update()

        self._problem = problem

def constantino(A, C, k, gap):
    """
    Polynomial-sized CCMcP Edge-Extended Model
    See Constantino et al. (2013)
    """
    t_0 = time.clock()
    _ = '*'
    m = Model()
    m.modelsense = GRB.MAXIMIZE
    m.params.mipgap = gap
    # m.params.timelimit = 60 * 60
    # m.params.nodefilestart = 1.0
    # m.params.nodefiledir = './.nodefiledir'
    # m.params.presparsify = 0
    # m.params.presolve = 0

    n = A.shape[0]
    vars = {}
    edges = tuplelist()

    print('[%.1f] Generating variables...' % (time.clock() - t_0))

    # Variables
    for l in range(n):
        for i in range(l, n):
            for j in range(l, n):
                if A[i, j] == 1:
                    e = (l, i, j)
                    edges.append(e)
                    w = 2 if j in C else 1
                    var = m.addVar(vtype=GRB.BINARY, obj=w)
                    vars[e] = var

        if l % 100 == 0 and l != 0:
            print('[%.1f] l = %d' % (time.clock() - t_0, l))

    m.update()

    print('[%.1f] Generated variables' % (time.clock() - t_0))
    print('[%.1f] Adding flow constraints...' % (time.clock() - t_0))

    # Constraint (2): Flow in = Flow out
    for l in range(n):
        for i in range(l, n):
            # Flow in
            lhs_vars = [vars[e] for e in edges.select(l, _, i)]
            ones = [1.0]*len(lhs_vars)
            lhs = LinExpr()
            lhs.addTerms(ones, lhs_vars)

            # Flow out
            rhs_vars = [vars[e] for e in edges.select(l, i, _)]
            ones = [1.0]*len(rhs_vars)
            rhs = LinExpr()
            rhs.addTerms(ones, rhs_vars)

            # Flow in = Flow out
            m.addConstr(lhs == rhs)

        if l % 100 == 0 and l != 0:
            print('[%.1f] l = %d' % (time.clock() - t_0, l))

    print('[%.1f] Added flow constraints' % (time.clock() - t_0))
    print('[%.1f] Adding cycle vertex constraints...' % (time.clock() - t_0))

    # Constraint (3): Use a vertex only once per cycle
    for i in range(n):
        c_vars = [vars[e] for e in edges.select(_, i, _)]
        ones = [1.0]*len(c_vars)
        expr = LinExpr()
        expr.addTerms(ones, c_vars)
        m.addConstr(expr <= 1.0)

        if i % 100 == 0 and i != 0:
            print('[%.1f] V_i = %d' % (time.clock() - t_0, i))

    print('[%.1f] Added cycle vertex constraints' % (time.clock() - t_0))
    print('[%.1f] Adding cycle cardinality constraints...' % (time.clock() - t_0))

    # Constraint (4): Limit cardinality of cycles to k
    for l in range(n):
        c_vars = [vars[e] for e in edges.select(l, _, _)]
        ones = [1.0]*len(c_vars)
        expr = LinExpr()
        expr.addTerms(ones, c_vars)
        m.addConstr(expr <= k)

        if l % 100 == 0 and l != 0:
            print('[%.1f] l = %d' % (time.clock() - t_0, l))

    print('[%.1f] Added cycle cardinality constraints' % (time.clock() - t_0))
    print('[%.1f] Adding cycle index constraints...' % (time.clock() - t_0))

    # Constraint (5): Cycle index is smallest vertex-index
    for l in range(n):
        rhs_vars = [vars[e] for e in edges.select(l, l, _)]
        ones = [1.0]*len(rhs_vars)
        rhs = LinExpr()
        rhs.addTerms(ones, rhs_vars)

#.........这里部分代码省略.........

def two_cycle(A, C, gap):
    """
    Solve high-vertex dense graphs by reduction to
    weighted matching ILP.
    """
    _ = '*'
    m = Model()
    m.modelsense = GRB.MAXIMIZE
    m.params.mipgap = gap
    m.params.timelimit = 60 * 60

    n = A.shape[0]
    vars = {}
    edges = tuplelist()

    # model as undirected graph
    for i in range(n):
        for j in range(i+1, n):
            if A[i, j] == 1 and A[j, i] == 1:
                e = (i, j)
                edges.append(e)
                w_i = 2 if i in C else 1
                w_j = 2 if j in C else 1
                w = w_i + w_j
                var = m.addVar(vtype=GRB.BINARY, obj=w)
                vars[e] = var

    m.update()

    # 2 cycle constraint <=> undirected flow <= 1
    for i in range(n):
        lhs = LinExpr()
        lhs_vars = [vars[e] for e in chain(edges.select(i, _), edges.select(_, i))]
        ones = [1.0]*len(lhs_vars)
        lhs.addTerms(ones, lhs_vars)

        m.addConstr(lhs <= 1)

    m.optimize()
    m.update()

    cycles = [list(e) for e in edges if vars[e].x == 1.0]
    return cycles, m.objval

class BFPBackupNetwork_Continuous(object):
    """ Class object for buffered failure probability-based model.

    Parameters
    ----------
    Gamma: importance sampling vector
    Nodes: set of nodes
    Links: set of links
    Capacity: capacities per link based based on random failures 
    Survivability: desired survivabiliy factor (epsilon)
    K: number of random scenarios
    
    Returns
    -------
    BackupCapacity: set of capacity per backup link. 
    BackupRoutes: set of backup links
    
    """

    # Private model object
    model = []

    # Private model variables
    BackupCapacity = {}
    bBackupLink = {}
    z0 = {}
    z = {}

    def __init__(self):
        """
        Constructor
        """

    def LoadModel(self, Gamma, Nodes, Links, Capacity, Survivability, NumSamples):
        """ Load model.
    
        Parameters
        ----------
        Gamma : importance sampling vector
        
        """
        self.Links = tuplelist(Links)
        self.Capacity = Capacity

        # Create optimization model
        self.model = Model("Backup")

        # Create variables
        for i, j in self.Links:
            self.BackupCapacity[i, j] = self.model.addVar(
                vtype=GRB.CONTINUOUS, lb=0, name="Backup_Capacity[%s,%s]" % (i, j)
            )
            # self.BackupCapacity[i,j] = self.model.addVar(lb=0, name='Backup_Capacity[%s,%s]' % (i, j))
        self.model.update()

        for i, j in self.Links:
            for s, d in self.Links:
                self.bBackupLink[i, j, s, d] = self.model.addVar(
                    vtype=GRB.BINARY, name="Backup_Link[%s,%s,%s,%s]" % (i, j, s, d)
                )
        self.model.update()

        for i, j in self.Links:
            for k in range(NumSamples):
                self.z[k, i, j] = self.model.addVar(lb=0, name="z[%s][%s][%s]" % (k, i, j))
        self.model.update()

        for i, j in self.Links:
            self.z0[i, j] = self.model.addVar(lb=-GRB.INFINITY, name="z0[%s][%s]" % (i, j))
        self.model.update()

        self.model.modelSense = GRB.MINIMIZE

        self.model.setObjective(quicksum(self.BackupCapacity[i, j] for i, j in self.Links))
        self.model.update()

        # ------------------------------------------------------------------------#
        #                    Constraints definition                              #
        #                                                                        #
        #                                                                        #
        # ------------------------------------------------------------------------#

        # Buffer probability I
        for i, j in self.Links:
            self.model.addConstr(
                self.z0[i, j]
                + 1 / (NumSamples * Survivability) * quicksum(self.z[k, i, j] for (k) in range(NumSamples))
                <= 0,
                "[CONST]Buffer_Prob_I[%s][%s]" % (i, j),
            )
        self.model.update()

        # Link capacity constraints
        for i, j in self.Links:
            for k in range(NumSamples):
                if Gamma == None:
                    self.model.addConstr(
                        (
                            quicksum(self.bBackupLink[i, j, s, d] * Capacity[k, s, d] for s, d in self.Links)
                            - self.BackupCapacity[i, j]
#.........这里部分代码省略.........