A Metaheuristic Algorithm for the Fleet
Size and Mix Vehicle Routing Problem
with Pickups and Deliveries
Developing an Approach for Finding Solutions to Realistic
Large-Scale Instances

Master’s Thesis in
Engineering Mathematics and Computational Science
& Computer Science - Algorithms, Languages and Logic

FATEMEH ALIBAKHSHI & PAWEL KASINSKI

DEPARTMENT OF ELECTRICAL ENGINEERING

CHALMERS UNIVERSITY OF TECHNOLOGY
Gothenburg, Sweden 2025
www.chalmers.se

www.chalmers.se




Master’s thesis 2025

A Metaheuristic Algorithm for the Fleet Size and
Mix Vehicle Routing Problem with Pickups and

Deliveries

Developing an Approach for Finding Solutions to Realistic
Large-Scale Instances

FATEMEH ALIBAKHSHI
& PAWEL KASINSKI

Department of Electrical Engineering
Division of Systems and Control

Chalmers University of Technology
Gothenburg, Sweden 2025



A Metaheuristic Algorithm for the Fleet Size and Mix Vehicle Routing Problem
with Pickups and Deliveries
Developing an Approach for Finding Solutions to Realistic Large-Scale Instances
FATEMEH ALIBAKHSHI & PAWEL KASINSKI

© FATEMEH ALIBAKHSHI & PAWEL KASINSKI, 2025.

Supervisors: Sunney Fotedar & Toheed Ghandriz, Volvo Group Trucks Technology
Academic Supervisor: Jiaming Wu, Dept. of Architecture and Civil Engineering
Examiner: Balázs Adam Kulcsár, Dept. of Electrical Engineering

Master’s Thesis 2025
Volvo Group Trucks Technology
Department of Electrical Engineering
Division of Systems and Control
Chalmers University of Technology
SE-412 96 Gothenburg
Telephone +46 31 772 1000

Cover: A visualization of a solution to the problem where a time series is shown
for each used vehicle. The dark gray intervals indicate traveling, the blue and red
intervals indicate loading and unloading of cargo respectively, and the green ones
indicate charging. The vehicles 1 and 2 are part of the fleet but have not been used.

Printed by Chalmers Reproservice
Gothenburg, Sweden 2025

iv



A Metaheuristic Algorithm for the Fleet Size and Mix Vehicle Routing Problem
with Pickups and Deliveries
Developing an Approach for Finding Solutions to Realistic Large-Scale Instances
FATEMEH ALIBAKHSHI & PAWEL KASINSKI
Department of Electrical Engineering
Chalmers University of Technology

Abstract
This thesis addresses an extended variant of the fleet size and mix vehicle routing
problem with aspects such as split pickups and deliveries, different loading and un-
loading methods, and the charging of electric vehicles. The aim of the work has been
to be able to efficiently solve large-scale instances of the problem based on real data.
To tackle this objective, a metaheuristic algorithm has been developed based on
the adaptive large neighborhood search framework. A procedure for destroying and
repairing a solution has been designed and implemented such that the challenging
aspects of the problem definition are handled appropriately. Moreover, the use of
simulated annealing and the technique of applying noise to diversify the search have
been explored. The algorithm was run on three realistic instances, in a set of shorter
runs to investigate the effect of adding noise and using simulated annealing. Then,
three relatively long runs were conducted to obtain definitive solutions to the in-
stances. The algorithm managed to find a promising solutions to all three instances,
and the objective value of the best solutions found were significantly lower than the
objective values of the initial solutions.

Keywords: vehicle routing problem, pickup and delivery, heterogeneous fleet, split
pickups and deliveries, charging, alns, adaptive large neighborhood search, heuristic,
operations research, logistics

v





Acknowledgements
First and foremost, we would like to thank Sunney Fotedar. Thank you for showing
great interest in the project, always finding time for discussions with us, and for
supporting us from the beginning to the end. Your input has truly helped us shape
this master thesis.

We would also like to thank Jiaming Wu; for sharing insights about his experience
with the ALNS framework, and for giving us advice regarding both the project and
how the master thesis should be written.

Furthermore, we would like to express our gratitude to Toheed Ghandriz for propos-
ing this project and introducing us to it, and finally, we would like to thank Balázs
Adam Kulcsár for keeping in touch with us throughout the project and for making
the examination process feel straightforward.

Fatemeh Alibakhshi & Pawel Kasinski, Gothenburg, June 2025

vii





List of Acronyms

Below is the list of acronyms used in this thesis, listed in alphabetical order.

Acronym Definition

ALNS Adaptive Large Neighborhood Search
AST Additional Semitrailer
BEV Battery Electric Vehicle
EV Electric Vehicle
FSMVRP Fleet Size and Mix Vehicle Routing Problem
GA Genetic Algorithm
LNS Large Neighborhood Search
MILP Mixed-Integer Linear Programming
OBL On-Board Lift
OBW On-Board Waiting
PDPTW Pickup and Delivery Problem with Time Windows
SA Simulated Annealing
SC Straddle Carrier
SoC State of Charge
TS Tabu Search
VNS Variable Neighborhood Search
VRP Vehicle Routing Problem

ix



x



Contents

List of Acronyms ix

List of Figures xv

List of Tables xvii

1 Introduction 1
1.1 The Vehicle Routing Problem . . . . . . . . . . . . . . . . . . . . . . 1
1.2 Overview of the Problem . . . . . . . . . . . . . . . . . . . . . . . . . 2
1.3 Context . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
1.4 Aim . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3
1.5 Delimitations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3

2 Literature Review 5
2.1 Classical and Metaheuristic Approaches to the Vehicle Routing Problem 5
2.2 Energy-Aware and Electric Vehicle Routing Problems . . . . . . . . . 5
2.3 Machine Learning for Vehicle Routing Problems . . . . . . . . . . . . 6
2.4 Comparative Analysis of Solution Methods . . . . . . . . . . . . . . . 7
2.5 Research Gaps and Summary . . . . . . . . . . . . . . . . . . . . . . 7

3 Theoretical Foundations 9
3.1 Vehicle Routing Problem . . . . . . . . . . . . . . . . . . . . . . . . . 9

3.1.1 The Complexity of the Vehicle Routing Problem . . . . . . . . 11
3.2 Large Neighborhood Search . . . . . . . . . . . . . . . . . . . . . . . 11

3.2.1 Fundamental Concepts . . . . . . . . . . . . . . . . . . . . . . 11
3.2.2 LNS Algorithmic Framework . . . . . . . . . . . . . . . . . . . 12
3.2.3 Acceptance Criterion . . . . . . . . . . . . . . . . . . . . . . . 13
3.2.4 Design of Destroy and Repair Operators . . . . . . . . . . . . 14

3.3 Adaptive Large Neighborhood Search . . . . . . . . . . . . . . . . . . 14
3.3.1 Motivation and Overview . . . . . . . . . . . . . . . . . . . . 15
3.3.2 ALNS Algorithmic Framework . . . . . . . . . . . . . . . . . . 15
3.3.3 Destroy and Repair Operator Selection . . . . . . . . . . . . . 17
3.3.4 Acceptance Criterion . . . . . . . . . . . . . . . . . . . . . . . 18
3.3.5 Design of Destroy and Repair Operators . . . . . . . . . . . . 18
3.3.6 Properties and Strengths of ALNS . . . . . . . . . . . . . . . 18

3.4 Dijkstra’s Algorithm . . . . . . . . . . . . . . . . . . . . . . . . . . . 20

xi



Contents

4 Mathematical Model 23
4.1 Terminology . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23
4.2 Problem Description . . . . . . . . . . . . . . . . . . . . . . . . . . . 24
4.3 Sets, Parameters, and Decision Variables . . . . . . . . . . . . . . . . 25

4.3.1 Sets . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
4.3.2 Decision Variables . . . . . . . . . . . . . . . . . . . . . . . . 26
4.3.3 Parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26

4.4 Objective Function of the Model . . . . . . . . . . . . . . . . . . . . . 28
4.5 Model Constraints . . . . . . . . . . . . . . . . . . . . . . . . . . . . 29

4.5.1 Routing Constraints . . . . . . . . . . . . . . . . . . . . . . . 29
4.5.2 Loading and Unloading Constraints . . . . . . . . . . . . . . . 29
4.5.3 Battery Energy Constraints . . . . . . . . . . . . . . . . . . . 31
4.5.4 Time Constraints . . . . . . . . . . . . . . . . . . . . . . . . . 32
4.5.5 Non-Negativity and Binary Constraints . . . . . . . . . . . . . 35

5 Implementation 37
5.1 Rationale for Choosing the ALNS Framework . . . . . . . . . . . . . 37
5.2 Solution Representation . . . . . . . . . . . . . . . . . . . . . . . . . 38
5.3 Removal Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 39

5.3.1 Random Removal . . . . . . . . . . . . . . . . . . . . . . . . . 39
5.3.2 Worst Removal . . . . . . . . . . . . . . . . . . . . . . . . . . 40
5.3.3 Single Route Removal . . . . . . . . . . . . . . . . . . . . . . 40

5.4 Insertion Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 40
5.4.1 Greedy Insertion . . . . . . . . . . . . . . . . . . . . . . . . . 44
5.4.2 Regret-2 Insertion . . . . . . . . . . . . . . . . . . . . . . . . . 44
5.4.3 Regret-3 Insertion . . . . . . . . . . . . . . . . . . . . . . . . . 44
5.4.4 Restricted Greedy Insertions . . . . . . . . . . . . . . . . . . . 45

5.5 Checking Feasibility . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
5.6 Repairing Infeasible Routes by Adding Charging . . . . . . . . . . . . 46
5.7 Finding a Charging Path . . . . . . . . . . . . . . . . . . . . . . . . . 46
5.8 The ALNS Procedure . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
5.9 Applying Noise to Diversify the Search . . . . . . . . . . . . . . . . . 48
5.10 Simulated Annealing . . . . . . . . . . . . . . . . . . . . . . . . . . . 48

6 Results 51
6.1 Setup . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 51

6.1.1 Problem Instances . . . . . . . . . . . . . . . . . . . . . . . . 51
6.1.2 Fleet Composition . . . . . . . . . . . . . . . . . . . . . . . . 52
6.1.3 Base Algorithm Configuration . . . . . . . . . . . . . . . . . . 53
6.1.4 Algorithm Variants . . . . . . . . . . . . . . . . . . . . . . . . 54

6.2 Initial Solutions . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54
6.3 Shorter Algorithm Runs: Comparing Variants . . . . . . . . . . . . . 55
6.4 Longer Algorithm Runs: Searching for Better Solutions . . . . . . . . 62

7 Discussion 71
7.1 The Effect of Noise and Simulated Annealing . . . . . . . . . . . . . . 71
7.2 Gauging the Quality of the Solutions . . . . . . . . . . . . . . . . . . 72

xii



Contents

7.3 Analyzing the Vehicles Used . . . . . . . . . . . . . . . . . . . . . . . 72
7.4 Limitations of the Algorithm . . . . . . . . . . . . . . . . . . . . . . . 73
7.5 Potential for Future Work . . . . . . . . . . . . . . . . . . . . . . . . 74

8 Conclusion 75

Bibliography 77

xiii



Contents

xiv



List of Figures

1.1 Illustration of a Solution to a VRP Instance . . . . . . . . . . . . . . 2

