Multi-Use Optimization of a Depot for Battery-Electric Heavy-Duty Trucks

By inergency On Feb 26, 2024

2.1. Optimization Model

The optimization model eFlame was primarily developed to optimize several use cases for bidirectional charging separately. In [14,26], the use cases arbitrage trading and self-consumption optimization are elaborated upon. The use case peak shaving is dealt with in [7]. The novelty of the present paper is the combination of the use cases in the context of a multi-use optimization that was not implemented before. Figure 1 shows all power flows relevant to the optimization. At this point, we describe the optimization problem covering decision variables, objective function, and constraints. Linear optimization is a method that can be used to solve problems where the objective function and constraints are linear functions of the decision variables. The constraints can be formulated as equalities or inequalities. A mixed-integer linear program (MILP) problem includes integer decision variables. A comprehensive introduction to linear optimization is given in [27]. Since the model was primarily developed for optimizing battery-electric cars, the vehicles are generally referred to as EVs in the following model description.

For all decision variables, the non-negativity constraint applies. The constraint is exemplarily defined in Equation (1) for the received power

P_{t}^{G C P, i n}

, and the feed-in power

P_{t}^{G C P, o u t}

at the Grid Connection Point (GCP), but it can be applied to the remaining decision variables. The total number of time steps t in the observation horizon is represented by n.

$P_{t}^{G C P, i n} ⩾ 0, P_{t}^{G C P, o u t} ⩾ 0 \forall t \in T = \{1, . . ., n\}$

(1)

The photovoltaic (PV) generation is not a decision variable, but it may be influenced during the optimization via the curtailment $P_{t}^{c u r t}$ . With this optimization variable, the generation of the PV system can be reduced, e.g., to prevent feed-in at negative prices. Using the decision variable $P_{t}^{G C P, p e a k}$ , the maximum power at the grid connection point is determined. The charging power $P_{t}^{c h a r g e}$ and discharging power $P_{t}^{d i s c h a r g e}$ and the energy capacity of the battery $E_{t}^{E V}$ are further decision variables that are related to the EVs. Furthermore, there is the decision variable $P_{t}^{v 2 g}$ , which is used to observe how much energy from the vehicles is fed back into the grid. The remaining decision variables $b_{t}^{c h a r g e}$ , $b_{t}^{d i s c h a r g e}$ , $b_{t}^{o u t}$ , and $b_{t}^{i n}$ are boolean variables, which are used to ensure that the power flow exchanges with the vehicles and the grid connection point are only in one direction at any time instant.

The objective of the optimization model is to maximize the revenue. The established objective function shown in Equation (2) consists of four terms: the cashflow from arbitrage trading at the spotmarket

{CF}^{s p o t}

, costs through levies

C^{l e v i e s}

, costs through grid fees

C^{g r i d f e e}

, and a term that evaluates the opportunity costs due to battery degradation

C^{b a t, d e g}

. The different terms are defined in the following.

$\begin{matrix} max ({CF}^{s p o t} - C^{l e v i e s} - C^{g r i d f e e} - C^{b a t, d e g}) \end{matrix}$

(2)

The cash flow, the difference between cash in- and outflows, from arbitrage trading at the spot market

C F^{s p o t}

is calculated in Equation (3). Different market data can be selected for the price time series

p_{t}^{i n}

and

p_{t}^{o u t}

, but constant values may also be used.

${CF}^{s p o t} = \sum_{t = 1}^{n} ({P_{t}}^{G C P, i n} \cdot p_{t}^{i n} \cdot Δ t - {P_{t}}^{G C P, o u t} \cdot p_{t}^{o u t} \cdot Δ t) \forall t \in T$

(3)

Consumers have to pay a gridfee

C^{g r i d f e e}

to the Distibution System Operator (DSO) for the use of the grid infrastructure. In Germany, the grid fee for commercial customers is divided into a usage price

p^{u s a g e}

and a capacity price

p^{c a p}

. The usage price depends on the energy consumed, whereas the capacity price depends on the annual peak power. The grid fee is included in the objective function of Equation (4).

$C^{g r i d f e e} = \sum_{t = 1}^{n} {P_{t}}^{GCP, i n} \cdot p^{u s a g e} \cdot Δ t + P_{t = n}^{GCP, p e a k} \cdot p^{c a p} \forall t \in T$

(4)

