# The Waste-to-Biomethane Logistic Problem:
A mathematical optimization approach

Víctor Blanco†, Yolanda Hinojosa‡ and Victor Zavala⋆
†Institute of Mathematics, Universidad de Granada
‡Institute of Mathematics, Universidad de Sevilla
⋆ Wisconsin Institute for Discovery, University of Wisconsin-Madison

###### Abstract.

In this paper, we propose a new mathematical optimization approach to make
decisions on the optimal design of a logistic system to produce biogas from
waste. We provide an integrated model that allows decision makers to optimally
determine the locations of different types of plants and pipelines involved in
the logistic process, as well as the most efficient distribution of the
products (from waste to biomethane) along the supply chain. We analyze the
mathematical model and reduce its size, being able to solve realistic
instances in reasonable CPU times. The results of our computational
experiments, both in synthetic and in a case study instance, prove the
validity of our proposal in practical applications.

###### Key words and phrases:

Logistics; Green Energy; Facility Location; Mathematical Optimization; Biogas;
Supply Chain

##  1\. Introduction

Decarbonization of energy has emerged as a critical global priority in recent
years due to the urgent need to combat climate change. Most World
Organizations have recognized the significance of transitioning to cleaner and
more sustainable energy sources in the next few years to reduce greenhouse gas
emissions. Different agreements have been signed to promote this change.
Specifically, the United Nations, through its Sustainable Development Goals,
has set different targets to ensure affordable, reliable, sustainable, and
modern energy for all by 2030 (United Nations, 2015). The International Energy
Agency (IEA) has been actively promoting decarbonization efforts by providing
policy recommendations, conducting research, and facilitating international
cooperation (International Energy Agency, 2021). Additionally, the
Intergovernmental Panel on Climate Change (IPCC) (Legg, 2021) has been
instrumental in assessing the impacts of climate change and highlighting the
importance of decarbonizing the energy sector. In 2019, the European
Commission presented the European Green Deal (Fetting, 2020), an initiative
that aimed to be a statement of intent on the road-map that Europe should
follow for the implementation of the United Nations’ Sustainable Development
Agendas for 2030 and 2050, designed to mitigate the effects of climate change.
Various measures have been taken by these organizations, including advocating
for renewable energy investments, promoting energy efficiency, encouraging the
use of electric vehicles, and supporting the development of innovative
technologies such as carbon capture and storage. These collective efforts by
world organizations are crucial in driving the global transition towards a
low-carbon future and mitigating the adverse effects of climate change.

One of challenges is to decarbonize the energy system by developing a new
power sector based on renewable sources. It is a known fact that one of the
main strategies to achieve the proposed goal is the use of biogas as an
alternative renewable energy source to carbon-based energies, since it
contributes the reduction of greenhouse gases but also to the development of
the circular economy through the anaerobic digestion of organic waste from
different sources and its transformation into fuel. Since biomethane is the
same molecule as natural gas, it can be distributed via the existing gas
distribution networks, facilitating the transition from natural gas to biogas
energy. Thus, the installation of anaerobic digestion plants for the
conversion of organic waste and livestock manure into biomethane has gained
prominence in recent years as a sustainable waste-to-energy solution and many
countries are making special interest in the assessment of the installation
and adaptation of these biogas plants. For instance, the U.S. Farm Bill
(National Sustainable Agriculture Coalition, 2023), periodically reauthorized
by Congress, includes provisions that support the development and utilization
of biogas from agricultural sources. The bill provides funding for biogas
projects, such as the installation of anaerobic digesters plants on farms or
near them, which capture methane emissions from manure and convert them into
usable biogas for electricity generation or transportation fuel.

Analyzing the logistic systems behind the implementation of different modes of
energy, particularly biogas, is of paramount importance in ensuring the
successful and sustainable integration of renewable energy sources into our
global energy landscape. Biogas, derived from organic waste materials through
anaerobic digestion, presents a promising alternative to conventional fossil
fuels. However, its widespread adoption hinges on addressing intricate
logistic challenges. Understanding the logistics involved in the collection,
transportation, and processing of organic feedstocks is crucial for optimizing
efficiency and minimizing environmental impact. Rigorous analyses on the
logistical aspects of biogas production in different countries reveal
opportunities for streamlining supply chains, enhancing resource utilization,
and reducing overall costs (see, e.g. Muradin and Foltynowicz, 2018; Kwaśny
and Balcerzak, 2017). Effective logistic planning also ensures reliable and
consistent feedstock availability, a key factor in the stable operation of
biogas plants. By delving into the logistics of biogas implementation, we can
develop strategies that not only bolster the economic viability of this
renewable energy source but also contribute to the broader goal of achieving a
more sustainable and resilient energy infrastructure.

Among the decisions to make in the design of an efficient logistic systems,
one of the most critical ones is selecting the optimal sites to locate the
different types of plants involved in the system (biogas plants, pre-treatment
plants, and liquefaction plants). This is a difficult task since it involves
different agents, different production and conversion technologies, different
types of demand centers, and it has to be coordinated with the different
distribution processes (Jafaru et al., 2018). At this point, mathematical
optimization plays a very important role since it is the main tool to
determine solutions to complex logistics system like the one we are dealing
with here. This fact has caused the publication of a large amount of
literature related to this topic. For instance, Mathematical optimization has
been used in Tampio et al. (2022) for the design of a cost-optimal processing
route for a biogas anaerobic digestion plant to produce fertilizer products
based on specified regional needs. In Balaman and Selim (2014) a mixed integer
linear programming model is developed to determine the appropriate locations
for the biogas plants and biomass storages. Scarlat et al. (2018) provides a
spatial analysis algorithm that uses data of manure production and collection
to assess the spatial distribution of the biogas potential in Europe, in order
to decide the location of bioenergy plants. A geographic information system
(GIS) based analysis is used in Valenti et al. (2018) to determine the size
and location of four biogas plants for the Catania province in Sicily. A
multi-objective optimization approach for the optimal location of biogas and
biofertilizer plants at the maximum proﬁt and the minimum environmental impact
is presented in Díaz-Trujillo and Nápoles-Rivera (2019) and it is applied for
a geographical region in Mexico. Mathematical methods have been used to study
the location of bioenergy plants in many other regions and countries (Park et
al. (2019) in North Dakota, Sultana and Kumar (2012) in the province of
Alberta, Egieya et al. (2020) in Slovenia, Amigun and von Blottnitz (2010) in
Africa, Silva et al. (2017) in Portugal, Soha et al. (2021) in Hungarian,
Igliński et al. (2012) in Poland, Delzeit and Kellner (2013) in Germany).
Furthermore, mathematical optimization tools have been already recognized as a
crucial tool to make decisions in the design of systems for biogas production.
Hernandez and Martin (2017), propose a non linear optimization model to
compute the optimal operating conditions of the reactors, the biogas
composition and the content of nutrients in the digestate in the production of
biodiesel from waste via anaerobic digestion. Hu et al. (2018) give a
mathematical optimization framework to analyze the economical and
environmental benefits of recovering biogas, caproic acid (C6), and caprylic
acid (C8) from different sources of organic waste. Ankathi et al. (2021)
consider a mixed integer linear optimization model to determine the location
and capacities of biogas plants based only in the location and production of
the farms. In terms of locational decisions, some of the classical
mathematical approaches have been already adapted to incorporate green goals
(see e.g. Martínez and Fransoo, 2017), although more advances are expected
within the next few year based on what our society needs.

On the other hand, mathematical optimization is recognized as a fundamental
tool for the design and modeling of multiple logistics, transportation, and
supply chain problems (see e.g. Hinojosa et al., 2000, 2008; Blanco et al.,
2010, 2020, 2022, 2023; Blanco and Martínez-Antón, 2024, among many others).
Specifically, is particularly useful to construct robust and effective supply
chains to integrate bioenergy in any economy (see e.g. De Meyer et al., 2015;
Sarker et al., 2018; Tominac et al., 2022; Cambero and Sowlati, 2014; Ghaderi
et al., 2016, and the references therein) . In particular, in order to
implement an efficient and sustainable biogas distribution system, it is
necessary to properly design a robust logistic plan for all the elements
involved in the process, such as manure, biomethane, waste, etc. In this
phase, it is crucial to decide not only where to locate the anaerobic
digestion plants but also, where to locate pre-treatment plants for cleaning
and drying the waste or transshipment plants to distribute the products, how
to collect the manures from the farms or fields and sending them to the
plants, how to link (if needed) the different types of plants, how to
distribute the final biomethane to the gas distribution network or to external
clients, how to dispatch possible fertilizers back to some of the farms, etc.
A few works have already proposed models to make decision in this setting.
Jensen et al. (2017) propose minimum cost flow based model for finding the
optimal production and investment plan for a biogas supply chain by means of
minimizing the transportation cost on an existing supply chain network. Three
layers of analysis for designing optimal animal waste supply for anaerobic
bio-digestion are detailed in Mayerle and Neiva de Figueiredo (2016): (1)
Identification of the optimal locations for the anaerobic digestion plants the
farms’ information; (2) specification of the optimal distribution system; and
(3) an operational layer that includes scheduling optimal biomass collection
from each farm to minimize biogas loss. Sarker et al. (2019) studies the
optimization of the supply chain cost for a biomethane gas production system
which is organized in four stages (collecting feedstock to hubs located
according to zip code areas, transporting feedstock from hubs to reactor(s),
transporting biomethane gas from reactor(s) to condenser(s), and shipping the
liqueﬁed biomethane gas from condensers to final demand points). The problem
is formulated as a mixed-integer non-linear mathematical optimization problem
and, due to the complexity of the non-lineal model, a genetic heuristic
algorithm is propose to solve it.

