As the core equipment of the horizontal transport system in ACTs, the stable and efficient operation of AGVs is paramount. The AGV scheduling problem, which focuses on optimally allocating limited vehicle resources to pending container transport tasks, has been extensively studied. Traditional research primarily targets the minimization of makespan or delay time through various methodological approaches, including mathematical programming and heuristic algorithms. For instance, scholars have conducted in depth research from the perspectives of job assignment and scheduling optimization (Che et al., 2024; Guan et al., 2025), as well as conflict-free path planning (Cao et al., 2023). These studies established the foundational framework for AGV operations but often treated batteries as perpetually available resources, overlooking the constraints imposed by energy replenishment.

With the maritime industry’s paradigm shift toward green logistics and decarbonization, energy management has become a critical area in AGV operations. Recent literature has moved beyond simple efficiency objectives to incorporate energy constraints and carbon emission reduction. For example, Zhang et al. (2025b) proposed low-carbon scheduling methods for electric transportation vehicles, emphasizing the trade-off between operational efficiency and energy consumption. Zhou et al. (2024) modeled the charging decision as a Markov decision process to manage AGV charging sequences. Furthermore, for the infrastructure level, Zhen et al. (2025) highlighted the operational advantages of the battery swapping mode over traditional plug-in charging, noting its ability to maximize vehicle availability through rapid energy replenishment. Zhao and Liang (2025) proposed a novel shallow charge and shallow discharge charging strategy based on idle time, which significantly moves beyond simple full charge and full discharge cycles to preserve battery life.

However, despite the rich literature on short-term operational scheduling, there is a disconnect between daily operations and long-term asset management. Most existing studies treat the capacity of BSS and the number of batteries as fixed parameters or infinite resources during the scheduling phase. Systematic research on the strategic level remains relatively scarce.

Regarding terminal equipment investment and asset management, significant work has focused on the electrification and retrofitting of yard cranes (YCs), which are major contributors to terminal emissions. To actively respond to the green port initiative, scholars have developed optimization models to balance the trade-off between capital expenditure and environmental benefits. For instance, Tan et al. (2022) proposed a strategic decision method for YC transformation, constructing a model to minimize the total investment and operation costs while meeting specific carbon emission reduction targets. Similarly, Ding et al. (2021) addressed the transition from diesel-powered to electric-driven rubber-tyred gantry cranes. They emphasized that retrofitting decisions must be integrated with operational planning, developing a model to determine the optimal investment timing and deployment to maintain yard handling capabilities during the transition period. Peng et al. (2016) explored the optimal allocation of resources for YC networks. Rather than focusing on individual equipment, they modeled the energy replacement problem at a network level to minimize total carbon dioxide emissions, demonstrating the importance of systemic resource planning.

In addition to YCs, the investment and configuration of horizontal transport equipment have also attracted attention, particularly under complex operational constraints. Zhang et al. (2024a) investigated the retrofitting and deployment of yard trucks for hazardous material transportation. Their study highlights the complexity of managing a mixed fleet—comprising diesel, electric, and LNG trucks—where investment decisions must satisfy strict safety regulations alongside environmental goals.

Despite these advancements in equipment investment, a critical gap remains regarding the energy infrastructure itself, specifically BSS. Existing research primarily treats equipment (like cranes and trucks) as independent assets to be purchased or retrofitted. However, BSS investment involves a set of multi-period challenges. It requires managing a pool of circulating assets subject to performance battery degradation.

From a methodological perspective, mathematical programming is the dominant approach in terminal equipment investment research. Deterministic optimization models, such as those employed by Tan et al. (2022) and Ding et al. (2021), effectively handle scenarios with fixed parameters, determining optimal investment timing and configuration schemes through exact or heuristic algorithms. However, given that investment problems often span long time horizons, handling uncertainty is crucial. Neglecting these stochastic factors may lead to rigid strategies that fail to adapt to actual demand variations, resulting in either substantial capital waste from over-investment or service interruptions due to supply shortages. To address this, stochastic programming has been widely applied. For example, Zhang et al. (2024b) constructed a two-stage stochastic MILP model to handle uncertainties in workloads and carbon quotas. Similarly, Peng et al. (2016) combined mathematical resource allocation models with simulation to cope with high uncertainties in the terminal transportation network. These stochastic programming methods often rely on precise assumptions about probability distributions of uncertain parameters, which are difficult to obtain accurately in practice. This reliance can lead to suboptimal decisions if the assumed distribution is different from reality.

