Optimization of sawmill residues collection for bioenergy production

Optimization of sawmill residues collection for bioenergy production

Applied Energy 202 (2017) 487–495 Contents lists available at ScienceDirect Applied Energy journal homepage: www.elsevier.com/locate/apenergy Optim...

1MB Sizes 0 Downloads 54 Views

Applied Energy 202 (2017) 487–495

Contents lists available at ScienceDirect

Applied Energy journal homepage: www.elsevier.com/locate/apenergy

Optimization of sawmill residues collection for bioenergy production David S. Zamar a,⇑, Bhushan Gopaluni a, Shahab Sokhansanj a,b a b

Department of Chemical and Biological Engineering, University of British Columbia, Vancouver, BC V6T 1Z3, Canada Resource and Engineering Systems Group, Environmental Sciences Division, Oak Ridge National Laboratory, Oak Ridge, TN 37831-6422, United States

h i g h l i g h t s  A real-life stochastic vehicle routing problem is presented and solved.  The goal is to optimize the collection of biomass residues from a set of sawmills.  A model that optimizes the energy returned on energy invested (EROEI) is proposed.  Accounting for the randomness in the input parameters makes the model more robust.  Our solution is shown to be more accurate and precise than a deterministic one.

a r t i c l e

i n f o

Article history: Received 1 March 2017 Received in revised form 5 May 2017 Accepted 22 May 2017

Keywords: Renewable energy systems Routing algorithms Robust estimation Uncertainty Stochastic approximation Logistics modeling

a b s t r a c t The collection of sawmill residues is an important logistic activity for the pulp and paper industry, which uses the biomass as a source of energy. We study a vehicle routing problem for a network composed of a single depot and several sawmills. The sawmills serve as potential suppliers of biomass residues to the depot, which in turn processes and distributes the residues to the pulp and paper mills. This problem consists of identifying the best daily routing schedule for a fixed number of trucks. The objective is to maximize the ratio of energy returned on energy invested, while satisfying a minimum daily amount of dried biomass residues. There are several random components in the problem, including the availability and moisture content of the biomass residues. We use a combination of scenario analysis and heuristics to solve this stochastic vehicle routing problem. A performance comparison of the proposed method reveals an estimated daily energy savings of 6 GJ over the benchmark method. Ó 2017 Elsevier Ltd. All rights reserved.

1. Introduction The deleterious impacts of climate change coupled with the ongoing urbanization of countries around the world have led to a global effort to reduce greenhouse gas emissions for energy production, especially in the transportation sector. The transport sector is responsible for approximately one quarter of greenhouse gas emissions in both Europe and America, making it the second largest emitting sector after energy [1,2]. A plethora of studies show that urban freight transport could be vastly more efficient. According to the European Commission, 24% of commercial trucks that operate in Europe are empty [2]. In the United States, commercial trucks drive an estimated 19 billion needless miles each year [3]. Thus, significant economic and environmental savings may be achieved by reducing truck transport.

⇑ Corresponding author. E-mail addresses: [email protected] (D.S. Zamar), bhushan.gopaluni@ ubc.ca (B. Gopaluni), [email protected] (S. Sokhansanj). http://dx.doi.org/10.1016/j.apenergy.2017.05.156 0306-2619/Ó 2017 Elsevier Ltd. All rights reserved.

This paper presents a modeling framework to optimize the collection of biomass residues from sawmills in the presence of uncertainty. Recognizing that the goal of procuring the biomass is for the production of energy, we cast the problem as an energy optimization task. Biomass feedstocks that are brought to a combined heat and power (CHP) plant for conversion into energy must meet certain specifications in order for the process to be safe and successful. The wet basis moisture content is used to describe the water content of biomass and is defined as the percentage equivalent of the ratio of the weight of water to the total weight of the biomass. The most important property of biomass feedstocks with regard to transport and energy conversion is its moisture content. This is for two reasons. First, the biomass feedstock might have to undergo drying if the moisture content exceeds that accepted by the facilities boiler, with an ensuing energy expense. Second, the unnecessary transportation of water embedded in the biomass feedstock increases the energy utilized for transportation. Hence, our objective is to schedule the collection activities and identify

488

D.S. Zamar et al. / Applied Energy 202 (2017) 487–495