In this paper, we introduce the Waste-to-Biomethane Logistic Problem and
provide a general and flexible mathematical optimization model to make
decisions of locating biogas, pre-treatment, and liquefaction plants and
distributing the different types of elements along a complex network that
integrates all the aspects mentioned above. Thus, we consider a supply chain
for the biomethane gas production combining six stages: (A) collecting waste
and sending it to pre-treatment plants (independently of the zip code); (B)
delivering contaminant-free organic waste from pre-treatment plants to
anaeorobic digestion plants; (C) constructing pipelines for transporting
biomethane from anaerobic digestion plants to injection points of an existing
gas pipeline network; (D) constructing pipelines for transporting biomethane
from anaerobic digestion plants to biomethane liquefaction plants; (E)
dispatching fertilizer from anaerobic digestion plants to waste sources and;
(F) shipping liquefied biomethane from biomethane liquefaction plants to
customer demand points. We jointly integrate all these complex stages into a
mixed integer linear programming (MILP) model seeking to minimize the overall
transportation cost of the system by assumming that a limited budget is given
to install all the plants and pipelines. Our MILP model can be efficiently
solved using off-the-shelf optimization solvers (as Gurobi, CPLEX, or FICO)
after proving some theoretical results for reducing the number of variables.
This model is applied and tested in a real-word dataset based on the region of
upper Yahara Watershed in the state of Wisconsin (see Sampat et al., 2019,
2021, for further information about this dataset).

The remainder of this paper is organized as follows. In Section 2 we introduce
the Waste-to-Biomethane Logistic Problem (W2BLP) and present our modeling
assumptions. Section 3 is devoted to detail the mathematical optimization
model that we develop for the problem. We also prove in this section a
theoretical result that allows us to significantly reduce the number of
variables in the model. In Section 4, we report the results of our
computational experiments on realistic synthetic instances. A case study is
presented in Section 5 where the model is applied to the real-world dataset
based on the region of upper Yahara Watershed. The paper closes in Section 6
with some conclusions and future research.

##  2\. The Waste-to-Biomethane Logistic Problem

The problem under analysis consists on efficiently use the agricultural waste,
as livestock manure, energy crops, municipal waste, etc, to be transformed
into biomethane. This transformation requires different phases. First, from a
given set of waste sources (WS), as farms or residual storages, each of them
producing an amount of waste, it is transported into pre-treatment (PT) plants
where the non-organic material is removed, and what remains is dried, pressed
and prepared. The location of the PT plants is to be determined among a finite
set of potential positions.

The dried contaminant-free organic waste (DOW) obtained after the pre-
treatment, is delivered to an anaerobic digestion (AD) plant, where it is
transformed into biomethane. The position of the AD plants is also to be
decided among a finite set of potential locations. Biomethane is then either
directly injected into an existing gas pipeline network (GPN) or processed
into liquified natural gas (LNG). We assume that we are given a finite set of
injection points in the GPN and that pipelines connecting the AD plants with a
selected set of injections point are to be built. A minimum percent of the
produced biomethane is to be injected in the GPN, and the remainder biomethane
is transformed into LNG and distributed to a given set of extra customers
(EC), each of them with a given demand. This minimum percent represents the
trade-off between the gas injection to the network with respect to the amount
of gas served to the external customers. The liquefaction of biomethane
requires the construction of Biomethane Liquefaction (BL) plants as well as
pipelines connecting the AD plants with them. Then, the liquified gas is
delivered to the customers using tank trucks. Finally, the remaining material
in the DOW to biomethane transformation process (digester solids) is then
returned back from the AD plants to some of the WS where it can be used as
fertilizer. In Figure 1 we illustrate the different elements that appear in
this logistic problem. The names are differentiated by color. In blue color we
indicate the input data (WS, GPN and EC), whereas in green color we highlight
the decisions to be made on: the (PT, AD and BL) plants and the (to GPN and
BL) pipelines to be installed. The transportation routes are highlighted with
red arrows in the plot.

Figure 1. The process modeled in the W2BLP.

The following hyphoteses are assumed in our model, based on the technological
requirements of the process:

  * •

A minimum percentage of the total amount of waste in the WS must be collected
and sent to the PT plants to satisfy the gas requirements of the GPN. Since a
given percent of the produced biomethane is to be injected in the GPN, this
hypothesis implies that a minimum amount of biogas must be produced to be
injected to the network.

  * •

Each PT plant sends its whole produced DOW to the AD plants. In the AD plants,
the DOW received is transformed into biomethane. The digester solids of this
transformation can be delivered to the WS as fertilizers.

  * •

We assume that each WS receives a given proportion of the total amount of
digester solid produced at the AD plants. In case this proportion is $0$, the
WS does not receive any digester solid (as in the case of residual storages
that may not be interested in the product). We do not force delivering the
whole digester solid produced in the AD plants to the WS.

  * •

The whole biomethane production obtained in the AD plants must be sent. It is
shared into the GPN and the BL plants but, as already mentioned, a minimum
percentage of the total biomethane production is required to be injected into
the GPN.

  * •

Each BL plant sends its whole LNG production to the EC.

  * •

The demand of the EC is not assumed to be fully satisfied. Instead, we
consider that every EC receives a pre-specified proportion of the total LNG
produced at the BL plants. This proportion may depends on its demand or other
preference criterion that the agents consider reasonable.

  * •

A percent of the amount (weight and volume) of the different products is
assumed to be lost at each phase of the whole process. Specifically, a given
percentage of the amount of waste, DOW, and biomethane received in PT, AD and
BL, respectively, is lost during the process. This technical hyphotesis is
derived from the chemical processes to obtain the different types of products
from their raw materials.

The efficient management of this system requires both deciding the location of
the different types of plants (PT, AD, and BL) and pipelines (connecting AD
with GPN and AD with BL) to be constructed, as well as the distribution of the
product along the different phases of this process.

There is a cost for opening each potential location of each different type of
plants, and a construction cost for unit length for the pipelines. We assume
that a budget is given to install all the plants and the pipelines.

For the distribution of the products we assume that a unit transportation cost
is provided for the different phases, namely, delivering from WS to the PT
plants (waste), from PT to the AD plants (DOW), from the BL plants to the EC
(LNG), and from AD to the WS (digester solid).

The Waste-to-Biomethane Logistic Problem (W2BLP) consists of deciding the
location of the different plants and pipelines, and the way the product is
distributed, minimizing the overall transportation cost with the given
installation budget.

##  3\. Mathematical Optimization Model

In this section we present the mathematical optimization model that we propose
for the W2BLP. First, we define the parameters and variables that we use in
our model.

###  3.1. Parameters

The following parameters are assumed to be given to our model.

####  3.1.1. Index Sets

$-$ $N=\\{1,\ldots,n\\}$: Index set for the set WS.
$-$ $P_{1}=\\{1,\ldots,m_{1}\\}$: Index set for the set of potential PT
plants.
$-$ $P_{2}=\\{1,\ldots,m_{2}\\}$: Index set for the set of potential AD
plants.
$-$ $P_{3}=\\{1,\ldots,m_{3}\\}$: Index set for the set of potential BL
plants.
$-$ $I=\\{1,\ldots,|I|\\}$: Index set for Injection Points in GPN.
$-$ $C=\\{1,\ldots,|C|\\}$: Index set for the set EC.

####  3.1.2. Production parameters

