**About the Editors**

#### **Leonardo Caggiani**

Leonardo Caggiani is a research fellow at the Department of Civil, Environmental, Land, Building Engineering, and Chemistry of the Polytechnic University of Bari, Italy. He was awarded a Ph.D. in Road and Transportation Systems, Land and Technological Innovation at Polytechnic University of Bari, Italy, in 2009. He received the best Ph.D. thesis award from the Italian Society of Transportation Academics in 2012 and has been rewarded, together with co-authors, with the paper award at the mobil.TUM 2016 conference. He participated in an Italian research project (PRIN2007) on container ports competitiveness in the Mediterranean; he also took part in the Green Intermodal Freight Transport European project (GIFT) of the SEE program. His research interests concern quantitative methods for sustainable mobility, transportation networks design, and equity issues in transportation. Among others, bike-sharing systems are the subject of the most important studies he has carried out.

#### **Rosalia Camporeale**

Rosalia Camporeale is an Associate Senior Lecturer at the Department of Technology and Society, Division of Transport and Roads, Lund University, Sweden. She has expertise in transportation network design problems, sustainable transport, accessibility, and equity issues. During her Ph.D., granted at the Polytechnic University of Bari in March 2017, she developed a methodology to plan and design public transport routes to meet the needs of communities fostering equitable accessibility. Rosalia is currently involved in two JPI Urban Europe projects with a focus on competing and complementary mobility in urban contexts, and she is Project Manager of two national research projects: one dealing with accessibility analysis and distributional effects connected to Swedish railway investments, and the other one concerning the vulnerability of public transport in the COVID-19 era.

### *Editorial* **Toward Sustainability: Bike-Sharing Systems Design, Simulation and Management**

**Leonardo Caggiani 1,\* and Rosalia Camporeale <sup>2</sup>**


#### **1. Introduction**

Bike-sharing systems (BSSs) are a mobility service of public bicycles available for shared use that is becoming increasingly popular in urban contexts. These shared systems provide city users with an alternative, low-carbon, and ecologically sustainable transportation mode (especially suited for short-distance trips) that can significantly reduce traffic congestion, air pollution and noise in city centers, and supports the growth of greener urban environments.

Different issues and challenges have been discussed in previous studies with regard to these systems [1]. Among them are BSS planning and design problems, especially concerning station locations, system simulation and operation problems, such as user demand forecasting and bicycle relocation [2,3]. In this framework, new possible solutions are constantly suggested, each one with its own strengths and weaknesses. Dockless systems (also known as free-floating BSSs) have started to become popular alongside station-based ones, both in big cities and smaller urban environments [4]. At the same time, together with regular bicycles, electric/pedal-assisted bicycles are also being used [5]: in the Vélib' BSS in Paris, for example, a mixed system with both traditional and electric bicycles has recently been implemented [6].

The goal of this Special Issue is to discuss new challenges in the simulation and management problems of both traditional and innovative BSSs, to ultimately encourage the competitiveness and attractiveness of BSSs and contribute to the further promotion of sustainable mobility. We have selected thirteen papers for publication in this Special Issue. Their contributions are summarized and discussed in the following section.

#### **2. Synopsis of the Contributions**

One of the common challenges facing all BSS operators is managing the practical problem of mismatch of bike supply and user BSS demand. To maintain the quality of service to a certain level, these systems need bicycle relocation operations to compensate for imbalances in the network.

Jia et al. (2021) (contribution 1) contribute in this sense by suggesting a new bikesharing rebalancing problem that considers multi-energy mixed fleets and traffic restrictions (aspects mostly neglected in previous studies), using a mixed-integer programming model with the objective of minimizing the total rebalancing cost of the fleet. Their results and sensitivity analysis seem to confirm the efficacy of the algorithm to reduce the total cost associated with BSS rebalancing operations.

A different approach is proposed by Lahoorpoor et al. (2019) (contribution 2), with their bottom-up cluster-based model. They start from an investigation of spatial and temporal patterns of bike-sharing trips, aiming at discovering groups of correlated stations with an agglomerative clustering method. Intra-cluster and inter-cluster rebalancing levels

**Citation:** Caggiani, L.; Camporeale, R. Toward Sustainability: Bike-Sharing Systems Design, Simulation and Management. *Sustainability* **2021**, *13*, 7519. https:// doi.org/10.3390/su13147519

Received: 1 July 2021 Accepted: 3 July 2021 Published: 6 July 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

are considered, and relocation tours are optimized using a single objective genetic algorithm that minimizes the tour length significantly, which ultimately corresponds to a direct cost to the operator and indirect cost to the sustainability of BSSs.

Another crucial issue not often discussed concerns the disorderly parking of freefloating shared bikes. In their study, Jiang et al. (2019) (contribution 3) try to collect, as comprehensively as possible, the causes of such parking behavior through a two-phase questionnaire survey followed by factor analysis. Their investigation, carried out in China (where the problem is particularly acute), aims at facilitating decision-making by governments and enterprises for reference.

Several studies have already attempted to explore the factors that may affect the willingness to use BSS: individual socio-demographic characteristics (gender, age, occupation, education level, monthly income, household bicycle ownership, etc.), individual travel patterns (trip mode, travel time, trip purpose, etc.), transportation infrastructure, land use and built environment characteristics, bike-sharing facilities, and environmental conditions. Different angles and perspectives are presented in the papers collected within this Special Issue.

For instance, the stated preference survey designed and conducted by Politis et al. (2020) (contribution 4) targets car and bus users as well as pedestrians. The results highlight a distinctive set of factors and patterns: the choice of preferred transport mode is most sensitive to travel time and cost of the competitive travel options. According to their findings from Thessaloniki, Greece, BSS seems to be a more attractive option for certain socio-demographic groups and seems to mainly attract bus users and pedestrians rather than car users.

