Preface to the “2024 v1 edition”

I’m writing this preface in March 2024, having originally wrote the bulk of this document in the COVID lockdown of 2020, inspired by several conversations at the time about South Africa’s electricity troubles. At the time I had built some basic capacity expansion models for SA to understand for myself what future investments should be made for the SA grid. My own modelling results supporting the (much more sophisticated modelling) results from the CSIR and others. I built these models using OSeMOSYS (the “software,” which is explained later).

Back then (and many times since) I found myself in tense discussions with friends over the value of modelling and the accuracy thereof in the real world. I found a glaring lack of appreciation for what the capacity expansion modelling (CEM) is good at, and what it is not good at. At the heart of this was a lack of understanding of energy modelling in general. I felt that a “beginners guide” to CEM was a gap that I would try to fill. This document is the result of that attempt.

I made the first, very draft version of this doc in 2020, shared with some friends and found it to be mildly useful. There were many gaps, but it was helpful in some ways. I never put the document online because of its very rough nature. After that first version, I made some half-hearted occasional attempts to update the document, and so ended up with a marignally improved “2022 version” and then even a “2023 version.” But these were effectively light touch updates.

Now, in 2024, I find myself with a renewed urge to update this doc (for many reasons, which I will discuss elsewhere). The doc is still extremely rough, but I am now of the view that I can put it in the world in its rough state, with various disclaimers and disclosures and warnings. My new intention is to update this document “in public” because I now think that even in rough format, it will be of more use to me in the public domain than remaining on some select whatapp groups.

For this early 2024 version, I’ve made some light tough edits, to assist with flow mostly.

I’m not sure how far I will get with this update though. Time will tell.

Without further rambling, here are the warnings:

Warnings!

• This document is very rough, and still very much a “work-in-progress” (WIP)
• Therfore there are comments in this doc which are directed at myself and only I understand
• There are also probably errors - this is not (yet) a peer-reviewed piece of work.
• The bulk of this was written several years ago in 2020, so some numbers might feel odd. This is OK, the numbers don’t need to be accurate and are rather illustrative.
• In general, read with an open mind and send me comments please.

Version control

• This is the “2024 v1 edition”
• I intend to use this section to link to previous versions as I update this doc

Introduction

When someone says South Africa needs more coal, or more solar, or more nuclear, how do they know this? Well, they’re usually referring to a “model” someone built somewhere11 In SA, most people making these statments are “picking up stompies” as we sometimes say, and don’t have a good view of the overall picture.

“Capacity Expansion Models” (CEMs) are the type of model that gives the best answers to what we need to build in the future to meet electricity demand under certain constraints and objectives. CEMs are basically large mathematical optimisation problems that solve for the lowest possible “system cost” under certain constraints.

Do you remember “linear programming” from high school maths? Yeah, it’s essentially that. Just bigger. And more complex. Way bigger and way more complex.

As I discussed in my preface (really, go read my preface above), my intention here to try to explain CEMs to a beginner. I assume that my audience22 I should do a proper “assumed audience” piece like Maggie describes here has some basic mathematical knowledge, and for now some basic knowledge of electricity (I will try to strip this away over time).

The case study is ridiculously simple and idealistic. The point is to get an understanding, not to derive useful insights at this point.

Okay, here we go:

Overall Process

2024 edit: I’ve left the section below almost untouched. The original idea was to give an overview of how this “exploration explainer” flows. In time, I will make this section into actual full sentences and give it some meat.

• Uses the style of “introduction to operations research textbooks,” because this style has had success (anectodally in explaining linear programming)
• Case study that solves “by inspection”33 insert wikipedia definition: “If an expression can be evaluated by straightforward application of simple techniques and without recourse to extended calculation or general theory, then it can be evaluated by inspection.”
• Then we’ll solve by basic linear optimisation
• Then we’ll solve using OSeMOSYS via MoManI
• Then we’ll explore the OSeMOSYS model to understand it better.

An introduction to Linear Programming

See the MIT book. Rework an example. Probably the production example. http://web.mit.edu/15.053/www/AppliedMathematicalProgramming.pdf

