动力电池逆向物流网络两阶段鲁棒优化研究

doi:10.12677/mse.2025.141031

期刊菜单

动力电池逆向物流网络两阶段鲁棒优化研究
Research on Two-Stage Robust Optimization of Reverse Logistics Network for Power Batteries

DOI: 10.12677/mse.2025.141031, PDF, HTML, XML,
作者: 胡立沛, 贾永基^*：东华大学旭日工商管理学院，上海
关键词: 动力电池；逆向物流；鲁棒优化；C&CG算法；Power Battery； Reverse Logistics； Robust Optimization； C&CG Algorithm

摘要: 国内动力电池即将迎来报废高峰期，但目前废旧动力电池回收网络尚未完善，对此有必要为其建立高效的逆向物流回收网络以提高回收效率和经济效益。此外存在许多缺乏专业回收能力的非正式组织，导致网络中废旧动力电池回收量具有不确定性。因此考虑各类动力电池回收量不确定性影响，以建设成本和运输成本最小为目标建立了逆向物流两阶段鲁棒优化模型。之后采用C&CG算法进行求解，通过随机规模算例结果验证了模型和方法有效性，为建立废旧动力电池回收网络提供了参考。

Abstract: The domestic power battery sector is approaching a peak in scrap disposal; however, the current recycling network for used batteries is not yet well-established. It is essential to create an efficient reverse logistics recycling network to enhance recovery efficiency and economic benefits. Additionally, many informal organizations lack specialized recycling capabilities, leading to uncertainty in the quantity of used power batteries recovered in the network. Therefore, considering the uncertainty in the recovery volume of various types of power batteries, a two-stage robust optimization model for reverse logistics has been developed, aiming to minimize construction and transportation costs. The C&CG algorithm was subsequently employed for solution, and results from randomly scaled examples validated the model and methodology, providing a reference for establishing a recycling network for used power batteries.

文章引用：胡立沛, 贾永基. 动力电池逆向物流网络两阶段鲁棒优化研究[J]. 管理科学与工程, 2025, 14(1): 308-323. https://doi.org/10.12677/mse.2025.141031

1. 引言

虽然我国新能源汽车产销量日益增长，已经连续九年保持全球第一，但动力电池作为新能源汽车的核心零部件，使用寿命通常为5到8年，当其续航能力不能支持乘用电动汽车的日常使用时便进入报废回收阶段。因此新能源汽车产业的快速发展预示了我国未来将出现动力电池大规模退役问题。回收处理动力电池主要有资源回收和梯次利用两种方式，而我国目前市面上流通的动力电池主要为三元锂电池和磷酸铁锂电池两种类型[1]，两者因制造材料和循环充放电次数存在差异，导致其回收方式倾向性不同。如废旧三元锂电池中含有锂、钴、镍等有价贵金属且电池容量衰减速度快，直接进行资源回收利用提取有价金属等原材料相比梯次利用具有经济价值。与此相反磷酸铁锂电池由于容量衰减速度缓慢且贵金属含量低所以更适合先进行梯次利用后再进行资源回收，其可在储能系统和电动车充电站等梯次利用，以充分利用能源。此外我国废旧动力电池回收利用体系尚未完善，存在大量废旧动力电池被私人小作坊等非正式组织私自拆解倒卖情况[2]，造成了废旧动力电池回收量存在不确定性。因此根据不同废旧动力电池的回收方式不同，考虑各类电池回收数量不确定性下研究废旧动力电池回收网络设计，对研究动力电池回收网络优化布局具有重要的现实意义。

2. 文献综述

目前关于对废旧动力电池回收问题的相关研究主要集中在回收模式[3]-[5]、回收技术[6]-[8]和网络设计[9] [10]，其中不确定环境下动力电池逆向物流网络设计方面的相关研究较少，大多采用模糊规划方法求解，较少学者运用鲁棒优化方法，如Lin等[11]构建了考虑回收数量和电池质量不确定性的混合整数线性规划模型，并采用模糊规划方法求解模型得到设施选址和流量分配方案。刘娟娟等[12]采用了三角模糊数描述回收量、回收技术水平等不确定性参数，同时在考虑经济和环境效益双目标下构建了动力电池逆向物流网络的多目标模糊规划模型，之后分析总结了关键参数对目标函数值的影响程度，为管理者提供决策支持。Wang等[13]利用模糊规划解决了中国西安市退役电池回收网络的优化问题，考虑了回收电池的数量和质量不确定性，并引入三角模糊函数将不确定条件转化为确定条件。Mu等[14]从经济、环境和社会的角度开发了可持续动态逆向物流网络模型，采用三角模糊数处理不确定参数，并通过多目标线性加权求解多目标组合优化模型。

以上文献对不确定环境下废旧动力电池逆向物流回收网络研究较少，且涉及不确定性的处理大多采用模糊规划方法，其受决策者主观影响较大。而鲁棒优化方法抗风险强，是处理不确定性的前沿方法之一。再加上国内动力电池主要为三元锂电池和磷酸铁锂电池这两种动力电池为主，因此本文考虑该两类电池的回收量不确定性下构建动力电池逆向物流网络两阶段鲁棒优化模型，并采用CCG算法进行求解。

3. 问题描述与模型构建

3.1. 问题描述及假设

(1) 问题描述

本文设计了两类废旧动力电池的逆向物流网络，主要由区域收集点、收集中心、检验拆解中心、梯次利用中心、资源回收中心、废弃物处理中心、二手市场、梯次利用市场、原材料市场组成，其中梯次利用市场指储能系统和电动车充电站等领域市场。如图1所示，区域收集点A是汽车4S店和动力电池租赁网点等组成，主要负责电池更换等功能；各回收主体企业合作或自营的收集中心K负责回收区域收集点的废旧动力电池，其具备长期存储动力电池的条件；检验拆解中心I负责废旧动力电池的检验和拆解，对于两种动力电池，根据物料清单拆卸为电池模组、金属外壳、冷却剂等各零部件，同时经检验的不同电池模组将根据电池类型发往不同的回收设施：三元锂电池模组将由资源回收中心R进行资源回收利用，以提取生产动力电池的有价金属等原材料，而磷酸铁锂电池模组将由梯次利用中心N负责进行重新集成和再加工等以生产梯次产品。其他可重复使用的部件如金属外壳等共同运送至二级市场SM；剩下不可回收物品如冷却剂则由废弃物处理中心L负责对进行无害化处理；最后原材料市场RM和梯次利用市场EM分别负责接收贵金属原材料和梯次电池产品。