When looking at the infrastructure, the provision of a connected bikeway network has been proven one of the main measures to motivate cycling, since it is directly connected to cyclists' safety. In this regard, Shui and Chan (2019) (contribution 5) propose a novel bikeway design problem that combines a genetic algorithm and an elimination heuristic, and that aims at covering all demand sources and minimizing the total travel time of all cyclists under budget constraints. Their model, tested in two Hong Kong new towns, is not only applicable to new system designs but can also capture the existence of built bikeways and bike stations for system expansion.

Wu and Chen (2019) (contribution 6) support improvement of the night visibility of cyclists by evaluating the differences between shared and private bikes with five types of visibility aids. Their goal is to help policymakers incorporate suitable visibility aids within bike-sharing programs, enhancing the overall traffic safety conditions.

It is also of great importance to understand the motivations and barriers underlying the usage of shared bicycles. The study by Xu et al. (2020) (contribution 7) focuses on free-floating BSSs, adopting an extended theory of planned behavior (TPB) to examine psychological determinants of intention and actual behavior of users. The results, based on an online survey in Beijing, show important implications for planners and lead to several suggestions proposed to support the policymaking of the system.

More specifically, Xiao and Wang (2020) (contribution 8) target as research object of their study the brand choice of bike-sharing in China (namely Hellobike, Mobike, and Ofo). Using a conditional Logic model calibrated on data from an online questionnaire survey, they explore the influence of socio-economic attributes of cyclists and their subjective evaluations, providing a basis for traffic management departments to quantitatively evaluate performances of bike-sharing companies, and assessing the distribution of the total volume among them.

Bardi et al. (2019) (contribution 9) focus on e-bike sharing programs for cruise tourists, an additional niche of operation for bike-sharing systems. They try to understand the major driving forces that lead to the development of these programs, and the major motivating factors for cruise tourists to participate in e-bike sharing services. An ordered probit model is specified to identify the relationship among the variables influencing e-bike sharing

usage and satisfaction of cruise tourists, and interesting interpretations are provided in terms of the relative importance of significant variables.

Most existing studies mainly discuss the relationship between BSSs and external environments, while studies from the perspective of the relationship between internal stations of BSSs are insufficient. Yao et al. (2019) (contribution 10) try to fill this gap with their research. They construct the public bicycle networks of different urban areas (based on real-time data of the Nanjing public BSS) using Gephi software. Using complex network theory and a geographic visualization method, they aim to analyze internal correlation characteristics of BSSs and better understand the station usage.

Considering the large spread of BSSs, it is crucial to gain a better comprehension of the differences between these systems, hence the search for strategies to classify and compare them. An interesting possible approach for clustering different bike-sharing systems around the world can be found in the article by Mátrai and Tóth (2020) (contribution 11). They have gathered data about existing BSSs, grouping them into four categories (public, private, mixed, and other) as the first step for further identification of their common features, which can help to find similar systems and identify problems and best practices in early stages.

Moreover, Caggiani et al. (2021) (contribution 12) suggest a method to evaluate the efficiency of BSSs based on data envelopment analysis (DEA) in order to assist the decisions regarding the performance evaluation of BSS stations. A pool of input and output variables supported by literature, reports, and BSS planning guides is considered, and application to the Malmöbybike system, in Sweden, shows how this approach can provide a reliable evaluation of BSS efficiency.

We conclude the synopsis of this Special Issue's contributions with the study by Nikitas (2019) (contribution 13), which has the ambition to reinvent the formula of longterm success for bike-sharing operations by developing policy and business lessons that will help policymakers and transport providers in establishing and managing these (and other micro-mobility) schemes more effectively. Their findings are supported by primary data from two survey-based studies in Sweden and Greece.

#### **3. What the Future of BSS Holds**

The future of bike-sharing systems is of course unknown, but some speculations are possible based on what has been observed, and what trends seem to be arising.

Technological innovations are definitely contributing to a considerable change in the way of using and owning all kinds of vehicles and having an impact on all transport systems. The idea of geofencing—that is, a virtual boundary around a predetermined area or building [7]—might represent a compromise between traditional station-based and free-floating BSSs, facilitating the benefits and alleviating the challenges associated with these systems [4,8,9]. Designated operating areas to pick up and drop off vehicles could help in overcoming some docked BSS limitations (i.e., insufficient racks or station malfunctions), retaining to a certain extent the parking flexibility provided by free-floating BSSs without hindering pedestrians and/or blocking cycle paths or traffic flows.

A larger differentiation among vehicles can be foreseen. Alongside traditional bike sharing, BSSs with alternative vehicles (mixed-fleet) can attract more users and help satisfy more necessities. One possible option is represented by BSSs using e-bikes, which are superior to conventional bicycles in the ability to traverse longer distances and reach higher speeds, and in greater ease of use, especially over hilly terrains [10]. Another option can be BSSs using traditional and cargo, or e-cargo, bikes. This type of bicycle has recently gained attention as a possible urban mode of transport, particularly for families with children, or to carry heavy shopping or goods [11,12].

The most critical key usage barrier to the future development of BSSs concerns the lack of adequate cycling infrastructure (e.g., bike lanes, cycle paths, parking racks) that, in turn, is directly related to better road safety for cyclists. Poor traffic safety and insufficient bikefriendly infrastructure are the main reasons that cause reluctance to use BSSs (contribution 13, [13,14]).