To overcome these limitations, DRO has recently received widespread attention due to its robustness against distribution deviation (Mohajerin Esfahani and Kuhn, 2018). The DRO method constructs an ambiguity set for uncertain parameters, possessing both the interpretability of mathematical models and the generalizability of statistical methods. It has been widely applied in fields such as transportation systems (Zhang et al., 2018; Ghosal and Wiesemann, 2020; Zhang et al., 2021b) and terminal operations. For instance, Agra and Rodrigues (2022) pioneered its application in berth allocation using a Wasserstein-based ambiguity set to address vessel handling time uncertainty. Subsequently, Rodrigues and Agra (2024) extended this methodology to quay crane scheduling, proving that robust models can effectively balance operational costs and reliability.

Motivated by these advancements, this paper adopts the DRO method to handle the uncertainty of battery swapping demand. This method summarizes the patterns of battery swapping demand from historical data, which can effectively prevent the overfitting of theoretical models and improve the generalizability of the model. The method proposed in this paper provides scientific support for the long-term investment decision of battery assets in BSSs, bridging the gap between theoretical robustness and practical application.

The mathematical symbols used in this paper are defined as follows, based on the problem description in the Section 3.

 Sets:

K        Set of battery age indices, K={0,1,2,…,Klife−1}, where k represents the cumulative usage cycles of the battery;

T        Set of full time indices, T={0,1,…,T,T+1}, where T is the length of the planning horizon;

Tdec    Set of decision period indices, Tdec={0,1,…,T};

Toper   Set of operational period indices, Toper={1,2,…,T};

Ω   Set of demand uncertainty scenarios, |Ω| denotes the number of scenarios.

 Parameters:

dbuy      Procurement cost of a single battery;

dhold      Maintenance cost of a single battery per period;

Klife      Battery lifecycle threshold;

ξt           Battery swapping demand in period t;

ξω          Battery swapping demand in scenario ω ∈ Ω;

Bt          Budget upper limit in period t;

B           Total budget upper limit for the entire planning horizon;

r             Annual discount rate (for NPV calculation);

ρ            End-of-period asset conversion rate;

N           Number of periods in a year;

F(θ)     Distribution ambiguity set based on Wasserstein distance, with radius parameter θ;

Є          Maximum violation probability for the distributionally robust chance constraint;

y¯k,0         Number of batteries with k cumulative usages at the initial time.

 Decision Variables:

xt         Number of new batteries purchased at the end of period t;

yk,t         Number of batteries with k cumulative usages owned at the beginning of period t;

uk,t         Number of batteries with k cumulative usages utilized in period t.

To address the uncertainty of battery swapping demand without requiring precise knowledge of its underlying probability distribution, we employ a data-driven approach using Wasserstein distance-based ambiguity sets. The Wasserstein ambiguity set encompasses all probability distributions that lie within a specified statistical distance from the empirical distribution constructed from historical samples, where the radius parameter θ controls the level of robustness.

Let ξ˜ denote the random variable representing the uncertain battery swapping demand. We construct the ambiguity set F(θ) using the following formulation:

Equations 1–5 establish the mathematical framework for our ambiguity set construction. Equation 1 defines the ambiguity set F(θ) by incorporating all probability distributions ℙ supported on U that reside within a Wasserstein distance θ of the empirical distribution ℙ^. The distance metric is formalized in Equation 2, which utilizes the absolute deviation |ξ˜−ξ^| to quantify the cost of transporting probability mass from the empirical distribution to the true distribution. Equations 3 and 4 specify that the uncertain demand ξ˜ follows the candidate distribution ℙ, while the historical samples ξ^ follow the empirical distribution ℙ^. Finally, Equation 5 defines the joint distribution ℙ¯ within the set P(ℙ,ℙ^), representing all valid joint distributions with marginals ℙ and ℙ^. The support set U ensures the physical meaningfulness of the demand values (e.g., non-negativity).

To accurately model this complex decision problem, the following reasonable assumptions are made:

Regarding Assumption 2, while battery degradation is inherently non-linear, linear simplification is adopted to ensure the tractability of the strategic optimization model. Furthermore, since the procurement decisions are primarily driven by the maximum service life rather than the intermediate degradation process, this approximation has a negligible impact on the optimal investment strategy, affecting only the minor book value in the end period.

The optimization problem for the battery replacement strategy, considering demand uncertainty, can be formulated as the following model with distributionally robust chance constraints:

The objective function (Equation 6) aims to minimize the total NPV cost of the BSS, including three parts. The first term represents the present value of new battery procurement costs. The second term is the present value of battery holding costs, calculated based on the total number of batteries owned in each period and discounted. The third term represents the book value of the remaining batteries at the end, due to asset recovery. Constraints (Equation 7) are the distributionally robust chance constraints, ensuring the reliability of the battery swapping service under demand uncertainty. This constraint is based on the Wasserstein ambiguity set F(θ), requiring that the total supply capacity of the BSS meets the uncertain demand with a probability of at least 1−Є under the worst-case distribution. Compared to traditional deterministic demand constraints, this chance constraint maintains an acceptable service level by controlling the violation probability Є. Constraints (Equation 8) limit the usage quantity of each age group to no more than the number of batteries owned in that group. Constraints (Equation 9) state that the new batteries in the next period originate from the batteries purchased in the current period. Constraints (Equation 10) are the core constraints for age evolution. At the beginning of period t, there are yk,t batteries with k cumulative usages, of which uk,t are used, and their cumulative usage increases to k+1; the unused yk,t−uk,t batteries retain k cumulative usages; meanwhile, uk−1,t batteries that had k−1 cumulative usages in the previous period become k cumulative usages after being used. Therefore, the number of batteries with k cumulative usages at the beginning of period t+1 equals the sum of k-aged batteries unused in this period and (k−1)-aged batteries used in the previous period. Constraints (Equation 11) are the per-period budget constraint, limiting the sum of procurement expenditure and holding costs in each operational period ( t≥1) to no more than the budget upper limit Bt. The initial procurement at t=0 is treated as a one-time construction investment and is excluded from this operational budget constraint. Constraints (Equation 12) are the total budget constraint, limiting the cumulative procurement expenditure and holding costs over the entire planning horizon. Constraints (Equation 13) set the initial age distribution based on the pre-existing inventory y¯k,0. Constraint (Equation 14) stipulates that no new batteries are purchased during the liquidation period. Constraints (Equations 15–17) define the domain of the decision variables, where all variables are non-negative integers.

To ensure the reliability of battery swapping services, we require that the probability of supply meeting demand is at least 1−Є under the worst-case distribution within the ambiguity set. This is formally expressed as the following distributionally robust chance constraint (7). Directly solving it is computationally intractable due to the non-convexity of the chance constraint and the infinite dimension of the ambiguity set. To address this, we introduce the CVaR as a convex conservative approximation. Based on the fundamental property that CVaR serves as a rigorous upper bound for the Value-at-Risk, we can transform the original chance constraint into a mathematically tractable form through the transformation shown in Equation 18:

The implication in step 3 holds because CVaRЄ(Z)≤0 is a sufficient condition for ℙ(Z>0)≤Є, providing a safe approximation for the original problem (Rockafellar and Uryasev, 2000).

To further simplify the constraint and decouple the decision variables from the uncertainty, we exploit the translation invariance property of coherent risk measures.

Equation 19 shows that the total battery supply in each period must exceed the worst-case CVaR of the random demand.

It’s difficult to compute the worst-case CVaR of the uncertain demand ξ_t over the infinite-dimensional Wasserstein ambiguity set. Based on strong duality, we derive a closed-form expression that depends only on the historical data and the robustness radius.

Theorem 1. Consider the uncertainty set Ωwith |Ω|historical scenarios. Let the observed demand values be denoted by ξ^ωfor each scenario ω ∈ Ω. Under a confidence level ∈and a Wasserstein ambiguity set F(θ)with radius θ, if the historical samples are sorted in descending order such that ξ^1≥ξ^2≥⋯≥ξ^|Ω|, the worst-case CVaR can be calculated as:

where ∑ω=10ξ^ω=0and ξ^|Ω|+1=0are defined to handle boundary cases. The proof is provided in Appendix 1.

Theorem 1 decomposes the worst-case risk into three interpretable components: 1) The average value of the worst scenarios, capturing the tail risk inherent in the historical data; 2) A fractional adjustment term for the boundary scenario, ensuring the continuity of the risk measure; 3) A robustness term, where a larger Wasserstein radius imposes a stricter penalty for distributional deviation, offering stronger protection.

Through Theorem 1, we obtain the tractable formulation. The original distributionally robust chance constraint (7) is replaced by the Equation 21.

This transformation eliminates the infinite-dimensional constraints. The resulting model consists solely of a linear objective function and linear constraints, allowing it to be efficiently solved to global optimality using standard commercial solvers such as Gurobi or CPLEX.