3.1 Illustration of the LNS Process . . . . . . . . . . . . . . . . . . . . . 13
3.2 Illustration of ALNS Neighborhoods . . . . . . . . . . . . . . . . . . . 15
3.3 Comparison of LNS and ALNS . . . . . . . . . . . . . . . . . . . . . 20

6.1 Progression of the Search (Instance 01, Basic Variant, Random Seed=3,
500 Iterations) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57

6.2 Best Solution to Instance 01 from the Shorter Runs . . . . . . . . . . 57
6.3 Progression of the Search (Instance 02, SA Variant, Random Seed=3,

500 Iterations) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
6.4 Best Solution to Instance 02 from the Shorter Runs . . . . . . . . . . 59
6.5 Progression of the Search (Instance 03, SA Variant, Random Seed=3,

1000 Iterations) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 61
6.6 Best Solution to Instance 03 from the Shorter Runs . . . . . . . . . . 62
6.7 Progression of the Search (Instance 01, 10 000 Iterations) . . . . . . . 63
6.8 Best Solution to Instance 01 Found in the Longer Run . . . . . . . . 64
6.9 Distribution of the Weights after the Longer Run on Instance 01 . . . 64
6.10 Progression of the Search (Instance 02, 10 000 Iterations) . . . . . . . 65
6.11 Best Solution to Instance 02 Found in the Longer Run . . . . . . . . 66
6.12 Distribution of the Weights after the Longer Run on Instance 02 . . . 66
6.13 Progression of the Search (Instance 03, 10 000 Iterations) . . . . . . . 67
6.14 Best Solution to Instance 03 Found in the Longer Run . . . . . . . . 68
6.15 Distribution of the Weights after the Longer Run on Instance 03 . . . 68

xv



List of Figures

xvi



List of Tables

2.1 Optimization Techniques for VRPs . . . . . . . . . . . . . . . . . . . 7

4.1 Sets in the Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
4.2 Decision Variables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
4.3 Model Parameters . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27

5.1 Example of the Solution Representation . . . . . . . . . . . . . . . . 39

6.1 Details Regarding the Problem Instances . . . . . . . . . . . . . . . . 52
6.2 Details Regarding the Fleet . . . . . . . . . . . . . . . . . . . . . . . 53
6.3 Details Regarding the Initial Solutions . . . . . . . . . . . . . . . . . 55
6.4 Summary of the Shorter Runs on Instance 01 . . . . . . . . . . . . . 56
6.5 Summary of the Shorter Runs on Instance 02 . . . . . . . . . . . . . 58
6.6 Summary of the Shorter Runs on Instance 03 . . . . . . . . . . . . . 60
6.7 Summary of the Longer Runs on All Instances . . . . . . . . . . . . . 63

xvii



List of Tables

xviii



1
Introduction

Heavy-duty vehicles account for approximately two percent of vehicles on the road,
but contribute to over 25 percent of greenhouse gas emissions from road transport
in the European Union [1]. To reduce the climate footprint of this group of vehi-
cles, a proposed solution is electrification. However, as electric heavy-duty vehicles
are being commercialized, new challenges in infrastructure, routing, and investment
planning emerge for customers of such vehicles. As customers want to make sure
that their businesses will be profitable after investing in, for example, battery electric
trucks, it is desirable to be able to estimate the economic effects that the purchase
will have on their operations. Moreover, once a customer has extended or replaced
their fleet of trucks with battery electric trucks, they want to plan their deliveries
in an optimal way to minimize their costs while meeting the demands of their cus-
tomers. To address these problems, it is necessary to be able to solve a new variant
of the vehicle routing problem in an efficient manner.

1.1 The Vehicle Routing Problem

The vehicle routing problem (VRP) is an established problem, and determining an
optimal solution is known to be NP-hard [2]. In general, the problem is to find routes
originating from a depot for vehicles in a fleet such that goods are delivered to a set
of customers while costs are kept to a minimum. Since the problem first appeared
in an academic paper by George Dantzig and John Ramser in 1959, a multitude of
variants have been introduced and studied [3], [4]. One of the variants is referred to
as the fleet-size and mix-vehicle routing problem (FSMVRP). The differing aspects
of FSMVRP compared to VRP is that the fleet size, i.e. the number of vehicles
used, is not predetermined, and that the costs and capacities associated with the
vehicles can vary [5]. In this thesis, a more complex vehicle routing problem that
incorporates the aspects of FSMVRP is to be solved. Below, in Figure 1.1, a simple
instance of the vehicle routing problem and an example solution to it are shown.

1



1. Introduction

C1

C2 C3

C4

C5

C6

C7
C8

C9

C10

C11

C12

C13

C14

C15

C16
C17

C18

C19

C20

Figure 1.1: An illustration of a solution with four vehicle routes to a vehicle routing
problem instance with a central depot.

1.2 Overview of the Problem
In the variant of the VRP that will be treated in this master thesis, the fleet of
vehicles can include both diesel vehicles and electric vehicles. The vehicles can have
different mass capacities, volume capacities, and in the case of electric vehicles,
battery capacities. Similarly to FSMVRP, not all vehicles from the fleet have to
be used in the solution. The addition of electric vehicles is significant as it comes
with the added challenge of planning when and how to charge when searching for a
solution.

An instance of the problem is represented as a graph where the locations are nodes
and energy consumption values are arcs. Since the vehicles have varying specifica-
tions, they have differing energy values for the arcs. A node can be a depot, customer
location, or a standalone charging station. At a customer node, a vehicle can either
load or unload cargo, but a necessary condition is that the vehicle and the customer
have matching loading types. While loading or unloading, an electric vehicle may
also charge in parallel if possible at the considered customer. Furthermore, each
vehicle start from a designated home depot location, and has to end its route at this
location.

The problem is a pickup and delivery VRP variant, meaning that cargo is loaded at
customer locations instead of at the depot and is then unloaded at other customer
locations. A pair consisting of a pickup location and the corresponding delivery
location will be denoted as a request, and each request will be associated with the
amount of mass to be transported between these nodes. Moreover, since the mass

2



1. Introduction

associated with a request can exceed the mass capacity of a vehicle, the request may
need to be served by multiple vehicles carrying part of the cargo, or by the same
vehicle going back and forth. This aspect of the problem will be referred to as split
pickups and deliveries. At a later point in the thesis, in chapter 4, the problem is
defined more rigorously using mathematical notation.

1.3 Context
Prior to this master thesis, extended versions of FSMVRP have been studied, and
recently an approach using column generation was proposed for a variation that
shares many similarities with the problem that will be dealt with in this thesis
[6]. Column generation is an exact approach, and exact approaches are one of
the main directions where research effort is put in the field. However, the vehicle
routing problem cannot be solved by a polynomial-time algorithm and the time to
find an optimal solution scales exponentially with regard to instance size. Thus,
it is not feasible to find optimal solutions for large instances. Therefore, another
prominent research direction in the field regards inexact approaches that utilize, for
example, metaheuristics, matheuristics, or deep reinforcement learning [7], [8], [9].
These so called heuristic approaches have been shown to be effective for large-scale
instances of vehicle routing problems and have performed well in international VRP
competitions [10].

1.4 Aim
The aim of this master thesis is to study how a heuristic method can be applied
to solve the previously described extended variant of FSMVRP. More specifically,
the focus is to adapt the Adaptive Large Neighborhood Search (ALNS) framework
to be able to find good solutions for relatively large instances of the variant of the
problem, in a time efficient manner. As such, the research question for this thesis
can be formulated as follows:

- How to develop an algorithm for solving the defined problem such that its
complexities are handled and large-scale instances can be solved efficiently?

1.5 Delimitations
To reduce the scope of the project to some extent, a few assumptions have been
made. First, it is assumed that the pickups and deliveries will always be one to one
relations. In other words, if it is known where certain cargo has been picked up, it
is known where it should be delivered.

Second, it has been assumed that a vehicle can only load or unload at a node, and
never perform both actions.

3



1. Introduction

Third, it is assumed that it is possible to serve all requests in the considered problem
instances using the vehicles from the fleet. Thus, solutions where some customer
are unable to be served will not be accepted and consequently there is no need for a
request bank for storing such requests like the one proposed by Røpke and Pisinger
in the academic paper "An adaptive large neighborhood search heuristic for the
pickup and delivery problem with time windows" [11].

4



2
Literature Review

The Vehicle Routing Problem (VRP) is one of the most studied problems in trans-
portation and logistics. Over the years, researchers have developed many versions
of this problem to reflect more realistic conditions, such as time windows, split pick-
ups and deliveries, and heterogeneous vehicles. Recently, there has been growing
interest in solving VRPs with electric vehicles (EVs), due to the increasing demand
for more sustainable transportation, especially in urban areas.

2.1 Classical and Metaheuristic Approaches to the
Vehicle Routing Problem

Traditional methods for solving VRPs include classical heuristics like the nearest
neighbor algorithm and Clarke-Wright savings algorithm [12], [13]. These methods
are fast and simple, but they often produce low quality solutions when the problem
includes complex constraints.

To improve solution quality, many researchers turned to metaheuristics, such as
Genetic Algorithms (GA) [14], Tabu Search (TS) [15], Simulated Annealing (SA)
[16], and Variable Neighborhood Search (VNS) [17]. These methods are more flexible
and better at exploring the solution space, especially in large and complex problems.

Among these, Adaptive Large Neighborhood Search (ALNS) has become one of the
most popular algorithms for solving rich VRPs. It was introduced by Røpke and
Pisinger [11] and improves on the original Large Neighborhood Search by using an
adaptive mechanism to select the best operators during the search. ALNS has been
successfully applied to many versions of VRP and is known for its flexibility and
good performance.

2.2 Energy-Aware and Electric Vehicle Routing
Problems

As more companies adopt electric vehicles, routing models must also include energy
constraints, such as battery capacity and the need to stop at charging stations.

5



2. Literature Review

Several studies have explored how to include these aspects in VRP models.

For example, Hiermann et al. [18] studied a version of the problem where the fleet
includes both electric and diesel vehicles. Their model considers battery limits and
charging needs, while also respecting time windows and vehicle capacities. They
used a hybrid approach that combines ALNS with an exact Set Partitioning model.
This allows them to explore many possible routes using ALNS and then choose the
best ones using optimization.

Keskin and Çatay [19] also used a hybrid method for the electric VRP with time
windows. They focused on fast charging stations and used ALNS to generate routes
and a mathematical model to fine-tune charging decisions. Their results show that
combining metaheuristics with exact models can lead to better solutions.

Zhou et al. [20] added even more realism by including waiting times at charging
stations. Their model uses a two-step approach: first, it finds possible routes using
heuristics, and then it adjusts charging times using a queueing model. This is
important because, in real life, EVs often have to wait at busy charging stations.
Even short delays can affect delivery times and energy usage.

These studies all support the idea that routing and charging decisions should not be
made separately. They show that hybrid models, especially those using ALNS, are
well-suited for solving modern VRPs that involve energy constraints and complex
delivery requirements.

2.3 Machine Learning for Vehicle Routing Prob-
lems

In recent years, machine learning (ML) has started to play a role in solving VRPs.
ML can help by learning patterns from past data or by improving decision-making
in real time [21], [22].

Some studies use supervised learning methods like decision trees or neural networks
to predict travel time, energy use, or customer demand. For example, supervised
models have been used to forecast customer order volumes, estimate travel dura-
tions under uncertainty, and predict vehicle energy consumption [23], [24]. These
predictive models can support traditional algorithms or help reduce the search space,
leading to more efficient solutions.