Also quite relevant is the need for strategic solutions and infrastructure investments that could help in reducing pollution on urban cycle paths. Cycling in downtown areas, especially during the commute, may expose cyclists to air pollutants harmful to human health in large quantities [15,16]. Moreover, because of their physical activity, cyclists often have much higher respiration rates than people who travel by car, and consequently inhale more air pollutants over the same time [17]. The choice of paths is very important to reduce cyclists' exposure to air pollution [18]. Longer cycling routes toward the preferred destination could sometimes significantly lower this exposure. For instance, a recent study done in Coimbra [19] has shown that a 6% increase in distance and time can reduce the exposure to particulate matter and carbon monoxide related to traffic emissions by almost one-third, without requiring any additional physical effort. Hence, it is essential to acquire proper knowledge of the parameters influencing air pollution and noise along cycling facilities to better inform the planning and design of urban bicycle networks [20].

Finally yet equally importantly, there is a need to resolve challenges related to the "new normal" after the COVID-19 pandemic. The role of sustainable transport has been, in a way, redefined, and although cycling per se seems to have had a positive surge [21], the shared use of equipment in BSSs may cause concerns. Micro-mobility systems (such as bike-sharing) can definitely provide a safe alternative transport mode, but it is important that operators expand their efforts in performing and communicating precautionary actions and policies to support community health, to maintain and promote BSSs' roles during the last stages of the pandemic and afterwards [22].

#### **4. List of Contributions**


**Author Contributions:** The authors' contribution is equal. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** Not applicable.

**Acknowledgments:** We would like to thank all the authors who have contributed to this Special Issue for their professional work. We are grateful to the reviewers whose critical support and valuable comments have contributed to improving significantly the quality of the collection. Special thanks go to the editorial assistance office of MDPI for their support throughout the review and publication process of this Special Issue.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


### *Article* **The Bike-Sharing Rebalancing Problem Considering Multi-Energy Mixed Fleets and Traffic Restrictions**

**Yongji Jia \* , Wang Zeng, Yanting Xing, Dong Yang and Jia Li \***

**Abstract:** Nowadays, as a low-carbon and sustainable transport mode bike-sharing systems are increasingly popular all over the world, as they can reduce road congestion and decrease greenhouse gas emissions. Aiming at the problem of the mismatch of bike supply and user demand, the operators have to transfer bikes from surplus stations to deficiency stations to redistribute them among stations by vehicles. In this paper, we consider a mixed fleet of electric vehicles and internal combustion vehicles as well as the traffic restrictions to the traditional vehicles in some metropolises. The mixed integer programming model is firstly established with the objective of minimizing the total rebalancing cost of the mixed fleet. Then, a simulated annealing algorithm enhanced with variable neighborhood structures is designed and applied to a set of randomly generated test instances. The computational results and sensitivity analysis indicate that the proposed algorithm can effectively reduce the total cost of rebalancing.

**Keywords:** sustainable transport; bike-sharing rebalancing problem; multi-energy mixed fleets; traffic restrictions; simulated annealing; variable neighborhood structures

**Citation:** Jia, Y.; Zeng, W.; Xing, Y.; Yang, D.; Li, J. The Bike-Sharing Rebalancing Problem Considering Multi-Energy Mixed Fleets and Traffic Restrictions. *Sustainability* **2021**, *13*, 270. https://doi.org/10.3390/ su13010270

Received: 16 December 2020 Accepted: 27 December 2020 Published: 30 December 2020

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2020 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

#### **1. Introduction**

Nowadays, bike-sharing systems (BSSs), as a low-carbon and sustainable transport mode, are becoming more and more popular across the global, as they can reduce road congestion and decrease greenhouse gas emissions caused by motorized transportation [1,2]. The first BSS was introduced in Amsterdam in 1965 [3] and there are now more than 1500 active BSSs [4] and this number is growing at an increasing rate [5,6]. Recently, some scholars begin to pay attention to the practical problem of mismatch of bike supply and user demands in the BSS, which is a common challenge to all BSS operators [7]. Some operators are trying to meet user demand by placing bikes in cities in large numbers, but this creates congestion on city streets and is not sustainable. In China, the government has introduced policies to restrict operators from placing too many bikes. Thus, these operators have to transfer bikes from the surplus stations to the deficiency stations by means of vehicles so that the BSS can be rebalanced. This problem is known as the Bike-sharing Rebalancing Problem (BRP) [8].

Originally, the BSS operators employ internal combustion vehicles (ICVs) to do the rebalancing operations. However, the development of electric vehicles (EVs) has made great progress over the past few years and many BSS operators are discovering that EVs bring more advantages than ICVs beyond environmental benefits, such as less maintenance, less noise pollution and reduced driving cost [9]. The Chinese government has strongly supported companies to develop sustainable EVs in recent years, providing many policy bonuses for the purchase and use of EVs. Furthermore, more and more cities have implemented traffic restrictions on ICVs to reduce carbon emissions and alleviate traffic congestions [10]. In these cities, ICVs cannot provide services for certain BSS stations in the restricted areas, while EVs can. Therefore, in many practical situations, a multi-energy mixed fleet, consisting of both EVs and ICVs, are used to perform the rebalancing operations.

In this study, the original version of the BRP is extended by assuming the fleet of vehicles to be multi-energy types. The BRP considering Multi-energy Mixed Fleets and Traffic Restrictions (BRP-MMFTR) is considered to be an NP-hard problem, for it originates from the classical Vehicle Routing Problem (VRP). We formulate the BRP-MMFTR as a mixed integer programming model based on one commodity pickup and delivery problem [11], with the objective of minimizing the total rebalancing cost of the mixed fleet. Our model has more complicate constraints than general BRP, such as battery capacity limits for EVs and traffic restrictions for ICVs. It is therefore more complex and more difficult to solve. As the size of bike-sharing stations increases, the calculation time will increase exponentially. Hence, the real-life BRP-MMFTR instances cannot be solved exactly within acceptable computation time.

