import pyomo.environ as pyo
def build_model(v, a=0):
model = pyo.ConcreteModel()
# creat indices and sets
model.I = pyo.Set(initialize=range(1,51))
model.J = pyo.Set(initialize=list(range(1,201)))
# #input data
# mu = 5
# sigma = 10
# def c_init(model, i, j):
# c_data = np.random.randn(1) * sigma +mu
# model.f = pyo.Param(model.I, initialize = 50)
# model.c = pyo.Param(model.I, model.J, initialize = c_init)
# input data
model.f = pyo.Param(model.I, initialize = 50)
model.c = pyo.Param(model.I, model.J, initialize=v)
# decision varibles
model.x = pyo.Var(model.I, model.J, within=pyo.Binary) # give a initial value is crucial
model.y = pyo.Var(model.I, within=pyo.Binary)
# objective function
def Obj_rule(m):
return pyo.summation(m.f, m.y) + sum(sum(m.c[i,j] * m.x[i,j] for j in m.J) for i in m.I)
model.Obj = pyo.Objective(rule=Obj_rule)
# Constraints
def must_be_served(m, j):
return sum(m.x[i,j] for i in m.I) == 1
model.must_be_served_constraints = pyo.Constraint(model.J, rule=must_be_served)
if a == 0:
def facility_open(m, i):
return sum(m.x[i,j] for j in m.J) <= 100*m.y[i]
model.facility_open_constraint = pyo.Constraint(model.I,rule=facility_open)
else:
def facility_open(m, i,j):
return m.x[i,j]<= m.y[i]
model.facility_open_constraint = pyo.Constraint(model.I, model.J, rule=facility_open)
import random
import pyomo.environ as pyo
from pyomo.opt import SolverFactory
import time
solver = SolverFactory('coptampl')
import Lab5_homework_model1_Q4 as Lab5_homework_model1
def v_init():
v = dict()
for i in list(range(1,51)):
for j in list(range(1,201)):
#v[i,j] = np.random.randn(1)* 5 + 10
v[i,j] = random.random()*5 + 10
return v
for itr in range(5):
v = v_init()
model1 = Lab5_homework_model1.build_model(v,a=0)
model2 = Lab5_homework_model1.build_model(v,a=1)
time_st1 = time.time()
results1 = solver.solve(model1, tee = True)
print(results1)
time_end1 = time.time()
print("The solver returned a status of:" + str(results1.solver.status))
results2 = solver.solve(model2, tee = True)
print(results2)
time_end2 = time.time()
print("The solver returned a status of:" + str(results2.solver.status))
with open('Q4.txt','a') as f:
f.write('ObjectVal of model1 for ' + str(itr+1) + ' times result:' +str(pyo.value(model1.Obj)) + ' used time(sec): ' + str(time_end1-time_st1))
f.write('\r')
f.write('ObjectVal of model2 for ' + str(itr+1) + ' times result:' +str(pyo.value(model2.Obj)) + ' used time(sec): ' + str(time_end2-time_end1))
f.write('\r')
import pyomo.environ as pyodef build_model(v, a=0): model = pyo.ConcreteModel() # creat indices and sets model.I = pyo.Set(initialize=range(1,51)) model.J = pyo.Set(initialize=list(range(1,201))) # #input data # mu = 5 # si.