Open-access A TWO-LEVEL MULTI-OBJECTIVE MODEL AND SIMULATED ANNEALING METAHEURISTIC FOR COST, CO2 EMISSION, AND CUSTOMER SATISFACTION OPTIMIZATION IN LAST-MILE DELIVERY

ABSTRACT

Last-mile delivery-the final and most critical stage of fulfilling customer orders-is typically the most costly and time-consuming component of the distribution process. Optimizing this stage is therefore essential for improving the overall efficiency and sustainability of logistics operations. This study begins with a comprehensive review of the principal contributions to last-mile delivery optimization. Building on insights from existing models and metaheuristics, we propose a new two-level framework in which customers are served through distribution centers, ensuring that each customer location lies within a designated coverage area. The framework consists of two interconnected levels: depot to distribution center routing and distribution center to customer delivery. To support this structure, we develop a new optimization model and a Simulated Annealing-based metaheuristic inspired by the vehicle routing problem with time windows and tailored specifically to last-mile customer service operations. Computational results demonstrate that the proposed approach provides both economic benefits-through reduced travel distance-and environmental advantages through lower CO2 emissions. Moreover, the routing plan generated at the second level allows for the quantification of customer satisfaction-averaging 72.5% for the optimization model and 100% for the metaheuristic-thereby providing a more comprehensive and robust evaluation of delivery.

Keywords:
combinatorial optimization; logistics; last-mile; vehicle routing; coverage area

1 INTRODUCTION

Urban population growth has progressively intensified the demand for high-quality services capable of meeting the increasingly diverse needs of citizens. Addressing these needs requires the application of analytical and theoretical frameworks to support a society that is continually becoming more demanding. Within this context, the e-commerce industry has experienced remarkable expansion, particularly across emerging economies. For instance, in India, e-commerce revenues rose from 14 billion USD in 2014 to 50 billion USD in 2019, with projections estimating a further increase to 200 billion USD by 2026.

The COVID-19 pandemic-triggered by the emergence of the novel severe acute respiratory syndrome coronavirus (SARS-CoV-2)-has also played a critical role in reshaping consumer behavior, fostering a strong preference for home delivery of goods (Tiwari & Sharma, 2023). As a result, logistics has become even more central to modern society, as it must continuously adapt to meet evolving customer expectations. Consequently, the volume and complexity of urban delivery operations have expanded substantially.

A review conducted by Jazemi et al. (2023) identifies various emerging approaches in last-mile delivery, such as crowdshipping, parcel lockers, sidekick delivery, and optional collection points. This review provides valuable insights into the strengths and limitations of existing strategies, offering a foundation for identifying future research opportunities. Complementarily, Giuffrida et al. (2023) examine the integration of optimization and machine learning techniques in last-mile logistics, underscoring the growing relevance of hybrid approaches.

Efficient last-mile delivery requires substantial resource allocation, making the reduction of operational costs through optimized routing and scheduling a priority in e-commerce logistics. Özarik et al. (2021) proposed a mixed-integer linear programming model combined with the Adaptive Large Neighborhood Search (ALNS) metaheuristic, demonstrating cost reductions of up to 40%. More recently, Tausif (2023) explored the use of autonomous vehicles and drones to support last-mile operations, developing a mixed-integer linear programming model that yielded promising results.

The evolution of transportation methods for last-mile delivery, alongside emerging innovations and operational challenges, has fueled a growing academic interest in the field. Zhu et al. (2023), through a systematic review of 150 articles following PRISMA guidelines, emphasize the rapid proliferation of proposed approaches and the resulting opportunities for identifying new research directions.

Transportation remains a critical element in service delivery, as it directly influences operational efficiency and cost minimization. Building on this premise, the present study proposes a novel two-level optimization model, complemented by a Simulated Annealing-based metaheuristic, to address last-mile delivery challenges in urban environments. For a more detailed illustration of the research framework and its relevance, refer to Figure 1.

Figure 1
Two-level framework to last-mile delivery : Depot → Distribution center → Customers.

The remainder of the article is structured as follows: Section 2 reviews the literature relevant to this study. Section 3 presents the proposed two-level framework. Section 4 discusses the computational results, and Section 5 outlines the main conclusions and references.

2 REVIEW OF SPECIALIZED LITERATURE

2.1 Combinatorial Optimization

Combinatorial optimization comprises a broad class of problems in which the goal is to identify an optimal or near-optimal solution from a discrete and finite, yet often exponentially large, set of feasible alternatives. As outlined by Peres & Castelli (2021), optimization methodologies can be generally classified into three major categories: enumerative, deterministic, and stochastic approaches. Enumerative methods explore all possible solutions exhaustively and guarantee the identification of a global optimum. However, their computational burden grows exponentially with problem size, making them impractical for medium- and large-scale real-world instances.

Deterministic and stochastic strategies provide more tractable alternatives. Deterministic approaches rely on mathematical formulations and structured search rules that systematically guide the exploration of the solution space. Stochastic approaches, in contrast, incorporate probabilistic mechanisms designed to enhance diversification and prevent premature convergence to suboptimal regions. Although neither guarantees optimality, both types of methods offer a more favorable balance between computational efficiency and solution quality-an essential trade-off for complex logistical and routing applications.