To handle BRP-MMFTR in a runtime acceptable to BSS operators, we propose the Simulated Annealing algorithm with Variable Neighborhood structures (SAVN) to obtain the optimal solution. To avoid getting into the trap of local optima and enhance the exploratory capability, several variable neighborhood structures are incorporated in our algorithm.

The contribution of this paper lies on the following:


The remainder of this paper is organized as follows. Section 2 presents the literature review on the related problems. Section 3 formally presents a mixed integer programming formulation for BRP-MMFTR. Section 4 describes the procedure and key components of SAVN. Section 5 contains the computational experiments and a sensitivity analysis. Finally, Section 6 presents some concluding remarks about this work.

#### **2. Literature Review**

The BRP has been receiving considerable attention from the literatures during the past decade. Most of studies have focused on optimization models that maximize profits or minimize costs. Dell' Amico et al. provide four mixed integer programming models with the objective of minimizing total costs, where a fleet of capacitated vehicles is employed to relocate the bikes [12]. Erdo ˘gan et al. and Cruz et al. solve the same model proposed by Dell' Amico with a single vehicle, where multi-visit strategy is considered at each station [11,13]. Duan et al. focus on multi-vehicle BRP with the objective of minimizing the total travel distance, and a greedy algorithm is proposed [14]. Casazza et al. and Bulhões et al. incorporate time into the constraints such that each route does not exceed service time limitation [5,15].

All the problems above need to fully meet the rebalancing demands in the BSS, which is difficult to achieve in reality. Hence, some researchers are trying to relax these constraints, and setting the objective to maximize the satisfied rebalancing demand. Papazek et al. set the primary goal to minimize the absolute deviation between target and final fill levels for all stations [16]. Gaspero et al. solve the problem with the aim of minimizing the weighted sum of the total travel time and the total absolute deviation from the target number of bikes [17], while Raviv and Kolka take it as a penalty cost [18]. Faulty bikes are considered by Wang and Szeto with the objective of minimizing the total carbon emission of all vehicles [19]. Usama et al. also consider replacing faulty bikes in the system with the following objectives: User dissatisfaction and vehicle routing costs [20].

Actually, BRP is a large-scale uncertainty problem, due to the existence of uncertainty of user demand [21]. Hence, some scholars have turned their attention to the BRP with uncertain user demand. There are roughly three types of methods for dealing with these uncertainty problems. First, prediction is a common method to resolve uncertainty. Alvarez-Valdes et al. estimate the unsatisfied demand to guide rebalancing algorithms [3]. Zhang et al. propose a dynamic bike rebalancing method that considers both bike rebalancing, vehicle routing and the prediction of inventory level and user arrivals [22]. Second, some scholars divide the uncertainty into multiple independent stochastic scenarios for research. Dell'Amico et al. develop the stochastic programming model with the objective of minimizing the travel cost and the penalty costs for unfulfilled demands [23]. Maggioni et al. propose the two-stage and multi-stage stochastic optimization models to determine the optimal number of bikes to assign to each station [24]. Yuan et al. present a mixed integer programming model with the objective of minimizing the daily costs including the capital cost and the expected operating cost [25]. Third, some studies deal uncertainty directly with dynamic factors. Legros focuses on dynamic rebalancing strategy in the BSS with the objective of minimizing the rate of arrival of unsatisfied users, and solve it by a Markov decision process approach [26]. You develops a constrained mathematical model to deal with a multi-vehicle BRP, aiming to develop dynamic decisions to minimize the sum of the travel costs and unmet costs under service level constraints over a planning horizon [27].

However, none of these literature works considers multi-energy mixed fleets. Goncalves et al. consider the vehicle routing problem with pickup and delivery, and propose a heterogeneous vehicle routing model based on ICVs and EVs [28]. Sassi et al. formulate the heterogeneous electric vehicle routing problem with time dependent charging costs, in which a set of customers have to be served by a mixed fleet of vehicles composed of ICVs and EVs [29]. A mixed fleet of ICVs and EVs are also considered in Goeke and Schneider [30]. The authors formulate the electric vehicle routing problem with time windows and mixed fleet, which is solved through adaptive large neighborhood search algorithm. Charging times vary according to the battery level when the EV arrives at the charging station, and charging is always done up to maximum battery capacity. Macrina et al. present a new variant of the green vehicle routing problem with time windows and propose an iterative local search heuristic to optimize the routing of a mixed vehicle fleet, composed of EVs and ICVs [31].

In summary, there are abundant studies about BSS, but no model considers multienergy mixed fleets and traffic restrictions. For real-life BSSs in big cities, it is more realistic and more useful to consider a fleet of EVs and ICVs. Therefore, in this paper, we focus on a new and practical variant of BRP under the background of multi-energy mixed fleets (composed of EVs and ICVs) and traffic restrictions to ICVs. Our contribution to the development and application of the BSS model is twofold. On the one hand, we formulate BRP-MMFTR as a mixed integer programming model with the objective of minimizing the total rebalancing cost. On the other hand, because this model is too complex to be solved accurately, we develop SAVN to solve it. Our work will expand the existing knowledge on modeling BSS.

#### **3. Model Formulation**

#### *3.1. Model Descripion and Notations*

As shown in Figure 1, BRP-MMFTR studied in this paper can be described as follows: There are two types of vehicles available, namely EVs and ICVs. Vehicles start from the depot with no inventory of bikes, then serve all stations by sequentially loading excess bikes or unloading insufficient bikes, and finally return to the depot after serving all stations. Each station can be served by a vehicle only once. During the rebalancing process, the number of bikes carried by a vehicle cannot exceed its maximal capacity. EVs need to visit the charging stations during the service process. To reduce the complexity of our model, the EV battery must be fully charged at any charging station and the depot. In addition, ICVs are not allowed to enter the traffic restricted area. The objective of BRP-MMFTR is to minimize the total rebalancing cost including the vehicles' fixed costs and traveling energy costs, recharging costs for EVs, and carbon emissions for ICVs.