collection routes that maximize the total energy returned on energy invested (EROEI) taking into consideration the transportation and drying of the biomass feedstock. The EROEI is defined as the ratio of the amount of usable energy acquired from a particular energy resource to the amount of energy expended to obtain that energy resource. Many processes in the biomass supply chain may randomly vary each time they are performed. For example, the availability and quality of the biomass residues as well as the time required to load, unload and transport will vary on a case by case basis. These disturbances can be caused by spatial and temporal variations in a seemingly unpredictable fashion that can be best represented in the model with the use of random variables. A main source of uncertainty driving these random variations can be attributed to weather, such as precipitation, wind, and exposure to the environment. In the context of biomass procurement for bioenergy, variability in the biomass quality (i.e. moisture content, calorific value, etc.) and availability can have a significant impact on the cost of energy production and the energy conversion efficiency. As such, the modeling of these activities should consider both the biomass quality and availability as random parameters with probability distributions. Modeling such variations is often done using specified stochastic distributions with mean and variance estimated from historical data. The sampled values of the random variables were obtained using random number generators. All other model input parameters (i.e. number and location of sawmills, travel distances, number of trucks, truck capacities, etc.) are considered fixed. The values assigned to the fixed parameters were established per the study criteria or derived from the literature. Our motivating case study deals with the planning of truck routes for the collection of sawmill residues (or waste) in the Lower Mainland region of British Columbia, Canada [4]. A symbiotic relationship exists between sawmills and pulp mills. Approximately 50% of a log, by volume, gets turned into lumber at a sawmill. The waste residues from the process are utilized by the pulp mills to produce both pulp and excess green energy. The majority of pulp mills rely on purchased sawmill residue chips for most, if not all, of their chip supply. Consequently, the sale of residue chips has become an essential revenue stream for the sawmills. The pulp mills have the necessary expertise, infrastructure and potential to be future large scale producers of biomass based transportation fuels [5]. The real-life sawmill residues collection problem under consideration may be described as follows. There are a total of 25 sawmills in the region and a single depot. The location of the sawmills and the depot are given along the streets of a defined road network. The biomass residues produced by the sawmills must be collected by a fleet of trucks with known capacities. The average daily amount of biomass residues produced by each sawmill are subject to variability. Each truck may collect biomass residues from several sawmills before returning to the depot to unload. Truck drivers work 8 h shifts that start at 9 am and end at 5 pm. The trucks leave the depot at the start of the day at 9 am and are allowed to make several, potentially different, routes in a single day. A truck must return to the depot to unload after completing a route. In addition, all trucks must return to the depot before the end of the day at 5 pm. The amount of biomass residues that should be collected on a daily basis is determined by the energy demand from the pulp mills that are being supplied by the depot. Information regarding the mass, measured in green tonnes (gt), and energy content, measured in gigajoules per tonne (GJ t1 ), of the biomass residues available at each sawmill are unknown and highly variable. The energy content of the biomass residues

depends on the moisture content, which is known to play a significant role in the optimization of biomass supply chains [6]. The biomass residues collection problem can be viewed as a stochastic vehicle routing problem (SVRP). The vehicle routing problem (VRP) is a very well known and widely studied problem in combinatorial optimization. The general objective is to route the vehicles, with each route starting and ending at the depot, so that all customer supply demands are met and the total travel distance is minimized. As this is a computationally very hard problem, which cannot be solved by exact methods, in practice heuristics are typically used for this purpose [7]. The SVRP arises when some of the elements of the VRP are not known exactly, such as the travel times, product availability and customer demand. Prevalent techniques used to solve the SVRP minimize the expected cost of the assigned routes, including those incurred by recourse actions, while satisfying a set of chance constraints. In most cases, simplifying assumptions are used in order to fit probabilistic models that rely on mathematical theorems and lemmas [8]. In our experience these techniques generally result in computationally demanding and unscalable models [9]. Moreover, the explicit use of the expected value in solving the SVRP implies that the decision maker is primarily interested in the average performance of the solution and is not concerned with its variance and other features of its random behavior. To avoid these drawbacks, the SVRP presented in this paper is solved using the quantilebased scenario analysis (QSA) approach [9]. The QSA approach analyzes the performance of solutions obtained from solving deterministic realizations of the stochastic problem and identifies the solution that optimizes chosen quantiles of the stochastic objective function, subject to satisfying conditions on given quantiles of the constraint distribution. An advantage of this approach is that it requires only that each deterministic version (i.e. scenario) of the stochastic problem be solvable. The contributions of this paper are threefold. First, it demonstrates the methodological and practical advantages of addressing the problem of procurement planning for bioenergy production in a stochastic manner in comparison to a purely deterministic approach. Second, it is shown that the procurement of biomass for energy production can be presented in a mathematical framework that accounts for both the energy utilized for procurement and the energy obtained during conversion. Third, it illustrates how randomness can be directly incorporated into any existing deterministic VRP, often encountered in energy production, and solved as an SVRP by a direct application of the QSA approach. The remainder of this paper is organized as follows. A review of the literature on SVRP studies is discussed in Section 2. The sawmill residue collection model and its input requirements are presented in Section 3.1. The heuristic routing and scheduling methodologies are explained in Section 3.2, followed by the presentation and discussion of results in Section 4. Main conclusions are provided in Section 5. 2. Literature review The deterministic VRP has received some attention in recent years in the area of energy efficiency, production and management. For example, Araujo et al. [10] studied the collection of waste frying oil for the production of biodiesel in Rio de Janeiro, Brazil. Juan et al. [11] investigated several open research problems related to the introduction of electric vehicles in logistics and transportation, including strategic and planning challenges associated with hydrogen-based electric vehicles. Shao et al. [12] studied an electric vehicle routing problem in Beijing to address some operational issues such as range limitation and charging demand. Jokinen et al. [13] present a mathematical model that minimizes the total costs

489

D.S. Zamar et al. / Applied Energy 202 (2017) 487–495

associated with transporting liquefied natural gas. Yagcitekin & Uzunoglu [14] present a smart charging management algorithm strategy to successfully route electric vehicles to the most suitable charging point, thereby decreasing the charging costs and prevent the overloading of transformers. These studies, however, do not address parameter uncertainty in their respective optimization problems. While stochastic optimization has been used in the context of energy, attention has been primarily focused on energy storage and flow characteristics, see for example the works of Moura et al. [15] and Kou et al. [16]. On the other hand, despite its importance and frequent appearance in energy related problems, SVRP has been mostly absent from this literature. Approaches for solving the SVRP may be found throughout the literature on research and development (R&D). Nevertheless, a recent review of this area reveals the lack of a comprehensive methodology that can be applied to the general SVRP [17]. In most instances surveyed, the SVRP is transformed into a deterministic equivalent that can be solved using mathematical programming solvers [18,19]. The predominant approach for solving the SVRP is to use ‘‘here-and-now-optimization” techniques, where decisions are made based on information available at the time they are taken. When some of the parameters are random, it is no longer possible to require that all constraints be satisfied for all realizations of the random parameters. Therefore, the decision maker may require that the constraints be satisfied with a given probability [20–22] or the incorporation of corrective actions be taken when a constraint is violated [23,24]. Published studies that investigate different variants of the SVRP include urban public bus transport systems [25], plant to customer