More advanced approaches use reinforcement learning (RL), where the model learns
routing strategies by interacting with the environment. Deep reinforcement learning,
including methods like Pointer Networks [25] and attention-based models trained via
policy gradients [26], has been applied to VRP and related problems. These neural
combinatorial optimization methods aim to build solutions directly, without the need
for handcrafted heuristics, and have shown promising results in routing, scheduling,

6



2. Literature Review

and other discrete optimization tasks.

However, ML-based methods still face challenges. They often require large amounts
of training data, can be difficult to generalize to new instances, and may be hard to
scale to very large or complex problem settings [27]. Also, they are not yet widely
used in practical logistics environments, partly due to their complexity and limited
interpretability.

Although ML is not the main focus of this thesis, it is an important research direction
and could be combined with metaheuristics in future work. For example, ML models
could be used to guide the selection of operators in ALNS or to estimate energy use
during routing.

2.4 Comparative Analysis of Solution Methods
Table 2.1 provides a summary of the most common approaches used to solve VRPs,
comparing their main strengths and weaknesses. Each method has its own advan-
tages, depending on the problem size, constraints, and computational resources.
Classical heuristics are fast but often not suitable for rich VRPs. Metaheuristics
offer better solution quality and are more flexible but need careful tuning. Ex-
act methods guarantee optimality, but are limited to small instances due to high
computational cost. ALNS strikes a balance between flexibility and performance,
especially in large and complex problems, whereas machine learning-based models
are promising but still developing and face challenges with generalization and data
requirements.

Table 2.1: Overview of common optimization techniques for VRPs

Method Type Strengths Weaknesses

Classical Heuristics Easy to implement, fast Poor performance on complex
or large-scale problems

Metaheuristics
(e.g., ALNS)

Flexible, good solution quality,
handles rich constraints

Requires parameter tuning, no
optimality guarantee

Exact Methods Optimal solutions, mathemat-
ically rigorous

High computation time, lim-
ited scalability

ML-based Models Learns from data, adapts to
real-time input

Requires large datasets, lim-
ited practical deployment

2.5 Research Gaps and Summary
To summarize, many different approaches have been used to solve VRPs, including
classical heuristics, metaheuristics, exact methods, and machine learning. Among

7



2. Literature Review

these, ALNS stands out as an effective framework for solving complex, real-world
problems. It offers a good balance between performance and flexibility, and it works
well when combined with other methods.

There are still important challenges in this area. One of them is how to properly
handle split pickups and deliveries. Some earlier approaches address this by model-
ing split loads explicitly, such as the work by Nowak et al. [28], which introduces a
method for the Pickup and Delivery Problem with Split Loads. While this provides
more flexibility than simple node duplication, it also adds significant modeling and
computational complexity. There remains a need for heuristic methods that can
handle split service in a more integrated and scalable way, especially in cases where
requests may be shared across multiple vehicles or revisited within the same route.

Another major challenge is how to integrate electric vehicle charging into the routing
process. In many cases, charging is planned ahead of time, but in reality, conditions
can change, such as stations might be full, or vehicles might use more energy than
expected. Dealing with this in real time, especially under uncertainty, is still very
difficult.

In this work, we have focused on addressing both of these challenges. We investigate
a more flexible way to handle split pickups and deliveries, and we also explore how
charging decisions can be better integrated into routing. Both are important steps
toward building more realistic and adaptable vehicle routing systems.

8



3
Theoretical Foundations

This chapter presents the theoretical foundations and algorithms underlying the
solution approach used in this thesis. It begins by presenting the Vehicle Routing
Problem (VRP) from a mathematical perspective, and explaining the computational
complexity of the problem and its variants. This motivates the use of heuristic and
metaheuristic methods, as exact algorithms become impractical for large or realistic
instances. The core of the chapter, however, focuses on two powerful metaheuristics:
Large Neighborhood Search (LNS) and its extension, Adaptive Large Neighborhood
Search (ALNS). These techniques are particularly effective for solving complex and
large-scale VRP instances that involve multiple constraints, such as capacity limits
or heterogeneous fleets. LNS and ALNS form the backbone of the solution strategy
adopted in this study.

3.1 Vehicle Routing Problem
The Vehicle Routing Problem (VRP) is a well-known problem in transportation and
logistics [29], [2] , [4]. It deals with planning the best routes for a fleet of vehicles
that deliver goods from one central location (called a depot) to several customers.
Each customer needs a specific amount of goods, and each vehicle can only carry a
limited load. The goal is to reduce the total delivery cost while making sure that
every customer is visited only once, no vehicle carries more than it should, and all
the planned routes are possible to follow. A comprehensive treatment of this model
and its variants can be found in Toth and Vigo [30].

Let G = (V ,A) be a complete directed graph, where V = {0, 1, 2, . . . , n} is the
set of nodes, with node 0 denoting the depot and nodes 1 through n representing
customers. Let K = {1, 2, . . . , m} denote the set of vehicles. The set of arcs is
defined as A = {(i, j) ∈ V × V | i ̸= j}, where each arc (i, j) is associated with a
cost cij ≥ 0.

Each customer i ∈ V \ {0} has a demand di > 0, and each vehicle k ∈ K has a
capacity limit Q.

The decision variables are defined as follows:

9



3. Theoretical Foundations

uk
i ∈ R≥0 represents the cumulative load of vehicle k upon arrival at customer i.

xk
ij =

1, if vehicle k travels directly from node i to node j,

0, otherwise.

yk
i =

1, if customer i is visited by vehicle k,

0, otherwise.

The VRP can then be formulated as:

min
∑
k∈K

∑
(i,j)∈A

cijx
k
ij (3.1)

s.t.
∑
k∈K

yk
i = 1 ∀i ∈ V \ {0} (3.2)

∑
j∈V\{i}

xk
ij −

∑
j∈V\{i}

xk
ji =

1, if i = 0
0, if i ∈ V \ {0}

∀i ∈ V , k ∈ K (3.3)

yk
i =

∑
j∈V

xk
ij ∀i ∈ V \ {0}, k ∈ K (3.4)

uk
i − uk

j + Qxk
ij ≤ Q− dj ∀(i, j) ∈ A, k ∈ K (3.5)

di ≤ uk
i ≤ Q ∀i ∈ V \ {0}, k ∈ K (3.6)

xk
ij ∈ {0, 1} ∀(i, j) ∈ A, k ∈ K (3.7)

yk
i ∈ {0, 1} ∀i ∈ V , k ∈ K (3.8)

Objective (3.1) minimizes the total cost incurred by all vehicle routes. Constraint (3.2)
ensures that each customer is visited exactly once by a single vehicle. Flow conser-
vation for each vehicle is maintained through constraint (3.3), which ensures that
the number of routes entering and leaving each node is balanced. Constraint (3.4)
links the routing variables to the customer assignment decisions. Subtour elimina-
tion and capacity feasibility are enforced by the Miller–Tucker–Zemlin (MTZ) con-
straints in (3.5), while Constraint (3.6) bounds the cumulative vehicle load. Finally,
Constraints (3.7) and (3.8) enforce the binary nature of the routing and customer
assignment variables, respectively.

This VRP formulation supports explicit modeling of individual vehicle routes and is
widely used in extensions involving heterogeneous fleets, time windows, or electric
vehicles [30].

10



3. Theoretical Foundations

3.1.1 The Complexity of the Vehicle Routing Problem
The Vehicle Routing Problem is known to be NP-hard [31], meaning that the com-
putational effort required to solve it exactly increases exponentially with problem
size. As a result, solving large instances optimally is often impractical, particularly
in real-world contexts where timely decisions are crucial.

The computational complexity of the VRP has been the subject of extensive re-
search. A seminal study by Lenstra and Rinnooy Kan [31] analyzed the complexity
of various vehicle routing and scheduling problems, showing that many of them are
computationally intractable when approached with exact methods, especially as the
number of customers, constraints, or vehicle types grows.

Due to these challenges, most practical applications rely on heuristic or metaheuristic
approaches. Although these methods do not guarantee optimality, they can produce
high-quality solutions within acceptable computation times. As a result, they are
widely used in real-world applications such as logistics, freight distribution, and
transportation planning.

3.2 Large Neighborhood Search
Large Neighborhood Search (LNS) is a heuristic method designed to solve complex
combinatorial problems with large solution spaces. It was first applied to vehicle
routing by Shaw [32], and later extended and improved by Pisinger and Røpke [11],
[33], who introduced adaptive mechanisms to enhance its performance.

The key idea in LNS is to change a large part of the current solution, instead of
making small changes like traditional local search methods. This helps avoid getting
stuck in a low-quality solution and gives the algorithm a better chance to find an
improved one. The method works in two steps: it destroys part of the solution, and
then repairs it to create a new one.

3.2.1 Fundamental Concepts
Consider a combinatorial optimization problem defined over the set of feasible solu-
tions X . The goal is to find a solution x ∈ X that minimizes a given cost function
c(x), formally expressed as:

minimize c(x)
subject to x ∈ X ,

where c(x) denotes the objective value (or cost) associated with solution x.

In LNS, we use two important steps:

Destroy step: This removes some parts of the current solution, creating a
partial or even broken solution.

11



3. Theoretical Foundations

Repair step: This tries to fix the broken solution and make it valid again.

Let’s call the set of broken or incomplete solutions X partial. Then, we define:

Destruction operator d : X → X partial selectively removes a subset of elements
from a solution, thereby creating a partial structure.

Repair operator r : X partial → X reconstructs a feasible solution by reinserting
the removed elements.

x
destroy−−−−→ xpartial repair−−−→ x′,

The neighborhood explored by LNS around a solution x is not defined directly, but
instead through a process. It includes all the possible solutions that can be created
by first using a destruction operator to partially break the current solution, and then
applying a repair operator to rebuild it. This process can be described formally as:

N (x) = {r(d(x)) | d and r are applied to x, and x ∈ X}.

Instead of listing all possible neighboring solutions, LNS selects a few candidates
from this large neighborhood in each iteration. This makes it possible to explore
the solution space efficiently without needing to examine every option.

3.2.2 LNS Algorithmic Framework
The LNS works by gradually improving solutions through repeated destruction and
repair steps. During its execution, it keeps track of three main variables: x which
is the current solution being worked on; xb, the best solution found so far; and xt,
a temporary solution created in each iteration.

The structure of the overall LNS process, as depicted in Algorithm 1, is based on
the methodology proposed by Pisinger and Røpke [33].

Algorithm 1 Large Neighborhood Search
Require: An initial feasible solution x
Ensure: Best solution xb

1: xb ← x
2: repeat
3: Generate a candidate solution: xt ← r(d(x))
4: if accept(xt, x) then
5: Update current solution: x← xt

6: end if
7: if c(xt) < c(xb) then
8: Update best solution: xb ← xt

9: end if
10: until a stopping criterion is satisfied
11: return xb

12



3. Theoretical Foundations

At each iteration, the destroy and repair operators are applied to the current solution
x to produce a new candidate solution xt. This candidate may replace the current
solution depending on an acceptance criterion, and if it improves the best-known
solution, xb is updated. The algorithm terminates once a stopping condition, such
as a maximum number of iterations or a time limit, is met.

An illustration of this process is shown in Figure 3.1, where a portion of the solution
is removed and later reintegrated in a new configuration, demonstrating how LNS
explores different neighborhoods in the search space.