**Figure 1.** A schematic example of BRP-MMFTR.

The assumptions related to BRP-MMFTR are given as follows:


The notations of sets, parameters, and decision variables used in our model are listed in Table 1.

*k*

*K* =

*K* = ∪

*K* =

=

⊆

⊆


**Table 1.** Description of notations.

#### *3.2. The Total Rebalancing Cost*

The total rebalancing cost of BRP-MMFTR include four parts: Vehicles' fixed costs, vehicles' traveling energy costs, EV recharging costs, and ICVs' carbon emission costs. The vehicles' fixed costs and traveling energy costs are the basic components of objectives in most traditional VRP models. In the objective of our model, the recharging costs reflect the additional charging time during the operation of EVs, while the carbon emission costs are related to the greenhouse gas emission of ICVs.

(1) Fixed costs (C1)

The fixed costs of vehicles are the costs of using vehicles for rebalancing operations. They are different for EVs and ICVs, and the equation of vehicles' fixed costs in our model is shown as follows:

$$\mathbf{C}\_1 = \sum\_{k \in K} \sum\_{j \in N, \ j \neq 0} x\_{0j}^k \cdot F^k \tag{1}$$

#### (2) Traveling energy costs (C2)

The traveling energy costs are composed of two types of costs, that is, fuel costs of ICVs and electricity costs of EVs. Both of them are only associated with the travel distance. And the equation of vehicles' traveling energy costs is shown as follows:

$$\mathbf{C}\_2 = \sum\_{k \in K} \sum\_{(i,j) \in A, i \neq j} \mathbf{C}\_{ij}^k \cdot d\_{ij} \cdot x\_{ij}^k \tag{2}$$

(3) Recharging costs (C3)

The travel distance of an EV is limited by its maximum battery level. Hence, it needs to be charged when its residual charge level lower than the safe residual charge level. The recharging cost is only related to the charging time. The equation of recharging costs is shown as follows:

$$\mathbf{C}\_{\mathfrak{I}} = \sum\_{k \in K\_{\mathfrak{I}}} \sum\_{\mathfrak{e} \in E} \mathbf{C}\_{\mathfrak{w}} \cdot t\_{\mathfrak{e}}^{k} \cdot y\_{\mathfrak{e}}^{k} \tag{3}$$

#### (4) Carbon emission costs (C4)

In the process of rebalancing, a large amount of CO<sup>2</sup> is generated by ICVs from their fuel consumptions, resulting in greenhouse effect. By reducing the costs of carbon emissions, the total rebalancing cost is reduced. The equation of carbon emission costs is shown as follows:

$$\mathbf{C}\_{4} = \sum\_{k \in K\_{1}} \sum\_{(i,j) \in A, i \neq j} L \cdot \mathbf{C}\_{L} \cdot d\_{ij} \cdot \mathbf{x}\_{ij}^{k} \tag{4}$$

#### *3.3. Model Establishment*

The mixed integer programming model for BRP-MMFTR can be written as follows:

$$\begin{aligned} \min Z &= \quad \mathsf{C}\_{1} + \mathsf{C}\_{2} + \mathsf{C}\_{3} + \mathsf{C}\_{4} = \sum\_{k \in K} \sum\_{j \in N, j \neq 0} \mathsf{x}\_{0j}^{k} \cdot F^{k} + \sum\_{k \in K} \sum\_{(i, j) \in A, i \neq j} \mathsf{C}\_{ij}^{k} \cdot d\_{ij} \cdot \mathsf{x}\_{ij}^{k} \\ &+ \sum\_{k \in K\_{2}} \sum\_{e \in E} \mathsf{C}\_{w} \cdot t\_{e}^{k} \cdot y\_{e}^{k} + \sum\_{k \in K\_{1}} \sum\_{(i, j) \in A, i \neq j} \mathsf{L} \cdot \mathsf{C}\_{L} \cdot d\_{ij} \cdot \mathsf{x}\_{ij}^{k} \end{aligned} \tag{5}$$

*Subject to*:

$$\sum\_{j \in N, j \neq 0} x\_{0j}^k = \sum\_{i \in N, i \neq 0} x\_{i0}^k = 1 \quad \forall \ k \in K \tag{6}$$

$$\sum\_{k \in \mathcal{K}} \sum\_{j \in \mathcal{N}, j \neq i} \mathfrak{x}\_{ij}^{k} = 1 \,\,\forall \, i \in \mathcal{N} \tag{7}$$

$$\sum\_{\substack{j \in \mathcal{N}, j \neq i}} \mathfrak{x}\_{ij}^{k} - \sum\_{j \in \mathcal{N}, j \neq i} \mathfrak{x}\_{ji}^{k} = \mathbf{0} \,\,\forall \, i \in \mathcal{N} \,\forall \, k \in \mathcal{K} \tag{8}$$

$$\sum\_{i \in S} \sum\_{j \in S} \mathbf{x}\_{ij}^k \le |\mathcal{S}| - 1 \,\,\forall \, \mathcal{S} \subseteq \mathcal{N}\_\prime \,\, |\mathcal{S}| \ge 2 \tag{9}$$

$$\sum\_{\substack{j \in D, j \neq i}} \mathfrak{x}\_{ij}^{k} = 0 \,\,\forall \, i \in N \,\forall \, j \in D \,\forall \, k \in K\_1 \tag{10}$$