Optimization problems can be categorized as unconstrained or constrained. Unconstrained optimization focuses solely on maximizing or minimizing an objective function. Constrained optimization, by contrast, introduces feasibility requirements-such as vehicle capacity limits, time windows, and geographic restrictions-that shape the set of allowable solutions. Last-mile delivery routing is inherently a constrained optimization problem due to its operational, spatial, and temporal conditions.

Formally, a combinatorial optimization problem is defined by a solution space S, an objective function f , and a set of constraints Ω. The task is to determine a globally optimal solution s ∗ ∈ S such that f(s ∗) ≤ f(s) for all s ∈ S in a minimization framework (or the reverse inequality for maximization). While this formulation is conceptually simple, the cardinality of the solution space typically grows exponentially with the number of decision variables, making exact algorithms computationally prohibitive for large instances.

For this reason, heuristics and metaheuristics-such as Simulated Annealing, Tabu Search, Genetic Algorithms, and Variable Neighborhood Search-have become indispensable tools for addressing large-scale combinatorial optimization problems. These methods provide high-quality approximate solutions within reasonable computational times, making them especially valuable in real-world last-mile delivery scenarios characterized by multiple interacting constraints and diverse performance objectives.

2.2 Logistics and Last-Mile Delivery

Logistics encompasses the coordinated management of material, informational, and financial flows across supply chains, with the goal of ensuring efficient, reliable, and cost-effective delivery of goods and services. As highlighted by Zijm et al. (2019), modern logistics incorporates a wide range of specialized domains, including autonomous transportation systems, blockchain-enabled supply chains, green and sustainable logistics, humanitarian logistics, and reverse logistics. Advances in digitalization, data analytics, and automation have further transformed logistical operations, enabling faster decision-making and greater adaptability to demand fluctuations.

Within this broader context, last-mile delivery represents the final stage of the distribution process, involving the transportation of goods from local depots or distribution centers to the end consumer-typically in densely populated urban environments. According to Boysen et al. (2021), last-mile delivery is characterized by four primary challenges: (a) the rapid growth of parcel demand driven by urbanization and the rise of e-commerce; (b) environmental and sustainability concerns associated with urban congestion, fuel consumption, and CO2 emissions; (c) high operational costs, particularly those related to labor and vehicle utilization; and (d) increasing customer expectations regarding delivery speed, reliability, and flexibility.

Recent research underscores the ongoing transformation of last-mile delivery systems. Jazemi et al. (2023) identify emerging delivery models-including crowdshipping, parcel lockers, sidekick delivery, and optional pickup points-that seek to alleviate urban congestion and improve service quality. Each of these models introduces distinct operational trade-offs, such as increased reliance on third-party actors, varying degrees of infrastructure investment, and different impacts on routing complexity.

Furthermore, the integration of machine learning and optimization techniques has opened new avenues for improving last-mile operations. See Figure 2. Giuffrida et al. (2023) highlight the growing use of predictive analytics, classification models, and hybrid OR-ML frameworks to enhance demand forecasting, routing efficiency, and resource allocation. Such data-driven approaches enable logistics providers to better anticipate customer behavior, adapt routing plans dynamically, and optimize operational costs.

Figure 2
Bibliometric analysis of last-mile delivery and its connections to optimization, metaheuristics, and machine learning.

Vehicle routing problems (VRPs) constitute the analytical backbone of last-mile delivery planning. Numerous VRP extensions incorporate practical constraints such as time windows, heterogeneous fleets, environmental objectives, energy consumption, and electric vehicle limitations. Examples include the Green Vehicle Routing Problem (Sabet & Farooq, 2022), which emphasizes environmental sustainability, and models involving autonomous or electric vehicles, as examined by Stamadianos et al. (2023). These variants capture the increasing complexity and multidimensionality of modern urban logistics systems.

Hybrid optimization approaches continue to gain prominence in the literature. Moradi et al. (2023) propose a mixed-integer programming formulation combined with a hybrid Variable Neighborhood Search-Simulated Annealing method to address autonomous electric vehicle routing. Similarly, Li & Yang (2024) develop a variant of the delivery-and-collection routing problem solved using Variable Neighborhood Search and Randomized Tabu Thresholding, demonstrating the effectiveness of hybrid heuristics for dealing with complex operational environments.

Machine learning has also been applied to performance assessment, demand forecasting, and service quality evaluation within last-mile logistics. Brochado et al. (2024) and Aljohani (2024) provide evidence that ML-based models can deliver valuable insights into customer satisfaction, operational performance, and network efficiency. Meanwhile, Boresta et al. (2024) explore the integration of OR and ML for equitable cost allocation and balanced resource sharing, showing promising results for sustainable and efficient logistics operations.

Finally, the inherent complexity of last-mile delivery-stemming from fluctuating demand, urban congestion, and operational constraints-necessitates the use of approximate methods such as heuristics, metaheuristics, and reinforcement learning. Recent studies, including Senna et al. (2024), highlight the growing relevance of these techniques for developing robust, scalable, and adaptive routing strategies suited to dynamic urban environments. See Figure 2.

2.3 Research Gap and Contributions

