Calliope math helper function for clipping values

[joeseakinsarup]

What can be improved?

Context:
I can't work out how to force, for example, heat generation technologies to act as either baseload or peaking without using cost based incentives. In my scenario I want a gas boiler to only be used if a heat pump is working at max capacity and there is still demand to be met, and I want to keep the gas boiler as cheaper to run than the heat pump.
My attempts so far have revolved around trying to set the gas boiler flow_out <= demand - heat pump max cap. This would work if it didn't introduce, when the demand < heat pump max out, a (less than) negative constraint for the gas boiler's flow out (which I believe is making the model infeasible as this conflicts with other min 0 constraint). I (crudely) checked my theory by setting the constraints to 0 if they were negative using the model.backend and then successfully solving the model.
I was also unable to introduce logic to only apply a 0 boiler flow out limit where demand < heat pump max cap.

Possible solution:
I think one option would be to able to apply the xarray.DataArray.where such that you can apply a default of 0 (i.e.a.where(a >0, 0)). Similarly, being able to create some kind of mask (e.g. demand <= tech max flow out) in a global expression or sub_expression that can be used in a where.

Alternatively, if this is already possible I would really appreciate an example in the docs for this kind of logic!

Version

v0.7.0.dev6

[brynpickering]

We've just expanded what you can define as math in the where string (#834), which will be available in the next release. With that, if you are OK linking boiler output to the max possible heat pump capacity (or flow_cap is just a parameter), you could do:

constraints:
  boiler_operation_limit:
    foreach: [nodes, timesteps, carriers]
    where: sink_use_equals[techs=<demand-tech>] < flow_cap_max[techs=<heat-pump-tech>]
    equations:
      - expression: flow_out[techs=<boiler-tech>] == 0

What you can also do already is give an additional cost to the boiler on a different cost class, which you remove later:

techs:
  <boiler-tech>:
    cost_flow_out: # or cost_flow_in
      data: [<boiler-cost>, <heat-pump-cost - boiler-cost + 1>] # make total impact on the objective slightly higher for the boiler
      index: [monetary, dummy]
      dims: costs
parameters: 
  objective_cost_weights:
      data: [1, 1]
      index: [monetary, dummy]
      dims: costs

This will penalise the boiler operation but then when you come to calculate total costs, you can just ignore the dummy dimension.

[joeseakinsarup]

We've just expanded what you can define as math in the where string (#834), which will be available in the next release. With that, if you are OK linking boiler output to the max possible heat pump capacity (or flow_cap is just a parameter), you could do:

constraints:
boiler_operation_limit:
foreach: [nodes, timesteps, carriers]
where: sink_use_equals[techs=] < flow_cap_max[techs=]
equations:
- expression: flow_out[techs=] == 0
What you can also do already is give an additional cost to the boiler on a different cost class, which you remove later:

techs:
:
cost_flow_out: # or cost_flow_in
data: [, <heat-pump-cost - boiler-cost + 1>] # make total impact on the objective slightly higher for the boiler
index: [monetary, dummy]
dims: costs
parameters:
objective_cost_weights:
data: [1, 1]
index: [monetary, dummy]
dims: costs
This will penalise the boiler operation but then when you come to calculate total costs, you can just ignore the dummy dimension.

Ah amazing, thank you for the response, #834 sounds like it would be the best solution. For now I am using a workaround of limiting gas availability to the node based on demand - max heat pump cap.

The dummy cost suggestion makes sense, however, I left out from my initial message that in our model there are a number of nodes with various combinations of allowed technologies (such as thermal storage) and pre-defined fixed (or flexible) capacities, and we are optimising the usage of these technologies via the monetary costs. I'd be worried that introducing the extra cost would therefore have unexpected consequences. For example, in a node with a heat demand, a heat pump, some thermal storage, and a gas boiler, we want the storage to be sized primarily to shift heat pump usage to utilise cheaper electricity, however, if we inflate the boiler running costs (through the dummy class) that could mean that storage might instead be oversized to reduce the boiler usage.

[brynpickering]

In such a circumstance you could also consider giving the heat pump a buffer storage (include_storage: true). Technically it is a buffer on the input (electricity) but if you move the COP to flow_in_eff then it would be effectively the same as representing a heat storage (since the energy flowing in has already been converted to heat). This would mean only the heat pump could use the storage. Or do you still want the option for the boiler to feed into the storage, but only if the heat pump is already being fully utilised?

[joeseakinsarup]

Interesting, I had forgot about that option! However, wouldn't other heat generation technologies still be able to indirectly utilise this storage by freeing up the heat pump to store heat, and the boiler being artificially more expensive could then still force the heat pump buffer to become oversized (as the heat pump could run well in advance of when gas would be needed to store heat so the gas doesn't need to be used)? We would also want non-heat pump techs to contribute to the storage directly (i.e. at their cost rather than at the heat pumps electricity cost in that timestep) as we require storage to meet peaks (which could be prolonged and coincide with when electricity is very costly).

[brynpickering]

Fair enough. If you have heat pump capacity fixed then you could pre-compute the boiler output maximum and then add a timeseries constraint like you're already doing by limiting gas supply, but a bit more explicit:

constraints:
  limit_boiler_output:
    foreach: [nodes, timesteps]
    where: lookup_zero_boiler_output
    equations:
      - where: lookup_zero_boiler_output
         expression: flow_out[techs=boiler, carriers=heat] == 0
      - where: NOT lookup_zero_boiler_output
         expression: flow_out[techs=boiler, carriers=heat] <= flow_out[techs=ashp, carriers=heat] - flow_in[techs=demand_heat, carriers=heat]

