A solar plant is supplying power to three water pumps, and other appliances and I want to find when should I run my water pumps to make the most of the solar plant (since its production through days is not constant).

I tried to implement such an optimization problem on Pyomo, but he provides me 0 solution, which is quite surprising since it doesn’t appear to be a very complex problem.

My modelisation follows these rules :

• power_m_pyomo contains the average yield power of the solar panel per month (indexed by m) and per hour between 8 and 17 (indexed by Id),
• base_demand contains the electricity demand of the appliances, which is equals to 90 kW every hour between 8 and 17,
• model.P contains the number of water pumps that should be running each hour. So it’s an integer smaller or equal to 3. I only need 9 “equivalent hour of water pump running per day”, which means that the sum of model.P[i] should be equal to 9.

My code is :

# Importation of the supplied power
df = pd.read_csv('test_power.csv', header=0, names=["Hour", "Power [kW]"])
power_supplied = df['Power [kW]'].tolist()

# Creation of the average monthly yield power 
power_m = monthly_profile_creation(power_supplied)

# Conversion of the power supplied into something that can be used by Pyomo
power_m_pyomo = {}
for m in range(1,13):
    for h in range(8, 18):
        power_m_pyomo[(m,h)] = power_m[m][h]
    

# Creation of the model (now empty)
model = pyomo.ConcreteModel()

# Creation of index for the peak demand
model.M = pyomo.Set(initialize = range(1, 13))
model.Id = pyomo.Set(initialize = range(8, 18))

# Parameters 
base_demand = list_to_dico([90]*10, 8)
model.base_demand = pyomo.Param(model.Id, initialize=base_demand)
model.power_supplied = pyomo.Param(model.M, model.Id, initialize=power_m_pyomo)

# Variables
model.P = pyomo.Var(model.Id, domain=pyomo.Integers)

# Objective function

def objective_func(model):
    power_stored = 0
    for m in model.M:
        for i in model.Id:
            power_stored += model.power_supplied[m,i] - model.base_demand[i] - 24*model.P[i]
    return power_stored
model.objective = pyomo.Objective(rule = objective_func, sense=pyomo.minimize)

# Constraints

def constraint_pump_1(model, i):
    return model.P[i] <= 3
model.constraint_pump_1 = pyomo.Constraint(model.Id, rule=constraint_pump_1)

def constraint_pump_2(model):
    return sum([model.P[i] for i in model.Id]) == 9
model.constraint_pump_2 = pyomo.Constraint(rule=constraint_pump_2)

solver = pyomo.SolverFactory("cbc")
results = solver.solve(model, keepfiles=False, logfile="solve.log")
print(results)

Can someone please help me ?