P max Jðx; dÞ ¼

P

ði;jÞ2As

P k2K

l2Rk djkl

h



how much the vehicle weighs on its own, without any cargo or passengers. Thus, for a given truck k and route l,

X

uijkl ¼

dykl þ fk ;

ð1Þ

yj y;j2V s

where y  j if and only if sawmill y is visited before sawmill j. Let h be the required moisture content of the biomass residues for conversion into heat and power at the pulp mills. Let g represent the calorific value of the biomass at a moisture content of h, and T be the temperature at which it is stored. Let k represent the energy required to transport a green tonne of biomass a distance of one kilometer (km). The time required to load a truck at node j 2 V s is given by aj þ bj djkl , where aj is the truck setup time at node j. Similarly, the time required to unload a truck at the depot is given by a0 þ b0 ck . The parameters tij ; aj ; bj ; wj and xj ; 8ði; jÞ 2 A, are unknown quantities, which are assumed to be normally distributed. Thus, t ij ; aj ; bj ; xj and wj are stochastic variables. A summary of the model parameters is provided in Table 1. For a specific realization of the random variables t ij ; aj ; bj ; xj , and wj ; 8ði; jÞ 2 A, the deterministic problem can be modeled using the following additional variables. Let uk ¼ 1; k 2 K, if and only if the driver of truck k has had a lunch break and uk ¼ 0 otherwise. Let c be the time alloted for a lunch break. Define the set of arcs between nodes that do not terminate at the depot as As ¼ fði; jÞ j i 2 V; j 2 V s ; i – jg. Let xijkl ¼ 1 if and only if truck k 2 K uses arc ði; jÞ on its lth assigned route and xijkl ¼ 0 otherwise. The optimization model for the SVRP is given in Eq. (2).



i

xj Þ g  xj  hð1 g þ 4:19  103 ð100  TÞ þ 2:26 xijkl  g 1h P P P ði;jÞ2A k2K l2Rk xijkl  dij  uijkl  k