$-$ $w_{i}$: waste produced at WS $i\in N$.
$-$ $\mathrm{W}=\displaystyle\sum_{i\in N}w_{i}$: Overall production of waste
in all the WS.
$-$ $\delta$: Proportion of waste transformed into DOW at the PT plants
($0<\delta<1$).
$-$ $\gamma_{1}$: Proportion of DOW transformed into biomethane at the AD
plants ($0<\gamma_{1}<1$).
$-$ $\gamma_{2}$: Proportion of DOW transformed into fertilizers at the AD
plants ($0<\gamma_{2}<1$). We assume that $\gamma_{1}+\gamma_{2}<1$.
$-$ $\beta$: Proportion of biomethane transformed into LNG at the BL plants
($0<\beta<1$).
$-$ $D_{\ell}$: Demand of EC $\ell\in C$.
$-$ $p$: Lower bound percent of the total waste produced at the WS that must
be collected and sent to the PT plants ($0<p\leq 1$).
$-$ $q$: Lower bound percent of the total biomethane produced at the AD plants
that must be injected into the GPN ($0<q\leq 1$).
$-$ $R_{i}$: Proportion of the total amount of fertilizer to be deliver to WS
$i\in N$ ($\displaystyle\sum_{i\in N}R_{i}\leq 1$).
$-$ $\alpha_{\ell}$: Proportion of the total amount of LNG at the BL plants
that will be received by demand point $\ell\in C$ ($\displaystyle\sum_{\ell\in
C}\alpha_{\ell}=1$).

####  3.1.3. Costs

Set-up Costs:

These costs represent the installation costs of the different type of plants
and links. They may include not only the construction cost of the different
plants and pipelines, but also other costs related with the maintenance and
use of them, as labor costs, land costs, etc.

$-$ $f^{1}_{j}$: Set-up cost for opening PT plant $j\in P_{1}$.
$-$ $f^{2}_{j}$: Set-up cost for opening AD plant $j\in P_{2}$.
$-$ $f^{3}_{j}$: Set-up cost for opening BL plant $j\in P_{3}$.
$-$ $h_{jk}$: Set-up cost for installing a pipeline linking AD plant $j\in
P_{2}$ with BL plant $k\in P_{3}$ .
$-$ $g_{j\ell}$: Set-up cost for installing a pipeline linking AD plant $j\in
P_{2}$ with injection point in GPN $\ell\in I$.
$-$ $\mathrm{B}$: Total budget for installing plants and pipelines.

Transportation Costs:

These costs represent the transportation costs of the different products
(waste, DOW, LNG, and digester solid) that are usually transported by trucks
from the different geographical locations of the elements in the system.

$-$ $c^{1}_{ij}\geq 0$: Unit transportation cost of waste from WS $i\in N$ to
PT plant $j\in P_{1}$.
$-$ $c^{2}_{jk}\geq 0$: Unit transportation cost of DOW from PT plant $j\in
P_{1}$ to AD plant $k\in P_{2}$.
$-$ $c^{C}_{j\ell}\leq 0$: Unit transportation cost (profit) of LNG from BL
plant $j\in P_{3}$ to EC $\ell\in C$. Note that this non positive cost
represent the benefit of delivering to the EC each unit of liquified gas for
their particular use.
$-$ $c^{N}_{ji}\geq 0$: Unit transportation cost of digester solid from AD
plant $j\in P_{2}$ to WS $i\in N$.

Note that the transportation costs for the biogas through the pipelines are
already included in the set-up costs described above.

###  3.2. Variables

Our model makes decision on the different plants and pipelines that are
opened, as well as the amount of product delivered at the different phases of
the process.

Binary variables:

The following variables determine whether a plant or a pipeline should be
installed.

$y_{j}^{1}=\begin{cases}1&\mbox{if PT $j$ is open,}\\\
0&\mbox{otherwise}\end{cases}\ \ \forall j\in P_{1},\ \ \quad
y_{j}^{2}=\begin{cases}1&\mbox{if AD $j$ is open,}\\\
0&\mbox{otherwise}\end{cases}\ \ \forall j\in P_{2},$

$y_{j}^{3}=\begin{cases}1&\mbox{if BL $j$ is open,}\\\
0&\mbox{otherwise}\end{cases}\ \ \forall j\in P_{3},\ \ \quad
z_{jk}=\begin{cases}1&\mbox{if a pipeline linking AD $j$ }\\\ &\mbox{ with BL
$k$ is built $\qquad\forall j\in P_{2},k\in P_{3}$,}\\\
0&\mbox{otherwise}\end{cases}\ \ $

$t_{j\ell}=\begin{cases}1&\mbox{if a pipeline linking AD plant $j$ with }\\\
&\mbox{injection point in GPN $\ell$ is built }\\\
0&\mbox{otherwise}\end{cases}\qquad\forall\,j\in P_{2},\quad\ell\in I,$

Continuous Variables:

The following variables of our model decide the amount of product delivered
between the different sites in the system.

$x_{ij}^{1}:$ amount of waste delivered from WS $i$ to the PT plant $j$,
$\forall$ $i\in N$, $j\in P_{1}$.

$x_{jk}^{2}:$ amount of DOW delivered from PT plant $j$ to AD plant $k$,
$\forall$ $j\in P_{1}$, $k\in P_{2}$.

$x_{jk}^{3}:$ amount of biomethane delivered from AD plant $j$ to BL plant
$k$, $\forall$ $j\in P_{2}$, $k\in P_{3}$.

$x_{j\ell}^{I}:$ amount of biomethane delivered from AD plant $j$ to GPN’s
injection point $\ell$, $\forall$ $j\in P_{2}$, $\ell\in I$.

$x_{ji}^{N}:$ amount of digester solid delivered from AD $j$ to WS $i$,
$\forall$ $j\in P_{2}$, $i\in N$.

$x_{j\ell}^{C}:$ amount of LGN delivered from BL plant $j$ to EC $\ell$,
$\forall$ $j\in P_{3}$, $\ell\in C$.

####  3.2.1. Objective Function

The goal of the W2BLP is to minimize the total transportation cost of the
system. For that, one must decide the optimal position of the different plants
and pipelines with the given budget $B$. These locations have a direct impact
in the transportation cost.

With the above set of variables, the overall transportation cost can be
written as follows:

$\displaystyle\displaystyle\sum_{i\in N}\displaystyle\sum_{j\in
P_{1}}c_{ij}^{1}x^{1}_{ij}+\displaystyle\sum_{j\in
P_{1}}\displaystyle\sum_{k\in
P_{2}}c_{jk}^{2}x^{2}_{jk}+\displaystyle\sum_{j\in
P_{2}}\displaystyle\sum_{i\in N}c^{N}_{ji}x^{N}_{ji}+\displaystyle\sum_{j\in
P_{3}}\displaystyle\sum_{\ell\in C}c^{C}_{j\ell}x^{C}_{j\ell}.$ (1)

###  3.3. Constraints

The assumptions of our problem are adequately established by the following
constraints:

  * •

The budget for installing the different plants and pipelines must be
satisfied:

$\displaystyle\sum_{j\in P_{1}}f^{1}_{j}y_{j}^{1}+\displaystyle\sum_{j\in
P_{2}}f_{j}^{2}y_{j}^{2}+\displaystyle\sum_{j\in
P_{3}}f^{3}_{j}y_{j}^{3}+\displaystyle\sum_{j\in P_{2}}\displaystyle\sum_{k\in
P_{3}}h_{jk}z_{jk}+\displaystyle\sum_{j\in P_{2}}\displaystyle\sum_{\ell\in
I}g_{j\ell}t_{j\ell}\leq\mathrm{B}$ (2)

  * •

Flow Conservation Constraints: The product must be adequately routed through
the intermediate stages of the supply chain and enforced by the following
constraints:

$\displaystyle\sum_{k\in P_{2}}x^{2}_{jk}=\delta\displaystyle\sum_{i\in
N}x^{1}_{ij},\forall j\in P_{1},$ (3) $\displaystyle\sum_{k\in
P_{3}}x^{3}_{jk}+\displaystyle\sum_{\ell\in
I}x^{I}_{j\ell}=\gamma_{1}\displaystyle\sum_{k\in P_{1}}x^{2}_{kj},\forall
j\in P_{2}.$ (4) $\displaystyle\sum_{i\in
N}x^{N}_{ji}=\gamma_{2}\displaystyle\sum_{k\in P_{1}}x^{2}_{kj},\forall j\in
P_{2},$ (5) $\displaystyle\sum_{\ell\in
C}x^{C}_{j\ell}=\beta\displaystyle\sum_{k\in P_{2}}x^{3}_{kj},\forall j\in
P_{3},$ (6)

Constraints (3) ensure that all the waste received at a PT plant is delivered
to the AD plants once it is processed ($100\delta\%$ of the waste is
transformed into DOW). $100\gamma_{1}\%$ of the amount of DOW received at each
AD plant is transformed into biomethane and fully delivered either to the
injection points or to the BL plants (Constraint (4)) and $100\gamma_{2}\%$ is
converted into digester solid and delivered to the WS (Constraint (5)).
Finally, Constraints (6) assure that the amount of biomethane received at the
BL plants is converted into $100\beta\%$ of LNG and fully ddelivered to the
EC.

  * •