Additionally, various taxes and levies are charged on electricity, and those are summarized through

C^{l e v i e s}

and shown in Equation (5). Stationary battery storage may be partially exempt from levies, and such an exemption is also being discussed for bidirectional vehicles. The problem is to determine how much energy is actually fed back into the grid. This is especially problematic in combination with PV systems. Via the subtrahend of Equation (5), a partial exemption from the levies on energy fed back into the grid is implemented. The decision variable

P_{t}^{v 2 g}

represents the power the vehicles feed into the grid and is introduced later in (17) and (18). It is an auxiliary variable calculated from the other power variables and is therefore not directly included in the power balance following in Equation (7). A partial exemption may be dynamically parameterized via the levies on V2G

p^{l e v i e s, v 2 g}

that are still charged even if the energy is fed back into the grid. If

p^{l e v i e s, v 2 g}

is set equal to

p^{l e v i e s}

, no exemption occurs. A full exemption can be achieved by setting

p^{l e v i e s, v 2 g}

equal to zero.

$\begin{matrix} C^{l e v i e s} = \sum_{t = 1}^{n} {P_{t}}^{GCP, i n} \cdot p^{l e v i e s} \cdot Δ t - \sum_{t = 1}^{n} {P_{t}}^{v 2 g} \cdot (p^{l e v i e s} - p^{l e v i e s, v 2 g}) \cdot Δ t \forall t \in T \end{matrix}$

(5)

The opportunity costs from battery degradation

C_{b a t, d e g}

are included in the optimization problem by using Equation (6) based on [28]. The calculation of the degradation costs

C_{b a t, d e g}

is based on the use of the battery and determined by the decrease of the available capacity

C^{l o s s}

from a cycling aging model. The costs result primarily from the total charge quantity throughput, which is defined by the charging and discharging power. The price of the battery is represented by

c_{b a t, b u y}

, and

E^{E V, m a x}

is the capacity of the battery. The used model assumes the end of life of the battery at a loss of 20% of the initial capacity.

$\begin{matrix} C_{b a t, d e g} = \frac{c_{b a t, b u y} \cdot E^{E V, m a x}}{20 %} C_{l o s s} (P_{t}^{E V, c h a r g e}, P_{t}^{E V, d i s c h a r g e}) \forall t \in T \end{matrix}$

(6)

The optimization model is restricted by several constraints concerning the GCP and the EVs. We start by introducing the boundary conditions of the GCP. According to the law of conservation of energy, the incoming power flows at the GCP must be equal to the outgoing power flows. This is ensured by Equation (7). The load profile of the building

P_{t}^{b u i l d}

is integrated into the optimization as a static time series.

$\begin{matrix} P_{t}^{GCP, i n} + \sum_{i = 1}^{n_{E V}} P_{t}^{E V, d i s c h a r g e} + P_{t}^{P V} = \\ P_{t}^{GCP, o u t} + \sum_{i = 1}^{n_{E V}} P_{t}^{E V, c h a r g e} + P_{t}^{c u r t} + P_{t}^{b u i l d} \forall t \in T \end{matrix}$

(7)

For the determination of the grid fee

C^{g r i d f e e}

in Equation (4), the annual peak power at the GCP

P_{t}^{GCP, p e a k}

is required. Using Equation (8), the power peak is updated continuously during the optimization. Thus, the last time step n contains the annual power peak.

$P_{t}^{GCP, p e a k} ⩾ P_{t}^{GCP, i n}, P_{t}^{GCP, p e a k} ⩾ P_{t - 1}^{GCP, p e a k} \forall t \in T$

(8)

Equations (9) and (10) are introduced to prevent energy from being purchased and fed in simultaneously at the GCP. In consequence, the boolean decision variables

b_{t}^{i n}

and

b_{t}^{o u t}

are used.

P^{GCP, m a x}

describes the maximum grid connection capacity, which results from the transformer and structural conditions at the grid connection point. The combination of Equations (8) and (9) ensures the grid connection capacity

P^{GCP, m a x}

is always greater than or equal to the annual power peak

P_{t}^{GCP, p e a k}

$\begin{matrix} P^{GCP, m a x} \cdot b_{t}^{i n} ⩾ P_{t}^{GCP, i n}, P^{G C P, m a x} \cdot b_{t}^{o u t} ⩾ P_{t}^{GCP, o u t} \forall t \in T \end{matrix}$