Despite these advances, most existing contributions address last-mile delivery either from a single optimization level or by focusing on a limited subset of performance criteria. In particular, the integration of multi-level distribution structures, environmental performance metrics, and customer satisfaction indicators within a unified optimization framework remains relatively under-explored. Furthermore, many models assume a direct connection between depots and customers, without explicitly modeling intermediate distribution centers and their associated coverage areas.

The present study contributes to the literature by addressing these limitations through a two-level multi-objective optimization framework that simultaneously considers routing efficiency, environmental impact, and service quality. Unlike traditional single-level vehicle routing formulations, the proposed approach explicitly models the hierarchical structure of urban logistics systems, distinguishing between upstream distribution from a central depot to distribution centers (Level I) and downstream delivery from distribution centers to customers located within predefined coverage areas (Level II). This structure allows for a more realistic representation of urban logistics networks commonly observed in practice.

Another key contribution lies in the integration of customer satisfaction as a quantitative performance indicator directly embedded within the optimization framework. While prior studies often evaluate service quality separately from routing decisions, this research incorporates time window compliance as a measurable objective, enabling simultaneous evaluation of operational efficiency and service reliability. Additionally, the proposed framework incorporates CO2 emission estimation based on fuel consumption models. This integrated perspective extends existing green logistics models by explicitly linking routing decisions to environmental performance within a multi-level structure.

From a methodological standpoint, this study develops a Simulated Annealing metaheuristic specifically adapted to the two-level routing structure, incorporating neighborhood search operators such as relocation, swap, and 2-opt improvements. While Simulated Annealing has been widely applied in combinatorial optimization, its tailored implementation for a hierarchical last-mile delivery structure with coverage areas represents a novel contribution. The metaheuristic is designed to maintain feasibility with respect to vehicle capacity and time window constraints while exploring the solution space efficiently.

Finally, the proposed approach contributes to the literature by providing a comprehensive performance evaluation framework that simultaneously measures: total traveled distance; operating costs; CO2 emissions; delivery completion rate and customer satisfaction levels.

This multi-criteria evaluation enables decision-makers to analyze trade-offs between economic efficiency, environmental sustainability, and service quality, which are critical dimensions in modern urban logistics systems.

3 VEHICLE ROUTING ASSOCIATED WITH LAST-MILE DELIVERY

This section introduces the proposed two-level methodological framework. The first level focuses on the upstream distribution process, in which a central depot located on the outskirts of the city supplies multiple distribution or sales centers positioned strategically within the urban area. The second level concerns the downstream stage, where customers located within predefined coverage areas are served by their associated distribution centers. At this level-corresponding to last-mile operations-an urban delivery vehicle follows an optimized route to satisfy customer demand. Both levels of the framework are formulated through a binary integer linear programming model and complemented by a Simulated Annealing metaheuristic designed to address the computational challenges of larger instances. See Figure 1.

3.1 Optimization Model

3.1.1 Level I: Multi-Product Routing from Depot to Distribution Centers

The first level handles product distribution from a single authorized depot to each distribution center, taking into account the ordered quantities and vehicle capacity constraints. The underlying formulation is based on a vehicle routing structure and incorporates several elements essential for modeling this stage.

Sets:

  1. N = CD ∪ D: Network composed of distribution centers (CD) and depot (D).

  2. K: Set of vehicles assigned to distribution-center replenishment.

  3. P: Set of products.

Parameters:

  1. distij: Distance between nodes i and j in the distribution network.

  2. CVp: Capacity of the delivery vehicle for product p ∈ P.

  3. djp: Demand distribution center j for product p.

Decision Variables:

  1. Xij: Binary variable equal to 1 if the route between nodes i and j is selected, and 0 otherwise.

  2. Yijp: Quantity of product p transported on route (i, j).

  3. ui: Auxiliary variable used for subtour elimination.

M i n : ∑ i , j ∈ N : i ≠ j d i s t i j X i j (1)

subject to ∑ i ∈ N : i ≠ j X i j = 1 , ∀ j ∈ C D (2)

∑ i ∈ N : i ≠ j X j i = 1 , ∀ j ∈ C D (3)

∑ i ∈ D X i j = ∣ K ∣ , ∀ j ∈ C D (4)

∑ i ∈ D X j i = ∣ K ∣ , ∀ j ∈ C D (5)

∑ i ∈ N : i ≠ j Y i j p = d j p , ∀ j ∈ C D ; p ∈ P (6)

Y i j p ≤ C V p X i j , ∀ i , j ∈ N ; p ∈ P (7)

u i - u j + ∣ C D ∣ X i j ≤ ∣ C D ∣ - 1 , ∀ i , j ∈ C D ; i ≠ j (8)

The objective function (1) minimizes the total traveled distance ensuring efficient deployment of the available vehicles. The constraints (2) and (3) guarantee that each distribution center is visited exactly once. Equations (4) and (5) guarantee that the number of vehicles departing from the depot is equal to the number returning, i.e., control the number of vehicles to be used.

Constraint (6) ensures that the total amount of products delivered to each distribution center matches its corresponding demand. Constraint (7) establishes that the quantity of products transported along each selected route (i, j) is feasible as long as vehicle capacity is not exceeded, thereby creating a linkage between the routing decision variable and the quantity of products transported. Finally, constraint (8) prevents the formation of subtours by forcing the variable u to maintain a consistent ordering of the visited nodes, ensuring that each route forms a single continuous tour that starts and ends at the depot.

