Development of a biogas distribution model for a wastewater treatment plant: a mixed integer linear programming approach

In 2015, Northumbrian Water Limited (NWL) was the only wastewater company in the UK to use all sludge after wastewater treatment to produce renewable electricity (Northumbrian Water 2015). Recently it has added a gas-to-grid injection plant at its sewage treatment works (STW) in Howdon, Tyneside, whereby biomethane/biogas produced from anaerobic digesters is purified prior to injection into the National Grid (Martins & Wilson 2016). The site typically digests up to 40,000 tonnes of sewage sludge (dry solids) annually producing renewable biogas. The STW can be considered to have two sections: effluent water treatment, and sludge treatment. Sludge treatment occurs at the advanced anaerobic digestion (AAD) plant. The relatively new processing techniques (CAMBI thermal hydrolysis) involved in the AAD plant results in significant operational challenges to maximise the economic performance.

Wastewater treatment plants (WWTPs) are continuously being driven towards increased efficiency of plant operation due to regulatory pressures (O'Brien et al. 2011; Mulas et al. 2015), and it is widely thought that there is considerable room for improving the efficiency of WWTPs (Longo et al. 2016). If companies are able to schedule their operations successfully around supply and demand of energy, there are opportunities for annual OPEX savings of 2–5% (Merkert et al. 2015), and the use of energy recovery techniques and biogas can reduce a WWTP carbon footprint by 10% (Singh et al. 2016). In addition, WWTPs have an abundance of low-grade thermal energy and organic substances in sludge (Zhang et al. 2017) which, if harnessed, can reduce carbon emissions and energy cost/consumption.

Typically, biogas produced on a wastewater anaerobic digestion (AD) plant in the UK is used to generate electricity only. The AAD plant at Howdon has the ability to use the biomethane in different ways. The AAD plant is rare in that it has three possible uses for biogas produced on site, shown in Figure 1: upgrade for injection into the national gas grid; burning it in combined heat and power (CHP) engines; or burning it in steam boilers for the thermal hydrolysis plant. There are three CHP engines, three steam boilers and one gas upgrade/injection plant on site (Figure 2). If required, the plant may flare excess biogas under emergency circumstances, though operators minimise this as much as possible. To ensure overall sludge processing remains unimpeded the plant may draw natural gas from the national gas grid to be used in the CHP engines or steam boilers. Nationally, many similar sites are looking to upgrade to produce biomethane for injection into the National Gas grid to take advantage of the government's renewable heat incentive (Hale 2017).

CHP electricity generation is an effective way to reduce energy consumption and carbon emissions (Mago & Chamra 2009; Bischi et al. 2014); these systems typically take a fuel source and efficiently convert it into usable heat and energy (Andersen & Lund 2007) and in the case of the Howdon WWTP the CHP units are ‘gas engines’ provided by MWM (Martins & Wilson 2016) that burn the non-purified biogas from the digesters to generate heat and electricity.

According to M. Katebi et al., ‘The barrier to the successful implementation of the control system is not the control software or hardware, but rather the problem of designing control systems that are integrated with the plant operation and management and have a high degree of local autonomy, flexibility and reliability’ (Katebi et al. 1998). Typically, in the water industry, advanced process control is uncommon and standard control schemes use on-off control via PLC (programmable logic controller) and SCADA (supervisory control and data acquisition) systems to control localised processes to within specified limits, without consideration of upstream plant behaviour (O'Brien et al. 2011). Widespread use of SCADA-type system technology permits the exploitation of more advanced supervisory concepts and system control (Ocampo-Martinez 2010).

SCADA based systems can store vast amounts of historical data regarding plant operations. Such large amounts of data can be used for retrospective analysis and learning to make improvements to future operations. An example of the application of retrospective learning is in the aviation industry (Dambier & Hinkelbein 2006), where analysis and learning after incidents creates improved safety procedures. Typically, retrospective learning techniques are not used in regard of operational aspects of a process (Cummings et al. 2017). However, in work conducted by T. Cummings et al. in 2017, retrospective optimisation (RO) and learning was used to develop improved scheduling procedures for multiple nitrogen liquefier units and the development of electricity spot pricing for a pricing predictor (Cummings et al. 2017). In the case of Howdon WWTP, site systems store several years of operational data which is not currently being used for any RO or analysis. To the authors' knowledge, there has been no development of models of WWTPs based on RO in literature, though it should be noted that previous work has been done developing predictors or forecasts of energy pricing using grey prediction models (Gu et al. 2015).

To perform RO for this site, we develop a model of the process and use mixed integer linear programming (MILP) to determine optimal historic performance that can be compared to actual historic operations. MILP takes a series of linear relationships made up of equalities, inequalities, integer defined parameters and upper and lower bound constraints to minimise a linear objective function. Some examples of various MILP applications include: scheduling of a polymerisation reactor (Kelley et al. 2018), multiproduct milk processing (Touil et al. 2016), multi-week scheduling for an ice cream processing facility (Wari & Zhu 2016) and scheduling of a multiple cryogenic air separation unit and compressor plant (Adamson et al. 2017a, 2017b), where the later was performed within a Microsoft Excel spreadsheet. This highlights the effectiveness and ease of application of being able to model a process in a MILP form, if the problem can be posed in such a way to take advantage of MILP software.

There is a clear gap in literature of designing control schemes, optimisation techniques and models for gas distribution of AAD sections of WWTP, as literature-discussed here focuses on developing control strategies or models of the effluent treatment side of a WWTP. It should be noted that whilst there are many studies on improving the yield of biomethane produced from anaerobic digesters, such as using calcium or enzymatic pre-treatment of sewage (Agabo-García et al. 2019; Salama et al. 2020), this paper focuses on the optimal way to use the biomethane produced, not improving yield or developing a model of the anaerobic digesters.

Currently, there is no model on site to advise the optimal distribution of biogas produced or how much natural gas must be purchased to ensure optimal cost and sludge processing performance, nor is there evidence of such a model in use at a WWTP internationally. Here, for the first time, we report the development and application of a realistic model and its solution using MILP, which is used to advise the optimal daily operational strategy that will minimise of gas distribution costs. Using the model, we perform RO using historic data to determine how optimal previous operational strategies were on a cost minimisation basis.

The remainder of this paper is structured as follows. In the methods section we outline the site constraints, formulation of the MILP model and develop the objective function for optimisation. In the results we show example optimisation results from the optimiser, show the results of the RO performed, consider the importance of maintenance and improved breakdown robustness, and consider current model limitations. Finally, our conclusions are drawn and presented alongside proposals for further work and improvements to this model.

To model the site and perform RO, processing limitations of gas flow for each unit must be known. The daily processing limits (Table 1) were determined through retrospective analysis of historical plant operational data and through discussions with operational managers on site. Flow limits are different for biogas and natural gas volumes. The minimum gas flow for an engine to operate is 50% the maximum flow. Currently, should there not be enough biogas produced to satisfy operation of a CHP engine, operators require that an engine run on natural gas and have the gas flow rate to that engine set to the maximum.

The gas holders on site that capture and intermediately store the biogas produced are typically able to store up to an hour's production of biogas, should the biogas flow downstream be interrupted. Therefore, for the purpose of this model the daily volume of biogas produced must all be utilised in the available units on site.

The cost parameters for each gas flow were taken from the OPEX reporting features available on the onsite SCADA system. Costs tend to be in units of p/kWh or equivalent, so the cost parameters in the model are combined with standardised conversion factors to enable a simple cost in £/Nm³ (GBP per normalised cubic metre) of gas; gas flows on site are reported in Nm³. For the purposes of confidentiality, the cost function parameters cannot be reported here.

For a realistic optimisation, both costs and revenues of gas streams on site are considered; injecting biomethane into the National Grid through the BUP creates a revenue as the gas is sold (⁠⁠), burning biogas in the boilers and CHP engines has no associated cost (⁠⁠) whilst purchasing natural gas for use anywhere has a cost associated (⁠⁠).

Should the optimiser advise that flaring gas should be done, it will affect the actual minimised cost; after optimisation, when calculating the total daily operational cost on site, if required the penalty term is subtracted from the minimised cost to reflect actual operational costs, and for more accurate RO.

Under these conditions, when takes the value of 1 then only biogas may flow to engine ‘i’, whereas a value of 0 denotes a natural gas flow. is a variable determined by the model to optimise gas distribution. The binary variable provides the gas selection functionality by allowing to alter the upper and lower constraints of the inequalities.

Retrospective analysis of the different fuel flows to the boilers also revealed no distinguishable difference in processed volume of biogas or natural gas flows, hence each gas type has equal weighting in this constraint.

In an ideal setting there would be no limit to the volume of biomethane that can be injected into the National Grid. However, the total volume that can be injected is subject to local demand and gas network pressures; if too much biomethane is injected too quickly the pressure in the grid could rise too high for continued injection, whereby Northern Gas Networks (the local gas network operator) would shut off grid injection from the site; this is known as going into a ‘reject’ state. Operators have discovered that the maximum daily volume of biogas that can be processed through the BUP is currently around 40,000 Nm³, which was validated through retrospective analysis of BUP processing volumes, grid injection volumes and ‘reject’ state instances.

The optimiser model was developed using MATLAB using function intlinprog, which is a MILP algorithm in the Optimisation Toolbox package. The optimiser model was implemented using this function and contains 24 constraints, 14 possible gas flows with 6 binary variables. The optimal gas flows and binary variable values are obtained through minimisation of the cost function to give optimal gas distribution on an economic basis, whilst maintaining site operability.

The total historic daily biogas volume produced on site is the only historic value passed to the optimiser. The model then calculates the volumes of biogas and natural gas that are required to be distributed to each unit to satisfy daily operation at the minimum cost, subject to the constraints outlined.

The function intlinprog will always aim to minimise the objective function. As such, all revenue parameters are negative and cost parameters are positive inside the optimiser; a negative optimal value provided by the optimiser indicates a potential revenue, whilst a positive value indicates a cost. After optimisation, the resulting ‘optimal value’ is multiplied by −1 so that a positive value represents revenue and negative a cost, as one would expect. When displaying results of optimisations visually, total daily revenue is used for a direct comparison to site operations.

To perform RO, historic plant operational data was used for the 12-month period. For a given date, the historic total daily production of biogas on site is passed to the optimiser which then provides operators with an optimised minimum cost and the optimal daily strategy for operating the plant for that date. The optimal strategy provided can then be compared to actual gas distribution on site, using historic data. For each day in the 12-month data set, the optimised minimum cost for RO was compared to historical site operation (and the subsequent cost of operating using that strategy) to compare how the site has performed in the past.

During RO, the model constraints are adjusted to reflect planned site maintenance; for example, if a CHP engine was offline on a specific date for annual safety inspections, then the upper and lower gas flow constraints for biogas and natural gas for that engine would be set to 0 Nm³. This allows the model to perform optimisation with the correctly available units on site.

For a given daily volume of biogas produced on site, the optimiser provides the operator with a visual operation strategy for the optimised minimum cost; Figure 3 (left graph) shows the results to maximise cost reductions for a daily biogas production of 40,000 Nm³ (a typical production level for the site). The rightmost graph in Figure 3 displays a version, provided by the optimiser, showing the percentage daily utilisation of each unit.

For a daily biogas production of 40,000 Nm³, the site should be operated according to the strategy in Figure 3 for optimal cost efficiency: to inject all biomethane into the National Grid, use natural gas in the CHP engines at 100% load to generate electricity on site and use natural gas in the steam boilers to create steam as required. For a given daily volume of biogas, the optimiser provides a fast and reliable result in a matter of seconds.

Performing RO of the plant over multiple days is difficult to achieve if the operators must sift through multiple graphs for each day that is analysed. For an easier visual indication of how the optimiser suggests the plant should be operated with regards to fuel types for each unit, a Gantt chart showing the gas type selected for is presented – Figure 4 shows an example date range from the full RO and compares this to historical operation.

The colour of each block represents a particular fuel type for each unit, typically biomethane or natural gas. However, in reality there may be times where both fuels were recorded historically, as shown in yellow. When switching an engine from biomethane to natural gas, historical operations will show both fuel types were used but for clarity each fuel would have been combusted separately, not as a blend. For example, Figure 4 suggests both fuel types were used historically in the CHP engines on some days, but for CHP engines this actually show days where a gas type switchover took place.

The historical fuel type shown in the boilers is inferred from historical gas flow data, and on occasion the flow meters indicate that two fuel types were used when this may not have been the case, hence the multiple ‘yellow’ days for boilers. Specifically, in the case of Boiler 1 it may be an indication of an issue with the gas flow meters on site. It may also be an indication that fuel switchovers in the boilers happen more often.

In this figure, it may be seen that the optimiser advised that, between days 6 and 11, the plant should have swapped some CHP engines and boilers to run on biomethane, which the plant mostly managed to achieve with the exception of the boilers.

Each daily optimisation has an associated minimised operational cost, which is used to evaluate RO of plant operations. The minimum cost found by the optimiser is compared to the cost associated with actual historical distribution on site by plotting one against the other. The daily optimal operational cost was compared to the actual historical cost over the 12-month period November 2017 to October 2018, the historical plant data available, and is shown in Figure 5.

In order that RO is representative of site operations, planned site maintenance is included in the RO. Using site maintenance logs, if a unit was taken offline for a full day for scheduled repairs or safety inspections then the optimiser will not have that unit available for the given day of RO.

It should be noted that, whilst the parameters of the model are static, for each 3-month moving period, the average maximum gas flows to each unit (biogas and natural gas) are recorded and compared to the gas flow parameters of the model. If the model parameters are more than 10% away from the 3-month average, a warning is displayed during RO such that the user can determine and select an improved parameter. As of yet, the warning parameters have not been triggered, indicating that the initial parameters are valid across the whole RO horizon. It is also why there is a 10% band on the RO graphs.

A band of 10% of the maximum revenue achieved was chosen as the tolerance for optimality as the model parameters do not deviate more than 10% from the historic use of the units on site. Of the 366 days analysed, based on a tolerance of ±10% of the highest revenue, RO indicates the plant operated to within tolerance 98% of the time (332 days). This provides some validation to operators that their historic operational strategy was optimal. However, the 2% non-optimal operational days have a potential lost revenue of almost £350,000 (£349,291); if the site was optimally operated 100% of the time, there is the potential for an increase in revenue of 7.8%. The non-optimal operational days are likely due to unplanned breakdowns and maintenance on site, which can be investigated by comparison with site maintenance logs.

It is paramount that planned shutdowns are logged and programmed into the optimiser to allow for the best comparison of RO and historic operation. However, unplanned shutdowns cannot be pre-configured in the model. RO will highlight non-optimal days and, typically, will highlight the severity (if any) of unplanned shutdowns across the site and will produce a monetary value attributed to each shutdown. However, should an operator wish to change the number of available units (e.g. CHP engines) on site for scenario modelling, they are able to do so.

Figures 6 and 7 demonstrate the importance of including planned shutdowns in the model; by not including them, the model would assume all units were available and thus it would appear the plant operated sub-optimally yet this was not the case historically. It should be noted here that the model has no preference of which engine should be in operation, hence on the example date shown in Figure 7 the optimiser selects a different engine used for biogas – the actual engine used is not important to the model, and is an operator decision to maintain the health of the engines.

Using the site maintenance log, planned and unplanned downtime of units on site can be shown on the annual RO graph in Figure 5. The boiler, CHP engine and the BUP unplanned outages are plotted on Figure 8. Based on the distribution of outages on Figure 8, boiler downtime has no major impact on the historic operability of the site. The CHP engines would appear to have more of an impact on site operations as more of the outliers are associated with unplanned outages.

However, almost every non-optimal outlier on Figure 5 is associated with an unplanned outage (full or part day) with the BUP shown on Figure 8. Gas injection into the grid is where the plant makes the majority of gas revenue, and the fluctuations of daily gas injection volumes account for most of the variation in the plant operation; whilst a maximum biogas processing volume of 40,000 Nm³ is typical, this limit is influenced by external factors and on occasion the site may be able to inject more biomethane (and subsequently process more biogas) than the model accounts. It also stands to reason that unplanned outages with this unit will likely cause the plant to operate sub-optimally.

Crucially, full day planned maintenance outages do not cause the site to operate sub-optimally – the plant operates within optimal tolerance on each of the Full Day planned maintenance plots on Figure 9. Currently, the model can only handle planned maintenance where a unit is due to be offline for the full day – planned maintenance taking less than a full day cannot accounted for in the current model. The potential shortfall in revenue of ∼£350 k is mostly due to due to unplanned breakdowns and planned maintenance that does not take a full day.

In an ideal world, process plants would always run optimally without any unexpected breakdowns or unforeseen circumstances preventing optimal operation. In reality, this is not the case and should be noted when considering the potential 7.8% increase in revenue available. If the plant operated to within optimal tolerance 98% of the time, then to what extent could this be practically increased further? Whilst undesirable, it is expected that process plants will ultimately have unexpected downtimes, especially with new technologies and installations. Acceptable levels of lost revenue due to unexpected maintenance is a decision for site managers; the RO presented here is supporting investment decision making – allowing operational managers to quantify the benefit of maintenance expenditure

After investigating the possible cause of outliers from RO (Figure 5), it is clear the main focus for improvement should be the BUP. Each of these outliers corresponds to a day where the BUP failed (full day or part day, Figure 8) or whether the BUP was planned to be offline for maintenance for part of the day (Figure 9). Failures in the BUP could be mechanical (such as a component breaking) or ‘gas-spec’ failures; if the gas composition does not meet the correct standard for injection it will enter a ‘failed’ state and reject gas until the specification target is met.

It is worth noting that the BUP part of the AAD plant is a recent addition to the site, with installation completing late 2015. It is a new operation for operators and for the local gas network to manage correctly; the data used for RO is still within the BUP's early years on site, so a higher number of unexpected downtime would be anticipated due to the complexity and age of the technology. It is believed by operators that the BUP has operated more robustly in recent years, but RO has not been performed over a later date yet so this cannot be confirmed in this paper.

The highlighted data on Figure 8 is also a good indicator of the number of unexpected instances on site, their origin and how impactful breakdowns of individual units are with regards to economic operations. As the number of planned part-day maintenance instances is low, RO of the plant using this model is a good indication of the potential improvement for plant robustness and could help provide justification for future investments.

Currently, the Biogas Distribution Model is a retrospective optimiser only, as a future prediction of biogas production is currently unavailable. Work will be carried out modelling the digesters to predict this biogas production, to provide operators with a full sludge processing optimisation model to allow for improved operational forecasting.

In addition, the model does not consider electrical costs associated with the CHP engines, electrical import tariffs or site load with regards to electricity and heat and assumes that CHP engines must be in operation provided they are not offline for maintenance. Electricity generated through the CHP engines is used on site to reduce electrical consumption from the National Grid and heat produced from them is reused to reduce gas requirements for the steam boilers on site (the steam is used as part of the thermal hydrolysis process). Clearly, shutting down a CHP engine will increase the electrical demand on the National Grid, with site electrical costs increasing due to reduced CHP electrical output, as well as increasing the heat load on the steam boilers.

To incorporate electricity generation, it will require splitting the model into smaller optimisation steps. The model currently uses daily values for processing limits, however electrical import/export costs on site are ‘fixed-variable’ (vary half hourly but are fixed and known) meaning that each half hour the CHP electricity generated on site will offset a different import cost, thus will require its own optimisation step. This would present unique and interesting challenges, as each optimisation step would require linking together, to prevent excessive ramping (turning on and off) of CHP engines, or fuel type switching every half hour – fuel types for CHP engines cannot be switched over throughout the day. Such improvements to the model could be developed further to include rotational scheduling to improve the longevity of equipment on site.

This paper proposes a MILP model for gas distribution and optimisation of an AAD plant with multiple options available for biogas use on site. The optimisation model takes a single input of daily biogas volume produced, in m³, by the anaerobic digesters to provide operators with an optimal daily operational strategy; the strategy is provided to the operators in visual form for each day, though weekly strategies are also possible if required. The model also includes a penalty term on the flare stack for preventing unnecessary use of the emergency flare.

RO techniques were performed with historical operational data to determine the economic effectiveness of previous operational strategies. The model proposed uses MILP optimisation techniques to perform RO across the historical data horizon (one year) in a matter of seconds, with individual daily optimisations taking far less computational time.

Comparing the summated daily gas revenues used for RO of the plant indicates a potential increase in revenue of 7.8% (∼£0.35 m) for the period optimised if the plant was run optimally all of the time. Upon discussion with operators and analysis of the maintenance and breakdown logs, the ‘lost revenue’ and non-optimal performing days (points on Figure 5 not within tolerance band) are mostly due to unit malfunctions, indicated on Figure 8, and are mostly due to the loss in revenue from lower biomethane injection due to BUP unplanned outages (Figure 8, graphs ‘e’ and ‘f’).

The potential increase in annual revenue from gas distribution could be achieved through improved operational strategies or (more importantly) improved plant robustness, such as reducing the number of incidents of unplanned outages, with a specific focus on the biomethane injection plant. The potential increase in revenue could be used to justify future plant improvements.

Finally, improving the proposed model to include electrical costs and heat loads, as discussed previously, would transform the model into a full energy management model for the AAD plant: a possible ‘digital twin’, or a decision support system (DSS), that could be used for RO and to validate future operational strategies.

This work has been collaboratively funded by Newcastle University and Northumbrian Water Ltd (under grant EPN5095281) and also the Engineering and Physical Science Research council (EPSRC) as part of a grant of Centre for Doctoral Training in Engineering for the Water Sector (STREAM, EP/L015412/1). The authors wish to express thanks to the staff at Northumbrian Water (with special thanks to Andrew Moore, Luke Dennis and Ken Black), who have been supportive through this ongoing work and without whom this project would not be possible.

Data cannot be made publicly available; readers should contact the corresponding author for details.