lookup_zero_boiler_output would be a binary array where it's True whenever demand is less than the heat pump capacity. The second constraint ensures the boiler is still limited to meet demand rather than purely flow into the storage. It could even be made more strict:

flow_out[techs=boiler, carriers=heat] <= (flow_out[techs=ashp, carriers=heat] + flow_out[techs=heat_storage, carriers=heat]) - flow_in[techs=demand_heat, carriers=heat]

This still allows the model to size the storage even if there's a fixed flow_cap for the ashp. It obviously breaks down if you want the ashp flow_cap to be a decision variable. At that point, you'd need a binary variable:

variables:
  can_operate:
    foreach: [nodes, timesteps]
    domain: integer
    bounds:
      min: 0
      max: 1
constraints:
  set_can_operate:
    foreach: [nodes, timesteps]
    equations:
      - expression: flow_out[techs=ashp, carriers=heat] >= flow_cap[techs=ashp, carriers=heat] - BigM * ( 1 - can_operate)
  limit_boiler_output:
    foreach: [nodes, timesteps]
    equations:
      - expression: flow_out[techs=boiler, carriers=heat] <= flow_cap_max[techs=boiler, carriers=heat] * can_operate

Note

These are all untested. They might have parsing errors in them.

[joeseakinsarup]

The lookup_zero_boiler_output timeseries makes sense, and am I right in thinking the main benefit to this over loading in a timeseries of the max capacities would be it is a lighter file and you can encode the logic in the math rather than a separate pre-processing script/step? For this, should the expression be:

expression: flow_out[techs=boiler, carriers=heat] <= flow_in[techs=demand_heat, carriers=heat] - flow_out[techs=ashp, carriers=heat]
rather than
expression: flow_out[techs=boiler, carriers=heat] <= flow_out[techs=ashp, carriers=heat] - flow_in[techs=demand_heat, carriers=heat]
as we want a positive difference between demand and the ashp output (as in this case we know the askp<demand)?

On the binary variable addition, the below is my semi-raw thought process on trying to understand it, and so I think I answered some of my own questions later in the text. I've left it semi-raw in case working through the logic as I did helps others (who like me haven't used binary variables before, and don't fully understand them). Saying that, I'd really appreciate any corrections to my logic if I don't actually answer my own questions correctly!

I don't understand why the set_can_operate constraint will make the can_operate equal 1 when the ashp is operating at full capacity: wouldn't can_operate == 0 satisfy all instances of the constraint? I can see that can_operate == 1 would work when the ashp is working at full capacity, but I don't see why it choose to be 1 instead of 0 (if both are valid)?

Ah, is it because then the limit_boiler_output constraint will force the can_operate variable to be 1 in cases when the boiler has to be used to meet the peak (or whenever it is optimal to use it)? Then, as the can_operate forces the ashp to be working at full capacity (in set_can_operate) the constraint will only 'want to' turn this can_operate on when there is enough (>ashp_max_cap) demand, otherwise you would be generating more heat then you can deal with (if no heat storage or export)?

If I'm correct with that, does that then also mean that the boiler can't operate in times when you wanting to partially shift the heat pump usage by using storage (as the ashp will not be working at max cap)? For example, if a peak == ashp + boiler, you wouldn't be able to also meet it by doing peak == 0.5ashp + storage + boiler as can_operate, and, therefore, the boiler cap, will be 0. I guess this might not be a problem if the ashp isn't a fixed size as you then would be able to have 0.5ashp as your chosen capacity, and this is how you then are optimising between the size of the storage and the heat pump.

In my scenario, if the heat pump is a fixed capacity would I still be able to use a binary variable and construct constraints to allow for the boiler to work along side a partially loaded heat pump (as long as the demand > heat pump max capacity)? Or does no longer forcing max heat pump usage meant that there wouldn't be the required pressure to not turn the boiler on whenever it wants? I'm thinking you would almost want to construct a fake heat pump that would implement your binary logic, and then not constrain the real heat pump...

[brynpickering]

The binary variable is a standard trick in linear programming to introduce a binary switch into a problem, to represent an if/else statement. In this case, for there to be any ASHP outflow, can_operate must equal 0. if can_operate = 1, then the first constraint becomes flow_out >= flow_cap, which effectively says flow_out == flow_cap. If can_operate = 0 then the first equation becomes flow_out >= flow_cap - bigM. BigM is a very large number so it's kinda like saying flow_out >= -inf, i.e., no constraint. bigM could also be flow_cap_max[techs=ashp] for the same effect.

The second constraint then follows from that, if can_operate = 0 then the constraint becomes flow_out <= 0, i.e. no operation allowed. If can_operate = 1 then flow_out <= flow_cap_max, so it's effectively not a binding constraint. In this case, the boiler can now operate and, simultaneously, the ASHP has to be running on full output.

You could update this to limit the storage device as well. For instance, not allowing the storage device to be charged when the boiler is operated:

flow_in[techs=heat_storage, carrier=heat] <= flow_cap_max[techs=heat_storage, carrier=heat] * (1 - can_operate)

It's quite a powerful tool but one that does add complexity to the problem that you'll likely need a commercial solver to get to a solution in a reasonable amount of time.

[joeseakinsarup]

Thanks Bryn, I also just noticed the link to Modelling Energy Systems shared in the docs - I've got some learning to do!