Figure 1. Design of reverse logistics recycling network for power batteries

图1. 动力电池逆向物流回收网络设计

当动力电池进入报废期，它们将被送到区域收集中心A进行更换或废弃；退役的各类动力电池将从分散的区域收集点A集中运往到收集中心K；之后，收集到的动力电池将从收集中心K转运到检查拆卸中心I，由检验拆解中心I对其按照清单进行拆卸并检测分类电池模组；根据电池类型，三元锂电池模组将送往资源回收中心R回收提取钴、镍、锂等原材料，并出售给原材料市场RM；而磷酸铁锂电池模组将运送至梯次利用中心N，由其生产出梯次产品并出售给梯次利用市场EM。同时，在检验拆解中心I由所有电池拆解得来的可重复使用部件将送往二手市场SM进行销售；其余不可回收物则送往废弃物处理中心L集中进行无害化处理。

该逆向物流回收网络问题需要进行决策的是：针对两类动力电池回收方式不同特点，确定收集中心K和检验拆解中心I的选址点和各设施之间的流量分配。该优化问题旨在不确定环境下最小化由建设成本和运输成本构成的总成本。同时该问题具有两阶段决策特征，第一阶段确定收集中心K和检验拆解中心I的选址和数量；第二阶段的决策确定流量分配。基于此，以最小化总成本为目标，建立废旧动力电池逆向物流网络两阶段鲁棒优化模型。

(2) 模型假设

1) 不考虑各个节点的运输方式差异，假设各节点间只有一种运输方式；

2) 所有潜在设施的设施容量都是有限的，且所有设施的潜在位置都已确定；

3) 网络中相邻物流节点的单位运输成本和距离已知；

4) 考虑“双碳”目标以及对梯次利用市场需求空间巨大，所以本文假设所有市场都有充足的需求，所有原材料和梯次利用产品都能够卖出去；

5) 各种动力电池拆解后的各部分重量比例已知。

3.2. 数学模型构建

(1) 符号说明

集合：
A	废旧动力电池区域收集点集合， $a \in A$
K	收集中心备选点集合， $k \in K$
I	检验拆解中心备选点集合， $i \in I$
N	梯次利用中心集合， $n \in N$
R	资源回收中心集合， $r \in R$
L	废弃物处理中心集合， $l \in L$
EM	梯次利用市场集合， $e \in E M$
SM	二手市场集合， $b \in S M$
RM	原材料市场集合， $m \in R M$
P	动力电池类型集合， $p \in P, p = 1, 2$
	其中1代表三元锂电池，2代表磷酸铁锂电池
参数：
$U_{k}$	收集中心k的容量上限(吨)
$U_{i}$	检验拆解中心i的容量上限(吨)
$U_{n}$	梯次利用中心n的容量上限(吨)
$U_{r}$	资源回收中心r的容量上限(吨)
$Q_{p a}$	p类型废旧动力电池在区域收集点a的回收总吨数
$C_{a k}$	从区域收集点a运送至收集中心k的单位运输成本
$C_{k i}$	从收集中心k运送至检验拆解中心i的单位运输成本
$C_{i l}$	电池中不可回收物从检验拆解中心i运送至废弃物处理中心l的单位运输成本
$C_{i r}$	从检验拆解中心i运送至资源回收中心r的单位运输成本
$C_{i n}$	从检验拆解中心i运送至梯次利用中心n的单位运输成本
$C_{i b}$	电池中可重复使用部件从检验拆解中心i运送至二手市场b的单位运输成本
$C_{r m}$	由电池得到的原材料从资源回收中心r运送至原材料市场m的单位运输成本

续表

$C_{n e}$	由电池得到的梯次产品从梯次利用中心n运送至梯次利用市场e的单位运输成本
$α_{p 1}$	p类型电池可重复利用物重量占整个动力电池重量的比例
$α_{p 2}$	p类型电池拆解后的电池模组重量占整个动力电池重量的比例
$ω_{p}$	p类型电池模组中可进行资源回收的比例
$f_{k}$	建立收集中心k的固定成本
$f_{i}$	建立检验拆解中心i的固定成本
决策变量：
$Y_{k}$	0~1变量，若建立候选点收集中心k则为1，否则为0
$Y_{i}$	0~1变量，若建立候选点检验拆解中心i则为1，否则为0
$X_{p a k}$	p类型电池从区域收集点 $a$ 运送至收集中心k的运输量(吨)
$X_{p k i}$	p类型电池从收集中心k运送至检验拆解中心i的运输量(吨)
$X_{p i l}$	p类型电池中不可回收物从检验拆解中心i运送至废弃物处理中心l的运输量(吨)
$X_{p i r}$	p类型电池从检验拆解中心i运送至资源回收中心r的运输量(吨)
$X_{p i n}$	p类型电池模组从检验拆解中心i运送至梯次利用中心n的运输量(吨)
$X_{p i b}$	p类型电池中可重复使用部件从检验拆解中心i运送至二手市场b的运输量(吨)
$X_{p r m}$	由p类型电池得到的原材料从回收中心r运送至原材料市场m的运输量(吨)
$X_{p n e}$	由p类型电池得到的梯次产品从梯次利用中心n运送至梯次利用市场e的运输量(吨)

(2) 动力电池逆向物流网络两阶段鲁棒优化模型

$\begin{matrix} \min Z = \sum_{k \in K} f_{k} Y_{k} + \sum_{i \in I} f_{i} Y_{i} + \max_{Q \in ℒ} \min_{F (Y, Q)} (\sum_{p \in P} \sum_{a \in A} \sum_{k \in K} C_{a k} X_{p a k} + \sum_{p \in P} \sum_{k \in K} \sum_{i \in I} C_{k i} X_{p k i} \\ + \sum_{p \in P} \sum_{i \in I} \sum_{l \in L} C_{i l} X_{p i l} + \sum_{p \in P} \sum_{i \in I} \sum_{b \in S M} C_{i b} X_{p i b} + \sum_{p \in P} \sum_{i \in I} \sum_{r \in R} C_{i r} X_{p i r} \\ + \sum_{p \in P} \sum_{r \in R} \sum_{m \in R M} C_{r m} X_{p r m} + \sum_{p \in P} \sum_{i \in I} \sum_{n \in N} C_{i n} X_{p i n} + \sum_{p \in P} \sum_{n \in N} \sum_{e \in E M} C_{n e} X_{p n e}) \end{matrix}$ (1)

