如何编写Pyomo优化以选择可再生能源的最佳容量?

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

how to write a pyomo optimization to select optimal volume of renewables?

问题

Here's the translation of the code part you provided:

  1. import os
  2. import pandas as pd
  3. from pyomo.environ import *
  4. import datetime
  5. def model_to_df(model, first_period, last_period):
  6. # Need to increase the first & last hour by 1 because of pyomo indexing
  7. periods = range(model.T[first_period + 1], model.T[last_period + 1] + 1)
  8. spot = [value(model.spot_v[i]) for i in periods]
  9. load = [value(model.load_v[i]) for i in periods]
  10. slr1 = [value(model.slr1_size[i]) for i in periods]
  11. slr2 = [value(model.slr2_size[i]) for i in periods]
  12. slr3 = [value(model.slr3_size[i]) for i in periods]
  13. wnd1 = [value(model.wnd1_size[i]) for i in periods]
  14. wnd2 = [value(model.wnd2_size[i]) for i in periods]
  15. wnd3 = [value(model.wnd3_size[i]) for i in periods]
  16. slr1_gen = [value(model.slr1_gen[i]) for i in periods]
  17. slr2_gen = [value(model.slr2_gen[i]) for i in periods]
  18. slr3_gen = [value(model.slr3_gen[i]) for i in periods]
  19. wnd1_gen = [value(model.wnd1_gen[i]) for i in periods]
  20. wnd2_gen = [value(model.wnd2_gen[i]) for i in periods]
  21. wnd3_gen = [value(model.wnd3_gen[i]) for i in periods]
  22. total_gen = [value(model.total_gen[i]) for i in periods]
  23. spill_gen = [value(model.spill_gen[i]) for i in periods]
  24. spill_gen_sum = [value(model.spill_gen_sum[i]) for i in periods]
  25. spill_binary = [value(model.spill_binary[i]) for i in periods]
  26. self_cons = [value(model.self_cons[i]) for i in periods]
  27. df_dict = {
  28. 'Period': periods,
  29. 'spot': spot,
  30. 'load': load,
  31. 'slr1': slr1,
  32. 'slr2': slr2,
  33. 'slr3': slr3,
  34. 'wnd1': wnd1,
  35. 'wnd2': wnd2,
  36. 'wnd3': wnd3,
  37. 'slr1_gen': slr1_gen,
  38. 'slr2_gen': slr2_gen,
  39. 'slr3_gen': slr3_gen,
  40. 'wnd1_gen': wnd1_gen,
  41. 'wnd2_gen': wnd2_gen,
  42. 'wnd3_gen': wnd3_gen,
  43. 'total_gen': total_gen,
  44. 'spill_gen': spill_gen,
  45. 'self_cons': self_cons,
  46. 'spill_gen_sum': spill_gen_sum,
  47. 'spill_binary': spill_binary
  48. }
  49. df = pd.DataFrame(df_dict)
  50. return df
  51. LOCATION = r"C:\cbc-win64"
  52. os.environ["PATH"] = LOCATION + ";" + os.environ["PATH"]
  53. df = pd.DataFrame({
  54. 'hour': [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23],
  55. 'load': [100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100],
  56. 'spot': [65.4, 62.7, 60.9, 60.3, 61.8, 64.5, 65.9, 57.9, 39.7, 28.3, 20.9, 16.3, 18.1, 23.9, 32.3, 43.2, 59.3, 76.3, 80.5, 72.5, 73.1, 69.0, 67.9, 67.7],
  57. 'slr1': [0.00, 0.00, 0.00, 0.00, 0.00, 0.04, 0.20, 0.44, 0.60, 0.69, 0.71, 0.99, 1.00, 0.66, 0.75, 0.63, 0.52, 0.34, 0.14, 0.02, 0.00, 0.00, 0.00, 0.00],
  58. <details>
  59. <summary>英文:</summary>
  60. **Background**
  61. I am trying to write a pyomo optimization which takes in a customer&#39;s electricity load and the generation data of several renewable projects, then optimally solves for the lowest cost selection of renewable projects to minimize electricity consumption, subject to a few constraints.
  62. **What I&#39;ve tried**
  63. Using the pyomo readthedocs and stackoverflow. I&#39;ve written my first attempt (below), but I am having two issues.
  64. **Problem**
  65. 1. ERROR: Rule failed for Expression &#39;d_spill_var&#39; with index 0: PyomoException:
  66. Cannot convert non-constant Pyomo expression
  67. I think this is because I am trying to return a max(expr, 0) value for one of my dependent Expresions. However even if I change this, I still get issue 2 below;
  68. 2. RuntimeError: Cannot write legal LP file. Objective &#39;objective&#39; has nonlinear terms that are not quadratic.
  69. **Help Requested**
  70. Could someone please point me in the right direction to solving the above two issues? Any help would be greatly appreciated!
  71. **Code OLD v1.0**
  72. import os
  73. import pandas as pd
  74. from pyomo.environ import *
  75. import datetime
  76. def model_to_df(model, first_period, last_period):
  77. # Need to increase the first &amp; last hour by 1 because of pyomo indexing
  78. periods = range(model.T[first_period + 1], model.T[last_period + 1] + 1)
  79. spot = [value(model.spot[i]) for i in periods]
  80. load = [value(model.load[i]) for i in periods]
  81. slr1 = [value(model.slr1_size[i]) for i in periods]
  82. slr2 = [value(model.slr2_size[i]) for i in periods]
  83. slr3 = [value(model.slr3_size[i]) for i in periods]
  84. wnd1 = [value(model.wnd1_size[i]) for i in periods]
  85. wnd2 = [value(model.wnd2_size[i]) for i in periods]
  86. wnd3 = [value(model.wnd3_size[i]) for i in periods]
  87. d_slrgen_var = [value(model.d_slrgen_var[i]) for i in periods]
  88. d_wndgen_var = [value(model.d_wndgen_var[i]) for i in periods]
  89. d_spill_var = [value(model.d_spill_var[i]) for i in periods]
  90. d_selfcons_var = [value(model.d_selfcons_var[i]) for i in periods]
  91. df_dict = {
  92. &#39;Period&#39;: periods,
  93. &#39;spot&#39;: spot,
  94. &#39;load&#39;: load,
  95. &#39;slr1&#39;: slr1,
  96. &#39;slr2&#39;: slr2,
  97. &#39;slr3&#39;: slr3,
  98. &#39;wnd1&#39;: wnd1,
  99. &#39;wnd2&#39;: wnd2,
  100. &#39;wnd3&#39;: wnd3,
  101. &#39;d_slrgen_var&#39;: d_slrgen_var,
  102. &#39;d_wndgen_var&#39;: d_wndgen_var,
  103. &#39;d_spill_var&#39;: d_spill_var,
  104. &#39;d_selfcons_var&#39;: d_selfcons_var
  105. }
  106. df = pd.DataFrame(df_dict)
  107. return df
  108. LOCATION = r&quot;C:\cbc-win64&quot;
  109. os.environ[&quot;PATH&quot;] = LOCATION + &quot;;&quot; + os.environ[&quot;PATH&quot;]
  110. df = pd.DataFrame({
  111. &#39;hour&#39;: [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23],
  112. &#39;load&#39;: [100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100],
  113. &#39;spot&#39;: [65.4, 62.7, 60.9, 60.3, 61.8, 64.5, 65.9, 57.9, 39.7, 28.3, 20.9, 16.3, 18.1, 23.9, 32.3, 43.2, 59.3, 76.3, 80.5, 72.5, 73.1, 69.0, 67.9, 67.7],
  114. &#39;slr1&#39;: [0.00, 0.00, 0.00, 0.00, 0.00, 0.04, 0.20, 0.44, 0.60, 0.69, 0.71, 0.99, 1.00, 0.66, 0.75, 0.63, 0.52, 0.34, 0.14, 0.02, 0.00, 0.00, 0.00, 0.00],
  115. &#39;slr2&#39;: [0.00, 0.00, 0.00, 0.00, 0.03, 0.19, 0.44, 0.68, 1.00, 0.83, 0.90, 0.88, 0.98, 0.94, 0.83, 0.70, 0.36, 0.11, 0.02, 0.00, 0.00, 0.00, 0.00, 0.00],
  116. &#39;slr3&#39;: [0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.03, 0.17, 0.39, 0.87, 0.91, 1.00, 0.89, 0.71, 0.71, 0.85, 0.63, 0.52, 0.32, 0.12, 0.02, 0.00, 0.00, 0.00],
  117. &#39;wnd1&#39;: [1.00, 0.72, 0.74, 0.94, 0.69, 0.90, 0.92, 0.76, 0.51, 0.35, 0.31, 0.34, 0.37, 0.28, 0.35, 0.40, 0.39, 0.32, 0.42, 0.48, 0.74, 0.63, 0.80, 0.97],
  118. &#39;wnd2&#39;: [0.95, 0.67, 0.82, 0.48, 0.51, 0.41, 0.33, 0.42, 0.34, 0.30, 0.39, 0.29, 0.34, 0.55, 0.67, 0.78, 0.84, 0.73, 0.77, 0.89, 0.76, 0.97, 1.00, 0.91],
  119. &#39;wnd3&#39;: [0.32, 0.35, 0.38, 0.28, 0.33, 0.38, 0.41, 0.38, 0.51, 0.65, 0.54, 0.88, 0.93, 0.89, 0.90, 1.00, 0.90, 0.76, 0.76, 0.92, 0.71, 0.56, 0.52, 0.40]
  120. })
  121. first_model_period = df[&#39;hour&#39;].iloc[0]
  122. last_model_period = df[&#39;hour&#39;].iloc[-1]
  123. # **********************
  124. # Build Model
  125. # **********************
  126. model = ConcreteModel()
  127. # Fixed Paramaters
  128. model.T = Set(initialize=df.index.tolist(), doc=&#39;hourly intervals&#39;, ordered=True)
  129. model.load_v = Param(model.T, initialize=df.load, doc=&#39;customers load&#39;, within=Any)
  130. model.spot_v = Param(model.T, initialize=df.spot, doc=&#39;spot price for each interval&#39;, within=Any)
  131. model.slr1 = Param(model.T, initialize=df.slr1, doc=&#39;1MW output solar farm 1&#39;, within=Any)
  132. model.slr2 = Param(model.T, initialize=df.slr2, doc=&#39;1MW output solar farm 2&#39;, within=Any)
  133. model.slr3 = Param(model.T, initialize=df.slr3, doc=&#39;1MW output solar farm 3&#39;, within=Any)
  134. model.wnd1 = Param(model.T, initialize=df.wnd1, doc=&#39;1MW output wind farm 1&#39;, within=Any)
  135. model.wnd2 = Param(model.T, initialize=df.wnd2, doc=&#39;1MW output wind farm 2&#39;, within=Any)
  136. model.wnd3 = Param(model.T, initialize=df.wnd3, doc=&#39;1MW output wind farm 3&#39;, within=Any)
  137. # Variable Parameters
  138. model.slr1_flag = Var(model.T, doc=&#39;slr 1 on / off&#39;, within=Binary, initialize=0)
  139. model.slr2_flag = Var(model.T, doc=&#39;slr 2 on / off&#39;, within=Binary, initialize=0)
  140. model.slr3_flag = Var(model.T, doc=&#39;slr 3 on / off&#39;, within=Binary, initialize=0)
  141. model.wnd1_flag = Var(model.T, doc=&#39;wnd 1 on / off&#39;, within=Binary, initialize=0)
  142. model.wnd2_flag = Var(model.T, doc=&#39;wnd 2 on / off&#39;, within=Binary, initialize=0)
  143. model.wnd3_flag = Var(model.T, doc=&#39;wnd 3 on / off&#39;, within=Binary, initialize=0)
  144. model.slr1_size = Var(model.T, bounds=(0, 1500), doc=&#39;selected size in MWs&#39;, initialize=0, within=NonNegativeIntegers)
  145. model.slr2_size = Var(model.T, bounds=(0, 1500), doc=&#39;selected size in MWs&#39;, initialize=0, within=NonNegativeIntegers)
  146. model.slr3_size = Var(model.T, bounds=(0, 1500), doc=&#39;selected size in MWs&#39;, initialize=0, within=NonNegativeIntegers)
  147. model.wnd1_size = Var(model.T, bounds=(0, 1500), doc=&#39;selected size in MWs&#39;, initialize=0, within=NonNegativeIntegers)
  148. model.wnd2_size = Var(model.T, bounds=(0, 1500), doc=&#39;selected size in MWs&#39;, initialize=0, within=NonNegativeIntegers)
  149. model.wnd3_size = Var(model.T, bounds=(0, 1500), doc=&#39;selected size in MWs&#39;, initialize=0, within=NonNegativeIntegers)
  150. model.total_gen = Var(model.T, initialize=0, within=NonNegativeReals)
  151. # Dependent Expression Parameters
  152. def dependent_solar_gen(model, t):
  153. &quot;Total selected solar Generation&quot;
  154. return (model.slr1[t] * model.slr1_flag[t] * model.slr1_size[t]) + \
  155. (model.slr2[t] * model.slr2_flag[t] * model.slr2_size[t]) + \
  156. (model.slr3[t] * model.slr3_flag[t] * model.slr3_size[t])
  157. model.d_slrgen_var = Expression(model.T, rule=dependent_solar_gen)
  158. def dependent_wind_gen(model, t):
  159. &quot;Total selected wind Generation&quot;
  160. return (model.wnd1[t] * model.wnd1_flag[t] * model.wnd1_size[t]) + \
  161. (model.wnd2[t] * model.wnd2_flag[t] * model.wnd2_size[t]) + \
  162. (model.wnd3[t] * model.wnd3_flag[t] * model.wnd3_size[t])
  163. model.d_wndgen_var = Expression(model.T, rule=dependent_wind_gen)
  164. def dependent_spill(model, t):
  165. &quot;Volume of energy not consumed by customer (spilled into grid)&quot;
  166. expr = (model.d_slrgen_var[t] + model.d_wndgen_var[t]) - model.load_v[t]
  167. return max(0, expr)
  168. model.d_spill_var = Expression(model.T, rule=dependent_spill)
  169. def dependent_self_cons(model, t):
  170. &quot;Volume of energy consumed by customer&quot;
  171. expr = (model.d_slrgen_var[t] + model.d_wndgen_var[t]) - model.d_spill_var[t]
  172. return expr
  173. model.d_selfcons_var = Expression(model.T, rule=dependent_self_cons)
  174. # -----------------------
  175. # Constraints
  176. # -----------------------
  177. def min_spill(model, t):
  178. &quot;Limit spill renewables to 10% of total&quot;
  179. return model.d_spill_var[t] &lt;= 0.1 * (model.d_slrgen_var[t] + model.d_wndgen_var[t])
  180. model.min_spill_c = Constraint(model.T, rule=min_spill)
  181. def load_match(model, t):
  182. &quot;contract enough renewables to offset 100% load, even if its not time matched&quot;
  183. return (model.d_slrgen_var[t] + model.d_wndgen_var[t]) &gt;= model.load_v[t]
  184. model.load_match_c = Constraint(model.T, rule=load_match)
  185. # **********************
  186. # Define the income, expenses, and profit
  187. # **********************
  188. green_income = sum(model.spot_v[t] * model.d_spill_var[t] for t in model.T)
  189. black_cost = sum(model.spot_v[t] * (model.load_v[t] - model.d_selfcons_var[t]) for t in model.T)
  190. slr_cost = sum(40 * model.d_slrgen_var[t] for t in model.T)
  191. wnd_cost = sum(70 * model.d_wndgen_var[t] for t in model.T)
  192. profit = green_income - black_cost - slr_cost - wnd_cost
  193. model.objective = Objective(expr=profit, sense=maximize)
  194. # Solve the model
  195. # solver = SolverFactory(&#39;glpk&#39;)
  196. solver = SolverFactory(&#39;cbc&#39;)
  197. solver.solve(model, timelimit=10)
  198. results_df = model_to_df(model, first_period=first_model_period, last_period=last_model_period)
  199. print(results_df)
  200. **Solved**
  201. Thanks @airliquid, I managed to solve the issue thanks to your advice.
  202. What I had to do was linearize the max constraint, redefine dependent expressions as constraints, remove some redundant variables, and change the last two constraints to summations
  203. It&#39;s not the prettiest answer, but it works!
  204. **UPDATED CODE (V2)**
  205. import os
  206. import pandas as pd
  207. from pyomo.environ import *
  208. import datetime
  209. def model_to_df(model, first_period, last_period):
  210. # Need to increase the first &amp; last hour by 1 because of pyomo indexing
  211. periods = range(model.T[first_period + 1], model.T[last_period + 1] + 1)
  212. spot = [value(model.spot_v[i]) for i in periods]
  213. load = [value(model.load_v[i]) for i in periods]
  214. slr1 = [value(model.slr1_size[i]) for i in periods]
  215. slr2 = [value(model.slr2_size[i]) for i in periods]
  216. slr3 = [value(model.slr3_size[i]) for i in periods]
  217. wnd1 = [value(model.wnd1_size[i]) for i in periods]
  218. wnd2 = [value(model.wnd2_size[i]) for i in periods]
  219. wnd3 = [value(model.wnd3_size[i]) for i in periods]
  220. slr1_gen = [value(model.slr1_gen[i]) for i in periods]
  221. slr2_gen = [value(model.slr2_gen[i]) for i in periods]
  222. slr3_gen = [value(model.slr3_gen[i]) for i in periods]
  223. wnd1_gen = [value(model.wnd1_gen[i]) for i in periods]
  224. wnd2_gen = [value(model.wnd2_gen[i]) for i in periods]
  225. wnd3_gen = [value(model.wnd3_gen[i]) for i in periods]
  226. total_gen = [value(model.total_gen[i]) for i in periods]
  227. spill_gen = [value(model.spill_gen[i]) for i in periods]
  228. spill_gen_sum = [value(model.spill_gen_sum[i]) for i in periods]
  229. spill_binary = [value(model.spill_binary[i]) for i in periods]
  230. self_cons = [value(model.self_cons[i]) for i in periods]
  231. df_dict = {
  232. &#39;Period&#39;: periods,
  233. &#39;spot&#39;: spot,
  234. &#39;load&#39;: load,
  235. &#39;slr1&#39;: slr1,
  236. &#39;slr2&#39;: slr2,
  237. &#39;slr3&#39;: slr3,
  238. &#39;wnd1&#39;: wnd1,
  239. &#39;wnd2&#39;: wnd2,
  240. &#39;wnd3&#39;: wnd3,
  241. &#39;slr1_gen&#39;: slr1_gen,
  242. &#39;slr2_gen&#39;: slr2_gen,
  243. &#39;slr3_gen&#39;: slr3_gen,
  244. &#39;wnd1_gen&#39;: wnd1_gen,
  245. &#39;wnd2_gen&#39;: wnd2_gen,
  246. &#39;wnd3_gen&#39;: wnd3_gen,
  247. &#39;total_gen&#39;: total_gen,
  248. &#39;spill_gen&#39;: spill_gen,
  249. &#39;self_cons&#39;: self_cons,
  250. &#39;spill_gen_sum&#39;: spill_gen_sum,
  251. &#39;spill_binary&#39;: spill_binary
  252. }
  253. df = pd.DataFrame(df_dict)
  254. return df
  255. LOCATION = r&quot;C:\cbc-win64&quot;
  256. os.environ[&quot;PATH&quot;] = LOCATION + &quot;;&quot; + os.environ[&quot;PATH&quot;]
  257. df = pd.DataFrame({
  258. &#39;hour&#39;: [0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23],
  259. &#39;load&#39;: [100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100, 100],
  260. &#39;spot&#39;: [65.4, 62.7, 60.9, 60.3, 61.8, 64.5, 65.9, 57.9, 39.7, 28.3, 20.9, 16.3, 18.1, 23.9, 32.3, 43.2, 59.3, 76.3, 80.5, 72.5, 73.1, 69.0, 67.9, 67.7],
  261. &#39;slr1&#39;: [0.00, 0.00, 0.00, 0.00, 0.00, 0.04, 0.20, 0.44, 0.60, 0.69, 0.71, 0.99, 1.00, 0.66, 0.75, 0.63, 0.52, 0.34, 0.14, 0.02, 0.00, 0.00, 0.00, 0.00],
  262. &#39;slr2&#39;: [0.00, 0.00, 0.00, 0.00, 0.03, 0.19, 0.44, 0.68, 1.00, 0.83, 0.90, 0.88, 0.98, 0.94, 0.83, 0.70, 0.36, 0.11, 0.02, 0.00, 0.00, 0.00, 0.00, 0.00],
  263. &#39;slr3&#39;: [0.00, 0.00, 0.00, 0.00, 0.00, 0.00, 0.03, 0.17, 0.39, 0.87, 0.91, 1.00, 0.89, 0.71, 0.71, 0.85, 0.63, 0.52, 0.32, 0.12, 0.02, 0.00, 0.00, 0.00],
  264. &#39;wnd1&#39;: [1.00, 0.72, 0.74, 0.94, 0.69, 0.90, 0.92, 0.76, 0.51, 0.35, 0.31, 0.34, 0.37, 0.28, 0.35, 0.40, 0.39, 0.32, 0.42, 0.48, 0.74, 0.63, 0.80, 0.97],
  265. &#39;wnd2&#39;: [0.95, 0.67, 0.82, 0.48, 0.51, 0.41, 0.33, 0.42, 0.34, 0.30, 0.39, 0.29, 0.34, 0.55, 0.67, 0.78, 0.84, 0.73, 0.77, 0.89, 0.76, 0.97, 1.00, 0.91],
  266. &#39;wnd3&#39;: [0.32, 0.35, 0.38, 0.28, 0.33, 0.38, 0.41, 0.38, 0.51, 0.65, 0.54, 0.88, 0.93, 0.89, 0.90, 1.00, 0.90, 0.76, 0.76, 0.92, 0.71, 0.56, 0.52, 0.40]
  267. })
  268. # df.to_csv(&#39;example.csv&#39;, index=False)
  269. first_model_period = df[&#39;hour&#39;].iloc[0]
  270. last_model_period = df[&#39;hour&#39;].iloc[-1]
  271. # **********************
  272. # Build Model
  273. # **********************
  274. model = ConcreteModel()
  275. # Fixed Paramaters
  276. model.T = Set(initialize=df.index.tolist(), doc=&#39;hourly intervals&#39;, ordered=True)
  277. model.load_v = Param(model.T, initialize=df.load, doc=&#39;customers load&#39;, within=Any)
  278. model.spot_v = Param(model.T, initialize=df.spot, doc=&#39;spot price for each interval&#39;, within=Any)
  279. model.slr1 = Param(model.T, initialize=df.slr1, doc=&#39;1MW output solar farm 1&#39;, within=Any)
  280. model.slr2 = Param(model.T, initialize=df.slr2, doc=&#39;1MW output solar farm 2&#39;, within=Any)
  281. model.slr3 = Param(model.T, initialize=df.slr3, doc=&#39;1MW output solar farm 3&#39;, within=Any)
  282. model.wnd1 = Param(model.T, initialize=df.wnd1, doc=&#39;1MW output wind farm 1&#39;, within=Any)
  283. model.wnd2 = Param(model.T, initialize=df.wnd2, doc=&#39;1MW output wind farm 2&#39;, within=Any)
  284. model.wnd3 = Param(model.T, initialize=df.wnd3, doc=&#39;1MW output wind farm 3&#39;, within=Any)
  285. # Variable Parameters
  286. model.slr1_size = Var(model.T, bounds=(0, 1000), doc=&#39;selected size in MWs&#39;, initialize=0, within=NonNegativeIntegers)
  287. model.slr2_size = Var(model.T, bounds=(0, 1000), doc=&#39;selected size in MWs&#39;, initialize=0, within=NonNegativeIntegers)
  288. model.slr3_size = Var(model.T, bounds=(0, 1000), doc=&#39;selected size in MWs&#39;, initialize=0, within=NonNegativeIntegers)
  289. model.wnd1_size = Var(model.T, bounds=(0, 1000), doc=&#39;selected size in MWs&#39;, initialize=0, within=NonNegativeIntegers)
  290. model.wnd2_size = Var(model.T, bounds=(0, 1000), doc=&#39;selected size in MWs&#39;, initialize=0, within=NonNegativeIntegers)
  291. model.wnd3_size = Var(model.T, bounds=(0, 1000), doc=&#39;selected size in MWs&#39;, initialize=0, within=NonNegativeIntegers)
  292. model.slr1_gen = Var(model.T, initialize=0, within=NonNegativeReals)
  293. model.slr2_gen = Var(model.T, initialize=0, within=NonNegativeReals)
  294. model.slr3_gen = Var(model.T, initialize=0, within=NonNegativeReals)
  295. model.wnd1_gen = Var(model.T, initialize=0, within=NonNegativeReals)
  296. model.wnd2_gen = Var(model.T, initialize=0, within=NonNegativeReals)
  297. model.wnd3_gen = Var(model.T, initialize=0, within=NonNegativeReals)
  298. model.total_gen = Var(model.T, initialize=0, within=NonNegativeReals)
  299. model.spill_gen_sum = Var(model.T, initialize=0, within=Reals)
  300. model.spill_binary = Var(model.T, doc=&#39;to get max&#39;, within=Binary, initialize=0)
  301. model.spill_gen = Var(model.T, initialize=0, within=NonNegativeReals)
  302. model.self_cons = Var(model.T, initialize=0, within=NonNegativeReals)
  303. # -----------------------
  304. # Constraints
  305. # -----------------------
  306. # SIZE CONSTRAINTS
  307. def slr1_size(model, t):
  308. &quot;slr1 size&quot;
  309. return model.slr1_size[t] == model.slr1_size[1]
  310. model.slr_size1_c = Constraint(model.T, rule=slr1_size)
  311. def slr2_size(model, t):
  312. &quot;slr2 size&quot;
  313. return model.slr2_size[t] == model.slr2_size[1]
  314. model.slr_size2_c = Constraint(model.T, rule=slr2_size)
  315. def slr3_size(model, t):
  316. &quot;slr3 size&quot;
  317. return model.slr3_size[t] == model.slr3_size[1]
  318. model.slr_size3_c = Constraint(model.T, rule=slr3_size)
  319. def wnd1_size(model, t):
  320. &quot;wnd1 size&quot;
  321. return model.wnd1_size[t] == model.wnd1_size[1]
  322. model.wnd_size1_c = Constraint(model.T, rule=wnd1_size)
  323. def wnd2_size(model, t):
  324. &quot;wnd2 size&quot;
  325. return model.wnd2_size[t] == model.wnd2_size[1]
  326. model.wnd_size2_c = Constraint(model.T, rule=wnd2_size)
  327. def wnd3_size(model, t):
  328. &quot;wnd3 size&quot;
  329. return model.wnd3_size[t] == model.wnd3_size[1]
  330. model.wnd_size3_c = Constraint(model.T, rule=wnd3_size)
  331. # GENERATION EXPRESSIONS / CONSTRAINTS
  332. def slr1_gen(model, t):
  333. &quot;solar 1 generation&quot;
  334. return model.slr1_gen[t] == model.slr1[t] * model.slr1_size[t]
  335. model.slr_gen1_c = Constraint(model.T, rule=slr1_gen)
  336. def slr2_gen(model, t):
  337. &quot;solar 2 generation&quot;
  338. return model.slr2_gen[t] == model.slr2[t] * model.slr2_size[t]
  339. model.slr_gen2_c = Constraint(model.T, rule=slr2_gen)
  340. def slr3_gen(model, t):
  341. &quot;solar 3 generation&quot;
  342. return model.slr3_gen[t] == model.slr3[t] * model.slr3_size[t]
  343. model.slr_gen3_c = Constraint(model.T, rule=slr3_gen)
  344. def wnd1_gen(model, t):
  345. &quot;wind 1 generation&quot;
  346. return model.wnd1_gen[t] == model.wnd1[t] * model.wnd1_size[t]
  347. model.wnd_gen1_c = Constraint(model.T, rule=wnd1_gen)
  348. def wnd2_gen(model, t):
  349. &quot;wind 2 generation&quot;
  350. return model.wnd2_gen[t] == model.wnd2[t] * model.wnd2_size[t]
  351. model.wnd_gen2_c = Constraint(model.T, rule=wnd2_gen)
  352. def wnd3_gen(model, t):
  353. &quot;wind 3 generation&quot;
  354. return model.wnd3_gen[t] == model.wnd3[t] * model.wnd3_size[t]
  355. model.wnd_gen3_c = Constraint(model.T, rule=wnd3_gen)
  356. # TOTAL GENERATION
  357. def total_gen(model, t):
  358. &quot;sum of generation&quot;
  359. return model.total_gen[t] == model.slr1_gen[t] + model.slr2_gen[t] + model.slr3_gen[t] + \
  360. model.wnd1_gen[t] + model.wnd2_gen[t] + model.wnd3_gen[t]
  361. model.total_gen_c = Constraint(model.T, rule=total_gen)
  362. # SPILL GENERATION
  363. def spill_gen_sum(model, t):
  364. &quot;X &gt;= x1&quot;
  365. return model.spill_gen_sum[t] == model.total_gen[t] - model.load_v[t]
  366. model.spill_gen_sum_c = Constraint(model.T, rule=spill_gen_sum)
  367. def spill_check_one(model, t):
  368. &quot;X &gt;= x1&quot;
  369. return model.spill_gen[t] &gt;= model.spill_gen_sum[t]
  370. model.spill_check_one_c = Constraint(model.T, rule=spill_check_one)
  371. def spill_check_two(model, t):
  372. &quot;X &gt;= x2&quot;
  373. return model.spill_gen[t] &gt;= 0
  374. model.spill_check_two_c = Constraint(model.T, rule=spill_check_two)
  375. def spill_binary_one(model, t):
  376. &quot;X &lt;= x1 + M(1-y)&quot;
  377. return model.spill_gen[t] &lt;= model.spill_gen_sum[t] + 9999999*(1-model.spill_binary[t])
  378. model.spill_binary_one_c = Constraint(model.T, rule=spill_binary_one)
  379. def spill_binary_two(model, t):
  380. &quot;X &lt;= x2 + My&quot;
  381. return model.spill_gen[t] &lt;= 9999999*model.spill_binary[t]
  382. model.spill_binary_two_c = Constraint(model.T, rule=spill_binary_two)
  383. # SELF CONS
  384. def self_cons(model, t):
  385. &quot;X &lt;= x2 + My&quot;
  386. return model.self_cons[t] == model.total_gen[t] - model.spill_gen[t]
  387. model.self_cons_c = Constraint(model.T, rule=self_cons)
  388. # ACTUAL CONSTRAINTS
  389. def min_spill(model, t):
  390. &quot;Limit spill renewables to 10% of total&quot;
  391. return sum(model.spill_gen[t] for t in model.T) &lt;= 0.2 * sum(model.total_gen[t] for t in model.T)
  392. model.min_spill_c = Constraint(model.T, rule=min_spill)
  393. def load_match(model, t):
  394. &quot;contract enough renewables to offset 100% load, even if its not time matched&quot;
  395. return sum(model.total_gen[t] for t in model.T) &gt;= sum(model.load_v[t] for t in model.T)
  396. model.load_match_c = Constraint(model.T, rule=load_match)
  397. # **********************
  398. # Define the battery income, expenses, and profit
  399. # **********************
  400. green_income = sum(model.spot_v[t] * model.spill_gen[t] for t in model.T)
  401. black_cost = sum(model.spot_v[t] * (model.load_v[t] - model.self_cons[t]) for t in model.T)
  402. slr_cost = sum(40 * (model.slr1_gen[t] + model.slr2_gen[t] + model.slr3_gen[t]) for t in model.T)
  403. wnd_cost = sum(70 * (model.wnd1_gen[t] + model.wnd2_gen[t] + model.wnd3_gen[t]) for t in model.T)
  404. cost = black_cost + slr_cost + wnd_cost - green_income
  405. model.objective = Objective(expr=cost, sense=minimize)
  406. # Solve the model
  407. # solver = SolverFactory(&#39;glpk&#39;)
  408. solver = SolverFactory(&#39;cbc&#39;)
  409. solver.solve(model)
  410. # , timelimit=10
  411. results_df = model_to_df(model, first_period=first_model_period, last_period=last_model_period)
  412. print(results_df)
  413. results_df.to_csv(&#39;temp.csv&#39;, index=False)
  414. </details>
  415. # 答案1
  416. **得分**: 2
  417. 以下是翻译的内容
  418. 1. 您的模型写得很好清晰而有条理),但存在一些重要的线性规划问题以下是一些建议以帮助您解决这些问题....
  419. 您通过将变量相乘来使此模型非线性化这会给求解器带来巨大的复杂性在这种情况下是不必要的您可以使用一些约束包括一个"big-M"约束以消除您在将标志乘以容量两个变量...因此非线性性部分的非线性性您实际上正在执行以下操作

output[t] = energy[t] * size[t] * running[t]

  1. 其中energy是一个参数,其他的是变量,而`output`则被合并到一个表达式中。您可以通过将output作为一个新的变量,并使用3个约束来线性化它:

output[t] >= energy[t] * size[t] - (1 - running[t]) * M
output[t] <= energy[t] * size[t]
output[t] <= running[t] * M

  1. 这里的M是能量*大小的逻辑上限。有时候上限和下限都是不需要的,但是由于超额约束,您可能需要这个来强制执行您想要的"相等"约束。我建议将`output`设为一个新变量,以减轻您的复杂表达式。
  2. 2. 在尝试创建"合并表达式"时,您实际上正在创建一组"索引表达式"。只需在模型的T上使用求和语句:

model.d_slrgen_var = Expression(model.T, rule=dependent_solar_gen)

  1. 创建了一组表达式,而不是一个求和。
  2. 3. `max()`是一个非线性运算符,您不能在LP表达式中使用它或返回它。因此,引入另一个类似于以下的变量:

spill[t] ∈ {非负实数 }
spill[t] >= ... 对于模型T中的每个t

  1. 当您解决了所有这些问题并使其正常工作后,您可以考虑重新组织模型,并通过`[source, time]`双重索引主要变量,而不是按源命名它们,这样会更清晰一些,但可以留到以后再考虑。如果您遇到困难,请回复评论...
  2. <details>
  3. <summary>英文:</summary>
  4. Your model is well written (clear and organized), but shows some significant Linear Programming issues. A couple pointers to help you along the way....
  5. 1. You are making this model non-linear by multiplying variables together, which adds HUGE complication for the solver and is not necessary in this case. You can use a couple of constraints, including a &quot;big-M&quot; constraint to remove the portion where you are multiplying your flag times capacity (both variables...hence the nonlinearity). You are essentially doing this:
  6. ###
  7. output[t] = energy[t] * size[t] * running[t]
  8. where energy is a parameter, the others are variables, and `output` is rolled up into an expression. You can linearize that by just making output a variable and using 3 constraints:
  9. output[t] &gt;= energy[t] * size[t] - (1 - running[t]) * M
  10. output[t] &lt;= energy[t] * size[t]
  11. output[t] &lt;= running[t] * M
  12. Where M is some logical upper bound on energy*size. Sometimes the upper + the lower are not needed, but you probably need this to enforce the &quot;equality&quot; constraint you want because of the overage constraint.
  13. I would make `output` a new variable to alleviate the complex expressions you have.
  14. 2. In your attempt to make the &quot;rolled up expressions&quot; you are instead making a set of &quot;indexed expressions&quot;. Just use a summation statement over model.T. This:
  15. ###
  16. model.d_slrgen_var = Expression(model.T, rule=dependent_solar_gen)
  17. makes a **set** of expressions, not a summation
  18. 3. `max()` is a nonlinear operator, and you cannot return or use that in a LP expression. So.... Introduce another variable similar to:
  19. ###
  20. spill[t] {non-negative REALS }
  21. spill[t] &gt;= ... for each t in model.T
  22. After you manage to get all that ironed out ( ;) ) and working, you might want to re-organize the model and double-index your main variables by `[source, time]` instead of naming them all by source, which is a bit cleaner, but window dressing for later.
  23. Comment back if you are stuck...
  24. </details>

huangapple
  • 本文由 发表于 2023年2月24日 09:48:12
  • 转载请务必保留本文链接:https://go.coder-hub.com/75551933.html
匿名

发表评论

匿名网友

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

确定