(9)

$\begin{matrix} b_{t}^{o u t} + b_{t}^{i n} ⩽ 1 \forall t \in T \end{matrix}$

(10)

The following constraints are related to the EVs and apply separately for each EV. The energy balance of the vehicle battery must be maintained to preserve the physical consistency of the EVs. The energy stored in the EV battery in the first time step is defined by the constraint Equation (11). For the first time step, this equation defines the stored energy as equal to the initial stored energy plus the charged energy at the GCP minus the discharged energy and the energy consumed during trips

E_{t}^{E V, t r i p}

plus the energy charged at public stations

E_{t = 1}^{E V, p u b l i c}

. Constant efficiencies for charging

η^{E V, c h a r g e}

and discharging

η^{E V, d i s c h a r g e}

are considered.

$\begin{matrix} E_{t = 1}^{E V} = {S O C}_{t = 1}^{E V} \cdot E^{E V, m a x} + P_{t = 1}^{E V, c h a r g e} \cdot η^{E V, c h a r g e} \cdot Δ t \\ - P_{t = 1}^{E V, d i s c h a r g e} \cdot η^{E V, d i s c h a r g e} \cdot Δ t - E_{t = 1}^{E V, t r i p} + E_{t = 1}^{E V, p u b l i c} \end{matrix}$

(11)

For the remaining time steps, Equation (12) applies, where the initially stored energy is replaced by the stored energy of the previous time step.

$\begin{matrix} E_{t}^{E V} = E_{t - 1}^{E V} + P_{t}^{E V, c h a r g e} \cdot η^{E V, c h a r g e} \cdot Δ t \\ - P_{t}^{E V, d i s c h a r g e} \cdot η^{E V, d i s c h a r g e} \cdot Δ t - E_{t}^{E V, t r i p} + E_{t}^{E V, p u b l i c} \forall t \in \{2, . . ., n\} \end{matrix}$

(12)

Equation (13) ensures that the vehicles are always charged with a minimum State of Charge

{S O C}^{E V, d e p, m i n}

at departure. The condition is only valid for the time steps in which a vehicle departs, as indicated by the boolean variable

b_{t}^{E V, d e p}

. This variable is determined before the optimization based on the driving profiles and is only equal to one if the vehicle departs. To ensure that the condition can also be met if the vehicle is only plugged in for a short time and thus the minimum SOC cannot be reached, a buffer

E_{t}^{b u f f e r}

is integrated into the condition. This buffer is also determined before the optimization.

$\begin{matrix} E_{t}^{E V} + E_{t}^{b u f f e r} = {S O C}^{E V, d e p, m i n} \cdot E^{E V, m a x} \cdot b_{t}^{E V, d e p} \forall t \in T \end{matrix}$

(13)

Apart from public charging, each EV can only be charged or discharged if it is connected to a charging point at the GCP, and this is ensured by Equations (14) and (15). The boolean variable

b_{t}^{E V}

is determined before the optimization based on the driving profiles. If the vehicle is plugged in, the variable is one, and otherwise it is zero. We assume that each vehicle has its own charging point. To prevent the EVs from charging and discharging at the same time, the decision variables

b_{t}^{c h a r g e}

and

b_{t}^{d i s c h a r g e}

are added to Equations (14) and (15). Equation (16) prevents both variables from being equal to one simultaneously. If only unidirectional charging is considered, the Equations (15) and (16) are omitted, and

P_{t}^{E V, d i s c h a r g e}

is set to zero via a further boundary condition.

$b_{t}^{E V} \cdot b_{t}^{c h a r g e} \cdot P^{E V, c h a r g e, m a x} ⩾ P_{t}^{E V, c h a r g e} \forall t \in T$

(14)

$b_{t}^{E V} \cdot b_{t}^{d i s c h a r g e} \cdot P^{E V, d i s c h a r g e, m a x} ⩾ P_{t}^{E V, d i s c h a r g e} \forall t \in T$

(15)

$\begin{matrix} b_{t}^{c h a r g e} + b_{t}^{d i s c h a r g e} ⩽ 1 \forall t \in T \end{matrix}$

(16)

Finally, boundary conditions are required to determine the power fed back from the EVs into the grid

P_{t}^{v 2 g}

. This variable is necessary to calculate the exemption from levies in Equation (5). Therefore, we choose a power balance based approach and rearrange Equation (7) according to the discharged energy. Since power can only be fed into the grid if no energy is purchased,