s.t.

$\sum_{k \in K} X_{p a k} \geq {\tilde{Q}}_{p a} \forall a \in A, p \in P$ (2)

$\sum_{a \in A} X_{p a k} = \sum_{i \in I} X_{p k i} \forall k \in K, p \in P$ (3)

$\sum_{k \in K} X_{p k i} = \sum_{l \in L} X_{p i l} + \sum_{b \in S M} X_{p i b} + \sum_{r \in R} X_{p i r} + \sum_{n \in N} X_{p i n} \forall i \in I, p \in P$ (4)

$\sum_{b \in S M} X_{p i b} = α_{p 1} \sum_{k \in K} X_{p k i} \forall i \in I, p \in P$ (5)

$\sum_{r \in R} X_{p i r} = ω_{p} α_{p 2} \sum_{k \in K} X_{p k i} \forall i \in I, p \in P$ (6)

$\sum_{n \in N} X_{p i n} = (1 - ω_{p}) α_{p 2} \sum_{k \in K} X_{p k i} \forall i \in I, p \in P$ (7)

$\sum_{l \in L} X_{p i l} = (1 - α_{p 1} - α_{p 2}) \sum_{k \in K} X_{p k i} \forall i \in I, p \in P$ (8)

$\sum_{i \in I} X_{p i r} = \sum_{m \in R M} X_{p r m} \forall r \in R, p \in P$ (9)

$\sum_{i \in I} X_{p i n} = \sum_{e \in E M} X_{p n e} \forall n \in N, p \in P$ (10)

$\sum_{p \in P} \sum_{a \in A} X_{p a k} \leq U_{k} Y_{k} \forall k \in K$ (11)

$\sum_{p \in P} \sum_{k \in K} X_{p k i} \leq U_{i} Y_{i} \forall i \in I$ (12)

$\sum_{p \in P} \sum_{i \in I} X_{p i r} \leq U_{r} \forall r \in R$ (13)

$\sum_{p \in P} \sum_{i \in I} X_{p i n} \leq U_{n} \forall n \in N$ (14)

$Y_{k}, Y_{i} \in {0, 1} \forall k \in K, \forall i \in I$ (15)

$\begin{array}{l} X_{p a k}, X_{p k i}, X_{p i l}, X_{p i b}, X_{p i r}, X_{p i n}, X_{p r m}, X_{p n e} \geq 0 \\ \forall p \in P, a \in A, k \in K, i \in I, l \in L, b \in S M, r \in R, n \in N, m \in R M, e \in E M \end{array}$ (16)

其中目标函数式(1)表示选址建设成本和最坏情况下的运输成本，包括第一阶段收集中心和检验拆解中心的选址建设成本和第二阶段的各设施节点之间运输成本。第二阶段中max表示不确定集合中的最坏场景，min表示在给定场景和第一阶段决策变量下的最小化运输成本。约束(2)表示每个区域收集点回收的废旧动力电池总吨数等于其运往收集中心的所有运输量之和；约束(3)表示每个收集中点接收的废旧动力电池总吨数等于其运往检验拆解中心的所有运输量之和；约束(4)表示每个检验拆解中心接收的废旧动力电池总吨数等于其运往废弃物处理中心、资源回收中心、梯次利用中心以及二手市场的所有运输量之和；约束(5)表示每个检验拆解中心拆解完废旧动力电池后得到的可重复使用部件运输量等于其运往二手市场的所有运输量之和；约束(6)表示每个检验拆解中心拆解废旧三元锂电池后得到的电池模组运输量等于其运往资源回收中心的所有运输量之和；约束(7)表示每个检验拆解中心拆解废旧磷酸铁锂电池后得到的电池模组运输量等于其运往梯次利用中心的所有运输量之和；约束(8)表示每个检验拆解中心拆解完废旧动力电池后得到的不可回收物运输量等于其运往废弃物处理中心的所有运输量之和；约束(9)表示每个资源回收中心处理完三元电池模组后得到的原材料运输量等于其运往原材料市场的所有运输量之和；约束(10)表示每个梯次利用中心处理完磷酸铁锂电池模组后得到的梯次产品运输量等于其运往梯次利用市场的所有运输量之和；约束(11)~(14)为各设施容量限制；约束(15)~(16)为决策变量约束，前者表示0~1变量约束，后者表示流量非负约束。

4. 列和约束生成算法C&CG

4.1. 不确定回收量鲁棒优化处理

为了刻画回收量不确定，设废旧动力电池p在区域收集点a的回收量 ${\tilde{Q}}_{p a}$ 的取值范围为 $[{\bar{Q}}_{p a} - {\hat{Q}}_{p a}, {\bar{Q}}_{p a} + {\hat{Q}}_{p a}]$ ，其中 ${\bar{Q}}_{p a}$ 表示回收量的名义值， ${\hat{Q}}_{p a} \geq 0$ 表示回收量的最大偏离程度。实际上，当实际回收量大于名义回收量时，该问题决策才会出现解的恶化，因此该模型实际只考虑 ${\tilde{Q}}_{p a}$ 取值范围为 $[{\bar{Q}}_{p a}, {\bar{Q}}_{p a} + {\hat{Q}}_{p a}]$ 。引入预算控制参数 $Γ_{p a}$ 来表示模型的保守程度，即回收量 $Q_{p a}$ 出现数值偏离的数量，当 $Γ_{p a}$ 取0值则表示模型不考虑某类动力电池回收量变化的影响，此时模型的抗风险性较差；当 $Γ_{p a}$ 取最大值时则考虑某类动力电池回收量全部发生偏离，此时模型最保守。因此决策者可根据实际情况通过调整 $Γ_{p a} \in [0, | A |]$ 的取值，以控制模型的鲁棒性。另外在模型实际应用中，所有参数都存在不确定性太过于保守，所以本文认为至多有 $Γ_{p a}$ 个某类动力电池在区域收集点的回收量发生变化，同时还有小数部分代表某类动力电池在区域收集点的回收量改变 $(Γ_{p a} - Γ_{p a}) {\hat{Q}}_{p a}$ 。因此刻画不确定性集合为 $ℒ = {{\tilde{Q}}_{p a} : {\tilde{Q}}_{p a} = {\bar{Q}}_{p a} + h_{p a} {\hat{Q}}_{p a}, h_{p a} \in [0, 1], \sum_{a \in A} h_{p a} \leq Γ_{p a}, a = 1, 2, \dots, A, \forall p \in P}$ 。