The computational complexity of the transformed model is primarily determined by the planning horizon T and the battery lifecycle resolution Klife. Specifically, the number of integer decision variables scales linearly with the product of T and Klife ( O(T·Klife)). Although the problem size naturally expands as the planning horizon extends or the age becomes finer, this growth is polynomial rather than exponential. Consequently, even for long-term scenarios, the problem dimension remains within the efficient handling capabilities of standard commercial solvers. Furthermore, unlike individual battery tracking models, the complexity of our formulation is independent of the total number of batteries, ensuring scalability in large-scale terminals.

The parameters for the numerical experiments are derived from the actual operational data of a large-scale ACT in China. Specifically, the dataset reflects the aggregate charging requirements of a fleet comprising approximately 90 AGVs. Detailed demand characteristics are presented in Figure 2 in Section 6.2. We designed two scenarios to reflect different stages of infrastructure development:

According to the actual battery life degradation curve, a single battery can undergo approximately 2500 full charge-discharge cycles. Considering that each battery is used a maximum of 50 times per month in actual terminal operations, we define one month as one decision period, thus obtaining a battery lifecycle threshold K_life of 50 usage cycles.

For cost parameters, the procurement cost of a new battery d_buy is set to 100,000 Yuan/unit. This value is based on quotes from the terminal’s actual equipment supplier and aligns with current average market prices for industrial battery packs. The single-period holding cost d_hold is 500 Yuan/month, covering maintenance, storage, and safety monitoring fees. The end-of-period asset conversion rate ρ is set to 0.3. The annual discount rate r is set to 8%, which is the general standard of China. For the DRO model parameters, the Wasserstein distance radius is set to θ = 0.1 and the confidence level to ϵ = 0.2.

To ensure the proposed model captures the operational realities of container terminals, we analyzed daily battery swapping demand data collected from a large-scale ACT in China over a two-month period (October 29 to December 26). This demand dataset serves as the empirical basis for constructing the ambiguity set F(θ).

This section conducts a granular analysis of the 5-year planning horizon for an operational BSS, simulating the investment decision process within a Chinese port enterprise’s strategic framework. The scenario targets a BSS with an initial stock of 120 batteries at varying stages of degradation, with cumulative usage cycles uniformly distributed in [Klife/2,Klife−1]. The solver reached the optimal solution with a Gap of 0.01% in 24.5 seconds.

This section extends the analysis to a 20-year planning horizon. Unlike the 5-year short-term plan, this case assumes a scenario with zero initial inventory. The objective is to determine the optimal investment roadmap from scratch, using the demand characteristics of the target terminal to validate the model’s long-term strategic performance. The solver reached the optimal solution with a gap of 0.01% in 1461 seconds.

To evaluate the effectiveness of the CVaR approximation, this paper introduces a comparative method based on the expectation with the Wasserstein distance ambiguity set.

This section analyzes the sensitivity of the investment strategy to key parameters: the Wasserstein radius θ, maximum violation probability Є, battery lifecycle K_life, and per-period budget B_t. The experimental results are illustrated in Figure 7.

The Wasserstein radius θ quantifies the level of distributional ambiguity. We varied θ within the range [0.01,0.5]. As illustrated in the results, the total NPV cost exhibits an increasing trend with θ, rising from 21.4 million Yuan to 21.7 million Yuan. A larger radius implies a higher degree of uncertainty regarding the demand distribution, forcing the robust model to procure more batteries to hedge against worst-case scenarios. Conversely, a smaller θ assumes the empirical data is highly reliable, leading to a lower cost strategy.

The parameter Є represents the terminal operator’s tolerance for service level violations. We tested Є values from 0.1 to 0.3. The results indicate an inverse relationship between Є and total cost. When Є is low (e.g., 0.1, implying a 90% service level guarantee), the model tends to build a high safety stock, significantly increasing the NPV to 24.2 million Yuan. As Є relaxes to 0.3, the system tolerates higher tail risks, allowing for a lower inventory and reduced expenditures (19.3 million Yuan).

Considering the battery technology enhancement, as K_life increases from 40 to 60 cycles, the total NPV cost decreases significantly from 26.4 million Yuan to 18.6 million Yuan. Longer life of batteries reduces the frequency of replacement cycles, thereby lowering the total procurement volume over the planning horizon (from 233 to 137 units).

Finally, we analyzed the sensitivity of the solution to the budget constraints. We found that when the budget is tightened to 80% of the baseline, the cost increases to 21.7 million Yuan. This increase is due to the liquidity constraint. Strict per-period budgets prevent the operator from purchasing in bulk when needed. Instead, the model is forced to pre-purchase batteries in earlier periods to accumulate stock gradually, which increases the NPV. As the budget relaxes to 120%, the objective function value decreases to 21.3 million Yuan, as the relaxed constraint allows for optimal, just-in-time procurement.