Before (Original Route)

Depot

C1
C2

C3

C4

C5

X

X

After (Repaired Route)

Depot

C1
C2

C3

C4

C5

Figure 3.1: Illustration of the LNS process. Left: The initial solution is shown
as a gray dashed route visiting all customers. Two customers, C2 and C3, are
selected for removal during the destruction phase. Right: The repaired solution is
constructed by reinserting the removed customers using a different route. The solid
blue arrows show the new improved sequence, demonstrating the flexibility of LNS
in exploring alternative neighborhoods.

3.2.3 Acceptance Criterion
To balance exploitation (exploiting high-quality regions of the solution space) and
diversification (exploring new areas), LNS often incorporates an acceptance criterion
inspired by the principles of Simulated Annealing. As proposed by Pisinger and
Røpke [33], the probability of accepting a candidate solution xt, given the current
solution x, is defined as:

P
(
A(xt, x)

)
=


1, if c(xt) ≤ c(x),

exp
(
−c(xt)− c(x)

T

)
, otherwise,

where A(xt, x) denotes the acceptance probability, c() is the cost function, and T > 0
is a temperature parameter controlling the acceptance probability.

This acceptance criterion ensures that an improved solution (c(xt) ≤ c(x)) is always

13



3. Theoretical Foundations

accepted. However, if the candidate solution is worse, it may still be accepted with
a probability that decreases as the cost difference increases and as the temperature
T decreases. This allows the algorithm to escape local optima by accepting worse
solutions.

At the beginning of the search, a high temperature encourages exploration by in-
creasing the likelihood of accepting worse solutions. As the search progresses, the
temperature is gradually reduced, shifting the focus toward exploitation. A common
cooling schedule is given by:

Tnew = α · Told,

where 0 < α < 1 is a cooling rate parameter that controls how quickly the temper-
ature decreases over time.

3.2.4 Design of Destroy and Repair Operators
The performance of the LNS algorithm depends on how its destruction and repair
operations are designed [33].

In the destroy phase, part of the current solution is removed to create a partial, and
often infeasible solution [33]. This removal can be random to encourage exploration.
Strategic destruction helps balance between exploring new solutions and improving
existing ones.

The repair phase then rebuilds a feasible solution [33]. This can be done using fast
greedy methods that reinsert elements with minimal negative impact, or through
more advanced strategies like random insertion or based on optimization techniques,
including mixed integer programming. These methods can lead to better-quality
solutions, although they may require more computational effort.

Choosing how much of the solution to destroy is important: removing too little limits
the search to nearby solutions, while removing too much can make the solution hard
to repair, which increases computation time and might lower overall quality.

3.3 Adaptive Large Neighborhood Search
Adaptive Large Neighborhood Search (ALNS) is an extension of the traditional
Large Neighborhood Search (LNS) method, first introduced by Røpke and Pisinger
[11], [33]. Unlike standard LNS, ALNS does not rely on a single destroy and repair
strategy. Instead, it chooses from several different heuristics during the search pro-
cess. This flexibility makes the algorithm more robust and better able to adjust to
different types of problems and stages in the search.

14



3. Theoretical Foundations

3.3.1 Motivation and Overview
Traditional LNS relies on a fixed pair of destroy and repair heuristics throughout
the search process. In contrast, Adaptive Large Neighborhood Search (ALNS) ex-
tends this approach by employing a pool of heuristics and adjusting their selection
probabilities based on performance feedback [33]. Instead of using the same pair of
heuristics, ALNS evaluates the effectiveness of each one and updates their chances
of being selected accordingly. This adaptive mechanism enhances the algorithm’s
ability to respond to varying problem characteristics.

From the perspective of neighborhood search, ALNS can be viewed as navigat-
ing among multiple large neighborhoods, each defined by a distinct destroy-repair
heuristic pair. The choice of the next neighborhood is influenced by a selection strat-
egy that favors heuristics with better performance. This adaptivity supports both
diversification by exploring a broader portion of the solution space and exploitation
by reinforcing the application of successful heuristics.

N ∗

N1

N2
N3

N4

N5

x

Figure 3.2: Illustration of neighborhoods used by ALNS. The current solution is
marked as x. ALNS operates on structurally diverse neighborhoods N1, . . . , Nk, each
defined by a specific destroy-and-repair heuristic. These neighborhoods are all subsets
of a larger neighborhood N ∗, defined by modifying up to q variables — the maximum
degree of destruction allowed.

3.3.2 ALNS Algorithmic Framework
Let:

W − = {d1, d2, . . . , dk} be the set of destroy heuristics,

W + = {r1, r2, . . . , rℓ} be the set of repair heuristics.

Each operator di ∈ W − and rj ∈ W + is associated with a weight ρ−
i and ρ+

j ,
respectively. These weights determine the selection probability of each heuristic,
which is updated throughout the search based on performance feedback.

15



3. Theoretical Foundations

The overall procedure is summarized in Algorithm 2, following the structure pro-
posed in [33].

Algorithm 2 Adaptive Large Neighborhood Search (ALNS)
Require: Initial feasible solution x ∈ X
Ensure: Best solution found xb

1: Initialize xb ← x, ρ−
i ← 1, ρ+

j ← 1 for all i, j
2: repeat
3: Select destroy di ∈ W − and repair rj ∈ W + using roulette wheel selection

based on ρ−, ρ+

4: Generate candidate: xt ← rj(di(x))
5: if accept(xt, x) then
6: x← xt

7: end if
8: if c(xt) < c(xb) then
9: xb ← xt

10: end if
11: Update weights ρ−

i , ρ+
j according to reward function

12: until stopping criterion is satisfied
13: return xb

Algorithm 2 begins with an initial feasible solution x, which is also recorded as the
best-known solution xb. For each destroy heuristic di ∈ W − and repair heuristic
rj ∈ W +, an associated weight ρ−

i or ρ+
j is initialized to 1. These weights are later

used to guide the heuristic selection process.

During each iteration of the algorithm, one destroy and one repair heuristic are
selected. The selected destroy heuristic di is applied to the current solution x,
producing a partial solution. This partial solution is then completed using the
selected repair heuristic rj, resulting in a new candidate solution xt.

The candidate solution xt is accepted based on a predefined acceptance criterion
that allows the algorithm to accept worse solutions in order to escape local optima
and promote diversification. If xt is accepted, it replaces the current solution x.

If the candidate solution is better than the best solution found so far (c(xt) < c(xb)),
it is stored as the new best solution xb. Regardless of acceptance, the weights ρ−

i and
ρ+

j of the selected heuristics are updated according to a reward mechanism, which
increases the weight of heuristics that lead to better outcomes.

The algorithm continues until a stopping criterion is met, such as a time limit or
a maximum number of iterations. At the end of the search, the best solution xb

encountered during the process is returned.

16



3. Theoretical Foundations

3.3.3 Destroy and Repair Operator Selection
ALNS selects destroy and repair heuristics at each iteration using an adaptive selec-
tion mechanism, as proposed by Pisinger and Ropke [33]. Each heuristic di ∈ W −

(destroy) and rj ∈ W + (repair) is assigned a weight, denoted by ρ−
i and ρ+

j , respec-
tively. The probability of selecting a destroy heuristic di is computed by normalizing
its weight against the sum of all destroy heuristic weights:

ϕ−
i = ρ−

i∑|W −|
k=1 ρ−

k

.

Similarly, repair heuristics are selected with probability:

ϕ+
j =

ρ+
j∑|W +|

l=1 ρ+
l

.

These probabilities form the basis for a roulette wheel selection process, where
heuristics with higher weights are more likely to be chosen, but lower-weight heuris-
tics still have a chance, allowing for a balance between exploitation and diversifica-
tion.

After applying the selected pair of heuristics (da, rb), the resulting candidate solution
xt is evaluated, and a reward score y is assigned to this heuristic pair based on the
quality of xt:

y =


w1, if xt is a new global best,
w2, if xt improves the current solution,

w3, if xt is accepted,

w4, if xt is rejected,

where the weights satisfy w1 > w2 > w3 > w4 ≥ 0, reflecting the relative importance
of each outcome.

The heuristic weights are then updated using an exponential smoothing update rule
controlled by a parameter λ ∈ [0, 1]:

ρ−
a ← λρ−

a + (1− λ)y, ρ+
b ← λρ+

b + (1− λ)y.

This update rewards heuristics that contribute to the search progress, while gradu-
ally discounting past performance, allowing the algorithm to adaptively learn which
heuristics are most effective over time.

In summary, this adaptive mechanism, combining weighted roulette wheel selection
based on performance scoring and updating, enables ALNS to adjust its search
strategy to better explore and exploit the solution space.

17



3. Theoretical Foundations

3.3.4 Acceptance Criterion
Similar to the LNS approach, ALNS employs an acceptance criterion inspired by
Simulated Annealing to decide whether to accept a candidate solution xt. The
acceptance probability is defined as:

P
(
A(xt, x)

)
=


1, c(xt) ≤ c(x),

exp
(
−c(xt)− c(x)

T

)
, otherwise,

where A(xt, x) denotes the acceptance probability, c(·) is the cost function, and
T > 0 is a temperature parameter controlling the acceptance probability.

3.3.5 Design of Destroy and Repair Operators
The way destroy and repair operators are designed plays a major role in how well
the ALNS algorithm performs. According to [33], different types of operators are
used to balance exploration and improvement during the search.

Destroy operators that focus on diversification such as random removals methods,
help the algorithm explore new parts of the solution space by breaking up the current
solution structure. On the other hand, destroy operators focused on exploitation
target specific elements, such as the most costly or critical components, to create
opportunities for improved local solutions.

One popular destroy operator strategy removes elements that are similar in terms of
node, timing, or type. This approach enables more focused and flexible exploration
of local adjustments.

Repair operators can vary in how simple or complex they are. Greedy repair heuris-
tics work quickly by adding elements back in places that cause the least extra cost.
Regret-based methods, on the other hand, focus on reinserting the more difficult el-
ements first, to avoid making worse decisions later. In some cases, it is even possible
to use exact optimization techniques to solve small subproblems during the repair
phase. However, since these methods can be time consuming, the algorithm may
penalize slow repair strategies to maintain overall efficiency.

3.3.6 Properties and Strengths of ALNS
ALNS is a powerful and flexible framework for solving large and complex optimiza-
tion problems. One of its key strengths is its ability to adapt during the search
process. By learning from past performance, ALNS adjusts how it chooses heuris-
tics, allowing it to balance between exploitation and diversification. This makes the
method highly effective across different types of problem instances.

Another advantage of ALNS is that it requires relatively little parameter tuning.
Most of the parameters are related to the adaptive selection process, not the specific

18



3. Theoretical Foundations

heuristics themselves. In addition, the framework is very flexible. new destroy and
repair operators can easily be added, and the algorithm will automatically learn
when to use them based on their performance.

Unlike methods such as Variable Neighborhood Search (VNS), which follow a fixed
sequence of neighborhood structures, ALNS uses a probability-based approach that
responds to how well different heuristics perform. This makes ALNS especially
useful in problems where neighborhoods are defined by heuristics rather than strict
mathematical rules.

Figure 3.3 illustrates the key difference between the traditional LNS and ALNS.
In LNS, a fixed pair of destroy and repair heuristics is used throughout the search
process. This approach can be effective but may limit flexibility across different
problem landscapes.