4.2. 列和约束生成算法

两阶段鲁棒优化模型中第二阶段是存在着max-min结构，涉及变量和约束较多，难以通过现有商业求解器高效求解。因此为解决计算时间过长和占用运行内存过多问题，本文采用由Zeng [15]提出了列和约束生成算法(C&CG算法)来求解两阶段鲁棒优化模型，并运用KKT条件转化第二阶段问题。

根据C&CG算法流程对两阶段鲁棒优化模型进行重构，引入表示备选点收集中心和检验拆解中心容量的连续变量 $Z_{k}$ 、 $Z_{i}$ ，于是原模型的主问题MP为：

$\min Z = \sum_{k \in K} f_{k} Y_{k} + \sum_{i \in I} f_{i} Y_{i} + t$ (17)

s.t.

$Z_{k} \leq U_{k} Y_{k} \forall k \in K$ (18)

$Z_{i} \leq U_{i} Y_{i} \forall i \in I$ (19)

$Z_{k} \geq \sum_{p \in P} \sum_{a \in A} X_{p a k} (o) \forall k \in K, o \in O$ (20)

$Z_{i} \geq \sum_{p \in P} \sum_{k \in K} X_{p k i} (o) \forall i \in I, o \in O$ (21)

$\begin{array}{l} \sum_{p \in P} \sum_{a \in A} \sum_{k \in K} C_{a k} X_{p a k} (o) + \sum_{p \in P} \sum_{k \in K} \sum_{i \in I} C_{k i} X_{p k i} (o) + \sum_{p \in P} \sum_{i \in I} \sum_{l \in L} C_{i l} X_{p i l} (o) \\ + \sum_{p \in P} \sum_{i \in I} \sum_{b \in S M} C_{i b} X_{p i b} (o) + \sum_{p \in P} \sum_{i \in I} \sum_{r \in R} C_{i r} X_{p i r} (o) + \sum_{p \in P} \sum_{r \in R} \sum_{m \in R M} C_{r m} X_{p r m} (o) \\ + \sum_{p \in P} \sum_{i \in I} \sum_{n \in N} C_{i n} X_{p i n} (o) + \sum_{p \in P} \sum_{n \in N} \sum_{e \in E M} C_{n e} X_{p n e} (o) \leq t \forall o \in O \end{array}$ (22)

$\sum_{k \in K} X_{p a k} (o) \geq Q_{p a} (o) \forall a \in A, p \in P$ (23)

$\sum_{a \in A} X_{p a k} (o) = \sum_{i \in I} X_{p k i} (o) \forall k \in K, p \in P$ (24)

$\sum_{k \in K} X_{p k i} (o) = \sum_{l \in L} X_{p i l} (o) + \sum_{b \in S M} X_{p i b} (o) + \sum_{r \in R} X_{p i r} (o) + \sum_{n \in N} X_{p i n} (o) \forall i \in I, p \in P, o \in O$ (25)

$\sum_{b \in S M} X_{p i b} (o) = α_{p 1} \sum_{k \in K} X_{p k i} (o) \forall i \in I, p \in P, o \in O$ (26)

$\sum_{r \in R} X_{p i r} (o) = ω_{p} α_{p 2} \sum_{k \in K} X_{p k i} (o) \forall i \in I, p \in P, o \in O$ (27)

$\sum_{n \in N} X_{p i n} (o) = (1 - ω_{p}) α_{p 2} \sum_{k \in K} X_{p k i} (o) \forall i \in I, p \in P, o \in O$ (28)

$\sum_{l \in L} X_{p i l} (o) = (1 - α_{p 1} - α_{p 2}) \sum_{k \in K} X_{p k i} (o) \forall i \in I, p \in P, o \in O$ (29)

$\sum_{i \in I} X_{p i r} (o) = \sum_{m \in R M} X_{p r m} (o) \forall r \in R, p \in P, o \in O$ (30)

$\sum_{i \in I} X_{p i n} (o) = \sum_{e \in E M} X_{p n e} (o) \forall n \in N, p \in P, o \in O$ (31)

$\sum_{p \in P} \sum_{i \in I} X_{p i r} (o) \leq U_{r} \forall r \in R, o \in O$ (32)

$\sum_{p \in P} \sum_{i \in I} X_{p i n} (o) \leq U_{n} \forall n \in N, o \in O$ (33)

$Y_{k}, Y_{i} \in {0, 1} \forall k \in K, \forall i \in I$ (34)

$\begin{array}{l} X_{p a k} (o), X_{p k i} (o), X_{p i l} (o), X_{p i b} (o), X_{p i r} (o), X_{p i n} (o), X_{p r m} (o), X_{p n e} (o) \geq 0 \\ \forall p \in P, a \in A, k \in K, i \in I, l \in L, b \in S M, r \in R, n \in N, m \in R M, e \in E M \end{array}$ (35)

其中，O表示场景集合， $o \in O$ ；式(18)至式(21)表示收集中心和检验拆解中心的容量限制，其他约束的释义和上节相同；同时式(22)是一种最优割平面，添加该割平面目的旨在算法迭代过程中能够提高目标函数值的下界，其中辅助变量t表示的是第二阶段决策问题的最优目标函数值，而子问题SP数学模型为：

$\begin{array}{l} \max_{Q \in ℒ} \min_{F (Y, Q)} \sum_{p \in P} \sum_{a \in A} \sum_{k \in K} C_{a k} X_{p a k} + \sum_{p \in P} \sum_{k \in K} \sum_{i \in I} C_{k i} X_{p k i} + \sum_{p \in P} \sum_{i \in I} \sum_{l \in L} C_{i l} X_{p i l} \\ + \sum_{p \in P} \sum_{i \in I} \sum_{b \in S M} C_{i b} X_{p i b} + \sum_{p \in P} \sum_{i \in I} \sum_{r \in R} C_{i r} X_{p i r} + \sum_{p \in P} \sum_{r \in R} \sum_{m \in R M} C_{r m} X_{p r m} \\ + \sum_{p \in P} \sum_{i \in I} \sum_{n \in N} C_{i n} X_{p i n} + \sum_{p \in P} \sum_{n \in N} \sum_{e \in E M} C_{n e} X_{p n e} \end{array}$ (36)