3.1.2 Level II: Routing from Distribution Centers to Customers within Coverage Areas

The second level deals with the delivery of products to customers located inside predetermined coverage regions. Each distribution center serves as the origin and endpoint of a route executed by a smaller urban delivery vehicle. The objective is to determine an efficient route that satisfies customer demands while respecting their assigned time windows.

Customer satisfaction is quantified as equation (9).

s a t i s f a c t i o n = s a t i s f i e d c u s t o m e r s t o t a l c u s t o m e r s × 100 (9)

where a customer is considered satisfied if service occurs within the time window [a, b] assigned by the distribution center, i.e,

a ≤ s e r v i c e t i m e ≤ b (10)

As customers j are served, the variables delayed, T j , and delivered, E j quantify the delay in supplying the product and the successful completion of the delivery, respectively. A penalty is applied when the delivery is not completed within the expected conditions. In addition, two other objectives are addressed: quantifying the cost associated with implementing the proposal approach and evaluating its overall operational performance.

t o t a l c o s t = o p e r a t i n g c o s t p e r k m × d i s t i , j X i , j , ∀ i , j ∈ N N (11)

The amount of CO2 generated by the urban vehicles responsible for delivering products to customers was also quantified. For this purpose, both the fuel consumption rate and the CO 2 emission per liter of fuel were considered. The corresponding calculation is given by the following expression:

f u e l c o n s u m p t i o n r a t e × C O 2 e m i s s i o n × d i s t i j × X i j q k (12)

where q k represents the carrying capacity of the urban vehicle.

The proposed model consists of the following notations:

Sets:

  1. CD: Set of Distribution or sales centers.

  2. C: Set of customers that generate a demand to be served.

  3. NN = CD ∪ C: Network composed of distribution centers and customers within a defined coverage area.

  4. KK: Urban vehicles responsible for transporting products to be delivered to customers.

Parameters:

  1. distij: Distance between distribution centers and customers.

  2. pc: Factor that weights the impact of operating costs, delivery delays, and CO2 emissions.

  3. ckm: Operating cost per kilometer traveled by the urban vehicle.

  4. ch: Operating cost per hour of work performed by urban vehicles.

  5. qk: Urban vehicle capacity.

  6. ts: Service time.

  7. pr: Penalty applied for each hour of delivery delay.

  8. pne: Penalty for failure to complete product delivery to customer j.

  9. jl: Maximum allowable working hours.

  10. em: CO2 emission per liter.

  11. tc: Fuel consumption rate.

Decision variables:

  1. Xij: Binary decision variable equal to 1 if distribution center i serves customer j, and 0 otherwise.

  2. ui: Integer variable that specifies the sequence in which customers are visited, ensuring the elimination of subtours.

  3. Tj: Variable representing the delivery delay for customer j ∈ C.

  4. Ej: Binary decision variable equal to 1 if the delivery to customer j is successfully completed, and 0 otherwise, triggering the corresponding penalty.

M i n : ∑ i , j ∈ N N : i ≠ j d i s t i j X i j + p c × ∑ i , j ∈ N N c k m × d i s t i j X i j + ∑ i , j ∈ N N c h × d i s t i j / q k + t s / 60 × X i j (13)

M i n : p r ∑ j ∈ C T j + p n e ∑ j ∈ C 1 - E j (14)

M i n : ∑ i , j ∈ N N , k ∈ K K e m × t c × d i s t i j × X i j / q k (15)

subject to ∑ i ∈ N N : i ≠ j X i j = 1 , ∀ j ∈ N N (16)

∑ i ∈ N N : i ≠ j X j i = 1 , ∀ j ∈ N N (17)

u i - u j + ∣ N N ∣ X i j ≤ ∣ N N ∣ - 1 , ∀ i , j ∈ C , i ≠ j (18)

∑ i , j ∈ N N : i ≠ j d i s t i , j / 30 + t s / 60 X i , j ≤ j l (19)

T j ≥ ∑ i , j ∈ N N : i ≠ j d i s t i , j / 30 + t s / 60 X i , j - b , ∀ j ∈ C (20)

E j b - a ≥ b - a - T j , ∀ j ∈ C (21)

The objective function (13) is composed of several summation blocks. The first block minimizes the traveled distance according to the route defined by the variable X ij . The second block optimizes the operating costs associated with the distance traveled in kilometers, while the third block accounts for the operating costs related to service time. In equation (14), the two summation blocks represent the penalties incurred for delivery delays and for incomplete deliveries, respectively. Finally, equation (15) quantifies CO2 emissions based on fuel consumption along the selected route, as determined by X ij , assuming that the capacity of the urban vehicle is known.

Constraints (16) and (17) ensure that each customer is visited exactly once. Constraint (18) eliminates the possibility of subtours in the construction of the optimal routing plan. Constraint (19) imposes the working-day limit, including service time; that is, it ensures that the total vehicle operating time does not exceed the maximum allowable workday. When multiplied by the routing variable X ij , this constraint guarantees that only the selected routes are considered in the calculation.