Demand and Capacity Constraints: The different products obtained in the
process are adequately delivered according to demands and capacities.
Specifically, the WS must send waste and receive digester solid, the EC must
receive LNG and the GPN biomethane. To this end, the following constraints are
incorporated to our model:

$\displaystyle\sum_{j\in P_{1}}x^{1}_{ij}\leq w_{i},\forall i\in N,$ (7)
$\displaystyle\sum_{i\in N}\displaystyle\sum_{j\in P_{1}}x^{1}_{ij}\geq pW$
(8) $\displaystyle\sum_{j\in P_{2}}\displaystyle\sum_{\ell\in
I}x^{I}_{j\ell}\geq q\gamma_{1}\displaystyle\sum_{k\in
P_{1}}\displaystyle\sum_{j\in P_{2}}x^{2}_{kj}$ (9) $\displaystyle\sum_{j\in
P_{3}}x^{C}_{j\ell}\leq D_{\ell},\forall\ell\in C.$ (10)

Constraints (7) assure that the amount delivered from each WS is at most the
waste it produces. Constraint (8) forces to send, from the WS’s, at least
$100p\%$ of the total amount of waste. Constraint (9) states that at least
$100q\%$ of the produced biomethane must be injected in the GPN. Each EC, by
Constraints (10), at most receive its required demand. Note that the
transportation cost for the product in this phase is non positive, and then,
if possible, this demand will be satisfied.

Note also that, by (3), Constraints (8) and (9) can be equivalently rewritten
as:

$\displaystyle\sum_{j\in P_{2}}\displaystyle\sum_{\ell\in I}x^{I}_{j\ell}\geq
pq\delta\gamma_{1}W$ (11)

It may happen that if $q$ is large, the overall demand for the EC that can be
satisfied is small. Thus, one can establish a sharing rule for the available
LNG, i.e., a vector $(\alpha_{1},\ldots,\alpha_{|C|})\in[0,1]^{|C|}$
representing the percent of the LNG at the BL plants that will be received by
each EC. The following constraints assure the adequate verification of this
share:

$\displaystyle\sum_{k\in
P_{3}}x^{C}_{k\ell}=\alpha_{\ell}\beta\displaystyle\sum_{j\in
P_{2}}\displaystyle\sum_{k\in P_{3}}x^{3}_{jk},\forall\ell\in C.$ (12)

Finally, Constraints (13) force that each WS receives its given percentage of
the total production of digester solid:

$\displaystyle\sum_{j\in
P_{2}}x^{N}_{ji}=R_{i}\gamma_{2}\displaystyle\sum_{j\in
P_{1}}\displaystyle\sum_{k\in P_{2}}x^{2}_{jk},\forall i\in N,$ (13)

  * •

Distribution through open Plants and Links: The transportation of the product
in the different phases forces the plants and/or pipelines to be installed.
This is assured by the following constraints:

$\displaystyle\displaystyle\sum_{i\in N}x^{1}_{ij}$ $\displaystyle\leq
pWy_{j}^{1},$ $\displaystyle\forall j\in P_{1},$ (14)
$\displaystyle\displaystyle\sum_{k\in P_{1}}x^{2}_{kj}$ $\displaystyle\leq
p\delta Wy_{j}^{2},$ $\displaystyle\forall j\in P_{2},$ (15) $\displaystyle
x^{3}_{jk}$ $\displaystyle\leq p\delta\gamma_{1}Wz_{jk},$
$\displaystyle\forall j\in P_{2},\forall k\in P_{3},$ (16) $\displaystyle
x^{I}_{j\ell}$ $\displaystyle\leq p\delta\gamma_{1}Wt_{j\ell},$
$\displaystyle\forall j\in P_{2},\forall\ell\in I.$ (17)

These constraints assure that in case a plant or a link is not open, the
product cannot be distributed through the corresponding plant or link.

  * •

Compatibility of open plants and links: The installation of a pipeline linking
a plant with other plant or with an injection point in the GPN, is subject to
the previous instalation of the plants that the pipelines are connecting.
These conditions are imposed in our model with the following linear
constraints:

$\displaystyle z_{jk}\leq y_{k}^{3},\forall j\in P_{2},k\in P_{3}$ (18)
$\displaystyle z_{jk}\leq y_{j}^{2},\forall j\in P_{2},k\in P_{3}$ (19)
$\displaystyle t_{j\ell}\leq y_{j}^{2},\forall j\in P_{2},\ell\in I$ (20)

Summarizing all the previous descriptions, the W2BLP can be model with the
following mathematical optimization problem that we denote by $(W2BLP)_{0}$:

$\displaystyle\min$ $\displaystyle\displaystyle\sum_{i\in
N}\displaystyle\sum_{j\in P_{1}}c_{ij}^{1}x^{1}_{ij}+\displaystyle\sum_{j\in
P_{1}}\displaystyle\sum_{k\in
P_{2}}c_{jk}^{2}x^{2}_{jk}+\displaystyle\sum_{j\in
P_{2}}\displaystyle\sum_{i\in N}c^{N}_{ji}x^{N}_{ji}+\displaystyle\sum_{j\in
P_{3}}\displaystyle\sum_{\ell\in C}c^{C}_{j\ell}x^{C}_{j\ell}$ s.t.
$\displaystyle\eqref{con:setup}-\eqref{flow:1},\eqref{con:demand}-\eqref{con:zp-y2},$
$\displaystyle y^{1}_{j}\in\\{0,1\\},\forall j\in P_{1},$ $\displaystyle
y^{2}_{j}\in\\{0,1\\},\forall j\in P_{2},$ $\displaystyle
y_{j}^{3}\in\\{0,1\\},\forall j\in P_{3},$ $\displaystyle
z_{jk}\in\\{0,1\\},\forall j\in P_{2},k\in P_{3},$ $\displaystyle
t_{j\ell}\in\\{0,1\\},\forall j\in P_{2},\ell\in I,$ $\displaystyle
x_{ij}^{1}\geq 0,\forall i\in N,j\in P_{1},$ $\displaystyle x_{jk}^{2}\geq
0,\forall j\in P_{1},k\in P_{2},$ $\displaystyle x_{jk}^{3}\geq 0,\forall j\in
P_{2},k\in P_{3},$ $\displaystyle x_{j\ell}^{I}\geq 0,\forall j\in
P_{2},\ell\in I,$ $\displaystyle x_{ji}^{N}\geq 0,\forall j\in P_{2},i\in N,$
$\displaystyle x_{j\ell}^{C}\geq 0,\forall j\in P_{3},\ell\in C.$

The above model is a mixed integer linear programming (MILP) problem, which
can be difficult to solve in commercial solvers for real-world instances. In
the following result we show that the size of the model can be considerably
reduced based on the fact that once decided the potential AD plant which is
open, it will optimally distribute the biogas to the GPN through its less
costly injection point.

######  Lemma 3.1.

Let
$\bar{S}=(\bar{y}^{1},\bar{y}^{2},\bar{y}^{3},\bar{z},\bar{t},\bar{x}^{1},\bar{x}^{2},\bar{x}^{3},\bar{x}^{I},\bar{x}^{N},\bar{x}^{C})$
be an optimal solution of the model $(W2BLP)_{0}$. If $\bar{t}_{j\ell}=1$ for
some $j\in P$ and $\ell\in I$ then, there exists an optimal solution
$\hat{S}=(\bar{y}^{1},\bar{y}^{2},\bar{y}^{3},\bar{z},\hat{t},\bar{x}^{1},\bar{x}^{2},\bar{x}^{3},\hat{x}^{I},\bar{x}^{N},\bar{x}^{C})$
with:

$\hat{t}_{j\ell(j)}=1,\quad\hat{x}^{I}_{j\ell(j)}=\sum_{\ell\in
I}\bar{x}^{I}_{j\ell},\qquad\hat{t}_{j\ell}=0,\quad\hat{x}^{I}_{j\ell}=0,\quad\forall\ell\neq\ell(j),$

where $\ell(j)=\arg\min_{\ell\in I}g_{j\ell}$.

Furthemore, the overall set-up costs of $\hat{S}$ is smaller or equal than the
ones for $\bar{S}$.

###### Proof.

The proof follows straightforward since there are no transportation costs
between the AD plants and the injections points but linking them is accounted
in the budget constraint. Thus, replacing the link $(j,\ell)$ by $(j,\ell(j))$
does not affect the objective function but reduces the set-up costs. ∎

The result above allows us to simplify the model by reducing the number of
variables modeling the links and the flows between the AD plants and the GPN.
Specifically, one can replace the two-index variables $t$ and $x^{I}$ by one-
index variables that, abusing of notation, we denote as:

$t_{j}=\begin{cases}1&\mbox{if a pipeline linking AD plant $j$ with }\\\
&\mbox{injection point $l(j)$ is built }\\\
0&\mbox{otherwise}\end{cases}\qquad\forall\,j\in P_{2},$

and

$x_{j}^{I}:$ amount of biomethane delivered from AD plant $j$ to GPN’s
injection point $l(j)$, $\forall\,j\in P_{2}$.

Thus, if we denote by $g_{j}=g_{j\ell(j)}$ for all $j\in P_{2}$, one can
replace the two-indices $t$ and $x^{I}$variables in model $(W2BLP)_{0}$ by the
one-index $t$ and $x^{I}$variables above. Consequently, Constraints (2), (4),
(11), and (17) are replaced by the following inequalities and equations:

$\displaystyle\displaystyle\sum_{j\in
P_{1}}f^{1}_{j}y_{j}^{1}+\displaystyle\sum_{j\in
P_{2}}f_{j}^{2}y_{j}^{2}+\displaystyle\sum_{j\in
P_{3}}f^{3}_{j}y_{j}^{3}+\displaystyle\sum_{j\in P_{2}}\displaystyle\sum_{k\in
P_{3}}h_{jk}z_{jk}+\displaystyle\sum_{j\in P_{2}}g_{j}t_{j}\leq\mathrm{B},$
(21) $\displaystyle\displaystyle\sum_{k\in
P_{3}}x^{3}_{jk}+x^{I}_{j}=\gamma_{1}\displaystyle\sum_{k\in
P_{1}}x^{2}_{kj},\forall j\in P_{2}.$ (22)
$\displaystyle\displaystyle\sum_{j\in P_{2}}x^{I}_{j}\geq
pq\delta\gamma_{1}W,$ (23) $\displaystyle x^{I}_{j}\leq
p\delta\gamma_{1}Wt_{j},\forall j\in P_{2}.$ (24)

As one can observe in the next section, the reduced model above is able to
solve real instances of the W2BLP in reasonable CPU time. Our approach is then
a useful tool for making decisions on the logistic of biogas production
system.

##  4\. Computational Experiments

We have run a series of experiments to analyze the computational performance
of our approach. The main goal of the experiments is to determine the
computational limitations of our model and its ability to obtain solutions for
real-world instances.

We have randomly generated different instances with different sizes and
parameters. We generate the coordinates for the WS, the potential location for
the different plants (PT, AD,and BL), the injection points in GPN and the EC
uniformly in $[0,1000]\times[0,1000]$. For the sake of simplification, we
assume that the number and the location of WS ($n$) is the same that the
number and the location of potential PT plants ($m_{1}$), and that the number
and the location of potential AD plants ($m_{2}$) and potential BL plants
($m_{3}$) are also the same. Aditionally, the number of EC ($|C|$) also
coincides with the number of injection points in GPN ($|I|$). The value of
$n=m_{1}$ ranges in $\\{25,50,100,200,500\\}$, $m:=m_{2}=m_{3}$ ranges in
$\\{5,10,20,50,100\\}$ with
$m\in[\lceil\frac{n}{10}\rceil,\lfloor\frac{n}{2}\rfloor)$, and $d:=|I|=|C|$
ranges in $\\{10,20,50,100\\}$ with $d\leq\lceil\frac{n}{2}\rceil$ (here
$\lceil\cdot\rceil$ and $\lfloor\cdot\rfloor$ stand for the ceiling and floor
integer rounding functions, respectively).

Let ${\rm dist}_{ij}$ denote the Euclidean distance between locations $i$ and
$j$. We considered as unit transportation costs the Euclidean distances
between the different points (${\rm dist}_{ij}$), except for the
transportation costs of LNG from BL plants to the EC where the unit
transportation cost is considered as a profit and it is defined as:

$c_{j\ell}^{C}=-\frac{\displaystyle\max_{k\in P_{3},l\in C}({\rm
dist}_{kl})^{2}}{1+{\rm dist}_{j\ell}}$

These costs are designed to represent that as closer the extra customer to the
open BL plant the larger the profit to satisfy its unit demand.

The set-up costs for the PT plants and the BL plants were chosen all equal to
one unit (in million of $) and the installation cost for the AD plants wer
fixed to $5$ units (in million of $). The set-up costs for the pipelines
linking the AD plants with the injection points in GPN and the BL plants are
0.1 times the normalized (by the maximum) Euclidean distance between the
corresponding points.

For determining adequate budgets for a given instance, we first solve the
problem without the budget constraint (2). After solving the problem, we
compute its effective set-up cost by adding up the costs of the open plants
and pipelines that route a positive flow. We denote by $\widehat{B}$ this set-
up cost and consider as budget $\mathrm{B}=0.2\widehat{B}$. This budget
indicates that one can only use $20\%$ of the cost of makig the best decision
with no budget.

The slurry production at each farm ($w_{i}$ for $i\in N$) has been uniformly
generated in $[1,100]\cap\mathbb{Z}_{+}$ (here $\mathbb{Z}_{+}$ stands for the
set of nonnegative integer numbers). The percent of the overall produced
biomethane that must be injected into the existing GPN, $q$, ranges in
$\\{50\%,70\%,90\%\\}$. We have considered parameters $p=1$, $\delta=0.8$,
$\gamma_{1}=0.8$, $\gamma_{2}=0.15$, $\beta=0.7$,
$\alpha_{\ell}=\frac{D_{\ell}}{\sum_{l\in C}D_{l}}$ for all $\ell\in C$, and
$R_{i}=\frac{w_{i}}{W}$ for all $i\in N$. The LNG demand of each EC
($D_{\ell}$ for $\ell\in C$) was randomly generated in
$[\delta\gamma^{1}\beta,100\delta\gamma^{1}\beta]$.

The values of the parameters that we consider in our experiments are
summarized in Table 1.

Coordinates $n=m_{1}$ $m:=m_{2}=m_{3}$ $d:=|I|=|C|$ Unif$[0,1000]^{2}$
$\\{25,50,100,200,500\\}$ $\\{5,10,20,50,100\\}$ with
$m\in[\lceil\frac{n}{10}\rceil,\lfloor\frac{n}{2}\rfloor)$
$\\{10,20,50,100\\}$ with $d\leq\lceil\frac{n}{2}\rceil$
$c^{1}_{ij},\;c^{2}_{ij},\;c^{N}_{ij}$ $c^{C}_{j\ell}$ $f^{2}$ $f^{1}$,
$f^{3}$ ${\rm dist}_{ij}$ $-\frac{\displaystyle\max_{k\in P_{3},l\in C}({\rm
dist}_{kl})^{2}}{1+{\rm dist}_{j\ell}}$ $5$ $1$ $h_{jk}$ $g_{j}$ $B$ $w_{i}$
$0.1\times\frac{{\rm dist}_{jk}}{\max_{j^{\prime},k^{\prime}}{\rm
dist}_{j^{\prime}k^{\prime}}}$ $0.1\times\min_{\ell}\frac{{\rm
dist}_{j\ell}}{\max_{k^{\prime},\ell^{\prime}}{\rm
dist}_{k^{\prime}\ell^{\prime}}}$ $0.2\widehat{B}$
Unif$[1,100]\cap\mathbb{Z}_{+}$ $q$ $\delta=\gamma_{1}$ $\gamma_{2}$ $\beta$
$\\{50\%,70\%,90\%\\}$ $0.8$ $0.15$ $0.7$ $p$ $\alpha_{\ell}$ $R_{i}$
$D_{\ell}$ 1 $\frac{D_{\ell}}{\sum_{l\in C}D_{l}}$ $\frac{w_{i}}{W}$
Unif$[\delta\gamma^{1}\beta,100\delta\gamma^{1}\beta]$

Table 1. Summary of the parameters used in our experiments.

The model has been coded in Python 3.7 in an iMac with 3.3GHz with an Intel
Core i7 with 4 cores and 16GB 1867 MHz DDR3 RAM. We used Gurobi 9.1.2 as
optimization solver. A time limit of 6 hours was fixed for all the instances.
All the instances with $n\leq 100$ were optimally solved. For the larger
instances we fixed a MIP Gap limit (the solver stops when reaching such a
limit and outputs the best feasible solution). For $n=200$ we fix a MIPGap
limit of $2\%$, for $n=500$ of $4\%$. For each combination of parameters
$n,m,d,q$ we solved five instances. Thus, we have solved a total of $450$
instances.

(a) CPU times (log scale).

(b) MIP Gaps.

Figure 2. Performance profile of our computational experience for the
different values of parameter $q$.