ALNS extends this structure by maintaining a pool of heuristics and selecting a pair
based on their historical performance. After each iteration, the algorithm updates
the weights associated with each heuristic according to how well they performed,
allowing it to gradually favor more effective strategies. This adaptive layer helps
ALNS to better balance exploration and exploitation, improving its performance on
a wide range of complex optimization problems [33].

19



3. Theoretical Foundations

Current Solution x

Destroy d(x)

Repair r(d(x))

Candidate xt

Accept or Reject

LNS (Fixed Heuristics) ALNS (Adaptive Heuristics)

Current Solution x

Select di, rj (via weights)

Apply di(x) and rj(x)

Candidate xt

Accept or Reject

Update Weights
W

ei
gh

t
up

da
te

Figure 3.3: Comparison of LNS and ALNS. LNS uses a fixed destroy/repair strat-
egy, while ALNS adaptively selects from multiple heuristics and updates their weights
based on performance.

3.4 Dijkstra’s Algorithm
Dijkstra’s algorithm was invented in 1956 and is used for finding the shortest paths
from one node, to all other nodes in a graph [34]. The algorithm starts at a given
node and repeatedly visits the nearest unvisited node while keeping track of the
shortest distance and the predecessor to each node.

In the algorithm, initial shortest distances are first set for all nodes. The value
associated with the starting node is set to zero, while the values of the remaining
nodes are initialized to infinity [34]. Then the following is done: the unvisited node
with the lowest shortest distance is selected, and for each of its neighboring nodes,
the distance to the neighboring node via the selected node is calculated by adding

20



3. Theoretical Foundations

the shortest distance to the selected node, and the distance between the selected
node and the neighboring node. If this sum is lower than the saved shortest distance
to the neighboring node, it becomes the new shortest distance to that node and the
selected node is saved as the predecessor of the neighboring node. Finally, the
currently selected nodes is marked as visited, and the procedure repeats, starting
with a new unvisited node being selected.

Once all nodes have been visited, the shortest paths originating from the starting
node can be extracted by starting at a certain node and going to the preceding node
that was saved the shortest distance was found, then going to the preceding node
of that preceding node, and so on until the starting node is reached [34]. At that
point, a sequence of nodes has been obtained that makes up the shortest path.

21



3. Theoretical Foundations

22



4
Mathematical Model

This chapter presents the mathematical model developed to address the extended
Fleet Size and Mix Vehicle Routing Problem (FSMVRP). The model is based on
a flow formulation and integrates several real-world considerations, including het-
erogeneous vehicle types (electric and diesel), detailed loading and unloading op-
erations, energy consumption, and the possibility for vehicles to revisit locations.
The transportation network is represented as a directed graph, and the formula-
tion incorporates practical constraints that reflect both logistical requirements and
physical limitations commonly encountered in real-world operations.

4.1 Terminology

Before presenting the formal mathematical formulation, we introduce the key terms
and concepts used throughout the model. These definitions provide the necessary
context for understanding the structure of the transportation network, vehicle op-
erations, and logistical constraints included in the problem at hand.

Node A location in the network, representing either a customer node
(for pickup or delivery) or a charging station node (for battery
charging).

Arc A directed connection between two nodes in the transporta-
tion network, representing a possible movement or traversal
by a vehicle.

Fleet The complete set of vehicles available for assignment, which
may differ in volume capacity, mass capacity, or propulsion
system (e.g., electric or diesel).

Home Depot A central node in the transportation network representing the
main facility. Each vehicle must complete its route at its cor-
responding home depot, although it may start from a different
node.

23



4. Mathematical Model

Step A discrete index n representing a movement or transition in
the route from one node to the next node, where n = 1 marks
the first move.

Route A planned sequence of steps (nodes and arcs) followed by a
vehicle from origin to destination.

Service Time The time spent at a node performing loading, unloading, or
charging operations.

Loading/Unloading Methods
Cargo handling techniques available at nodes, with multiple
types possibly supported per vehicle and location.

A more explanation of the loading and unloading methods considered in this study is
provided below. These methods are based on practical handling processes described
by Toheed Ghandriz in [35] and are incorporated with specific service time rates.

The model considers four types of loading and unloading operations:

(i) Additional Semitrailer (AST): Attaching or detaching an additional semitrailer
to or from the vehicle.

(ii) On-Board Lift (OBL): Utilizing a lift mechanism mounted on the vehicle to
handle loading or unloading.

(iii) On-Board Waiting (OBW): Waiting at the node while external lift vehicles
perform cargo handling.

(iv) Straddle Carrier (SC): Using a straddle carrier at the node to transfer con-
tainers to or from the vehicle.

4.2 Problem Description
This section presents the description of the problem under consideration. The un-
derlying structure is based on a directed graph G = (V ,A), where V is the set of
nodes and A is the set of directed arcs representing feasible movements between
nodes.

The node set V is partitioned as:

V = V0 ∪ Vc ∪ Vs,

where V0 includes the home depot, Vc represents customer locations, and Vs contains
charging stations. Notably, a node may simultaneously belong to both Vc and Vs if
it serves both functions.

24



4. Mathematical Model

The arc set A consists of all valid transitions between nodes:

A = {(i, j) | i, j ∈ V , i ̸= j, i ̸= 0+, j ̸= 0−}.

Each arc (i, j) ∈ A represents a direct route from node i to node j and, the home
depot is modeled as two separate nodes: 0− (the starting home depot) and 0+ (the
ending home depot).

4.3 Sets, Parameters, and Decision Variables
This section defines the fundamental components used in the mathematical formu-
lation of the problem. The model relies on a set of indices and elements to represent
the network.

4.3.1 Sets
The model is built using several sets that help define the structure of the problem.
These sets are used as the basic indexing system throughout the mathematical
model.

Table 4.1: Sets used in the optimization model

Set Description

Kel Set of electric vehicles

Kd Set of diesel vehicles

K Set of all vehicles: K = Kel ∪ Kd

Vc Customer nodes

Vs Charging station nodes

Vo The home depot

V All nodes in the network: V = Vs ∪ Vc ∪ Vo

VL Loading nodes

VU Unloading nodes

N Set of steps: N = {1, ..., N1}

A Set of arcs between nodes

H Loading and unloading methods

1N denotes the maximum number of steps across all vehicles, i.e., N = maxk∈K Nk, where Nk

is the last step for vehicle k.

25



4. Mathematical Model

4.3.2 Decision Variables
The decision variables capture the routing behavior, product flows, energy con-
sumption, and service activities of vehicles across the network. These are the key
quantities optimized in the model.

Table 4.2: Decision variables used in the optimization model

Variable Description

yk,n
ij

Load carried by vehicle k ∈ K on arc (i, j) ∈ A at step n ∈ N
[kg].

mk,n,h
Load,i

Quantity loaded on vehicle k using method h ∈ H at node
i ∈ VL during step n [kg].

mk,n,h
Unload,i

Quantity unloaded from vehicle k using method h at node i ∈
VU during step n [kg].

zk,n
ij Energy consumed by vehicle k on arc (i, j) at step n [kWh].

Ek,n
ij

Energy needed for vehicle k to traverse arc (i, j) at step n
[kWh].

qk,n
i Battery level of vehicle k upon arrival at node i at step n [kWh].

pk,n
i

Amount of electricity charged by vehicle k at node i during
step n [kWh].

sk,n
i Service time of vehicle k at node i during step n [h].

τ k,n
i Arrival time at node i for vehicle k at step n [h].

βk
Total working time (service time and travel time) of vehicle k
[h].

αk Effective wage-incurring time of vehicle k [h].

xk,n
ij

Binary variable indicating whether vehicle k travels on arc (i, j)
at step n (1 if traveled, 0 otherwise).

δfix
k

Binary variable indicating whether vehicle k is used, i.e., per-
forms at least one loading operation (1 if used, 0 otherwise).

4.3.3 Parameters
This subsection lists the model’s fixed input data. These parameters define the
operational environment and constraints under which optimization occurs.

26



4. Mathematical Model

Table 4.3: Model parameters used in the optimization formulation

Parameter Description

cd Cost of diesel per liter (€/liter).

cel
i Electricity cost for charging at node i ∈ V [€/kWh].

cw Driver wage cost per hour [€/h].

cf
k Deployment cost associated with using vehicle k ∈ K [€].

ρ Density of cargo [kg/m3].

CV
k Maximum cargo volume capacity for vehicle k [m3].

CM
k Maximum mass capacity for vehicle k [kg].

DUnload,i Delivery demand at unloading node i [kg].

DLoad,i Pick-up demand of at loading node i [kg].

Bk Battery capacity for electric vehicle k [kWh].

Γk
ij

Coefficient for calculating energy consumption per unit load on
arc (i, j) for vehicle k [kWh/kg].2

Πk
ij Fixed energy consumption on arc (i, j) for vehicle k [kWh].2

tij Travel time between node i and node j [h].

tLU,h Time per unit mass for loading/unloading using method h [h/kg].

tmax
k Maximum available working time for vehicle k [h].

rel
g Charging rate per electricity unit for charging type g [h/kWh].

M
A sufficiently large constant related to energy flow constraints
[kWh].

Q A sufficiently large constant used for constraint activation [kg].

Ω Conversion constant for balancing diesel and electric energy units
(liter/kWh).

Lk,h
i

Binary parameter indicating whether vehicle k is allowed to load
or unload using method h at node i (1 if allowed, 0 otherwise).

2Both Γk
ij and Πk

ij are derived from detailed vehicle dynamics and propulsion parameters not
included in this report.

27



4. Mathematical Model

4.4 Objective Function of the Model
Using the sets, parameters, and decision variables introduced earlier, this section
now presents the full mathematical model for the VRP, including all the necessary
constraints.

minimize
∑
k∈K

cw · αk︸ ︷︷ ︸
Wage cost

+
∑
k∈K

δfix
k · c

f
k︸ ︷︷ ︸

Vehicle deployment cost

+
∑

k∈Kd

∑
n∈N

∑
(i,j)∈A

Ω · cd · zk,n
ij︸ ︷︷ ︸

Diesel energy cost

+
∑

k∈Kel

∑
n∈N

∑
i∈Vs

cel
i · p

k,n
i︸ ︷︷ ︸

Electric charging cost

(4.1)

The objective function (4.1) minimizes the total cost of routing a heterogeneous
fleet consisting of both electric vehicles (EVs) and diesel vehicles. It integrates
cost elements related to driver wages, fixed vehicle usage, fuel consumption, and
electricity charging:

Wage cost: This term reflects the cost of active driver working time. It
includes both service time and travel time, but is only counted if the vehicle
performs at least one loading operation (i.e., when δfix

k = 1). The variable αk

is activated through linear constraints mentioned in (4.5.4) that link it to the
actual working time βk and the binary service indicator δfix

k .

Vehicle deployment cost: Each vehicle that serves at least one customer
(by performing a loading operation) is assigned a deployment cost. This cost
is included in the model to encourage using fewer vehicles when possible and
to avoid unnecessary vehicle usage.

Diesel fuel cost: This term models the energy cost of diesel vehicles based on
energy consumption zk,n

ij over arc (i, j) at step n, scaled by a unit conversion
factor Ω to match fuel and energy units. It is only applied to vehicles in the
diesel subset Kd.

Electricity charging cost: Electric vehicles incur a cost based on the amount
of energy pk,n

i charged at node i, with cost rates cel
i .

The structure ensures that costs are only applied when the vehicle is actually doing

28



4. Mathematical Model