Condition (20) represents the delivery delay experienced by the customer; in other words, it measures the extent to which the delivery occurs after the upper bound b of the assigned time window. Finally, condition (21) ensures that a delivery is considered complete only if it is performed within the customer’s specified time window.

3.2 Simulated Annealing Metaheuristic

3.2.1 Level I: Multi-Product Routing from Depot to Distribution Centers

Algorithm 1 applies a Simulated Annealing (SA) metaheuristic to solve the Level I routing problem, where a depot must deliver multiple products to several distribution centers. The goal is to minimize the total travel distance while ensuring that vehicle capacity is not exceeded. The algorithm starts by generating an initial feasible solution, assigning each distribution center to a vehicle and forming routes that begin and end at the depot. The cost of this initial solution is calculated by summing all route distances and adding a penalty whenever a route exceeds the vehicle’s capacity.

Algorithm 1
Level I: Depot to Distribution Centers based on Simulated Annealing

Next, the algorithm improves the solution using SA. It creates neighboring solutions by applying simple modification operators such as relocating a center to another route, swapping centers between routes, or improving the internal order of visits through a 2-opt move. Each modified solution is evaluated, and it replaces the current solution if it has a lower cost. If the new solution is worse, it may still be accepted with a probability that decreases as the temperature parameter cools. This mechanism helps the search escape local minima.

The temperature gradually decreases according to a cooling schedule until it reaches a minimum threshold. At the end of the SA process, a final 2-opt local improvement is applied to each route to remove unnecessary detours. A detailed description can be found in Table 1.

Table 1
Simulated Annealing Parameter Configuration for Level I

3.2.2 Level II: Routing from Distribution Centers to Customers within Coverage Areas

Algorithm 2 applies a SA metaheuristic to the Level II routing problem, where each distribution center must deliver products to customers located within its coverage area. The algorithm seeks to generate efficient delivery routes that minimize travel distance while also measuring CO2 emissions, delivered products, and customer satisfaction.

Algorithm 2
Level II: Distribution Centers to Customers based on Simulated Annealing

The procedure begins by creating the customer set for each distribution center, assigning random coordinates, time windows, product demands, and service times. Using these data, all pairwise distances are computed through Euclidean metrics. For each distribution center, an initial route is constructed by arranging its assigned customers in a random visiting order.

Simulated Annealing is then used to improve this route. In each iteration, a neighbor route is generated by swapping the positions of two customers. The new route is accepted if it shortens the total travel distance or, if worse, with a probability determined by the current temperature, allowing the search to escape local minima. The temperature decreases gradually according to a cooling schedule until the algorithm completes the specified number of iterations. For each final route, the algorithm computes key performance indicators: total distance traveled, CO2 emissions based on fuel consumption models, total delivered products, and the percentage of customers served within their time windows, which represents the satisfaction level. A detailed description of the parameters of the SA can be found in Table 2.

Table 2
Simulated Annealing Parameter Configuration for Level II.

4 RESULTS AND DISCUSSION

This section presents and analyzes the computational results obtained from implementing both the optimization model and the Simulated Annealing metaheuristic. The optimization model was formulated using the PuLP library, versions 2.8.0, and solved with solver GLPK versions 5.0 in Python, while the metaheuristic was likewise implemented in Python 3.12.3, incorporating the required numerical and algorithmic libraries such as numpy, networkx, matplotlib, random, math and time. All experiments were conducted on a 64-bit HP laptop equipped with an Intel(R) Core(TM) i7-8550U CPU (1.80-2.00 GHz) and 8 GB of RAM, operating under Windows 10 Pro.

4.1 Results Using the Optimization Model

This Case evaluates the performance of the proposed model within a two-level distribution structure composed of a single depot, five distribution centers, three product categories, and three vehicles at Level I. See Table 3. Level II incorporates 26 customers, each assigned a time window defined by the corresponding distribution center. The time windows are random, because after the demand generated by the customer, the distribution center will define the best time to serve the customer through delivery. Demands are randomly generated because it varies over time; it is the customer who decides which products to buy, based on their daily interests. The deliveries are carried out by five urban vehicles, each with a capacity of 15 weight units.

Table 3
Depot location and distribution centers for level I.

The results obtained from Level I, summarized in Table 4 and Figure 3, show the complete distribution of the three product types to the designated distribution centers. The total travel distance required to supply all centers amounts to 319.51 kilometers.

Table 4
Level I: Depot → Distribution Centers using optimization model.

Figure 3
Level I: Vehicle routing from depot to distribution centers A, B, C, D and E identified with blue nodes including the distance between them.

Once the distribution centers are replenished, Level II proceeds with routing to customers. The results for each distribution center -including travel distances (in kilometers), delivered product quantities, service time compliance, (Hour: Minutes), CO2 emissions (in kilograms), and total costs-are shown in Tables 5 through 9 and Figures 4, 5, and 6. These findings reflect the vehicle routes generated for each distribution center, the percentage of satisfied customers, and the environmental impact associated with each route. In all cases, the model successfully assigns customers to their nearest center and constructs feasible delivery tours that adhere to vehicle capacities and time window constraints. Example 11 : 00 ∈ [11 : 00; 18 : 00]. In addition, the CO2 emissions generated by the vehicle assigned to each coverage area.