This book is out of print (as far as I can tell, and this is verified by an MIT OCW course here). This MIT site also notes this book as “available free online” so based on all of this, I think it’s safe to say this material can be used by us for this tutorial.

The book itself was published in 1977 - here’s the Google Book’s Link.

Need to choose an example with two constraint sets, at least. (so that it matches later with the energyburg model)

Things to know before exploring the case study:

• Difference between Power and energy (GW vs GWh). (Capacity/size vs production/generation/output)
• Basics of the electricity system, i.e. generation, transmission, distribution (i.e. there are “plants” which make electricity and “wires” which move the electricity to consumers)
• Capital costs (cost of building a plant/asset, once off cost), fixed operating costs (cost of operating an asset/plant per year, that does not change no matter you use the plant), variable operating costs(cost of operating an asset/plant per year, based on how much you use that plant), etc.
• Some other stuff I’ll add as I go along

The “Energyburg” Case Study

Energyburg is a fictional island named after my home town of Johannesburg. Here is the context of the case study:

Energyburg was recently subject to a storm that destroyed all it’s power generation, and now it needs to rebuild. The transmission and distribution network are however still in place.

The island has a consistent electricity demand of 1 GW. This never changes (so in energy terms the annual demand is 8760 GWh_el)

Additionally, there are no seasonal variations on the island (weather is the same all year round). Every day consists of 12 hours of daylight and 12 hours of night. It is assumed that the sun rises instantly and provides the same output for every moment of the 12 hours. We have idealised the sun as a switch here, to make modelling easier.

Due to various reasons, the island has only two options for power generation. Diesel generators and solar PV. Both of these are assumed to be able to be built instantaneously.

Furthermore, everyone from the island will be leaving in exactly one year from tomorrow, so the energy planners have to only provide power for one more year. We will magically build everything we need today, so tomorrow Energyburg will have electricity, and needs to have electricity for 8760 hours (~1 year) from tomorrow.

The costs associated with the two options are as follows. The currency used is the South African Rand (R), a random choice to be more consistent with the naming convention of the location.

The capital costs (or the costs of building the plants instantaneously, sometimes called the “overnight” costs) of the diesel generator is R1 Million for a 1 GW plant. The capital cost for the PV plant is R2 million for a 1 GW plant.

The fixed operating costs are R100,000-00 per GW for the diesel plant, and R10,000-00 for per GW of the PV plant.

The variable operating costs: R10,000/GWh for the diesel and R200/GWh for PV.

Energyburg wants to minimise the cost of meeting demand for this year.

Other assumptions:

• Copperplate model (ignoring transmission and distribution)
• Efficiencies - we do not account for any efficiency losses in the generating plant

Solving by inspection

A page of mathematical jargon on wikipedia says that “by inspection” means something can be solved through simple analyse without extended analysis. This problem is solvable by inspection and deliberately so - we are able to test our model against the solved analysis.

We want the lowest total cost, or more precisely, the lowest system cost. (IV to expand on what this means)

We have two technologies available to Energyburg, so there are three options:

• Option 1: Diesel Only
• Option 2: Diesel and PV
• Option 3: PV Only

Since we have two time periods, day and night, and because the sun does not shine at night, option 3 (the PV only option) can be easily eliminated. PV alone cannot meet demand.

So let’s start with option 1 and work out what it would cost. If we used only diesel, the capital cost would be R1 million (we would build a 1 GW diesel plant to meet the constant demand). The fixed operating cost would be R100,000-00. To get the variable cost, we know that we annual demand (in energy terms) based on a constant power demand of 1 GW is 8760 GWh_el. Multiplying 8760 by the per unit of energy cost for diesel (R10,000-00) gives us R87.6m for the variable costs. Adding up the three cost components yields a total cost to supply Energyburg with electricity for the year of R88.7mn.

For option 2, here we build both plants, the PV plant and the diesel plant. Since the sun does not shine at night, we will need a 1 GW diesel plant to meet the 1 GW constant demand at night. For the diesel plant, we have capex cost of R1m, and a fixed operating cost of R100k. However, for variable costs, only the night time demand is satisfied by the diesel plant. The night time demand is 4380 GWh_el (due to the perfect 12 hour split). So multiplying the 4380 GWh_el by the per unit of energy cost for diesel (R10,000-00/ GWh_el) gives us R43.8m for variable costs. Adding them all up gives R44.9mn for the diesel plant.