P_{t}^{GCP, i n}

is set to zero. The discharged power is replaced by the introduced decision variable

P_{t}^{v 2 g}

, resulting in Equation (17). The boundary condition in Equation (18) ensures that

P_{t}^{v 2 g}

cannot become greater than the feed-in power.

$\begin{matrix} P_{t}^{v 2 g} ⩽ P_{t}^{GCP, o u t} - P_{t}^{P V} + P_{t}^{c u r t} + P_{t}^{b u i l d} + \sum_{i = 1}^{n_{E V}} P_{t}^{E V, c h a r g e} \forall t \in T \end{matrix}$

(17)

$\begin{matrix} P_{t}^{v 2 g} ⩽ P_{t}^{G C P, o u t} \forall t \in T \end{matrix}$

(18)

Since the model is intended to examine entire years and since the use of boolean variables makes it a mixed-integer optimization problem, the computational effort required to solve the problem is rather high. In order to be able to solve it with a reasonable computational effort, the model is computed as a rolling optimization. The determination of the annual power peak is a special aspect of the rolling optimization, which will be explained in the following using the schematic diagram in Figure 2. For rolling optimization, the whole optimization period is divided into m smaller optimization time periods of uniform size. In individual optimization steps, each of the smaller optimization periods is optimized one after the other. The results of an optimization step are passed as start values to the next step. By using an overlapping period, we increase the prediction horizon for the optimization. After the m-th step, the first run of the optimization is finished. According to Equation (8), the power peak is continuously updated as shown in Figure 2 below. As can be seen in the figure, the first optimization steps are limited by a lower power peak compared to the later steps. Therefore, in a second optimization run, the affected steps before the occurrence of the annual power peak are optimized again with the updated power peak.

The sequence of the used optimization model eFlame is illustrated schematically in Figure 3. After importing the input parameters and input data described below, the optimization problem is set up. The optimization problem is solved sequentially considering the charging strategies: uncontrolled charging (ref), unidirectional charging (uni), and bidirectional charging (bidi). The results are examined separately for each charging strategy.

2.2. Input Data

As mentioned in Section 1, prior research on the topic of BETs has relied on assumptions regarding driving profiles. In this paper, we had the opportunity to use real-life data from a depot of a freight forwarding company in Germany. The company primarily operates in the short-haul segment. The data were provided within the framework of the project NEFTON in which partners from industry and science jointly develop a Megawatt Charging System (MCS) for BETs. Mobility data of the company’s trucks, historical load profiles of its buildings, and information about the PV system are included in the data. The selected depot can serve as a real-life example.

In the project NEFTON, driving data from several fleets of German fleet operators were recorded using high-resolution GPS data loggers. The recorded dataset includes 1.26 million km of driving data and is openly available in anonymized form in [29]. Only the driving data of the depot under consideration were extracted from this dataset. Since the data were recorded for trucks using diesel fuel, our investigation builds on the observation that the company desires to keep its services in the same way with electric trucks. The data are available for different lengths of time and were extended to uniform periods using a Markov process. To avoid oversizing the vehicle batteries, the missions in the dataset are divided into two clusters depending on the distance traveled. Missions with a distance of more than 200 km are grouped into the cluster regional transport and those with less than 200 km into the cluster local transport, which is similar to the classification of [30]. The annual driving profiles are taken as given and are presented in the following. Figure 4 shows the average percentage of vehicles in different locations for the two clusters. It can be seen that especially the mobility profiles from the Local Transport cluster have very high idle times at the depot and that at least 50% of the BETs are always present at the depot. On weekends and at night, most of the vehicles are located at the depot. The driving profiles of the cluster Regional Transport show significantly lower idle times at the depot. During daytime on weekdays, 80% of the vehicles are absent. On weekends, almost 40% are not at the depot. In addition, the driving profiles of the Regional Transport cluster show high parking durations in industrial areas and other locations. The difference between the two clusters is also evident from the characteristic values included in Table 1. The annual kilometrage of the Local Transport cluster is about 14,000 km. This is significantly lower than the kilometrage of the Regional Transport cluster of about 66,000 km. The electrical energy consumption for the driving profiles is determined using the model from [31]. The average annual consumption determined in this way is also included in Table 1. The variables