Table 5
Level II: Distribution Center A → Customers.

Table 6
Level II: Distribution Center B → Customers.

Table 7
Level II: Distribution Center C → Customers.

Table 8
Level II: Distribution Center D → Customers.

Table 9
Level II: Distribution Center E → Customers.

Figure 4
Product delivery routes to customers (green nodes) from Centers A and, including the distance between them.

Figure 5
Delivering products: Center C to customers(green nodes), including the distance between them.

Figure 6
Product delivery routes to customers(green nodes) from Centers D and E, including the distance between them.

Distribution Center A, see Table 5, exhibits a balanced routing structure characterized by moderate travel distances between consecutive customers and stable service times within the specified time windows. The resulting 100% satisfaction level indicates that the geographical clustering of customers around this center enables efficient route sequencing and minimal service delays. Furthermore, the relatively low CO2 emissions suggest that the spatial density of customers contributes positively to environmental efficiency, since shorter inter-customer distances reduce fuel consumption. The routing for serving customers is shown in Figure 4 (a).

Distribution Center B demonstrates one of the most efficient configurations in terms of environmental performance, producing the lowest CO2 emissions among all centers. This outcome can be attributed to the compact spatial arrangement of customers assigned to this center, allowing the vehicle to complete deliveries within short travel intervals. The reduced total distance traveled also contributes to a lower total operating cost, indicating that this center plays a key role in enhancing both economic and environmental sustainability within the proposed framework. See Table 6 and the routing for serving customers in the Figure 4 (b).

Distribution Center C presents a different routing pattern characterized by longer internode distances compared to Centers A and B. See Table 7. Although customer satisfaction remains at 100%, the increased travel distances lead to higher CO2 emissions and operational costs. This suggests that the spatial dispersion of customers assigned to Center C generates less efficient routing conditions. Nevertheless, the model successfully constructs feasible routes that maintain service reliability despite the greater geographic spread of demand points. The routing for serving customers is shown in Figure 5.

Distribution Center D shows relatively higher travel distances between certain customer pairs, which increases total route length and associated emissions. See Table 8. This behavior reflects the sensitivity of last-mile performance to the spatial configuration of customers within the coverage area. When customers are less geographically clustered, the routing algorithm must allocate longer travel segments, which directly impacts cost and environmental indicators. Despite these challenges, the model ensures compliance with service time constraints, maintaining acceptable satisfaction levels. The routing for serving customers is shown in Figure 6 (a).

Distribution Center E exhibits a routing structure comparable to that of Centers A and B, where the distribution of customers within the coverage radius allows efficient sequencing of deliveries. The resulting performance indicators show a favorable balance between cost, emissions, and service quality, reinforcing the importance of appropriately defining coverage areas when designing urban logistics systems. See Table 9. The routing for serving customers is shown in Figure 6 (b).

Overall, the individual analysis of distribution centers highlights the importance of spatial demand distribution in determining last-mile performance. Centers serving geographically compact customer clusters tend to generate lower travel distances, reduced CO2 emissions, and lower operating costs, while dispersed customer configurations increase route complexity and environmental impact. These findings confirm that the effectiveness of the proposed two-level optimization framework depends not only on the routing algorithm itself but also on the strategic positioning of distribution centers and the definition of their coverage areas.

This analysis also illustrates how the proposed model can support decision-makers in evaluating alternative configurations of distribution networks. By identifying centers associated with higher operational costs or environmental impact, planners can redesign coverage areas, relocate distribution facilities, or adjust fleet allocation strategies in order to improve system-wide efficiency and sustainability.

Table 10 presents a comprehensive summary of the results obtained using the optimization model for Level II, which were presented in Tables 5 through 9.

Table 10
Level II: Distribution Center → Customers using Optimization model.

4.2 Results Using the Simulated Annealing Metaheuristic

To evaluate the performance of the proposed Simulated Annealing algorithm (Algorithms 1 and 2), the computational experiments, for level I and II, were performed using the same information from subsection 4.1, including Table 3. The results, summarized in Table 11, where total travel distance was 319.52 kilometers and Table 12, demonstrate the metaheuristic’s ability to generate high-quality solutions for each distribution center efficiently.

Table 11
Level I: Depot → Distribution Centers using Simulated Annealing.

Table 12
Level II: Distribution Center → Customers using Simulated Annealing.

Describing in detail, Table 12, a 100% customer satisfaction level was achieved, as all customers were served within their assigned time windows. Moreover, the total distance traveled remained competitive, indicating that the applied metaheuristic was capable of identifying efficient routes despite the lack of guarantees of global optimality. From an operational perspective, all delivered products fully met customer demands. Finally, CO2 emissions remained within reasonable ranges and were consistent with the distances traveled, reflecting an environmentally aligned performance in accordance with the observed logistical efficiency.

4.3 Comparative analysis between optimization model and Simulated Annealing

The research proposals are evaluated at Level II, as this stage incorporates a greater volume of data into the computational process associated with last-mile delivery operations, including five distribution centers and up to 118 geographically dispersed customers. Accordingly, a comparative analysis between the results derived from the exact optimization model (Table 10) and those obtained through the simulated annealing (SA) metaheuristic (Table 12) yields significant insights into the trade-offs among solution optimality, service quality, and practical applicability within last-mile distribution systems.