s.t.

式(2)至式(10)

$\sum_{p \in P} \sum_{a \in A} X_{p a k} \leq Z_{k}^{*} \forall k \in K$ (37)

$\sum_{p \in P} \sum_{k \in K} X_{p k i} \leq Z_{i}^{*} \forall i \in I$ (38)

式(13)至式(16)

式(37)和式(38)的 $Z_{k}^{*}$ 、 $Z_{i}^{*}$ 为主问题求解得到的最优解取值，将其传入子问题进行求解以获得最坏的场景，以及对应该场景的新变量和新约束条件。同时注意到子问题SP是一个max-min双层优化问题，难以直接求解。因此本文将采用KKT条件将内层min问题做等价转换，进而将子问题转为易求解的单层优化问题。

首先令 $β_{p a}$ ， $β_{p k}$ ， $θ_{p i}^{1}$ ， $θ_{p i}^{2}$ ， $μ_{p i}^{1}$ ， $μ_{p i}^{2}$ ， $μ_{p i}^{3}$ ， $δ_{p r}$ ， $ε_{p n}$ ， $τ_{k}$ ， $λ_{i}$ ， $ν_{r}$ ， $π_{n}$ 为内层min问题中约束条件式(2)~(10)、式(37)~(38)、式(13)~式(14)的对偶变量，根据对偶理论得到内层问题的对偶问题模型如下：

$\max \sum_{p \in P} \sum_{a \in A} {\tilde{Q}}_{p a} β_{p a} - \sum_{k \in K} U_{k} Y_{k} τ_{k} - \sum_{i \in I} U_{i} Y_{i} λ_{i} - \sum_{r \in R} U_{r} ν_{r} - \sum_{n \in N} U_{n} π_{n}$ (39)

s.t.

$β_{p a} - β_{p k} - τ_{k} \leq C_{a k}, \forall a \in A, k \in K$ (40)

$\begin{array}{l} β_{p a} - θ_{p i}^{1} + α_{p 1} θ_{p i}^{2} + α_{p 2} ω_{p} μ_{p i}^{1} + α_{p 2} (1 - ω_{p}) μ_{p i}^{2} \\ + μ_{p i}^{3} (1 - α_{p 1} - α_{p 2}) - λ_{i} \leq C_{k i}, \forall p \in P, k \in K, i \in I \end{array}$ (41)

$θ_{p i}^{1} - μ_{p i}^{3} \leq C_{i l}, \forall p \in P, i \in I, l \in L$ (42)

$θ_{p i}^{1} - θ_{p i}^{2} \leq C_{i b}, \forall p \in P, i \in I, b \in S M$ (43)

$θ_{p i}^{1} - μ_{p i}^{1} - δ_{p r} - ν_{r} \leq C_{i r}, \forall p \in P, i \in I, r \in R$ (44)

$δ_{p r} \leq C_{r m}, \forall p \in P, r \in R, m \in R M$ (45)

$θ_{p i}^{1} - μ_{p i}^{2} - ε_{p n} - π_{n} \leq C_{i n}, \forall p \in P, i \in I, n \in N$ (46)

$ε_{p n} \leq C_{n e}, \forall p \in P, n \in N, e \in E M$ (47)

$β_{p a}, τ_{k}, λ_{i}, ν_{r}, π_{n} \geq 0 \forall p \in P, k \in K, i \in I, r \in R, n \in N$ (48)

$β_{p k}, θ_{p i}^{1}, θ_{p i}^{2}, μ_{p i}^{1}, μ_{p i}^{2}, μ_{p i}^{3}, δ_{p r}, ε_{p n} 无约束 \forall p \in P, a \in A, k \in K, i \in I, r \in R, n \in N$ (49)

其次根据互补松弛定理得到以下互补松弛约束：

$(C_{a k} - β_{p a} + β_{p k} + τ_{k}) X_{p a k} = 0, \forall p \in P, a \in A, k \in K$ (50)

$\begin{array}{l} (C_{k i} - β_{p a} + θ_{p i}^{1} - α_{p 1} θ_{p i}^{2} - α_{p 2} ω_{p} μ_{p i}^{1} - α_{p 2} (1 - ω_{p}) μ_{p i}^{2} \\ - μ_{p i}^{3} (1 - α_{p 1} - α_{p 2}) + λ_{i}) X_{p k i} = 0, \forall p \in P, k \in K, i \in I \end{array}$ (51)

$(C_{i l} - θ_{p i}^{1} + μ_{p i}^{3}) X_{p i l} = 0, \forall p \in P, i \in I, l \in L$ (52)

$(C_{i b} - θ_{p i}^{1} + θ_{p i}^{2}) X_{p i b} = 0, \forall p \in P, i \in I, b \in S M$ (53)

$(C_{i r} - θ_{p i}^{1} + μ_{p i}^{1} + δ_{p r} + ν_{r}) X_{p i r} = 0, \forall p \in P, i \in I, r \in R$ (54)

$(C_{r m} - δ_{p r}) X_{p r m} = 0, \forall p \in P, r \in R, m \in R M$ (55)

$(C_{i n} - θ_{p i}^{1} + μ_{p i}^{2} + ε_{p n} + π_{n}) X_{p i n} = 0, \forall p \in P, i \in I, n \in N$ (56)

$(C_{n e} - ε_{p n}) X_{p n e} = 0, \forall p \in P, n \in N, e \in E M$ (57)

$(Z_{k}^{*} - \sum_{p \in P} \sum_{a \in A} X_{p a k}) τ_{k} = 0 \forall k \in K$ (58)

$(Z_{i}^{*} - \sum_{p \in P} \sum_{k \in K} X_{p k i}) λ_{i} = 0 \forall i \in I$ (59)

$(U_{r} - \sum_{p \in P} \sum_{i \in I} X_{p i r}) ν_{r} = 0 \forall r \in R$ (60)