$$\max \{ \mathbf{0}, \mathbf{G}\_{\dot{\mathbf{i}}}, \mathbf{G}\_{\dot{\mathbf{j}}} \} \quad \cdot \ x\_{\dot{\mathbf{i}}\dot{\mathbf{j}}}^{k} \le g\_{\dot{\mathbf{i}}\dot{\mathbf{j}}}^{k} \le \min \left\{ \mathbf{Q}^{k}, \mathbf{Q}^{k} + \mathbf{G}\_{\dot{\mathbf{i}}}, \mathbf{Q}^{k} - \mathbf{G}\_{\dot{\mathbf{j}}} \right\} \quad \cdot x\_{\dot{\mathbf{i}}\dot{\mathbf{j}}}^{k} \forall \mathbf{i}, \mathbf{j} \in \mathbf{N}, \, \mathbf{i} \ne \mathbf{j} \tag{11}$$

$$\sum\_{k \in K} \sum\_{\substack{j \in N, j \neq i}} g\_{ij}^k - \sum\_{k \in K} \sum\_{j \in N, j \neq i} g\_{ji}^k = G\_i \,\forall \, i \in N \tag{12}$$

$$0 \le \mathcal{g}\_{ij}^k \le \mathcal{Q}^k \cdot \mathfrak{x}\_{ij}^k \quad \forall \, i, j \in \mathcal{N}, \, i \ne j, \forall \, k \in \mathcal{K} \tag{13}$$

$$t\_e^k = y\_e^k \cdot \left[ \left( \mathbb{R} - u\_{ek}^1 \right) / \lambda \right] \quad \forall \, k \in \mathbb{K}\_2, \forall \, e \in E \tag{14}$$

$$0 \le u\_{\vec{j}k}^1 \le u\_{i\vec{k}}^2 - r \cdot d\_{i\vec{j}} \cdot x\_{i\vec{j}}^k + \mathcal{R} \cdot \left(1 - x\_{i\vec{j}}^k\right) \quad \forall \, i, j \in \mathcal{N}, \, i \ne j \,\forall \, k \in \mathcal{K}\_2 \tag{15}$$

$$\mathbf{r}\_{ik}^{1} \ge \boldsymbol{\mu} \cdot \boldsymbol{\mathcal{R}} \quad \forall \, i \in \mathcal{N}\_{\prime} \forall \, k \in \mathcal{K}\_{2} \tag{16}$$

$$
\mu\_{0k}^2 = R \quad \forall \ k \in K\_2 \tag{17}
$$

$$
\mu\_{\rm ek}^2 = y\_{\rm e}^k \cdot \mathbb{R} \quad \forall \, k \in \mathbb{K}\_2 \,\forall \, e \in E \tag{18}
$$

$$
\mu\_{ik}^1 = \mu\_{ik}^2 \quad \forall \, i \in N \/ \forall \, k \in K\_2 \tag{19}
$$

$$x\_{ij\prime}^k \ y\_e^k = [0, 1] \quad \forall i \in N \forall j \in N \forall k \in K \forall e \in E \tag{20}$$

Objective function (5) minimizes the total rebalancing cost, including vehicles' fixed costs, vehicles' traveling energy costs, EVs' recharging costs, and ICVs' carbon emission

*u*

costs. Constraints (6) ensure that each vehicle starts at the depot and returns to the depot at the end of its route. Constraints (7) guarantee that each station is served exactly once. Constraints (8) refer to the usual flow conservation. Constraints (9) can avoid subtours and thus guarantee route-connectivity. Constraints (10) indicate that ICVs are restricted from entering the restricted area. Constraints (11) give the upper and lower bounds of the number of bikes loaded by vehicles. Constraints (12) and (13) ensure that the vehicle's maximum capacity is not exceeded. Constraints (14) and (15) are the EVs' charging functions and power consumption functions, respectively. Constraints (16) indicate the safe residual charge constraints of EVs. Constraints (17) and (18) show that EVs are fully charged when leaving the depot or a charging station. Constraints (19) guarantee the residual charges of EVs are the same after serving a BSS station. Constraints (20) are the binary decision variables.

#### **4. Simulated Annealing Algorithm with Variable Neighborhood Structures**

Simulated Annealing (SA) is a heuristic method to solve various NP-hard optimization problems. It can expand the exploration capability by accepting worsen solutions with some probability. This has benefit to reduce the probability of getting trapped in local optima. In order to improve the search efficiency and get solutions with higher quality, we introduce the variable neighborhood structures [32] into the framework of SA, and propose SAVN to deal with the real-life BRP-MMFTR instances. The procedure and key components of our algorithm are described in detail below.

#### *4.1. The Procedure of Our Algorithm*

Given an initial solution, SAVN starts from the initial temperature *T*0. During the search process, the algorithm randomly selects a neighborhood structure to transform the current solution *S* into a randomly generated feasible neighbor *S* ′ . Note that the cost of a feasible solution *S*, namely *Z*(*S*), is evaluated with the Equation (5). If the cost of *S* ′ is less than the cost of *S*, *S* ′ will definitely be accepted. Otherwise, the acceptance probability of *S* ′ is *p* = *exp* − *Z*(*S* ′ )−*Z*(*S*) *T* , where *T* is the current temperature. For each *T*, this process is performed *Len* times. Then, *T* decreases by multiplying the cooling rate *α*.

Repeating the above processes until the stop criterion is met, that is, the unimproved number of the best solution so far reaches the pre-specified MaxUN. Furthermore, when *T* is less than 0.01, *T* will increase to help the algorithm escape from the local optima [33]. We double *T<sup>b</sup>* first (but no more than *T*0) and then set *T* = *T<sup>b</sup>* . The pseudocode of SAVN is shown in Algorithm 1.

#### *4.2. Initial Solution Generation*

Considering the particularity of BRP-MMFTR, we firstly classify all BSS stations into two types by region: stations in restricted areas and stations in non-restricted areas. Then a three-step algorithm is proposed to generate the high-quality initial solution.