Now, for the PV plant, the capex is R2mn, the fixed operating cost is R10k. The variable cost is worked out in a simliar fashion to the diesel plant: 4380 multiplied by R200, which gives R876,000-00. Summed up, we have a total for solar of R2.886m.

Adding the solar and diesel costs together, we have a grand total of R47.786m for option 2.

Now to compare the two options. Option 2 turns outs to be cheaper than Option by R40.914m! (Option 1 = R88.7 and Option 2 = R47.786m).

So from this it’s clear that Option 2 should be pursued. That is, both a PV and diesel plant should be build. And the PV plant will supply day time electricity while diesel will supply the night time demand. This will supply electrcity to Energyburg for a total system cost of R47.786m.

Solving by linear optimisation

(Note: most of this remains unchanged from the original 2020 version)

Before we start problem solving, let’s generate a table of the data we know:

(The tables above will need to be converted from an image to an actual table at some point)

Using the approach from Intro to Linear Programming, we know that we need to select the decision variables.

In this case, we let the decision variables be the rate at which each plant operates. I.e. the electricity (i.e. energy) produced by each plant in each time period. We have two time periods here: day and night. And we have two technologies. So we have:

\(x_{diesel,day}\) : electricity produced by the diesel plant during the day (in energy/“day period44 While the units of this variable are in GWh/day, this refers to all the ”day time” in the year. So not daily, but rather 12hrs *365 effectively, which is 4380 hours. If that makes sense”)
\(x_{diesel,night}\) : electricity produced by the diesel plant during the night
\(x_{PV,day}\) : electricity produced by the PV plant during the day
\(x_{PV,night}\) : electricity produced by the PV plant during the night

We also need to know the size/capacity of what we need to build,so we add two more decision variables:
\(x_{DieselSize}\): the size/capacity of the diesel plant to be built
\(x_{PVSize}\): the size/capacity of the solar plant to be built

Now the second step is to determine the objective function of this problem. We want to minimise costs. We know that:

\(Costs = Capital Costs + Fixed Operating Costs + Variable Operating Costs\)

But now, using the table of data costs, capital costs are equal to:

\(Capital\;Costs = Diesel\;Capital\;Costs\; +\; PV\;Capital\;Costs\) \(Capital\;Costs = x_{DieselSize}*CapCost_{diesel}+x_{PVSize}*CapCost_{PV}\)

\(Fixed Operating Costs = x_{DieselSize}*FixOpCost_{diesel}+x_{PVSize}*FixOpCost_{PV}\)

\(Variable Operating Costs = ((x_{Diesel,day}+x_{Diesel,night})*VarOpCost_{diesel})+((x_{PV,day}+x_{PV,night})*VarOpCost_{PV})\)

Putting it all together, we get an objective function as follows:

\(Min\;Costs = x_{DieselSize}*CapCost_{diesel}+x_{PVSize}*CapCost_{PV} + x_{DieselSize}*FixOpCost_{diesel}+x_{PVSize}*FixOpCost_{PV} + ((x_{Diesel,day}+x_{Diesel,night})*VarOpCost_{diesel})+((x_{PV,day}+x_{PV,night})*VarOpCost_{PV})\)

The third step is to generate the constraints.

The first set of constraints is related to ensuring that the demand in each time period (or timeslice) is met. To do that, we firstly need to calculate this demand.

So let \(demand_{night}\) = 4380 (8760 hours in a year, *0.5 since night is half the time). It’s the same with \(demand_{day}\). Units are energy, so GWh in our case.

Now for the constraints:

\(C1day:x_{diesel,day}+x_{PV,day}\geq demand_{day}\) \(C2night:x_{diesel,night}+x_{PV,night}\geq demand_{night}\)

So the energy produced by diesel and PV during the day summed together must equal or exceed the energy demanded during the day. Same for the night.

The second set of constraints is related to capacity. We need to ensure that the size/capacity of the plant that’s built by the model (\(x_{DieselSize}\) and \(x_{PVSize}\)) is sufficient in order to produce the energy required in each timeslice.

