1

I have a problem in finding the soulution for the status of my solver due to the fact that I want to find an optimal solution for my LP. Below you can find the code for the problem I want to solve. I have no idea why the optimizer is not able to find an optimal solution. Since I am a beginner in programming I deeply appreciate any help. Thank you in advance!!

`

import pyomo.environ as pe
from Timer import Timer
from Pvsystem import Pvsystem
from Loadprofile import Loadprofile

#model object
LPopt = pe.ConcreteModel()

#objects:
my_timer = Timer()
my_pv = Pvsystem()
my_load = Loadprofile()

#define 'index sets'
LPopt.T = pe.RangeSet(0, 8759, 1)

#declare pv_power as a variable
LPopt.pv_power = pe.Var(LPopt.T, domain=pe.NonNegativeReals) 
#define constraint list for pv output
LPopt.pv_constraint = pe.ConstraintList()
pv_value = (my_pv.get_total_pv_output()) #in kW
#add constraints to the pv constraint list
for t in range(my_timer.get_total_time_steps()):
    LPopt.pv_constraint.add(float(pv_value[t]) == LPopt.pv_power[t])

#declare load_power as a variable
LPopt.load_power = pe.Var(LPopt.T, domain=pe.NonNegativeReals)
#define constraint list for load
LPopt.load_constraint = pe.ConstraintList()
load_value = (my_load.get_load_profile()) #in kW

#add constraints to the load constraint list
for t in range(my_timer.get_total_time_steps()):
    LPopt.load_constraint.add(float(load_value[t]) == LPopt.load_power[t])

#parameters
#maximum battery capacity of the storage device
LPopt.max_bat_capacity = pe.Param(default=750, within=pe.NonNegativeReals) #in kWh

#battery charge power efficiency
LPopt.eta_bat_ch = pe.Param(default=0.975, within=pe.NonNegativeReals) #unitless

#limited grid import power
LPopt.max_grid_import_power = pe.Param(default=4500.0, within=pe.NonNegativeReals)
#limited grid export power
LPopt.max_grid_export_power = pe.Param(default=1000.0, within=pe.NonNegativeReals)

#variables
#battery input power as a variable indexed by the set T
LPopt.bat_ch_power = pe.Var(LPopt.T, domain=pe.NonNegativeReals, initialize=0.0)
#battery charge power as a variable indexed by the set T
LPopt.bat_dch_power = pe.Var(LPopt.T, domain=pe.NonNegativeReals, initialize=0.0)
#battery state energy as a variable indexed by the set T
LPopt.bat_energy = pe.Var(LPopt.T, domain=pe.NonNegativeReals, initialize=0.0)
#previous battery content
LPopt.bat_energy_prev = pe.Var(LPopt.T, domain=pe.NonNegativeReals, initialize=0.0)

#grid import power as a variable indexed by the set T
LPopt.grid_import_power = pe.Var(LPopt.T, domain=pe.NonNegativeReals, initialize=0.0)
#grid export power as a variable indexed by the set T
LPopt.grid_export_power = pe.Var(LPopt.T, domain=pe.NonNegativeReals, initialize=0.0)

#create the objective function
LPopt.value = pe.Objective(expr=sum(LPopt.grid_import_power[t] + LPopt.grid_export_power[t] for t in LPopt.T), sense=pe.minimize)

#constraints
#grid constraint
LPopt.grid_import_constraint = pe.ConstraintList()
for t in range(my_timer.get_total_time_steps()):
    LPopt.grid_import_constraint.add(LPopt.grid_import_power[t] <= LPopt.max_grid_import_power)

LPopt.grid_export_constraint = pe.ConstraintList()
for t in range(my_timer.get_total_time_steps()):
    LPopt.grid_export_constraint.add(LPopt.grid_export_power[t] <= LPopt.max_grid_export_power)

LPopt.grid_constraint = pe.ConstraintList()
for t in range(my_timer.get_total_time_steps()):
    LPopt.grid_constraint.add(LPopt.grid_import_power[t] - LPopt.grid_export_power[t] == LPopt.load_power[t] + LPopt.bat_ch_power[t] - LPopt.bat_dch_power[t] - LPopt.pv_power[t])

#battery energy constraint
LPopt.bat_energy_constraint = pe.ConstraintList()
for t in range(my_timer.get_total_time_steps()):
    if t == 0:
        LPopt.bat_energy_constraint.add(LPopt.bat_energy_prev[t] == 0.0)
    else:
        LPopt.bat_energy_constraint.add(LPopt.bat_energy_prev[t] == LPopt.bat_energy[t-1])
    LPopt.bat_energy_constraint.add(LPopt.bat_energy[t] == LPopt.bat_energy_prev[t] + my_timer.get_deltaT() * (LPopt.eta_bat_ch * LPopt.bat_ch_power[t] - LPopt.bat_dch_power[t]))

#battery capacity constraint        
LPopt.battery_capacity_constraint = pe.ConstraintList()
for t in range(my_timer.get_total_time_steps()):
    LPopt.battery_capacity_constraint.add(LPopt.bat_energy[t] <= LPopt.max_bat_capacity)

# Optimize
solver = pe.SolverFactory('cbc')
status = solver.solve(LPopt)

# Print the status of the solved LP
print("Status = %s" % status.solver.termination_condition)

`

I guess that one of my constraints/ or the objective itself might be the problem. Hope to get some help from you!

POLDI1912
  • 21
  • 4
  • 2
    The objective is not the cause of the problem being infeasible. Some solvers can give good messages if the infeasibility is detected in the presolver. Some solvers also have IIS capabilities: report the smallest set of constraints that make things infeasible. Otherwise you can try to relax some constraints in the hope that reveal some useful info. Finally, a fully elastic formulation can help. – Erwin Kalvelagen Jul 05 '18 at 10:13
  • Thanks for your help! Really appreciate it. I relaxed one of my constraints and then it worked out fine! – POLDI1912 Jul 06 '18 at 08:59
  • Check this: https://stackoverflow.com/questions/51044262/finding-out-reason-of-pyomo-model-infeasibility – oakca Aug 04 '18 at 15:10

0 Answers0