something important; like loading products or using energy.

4.5 Model Constraints
This section introduces the mathematical constraints that guide how the vehicle
routing problem works in this model. These rules make sure everything is feasible
and realistic in terms of vehicle movements, product deliveries, energy use, and time
limits. To make the explanation clearer, the constraints are organized into groups
based on the main parts of the problem: routing, product flow, battery energy, time
management, and how variables are defined.

4.5.1 Routing Constraints
The routing constraints control how vehicles move through the network, making
sure the routes are logical and continuous. These rules ensure that each vehicle
starts from the home depot, follows a step-by-step path through the network, and
eventually arrives at its final destination.

∑
j∈V:

(j,i)∈A

xk,n
ji −

∑
j∈V:

(i,j)∈A

xk,n+1
ij = 0, i ∈ V \ Vo, k ∈ K, n ∈ N \ {Nk}

(4.2)∑
j∈V:

(0−,j)∈A

xk,1
0−j = 1, k ∈ K (4.3)

∑
i∈V:

(i,0+)∈A

∑
n∈N

xk,n
i0+ = 1, k ∈ K. (4.4)

Constraint (4.2) makes sure that the vehicle routes are continuous. For each vehicle
k, and each step n, the number of incoming arcs to node i must equal the number
of outgoing arcs, except at the home depot.

Constraint (4.3) ensures that every vehicle begins its route at the home depot by
enforcing exactly one outgoing arc from that at the first step.

Constraint (4.4) guarantees that each vehicle terminates its route at the end depot
0+, by requiring exactly one incoming arc into this node across all steps.

4.5.2 Loading and Unloading Constraints
This group of constraints governs the flow of multiple product types, ensuring that
vehicle capacity limits are respected, and that loading and unloading activities meet
customer demands while preserving consistency throughout the network.

29



4. Mathematical Model

yk,n
ij ≤ CM

k xk,n
ij , (i, j) ∈ A, k ∈ K, n ∈ N , (4.5)

yk,n
ij

ρ
≤ CV

k xk,n
ij , (i, j) ∈ A, k ∈ K, n ∈ N , (4.6)

∑
i∈VL

∑
n∈N

∑
h∈H

mk,n,h
Load,i ≤ Q · δfix

k , k ∈ K, (4.7)

∑
n∈N

∑
h∈H

mk,n,h
Load,i ≤ CM

k , i ∈ VL, k ∈ K, (4.8)

∑
n∈N

∑
h∈H

mk,n,h
Load,i

ρ
≤ CV

k , i ∈ VL, k ∈ K, (4.9)

∑
k∈K

∑
n∈N

∑
h∈H

mk,n,h
Load,i = DLoad,i, i ∈ VL, (4.10)

∑
k∈K

∑
n∈N

∑
h∈H

mk,n,h
Unload,i = DUnload,i, i ∈ VU , (4.11)

∑
j∈V

(i,j)∈A

yk,n+1
ij −

∑
j∈V:

(j,i)∈A

yk,n
ji =

∑
h∈H

mk,n,h
Load,i, i ∈ VL, k ∈ K, n ∈ N , (4.12)

∑
j∈V:

(j,i)∈A

yk,n
ji −

∑
j∈V:

(i,j)∈A

yk,n+1
ij =

∑
h∈H

mk,n,h
Unload,i, i ∈ VU , k ∈ K, n ∈ N , (4.13)

mk,n,h
Load,i ≤ DLoad,i . Lk,h

i , i ∈ VL, h ∈ H, k ∈ K, n ∈ N , (4.14)

mk,n,h
Unload,i ≤ DUnload,i . Lk,h

i , i ∈ VU , h ∈ H, k ∈ K, n ∈ N . (4.15)

Constraints (4.5) and (4.6) ensure that the load carried by each vehicle respects
its mass and volume capacities, respectively. They enforce that the total quantity
transported does not exceed the vehicle’s physical constraints on any arc.

Constraint (4.7) ensures that the binary indicator δfix
k is set to one if vehicle k per-

forms at least one loading activity. In cases where no loading occurs, the summation
on the left hand side evaluates to zero, and the minimization objective naturally
favors δfix

k = 0, thereby avoiding unnecessary fixed vehicle and wage costs. The con-
stant Q is set as a sufficiently large upper bound representing the maximum load
a vehicle could reasonably carry, ensuring that the constraint activates only when
actual service is performed.

Constraints (4.8) and (4.9) ensure that at every loading node, the total quantity of

30



4. Mathematical Model

loaded cargo must not exceed the vehicle’s physical carrying capabilities in terms of
weight and space.

Constraints (4.10) and (4.11) guarantee that each product is properly picked up
from home depot and fully delivered to its destination. These rules ensure that all
product demands are completely met throughout the network.

Constraints (4.12) and (4.13) ensure flow balance at each node in the network.
Specifically, constraint (4.12) applies to loading nodes (VL) and ensures that the
difference between outgoing and incoming flows equals the quantity of goods loaded
at the node. Conversely, constraint (4.13) applies to unloading nodes (VU) and
ensures that the difference between incoming and outgoing flows matches the quan-
tity of goods unloaded. Together, these constraints preserve the continuity of the
transported load across the vehicle routes.

Finally, constraints (4.14) and (4.15) make sure that loading and unloading op-
erations are allowed. They check that products are only handled where they are
actually available and that each vehicle is permitted to load or unload product at
certain nodes, based on Lk,h

i .

4.5.3 Battery Energy Constraints
These constraints control how battery energy is used and managed for electric vehi-
cles along their routes. They make sure the energy levels stay consistent with how
the vehicles move, respect the battery capacity limits, and correctly account for any
recharging that happens during the trip.

qk,1
i = Bk, k ∈ Kel, i ∈ Vo \ 0+, (4.16)

Ek,n
ij = Γk

ij · y
k,n
ij + Πk

ij, (i, j) ∈ A, k ∈ K, n ∈ N ,

(4.17)
Ek,n

ij −M(1− xk,n
ij ) ≤ zk,n

ij , (i, j) ∈ A, k ∈ K, n ∈ N ,
(4.18)

zk,n+1
ij −Bk(1− xk,n+1

ij ) ≤ qk,n
i + pk,n

i − qk,n+1
j , (i, j) ∈ A, k ∈ Kel, n ∈ N \ {Nk},

(4.19)
qk,n

i + pk,n
i ≤ Bk, i ∈ Vs \ 0+, k ∈ Kel, n ∈ N ,

(4.20)
qk,Nk

0+ + pk,Nk
0+ = Bk, k ∈ Kel. (4.21)

Constraint (4.16) initializes the battery level of each electric vehicle k to its full
capacity Bk at the beginning of the route.

Constraint (4.17) computes the total energy requirement Ek,n
ij to traverse arc (i, j) at

31



4. Mathematical Model

step n. This energy includes two parts: one that depends on how much the vehicle
is carrying, and another fixed part, Πk

ij,

Constraint (4.18) ensures consistency between routing decisions and energy con-
sumption. When arc (i, j) is used (xk,n

ij = 1), the consumed energy zk,n
ij must equal

or exceed the required energy Ek,n
ij . If the arc is not traversed (xk,n

ij = 0), the con-
straint is relaxed via a big-M to avoid enforcing energy requirements on inactive
arcs.

Constraint (4.19) ensures that the battery level at the next visited node j accounts
for the previous battery level, any energy recharged at node i, and the energy con-
sumed along arc (i, j). The formulation applies only to electric vehicles and excludes
the final step Nk, since there’s no next step after that point.

Constraint (4.20) ensures that the total energy in the battery, does not exceed the
vehicle’s battery capacity Bk. This constraint applies to all charging nodes except
the end depot 0+, which is treated separately.

Earlier formulation where recharging at 0+ was disallowed, constraint (4.21) allows
the vehicle to recharge at the depot before completing its operations.

4.5.4 Time Constraints
Charging Rate and Charging Time Definition

Let rel denote the charging rate (in hours) required to supply 1 kWh, defined as:

rel = 1
P

,

where P is the charging power in kW.

The charging time at a node depends on the quantity of energy charged. Given this,
the charging time incurred by vehicle k at node i ∈ Vs during step n ∈ N is:

sk,n
ch,i = rel · pk,n

i ,

where pk,n
i is the amount of energy (in kWh) charged.

Service Time Definition

At each visited node i ∈ V , the service time, sk,n
i , must account for both charging and

product handling activities. Since, charging, loading, and unloading can happen at
the same time, the service time is defined as the maximum of these three durations:

sk,n
i = max

{
sk,n

ch,i, sk,n
L,i , sk,n

U,i

}
, ∀i ∈ V , k ∈ K, n ∈ N ,

32



4. Mathematical Model

with loading and unloading times given by:

sk,n
L,i =

∑
h∈H

tLU,h ·mk,n,h
Load,i,

sk,n
U,i =

∑
h∈H

tLU,h ·mk,n,h
Unload,i.

Here, tLU,h is the time per unit (e.g., per kg) for loading/unloading using method
h ∈ H.

Therefore, time constraints are given by:

∑
i∈V

∑
n∈N

sk,n
i +

∑
(i,j)∈A

∑
n∈N

tij · xk,n
ij ≤ tmax

k , k ∈ K, (4.22)

∑
i∈V\{0+}

∑
n∈N \{Nk}

sk,n
i +

∑
(i,j)∈A

∑
n∈N

tij · xk,n
ij ≤ βk, k ∈ K, (4.23)

βk −W (1− δfix
k ) ≤ αk, k ∈ K, (4.24)

sk,n
i + tij − tmax(1− xk,n+1

ij ) ≤ τ k,n+1
j − τ k,n

i , (i, j) ∈ A, k ∈ K, n ∈ N \ {Nk},
(4.25)

τ k,Nk
0+ + sk,Nk

0+ ≤ tmax, k ∈ K, (4.26)

αk ≤ tmax
k · δfix

k , k ∈ K, (4.27)

sk,n
ch,i ≤ sk,n

i , i ∈ Vs, k ∈ Kel, n ∈ N (4.28)

sk,n
L,i ≤ sk,n

i , i ∈ VL, k ∈ K, n ∈ N (4.29)

sk,n
U,i ≤ sk,n

i , i ∈ VU , k ∈ K, n ∈ N (4.30)

Constraints (4.22) represent the total working time of each vehicle k, composed of
both the service time sk,n

i at visited nodes and the travel time tij over traversed arcs,
and ensures that this total does not exceed the vehicle’s allowable working duration
tmax
k , thus maintaining route feasibility with respect to legal or contractual working

limits.

In contrast, constraint (4.23) introduces an auxiliary variable βk that upper-bounds
the same working time expression, again excluding activities at the end node 0+
in order to decouple time feasibility from cost modeling. This is necessary because
the final service at 0+ typically represents a terminal recharge operation, which
does not contribute to wage cost under the operational assumptions of this model.

33



4. Mathematical Model

Therefore, any service time at the end depot is deliberately excluded from wage cost
calculations while still being allowed in the time feasibility constraint.

Constraint (4.24) links the working time associated with βk to another variable αk,
that should incur wage costs. To control whether this cost is applied, the binary
variable δfix

k is used. If the vehicle is active (δfix
k = 1), the constraint forces αk to be

at least as large as βk, meaning all of the working time is counted for wage purposes.
If the vehicle is not used (δfix

k = 0), the constraint becomes:

βk −W ≤ αk,