4.3.1 Customer satisfaction: robustness vs optimality trade-off

One of the most significant differences between both approaches lies in customer satisfaction. While the optimization model achieves high satisfaction levels in most distribution centers, it fails to guarantee full compliance in all cases. In particular, Distribution Center D reaches only 80% satisfaction, indicating that some deliveries occur outside the assigned time windows. In contrast, the Simulated Annealing approach achieves 100% customer satisfaction across all distribution centers, ensuring that all deliveries are completed within their respective time windows.

This result highlights a critical practical advantage of the metaheuristic: The exact model prioritizes global optimality but may sacrifice service feasibility under tight constraints. The meta-heuristic, by contrast, implicitly prioritizes feasibility and adaptability, making it more suitable for real-world operations where service reliability is critical.

4.3.2 Total distance and operational efficiency

From the perspective of total traveled distance, the optimization model generally produces more efficient (shorter) routes, particularly in some centers:

  • Center B: 32.38 km (Optimization model) vs ∼ 42 - 51 km (SA)

  • Center C: 44.59 km (Optimization model) vs ∼ 46 - 52 km (SA)

This confirms that the exact model is more effective at minimizing distance, as expected from a mathematical programming formulation. However, the differences are not extreme, and in some cases (e.g., Center A and E), distances are nearly identical between both methods. This suggests that the metaheuristic is capable of generating near-optimal solutions in terms of routing efficiency.

4.3.3 Operational cost comparison

A similar pattern is observed in total costs: The optimization model yields lower operational costs, especially in centers with compact routing structures (e.g., B and C). The metaheuristic produces slightly higher costs due to longer routes and its stochastic nature.

Nevertheless, the increase in cost is moderate and can be interpreted as the “price” of achieving full customer satisfaction and improved feasibility under time constraints.

4.3.4 Environmental impact

CO2 emissions follow the same trend as distance: The optimization model achieves lower emissions in most centers (e.g., B: 0.25 vs ∼ 0.33 - 0.40). The metaheuristic results in slightly higher emissions due to increased travel distance.

However, emissions remain within comparable ranges, indicating that the environmental degradation introduced by the metaheuristic is limited and controlled.

4.3.5 Delivered products and service completeness

Both approaches successfully meet customer demand: Delivered products are consistent across both methods for most centers. This confirms that both approaches maintain high service completeness, ensuring operational reliability.

4.3.6 Practical implications for real-world implementation

From a practical standpoint, the SA approach demonstrates greater applicability in real logistics systems, due to: Its ability to guarantee full service compliance (100% satisfaction); its flexibility in handling complex and constrained environments and its capacity to generate high-quality solutions without excessive computational burden.

In contrast, the optimization model is more suitable as a benchmark or planning tool, providing: lower bounds for cost and distance.

4.4 Results with other scenarios

The critical phase of the proposals occurs at Level II due to the use of vehicle routing for last-mile delivery; therefore, to identify the strengths and weaknesses in large-scale scenarios, the proposals were tested across 13 scenarios, each of which specifies the number of customers served by each distribution center.

The comparative analysis between the results of the optimization model (Table 13) and the Simulated Annealing metaheuristic (Table 14), under identical scenarios and customer datasets, reveals clear differences in terms of efficiency, service quality, and scalability. First, the optimization model consistently achieves shorter travel distances, lower operational costs, and reduced CO2 emissions across all scenarios, with gaps that widen significantly as the problem size increases; for instance, in Scenario 13, the distance obtained by the metaheuristic is more than double that of the exact model. However, this superior efficiency is achieved at the expense of service quality, as customer satisfaction in the optimization model progressively declines from 96.15% to 26.27%, reflecting increasing difficulty in meeting time window constraints in larger instances. In contrast, Simulated Annealing ensures 100% customer satisfaction in all scenarios, prioritizing feasibility and service compliance. Furthermore, although the metaheuristic incurs higher costs and emissions due to less efficient routing, its main advantage lies in its computational stability, as processing time remains nearly constant even for large-scale scenarios, whereas the optimization model exhibits exponential growth, becoming impractical (e.g., 127.347 minutes in Scenario 13). Overall, these results highlight a fundamental trade-off between efficiency and robustness: the exact model is more suitable as a theoretical benchmark and for smaller problem instances, while the metaheuristic provides a more practical and scalable solution for real-world applications where customer satisfaction and response time are critical.

Table 13
Level II: Results for customers grouped in 13 scenarios using Model Optimization.

Table 14
Level II: Results for customers grouped in 13 scenarios using SA.

5 CONCLUSIONS

This study proposed a two-level optimization framework to address the complexities of last-mile delivery in urban environments. The first level focused on the replenishment of distribution centers from a central depot, while the second level modeled the delivery of products to customers located within predefined coverage areas. Both levels were formulated through binary integer linear programming and complemented by a Simulated Annealing metaheuristic designed to efficiently handle larger problem instances.