Step 1: For stations in the restricted areas, EVs must be employed. We use the insert algorithm to construct the initial routes for these stations. First, the station with the largest distance from the depot is selected as the first customer of an EV. Then, if the capacity constraints are satisfied, the other stations are inserted into the current EV route in turn. Otherwise, a new EV route will be generated and serve the remainder stations. Repeat this insertion process until all stations in the restricted areas are serviced.

Step 2: For the stations in non-restricted areas, we first judge whether a station can be inserted into the existing EV routes. For the stations that can be inserted, we insert them to the positions with the lowest incremental cost in the existing EV routes. For the other stations, we employ the same insert algorithm to generate some new ICV routes.

Step 3: Charging stations are allocated to EV routes. If the residual charge level of an EV is enough, this step can be skipped. Otherwise, when the residual charge level of an EV is lower than the designated safe-residual-charge level, the nearest charging station will be inserted into the EV route to increase its mileage.


#### *4.3. Neighborhood Structures*

Four neighborhood structures are employed in our algorithm. In the procedure of SAVN, when a new solution is needed, one of these four structures is randomly selected with equal probability. As shown in Figure 2, all structures are described by an instance with 5 stations.

Neighborhood structure (a): A randomly selected station is relocated to another position in the current solution.

Neighborhood structure (b): Two randomly selected stations exchange their positions.

Neighborhood structure (c): Two randomly selected segments exchange their positions, and their lengths are at most 3.

Neighborhood structure (d): A pair of stations is randomly selected and all the stations between them, including themselves, are reversed.

#### **5. Computational Results**

#### *5.1. Parameters and Experimental Data*

We have implemented SAVN using MATLAB R2018a, and all computational experiments have been carried on an Intel Core i5-8259U with 2.3 GHz CPU and 8 GB RAM running the macOS Catalina operating system.

The test instances of BRP-MMFTR are generated randomly. The BSS stations' coordinates are randomly generated in the range of abscissa [0,200] and ordinate [0,150]. The coordinate range of the traffic restricted area is the abscissa [80,160] and the ordinate [30,60]. The rebalancing demand of each BSS stations is randomly generated between [−20,20]. The depot coordinate is (100,75).

Table 2 shows an example instance with one depot, 18 BSS stations and five charging stations. And Table 3 displays the parameters of EVs and ICVs. The other parameters used in BRP-MMFTR instances are as follows: R = 75 kW·h, *r* = 0.5 kW·h/km, *µ* = 0.3, *λ* = 0.5 kW·h/min, *C<sup>w</sup>* = 0.4 CNY/min, *C<sup>L</sup>* = 0.06 CNY/kg [34], and *L* = 5 kg/km. Through preliminary tests, the parameters of SAVN are set as follows: *T*<sup>0</sup> = 50, *α* = 0.96, MaxUN = 50, and *Len* = *n* 2 (where, *n* is the number of BSS stations).

#### *5.2. Comparisons of Computational Results*

To verify the effectiveness of our SAVN, we compare it with SA and VNS (Variable Neighborhood Search) on the instance depicted above. Note that both SA and VNS use the same parameters of our SAVN. All the three algorithms were executed 10 times to weaken the randomness of the heuristic algorithm.

Figure 3 illustrates the vehicle routes obtained by SA, VNS and SAVN, respectively. The red lines represent the route of the ICV, while the blue lines represent the route of the EV. It can be observed that the solution obtained by SAVN has fewer detours, which make the vehicle's travel distance smaller than those obtained by SA and VNS. Furthermore, its traveling energy costs, recharging costs and carbon emission costs can be also decreased.

*X***-Axis** *Y***-Axis Demand Index** Depot 100 75 0 0 BSS stations 86 143 −6 1 75 35 3 2 86 57 6 3 99 36 11 4 132 25 6 5 138 44 −7 6 108 18 −9 7 100 98 −15 8 172 60 6 9 123 55 3 10 110 128 16 11 34 35 8 12 161 112 −14 13 140 103 7 14 191 71 −12 15 76 101 7 16 48 118 5 17 19 65 −4 18 Charging stations 30 66 C1 48 35 C2 126 73 C3 52 101 C4 145 27 C5

**Table 2.** Information of the depot, BSS stations and charging stations.



**Figure 3.** *Cont*.

**Figure 3.** Optimal vehicle routes obtained by three algorithms. (**a**) SA (**b**) VNS (**c**) SAVN.

Figure 4 displays the comparison results of these algorithms. The X-coordinate is the number of iterations and the Y-coordinate is the total rebalancing cost. It can be seen that the final solution obtained by SAVN is better than SA and VNS, although in the early stage of the search process, SAVN may be inferior to VNS. This is due to the better optimization performance of SAVN, which can better escape from the trap of local optima. Hence, our SAVN can be regarded as an effective algorithm to solve BRP-MMFTR.

**Figure 4.** Comparison of total rebalancing cost.

Table 4 lists the components of the total rebalancing cost in the BRP-MMFTR optimal solution obtained by SA, VNS, and SAVN. Obviously, for these three algorithms, SAVN is the best and SA is the worst for the total rebalancing cost. Although the electricity costs of EV obtained by SAVN is slightly worse than that of VNS, the fuel cost of ICV has been drastically reduced. In this way, SAVN's traveling energy cost is less than VNS, because it is the sum of the electricity cost of EV and the fuel cost of ICV. Similarly, although SAVN's recharging costs is higher than VNS for its EV has a longer travel distance, its carbon emission costs is lower. For practical BSS operators, SAVN's solution is significantly better than VNS, because it employs EVs to perform more rebalancing operations. As we know, EVs are environmentally friendly and they will replace ICVs in the near future.