where W is a large constant that serves as a “Big-M” parameter in mixed-integer
programming. Mathematically, this constraint has the following logic:

• If δfix
k = 1, then the constraint becomes:

βk ≤ αk,

ensuring that the wage-incurring time αk is at least equal to the actual working
time βk.

• If δfix
k = 0, then the constraint becomes:

βk −W ≤ αk.

Since W is chosen such that W ≥ maxk tmax
k , the left hand side becomes a

large negative number (i.e., βk −W ≪ 0). Therefore, the inequality is always
satisfied for any value of αk.

Thus, W enables the conditional enforcement of wage cost: αk must capture working
time only if the vehicle is used. To ensure correctness, W should be chosen as:

W ≥ max
k∈K

tmax
k ,

so that the logic holds in all cases.

Constraint (4.25) guarantees that the arrival time at the subsequent node reflects
both the service time at the current node and the travel time required to reach
the next node. Additionally, the term involving tmax(1 − xk,n+1

ij ) acts as a big-M
style that effectively deactivates the constraint when arc (i, j) is not traversed. This
formulation is also instrumental in preventing subtours, as it links the sequence of
time steps to the selected arcs, ensuring a coherent progression of time across the
route.

The constraint (4.26) ensures that the total duration of the vehicle’s route, from
initial departure to the return to the depot, does not exceed the maximum allowed
operational time.

Constraint (4.27) define another auxiliary variable αk. The logic is as follows:

34



4. Mathematical Model

• If the vehicle performs at least one loading operation (δfix
k = 1), then αk = βk,

meaning all working time contributes to wage cost.

• If the vehicle does not perform any customer service (δfix
k = 0), then αk = 0,

and no wage cost is accrued in the objective function.

This approach makes sure that wage costs are only applied to vehicles that actu-
ally perform deliveries or pickups. Vehicles that are only used for repositioning or
charging do not add to the wage cost.

Finally, constraints (4.28)–(4.30) ensure that the charging time sk,n
ch,i, loading time

sk,n
L,i , and unloading time sk,n

U,i at any visited node are each bounded by the service
time sk,n

i .

4.5.5 Non-Negativity and Binary Constraints
To keep the optimization model consistent and feasible, certain rules are applied to
the decision variables. All continuous variables must have non-negative values, and
all variables that represent logical choices must be binary.

yk,n
ij ≥ 0, ∀(i, j) ∈ A, k ∈ K, n ∈ N ,

Ek,n
ij ≥ 0, ∀(i, j) ∈ A, k ∈ K, n ∈ N ,

zk,n
ij ≥ 0, ∀(i, j) ∈ A, k ∈ K, n ∈ N ,

mk,n,h
Load,i ≥ 0, ∀i ∈ VL, k ∈ K, n ∈ N , h ∈ H,

mk,n,h
Unload,i ≥ 0, ∀i ∈ VU , k ∈ K, n ∈ N , h ∈ H,

qk,n
i ≥ 0, ∀i ∈ V , k ∈ K, n ∈ N ,

pk,n
i ≥ 0, ∀i ∈ V , k ∈ K, n ∈ N , g ∈ G,

sk,n
i ≥ 0, ∀i ∈ V , k ∈ K, n ∈ N ,

τ k,n
i ≥ 0, ∀i ∈ V , k ∈ K, n ∈ N ,

αk ≥ 0, ∀k ∈ K,

βk ≥ 0, ∀k ∈ K,

xk,n
ij ∈ {0, 1}, ∀(i, j) ∈ A, k ∈ K, n ∈ N ,

δfix
k ∈ {0, 1}, ∀k ∈ K.

35



4. Mathematical Model

36



5
Implementation

To solve the problem at hand within a reasonable amount of time, the ALNS frame-
work was adapted and implemented. The implemented algorithm uses the same
objective function and constraints as presented in section 4, but is not related to
the concept of linear programming. The basic idea of the ALNS framework is that
many iterations are performed and in each iteration the currently accepted solution
is destroyed and repaired in an attempt to find a better solution. For this problem,
destruction is done by removing visits where certain requests are served from the
current vehicle routes, and repairing is done by inserting visits so that the requests
are once more served. There are multiple ways in which the visits associated with
the requests can be removed and inserted, and in each iteration a removal method
and an insertion method is chosen randomly. As the search progresses, the methods
that contribute to new best solutions being found are rewarded and the probabilities
of these methods being chosen in future iterations increase - this is what adaptive
refers to in ALNS.

5.1 Rationale for Choosing the ALNS Framework

Among the many algorithms available for solving vehicle routing problems, ALNS
was selected for this thesis because it offers a strong balance between performance,
flexibility, and practicality. One of the key reasons for choosing ALNS is its ability to
adapt during the search process, which is especially valuable for complex problems
that include many different types of constraints.

The problem definition in this study is not a simple one, as it includes electric and
diesel vehicles, split pickups and deliveries, mass and volume constraints, as well as
charging decisions. Most traditional heuristics or metaheuristics struggle to handle
this level of detail. In contrast, ALNS allows the user to define and combine custom
destroy and repair operators, making it possible to guide the search in ways that
respect the specific rules of the problem. As the algorithm runs, it keeps track of
which operators are working best and adjusts their usage. This adaptive feature
helps the algorithm focus on the most promising parts of the solution space without
getting stuck too early.

37



5. Implementation

Further, previous research supports the use of ALNS in energy-aware VRPs, as
shown in the work by Røpke and Pisinger [11], Hiermann et al. [18], and Keskin
and Çatay [19].

5.2 Solution Representation

Due to aspects of the problem such as split pickups and deliveries, and charging, a
solution cannot simply consist of routes that are represented as sequences of nodes
like in many other VRP variations. In this thesis, a solution is a collection of
vehicle routes such that each vehicle has a route, and each route is a sequence of
visits. A visit contains information about which node was visited, the amount of
mass loaded/unloaded, and the amount of time that the vehicle stayed at the node.
Additionally, it holds information about how much energy was charged during the
visit. The visit objects contain relative values, e.g. mass change instead of total
mass, which was a deliberate choice as it removes the need for updating the data in
the objects after removals and insertions of visits are done. It is worth noting that
some vehicles may not be used, but the solution will still contain routes for these
vehicles. The reason for this is that a vehicle has to at least travel from the start
depot node to the end depot node, which makes up a route. These two nodes have
been created artificially and share the same geographical location, the location of
the home depot in the instance.

In addition to the solution, i.e. the vehicle routes, requests objects are used to keep
track of which vehicles have served which requests. A request contains information
about the pickup node and the delivery node for the request, the total amount of
cargo mass that needs to be transported to serve the request, how much cargo mass
is left to serve, and which vehicles have served the request. From here on out, when
referring to the removal or insertion of requests, what is meant is that visits in which
the requests are served are removed or inserted from the routes in the solution.

To give a clearer idea of what the solution representation looks like, a simple example
is presented below, in Table 5.1.

38



5. Implementation

Vehicle
ID

Node Load
Change
(kg)

Loading
Time (s)

Energy
Charged
(kWh)

Charging
Power
(kW)

1 1 0 0 0 -
1 2 +25 000 9000 100 250
1 3 -25 000 9000 0 -
1 4 0 0 200 250
2 1 0 0 0 -
2 2 +5000 1800 100 250
2 3 -5000 1800 0 -
2 4 0 0 200 250

Request
ID

Pickup
Node

Delivery
Node

Initial
Mass
(kg)

Remain-
ing Mass
(kg)

Served
By

1 2 3 30 000 0 1, 2

Table 5.1: An example of a solution in the form of a table where every row is a
visit, together with a table with the request that is served in the solution. In this
case, there is a split pickup and delivery as one request is served by two vehicles.

5.3 Removal Methods
In this subsection, the three removal methods used in this thesis are introduced.
The random removal and worst removal methods have previously been presented in
[11], and single route removal is based on the route elimination method from [20].

Each removal method selects which requests to remove using a distinct strategy,
and once the requests have been selected, the visits in which they are served are
removed from the solution. Moreover, visits done strictly to charge are removed, and
charging at customer nodes is removed. This is done to avoid the situation where
customer visits are removed from a route, but the charging stops coupled with those
visits remain even though they have become unsuitable. When a request is removed,
it is always removed from all the routes in which it is served. For instance, if the
request would be removed from the solution shown in table 5.1, it would be removed
completely, i.e. the visits to nodes two and three would be removed from the route
of vehicle one and the route of vehicle two.

5.3.1 Random Removal
Given a removal percentage, the number of requests to remove is obtained by mul-
tiplying the removal percentage with the total number of requests and rounding.
Then, this many requests are selected randomly and visits where these requests are

39



5. Implementation

served are removed from the solution.

5.3.2 Worst Removal
Using this method, the requests that are placed in the worst positions are removed.
Given a removal percentage, the number of requests to remove is decided. Subse-
quently, each request is removed from the current solution separately and the differ-
ence between the initial objective value and the new objective value is calculated.
In each such iteration, the request is saved along with the calculated difference.
Finally, the requests with the highest difference are removed. To reduce the risk
of the algorithm removing the same requests each time the method is applied to a
certain solution, there is a stochastic aspect when doing the removal: with some
probability, that can be adjusted using the parameter worst removal determinism
rate, the currently worst positioned request will be kept and the next in order will
be considered for removal instead.

5.3.3 Single Route Removal
In this method, a random vehicle that contributes to serving the requests in the
current solution is selected. Then, each request that is partially served by this
vehicle is removed from the solution. This means that every request from a single
vehicle route is removed, and that the routes of other vehicles that partially served
the same requests are affected. In theory, this method can reduce the number of
used vehicles as it may be possible to reinsert the removed requests into other vehicle
routes, and thereby reduce the objective value by avoiding some deployment cost.

5.4 Insertion Methods
After requests have been removed during an iteration, they are reinserted using one
of the insertion methods. The insertions methods include greedy insertion, and two
types of regret insertions which have previously been used in [11], and restricted
greedy insertions, a method introduced in this thesis.

Inserting a request means that a visit where a product is picked up, and a visit where
this product is delivered are inserted into a route such that the delivery happens
after the pickup. Each insertion method compares the best insertion options for the
different requests, and then the request insertion that is the best overall is selected.
This is done repeatedly until all requests have been inserted. The aforementioned
procedure is explained using pseudocode in Listing 5.1 and Listing 5.2. It is worth
noting that the different insertion methods determine which insertion is the best
overall in different manners.

To determine the best insertion option for each request that needs to be inserted,
a certain set of steps is executed (described using pseudocode in Listing 5.3). In
this procedure, every option of adding a pickup visit and delivery visit for a spe-
cific request into the current solution is checked. When considering the option of

40



5. Implementation

inserting pickup and delivery visits at specific positions in a route, the maximum
amount of mass that can be served is calculated so that the mass or volume capacity
of the vehicle is not exceeded at any point in the route. Then, the available loading
and unloading types of the pickup and delivery nodes are compared with the load-
ing/unloading types of the currently considered vehicle, and the shortest loading
and unloading times are determined and selected. Following this, it is checked if
the route created after the insertion would be feasible. If the route would not be
feasible due to an electric vehicle not having enough energy to reach a node, there
is an attempt to repair the route by adding charging, and then feasibility is checked
again. If the route would not be feasible for any other reason or if it cannot be
repaired, the insertion option is disregarded and the next option is tried.