Figure 2 depicts performance profiles with respect to CPU times and MIP Gaps.
In these pictures, we plots the percentage of instances solved up to each
value of the CPU or MIP Gap, respectively. We analyze the computational
performance of the model for the different values of the parameter $q$. As one
can observe, the model seems to perform slightly better for the the smaller
values of parameter $q$, although the differences are tiny. For this reason,
in what follows we will not differentiate between the different values of the
parameter $q$.

In Table 2 we summarize the results of our computational experiments. The
first three columns indicate the size of the instances. Each of the rows
summarize the results of a total of 15 instances (five random instances and
three values of the parameter $q$). Column CPUTime (secs) is the average CPU
time, in seconds, that the solver required to solve the instances. Column
MIPGap indicates the average percent MIP Gap. In case the instance is
optimally solved, the MIP Gap is $0\%$. Otherwise, this number reports either
the MIPGap obtained within the time limit (if it is greater than the MIP Gap
limit) or the MIP Gap when the MIP Gap limit was reached (that can be slightly
smaller than the MIP Gap limit). Finally, column UnSolved indicates the
percent instances summarized in the row that were not optimally solved or
reaching the MIP Gap limit within the time limit.

$n$ $m$ $d$ CPUTime (secs) MIPGap %UnSolved 25 5 10 0.67 0% 0% 10 10 1.71 0%
0% 25 Total 1.19 0% 0% 50 5 10 1.58 0% 0% 20 2.24 0% 0% 10 10 4.55 0% 0% 20
4.70 0% 0% 20 10 24.78 0% 0% 20 28.83 0% 0% 50 Total 11.11 0% 0% 100 10 10
32.08 0% 0% 20 48.22 0% 0% 50 28.14 0% 0% 20 10 109.59 0% 0% 20 181.59 0% 0%
50 179.04 0% 0% 100 Total 96.44 0% 0% 200 20 10 444.32 1.47% 0% 20 137.37
1.69% 0% 50 70.65 1.44% 0% 100 51.63 1.40% 0% 50 10 1195.32 1.77% 0% 20 413.42
1.81% 0% 50 46.63 1.54% 0% 100 40.78 1.45% 0% 200 Total 300.01 1.57% 0% 500 50
10 6946.33 4.45% 20.00% 20 1041.76 2.90% 0% 50 912.06 1.93% 0% 100 713.42
1.74% 0% 100 10 5656.14 3.10% 6.67% 20 3625.32 2.50% 0% 50 2366.02 2.47% 0%
100 3195.16 1.32% 0% 500 Total 3057.03 2.55% 3.33% Total 916.80 1.10% 0.89%

Table 2. Summary of our computational experience.

As mentioned before, we observe that all the instances with up to $n=100$ have
been optimally solved within the time limit. As expected, the computing time
increases with the value of $n$ and $m$, while it does not seem to depend on
the value of $d$. The average computing time over all these instances is
approximately $46.26$ seconds, being $1.19$ seconds for the instances with
$n=25$ , $11.11$ seconds for $n=50$, and $96.44$ for $n=100$. The average
computing time needed to obtain an optimal solution depending on the value of
$m$ is $1.5$ seconds for $m=5$, $19.9$ seconds for $m=10$, and $104.8$ seconds
for $m=20$.

Regarding the instances with $n=200$ we observe that all the instances have
been solved, which in this case means that all the instances have reached the
MIP gap limit within the time limit, being the average consuming time over all
the instances $300$ seconds and the average MIP gap $1.57\%$. One can observe
that, for these instances, the computing time is greater as the lower is the
value of the parameter $d$. This may be because the value of parameter $d$
does not affect to the number of binary variables in our improved formulation
and it might happen that the smaller the number of EC and injection points,
the more difficult it is to decide where to install the different plants to
satisfy them at optimal cost. This fact can also be observed for the largest
instances with $n=500$. For $d=10$, $20\%$ of the instances with $m=50$ and
$6.67\%$ of the instances with $m=100$ have not been solved within the time
limit, that is, these instances have not reached the MIP gap of $4\%$ in 6
hours. The rest of the instances have been solved within the time limit. The
average computing time for all the instances with $n=500$ was $3057$ seconds
and the average MIP gap, $2.55\%$.

The results shown in Table 2, together with Figure 2, allow us to get a clear
empirical evidence that using our improved formulation, the problem can be
solved (with the MIP gap limit) within the time limit using a commercial
solver. In fact, only 4 of the 450 instances (all of them with $n=500$) have
reached the time limit with a MIP Gap greater than the MIP Gap limit. This
number represents less than $1\%$ of the total instances and only $3.3\%$ of
the instances with $n=500$.

From our experiments, we conclude that our approach provides a useful tool to
make decisions about on the complex logistic problem behind the energy
transformation of waste into biogas in realistic instances using reasonable
computational resources. Thus, our model can also be used to evaluate
different alternatives based on different values for the set-up costs for the
plants or different transportation costs.

##  5\. Case Study

Figure 3. Lake Mendota in the Upper Yahara watershed region in Dane County, WI
(source: (Ma and Zavala, 2022)).

We tested our model in a real-world dataset based on the region of upper
Yahara watershed region in the State of Wisconsin (see Figure 3) whose data is
available at https://github.com/zavalab/JuliaBox/tree/master/Graph_S%26C. A
detailed description of this dataset can be found in (Ma and Zavala, 2022;
Sampat et al., 2019, 2021; Tominac et al., 2022). This region consists of 203
dairy farms whose waste production is known and whose locations are
predefined. These locations determine our set of WS and the potential
positions for the PT plants. The potential positions for the AD plants and the
BL were randomly generated in that region. The positions of the injections
points on the GPN were also randomly generated in this region and a network
connecting this points was constructed. The position of the EC where randomly
generated outside that region but in the minimum rectangular box containing
it.

In Figure 4 we plot the coordinates of the different agents involved in the
W2BLP. Red points in the plot indicate the WS locations (as larger the size of
the dot, larger the production of waste in the farm) and the potential
locations for PT, blue points are the potential locations for AD plants and
for BL plants, green points indicate the EC locations, grey dotted lines
represent the GPN and grey points the injection points of the GPN.

Figure 4. Input Data of the Case Study.

The parameters not provided in the above repository were fixed as follows:

  * •

$c^{1}_{ij}=0.3\times{\rm dist}_{ij}$, for $i\in N$ and $j\in P_{1}$.

  * •

$c^{2}_{jk}=0.15\times{\rm dist}_{jk}$, for $j\in P_{1}$ and $k\in P_{2}$.

  * •

$c^{N}_{ji}=0.15\times{\rm dist}_{ji}$, for $j\in P_{2}$ and $i\in N$.

  * •

$D_{\ell}$ was randomly generated in $[\delta\gamma^{1}\beta\times\min
w_{i},\delta\gamma^{1}\beta\times\max w_{i}]$.

  * •

The set-up budget was fix to $\kappa\times\widehat{B}$, for
$\kappa\in\\{0.1,0.2\\}$, where $\widehat{B}$ was computed as described in
Section 4:

$\widehat{B}=\begin{cases}472.33&\mbox{if $q=50\%$,}\\\ 472.45&\mbox{if
$q=70\%$,}\\\ 472.32&\mbox{if $q=90\%$,}\\\ \end{cases}$

  * •

The rest of parameters have been considered as described in Section 4 (q
ranges in {50%, 70%, 90%}, $p=1$, $\delta=0.8$, $\gamma_{1}=0.8$,
$\gamma_{2}=0.15$, $\beta=0.7$, $\alpha_{\ell}=\frac{D_{\ell}}{\sum_{l\in
C}D_{l}}$ for all $\ell\in C$, and $R_{i}=\frac{w_{i}}{W}$ for all $i\in N$).

The parameters required to reproduce the obtained results are available at
https://github.com/vblancoOR/w2blp.

$B$ | $q$ | CPU Time | MIPGap
---|---|---|---
$0.1\widehat{B}$ | 0.5 | 2595.82 | 0%
0.7 | 18093.03 | 0%
0.9 | 21602.8 | 0.1%
$0.2\widehat{B}$ | 0.5 | 84.67 | 0%
0.7 | 473.81 | 0%
0.9 | 3106.91 | 0%
Table 3. CPU times (in seconds) and MIP gaps for solving the instances in the
case study.

As in the previous section, we set a time limit of 6 hours for solving the
problem. In Table 3 we show the CPU times (in seconds) and MIP Gaps obtained
after running our model to the six different configurations of budget and $q$.
In that table, one can observe that except the case $B=0.1\widehat{B}$ with
$q=0.9$, all the instances were optimally solved within the time limit. The
case were the optimality is not guaranteed obtained a MIP gap of $0.1\%$ which
is insignificant. Furthermore, as expected, a more limited budget has a
significative impact in the CPU time required to solve the problem, being the
problems with budget $0.1\widehat{B}$ more challenging than those with budget
$0.2\widehat{B}$.