$(U_{n} - \sum_{p \in P} \sum_{i \in I} X_{p i n}) π_{n} = 0 \forall n \in N$ (61)

$(\sum_{k \in K} X_{p a k} - {\tilde{Q}}_{p a}) β_{p a} = 0 \forall a \in A, p \in P$ (62)

注意到以上互补松弛约束均为非线性化的，为求解模型需将其线性化，以式(50)为例，其他互补松弛约束同理，引入0~1变量 $φ_{p a k}$ ，式(50)等价于以下形式，其中 $ℳ$ 表示一个很大的常数：

$X_{p a k} \leq ℳ \cdot φ_{p a k}, \forall p \in P, a \in A, k \in K$ (63)

$C_{a k} - β_{p a} + β_{p k} + τ_{k} \leq ℳ \cdot (1 - φ_{p a k}), \forall p \in P, a \in A, k \in K$ (64)

对此采用KKT条件，将子问题SP的内层min问题转化为了“原问题约束”“对偶问题约束”以及“互补松弛约束”，结合外层max问题得到的子问题SP-CCG与SP等价。子问题SP-CCG数学模型如下：

$\begin{array}{l} \max \sum_{p \in P} \sum_{a \in A} \sum_{k \in K} C_{a k} X_{p a k} + \sum_{p \in P} \sum_{k \in K} \sum_{i \in I} C_{k i} X_{p k i} + \sum_{p \in P} \sum_{i \in I} \sum_{l \in L} C_{i l} X_{p i l} \\ + \sum_{p \in P} \sum_{i \in I} \sum_{b \in S M} C_{i b} X_{p i b} + \sum_{p \in P} \sum_{i \in I} \sum_{r \in R} C_{i r} X_{p i r} + \sum_{p \in P} \sum_{r \in R} \sum_{m \in R M} C_{r m} X_{p r m} \\ + \sum_{p \in P} \sum_{i \in I} \sum_{n \in N} C_{i n} X_{p i n} + \sum_{p \in P} \sum_{n \in N} \sum_{e \in E M} C_{n e} X_{p n e} \end{array}$ (65)

s.t.

式(2)至式(10)

式(37) (38)

式(13)至式(16)

式(40)至式(49)

式(63) (64)

$X_{p k i} \leq ℳ \cdot φ_{p k i}, \forall p \in P, k \in K, i \in I$ (66)

$\begin{array}{l} C_{k i} - β_{p a} + θ_{p i}^{1} - α_{p 1} θ_{p i}^{2} - α_{p 2} ω_{p} μ_{p i}^{1} - α_{p 2} (1 - ω_{p}) μ_{p i}^{2} \\ - μ_{p i}^{3} (1 - α_{p 1} - α_{p 2}) + λ_{i} \leq ℳ \cdot (1 - φ_{p k i}), \forall p \in P, k \in K, i \in I \end{array}$ (67)

$X_{p i l} \leq ℳ \cdot φ_{p i l}, \forall p \in P, i \in I, l \in L$ (68)

$C_{i l} - θ_{p i}^{1} + μ_{p i}^{3} \leq ℳ \cdot (1 - φ_{p i l}), \forall p \in P, i \in I, l \in L$ (69)

$X_{p i b} \leq ℳ \cdot φ_{p i b}, \forall p \in P, i \in I, b \in S M$ (70)

$C_{i b} - θ_{p i}^{1} + θ_{p i}^{2} \leq ℳ \cdot (1 - φ_{p i b}), \forall p \in P, i \in I, b \in S M$ (71)

$X_{p i r} \leq ℳ \cdot φ_{p i r}, \forall p \in P, i \in I, r \in R$ (72)

$C_{i r} - θ_{p i}^{1} + μ_{p i}^{1} + δ_{p r} + ν_{r} \leq ℳ \cdot (1 - φ_{p i r}), \forall p \in P, i \in I, r \in R$ (73)

$X_{p r m} \leq ℳ \cdot φ_{p r m}, \forall p \in P, r \in R, m \in R M$ (74)

$C_{r m} - δ_{p r} \leq ℳ \cdot (1 - φ_{p r m}), \forall p \in P, r \in R, m \in R M$ (75)

$X_{p i n} \leq ℳ \cdot φ_{p i n}, \forall p \in P, i \in I, n \in N$ (76)

$C_{i n} - θ_{p i}^{1} + μ_{p i}^{2} + ε_{p n} + π_{n} \leq ℳ \cdot (1 - φ_{p i n}), \forall p \in P, i \in I, n \in N$ (77)

$X_{p n e} \leq ℳ \cdot φ_{p n e}, \forall p \in P, n \in N, e \in E M$ (78)

$C_{n e} - ε_{p n} \leq ℳ \cdot (1 - φ_{p n e}), \forall p \in P, n \in N, e \in E M$ (79)

$τ_{k} \leq ℳ \cdot φ_{k} \forall k \in K$ (80)

$Z_{k}^{*} - \sum_{p \in P} \sum_{a \in A} X_{p a k} \leq ℳ \cdot (1 - φ_{k}) \forall k \in K$ (81)

$λ_{i} \leq ℳ \cdot φ_{i} \forall i \in I$ (82)

$Z_{i}^{*} - \sum_{p \in P} \sum_{k \in K} X_{p k i} \leq ℳ \cdot (1 - φ_{i}) \forall i \in I$ (83)

$ν_{r} \leq ℳ \cdot φ_{r} \forall r \in R$ (84)

$U_{r} - \sum_{p \in P} \sum_{i \in I} X_{p i r} \leq ℳ \cdot (1 - φ_{r}) \forall r \in R$ (85)

$π_{n} \leq ℳ \cdot φ_{n} \forall n \in N$ (86)

$U_{n} - \sum_{p \in P} \sum_{i \in I} X_{p i n} \leq ℳ \cdot (1 - φ_{r}) \forall n \in N$ (87)

$β_{p a} \leq ℳ \cdot φ_{p a} \forall a \in A, p \in P$ (88)

$\sum_{k \in K} X_{p a k} - {\tilde{Q}}_{p a} \leq ℳ \cdot (1 - φ_{p a}) \forall a \in A, p \in P$ (89)