distribution [26], waste collection [27,28,8], collection and delivery of goods [29], in-field scheduling for machinery fleets in biomass operations [30,31], and public transport in adverse weather conditions [32]. However, SVRP studies that analyze the collection of biomass feedstock in an urban setting or that focus on the optimization of energy production while minimizing energy utilization during the procurement and preprocessing stages remain mostly untouched areas in the R&D literature. 3. Methods 3.1. Optimization model In this section we formally define the SVRP model for this study. The problem is defined on a graph G ¼ fV; Ag, where the set of nodes V ¼ V d [ V s consists of a single depot V d ¼ f0g and m sawmills, V s ¼ f1; . . . ; mg, and A ¼ fði; jÞ j i; j 2 V; i – jg is a set of arcs. Let K ¼ f1; . . . ; ng be the set of trucks with payload capacities ck ; k 2 K: Let dij and t ij be the travel distance and time associated with arc ði; jÞ, respectively. In addition, each node, j 2 V s , has an inventory wj of biomass residues with an average moisture content xj and a corresponding mass flow rate bj for loading a truck. Define the set Rk ¼ f1; . . . ; rk g; k 2 K, where rk denotes the number of routes assigned to truck k. Let djkl represent the mass of the biomass residues picked up at node j 2 V s by truck k 2 K on its l 2 Rk assigned route. Let uijkl represent the sum of both the payload and curb weight, fk , of truck k when traveling from node i to node j on its lth assigned route. The curb weight is simply

s:t:

XX

x0jkl ¼ 1;

ð2aÞ

8k 2 K

ð2bÞ

j2V l2Rk

XX xi0kl ¼ 1; i2V l2Rk

XX

xijkl ¼

i2V l2Rk

XX

j  ck 6

xjikl ;

ð2cÞ

8j 2 V; k 2 K

ð2dÞ

i2V l2Rk

XX djkl 6 wj ; k2K l2Rk

8k 2 K

X

8j 2 V s

djkl 6 ck ;

ð2eÞ

8k 2 K; l 2 Rk

ð2fÞ

j2V s

8j 2 V s ; k 2 K; l 2 Rk 8k 2 K; l 2 Rk 2 f0; 1g; 8ði; jÞ 2 A; k 2 K

djkl P 0;

ð2gÞ

d0kl ¼ 0;

ð2hÞ ð2iÞ

xijkl XXX   djkl 1  xj þ h P d

ð2jÞ

j2V s k2K l2Rk

X ukl ¼ 1;

8k 2 K

ð2kÞ

l2Rk

XX xijkl ðtij þ aj þ bj djkl Þ

ði;jÞ2As l2Rk

þ

XX xi0kl ðt i0 þ a0 þ b0 ck Þ þ c 6 s;

8k 2 K

ð2lÞ

i2V l2Rk

To simplify the notation, the decision variables are compactly expressed as the route schedule, x, and the amount of residues, d, to pick up from each sawmill visited in each route. We assume that the energy obtained from a dry tonne of biomass residues at a moisture content of h is constant. The goal is to maximize the

490

D.S. Zamar et al. / Applied Energy 202 (2017) 487–495

Table 1 Optimization model notation. Symbol

Description

Units

Set notation/indices

V0: Vs: A: K: Rk

depot node sawmill nodes arcs connecting the nodes trucks routes assigned to truck k

– – – – –

V 0 ¼ f0g V s ¼ f1; 2; . . . ; mg A ¼ fði; jÞ j i; j 2 V s [ V 0 ; i – jg K ¼ f1; 2; . . . ; ng Rk ¼ f1; 2; . . . ; r k g; k 2 K

ck : fk : j: c: s: h wj : xj : k: g:

truck payload capacity truck curb weight minimum payload fill percentage per route lunch break duration maximum daily operating hours for a truck driver required moisture content of dry biomass inventory at node j wet basis moisture content of residues at sawmill j transport energy net calorific value of biomass at a moisture content of h setup time for loading a truck at sawmill j mass flow rate for loading a truck at sawmill j

t t % h h % gt %

k2K k2K – – –

Sets

Parameters

aj bj

a0

GJ t1 GJ t1 h gt h h

b0

setup time for unloading a truck at the depot mass flow rate for unloading a truck at the depot

dij tij uk :

travel distance between nodes i and j travel time between nodes i and j logical value indicating if truck k had lunch

gt h km h –

rk : xijkl : dikl :

number of routes assigned to truck k logical value indicating if truck k uses arc ði; jÞ on its lth assigned route amount of residues picked up by truck k at node i on its lth assigned route

– – gt

1

1

j 2 Vs j 2 Vs – – j 2 Vs j 2 Vs – – i; j 2 V 0 [ V s i; j 2 V 0 [ V s k2K

Variables

EROEI, denoted by the function Jðx; dÞ in Eq. (2a), which is the ratio of the amount of usable energy obtained from the biomass to the amount of energy used in transporting it. As the fuel economy is dependent upon the load of the truck [33], each kilometer traveled is multiplied (weighted) by the load transported. Consequently, the EROEI described in Eq. (2a) is a unitless ratio. The energy required for drying the biomass is also reflected in Jðx; dÞ. This is achieved by subtracting, from the numerator, the quantity of dried biomass that would be required to generate the energy needed to dry the delivered residues to the desired moisture content, h. Lemma 1 shows how to calculate this quantity. The optimization is done subject to the following constraints. Eqs. (2b) and (2c) ensure that all k trucks must leave and return to the depot at the end of each route. Eq. (2d) ensures that the inflow and outflow must be equal except for the depot node. Eq. (2e) ensures that the trucks cannot pick up more residues than are available at each sawmill node. Eq. (2f) ensures that each truck must come back at least j  100% full at the end of each route, where 0 < j 6 1. Eq. (2f) also ensures that the truck payload capacity ck is never exceeded. Eqs. (2g) and (2h) enforce nonnegativity and that the trucks are empty when leaving the depot at the start of each route. Eq. (2i) enforces binary variables. Eq. (2j) ensures that the minimum daily amount of required dry tonnes of sawmill residues, d , at a moisture content of h is met. A dry tonne refers the mass of the biomass after it is dried to a relatively low, consistent moisture level. Eq. (2k) ensures that each truck driver takes a lunch break. Finally, Eq. (2l) ensures that each truck operates a maximum of s hours. It is convenient to represent the solution to this problem as a set of collection runs. The lth collection run represents the ensemble of routes assigned to each truck on their lth route. As such, solutions to the SVRP are presented in this form. 3.2. Simulation model Scenario analysis approaches to solving complex stochastic problems has become more prominent in the literature Sharma

k2K k 2 K; i; j 2 V 0 [ V s ; l 2 Rk k 2 K; i 2 V s ; l 2 Rk

et al. [34]. We solve the SVRP presented earlier using the quantile-based scenario analysis (QSA) approach described in Zamar et al. [9]. A discrete event simulation model based on the heuristic illustrated in Fig. 1 was implemented to solve each scenario. To make the problem more manageable, we first cluster the sawmills based on their distance matrix and restrict attention to routes within each cluster. Thus, our approach does not consider routes that span across clusters, which are unlikely to be selected due to their Hamiltonian distance. The Hamiltonian distance of a set of nodes is the shortest path, made up from the edges connecting the nodes, which passes through each node. The QSA method samples scenarios from their underlying distribution and solves each scenario separately. The solution of each deterministic scenario problem is evaluated across the sampled scenarios to provide an estimate of the objective and constraint satisfaction distribution of each solution. Solutions are then ranked based on their corresponding objective and constraint distribution quantiles. For this application we maximize the 0.5 quantile (median) of the stochastic objective function subject to satisfying a quantile constraint on the procurement distribution. The objective distribution corresponds to the EROEI while the constraint distribution represents the fulfilled portion of the demand as described in Eqs. (2a) and (2j), respectively. We restrict attention to solutions that have a 90% probability of satisfying at least 90% of the demand across scenarios and apply our modeling framework to the previously described SVRP and solve it using the QSA approach. A map of the study region is shown in Fig. 2, which identifies the locations of 25 sawmills and a single depot. The depot requires a minimum of 180 dry tonnes of biomass at a moisture content of h ¼ 0:3  100%. The net calorific value of one dry tonne of biomass, g, is assumed to be fixed at 12.2 GJ t1 [35]. Similarly, the energy spent for transport is fixed at 0.0031 GJ per tonne-km [36]. The temperature, T, at which the biomass is stored at the depot is assumed to be fixed at 20 C. There are a limited number of identical trucks with a curb weight of fk ¼ 13:5 tonnes and a payload

491

D.S. Zamar et al. / Applied Energy 202 (2017) 487–495

Fig. 1. The logic flowchart for truck operations in the simulation model.

S5 S10

S6 S4 S11

S13 S3

S2

S9

S23

S7 S24

S12

S1

S8

S17

S20

S25 S21

S14

S19

S15 S22

DEP

S18

S16

Fig. 2. The study region comprises 25 sawmills (S1–S25) across the Lower Mainland of BC, and a single depot (DEP) located in Mission. The sawmills have been spatially clustered into four distinct color-coded groups. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

capacity of ck ¼ 30 tonnes, k 2 K, that are used to collect residues in the considered region. The random parameters that distinguish the scenarios are categorized as follows: (1) the biomass residues available at each sawmill, wj ; (2) the moisture content of the residue available at each

sawmill, xj ; (3) the time required to travel between the nodes (sawmills and the depot), t ij ; (4) the truck setup time at the depot and the sawmills, a0 and aj ; and (5) the load and unload flow rates at the sawmills and the depot, b0 and bj . The random parameters wj and xj are modeled based on published summary data [37] on

492

D.S. Zamar et al. / Applied Energy 202 (2017) 487–495

Table 2 Stochastic input data. Average values are shown for parameters with subscripts.

Total available residue (gt) Residue per Sawmill (gt)

Parameter P wj wj

Mean

Standard deviation

1010.2 40.4

31.6 10.3 3.8 0.01 0.01 0.001 0.002 0.6

Moisture content (%) Unload setup time (h) Load setup time (h) Unload time rate (h/gt) Load time rate (h/gt)

xj a0 aj bj

45.0 0.1 0.1 0.01 0.02

Travel time (h)

t ij

1.5

b0

biomass residue availability and moisture content at each of the sawmills. Established conversion factors were used to convert from wet to dry weight and energy content [38]. The random parameters aj and bj ; i 2 V, describing the truck setup, load and unload times, are modeled based on the values calculated by FPInnovations in their assessment of economically viable biomass [39]. The expected travel time, t ij , and the exact travel distance, dij ; 8ði; jÞ 2 A, between nodes were calculated using the RgoogleMaps package in R [40]. A summary of the distributions of the random parameters included in the model are shown in Table 2. The average available biomass residue at each sawmill was modeled as a multivariate normal distribution based on their average annual production. Spatial variogram models were fit using the sp package in R to represent the wet basis moisture content of the biomass residue available at each sawmill [41]. A variogram is a mathematical description of the spatial continuity of the data and is calculated using a measure of variability between pairs of points at various distances [42]. The travel time between a pair of nodes was modeled as an exponential random variable with a mean equal to the expected travel time. The mean and standard deviation of the travel time, averaged across all pairs of nodes, are included in Table 2. Similarly, the mean and standard deviation of the biomass residues available at each sawmill, averaged across all the sawmills, are provided in Table 2. For completeness, the mean and standard deviation of the total available biomass residues across all the sawmills are also given in Table 2. The number of sawmills, m ¼ 25, the number of available trucks, n ¼ 3, and the distance between nodes, dij ; 8ði; jÞ 2 A, are all known and assumed to be fixed. The number of trucks needed was established by adding one truck at a time and repeating the simulation until the required dry tonnes of daily biomass residues is achieved. The simulation model was implemented using the R system for statistical computing [43]. A total of 1000 scenarios were simulated by sampling from the appropriate distribution of each random parameter included in the model. The parameters were considered to be statistically independent. To control the amount of biomass procured on a daily basis, we restrict attention to solutions that have a 90% chance of satisfying 90% of the daily demand of 180 dry tonnes of biomass. Moreover we only consider solutions with a probability of 10% of exceeding the demand. This will compel the QSA solution to seldom exceed the required daily demand, but at the same time consistently meet at least 90% of the demand. Exceeding the daily demand cannot always be avoided as the trucks are constrained to return to the depot at least 70% full after each route in order to improve the transport efficiency. Thus,

j ¼ 0:70 in Eq. (2f). Moreover, as the availability of residues at each of the sawmills is unknown, the trucks gather as much residues as they can load at each sawmill they visit. The sawmills were clustered into four groups based on their travel distance matrix using the k-medoids method implemented in the package cluster [44]. The number of clusters was estimated by the method of optimum average silhouette width [45]. The four clusters are depicted in Fig. 2 and the cluster sizes are 3, 5, 7 and 10, respectively. Given the time constraints imposed by the sawmills and the truck drivers operating hours, each truck was able to perform a maximum of 3 routes per day. This fact was learned from the simulation results. The mean scenario (MS) approach is an alternative method that may serve as a benchmark for assessing the performance of the QSA method. The MS solution is obtained by substituting the random parameters /j ; xj ; a0 ; aj ; b0 ; bj , and tij in Eq. (2) by their corresponding expected values and solving the resulting deterministic constrained optimization problem. This is equivalent to solving the mean scenario, hence the name of the approach. While the MS approach is helpful for making decisions, its applicability to real world problems is sometimes limited due to its inability to account for uncertainties in the model parameters. Synthetic data are commonly generated in order to validate mathematical models. Model feature calibration for both the QSA and MS methods was done using a synthetic training dataset comprised of 1000 scenarios (i.e. realizations of the stochastic model parameters). A separate and independent synthetic dataset of the same size as the training dataset was used to perform an out-ofsample validation of the QSA and MS solutions. 4. Results and discussion The optimal routes selected for each truck by the QSA method are shown in Table 3. Each row in Table 3 corresponds to a collection run. From the first row, we see that the first collection run consists of sending the first truck to sawmill S22, the second truck to sawmills S15 and S18, and the third truck to sawmills S20, S23 and S25, in this order. We denote this route combination as [S22, (S15, S18), (S20,S23,S25)]. The total distance traveled in the first collection run is 52.0 km and an accumulated time of 5.6 h is required on average. The average amount of biomass obtained in the first collection run is 54:9 dry tonnes. Furthermore, the first collection run yields a performance ratio (EROEI) of 96:7. Recall, that the EROEI has been adjusted for the energy spent for transportation and drying of the biomass residues. Similarly, the performance of collection runs 2, and 3 are shown in the corresponding rows of Table 3. Notice the drop in the EROEI between collection runs. This is attributed to the fact that the sawmills closest to the depot are the first to become depleted, compelling the trucks to travel longer distances to collect the biomass residues. Sawmills S22, S2 and S25 are repeatedly visited as they are relatively close to the depot and have a large production of residues with a low moisture content. The average accumulated operating hours for the collection runs are shown in the last column of Table 3. The mean of this column gives the average operating time in hours for each of the three trucks (6.5 h). Evidently, two trucks would not be sufficient to collect the required daily amount of residues in an 8-h work day and employing four trucks would be wasteful as they would be

Table 3 QSA optimal route schedule. Truck 1

Truck 2

Truck 3

Distance (km)

Weight (t)

EROEI

Time (h)

S22 S14,S19 S2

S15,S18 S22,S19,S25 S17,S25,S16

S20,S23,S25 S2 S2,S24

52.0 118.3 218.3

54.9 57.1 59.8

96.7 43.0 22.3

5.6 6.5 7.5

493

1.00

1.00

0.75

0.75

Cumulative probability

Cumulative probability

D.S. Zamar et al. / Applied Energy 202 (2017) 487–495

0.50

0.25

Method MS QSA

0.50

0.25

0.00

0.00 40

45

130

50

150

EROEI

170

190

210

Procurement (dry tonnes)

Fig. 3. Performance distributions of the QSA and MS solutions. (a) Cumulative distribution function of the energy returned on energy invested (EROEI). The QSA solution is more efficient as it has a consistently larger EROEI ratio than that of the MS solution. (b) Cumulative distribution function of the total amount of biomass procured. The MS solution consistently delivers more dry tonnes of biomass than required thereby producing an undesirable surplus at the depot. The horizontal dotted line indicates where 90% of the demand (i.e. 162 dry tonnes) is achieved.

underutilized. Therefore, in this case, using three trucks may be considered optimal. The QSA solution collects an average of 171:8 dry tonnes of biomass residues per day, and achieves a median EROEI of 47:3. A computation time of 3.91 min was required by the QSA method to solve this problem running on an Intel Xeon 3.6 GHz processor. A comparison of the cumulative distribution function for the EROEI and amount of biomass residues collected by the QSA and MS solutions are shown in Fig. 3a and b, respectively. As can be

seen in Fig. 3a and b, the QSA solution is able to consistently obtain a greater EROEI than that of the MS solution while meeting 90% of the demand (marked by the vertical dotted line at 162 tonnes in Fig. 3b) at least 90% of the time. To validate the QSA solution shown in Table 3, we applied it to a new set of 1000 randomly generated scenarios. The newly obtained scenario data may be considered test data because they were not originally used to train (or fit) the QSA solution. In outof-sample evaluation, the QSA solution obtained a 90% probability 0.06

0.2

Density

Density

0.04

Method MS QSA

0.1 0.02

0.0

0.00 40

45

EROEI

50

125

150

175

200

225

Procurement (dry tonnes)

Fig. 4. Out-of-sample performance of the QSA and MS solutions. (a) Probability density function of the energy returned on energy invested (EROEI). The QSA solution is more efficient as it has a consistently larger EROEI ratio than that of the MS solution. (b) Probability density function of the total amount of biomass procured. The dotted lines show where 90% and 100% of the daily demand is met. The MS solution consistently delivers a great deal more than the required 180 dry tonnes of biomass thereby producing an undesirable surplus at the depot. In contrast, the QSA solution is more accurate and precise in delivering the required amount of biomass.

494

D.S. Zamar et al. / Applied Energy 202 (2017) 487–495

of collecting at least 89:7% of the demand and a median EROEI of 46:8 with a corresponding interquartile range (IQR) of 44.6–49.0. On the other hand, the MS solution obtained a 90% probability of collecting at least 94% of the demand and has a median EROEI of 44:1 with a corresponding IQR of 42.0–46.2. A Kruskal-Wallis test [46] of the difference in medians between the QSA and MS solu-

y ¼ h; dð1  xÞ þ y

ðA:3Þ

where h is the required moisture content. From Eq. (A.3) we have:

y ¼ dhð1  xÞ þ hy

ðA:4Þ

)yð1  hÞ ¼ dhð1  xÞ

ðA:5Þ

dhð1  xÞ ; 1h

ðA:6Þ

tions EROEI values resulted in a p-value of 2:2  1016 . Thus, the median EROEI of the QSA method (46.8) is statistically significantly greater than that of the MS approach (44.1) and results in an energy savings of approximately 6 GJ per day with a corresponding 95% confidence interval of 5.96–6.11 GJ. As can be seen in Fig. 4b, the MS solution tends to over collect and has an estimated 69:6% probability of exceeding the required demand of 180 dry tonnes per day. In contrast the QSA solution collects between 162 and 180 dry tonnes 75:2% of the time (area enclosed between the dotted lines in Fig. 4b) and only exceeds the demand 13:6% of the time. This is a convenient feature of the QSA solution because the depot has a limited storage capacity and throughput.

dx  y ¼ dx 

5. Conclusion



We developed a simulation approach to solve a stochastic vehicle routing problem, where the goal is to efficiently collect biomass residues from a set of sawmills. Several factors, such as the moisture content and the amount of available biomass residues are uncertain and treated as random. The optimal routing schedule was obtained by maximizing the energy returned on energy invested (EROEI). Our approach was found to be more efficient in terms of minimizing the EROEI and more accurate and precise in terms of meeting the daily demand than by simply solving the expected or mean scenario.

where the first term in the product corresponds to the energy required to heat one tonne of water to its boiling point plus the energy to change its state all divided by the energy stored in one tonne of biomass at a moisture content of h. From Eqs. (A.7) and (A.8) we have that the net energy gained (assuming 100% burner efficiency) is equivalent to that stored in

)y ¼

thus the amount of water, in tonnes, that must be evaporated from the delivered biomass is

ð100  TÞ4:19  103 þ 2:26



 1þ

ð100  TÞ4:19  103 þ 2:26

!#

g ðA:1Þ 3

tonnes of biomass at a moisture content of h, where 4:19  10 and 2.26 represent the specific heat and latent heat of vaporization of water in GJ t1 , respectively [47]. Proof. As the wet basis moisture content, x, of biomass is the percentage equivalent of the ratio of the weight of water to the total weight of the biomass, we can write



dx ; dx þ dð1  xÞ

ð100  TÞ4:19  103 þ 2:26

g

 dx 

 hð1  xÞ ¼d 1 x 1h 1þ

g

 dhð1  xÞ 1h ðA:9Þ



ð100  TÞ4:19  103 þ 2:26

ðA:8Þ

!# ðA:10Þ

tonnes of biomass at a moisture content of h.

Lemma 1. Let d be a delivered amount of biomass in green tonnes and x be its corresponding wet basis moisture content. Let h 6 x be the required moisture content of the biomass for conversion into energy and T be the temperature at which the biomass is stored. If the net calorific value of the biomass at a moisture content of h is g GJ t1 , then the net energy gained after drying (assuming 100% burner efficiency) is equivalent to that stored in

hð1  xÞ d 1 x 1h

 dhð1  xÞ ; dx  1h

  dhð1  xÞ d  ðdx  yÞ  z ¼ d  dx  1h

We would like to acknowledge the financial support from Mitacs, NSERC, and BioFuelNet NCE (BFN).



g



Appendix A. Accounting for the energy spent in drying

ðA:7Þ

If T is the temperature (in degrees Celsius) at which the biomass is stored and g is the net calorific value (in GJ t1 ) of the biomass at a moisture content of h, then the number of tonnes of biomass that would be required to generate the energy needed to evaporate dx  y tonnes of water is given by:

Acknowledgments

"

dhð1  xÞ : 1h

ðA:2Þ

where dx is the weight of water in the delivered biomass, in tonnes. Let y be the allowed weight, in tonnes, of the dried biomass, then y satisfies

References [1] European Commission. Commission takes action to make urban travel greener, better organised and more userfriendly; 2009. [Online; accessed 04October-2015]. [2] United States Environmental Protection Agency. Sources of greenhouse gas emissions; 2015. [Online; accessed 04-October-2015]. [3] Jaffe E. How the trucking industry could be vastly more efficient; 2015. [Online; accessed 04-October-2015]. [4] Zamar DS, Gopaluni B, Sokhansanj S, Ebadian M. Economic optimization of sawmill residues collection for bioenergy conversion. IFAC-PapersOnLine 2016;49:857–62. [5] Mercer International Group. BC’s Bioenergy Strategy: Examining Progress and Potential in the Province; 2010. . [6] Sosa A, Acuna M, McDonnell K, Devlin G. Controlling moisture content and truck configurations to model and optimise biomass supply chain logistics in Ireland. Appl Energy 2015;137:338–51. [7] Nuortio T, Kytojoki J, Niska H, Braysy O. Improved route planning and scheduling of waste collection and transport. Exp Syst Appl 2006;30:223–32. [8] Buhrkal K, Larsen A, Ropke S. The waste collection vehicle routing problem with time windows in a city logistics context. Proc – Soc Behav Sci 2012;39:241–54. [9] Zamar DS, Gopaluni B, Sokhansanj S, Newlands NK. A quantile-based scenario analysis approach to biomass supply chain optimization under uncertainty. Comput Chem Eng 2017;97:114–23. [10] Araujo VKWS, Hamacher S, Scavarda LF. Economic assessment of biodiesel production from waste frying oils. Bioresour Technol 2010;101:4415–22.

D.S. Zamar et al. / Applied Energy 202 (2017) 487–495 [11] Juan AA, Mendez CA, Faulin J, De Armas J, Grasman SE. Electric vehicles in logistics and transportation: a survey on emerging environmental, strategic, and operational challenges. Energies 2016;9:1–21. [12] Shao S, Guan W, Ran B, He Z, Bi J. Electric vehicle routing problem with charging time and variable travel time. Math Probl Eng 2017;2017:1–13. [13] Jokinen R, Pettersson F, Saxén H. An MILP model for optimization of a smallscale LNG supply chain along a coastline. Appl Energy 2015;138:423–31. [14] Yagcitekin B, Uzunoglu M. A double-layer smart charging strategy of electric vehicles taking routing and charge scheduling into account. Appl Energy 2016;167:407–19. [15] Moura SJ, Fathy HK, Callaway DS, Stein JL. A stochastic optimal control approach for power management in plug-in hybrid electric vehicles. IEEE Trans Control Syst Technol 2011;19:545–55. [16] Kou P, Gao F, Guan X. Stochastic predictive control of battery energy storage for wind farm dispatching: using probabilistic wind power forecasts. Renew Energy 2015;80:286–300. [17] Berhan E, Beshah B, Kitaw D, Abraham A. Stochastic vehicle routing problem: a literature survey. J Inform Knowl Manage 2014;13:1450022. [18] Smith SL, Pavone M, Bullo F, Frazzoli E. Dynamic vehicle routing with priority classes of stochastic demands. SIAM J Control Optim 2010;48:3224–45. [19] Zhang J, Lam WHK, Chen BY. A stochastic vehicle routing problem with travel time uncertainty: trade-off between cost and customer service. Networks Spatial Econ 2013;13:471–96. [20] Stewart WR, Golden BL. Stochastic vehicle routing: a comprehensive approach. Eur J Oper Res 1983;14:371–85. [21] Laporte G, Louveaux F, Mercure H. Models and exact solutions for a class of stochastic location-routing problems. Eur J Oper Res 1989;39:71–8. [22] Bastian C, Rinnooy Kan AH. The stochastic vehicle routing problem revisited. Eur J Oper Res 1992;56:407–12. [23] Laporte G, Louveaux FV, Mercure H. A priori optimization of the probabilistic traveling salesman problem. Oper Res 1994;42:543–9. [24] Gendreau M, Laporte G, Séguin R. An exact algorithm for the vehicle routing problem with stochastic demands and customers. Transport Sci 1995;29:143–55. [25] Osorio C, Bierlaire M. A simulation-based optimization framework for urban transportation problems. Oper Res 2013;61:1333–45. [26] Chepuri K, Homem-de Mello T. Solving the vehicle routing problem with stochastic demands using the cross-entropy method. Ann Oper Res 2005;134:153. [27] Pedrag M, Miomir J. The advanced system for dynamic vehicle routing in the process of waste collection. Mech Eng 2011;9:127–36. [28] Benjamin A, Beasley J. Metaheuristics for the waste collection vehicle routing problem with time windows, driver rest period and multiple disposal facilities. Comput Oper Res 2010;37:2270–80. [29] Mata-Pérez M, Saucedo-Martínez JA. Optimization of the distribution of steel pipes using a mathematical model. Dyna 2015;82:69–75.

495

[30] Sørensen CG, Bochtis DD. Conceptual model of fleet management in agriculture. Biosyst Eng 2010;105:41–50. [31] Orfanou A, Busato P, Bochtis DD, Edwards G, Pavlou D, Sørensen CG, et al. Scheduling for machinery fleets in biomass multiple-field operations. Computers and Electronics in Agriculture 2013;94:12–9. [32] Sumalee A, Uchida K, Lam W. Stochastic multi-modal transport network under demand uncertainties and adverse weather condition. Transport Res Part C: Emerg Technol 2011;19:338–50. [33] Fontaras G, Zacharof N-G, Ciuffo B. Fuel consumption and CO2 emissions from passenger cars in Europe – laboratory versus real-world emissions. Prog Energy Combust Sci 2017;60:97–131. [34] Sharma B, Ingalls RG, Jones CL, Huhnke RL, Khanchi A. Scenario optimization modeling approach for design and management of biomass-to-biorefinery supply chain system. Bioresour Technol 2013;150:163–71. [35] Krajnc N. Wood fuels handbook. Food and agricultural organization of the United Nations; 2015. . [36] Eom J, Schipper L, Thompson L. We keep on truckin’: trends in freight energy use and carbon emissions in 11 IEA countries. Energy Policy 2012;45:327–41. [37] BC Bioenergy Network and Biomass Availability Study Working Group. Biomass availability study for district heating systems; 2012. [Online; accessed 04-October2015]. [38] Briggs D. Forest products measurements and conversion factors: with special emphasis on the U.S. Pacific Northwest. Seattle Wash.: College of Forest Resources University of Washington; 1994. [39] MacDonald A. Assessment of economically accessible biomass. Technical Report Quebec, Canada: FP Innovations; 2009. [40] Loecher M, Ropkins K. RgoogleMaps and loa: unleashing R graphics power on map tiles. J Stat Software 2015;63:1–18. [41] Bivand RS, Pebesma EJ, Gómez-Rubio V. Applied spatial data analysis with R. 2nd ed. Springer; 2008. [42] Shekhar S, Xiong H. Variogram modeling. In: Encyclopedia of GIS. Boston, MA: Springer US; 2008. p. 1217. [43] R Core Team. R: a language and environment for statistical computing; 2016. . [44] Maechler M, Rousseeuw P, Struyf A, Hubert M, Hornik K. Cluster: cluster analysis basics and extensions. R package version 2.0.3 — For new features, see the ’Changelog’ file (in the package source); 2015. [45] Kaufman L, Rousseeuw PJ. Finding groups in data: an introduction to cluster analysis. John Wiley; 1990. [46] Hollander M. Nonparametric statistical methods. Hoboken, New Jersey: John Wiley & Sons, Inc.; 2013. [47] Perrot P. A to Z of thermodynamics. Supplementary series, vol. 27. Oxford University Press; 1998. , https://books.google.ca/books?id=n3hKx7u71_MC.