In Figures 5, 6, and 7 we show the results of our experiment for the Yahara
watershed dataset for the different choices of budget, $B$, and $q$, and the
different decisions that are made by our model.

In Figure 5, we show the results about the optimal locations for the PT plants
(thick red points), for the AD plants (blue squares), for the BL plants (blue
triangles) and for the links connecting AD plants with BL plants (blue lines)
and AD plants with injection points of the GPN (black lines). As can be
observed, the budget and the percent of production that must be injected in
the GPN has a direct impact in the structure of the obtained solutions. On the
one hand, as expected, as smaller the budget, smaller the number of plants and
pipelines that are open. In particular, the budget mostly affect to the
construction of pipelines to either inject to the GPN or transport to the BL,
but also in the number of installed plants. We also observed that the network
of the different installed plants and pipelines has more connected components
for the larger budget, whereas is almost connected for the small budget.
Regarding the value of $q$, it seems that as larger the value of $q$ closer
the new facilities to the GPN.

(a) $q=50\%$ and $\mathrm{B}=0.2\widehat{B}$.

(b) $q=50\%$ and $\mathrm{B}=0.1\widehat{B}$.

(c) $q=70\%$ and $\mathrm{B}=0.2\widehat{B}$.

(d) $q=70\%$ and $\mathrm{B}=0.1\widehat{B}$.

(e) $q=90\%$ and $\mathrm{B}=0.2\widehat{B}$.

(f) $q=90\%$ and $\mathrm{B}=0.1\widehat{B}$.

Figure 5. Optimal location of the different types of plants and pipelines for
different values of $B$ and $q$.

In Figure 6 we show the distribution network between the WS plants, the PT
plants, and the AD plants. Specifically, the links for which there is a
positive waste flow from WS to PT plants (red lines), a DOW flow from PT
plants to AD plants (purple lines) and digester solid flow from AD plants to
WS (yellow lines). The main observation that can be drawn from these plots is
that in most of the cases there are differentiated clusters of WS plants
(farms) that share the same PT plants and AD plants. This observations implies
that for larger datasets where the model could not be able to solve the
problem, the optimal solution would be adequately approximated by clustering
the WS (with an adequate criterion) and solve the problem separately for each
cluster.

(a) $q=50\%$ and $\mathrm{B}=0.2\widehat{B}$.

(b) $q=50\%$ and $\mathrm{B}=0.1\widehat{B}$.

(c) $q=70\%$ and $\mathrm{B}=0.2\widehat{B}$.

(d) $q=70\%$ and $\mathrm{B}=0.1\widehat{B}$.

(e) $q=90\%$ and $\mathrm{B}=0.2\widehat{B}$.

(f) $q=90\%$ and $\mathrm{B}=0.1\widehat{B}$.

Figure 6. Distribution network between WS, PT plants, and AD plants.

Finally, in Figure 7, we show the biomethane flow from AD plants to BL plants
(blue lines) and to injection points of the GPN (black lines), and LGN flow
from BL plants to EC (green lines). Note that some of the AD plants are
devoted only to give service to the GPN, others that serve exclusively the BL
plants (and then the EC), but one can also find AD plants that send part of
the production to the GPN and what remains, to the BL plants. This behaviour
would have never happen if an integrated model as the one we propose would not
have been considered. Note also that a single AD plant is allowed to send
biomethane to different BL plants.

(a) $q=50\%$ and $\mathrm{B}=0.2\widehat{B}$.

(b) $q=50\%$ and $\mathrm{B}=0.1\widehat{B}$.

(c) $q=70\%$ and $\mathrm{B}=0.2\widehat{B}$.

(d) $q=70\%$ and $\mathrm{B}=0.1\widehat{B}$.

(e) $q=90\%$ and $\mathrm{B}=0.2\widehat{B}$.

(f) $q=90\%$ and $\mathrm{B}=0.1\widehat{B}$.

Figure 7. Distribution network between AD plants, BL plants, injection points
of the GPN, and EC.

##  6\. Conclusions

In this paper, we propose a mathematical optimization approach to design an
optimal logistic system to distribute the different products involved in the
generation of biogas from waste. Specifically, given a set of waste storage
centers, an existing gas pipeline network, and a set of external customers, we
provide a decision aid tool to determine the number and optimal locations of
pre-treatment plants, anaerobic digestion plants, and biomethane liquefaction
plants, as well as the pipelines linking some of these plants. Additionally,
we provide a distribution plan to send the different wastes, to serve either
the gas pipeline network or the external customers with the produced biogas,
and to serve the waste storage centers with the generated fertilizer. All the
decisions are made by minimizing the transportation costs of the different
products and restricting the installation of plants and pipelines to a given
budget.

We develop a new mixed integer linear programming formulation for the problem
and prove some results that allow us to reduce the size of the model. With
this simplification, we were able to obtain optimal solutions for the problem
in real-world intances.

We report the results of an extensive battery of synthetic computational
experiments, concluding that our approach is suitable to be applied to
different settings and sizes. We also analyze the case study of the Yahara
watershed, and study the obtained solutions based on different parameters.

Our future research on this topic includes the incorporation of uncertainty in
some parameters of the model. Concretely, the production of waste in the farms
or waste storages is known to vary in the different seasons of the year, and
this production may has an impact in the solution of the problem. We will
design stochastic optimization models that take into account this uncertainty
to construct robust solution of the W2BLP. Additionally, we will consider
capacity constraints for the different plants and multiple periods for the
decisions made in our model. Instead of assuming that the plants and pipelines
are installed here-and-now, we will decide in which period of a time horizon
each of the installations is constructed and used, assuming that they have a
limited capacity and that a budget is available for each period. The
integrated model that consider uncertainty, capacity constraints and
multiperiod decisions will be closer to reality, but the mathematical
programming model for it would be prohibitive even for small instances. Thus,
a different solution approach will have to be designed to solve it, that would
not be exact, but heuristic.

## Acknowledgement

This research has been partially supported by Spanish Ministerio de Ciencia e
Innovación, AEI/FEDER grant number PID 2020 - 114594GBC21, RED2022-134149-T
(Thematic Network on Location Science and Related Problems); FEDER+Junta de
Andalucía projects B-FQM-322-UGR20, C‐EXP‐139‐UGR23, and AT 21_00032; VII
PPIT-US (Ayudas Estancias Breves, Modalidad III.2A); and the IMAG-Maria de
Maeztu grant CEX2020-001105-M /AEI /10.13039/501100011033. The authors also
acknowledge the support of the Wisconsin Institute of Discovery (University of
Wisconsin-Madison) for using its computational resources.