b_{t}^{E V, d e p}

E_{t}^{b u f f e r}

b_{t}^{E V}

E_{t}^{E V, t r i p}

and

E_{t}^{E V, p u b l i c}

are determined based on the driving profiles and serve as inputs for the optimization model.

In addition to the driving profiles, the load profile of the building of the depot

P_{t}^{b u i l d}

is another important input for the optimization. The used load profile shown in Figure 5b for an average week relies on real data of the depot. From the annual time series, the average was determined for each quarter-hour of the week as well as the ranges in which 80% and 100% of the values lie. The plot shows that there are significant load peaks in the evening hours on weekdays, indicating suitability for peak shaving. The load is significantly lower at weekends and at night than it is during the day on weekdays.

We assume that the depot pays variable electricity prices based on the prices of the electricity exchange. Therefore, we used the intraday auction prices as electricity prices

p_{t}^{i n}

and

p_{t}^{o u t}

for the optimization. In Europe, there are various short-term markets on the power exchange. One of those markets is the intraday auction. Due to the shorter time slices of quarter hours compared to the day-ahead market, in which hourly products are traded, this market offers higher price spreads. Thus, the revenue opportunities for flexibilities like bidirectional EVs are increased. The development of the prices of the intraday auction from the beginning of 2019 to the end of 2022 is shown in Figure 5a. As a consequence of the energy crisis, the price has risen from around 4 ct/kWh to a maximum of over 70 ct/kWh, and also the price spreads increased significantly.

The PV generation is determined as a time series depending on the historical irradiation data on CAMS level as a function of the orientation of the PV plant and its peak power [32]. The irradiation data are used for the location of the depot for the weather year 2012. The weather year is chosen based on the recommendation in [33].

2.3. Input Parameters

After introducing the data source and the model in the previous sections, the input parameters are presented in the following. For this purpose, we define a base scenario for which the input parameters are listed in Table 2. By varying different parameters of this base scenario, various sensitivities are examined. For the sensitivity analysis, one parameter of the base scenario is changed, while the rest of the parameters are left unchanged. The varied parameters of the sensitivity analysis are also included in the table. The base year is 2021, and the optimization is performed at a time step size of 15 min. As Figure 2 illustrates, we use a rolling approach and divide the examined years into 61 optimization steps. The observation period of each step is seven days, consisting of the optimization period of six days and one day of overlap. In contrast to real-world charging management systems that apply forecasts, we assume perfect foresight for each optimization step. In the base scenario, no exemption of levies on energy fed back into the grid is assumed. Therefore,

p^{l e v i e s}

is set to be equal to

p^{l e v i e s, v 2 g}

. However, the exemption is considered in the sensitivity analysis. In the base scenario, no limitation of the grid connection capacity is considered. Thus,

P^{GCP, m a x}

is set to the oversized value of 5 MW. A limitation of

P^{GCP, m a x}

is examined in the sensitivity analysis. The grid connection capacity is minimally limited to 700 kW, since a lower capacity would result in the curtailment of the PV system in times of high irradiation. The feed-in tariff of 0.06 EUR/kWh is an assumed value suitable for Germany and is only used in the reference simulation as

p_{t}^{o u t}

. It is also assumed that 30 BETs of the depot are electrified. The number of electric vehicles is one of the sensitivities examined. According to the distribution from the dataset, 30% of the vehicles are used for regional traffic and 70% are used for local traffic. The appropriate driving profiles are divided among the BETs according to the distribution, and a battery capacity of 250 kWh for local and 500 kWh for regional traffic is assumed. Based on [34], the price of the vehicle battery

c_{b a t, b u y}

is set to 139 EUR/kWh. The parameters of the PV system are selected according to the system of the real depot.

In addition to the year 2021 of the base scenario, the years 2019, 2020, and 2022 are also examined. For the optimization of the different years, several parameters have to be varied. In contrast, only one parameter is changed at a time in the sensitivity analysis presented below. The other parameters remain unchanged. In consequence, these year-dependent parameters are separated in Table 3. For the reference simulation, a constant price based on the average day-ahead price is assumed for

p_{t}^{i n}

[35]. For the levies, the real historical values for Germany from [36] are used. The prices for the grid fees are also based on historical values of the grid operator Netze BW, where the depot under consideration is located [37]. We use the prices for medium voltage networks and consider an annual usage time of less than 2500 h.