After trying the insertion options for a specific request, feasible insertion options
have been found if they exist. Out of these options, the one where the value ob-
tained by dividing the increase in objective value by the mass carried is the lowest
is saved. The division is done to avoid unfairly favoring insertions where less cargo
mass is transported. The previously described steps are repeated for all the requests
that need to be inserted (as shown in Listing 5.2), and then the insertion options
are compared to determine which insertion should be committed to (which is done
differently depending on the insertion method). Following this, the solution is mod-
ified by committing to the the overall best request insertion, and the cargo mass left
to serve for the selected request is updated. Next, the procedure is done again, and
this continues until all requests that were removed in this iteration are once again
fully served (as shown in Listing 5.1), or until the algorithm cannot proceed due to
there not being any feasible insertion options.

Below, the complete insertion process is described using pseudocode:
1 insertions(solution , requests_to_insert):
2 while requests_to_insert is not empty:
3 single_insertion(solution , requests_to_insert)
4

5 if the request associated with the insertion is fully served:
6 remove the request from requests_to_insert

Listing 5.1: Pseudocode for the entry point of the procedure of inserting all removed
requests.

1 single_insertion(solution , requests_to_insert):
2 overall_best_insertion = None
3

4 for request in requests_to_insert:
5 insertion = get_best_request_insertion(solution , request)
6

7 if insertion better than overall_best_insertion:
8 overall_best_insertion = insertion
9

10 if there is no overall_best_insertion:
11 raise exception
12 else:
13 commit to overall_best_insertion and update the solution

41



5. Implementation

14 update the request served in overall_best_insertion

Listing 5.2: Pseudocode for a function that selects the overall best feasible insertion,
commits to it, and updates the solution.

42



5. Implementation

1 get_best_request_insertion(solution , request):
2 feasible_insertions = []
3

4 for vehicle , route in solution:
5 for i from 1 to route.length ():
6 for j from i + 1 to route.length () + 1:
7 new_route = copy of route
8 mass_change = compute the max amount of mass that can be
9 carried between position i and j in the route

10 loading_time = compute loading time based on the loading
11 options of the vehicle and the pickup node
12 unloading_time = compute unloading time based on the loading
13 options of the vehicle and the delivery node
14

15 pickup_visit = create visit using the pickup node of
16 request , mass_change , and loading_time
17 delivery_visit = create visit using the delivery node of
18 request , -mass_change , and unloading_time
19

20 insert pickup_visit into new_route so that it is placed before
21 the visit currently at index i
22

23 insert delivery_visit into new_route so that it is placed
24 before the visit currently at index j
25

26 check the feasibility of new_route
27 if new_route not feasible due to violated energy constraints:
28 try repairing the route by adding or adjusting charging
29 check the feasibility of new_route
30 if new_route not feasible:
31 continue by checking the insertion for the next i and j
32

33 route_cost_delta = new_route.cost() - route.cost()
34

35 add information about the current insertion to
36 feasible_insertions , including route_cost_delta / mass_change
37

38 if feasible_insertions is empty:
39 raise exception
40 else:
41 sort feasible_insertions ascendingly by route_cost_delta /
42 mass_change
43

44 best_insertion = the first item in feasible_insertions
45 best_insertion.regret_value = compute the regret value using
46 feasible_insertions
47

48 return best_insertion

Listing 5.3: Pseudocode for a helper function that returns the best insertion option
for a specific request, out of the feasible insertion options that were found.

43



5. Implementation

5.4.1 Greedy Insertion
When greedy insertion is used, the feasible insertion with the lowest adjusted route
cost delta is selected as the best insertion during the procedure. In other words, the
insertion where the increase in objective value divided by the mass served in that
insertion is the lowest. Dividing by the mass served in the insertion is significant
because the mass carried affects the total weight of the vehicle, and thereby its
energy consumption. If the route cost delta was not adjusted, there could be a bias
in favor of insertions where less mass is served as these insertions would increase
the energy usage less, so the route cost delta would be lower. This behaviour is
unwanted as it could lead to more frequent transportations of small amounts of
cargo, instead of fewer transportations with greater amounts of cargo.

The adjusted route cost delta for a feasible insertion i of request r can be defined
mathematically as follows: ∆r,i = cr,i − c

mr,i

where cr,i is the route cost that would
be obtained after performing insertion i of request r, c is the route cost prior to
insertion, and mr,i is the mass transported between the pickup node and delivery
node of request r in insertion i. The values ∆r,i that exist for every feasible insertion
i of a request r can be ordered from lowest to highest. Let this ordering of the
values be denoted as ∆1

r, ∆2
r, ∆3

r, ..., ∆k
r . With such notation, greedy insertion can

be defined as a method where the insertion i with the lowest ∆r,i is committed to
for the request r with the lowest ∆1

r.

5.4.2 Regret-2 Insertion
When regret-2 insertion is used, the difference between the adjusted route cost deltas
of the two cheapest feasible insertions is calculated for each request, i.e vr = ∆2

r−∆1
r.

The resulting value for each request r is denoted as its regret value, vr. The request
insertion with the highest regret value is considered the best. If the regret values
of two different request insertions are the same, the insertion with the lowest ∆1

r is
chosen, like in greedy insertion, as a tie breaker. Thus, regret-2 insertion can be
defined as a method where the insertion i with the lowest ∆r,i is committed to for
the request r with the highest vr value.

5.4.3 Regret-3 Insertion
Regret-3 is an extension of the regret-2 insertion method. The difference is that
the regret value vr for each request r is calculated in the following way instead:
vr = ∆2

r−∆1
r +∆3

r−∆1
r. When this method is used, the best request insertion is the

one with the highest vr, and in the case of a tie the greedy insertion criteria apply
similarly as in the regret-2 insertion method. Formally, regret-3 insertion can be
defined as a method where the insertion i with the lowest ∆r,i is committed to for
the request r with the highest vr value. It should be noted that the regret value for
a request can only be calculated in this method if there exist at least three feasible
insertion options for that request. Compared to greedy insertion, the advantage
of regret-2 and regret-3 insertion is that these methods avoid deferring expensive

44



5. Implementation

insertions.

5.4.4 Restricted Greedy Insertions
This method is different from the others as it does not consider all requests when
determining the best insertion. Instead, the requests to insert are shuffled, and each
request is then inserted greedily in order. This approach is nondeterministic and
more time efficient as it considers fewer options. However, it is in theory less likely
to find as good of a solution as the unrestricted greedy approach. The method is
explained using pseudocode in Listing 5.4.

1 restricted_greedy_insertions(solution , requests_to_insert):
2 shuffle the requests in requests_to_insert
3

4 for request in requests_to_insert:
5 while request is not fully served:
6 single_insertion(solution , request)

Listing 5.4: When doing restricted greedy insertions, this procedure is done instead
of the one in the function named insertions. Here, the order in which the requests
are inserted is predetermined and random.

5.5 Checking Feasibility

In the process of checking the feasibility of a route, every adjacent pair of visits
in the route are iterated over. For each such pair, the following is done: First, it
is checked that there exists an arc between the nodes associated with the visits.
Then the elapsed time is updated by adding the travel time and service time at the
second visit, and then it is checked if the time horizon has been exceeded. Next,
the total cargo mass and total cargo volume is updated based on the mass change
in the second visit, and it is checked whether the mass capacity or the volume
capacity of the vehicle is exceeded. Furthermore, the energy balance of the vehicle
is calculated by subtracting the energy consumption of traveling between the nodes
of the currently considered visits, and it is checked if the energy balance goes below
zero. Additionally, it is checked if the vehicle has arrived at a node where charging
cannot be done with a SoC level below the value of a user-defined parameter named
min arrival soc, which can for example be set to 20%. This last condition is not
specified in the mathematical model, but helps avoiding situations where an electric
vehicle gets stuck because it does not have enough energy to reach a nearby charging
station when running the algorithm. Lastly, the charging done during the second
visit is taken into account, and it is checked if the energy balance after charging
exceeds the battery capacity of the vehicle. If any of the checks fails, the route is
considered infeasible, but if it is infeasible due to constraints related to energy, there
may be a possibility of repairing it.

45



5. Implementation

5.6 Repairing Infeasible Routes by Adding Charg-
ing

When a route has been determined to be infeasible due to violating constraints
related to energy, an algorithm is run that tries to repair the route. The algorithm
initiates the energy balance to the battery capacity of the vehicle, assuming that
the vehicle starts fully charged, and iterates over every pair of adjacent visits in the
route. For each pair of visits, where the first visit is referred to as the current visit,
and the second visit is referred to as the next visit, the following is done: First,
charging at the current visit is added to the energy balance, and if the energy of
the vehicle exceeds its battery capacity, the amount of energy charged is adjusted.
Next, the energy consumption for traveling between the nodes is calculated and
subtracted from the energy balance. If the energy consumption was negative, i.e.
there was energy recuperation, the charging done before traveling to the next node
is reduced if possible. Following this, it is checked if the next node was arrived at
with a negative SoC, or if the SoC is below the value of the parameter min arrival
soc and charging cannot be done at the next node. If this is the case, the energy
balance is updated as if the current arc was never traveled, and then a charging
path is searched for that starts at the current node and ends up at the next node. If
such a path is found, the path is inserted into the current route, and the procedure
continues. If at any point the next node could not be reached with an appropriate
SoC, and no charging path could be found, the route could not be repaired and is
deemed to be infeasible.

5.7 Finding a Charging Path
Before the ALNS procedure begins, charging paths are precomputed for each electric
vehicle type. This is done by running an algorithm based on Dijkstra’s algorithm
once for each node in the graph where that node is the starting node. In this
algorithm, the "shortest" paths from the given node to all other nodes in the graph
are searched for, and these paths are referred to as charging paths as they go through
nodes where the vehicle charges fully. The algorithm uses energy consumption as
the "distance", and throughout the search, the battery energy of the vehicle upon
arriving at a node is kept track off, making it possible to decide how much charging
needs to be done to charge fully. During the search, only nodes where charging can
be done are considered, and if a node cannot be reached with a non-negative SoC,
that option is disregarded. Moreover, if charging is possible at the start node, the
vehicle is charged there, which can reduce the number of visits in the charging path.

Once the search concludes, the shortest paths between the provided start node and
the other nodes are extracted. More specifically, these are paths that go between
the nodes such that the energy of the vehicle never becomes negative, charging to
the max is done at every intermediate visit, and the energy consumption is lower
than the other considered options. It should be noted that there may not exist a
charging path between two nodes, for example if charging stations are sparse in the

46



5. Implementation

problem instance.

After the charging paths have been precomputed, they are stored so that they can
be used when trying to repair a route in the main algorithm. When such a situ-
ation arises, and a charging path between a node A and a node B is needed, the
sequence of nodes from the precomputed charging path from A to B is utilized. As
the vehicle may have a different SoC compared to the starting SoC that was used
when precomputing the path, the charging amounts cannot be precomputed and
are therefore calculated when the charging path is needed. This is done by iterating
over the node sequence in the charging path, keeping track of the energy balance,
and adding charging such that the vehicle charges fully at every stop. Once a proper
path consisting of visits where charging amounts are specified has been obtained, it
can be inserted into a route of an electric vehicle in the reparation process.

5.8 The ALNS Procedure
To be able to run ALNS, initial feasible routes are needed. These are obtained by
starting with placeholder route for each vehicle that goes from the start home depot
node to the end home depot node, and then using greedy insertion to insert all
requests. ALNS starts with the solution that consists of these routes and is run for
the selected number of iterations.

In each iterat