Traveling energy cost is the sum of Electricity cost of EV<sup>1</sup> and Fuel cost of ICV<sup>2</sup> .

#### *5.3. Computational Results of More Instances*

To further verify the effectiveness and universality of SAVN, we tested the three algorithms in further instances. Theywere evaluated by nine test instances (three types, three instances per type), displayed in Table 5. The three types are small-size (*n* = 20), medium-size (*n* = 50) and large-size (*n* = 100). The instances are also generated randomly according to the method introduced in Section 5.1.

**Table 5.** Nine instances of three types.


To provide reliable statistics, each algorithm is executed 10 times for each instance. And the average CPU time (t/s), the average total rebalancing cost (TRC) and the standard deviation (Sd) are displayed in Table 6. The best values are marked with bold fonts.

**Table 6.** Comparative results on nine instances.


The best values are marked with bold.

From Table 6, we can draw some conclusions as follows: For the average CPU time, SA has obvious advantages, but its solution quality is the worst. In addition, the calculation time of SAVN and VNS is difficult to distinguish the pros and cons. For the average

total rebalancing cost, the solution obtained by SAVN can reduce 9.2% compared to SA or 3.4% compared to VNS in the small-size instances, 10.7% or 2.7% in the medium-size instances, and 15% or 3% in the large-size instances. With the number of stations increases, the solution quality of SAVN is getting better and better. Hence, from the perspective of solution quality and CPU time, SAVN can obtain a better solution than SA and VNS.

Furthermore, the standard deviation can reflect the stability of the algorithm. And the smaller the standard deviation is, the better the stability of the algorithm is. From Table 6, it is observed that the standard deviation of SAVN is the smallest in all the instances. Therefore, SAVN is more stable than the other two algorithms.

#### *5.4. Discussion on the Value of µ*

Parameter µ (the coefficient of safe-residual-charge level of EV) directly affects the charging time and then the recharging costs and the total travel distances. Hence, it is important to set a suitable value for µ. If the value of µ is set too low, there will be a probability of failing to drive to the nearest charging station. If the value of µ is set too high, the EV may frequently go to the charging station, which will cause a waste of energy and charging time. In this section, the sensitivity analysis of µ is analyzed using the 18-station instance in Section 5.2. The different values of µ are from 0.15 to 0.5 in increment of 0.05.

The SAVN is run 10 times for each µ and the trends of change with the increase of µ are illustrated in Figure 5. Obviously, the total rebalancing cost shows a clear trend of decreasing and then increasing. µ = 0.3 is the point with the lowest total rebalancing cost. In summary, parameter µ plays an important role in BRP-MMFTR and its value should be reasonably determined according to the layout of charging station system. If the quantity of charging stations is sufficient, the value of µ can be appropriately set lower. Otherwise, it should be set higher. Therefore, setting an appropriate value for µ is helpful to reduce the total rebalancing cost.

**Figure 5.** The trend of TRC changes with the increase of *µ*.

#### **6. Conclusions**

As a low-carbon and ecologically sustainable transportation mode, BSS has become a way to deal with the growing menace of global warming. In this study, we have proposed, modeled and solved BRP-MMFTR, which is a variant of BRP considering multi-energy mixed fleets and traffic restrictions. We first formulate BRP-MMFTR as a mixed integer programming model with the objective of minimizing the total rebalancing cost composed of vehicles' fixed costs, vehicles' traveling energy costs, EVs' recharging costs, and ICVs' carbon emission costs. Then SAVN is designed to solve this model, which is the simulated annealing algorithm enhanced with variable neighborhood structures. To illustrate the efficiency and efficacy of our algorithm, some test instances of BRP-MMFTR are generated randomly. The computational results reveal the huge advantage of SAVN, compared with SA and VNS. SAVN can achieve better solution in terms of solution quality and CPU time, outperforming those obtained by SA and VNS. In addition, SAVN is more stable than SA

and VNS. Finally, the sensitivity analysis results of parameter µ indicate that as µ increases, the total rebalancing cost shows a clear trend of decreasing and then increasing. Therefore, setting an appropriate value for µ is helpful to reduce the total rebalancing cost. In addition, the value of µ is not necessarily constant in the practical rebalancing operations, and it can be dynamically adjusted according to the real-time conditions.

For BSS operators, we provide the optimal vehicle scheduling suggestions for multienergy mixed fleets to minimize the total rebalancing cost when bike supply and user demand are not matched. In addition, we also focus on carbon emissions during the rebalancing process. BSS is originally a sustainable transportation mode, but if their operators use ICVs during rebalancing operations, the low-carbon benefits of BSS will be partially offset by the carbon emissions of ICVs. Therefore, it is obviously beneficial to create a green and low-carbon sustainable transportation system by using EVs instead of ICVs. The government should vigorously promote the sustainable development of EVs. Some policies, such as traffic restrictions on ICVs, can be introduced to encourage the BSS operators to purchase more EVs to minimize the carbon emissions.

Future research can extend our model and algorithm to solve more complex variants of BRP-MMFTR, such as considering the impact of average speed and speed variations to the traveling energy costs. Furthermore, dynamic or stochastic BRP-MMFTR is also a good research direction.

**Author Contributions:** Conceptualization, Y.J.; methodology, Y.J. and J.L.; software, W.Z. and Y.X.; validation, Y.J. and J.L.; writing—original draft preparation, Y.J. and W.Z.; writing—review and editing, D.Y. and J.L.; visualization, W.Z. and Y.X. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research was funded by a grant from Shanghai Planning Office of Philosophy and Social Science, grant No. 2018BGL018.

**Institutional Review Board Statement:** Not applicable.

**Informed Con sent Statement:** Not applicable.

**Data Availability Statement:** Data is contained within the article.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**