To do this, we firstly define a parameter we’ll call the “CapacityRatio”. For each technology, this represents the amount of energy that would be generated by that technology if one unit of that technology ran all year with no stop. So for diesel plant, if 1 GW (1 unit) of the plant ran all year, you would get 8760 GWh from that plant. So \(CapacityRatio_{diesel}\) = 8760. Same with the PV plant.

\(C3:x_{DieselSize}*CapacityRatio_{diesel}*day_{ratio} \geq x_{diesel,day}\) \(C4:x_{DieselSize}*CapacityRatio_{diesel}*night_{ratio} \geq x_{diesel,night}\) \(C5:x_{PVSize}*CapacityRatio_{PV}*day_{ratio} \geq x_{PV,day}\) \(C6:x_{PVSize}*CapacityRatio_{PV}*night_{ratio} \geq x_{PV,night}\)

So the LHS on the above inequality is basically the energy that could be produced by the diesel plant during the day (\(day_{ratio}=0.5\)), for a plant size \(x_{DieselSize}\). This (the LHS) has to be bigger or equal to than the energy produced by the diesel generator during the day (the RHS). This forces the model to appropriately size \(x_{DieselSize}\), while reducing it as much as possible because it appears in the objective function.

Put differently, the LHS of the above constraints basically calculates the energy (in this case electricity) that that specific size of diesel plant could generate in the entire year (by multiplying the plant size by it’s capacity ratio), and then scales it down for that specific timeslice (the scalling down is done by day_ratio and night_ratio, which are both 0.5 since they both last 0.5 times the year - ie 6 months each based on 12 hours ago each a day.)

2024 v1 edit: put even another way, this set of equations makes sure that the plant size (e.g. \(x_{DieselSize}\)) is big enough to produce the energy required of that plant in a particular time slice (e.g \(x_{diesel,day}\))

We have to do the same constraints for the diesel gen at night, and then for PV during the day and night. We don’t actually need constraint C6 above, but I left it in there because it’s simpler (doesn’t hurt either).

Now that we have the objective function and all the constraints, we can pick some software to solve this linear programming problem.

In my view the easiest way to solve this small LP is through spreadsheets. It’s the most easily accesible way (imo). To that end, I built two version - one in Excel and another in google sheets. The Google sheets one can be viewed here and illustrates how the problem is solved. I followed a tutorial that is linked on the sheet itself.

The Excel version is slightly more clearly setout (but not by much), and for now remains on my desktop55 5 IV Projects\6 Year 2020\1 OSeMOSYS 2020\7 Basic Excel LP.

In both models, once solving the above LP we see a lowest system cost of R47.786m, which matches perfectly with our “solved by inspection” version. So we have confidence that our LP is correct.

Solving using OSeMOSYS and MoManI

Conclusion to Part 1

While I didn’t exactly set out to have a part 1 and a part 2, it seems now that we’re at the conclusion of a natural part 1.

What we’ve seen is how a simple toy problem (Energyburg, with all it’s assumptions) can be solved by inspection, modelled as a Linear Progamme (LP) and solved in a spreadsheet, and modelled as an LP and solved in big boy software OSeMOSYS.

If you’ve spent enough time on this document, it should, I think, have been possible to go from beginner to at least understanding the basics of how CEM is done.

Of course, real problems have many, many more parts. But the principles and the approach stay the same.

(IV to expand on this)

Part 2: What’s happening underneath the CEM software bonnet?

So, we put all the data into MoManI and it sent the data to OSeMOSYS (on the back end) which gave us a nice LP which we solved using CPLEX. All very well. But what actually happened under the bonnet here?

Well, to answer this, we can use the original OSeMOSYS paper and the OSeMOSYS documentation (PDF) (online link). A careful reading of these docs, particularly the OSeMOSYS docs gives a good understanding of what’s happening in the equations of a typical model.