## References

  * Amigun and von Blottnitz (2010) B. Amigun and H. von Blottnitz.  Capacity-cost and location-cost analyses for biogas plants in Africa.  _Resources, Conservation and Recycling_ , 55(1):63–73, 2010.  ISSN 09213449.
  * Ankathi et al. (2021) S. Ankathi, D. Watkins, P. Sreedhara, J. Zuhlke, and D. R. Shonnard.  Gis-integrated optimization for locating food waste and manure anaerobic co-digestion facilities.  _ACS Sustainable Chemistry & Engineering_, 9(11):4024–4032, 2021.
  * Balaman and Selim (2014) Ş. Y. Balaman and H. Selim.  A network design model for biomass to energy supply chains with anaerobic digestion systems.  _Applied Energy_ , 130:289–304, 2014.
  * Blanco and Martínez-Antón (2024) V. Blanco and M. Martínez-Antón.  Optimal coverage-based placement of static leak detection devices for pipeline water supply networks.  _Omega_ , 122:102956, 2024.
  * Blanco et al. (2010) V. Blanco, L. Carpente, Y. Hinojosa, and J. Puerto.  Planning for agricultural forage harvesters and trucks: model, heuristics, and case study.  _Networks and Spatial Economics_ , 10:321–343, 2010.
  * Blanco et al. (2020) V. Blanco, E. Conde, Y. Hinojosa, and J. Puerto.  An optimization model for line planning and timetabling in automated urban metro subway networks. a case study.  _Omega_ , 92:102165, 2020.
  * Blanco et al. (2022) V. Blanco, G. González, Y. Hinojosa, D. Ponce, M. A. Pozo, and J. Puerto.  Network flow based approaches for the pipelines routing problem in naval design.  _Omega_ , 111:102659, 2022.
  * Blanco et al. (2023) V. Blanco, G. González, Y. Hinojosa, D. Ponce, M. A. Pozo, and J. Puerto.  The pipelines and cable trays location problem in naval design.  _Ocean Engineering_ , 286:115525, 2023.
  * Cambero and Sowlati (2014) C. Cambero and T. Sowlati.  Assessment and optimization of forest biomass supply chains from economic, social and environmental perspectives – a review of literature.  _Renewable and Sustainable Energy Reviews_ , 36:62–73, 2014\.
  * De Meyer et al. (2015) A. De Meyer, D. Cattrysse, and J. Van Orshoven.  A generic mathematical model to optimise strategic and tactical decisions in biomass-based supply chains (optimass).  _European Journal of Operational Research_ , 245(1):247–264, 2015.
  * Delzeit and Kellner (2013) R. Delzeit and U. Kellner.  The impact of plant size and location on profitability of biogas plants in Germany under consideration of processing digestates.  _Biomass and Bioenergy_ , 52:43–53, May 2013.
  * Díaz-Trujillo and Nápoles-Rivera (2019) L. A. Díaz-Trujillo and F. Nápoles-Rivera.  Optimization of biogas supply chain in Mexico considering economic and environmental aspects.  _Renewable Energy_ , 139:1227–1240, Aug. 2019.  ISSN 09601481.
  * Egieya et al. (2020) J. M. Egieya, L. Čuček, K. Zirngast, A. J. Isafiade, and Z. Kravanja.  Optimization of biogas supply networks considering multiple objectives and auction trading prices of electricity.  _BMC Chemical Engineering_ , 2(1):3, 2020.
  * Fetting (2020) C. Fetting.  The european green deal.  _ESDN report_ , 53, 2020.
  * Ghaderi et al. (2016) H. Ghaderi, M. S. Pishvaee, and A. Moini.  Biomass supply chain network design: An optimization-oriented review and analysis.  _Industrial Crops and Products_ , 94:972–1000, 2016.
  * Hernandez and Martin (2017) B. Hernandez and M. Martin.  Optimal integrated plant for production of biodiesel from waste.  _ACS Sustainable Chemistry & Engineering_, 5(8):6756–6767, 2017.
  * Hinojosa et al. (2000) Y. Hinojosa, J. Puerto, and F. R. Fernández.  A multiperiod two-echelon multicommodity capacitated plant location problem.  _European Journal of Operational Research_ , 123(2):271–291, 2000.
  * Hinojosa et al. (2008) Y. Hinojosa, J. Kalcsics, S. Nickel, J. Puerto, and S. Velten.  Dynamic supply chain design with inventory.  _Computers & operations research_, 35(2):373–391, 2008.
  * Hu et al. (2018) Y. Hu, M. Scarborough, H. Aguirre-Villegas, R. A. Larson, D. R. Noguera, and V. M. Zavala.  A supply chain framework for the analysis of the recovery of biogas and fatty acids from organic waste.  _ACS Sustainable Chemistry & Engineering_, 6(5):6211–6222, 2018.
  * Igliński et al. (2012) B. Igliński, R. Buczkowski, A. Iglińska, M. Cichosz, G. Piechota, and W. Kujawski.  Agricultural biogas plants in Poland: Investment process, economical and environmental aspects, biogas potential.  _Renewable and Sustainable Energy Reviews_ , 16(7):4890–4900, Sept. 2012.
  * International Energy Agency (2021) International Energy Agency.  Net zero by 2050: A roadmap for the global energy sector., 2021.  URL https://www.iea.org/reports/net-zero-by-2050.
  * Jafaru et al. (2018) E. Jafaru, Čuček Lidija, Z. Klavdija, I. Adeniyi, P. Bojan, and K. Zdravko.  Biogas supply chain optimization considering different multi-period scenarios.  _Chemical Engineering Transactions_ , 70:985–990, 2018\.
  * Jensen et al. (2017) I. G. Jensen, M. Münster, and D. Pisinger.  Optimizing the supply chain of biomass and biogas for a single plant considering mass and energy losses.  _European Journal of Operational Research_ , 262:744–758, 2017.
  * Kwaśny and Balcerzak (2017) J. Kwaśny and W. Balcerzak.  Production logistics and participation of biogas in obtaining primary energy in poland.  _Energy & Environment_, 28(4):425–436, 2017\.
  * Legg (2021) S. Legg.  Ipcc, 2021: Climate change 2021-the physical science basis.  _Interaction_ , 49(4):44–45, 2021.
  * Ma and Zavala (2022) J. Ma and V. M. Zavala.  Solution of large-scale supply chain models using graph sampling & coarsening.  _Computers & Chemical Engineering_, 163:107832, 2022.  ISSN 0098-1354.
  * Martínez and Fransoo (2017) J. C. V. Martínez and J. C. Fransoo.  Green facility location.  _Sustainable supply chains: a research-based textbook on operations and strategy_ , pages 219–234, 2017.
  * Mayerle and Neiva de Figueiredo (2016) S. F. Mayerle and J. Neiva de Figueiredo.  Designing optimal supply chains for anaerobic bio-digestion/energy generation complexes with distributed small farm feedstock sourcing.  _Renewable Energy_ , 90:46–54, 2016.
  * Muradin and Foltynowicz (2018) M. Muradin and Z. Foltynowicz.  Logistic aspects of the ecological impact indicators of an agricultural biogas plant.  _LogForum_ , 14(4), 2018.
  * National Sustainable Agriculture Coalition (2023) National Sustainable Agriculture Coalition.  What is the farm bill?, 2023.  URL https://sustainableagriculture.net/our-work/campaigns/fbcampaign/what-is-the-farm-bill/.
  * Park et al. (2019) Y. S. Park, J. Szmerekovsky, and A. Dybing.  Optimal location of biogas plants in supply chains under carbon effects: Insight from a case study on animal manure in north dakota.  _Journal of Advanced Transportation_ , 2019.
  * Sampat et al. (2019) A. M. Sampat, Y. Hu, M. Sharara, H. Aguirre-Villegas, G. Ruiz-Mercado, R. A. Larson, and V. M. Zavala.  Coordinated management of organic waste and derived products.  _Computers & Chemical Engineering_, 128:352–363, 2019\.
  * Sampat et al. (2021) A. M. Sampat, A. Hicks, G. J. Ruiz-Mercado, and V. M. Zavala.  Valuing economic impact reductions of nutrient pollution from livestock waste.  _Resources, Conservation and Recycling_ , 164:105199, 2021\.
  * Sarker et al. (2018) B. R. Sarker, B. Wu, and K. P. Paudel.  Optimal number and location of storage hubs and biogas production reactors in farmlands with allocation of multiple feedstocks.  _Applied Mathematical Modelling_ , 55:447–465, Mar. 2018\.
  * Sarker et al. (2019) B. R. Sarker, B. Wu, and K. P. Paudel.  Modeling and optimization of a supply chain of renewable biomass and biogas: Processing plant location.  _Applied Energy_ , 239:343–355, Apr. 2019.
  * Scarlat et al. (2018) N. Scarlat, F. Fahl, J.-F. Dallemand, F. Monforti, and V. Motola.  A spatial analysis of biogas potential from manure in europe.  _Renewable and Sustainable Energy Reviews_ , 94:915–930, 2018.
  * Silva et al. (2017) S. Silva, L. Alçada-Almeida, and L. C. Dias.  Multiobjective programming for sizing and locating biogas plants: A model and an application in a region of Portugal.  _Computers & Operations Research_, 83:189–198, July 2017\.  ISSN 03050548.
  * Soha et al. (2021) T. Soha, L. Papp, C. Csontos, and B. Munkácsy.  The importance of high crop residue demand on biogas plant site selection, scaling and feedstock allocation – A regional scale concept in a Hungarian study area.  _Renewable and Sustainable Energy Reviews_ , 141:110822, May 2021.
  * Sultana and Kumar (2012) A. Sultana and A. Kumar.  Optimal siting and size of bioenergy facilities using geographic information system.  _Applied Energy_ , 94:192–201, 2012.
  * Tampio et al. (2022) E. Tampio, F. Pettersson, S. Rasi, and M. Tuomaala.  Application of mathematical optimization to exploit regional nutrient recycling potential of biogas plant digestate.  _Waste Management_ , 149:105–113, 2022.
  * Tominac et al. (2022) P. A. Tominac, W. Zhang, and V. M. Zavala.  Spatio-temporal economic properties of multi-product supply chains.  _Computers & Chemical Engineering_, 159:107666, 2022.
  * United Nations (2015) United Nations.  Transforming our world: The 2030 agenda for sustainable development., 2015\.  URL https://sdgs.un.org/goals/goal7.
  * Valenti et al. (2018) F. Valenti, S. M. Porto, B. E. Dale, and W. Liao.  Spatial analysis of feedstock supply and logistics to establish regional biogas power generation: A case study in the region of sicily.  _Renewable and Sustainable Energy Reviews_ , 97:50–63, 2018\.