其中式(2)至(10)、式(37) (38)以及式(13)至式(16)为原问题约束条件，式(40)至式(49)为对偶问题约束条件，式(63)至式(89)为互补松弛约束。通过求解子问题SP-CCG得到最坏场景 $o^{*}$ ，并在主问题中加入该最坏场景对应的列和约束条件，以进行下一步迭代求解，直至满足收敛条件。C&CG算法框架如图2所示：

Figure 2. C&CG algorithm framework diagram

图2. C&CG算法框架图

5. 算例分析与数值实验

5.1. 参数设置

由于废旧动力电池回收产业处于发展阶段，缺乏详细分类别的报废动力电池数据，因此本文建立了6个不同规模算例，报废电池回收总吨数在区间[90~120]服从正态分布随机生成，同样在400 × 400 km区域里随机生成各设施点坐标。设置三元锂电池和磷酸铁锂电池回收量的最大偏离值 ${\hat{Q}}_{p a}$ 分别为40和20吨，同时鲁棒控制参数变量 $Γ_{p a}$ 均为 $| A | \times 0.7$ 。其他参数设置和各规模算例情况分别如表1和表2所示。

Table 1. Parameter settings

表1. 参数设置

参数	设置
$U_{k}$	1500
$U_{i}$	2500
$U_{n}$	1500
$U_{r}$	1500
$f_{k}$	650,000
$f_{i}$	2,050,000
$α_{p 1}$	三元锂电池：0.3 磷酸铁锂电池：0.35
$α_{p 2}$	三元锂电池：0.6 磷酸铁锂电池：0.55
$ω_{p}$	三元锂电池：1 磷酸铁锂电池：0
$C_{i j}$	$d_{i j} \times t$ ， $d_{i j}$ 为设施之间的距离， $t = 0.4$ 为单位运输成本

Table 2. Configuration table for different scale examples

表2. 不同规模算例设置表

算例	A	K	I	N	R	L	SM	RM	EM
Intance-1	8	4	2	2	2	2	2	2	2
Intance-2	12	5	2	2	2	2	2	2	2
Intance-3	15	6	3	3	3	2	2	2	2
Intance-4	18	7	3	3	3	2	3	3	3
Intance-5	21	8	5	3	3	2	3	3	3
Intance-6	25	10	5	3	3	2	3	3	3

5.2. 算例结果

算例测试平台为Intel Core i9-13900HX 2.20GHz处理器，内存16 GB的笔记本计算机，编程语言为python3.9，商业求解器为Gurobi11.0。各算例结果如表3所示，随着问题的复杂性增加，回收网络的总成本也随之上升。尽管问题规模扩大导致求解时间增加，各算例仍能在有限的迭代次数内得出最终结果，并且Gap均趋近于0。这表明C&CG算法适用于求解本模型，求解效率且效果良好。

Table 3. Results of C&CG algorithm examples

表3. C&CG算法算例结果

算例	目标值下界	目标值上界	Gap	迭代次数	求解时间
Intance-1	3729872.5	3729872.5	0.00%	1	0.02 s
Intance-2	5979890.3	5979890.3	0.00%	2	1.30 s
Intance-3	6640842.5	6640842.5	0.00%	1	3.32 s
Intance-4	7395375.7	7395375.7	0.00%	1	11.00 s
Intance-5	9527450.4	9527450.4	0.00%	2	45.66 s
Intance-6	10264675.7	10264675.7	0.00%	3	170.00 s

5.3. 决策结果分析

以算例6为例，不确定性场景下备选点选址、回收需求分配结果及其总成本见表4，可知，收集中心和检验拆解中心的建设数量与设施容量上限、识别出的最坏场景下两类废旧动力电池回收数量以及对应回收方式有关，5个收集中心和3个检验拆解中心即可满足该算例在最坏场景下全部废旧动力电池回收需求，且每个中心负责存放的各类动力电池总数量因动力电池回收方式有所差异，在收集中心和检验拆解中心的实际容量接近于总容量上限时，将会充分考虑各类动力电池回收数量实际占比以使得总成本最小化。其中，以设施建设成本和运输成本之和最小为优化目标，求解得到收集中心K5因为其与检验拆解中心I1更靠近资源回收中心R3，所以收集中心K5在向区域收集点回收废旧动力电池时会更倾向于回收三元锂电池，以缓解其他收集中心回收三元锂电池造成的运输成本上涨情况。从图3和图4可以得出收集中心K5的实际设施容量已达上限，且三元锂电池回收量占比较高，超过六成。同理，收集中心K3因由检验拆解中心I1 (靠近资源回收中心R3)服务，所以其回收三元锂电池回收量占比较高。为了更好观察决策结果的合理性，对所有收集点和经鲁棒模型求解后的设施具体位置和收集中心服务范围进行可视化处理，如图5所示。可以观察到A3、A4等区域收集点(图中黄色标记)同时由两个收集中心所共同服务，在回收废旧动力电池时会根据与梯次利用中心或资源利用中心距离远近倾向回收某类电池，如刚才所述收集中心K5对共同服务的收集点倾向运输三元锂电池，而收集中心K10和其附近的检验拆解中心I4因靠近梯次利用中心N1进而回收更多磷酸铁锂电池，与收集中心K5形成互补关系，以尽可能减少因电池回收方式不同造成的运输成本上涨情况。

Table 4. Example 6 decision results under uncertainty conditions

表4. 算例6不确定性下决策结果

检验拆解中心编号	收集中心编号	收集中心实际容量(吨)	区域收集点编号	总成本(元)
1	3	1000	4, 9, 16, 19, 23	10264675.7
1	5	1500	1, 3, 7, 12, 14, 15, 18, 20, 22
4	10	1500	3, 10, 12, 13, 14, 22, 24, 25
5	7	849	2, 5, 9, 10, 11, 17, 23
5	8	1500	4, 6, 7, 8, 18, 21

Figure 3. Capacity structure distribution diagram of collection center K under uncertainty

图3. 不确定性下收集中心K容量结构分布图

Figure 4. Capacity structure distribution diagram of decoupling center I under uncertainty

图4. 不确定性下检验拆解中心I容量结构分布图

Figure 5. Distribution map of service areas for collection points and facilities under uncertainty

图5. 不确定性下收集点与设施服务范围分布图

Figure 6. Distribution map of service areas for collection points and facilities under deterministic conditions

图6. 确定性下收集点与设施服务范围分布图