But reading these docs alone don’t give us a view of what’s happening under the bonnet in the case of Energyburg. To do this, we need to look at the underlying model. I’ve extracted the underlying model and put into this spreadsheet here (the process to do this will be covered in an Appendix, at some point). If you open the spreadsheet, you will find the model described by a matrix where all the equations are rows and the the columns represent all the variables (IV: note to self, we should make a matrix of our model above so we can easily relate the two. #FutureWork)

However, the matrix model is quiet large with many, many zeros (despite Energyburg being a very simple model, when I built the OSeMOSYS model via MoManI, I used all equations so we can get an appreciation of some of the underlying complexity). To make it easier to navigate the matrix model, I built some simple R scripts to help us view constraints (rows) and variables (columns) without the zeros (details of this will be (eventually) provided in an Appendix #FutureWork).

So between the docs and the tools, we can start exploring the model.

The most obvious place to start is the objective function, which we learn from the papers is the minimisation of “TotalDiscountedCost,” which in turn is the sum of TotalDiscountedCostByTechnology and TotalDiscountedStorageCost. We can see this from equation prefixed as “TDC2” in the OSeMOSYS docs (in the original paper there is a slightly different structure, with it being TDC1). We can look for “TDC2” in the spreadsheet (search) and we find equation 194, the very last equation.

Using the Rscript tool, we can pull out that equation and we have:

From here, we can see that TotalDiscountedCost is a summation of the the total discounted costs of each of the technologies, i.e. Diesel Gen, Solar PV and Diesel Import.

We can now “follow” one of these variables in the model. If we follow the diesel gen costs, we need to extract column x133 from the model:

Here we can see that x133 is linked to equation 191. So looking at equation 191:

We can see that TotalDiscountedCostByTechnology for Diesel Gen is the summation of Discounted Capital Investment of Diesel Gen (x22), Discounted Operating Costs (x25) and Discounted Salvage (x28) and Emission Penalties (x31). This makes sense. In our energyburg model we only have capital costs and ops costs, so x28 and x31 can be ignored for our purposes.

Let’s have a deeper look at x22 (X….22.DiscountedCapitalInvestment.Energyburg.Diesel_Gen.2021.):

So we see that X22 also shows up in equation 114, which has prefix CC2. Looking at that equation:

So this equation is “doing the discounting” but in our Energyburg model we have ignored discount rates. So this is a negligible piece for us (for this purpose).

Let’s have a look at variable x13.

x13: X….13.CapitalInvestment.Energyburg.Diesel_Gen.2021.:

We are now out of discounted land (finally, lol). Let’s have a look at Eq35:

Ah, here is something interesting: Remember how a GW of new diesel capacity costs R1mn to build? Here is it showing up in the model now. So every new GW of capacity for diesel gen (x35) adds 1e6 (R1mn) to the capital investment.

More importantly, this is effectively the link between costs (capital investment) and capacity (GW)

BTW, cross referencing this agains the OSeMOSYS doc is helpful here: on page 13 (17/47), under variables, we see that “NewCapacity[r,t,y]” is defined as “Newly installed capacity of technology t in year y” and importantly, is expressed in power units (which in our model is GW)

Let’s have a look at x35: X….35.NewCapacity.Energyburg.Diesel_Gen.2021.:

We can see that X35 is linked to eq29, which has an a prefix “CAa1.” Looking this up in the OSeMOSYS docs, we see that is part of the Capacity Adequacy A set of equations:

Let’s have a look at Eq29:29 CAa1_TotalNewCapacity[Energyburg,Diesel_Gen,2021]

Because we only have one year in this analysis, this equation is simply setting all the accumulated capacity equal to the new capacity in 2021 (the only year of our analysis)

Let’s have a look at X1: X…..1.AccumulatedNewCapacity.Energyburg.Diesel_Gen.2021.:

This takes us into “CAa2” equation, which is eq32: 32 CAa2_TotalAnnualCapacity[Energyburg,Diesel_Gen,2021]

So here we can see the model is linking new accumulated capacity with total annual capacity (x128). For completeness sake, “AccumulatedNewCapacity” is defined as “Cumulative newly installed capacity of technology t from the beginning of the time domain to year y” and has units of power. “TotalCapacityAnnual” is defined as “Total existing capacity of technology t in year y (sum of cumulative newly installed and pre-existing capacity)”

Remember that because we only have one year here in Energyburg, these all essentially collapse to the same thing in our case.

Let’s have a look at x128: X…128.TotalCapacityAnnual.Energyburg.Diesel_Gen.2021.

So here we clearly have a very interlinked variable, which makes sense. The Total Capacity of Diesel is important!

It’s worth pausing here and recapping our exploration: we started at the objective function of minimising costs, and have traced a path through the model which took us to capital cost and therefore new capacity (we went down the diesel path).

So using the OSeMOSYS docs, we can see that equation 105 (and 106), prefixed with CAa4, links capacity with activity. So let’s have a look at that equation.

And here we see how exactly capacity is linked to activity. RateOfActivity is defined as an intermediate variable which “represents the activity of technology t in one mode of operation and in time slice l, if the latter lasted the whole year” and has units of energy. Remember capacity has units of power. So this is basically ensuring (together with other equations) that the capacity is big enough to deliver the energy required.

However, here in Eq105, it is actually RateOfTotalActivity (x98) that is used as a variable. and not RateOfActivity itself. However, we should find the link between the two by following variable x98, and I expect to see it show up in an equation prefixed with “CAa3” (from the OSeMOSYS doc). Let’s check:

Sure enough, eq99 shows up. Let’s follow that.

And there we have the link between RateOfTotalActivity (x98) and RateOfActivity (x66). Now, the difference is basically due to different modes of operation, which we don’t have here in Energyburg, so these values will be the same.

RateOfActivity (x66) is an interesting variable. Let’s have a look.

Here we can see how it is linked to a bunch of equations, including equations prefixed with “EB.” According to the docs, this stands for “Energy Balance” and these equations “Ensures that demand for each commodity is met in each TimeSlice.” Aha.

So here is the link, presumably, between activity and demand. Let’s have a look at EBa1, which is equation 42:

Yep there it is. Recall that RateOfActivity is “the activity of technology t in one mode of operation and in time slice l, if the latter lasted the whole year.” Now RateOfProductionByTechnologyByMode is “the quantity of fuel f that technology t would produce in one mode of operation and in time slice l, if the latter lasted the whole year”

So given that we have only one mode of operation, variable x92 is the electricity produced by the diesel generator in the “day” timeslice.

Let’s have a closer look at x92:X….92.RateOfProductionByTechnologyByMode.Energyburg.AllDay.Diesel_Gen.1.Electricity.2021.

It links to another Energy Balance equation in equation 50, prefixed EBa2, so let’s have a look at that: 50 EBa2_RateOfFuelProduction2[Energyburg,Diesel_Gen,Electricity,2021,AllDay]

More mode related stuff, not revelant to us (its a sum over all modes, but we have only 1 mode). So let’s see what variable x81 is linked to: X….81.RateOfProductionByTechnology.Energyburg.AllDay.Diesel_Gen.Electricity.2021.

Continuing to follow the Energy Balance “path”, we can explore equation 128, prefixed with “EBa3”: 128 EBa3_RateOfFuelProduction3[Energyburg,AllDay,Electricity,2021]

So this is interesting and cool. This is essentially a summation of all the rate of production of electricity across all the technologies during the day timeslice.

Let’s follow x77:

And then lets look at eq136, prefixed “EBa7” :

This equation accounts for the year split (see OSeMOSYS docs), so basically in this case, it adjusts the RateOfProduction for the day timeslice by the ratio of time the day timeslice actually takes (0.5) Let’s follow x42:

And now looking at eq124, prefix “EBa11”:

And here at last we find the Energy Balance link between production, demand and usage! This is the big link.

This equation is essentially ensuring that the production of electricity in the day time slice is greater or equal to the sum it’s demand and usage.

“Use” or more accurately “UseByTech” is “Use of fuel f by technology t in time slice l.” In our case, electricity is not used by any tech, so this is going to be 0 in this equation. From the OSeMOSYS docs, we can see that there are several equations which produce this value (see EBa4, EBa5, EBa6 etc.).

Demand is also generated by other equations (see EQ_SpecifiedDemand etc.)

So that is essentially the end of the path we’ve traced. From costs, to capital costs, to capacity, to linking capacity to production, and production to demand.

The same process can be traced for op costs for diesel, and for all costs of solar. The same pattern plays itself out. So by now, we should have been able to see what’s happening underneath the bonnet.

(Now to just explain it more thoroughly).