为什么glpk将线性规划公式视为非线性规划？

huangapple
go评论124阅读模式
英文:

Why glpk sees a linear programming formulation as a nonlinear programming?

问题

I'm new to both Python and Pyomo. Nonetheless, I've made the below code following my previous question. However, if I run it using glpk solver, it says

> RuntimeError: Cannot write legal LP file. Objective 'cost' has nonlinear terms that are not quadratic.

I think the formulation is an LP. I'm guessing it is the calculation of parameter "demandsperyear"? But, it shouldn't be the cause as it is a parameter, not a variable, right?

If i use cbc, I've got errors that I cannot understand.

  1. from pyomo.environ import *
  2. m = AbstractModel() 
  4. m.techs = Set()
  5. m.products = Set()
  6. m.suppliers = Set()
  7. m.demands = Set()
  8. m.years = Set()
  10. m.initial_demands = Param(m.demands,m.products)
  11. m.growths = Param(m.demands,m.products)
  12. m.yields = Param(m.techs,m.products,m.suppliers)
  13. m.CO2emissions = Param(m.techs,m.suppliers)
  14. m.TAC = Param(m.techs,m.suppliers)
  15. m.max_installed = Param(m.techs,m.suppliers)
  17. m.installed_techs = Var(m.techs, m.suppliers, m.years, within=Binary)
  18. m.number_installed_techs= Var(m.techs, m.suppliers, m.years, bounds=(0,None))
  19. m.totalTAC = Var(bounds=(None,None))
  20. m.totalCO2 = Var(bounds=(None,None))
  21. m.totalTACperyear = Var(m.years, bounds=(None, None))
  22. m.totalCO2peryear = Var(m.years, bounds=(None, None))
  24. def total_TAC(m):
  25.     return m.totalTAC == sum(m.totalTACperyear[year] for year in m.years)
  26. m.cost = Objective(rule=total_TAC, sense=minimize)
  28. def total_CO2(m):
  29.     return m.totalCO2 == sum(m.totalCO2peryear[year] for year in m.years)
  30. m.totalCO2_cons = Constraint(rule=total_CO2)
  32. def total_TACperyear(m,year):
  33.     return m.totalTACperyear[year] == sum(m.number_installed_techs[tech,supplier,year]*m.TAC[tech,supplier] for tech in m.techs for supplier in m.suppliers)
  35. m.TACperyear_cons = Constraint(m.years, rule=total_TACperyear)
  37. def total_CO2peryear(m,year):
  38.     return m.totalCO2peryear[year] == sum(m.number_installed_techs[tech,supplier,year]*m.CO2emissions[tech,supplier] for tech in m.techs for supplier in m.suppliers) 
  40. m.CO2peryear_cons = Constraint(m.years, rule=total_CO2peryear)
  42. def demandperyear(m,demand,product,year):
  43.     if year == 1:
  44.         expression = m.initial_demands[demand,product]
  45.         return expression
  46.     expression = m.initial_demands[demand,product]*(1 + m.growths[demand,product]*year)
  47.     return expression
  49. m.demandsperyear = Param(m.demands, m.products,m.years, initialize=demandperyear)
  51. def fulfill_demands(m,supplier,demand,product,year):
  52.     if supplier == demand:
  53.         expression = sum(m.number_installed_techs[tech,supplier,year]*m.yields[tech,product,supplier] for tech in m.techs) >= m.demandsperyear[demand,product,year]
  54.         return expression
  55.     return Constraint.Skip
  57. m.fulfill_demands_cons = Constraint(m.suppliers,m.demands,m.products,m.years, rule=fulfill_demands)
  59. def max_installed(m,tech,supplier):
  60.     expression = sum(m.number_installed_techs[tech,supplier,year] for year in m.years) <= m.max_installed[tech,supplier]
  61.     return expression
  63. m.max_installed_cons = Constraint(m.techs,m.suppliers, rule=max_installed)
  65. def accum_installed(m,tech,supplier,year):
  66.     if year == 1:
  67.         expression = m.number_installed_techs[tech,supplier,year] == m.installed_techs[tech,supplier,year]
  68.         return expression
  69.     expression = m.number_installed_techs[tech,supplier,year] == m.number_installed_techs[tech,supplier,year-1] + m.installed_techs[tech,supplier,year]
  70.     return expression
  72. m.accum_cons = Constraint(m.techs,m.suppliers,m.years, rule=accum_installed)
  74. # parameter
  75. data = {None: {
  76.     'techs' : {None: ['tek_A', 'tek_B', 'tek_C']},
  77.     'products' : {None: ['prod_A', 'prod_B', 'prod_C']},
  78.     'suppliers' : {None: ['daerah_A', 'daerah_B', 'daerah_C']},
  79.     'demands' : {None: ['daerah_A', 'daerah_B', 'daerah_C']},
  80.     'years' : {None: [1, 2, 3, 4, 5,6,7,8,9,10]},
  81.     'initial_demands': {('daerah_A','prod_A'): 10.0, ('daerah_A','prod_B'): 10.0, ('daerah_A','prod_C'): 10.0,
  82.                         ('daerah_B','prod_A'): 20.0, ('daerah_B','prod_B'): 20.0, ('daerah_B','prod_C'): 20.0,
  83.                         ('daerah_C','prod_A'): 12.0, ('daerah_C','prod_B'): 12.0, ('daerah_C','prod_C'): 12.0},
  84.     'growths': {('daerah_A','prod_A'): 0.02, ('daerah_A','prod_B'): 0.02, ('daerah_A','prod_C'): 0.02,
  85.                         ('daerah_B','prod_A'): 
  87. <details>
  88. <summary>英文:</summary>
  90. I&#39;m new to both Python and Pyomo. Nonetheless, I&#39;ve made the below code following my previous question. However, if I run it using glpk solver, it says 
  91. &gt; RuntimeError: Cannot write legal LP file.  Objective &#39;cost&#39; has nonlinear terms that are not quadratic.
  93. I think the formulation is an LP. I&#39;m guessing it is the calculation of parameter &quot;demandsperyear&quot;? But, it shouldn&#39;t be the cause as it is a parameter, not a variable, right? 
  95. If i use cbc, I&#39;ve got errors that I cannot understand.

from pyomo.environ import *
m = AbstractModel()

m.techs = Set()
m.products = Set()
m.suppliers = Set()
m.demands = Set()
m.years = Set()

m.initial_demands = Param(m.demands,m.products)
m.growths = Param(m.demands,m.products)
m.yields = Param(m.techs,m.products,m.suppliers)
m.CO2emissions = Param(m.techs,m.suppliers)
m.TAC = Param(m.techs,m.suppliers)
m.max_installed = Param(m.techs,m.suppliers)

m.installed_techs = Var(m.techs, m.suppliers, m.years, within=Binary)
m.number_installed_techs= Var(m.techs, m.suppliers, m.years, bounds=(0,None))
m.totalTAC = Var(bounds=(None,None))
m.totalCO2 = Var(bounds=(None,None))
m.totalTACperyear = Var(m.years, bounds=(None, None))
m.totalCO2peryear = Var(m.years, bounds=(None, None))

def total_TAC(m):
return m.totalTAC == sum(m.totalTACperyear[year] for year in m.years)
m.cost = Objective(rule=total_TAC, sense=minimize)

def total_CO2(m):
return m.totalCO2 == sum(m.totalCO2peryear[year] for year in m.years)
m.totalCO2_cons = Constraint(rule=total_CO2)

def total_TACperyear(m,year):
return m.totalTACperyear[year] == sum(m.number_installed_techs[tech,supplier,year]*m.TAC[tech,supplier] for tech in m.techs for supplier in m.suppliers)

m.TACperyear_cons = Constraint(m.years, rule=total_TACperyear)

def total_CO2peryear(m,year):
return m.totalCO2peryear[year] == sum(m.number_installed_techs[tech,supplier,year]*m.CO2emissions[tech,supplier] for tech in m.techs for supplier in m.suppliers)

m.CO2peryear_cons = Constraint(m.years, rule=total_CO2peryear)

def demandperyear(m,demand,product,year):
if year == 1:
expression = m.initial_demands[demand,product]
return expression
expression = m.initial_demands[demand,product]*(1 + m.growths[demand,product]*year)
return expression

m.demandsperyear = Param(m.demands, m.products,m.years, initialize=demandperyear)

def fulfill_demands(m,supplier,demand,product,year):
if supplier == demand:
expression = sum(m.number_installed_techs[tech,supplier,year]*m.yields[tech,product,supplier] for tech in m.techs) >= m.demandsperyear[demand,product,year]
return expression
return Constraint.Skip

m.fulfill_demands_cons = Constraint(m.suppliers,m.demands,m.products,m.years, rule=fulfill_demands)

def max_installed(m,tech,supplier):
expression = sum(m.number_installed_techs[tech,supplier,year] for year in m.years) <= m.max_installed[tech,supplier]
return expression

m.max_installed_cons = Constraint(m.techs,m.suppliers, rule=max_installed)

def accum_installed(m,tech,supplier,year):
if year == 1:
expression = m.number_installed_techs[tech,supplier,year] == m.installed_techs[tech,supplier,year]
return expression
expression = m.number_installed_techs[tech,supplier,year] == m.number_installed_techs[tech,supplier,year-1] + m.installed_techs[tech,supplier,year]
return expression

m.accum_cons = Constraint(m.techs,m.suppliers,m.years, rule=accum_installed)

parameter

data = {None: {
'techs' : {None: ['tek_A', 'tek_B', 'tek_C']},
'products' : {None: ['prod_A', 'prod_B', 'prod_C']},
'suppliers' : {None: ['daerah_A', 'daerah_B', 'daerah_C']},
'demands' : {None: ['daerah_A', 'daerah_B', 'daerah_C']},
'years' : {None: [1, 2, 3, 4, 5,6,7,8,9,10]},
'initial_demands': {('daerah_A','prod_A'): 10.0, ('daerah_A','prod_B'): 10.0, ('daerah_A','prod_C'): 10.0,
('daerah_B','prod_A'): 20.0, ('daerah_B','prod_B'): 20.0, ('daerah_B','prod_C'): 20.0,
('daerah_C','prod_A'): 12.0, ('daerah_C','prod_B'): 12.0, ('daerah_C','prod_C'): 12.0},
'growths': {('daerah_A','prod_A'): 0.02, ('daerah_A','prod_B'): 0.02, ('daerah_A','prod_C'): 0.02,
('daerah_B','prod_A'): 0.02, ('daerah_B','prod_B'): 0.02, ('daerah_B','prod_C'): 0.02,
('daerah_C','prod_A'): 0.02, ('daerah_C','prod_B'): 0.02, ('daerah_C','prod_C'): 0.02},
'max_installed' : {('tek_A','daerah_A'): 100.0, ('tek_A','daerah_B'): 100.0, ('tek_A','daerah_C'): 100.0,
('tek_B','daerah_A'): 100.0, ('tek_B','daerah_B'): 100.0, ('tek_B','daerah_C'): 100.0,
('tek_C','daerah_A'): 100.0, ('tek_C','daerah_B'): 100.0, ('tek_C','daerah_C'): 100.0},
'TAC' : {('tek_A','daerah_A'): 100.0, ('tek_A','daerah_B'): 100.0, ('tek_A','daerah_C'): 100.0,
('tek_B','daerah_A'): 100.0, ('tek_B','daerah_B'): 100.0, ('tek_B','daerah_C'): 100.0,
('tek_C','daerah_A'): 100.0, ('tek_C','daerah_B'): 100.0, ('tek_C','daerah_C'): 100.0},
'CO2emissions' : {('tek_A','daerah_A'): 20.0, ('tek_A','daerah_B'): 20.0, ('tek_A','daerah_C'): 20.0,
('tek_B','daerah_A'): 20.0, ('tek_B','daerah_B'): 20.0, ('tek_B','daerah_C'): 20.0,
('tek_C','daerah_A'): 20.0, ('tek_C','daerah_B'): 20.0, ('tek_C','daerah_C'): 20.0},
'yields': {('tek_A','prod_A','daerah_A'):8.0, ('tek_A','prod_A','daerah_B'):8.0, ('tek_A','prod_A','daerah_C'):8.0,
('tek_A','prod_B','daerah_A'):8.0, ('tek_A','prod_B','daerah_B'):8.0, ('tek_A','prod_B','daerah_C'):8.0,
('tek_A','prod_C','daerah_A'):8.0, ('tek_A','prod_C','daerah_B'):8.0, ('tek_A','prod_C','daerah_C'):8.0,
('tek_B','prod_A','daerah_A'):8.0, ('tek_B','prod_A','daerah_B'):8.0, ('tek_B','prod_A','daerah_C'):8.0,
('tek_B','prod_B','daerah_A'):8.0, ('tek_B','prod_B','daerah_B'):8.0, ('tek_B','prod_B','daerah_C'):8.0,
('tek_B','prod_C','daerah_A'):8.0, ('tek_B','prod_C','daerah_B'):8.0, ('tek_B','prod_C','daerah_C'):8.0,
('tek_C','prod_A','daerah_A'):8.0, ('tek_C','prod_A','daerah_B'):8.0, ('tek_C','prod_A','daerah_C'):8.0,
('tek_C','prod_B','daerah_A'):8.0, ('tek_C','prod_B','daerah_B'):8.0, ('tek_C','prod_B','daerah_C'):8.0,
('tek_C','prod_C','daerah_A'):8.0, ('tek_C','prod_C','daerah_B'):8.0, ('tek_C','prod_C','daerah_C'):8.0}
}}

instance = m.create_instance(data)
instance.pprint()

solver = SolverFactory('glpk')
hasil = solver.solve(instance)

for tech in instance.installed_techs.values():
print(tech.value)

for number_of_techs_installed in instance.number_installed_techs.values():
print(number_of_techs_installed.value)

for TACperyear in instance.totalTACperyear.values():
print(TACperyear.value)

for CO2peryear in instance.totalCO2peryear.values():
print(CO2peryear.value)

print('Total TAC', instance.totalTAC())
print('Total CO2', instance.totalCO2())

  2. If I have to calculate the &quot;demandsperyear&quot; outside of the Pyomo, how should i do that?
  3. Thanks
  4. </details>
  5. # 答案1
  6. **得分**: 0
  7. Here is the translation of the code portions you provided:
  8. ```python
  9. 要获得某些东西，如@user2357112建议的，我需要解决`total_TAC`上的问题：
  10. def total_TAC(m):
  11. # return m.totalTAC == sum(m.totalTACperyear[year] for year in m.years)
  12. return sum(m.totalTACperyear[year] for year in m.years)
  13. 然后，在创建实例之前，我计算了`demandsperyear`：
  14. def demandperyear_outside(demand, product, year, initial_demands, growths):
  15. if year == 1:
  16. return initial_demands[demand, product]
  17. # 我假设接下来的年份采用复利增长，也许您想要更改这个
  18. return initial_demands[demand, product] * (1 + growths[demand, product] * year)
  19. demandsperyear_data = {
  20. (demand, product, year): demandperyear_outside(
  21. demand, product, year, data[None]['initial_demands'], data[None]['growths']) 
  22. for demand in data[None]['demands'][None] 
  23. for product in data[None]['products'][None] 
  24. for year in data[None]['years'][None]
  25. }
  26. 然后，我将其添加到`data`中：
  27. data[None]['demandsperyear'] = demandsperyear_data
  28. 通过这个操作，我将您的`m.demandsperyear`更改为：
  29. m.demandsperyear = Param(m.demands, m.products, m.years, initialize=data[None]['demandsperyear'])
  30. 我得到了一个**总CO2**（1400），但很遗憾**总TAC**是`None`。
  31. **更新我的第二次尝试**：我添加了一个`totalTAC`的约束条件，如下所示：
  32. def total_TAC_constraint(m):
  33. return m.totalTAC == sum(m.totalTACperyear[year] for year in m.years)
  34. m.totalTAC_cons = Constraint(rule=total_TAC_constraint)
  35. 现在我得到了结果：
  36. 总TAC 7000.0
  37. 总CO2 1400.0
英文:

To get something I had to solve, as @user2357112 suggested, the problem on total_TAC:

  1. def total_TAC(m):
  2. # return m.totalTAC == sum(m.totalTACperyear[year] for year in m.years)
  3. return sum(m.totalTACperyear[year] for year in m.years)

Then I calculated demandsperyear before creating the instance:

  1. def demandperyear_outside(demand, product, year, initial_demands, growths):
  2. if year == 1:
  3. return initial_demands[demand, product]
  4. # I assume compound growth for the next years, maybe you want to change this
  5. return initial_demands[demand, product] * (1 + growths[demand, product] * year)
  6. demandsperyear_data = {
  7. (demand, product, year): demandperyear_outside(
  8. demand, product, year, data[None][&#39;initial_demands&#39;], data[None][&#39;growths&#39;]) 
  9. for demand in data[None][&#39;demands&#39;][None] 
  10. for product in data[None][&#39;products&#39;][None] 
  11. for year in data[None][&#39;years&#39;][None]
  12. }

which I add to the data:

  1. data[None][&#39;demandsperyear&#39;] = demandsperyear_data

With this I change your m.demandsperyear to:

  1. m.demandsperyear = Param(m.demands, m.products, m.years, initialize=data[None][&#39;demandsperyear&#39;])

I'm getting a Total CO2 (1400) but sadly a None for Total TAC.

Update for my second attempt: I added a constraint for totalTAC as follows:

  1. def total_TAC_constraint(m):
  2. return m.totalTAC == sum(m.totalTACperyear[year] for year in m.years)
  3. m.totalTAC_cons = Constraint(rule=total_TAC_constraint)

and I'm now getting results:

  1. Total TAC 7000.0
  2. Total CO2 1400.0

答案2

得分: 0

The problem with your code as posted has nothing to do with linear/nonlinear. The problem is in your objective function. Your objective function (below) is written as an equality, as you would in a constraint.

你的问题与你发布的代码无关于线性/非线性。问题在于你的目标函数。你的目标函数（如下）被写成了一个等式，就像在一个约束中一样。

  1. def total_TAC(m):
  2.     return m.totalTAC == sum(m.totalTACperyear[year] for year in m.years)
  3. m.cost = Objective(rule=total_TAC, sense=minimize)

The objective function should be an expression of the variables that is, well, variable. If you substitute in this, all works. Of course, I don't know if this is the mathematical formulation you are looking for in your model, over to you on that. This just minimizes the sum as the objective.

目标函数应该是一个表达变量的表达式，它应该是可变的。如果你用以下方式替换，一切都能正常工作。当然，我不知道这是否是你在模型中寻找的数学表达方式，这取决于你。这只是将总和最小化作为目标。

  1. def total_TAC(m):
  2.     return sum(m.totalTACperyear[year] for year in m.years)
  3. m.cost = Objective(rule=total_TAC, sense=minimize)
英文:

The problem with your code as posted has nothing to do with linear/nonlinear. The problem is in your objective function. Your objective function (below) is written as an equality, as you would in a constraint.

  1. def total_TAC(m):
  2. return m.totalTAC == sum(m.totalTACperyear[year] for year in m.years)
  3. m.cost = Objective(rule=total_TAC, sense=minimize)

The objective function should be an expression of the variables that is, well, variable. If you substitute in this, all works. Of course, I don't know if this is the mathematical formulation you are looking for in your model, over to you on that. This just minimizes the sum as the objective.

  1. def total_TAC(m):
  2. return sum(m.totalTACperyear[year] for year in m.years)
  3. m.cost = Objective(rule=total_TAC, sense=minimize)

huangapple
  • 本文由 发表于 2023年5月7日 17:08:49
  • 转载请务必保留本文链接：https://go.coder-hub.com/76193018.html
  • glpk
  • pyomo
  • python
匿名

发表评论

匿名网友

:?: :razz: :sad: :evil: :!: :smile: :oops: :grin: :eek: :shock: :???: :cool: :lol: :mad: :twisted: :roll: :wink: :idea: :arrow: :neutral: :cry: :mrgreen:

确定