I have a convex programming problem in which I am constrained to several periods, each of these periods represents different times of a day in minutes. Assume we are constrained to 7 periods in the day, these periods consist [480, 360, 120, 180, 90, 120, 90].
Update to my thoughts on this:
Can the 7 intervals variable be transferred to a binary variable of 1440? This would mean we can calculated the level as needed.
I would assume to use these periods as a max for my integer variable which can be defined as X. X is a CVXPY variable X = cp.Variable(7). Performing my solution I create and define the problem by creating constraints, the constraints I want to work with are as follows:
- Target >= min target
- Reservoir level >= min level
- Reservoir level <= max level
I understand than in order to calculate the reservoir levels I must feed correct data to ensure calculation such as surface area, what is expected to leave the reservoir. The problem I am struggling with is that due to the shape of X I feel like I should be ensuring the reservoir isn't overfilling between periods, at the moment my calculation just checks at point 0, point 1 ... point 7 and this does satisfy the constraint, however in real world the issue I am facing is that we are exceeding the max level in between these points at stages and would need to factor this into account, how could we refactor the code below to account for this given the variable running times of pumps as set by X
Please see below the code that we currently are working with.
# Imports
import numpy as np
import cvxpy as cp
# Initial variables required to reproduce the problem
periods_minutes = [480, 360, 120, 180, 90, 120, 90]
energy_costs = np.array([0.19, 0.22, 0.22, 0.22, 0.22, 0.22, 0.19])
flow_out = np.array([84254.39998627, 106037.09985495, 35269.19992447, 47066.40017509, 26121.59963608, 33451.20002747, 20865.5999279])
# Constant variables
pump_flow = 9.6
pump_energy = 9.29
pump_flow_per_minute = pump_flow * 60
# CVXPY section
N = len(periods_minutes)
MIN_RUN_TIME = 10
running_time = cp.Variable(N, integer=True)
mins = np.ones(N) * MIN_RUN_TIME
maxs = np.ones(N) * periods_minutes
k = cp.Variable(N, boolean=True)
# Optimiation calculations
running_time_hours = running_time / 60
cost_of_running = cp.multiply(running_time_hours, energy_costs) * pump_energy
sum_of_energy = cp.sum(cost_of_running)
volume_cp = cp.sum(running_time*pump_flow_per_minute)
period_volume = running_time * pump_flow_per_minute
# Create a variable that will represent 1 if the period is runnning and 0 for the remeainder, 1440 total
# Example running_time[0] = 231, then this is 231 trues in the variable
# test = np.zeros((1, 1440))
# for i in range(N):
# for j in range(running_time[i]):
# test[0][j] = 1
# Reservoir information and calculations
FACTOR = 1/160.6
flow_in = running_time * pump_flow_per_minute
flow_diff = (flow_in - flow_out) / 1000
res_level = cp.cumsum(flow_diff) * FACTOR 2.01
# Constant constraints
min_level_constraint = res_level >= 1.8
max_level_constraint = res_level <= 2.4
volume_constraint = volume_cp >= 353065.5
# Build constraints
constraints = []
# Convert the integer variables to binary variables
# constraints = [test_cp[0] == 1]
# Append common constraints
constraints = [min_level_constraint]
constraints = [max_level_constraint]
constraints = [volume_constraint]
constraints = [running_time >= cp.multiply(k, mins)]
constraints = [running_time <= cp.multiply(k, maxs)]
# Objective definition
objective = cp.Minimize(cp.sum(sum_of_energy))
# Problem declaration
prob = cp.Problem(objective, constraints)
prob.solve(solver=cp.CPLEX, verbose=False)
# Each ith element of the array represents running time in minutes
running_time.value
Note that some variables are part of our external class: Surface Area: 160.6m² Min Level: 1.85m Max Level: 2.4m Pump Flow: 9.6l/s Pump Energy: 9kW
At the moment our outflow data for the reservoir is in 30 minute intervals, Ideally if we could develop a solution to allow for this in the sense of say a inflow matrix that accounted for the various running times over a period of time and accounted the volume such as imagine a variable output for X being [231 100 0 0 30 90 99] If we look at the first element being 231 I would expect something like this in our matrix given 480 minutes for the first element as the max running time for period 1, this yields 16 elements as 480/30.
The expected outcome given this would be something like
[17280 17280 17280 17280 17280 17280 17280 12096 0 0 0 0 0 0 0 0]
Figures shown above are in volumes 17280 being full 30 minute interval running and 12096 being 21 minutes of the period, 0 being not running. I hope to have provided enough information to entice people into looking at this problem and look forward to answering and queries you may have. Thanks for taking the time to read through my post.
CodePudding user response：
Problem
I assume that the pump starts running at the beginning of each time period, and stops after running_time seconds, until the next time period starts. We are checking the level at the end of each period but within the periods the level may get higher when the pump is working. I hope I've understood the problem correctly.
Solution
The constraint is:
res_level(t) < 2.4
The function is piecewise smooth, the pieces being separated by time period boundaries and the event of pump shutdown within each time period.
Mathematically we know that the constraint is satisfied if the value of res_level(t) is smaller than 2.4 at all critical points—i.e. piece boundaries and interior extrema.
We also know that res_level(t) is linear in the piece intervals. So there are no interior extrema—except in case of constant function, but in that case the value is already checked at the boundaries.
So your approach of checking res_level at ends of each period is correct, except that you also need to check the level at the times of pump shutdown.
From simple mathematics:
res_level(t_shutdown) = res_level(period_start) flow_in - (t_shutdown/t_period) * flow_out
In CVXPY this can be implemtend as:
res_level_at_beginning_of_period = res_level - flow_diff * FACTOR
flow_diff_until_pump_shutdown = (flow_in - cp.multiply(flow_out, (running_time / periods_minutes))) / 1000
res_level_at_pump_shutdown = res_level_at_beginning_of_period flow_diff_until_pump_shutdown * FACTOR
max_level_constraint_at_pump_shutdown = res_level_at_pump_shutdown <= 2.4
constraints = [max_level_constraint_at_pump_shutdown]
Running the code with this additional constraint gave me the following res_levels (levels at end of periods):
[2.0448792 1.8831538 2.09393089 1.80086488 1.96100436 1.81727335
2.0101401 ]