The computational results demonstrate that the proposed approach effectively reduces operational distances and associated costs while improving environmental performance through lower CO2 emissions. The routing plan generated for the second level makes it possible to quantify service quality through customer satisfaction, enabling a more comprehensive assessment of delivery performance. The metaheuristic proved particularly efficient for larger customer sets, consistently achieving high-quality solutions-including full satisfaction of customer time windows-within very short computational times.

Overall, the findings highlight the potential of optimization models and metaheuristic strategies to support decision-making in last-mile logistics. The proposed framework offers a balanced perspective on three key dimensions of urban distribution systems: economic efficiency, environmental impact, and customer satisfaction.

Future research may extend this work by incorporating additional real-world features such as heterogeneous fleets, electric or autonomous vehicles, dynamic travel times, or stochastic customer demand. Another promising direction involves integrating machine-learning-based forecasting models to further enhance routing accuracy and responsiveness. Finally, applying more advanced or hybrid metaheuristics to large-scale urban delivery networks may lead to improved scalability and robustness, strengthening the practical relevance of this framework.

Data Availability

All the data investigated are provided in the article.

References

  • ALJOHANI K. 2024. The role of Last-Mile Delivery Quality and Satisfaction in Online Retail Experience: An Empirical Analysis. Sustainability, MDPI, 16(11): 4743.
  • BORESTA M, PINTO DM & STECCA G. 2024. Bridging operations research and machine learning for service cost prediction in logistics and service industries. Annals of operations research, Springer link, 342: 113-139.
  • BOYSEN N, FEDTKE S & SCHWERDFEGER S. 2021. Last-mile delivery concepts: a survey from an operational research perspective. OR Spectrum, 43: 1-58.
  • BROCHADO AF, ROCHA EM, ADDO E & SILVA S. 2024. Performance evaluation and explainability of last-mile delivery. Procedia computer science, ScienceDirect, Elsevier, 232: 2478-2487.
  • GIUFFRIDA N, FAJARDO-CALDERIN J, MASEGOSA AD, WERNER F, STEUDER M & PILLA F. 2023. Optimization and machine learning applied to last-mile logistics: A review. Sustainability, MDPI , 14(9): 5329.
  • JAZEMI R, ALIDADIANI E, AHN K & JANG J. 2023. A review of literature on vehicle routing problems of last-mile delivery in urban areas. Applied sciences, MDPI, 13(24): 13015.
  • LI Y & YANG J. 2024. The last-mile delivery vehicle routing problem with handling cost in the front warehouse mode. Computers & Industrial Engineering, ScienceDirect, Elsevier, 190: 110076.
  • MORADI M, SADATI I & CATAY B. 2023. Last mile delivery routing problem using autonomous electric vehicle. Computers & Industrial Engineering, ScienceDirect , Elsevier, 184: 109552.
  • ÖZARIK SS, VEELENTURF LP, WOENSEL TV & LAPORTE G. 2021. Optimizing e-commerce last-mile vehicle routing and scheduling under uncertain customer presence. Transportation Research Part E: Logistics and Transportation Review, ScienceDirect, Elsevier, 148: 102263.
  • PERES F & CASTELLI M. 2021. Combinatorial optimization problems and metaheur´ıstics: Review, challenges, design, and development. Applied sciences, MDPI , 11(14): 6449.
  • SABET S & FAROOQ B. 2022. Green vehicle routing problem: State of the art and future directions. arXiv, Optimization and Control (math.OC), .
  • SENNA F, COELHO LC, MORABITO R & MUNARI P. 2024. An exact method for a last-mile delivery routing problem with multiple deliverymen. European Journal of Operational Research, Elsevier, 317(2024): 550-562.
  • STAMADIANOS T, KYRIANKAKIS NN, MARINAKI M & MARINAKIS Y. 2023. Routing problems with electric and autonomous vehicles: Review and potential for future research. Operatioms research forum, Springer link, 4(46).
  • TAUSIF IH. 2023. Las mile delivery: optimization model for drone-enabled vehicle routing problem. Emerging Minds Journal for Student Research, 1: 39-73.
  • TIWARI KV & SHARMA SK. 2023. An optimization model for vehicle routing problem in last-mile delivery. Expert Systems with Applications, ScienceDirect, Elsevier, 222: 119789.
  • ZHU X, CAI L, LAI PL, WANG X & MA F. 2023. Evolution, challenges, and opportunities of transportation methods in the last-mile delivery process. Systems, MDPI, 11(10): 509.
  • ZIJM WH, KLUMPP M, REGATTIERI A & HERAGU SSE. 2019. Operations, Logistics and Supply Chain Management. Springer International Publishing AG, part of Springer Nature 2019.

Funding

There is no funding for this research.

Editor responsible for the review

Editor-in-Chief: Antônio Augusto Chaves.

Publication Dates

  • Publication in this collection
    24 Aug 2026
  • Date of issue
    2026

History

  • Received
    13 Dec 2025
  • Accepted
    24 June 2026
location_on
Sociedade Brasileira de Pesquisa Operacional Rua Mayrink Veiga, 32 - sala 601 - Centro, 20090-050 , Tel.: +55 21 2263-0499 - Rio de Janeiro - RJ - Brazil
E-mail: sobrapo@sobrapo.org.br
rss_feed Acompañe los números de esta revista en su lector de RSS
Ir para arriba Notificar error