Table 5. Example 6 decision results under deterministic conditions

表5. 算例6确定性下决策结果

检验拆解中心编号	收集中心编号	收集中心实际容量(吨)	区域收集点编号	总成本(元)
1	3	1000	4, 8, 9, 16, 19, 21, 23	9449996.6
1	5	1500	1, 3, 6, 7, 12, 14, 15, 18, 20
4	10	1500	3, 7, 10, 12, 13, 14, 17, 22, 24, 25
5	7	1299	2, 5, 6, 8, 9, 11, 17, 21, 23

Figure 7. Capacity structure distribution diagram of collection center K under deterministic conditions

图7. 确定性下收集中心K容量结构分布图

Figure 8. Capacity structure distribution diagram of decoupling center I under deterministic conditions

图8. 确定性下检验拆解中心I容量结构分布图

此外，与该算例确定性情况下的求解决策结果相比，结合图6~8和表5所示，可知鲁棒优化模型识别出了“重要”的最坏场景，使其最优解在所有不确定场景中具有较强的抵抗不确定风险的能力，在面对各收集点不同类型回收需求上涨的潜在风险时，结合不同类型电池回收流向差异和设施容量上限考虑增设收集中心，以协调其他收集中心各类型回收需求。所以，增设了收集中心K8是因为收集中心K7已达到容量上限，需要另一个收集中心填补收集中心K3和K5因回收更多三元锂电池而造成磷酸铁锂电池回收缺口需求，以减少运输成本。

6. 结论

基于国内动力电池即将迎来报废高峰期的背景下，针对动力电池逆向物流网络设计问题，本文将两类动力电池回收量作为不确定性参数，以建设成本和运输成本最小化为目标，建立了两阶段鲁棒优化模型，同时引入预算不确定集处理回收量不确定性。之后采用C&CG算法重构原模型，并运用KKT条件将双层规划子问题转化成易求解的单层规划问题。随后通过随机生成的不同规模算例进行数值实验，最后验证了模型和算法的有效性。同时决策者应根据不同类型电池的回收数量和方式，确定回收设施的最终选址，以降低运输和建设成本。最后，本文模型设计和求解算法上仍存在一些缺陷，实际情况更为复杂，对应的算法复杂性更高。对此为使研究更加贴合实际和提高求解效率，在未来研究中可以考虑设施容量扩容和多周期库存管理问题，同时随着问题规模增大，探索和改进现有算法，找到合适有效的算法以提高模型求解效率也是未来研究的重要方向。

NOTES

^*通讯作者。

参考文献

[1]	黎华玲, 丁志英, 梁昌铖, 等. 废旧磷酸铁锂动力电池正极材料回收再生方法分析[J]. 科技资讯, 2021, 19(29): 69-71.
[2]	高艳红, 黎振东. 考虑非正式组织的动力电池生产者责任延伸回收模式研究[J]. 管理评论, 2024, 36(3): 225-236.
[3]	郝硕硕, 董庆银, 李金惠. 基于成本核算的废旧动力电池回收模式分析与趋势研究[J]. 中国环境科学, 2021, 41(10): 4745-4755.
[4]	姚海琳, 王昶, 黄健柏. EPR下我国新能源汽车动力电池回收利用模式研究[J]. 科技管理研究, 2015, 35(18): 84-89.
[5]	徐亦扬, 马宁, 冯宇, 等. 基于EPR制度的废旧动力电池问题探析[J]. 中国储运, 2023(10): 199-200.
[6]	王琨强, 周玲, 彭琦, 等. 退役锂离子动力电池梯次利用回收技术综述[J]. 电子质量, 2023(6): 109-113.
[7]	Wang, T., Liu, T. and Sun, H. (2023) Direct Recycling for Advancing Sustainable Battery Solutions. Materials Today Energy, 38, Article 101434. https://doi.org/10.1016/j.mtener.2023.101434
[8]	Tao, Y., Wang, Z., Wu, B., Tang, Y. and Evans, S. (2023) Environmental Life Cycle Assessment of Recycling Technologies for Ternary Lithium-Ion Batteries. Journal of Cleaner Production, 389, Article 136008. https://doi.org/10.1016/j.jclepro.2023.136008
[9]	Rosenberg, S., Glöser-Chahoud, S., Huster, S. and Schultmann, F. (2023) A Dynamic Network Design Model with Capacity Expansions for EOL Traction Battery Recycling—A Case Study of an OEM in Germany. Waste Management, 160, 12-22. https://doi.org/10.1016/j.wasman.2023.01.029
[10]	Wang, L., Wang, X. and Yang, W. (2020) Optimal Design of Electric Vehicle Battery Recycling Network—From the Perspective of Electric Vehicle Manufacturers. Applied Energy, 275, Article 115328. https://doi.org/10.1016/j.apenergy.2020.115328
[11]	Lin, J., Li, X., Zhao, Y., Chen, W. and Wang, M. (2023) Design a Reverse Logistics Network for End-of-Life Power Batteries: A Case Study of Chengdu in China. Sustainable Cities and Society, 98, Article 104807. https://doi.org/10.1016/j.scs.2023.104807
[12]	刘娟娟, 郭炎可. 考虑不确定性的电动汽车动力电池逆向物流网络设计[J]. 上海海事大学学报, 2021, 42(2): 96-102.
[13]	Wang, C., Feng, X., Woo, S., Wood, J. and Yu, S. (2023) The Optimization of an EV Decommissioned Battery Recycling Network: A Third-Party Approach. Journal of Environmental Management, 348, Article 119299. https://doi.org/10.1016/j.jenvman.2023.119299
[14]	Mu, N., Wang, Y., Chen, Z., Xin, P., Deveci, M. and Pedrycz, W. (2023) Multi-Objective Combinatorial Optimization Analysis of the Recycling of Retired New Energy Electric Vehicle Power Batteries in a Sustainable Dynamic Reverse Logistics Network. Environmental Science and Pollution Research, 30, 47580-47601. https://doi.org/10.1007/s11356-023-25573-w
[15]	Zeng, B. and Zhao, L. (2013) Solving Two-Stage Robust Optimization Problems Using a Column-and-Constraint Generation Method. Operations Research Letters, 41, 457-461. https://doi.org/10.1016/j.orl.2013.05.003

为你推荐

友情链接