Robust design of a strategic network planning for photovoltaic module recycling considering reclaimed resource price uncertainty

Abstract PhotoVoltaic (PV) power is one of the rapidly growing solar energy technologies worldwide. The installed PV power capacity has considerably increased over the past decades. Consequently, the End-of-Life management of used PV modules is becoming very urgent. This article investigates strategic recycling network planning for recycling PV module considering reclaimed resource price uncertainty. Based on a real case setting, the problem is first formulated as a risk-neutral model and then extended to a risk-averse model that considers the risk preferences of investors. Moreover, robust reformulations of the risk-neutral/risk-averse models are proposed to hedge against ambiguity in the probability distribution of the uncertain resource price. An outer approximation based solution approach is proposed to solve the robust reformulations. Numerical experiments with test data generated from the real case are carried out to demonstrate the benefits of the resulting robust model and managerial insights are explored. Lastly, conclusions are drawn and future research directions are outlined.


Introduction
Renewable energy sources such as solar energy, wind power, hydraulic power, and geothermal energy are playing an increasingly significant role in modern society, due to increasing environmental and energy issues. Solar energy is radiant light and heat from the sun that is harnessed using various continuously evolving technologies, such as solar heating, photovoltaics, solar thermal energy, solar architecture, molten salt power plants, and artificial photosynthesis. Compared with other energies, solar energy has significant advantages, such as providing reliable, secure, and independent energy (https://cleantechnica. com/2013/10/08/advantages-disadvantages-solar-power/). In 2011, the International Energy Agency (IEA) pointed out that: the development of affordable, inexhaustible and clean solar energy technologies will have huge longer-term benefits. It will increase countries' energy security through reliance on an indigenous, inexhaustible and mostly import-independent resource, enhance sustainability, reduce pollution, lower the costs of mitigating global warming, and keep fossil fuel prices lower than otherwise. These advantages are global.
The majority of solar energy technologies that are currently on the market are based on PhotoVoltaic (PV) effects, whereby an electric current is produced in materials when exposed to light (http://www.rsc.org/campaigning-outreach/ global-challenges/energy). PV power is currently one of the rapidly growing power-generation technologies worldwide (Cucchiella et al., 2015). The PV market has continuously grown since the late 1980s (Fthenakis, 2000), and the market for PV electricity generation has boomed over the past decade. According to statistical data from the IEA Photovoltaic Power Systems Program (Dudley and Dale, 2017), the cumulative installed PV power capacity reached 301 473 MW in 2016, and the growth rate per year reached 48:9% between 2005 and 2015. Figure 1 depicts the growth trends of the installed PV capacities of seven major countries. Among these countries, China has surpassed Germany, thereby becoming the country with the largest annual installed PV capacity in the world post 2015. Additionally, the anticipated total global installed solar capacity could grow to over 600 gigawatts by 2020 (Sch€ afer and Schmela, 2016).
With the growing installation of PV systems, the End-of-Life (EoL) management of used PV modules is becoming increasingly important (Sener and Fthenakis, 2014;Perez-Gallardo et al., 2018). First, from a legal perspective, the recent update of the European Waste Electrical and Electronic Equipment Directive classified EoL PV modules as electrical/electronic waste. Hence, ensuring the appropriate collection and recovery of EoL PV modules has become mandatory (Cucchiella et al., 2015). Second, from an environmental perspective, a few types of solar modules contain hazardous materials, such as cadmium, tellurium, lead, and selenium, whose toxicity can cause severe water, soil, and wildlife contamination (Mcdonald and Pearce, 2010). Third, from an economic perspective, PV modules have a long lag time from production to decommission (! 25 years). Based 2016; Perez-Gallardo et al., 2018), and deterministic problems are extensively considered. Most recently, Choi and Fthenakis (2014) developed a mathematical modeling framework to evaluate the economic feasibility of the macro-level RL planning (mixed-integer programming) and the micro-level recycling process planning (linear programming schemes) of PV modules. Hsueh and Lin (2015) proposed a network model to rank the alternatives for implementing the sorting process of RLs in the downstream PV industry. Perez-Gallardo et al. (2018) proposed a multi-objective framework to assess the recycling strategies of PV modules, such as closed-loop recycling, open-loop recycling, and semi-closed-loop recycling.
Despite the inherent uncertainties in the PV recycling network, related studies that consider uncertainty are scarce. Degel and Walther (2016) presented a robust strategic planning approach to analyze the early installation of appropriate collection and recycling infrastructures with a focus on resource criticalities (Achzet and Helbig, 2013;Gleich et al., 2013). Their approach was a scenario-based approach using the Hurwicz criterion, which calculates a linear combination of the best and worst objective values.

CLSC and RLs network design
The research on CLSC and RLs was triggered by a few early pioneering works (Fleischmann et al., 1997) and has attracted increasing attention from academia and practitioners over the past decade. Comprehensive reviews can be found in Guide and Wassenhove (2009), Ghiani et al. (2014), Govindan et al. (2015), Govindan and Soleimani (2017). To date, considerable efforts have been devoted to researching the development of different strategies/models for mitigating various types of inherent supply chain uncertainties (Tang, 2006;Tang and Musa, 2011;Heckmann et al., 2015) because these uncertainties, including uncertain demand, return, price, and cost, have significant impacts on supply chain network design.
Most of the related studies propose Mixed-Integer Linear Programming (MILP) models with the objective of expected costs/profits using two-stage stochastic programming. In addition, scenario-based approaches have been extensively adopted to solve the resulting models. Realff et al. (2004) designed a robust reverse production network for a carpet recycling problem, which was formulated as a scenario-based MILP to maximize net profits. Listeş and Dekker (2005) and   Listeş (2007) presented stochastic mixed-integer programming models with several scenarios to design product recovery networks considering demand and supply uncertainties. Salema et al. (2007) proposed a capacitated multi-product RL network model with product demand and return uncertainty using a multi-scenario approach. Lee and Dong (2009) proposed a two-stage stochastic programming model for modeling dynamic location and allocation problems with uncertain demand for forwarded products and the supply of returned products. El-Sayed et al. (2010) developed a stochastic MILP model to study a multi-period multi-echelon forward/reverse logistics network design under an uncertain demand problem. Amin and Zhang (2013) presented a multi-objective stochastic model (scenario-based) for a CLSC network under uncertain demand and returned products. Fattahi and Govindan (2017) addressed the design and planning of a stochastic forward/reverse logistics network over a planning horizon with multiple tactical periods. The collected amounts of used products with different quality levels were assumed to be dependent on the offered acquisition prices. Badri et al. (2017) developed a two-stage stochastic programming model for the value-based CLSC network design with uncertain demands and returns. The risk-averse modeling approaches have been investigated to obtain "immunized" solutions (Ahmed, 2006), because the aforementioned models, which are called riskneutral approaches in the literature, are not the adequate form for addressing CLSC/RL problems (Zeballos et al., 2016). Among Value at Risk (VaR), Conditional Value at Risk (CVaR), and Mean-Absolute Deviation (MAD), which are three of the well-known risk measures in financial and management optimization problems, CVaR is extensively adopted as a risk measure, due to its mathematical advantages such as subadditivity and convexity. Studies considering risk aversion incorporated CVaR into the framework of two-stage stochastic programming and scenario-based approaches were developed to address the resulting models.  and Soleimani and Govindan (2014) proposed two-stage stochastic programming models with mean-CVaR objective functions to address CLSC and RL network design problems under uncertainty, respectively. In addition, scenario-based approaches were adopted to linearize the CVaR functions. Zeballos et al. (2016) presented a comprehensive risk-averse multi-stage model for the design and planning problem of a CLSC with uncertain demands. Five objective functions that maximize the profit considering risk measures (CVaR) to manage the volatile conditions of the market were examined. Rahimi and Ghezavati (2018) proposed a risk-averse two-stage stochastic programming model to study a multi-period multi-objective RL network design problem under uncertainty.
Stochastic programming techniques generally require perfect information of random variable probability distributions (Klibi and Martel, 2010). Estimating the actual distribution of uncertain parameters is impossible, due to the absence of historical data (Zokaee et al., 2017). The Robust Optimization (RO) approach is an alterative approach for obtaining solutions immunized against uncertainty. RO was originally developed by El- Ghaoui and Lebret (1997) and Ben-Tal and Nemirovski (1998), and it has been successfully applied in many application fields, such as facility location (Baron et al., 2011) and portfolio optimization (Bonami and Lejeune, 2009;Hu and Mehrotra, 2015). Compared with the stochastic programming frameworks and scenario-based approaches, the major benefit of the RO approach is its tractability (Baron et al., 2011). Pishvaee et al. (2011) applied the RO approach for a CLSC network design problem and presented a tractable robust counterpart of a proposed deterministic MILP model by assuming stochastic demands and costs in an uncertain closed box set. Ramezani et al. (2013) presented a robust design for a multi-product, multi-echelon, closed-loop logistic network model by minimizing the maximum regret between the optimal objective function value and the resulting objective function for all possible scenarios. Hatefi and Jolai (2014) utilized the RO approach proposed by Bertsimas and Sim (2003) to cope with uncertain customer demands and the quantities and qualities of returned products for designing a CLSC network. Keyvanshokooh et al. (2016) presented a robust CLSC network design problem based on uncertain demands and returns by following the RO approach proposed by Bertsimas and Sim (2003). Safaei et al. (2017) considered demand uncertainty and presented a multi-period MILP for a cardboard CLSC by using the RO approach proposed by Mulvey et al. (1995).
Overall, although numerous approaches exist for dealing with the inherent uncertainty in the research field of CLSC/RL networks, most of the existing studies on PV recycling network design mainly focus on deterministic modeling, and studies on uncertainty remain scarce. Therefore, the current research fills this gap. Additionally, we summarize the literature on CLSC/RL network design problems incorporating risk measures in Table 1. In our research, we establish a distributional robust model for an RL network design problem considering technology selection decisions and a risk measure (CVaR ). Then, the model is reformulated as a second-order cone program. This work enriches the modeling methodology in the research on CLSC/RL network design problems.

Problem statement
We consider a two-layer recycling network with multiperiod recovery of EoL PV modules (PV modules in this study refer to EoL PV modules). According to a real case in Germany (detailed introduction is provided in Section 7), which has a representative recycling network structure, the network consists of collection and recycling centers (Figure 2), and the locations of the facilities are provided. In each period, the EoL PV modules in electricity stations are collected at the collection centers. Then, PV modules are sorted and delivered to the corresponding recycling centers for further processing. PV modules can be treated by different recycling technologies in the recycling centers, which include delamination, material separation and metal extraction or purification. The resources (scarce metals or other materials) acquired from the PV modules are sold for profit.
The recycling treatments of PV modules in the recycling centers are subject to the installed recycling technologies and their processing capacities. Thus, investors should determine the technology and its capacity to be invested in the recycling centers at each period. Apart from the investment decisions, investors must determine the type and amount of PV modules delivered to the corresponding recycling centers and the processing period of each technology. The objective of the planning problem is to maximize the revenue from selling the reclaimed resources from the collected PV modules to the varying price markets and minimize the costs associated with investment, operation, processing, and transportation.
Based on the real case, several assumptions are provided.
1. The prices of reclaimed resources are stochastic and assumed to be mutually independent over different resource types and varying planning periods. The market size of the resources is unlimited and all the reclaimed resources are timely sold. 2. The technologies and the corresponding equipment are available for all subsequent periods after they are invested. Several processing capacities associated with each technology can be installed. Each capacity has different investment costs. 3. The residual value of the invested equipment and technology is neglectable after the planning horizon. Future sales revenue and costs are discounted as a net present value with a discount rate. 4. All the collected PV modules must be processed within a corresponding period due to environmental and policy causes; thus, the inventory costs of PV modules at the end of each period are not considered.

Model formulation
The notation used throughout this article is summarized as follows: Sets and indices: Parameters: Auxiliary decision variables: Decision variables: Notably, a random variable is denoted byã, and boldface lowercase and uppercase letters are used to denote vectors (y) and matrices (Y), respectively. The spaces of vectors and symmetric matrices of dimension n are denoted by R n and S n , respectively.
The expected profit is equal to the expected sales revenue minus the total cost, which is the summation of the total transportation, total investment, total operational, and total processing costs. The incomes and costs are formulated as follows: Total transportation cost in period t: Total investment cost in period t: Total operational cost in period t: Total processing cost in period t:

Discount rate a it
Price of resource i in period t p mpi Yield of technology m reclaimed resource i from PV module p Q npt Amount of PV module p collected by collection center n in period t u jmt Unit processing cost of technology m with capacity j in period t q nlt Unit transportation cost from collection center n to recycling center l in period t d nl Distance between collection center n and recycling center l g jmt Investment cost of technology m with capacity j in period t f jmt The operational cost of technology m with capacity j in period t K jm The size of capacity j for technology m k mp 1 if technology m can process PV module p; 0 otherwise B A sufficiently large positive number h A non-negative weighted coefficient associated with risk preference Confidence level of CVaRr ; 2 ð0; 1 c 0 Control factor of price expectation l The mean of random vectorã (column vector converted from a it ) C The covariance matrix of random vectorã Continuous variable, VaRr of random revenueP V mlpt Continuous variable, amount of PV module p processed by technology m in recycling center l in period t W nlpt Continuous variable, amount of PV module p delivered to recycling center l from collection center n in period t X jmlt Binary variable, 1 if technology m with capacity j in recycling center l in period t is invested; 0 otherwise Y jmlt Continuous variable, available amount of technology m with capacity j in recycling center l in period t Overall, the discounted total cost in period t is In addition, the expected discounted total revenue (P) is In terms of the notation and formulation, the problem is formulated as follows: s:t: X jmlt 2 0; 1 f g; 8j 2 J; m 2 M; l 2 L; t 2 T: Objective function (7) is used to maximize the total expected profits of the recycling network of PV modules, which are equal to the discounted total revenues minus the discounted total costs. Constraints (8) and (9) are material balancing constraints associated with recycling and collection centers, respectively. The former states that the number of PV modules processed in a recovery center is equal to the amount of PV modules delivered to that recovery center. The latter stipulates that the total amount of PV modules delivered from a collection center is equal to the total amount of PV modules collected at that collection center. Constraint (10) indicates that the total amount of PV modules processed does not exceed its processing capacity, which is determined by the invested recycling techniques. Constraint (11) ensures that one PV module can only be processed by its eligible technology. Constraint (12) states that a recycling technology is available during all the forthcoming periods once it is invested. Constraints (13) and (14) are domain constraints.
Problem (7) can be rewritten as a risk-neutral two-stage stochastic programming problem based on the two-stage stochastic programming framework proposed by Birge and Louveaux (1997). The framework in its general form is formulated as follows (Noyan, 2012): where f ðX; xÞ ¼ c T X þ QðX; nðxÞÞ is the objective function of the first-stage problem and QðX; nðxÞÞ is the secondstage objective corresponding to the realization of the random data nðxÞ for the elementary event x.
We define an auxiliary decision variable z it as follows: which represents the amount of reclaimed resource i obtained in period t. Letã t ¼ fa it g i¼1;:::;I 2 R I andã ¼ fã t g t¼1;:::;T 2 R IÂT be vectors of the resource price. In addition, z (I Â T dimensions) is a column vector that is converted from z it by vectorization transformation. Thus, problem (7) is rewritten as the following risk-neutral twostage stochastic programming model: The objective function of the second stage problem (18), Pðz;ãÞ, is the revenue from the sale of reclaimed resources and is dependent on the realization of the uncertain priceã.

Risk-averse model incorporating risk preference
A risk-averse model that incorporates risk measures is proposed because the risk-neutral approach is inappropriate for a risk-averse investor (Ling et al., 2017). Although many variabilities exist, such as technology updates and cost fluctuations in the PV industry, the investment risk of recycling technology due to price volatility of reclaimed resources receives more attention from the perspective of investors (Degel and Walther, 2016). In addition, the construction of recycling centers and the development of recycling technologies are longterm strategic decisions. Thus, investments in PV module recycling technology still present high risks despite its increasing requirements for solving environmental and social issues. Incorporating risk preference into the modeling process is essential for mitigating the investment risks.
CVaR is a well-known risk measure that is defined in the context of loss random variables (Rockafellar and Uryasev, 2000). Compared with most CVaR models, we deal with a profit maximization problem. Thus, we conveniently define a CVaR-like risk measure, which is denoted as CVaRr , as a measure of satisfaction or utility function with the random revenueP in a manner similar to Zeballos et al. (2016), as follows: Definition 1. Let FP ðÁÞ represent the cumulative distribution function of the random revenue variableP. For any confidence level 2 ð0; 1, the VaRr ðPÞ is defined as the following -quantile: where VaRr denotes the minimum possible revenue during the time excluding the worse outcomes whose probability is less than . CVaRr ðPÞ is defined as the mean of the -quantile distribution ofP, that is, Notably, CVaRr considers the variability of revenues and provides a viable risk measure to bound the probability that the design and planning of a given RL network incur a large decrease in revenues under extreme scenarios. In the profit maximization context, a large CVaRr indicates a high capability of the logistics network to hedge against the risk of revenue variability.
Lemma 1. The CVaRr of revenueP can be approximated as the following one-variable convex optimization problem (21): where ðxÞ þ ¼ maxf0; xg, is a given confidence level and v is the VaRr of random revenueP.
Proof. The convex approximation of CVaRr is similar to that of the classic CVaR, and the technical details can be found in Rockafellar and Uryasev (2000).
w Furthermore, we propose a mean-risk model to mitigate the investment risks using CVaRr as the risk measure.
Proposition 1. The risk-averse two-stage stochastic programming model is formulated as follows: s:t: 8 ð Þ; 10 ð Þ; 11 ð Þ; 13 ð Þ; 16 ð Þ: Proof. Based on a well-known two-stage stochastic framework (Noyan, 2012;Soleimani and Govindan, 2014): where h is a non-negative weighted coefficient of the risk part and the large h stands for additional risk-averse investors, we obtain the mean-risk model by replacing qðÁÞ in Equation (24) with CVaRr : Based on the translation invariance of coherent risk measures (Noyan, 2012), CVaRr ða þ ZÞ ¼ a þ CVaRr ðZÞ holds, and we obtain: Rðz; v;ãÞ ¼ã T z þ hvÀ h ðvÀã T zÞ þ is defined, and then the proof is completed. w

Robust reformulation
A robust reformulation of the risk-averse two-stage stochastic programming model is proposed due to the existence of market volatility, which implies that the probability distribution A of the random priceã cannot be precisely identified. Hedging against ambiguity in probability distributions by using the minimal mean profit or utility over a set of possible probability distributions is prudent. We reformulate problem (22) as a max-min formulation: Furthermore, we decompose inf A E A ½Rðz; v;ãÞ into two subproblems: Subproblem 1 is the worst-case expected revenue, and subproblem 2 is the worst-case CVaRr . We utilize two different probability distributional sets to derive the robust formulations of the two subproblems, because utilizing identical probability distributional sets to derive robust formulations for different stochastic terms may lead to overconservatism or over-optimism (Goh and Sim, 2010).

Robust formulation of expected revenue: E A ½Pðz;ãÞ
The set of probability distributions is generally described by a set of distributions with known moments or bounded moments, such as first-and second-order moment information (Bertsimas et al., 2010;Delage and Ye, 2010). As subproblem 1 is the minimization of expected revenue E A ½Pðz;ãÞ, the general distributional sets with known moments are unsuitable for this subproblem, due to the existence of first moment ambiguity. We introduce a distributional set A 1 to describe the following uncertain first moment and avoid unduly optimistic risk assessments : where U is the set of all probability distributions on the measurable space and M is any closed convex set known to contain the support of the first moment E A ½ã. Suppose that the first moment E A ½ã belongs to an ellipsoid centered at an estimated mean l and C is an estimation of the covariance matrix ofã; l 2 R NÂT ; C 2 S NTÂNT . In addition, c 0 ! 0 controls the size of A 1 . If c 0 ¼ 0, then the volatility of the first moment E A ½ã no longer exists. Then, we reformulate subproblem 1 over distributional set A 1 : Constraint (31) can be expressed as an affine ball transformation of radius c 0 : Under ellipsoidal uncertainty: Thus, subproblem 1 is rewritten as

Robust formulation of
CVaRr : E A ½hv2 h ðv2ã T zÞ þ We employ a probability distributional set defined by Ghaoui et al. (2003) and Bertsimas et al. (2010) to obtain the robust formulation of subproblem 2: U and M have similar definitions as distributional set A 1 in Equation (29); l 2 R NÂT and Q 2 S NTÂNT are the first and second moments of random vectorã, respectively; and C 2 S NTÂNT is a covariance matrix. A 2 assumes that the first and second moments of random variables are exactly known and contains all the distributions with those moments. Since it helps to derive the analytical approximation with certain means and variances, A 2 is often used to formulate the worst-case CVaR (Ghaoui et al., 2003).
To obtain the robust formulation of subproblem 2 with a 2 A 2 , we first derive an upper bound of E A ½ðv Àã T zÞ þ . By the property of we have: Then, based on Jensen's inequality: Thus, we obtain: Then, the lower bound of E A ½hvÀ h ðvÀã T zÞ þ can be written as Combining the two robust formulations, we derive Proposition 2 as follows.
Proposition 2. Problem (28) is formulated as the Second-Order Cone Programming (SOCP) model: Proof. A ccording to Equations (34) and (39), problem (28) can be reformulated as: Defining which is a monotropic function with independent variable v, we obtain the optimal value of v by taking the first-and second-order derivatives with respect to v as follows: Therefore, HðvÞ is concave. Let H 0 ðvÞ¼0. Then, we obtain the maximum value s Furthermore, the following auxiliary decision variable c that is subject to a second-order cone constraint is introduced: c ! H Ã ðvÞ. In addition, due to the convex square-root function in Equation (42), the proof is completed by substituting c into Equation (42).

Solution approach
We propose an OA algorithm to solve the resulting model in this section. OA is a cutting plane algorithm originally proposed by Duran and Grossmann (1986). It iteratively decomposes a Mixed-Integer Nonlinear Program (MINLP) into a series of mixed-integer linear Master Problems (MPs) and nonlinear SubProblems (SPs). For a program with a maximized objective, the SP is able to obtain the values of continuous decision variables given fixed integer assignments and provides a lower bound of the MINLP. The MP is enhanced by adding valid cuts associated with nonlinear terms, and it provides an upper bound. The OA algorithm has been successfully applied in various areas, such as in the chemical industry (D ıaz and Bandoni, 1996), facility location problems (Shahabi et al., 2014;Ljubi c et al., 2018) and CLSCs (Zhang and Unnikrishnan, 2016;Liu and Zhang, 2017). Interested readers can refer to Duran and Grossmann (1986) and Bonami et al. (2008) for more details about the OA algorithm. Bonami et al. (2008) proved that the OA algorithm can converge to a global optimal solution for a convex MINLP. Therefore, the convexity of the proposed SCOP model shown in Proposition 2 ensures the optimality and applicability of the OA algorithm.

Initialization
The most conservative strategy for the considered problem is to invest all technologies with all types of capacities in every recycling center and every period, i.e., suppose that X 0 jlmt ¼ 1; 8j 2 J; l 2 L; m 2 M; t 2 T. Thus, we obtain an initial assignment solution of the integer decision variables.

OA subproblem
The subproblem is a nonlinear program that provides a feasible solution for the SOCP. After optimizing the SP, a corresponding lower bound can be obtained. The SP at the hth iteration is formulated with fixed-integer decision variables X h and is described as follows: where r is an auxiliary decision variable and Þð1 þ bÞ Àt :

OA master problem
The OA master problem is a linear approximation of the SOCP model, which provides an upper bound for the maximization program. Lemma 2 shows how the linear approximation of nonlinear constraint (47) Through simple algebra, inequality (49) holds, and this completes the proof. w In addition, define Q h ¼ fðj; m; l; tÞ :X h jmlt ¼ 1g and Q 0h ¼ fðj; m; l; tÞ :X h jmlt ¼ 0g at the hth iteration of the OA algorithm, respectively, and inequality (52) Thus, the MP at the h th iteration is summarized as follows: MP : maxd; Constraint (53) is the definition of the objective function. Constraint (54) states that the solution of the MP is larger than the lower bound at every iteration h. Additionally, this constraint ensures that the OA algorithm can be terminated and that a eÀoptimal solution can be obtained once the MP is infeasible (Fletcher and Leyffer, 1994). The procedure of the OA algorithm is summarized in Algorithm 1.

Numerical results
We conduct numerical experiments with test data generated based on real data in Germany to evaluate the performance of the resulting models and to explore managerial insights associated with the recycling network of PV modules. We evaluate the computational performance of the OA algorithm by comparing with CPLEX and the robustness of the solutions by comparing them with the SAA. The impacts of investor risk preference and technological advancement on the recycling network are also considered. A topology structure analysis based on North Rhine-Westphalia (NRW) is finally presented. We describe the test environment in detail before presenting the results.

Test environment and parameter estimation
Since the resulting model is a SOCP model, we solve it directly by CPLEX 12.7 and perform a comparison with the proposed OA algorithm in Section 7.2. The proposed OA algorithm is coded in Cþþ. The MP is solved by the MILP solver of CPLEX with the Benders strategy, which indicates that CPLEX applies a Benders decomposition technique to solve the linear MP (IBM, 2017), and the value of e is set to one. The SP is solved by the MIQCP solver of CPLEX 12.7. Moreover, we employ the SAA method, which is coded in Cþþ and solved by CPLEX 12.7, as a benchmark approach to evaluate the robustness of the SOCP in Section 7.3. Each problem instance is conducted on a computer with an Intel Core i9 Duo 3.40 GHz processor and 32 GB of RAM running the Windows 10 operating system. All reported CPU times are in seconds.
The test data in this research are generated based on real data collected from NRW, Germany. The considered recycling network of PV modules consists of 23 fixed collection centers (jNj¼23) and 10 potential recycling center locations (jLj¼10), which are illustrated in Figure 3 The pentagram in the figure represents the recycling center, and the black spot denotes the collection center. The planning horizon of the recycling network is between 2016 and 2030 (that is, jTj¼15).
Four PV modules (c-Si, a-Si, CdTe, and CIGS, jPj¼4) are considered, and the amount of modules at each period is set to their yearly installed amount between 1991 and 2006 (assuming that the life cycle of the PV modules is 25 years). Correspondingly, four recycling technologies (jMj¼4) are considered. Each recycling technology and the corresponding PV modules are listed in Table 2 Notably, Technology 4 can process all four types of PV modules, which is referred to as flexible technology in the remainder of this section. Technologies 1 and 2 can process c-Si and a-Si with different investments, operational and processing costs, and yields. Seven types of reclaimed resources (aluminum, silicon, cadmium, silver, copper, tellurium, and indium, jIj¼7)  can be recycled from the four types of PV modules by the given yields (p). Based on price panel data over the past 40 years obtained from the U.S. Geological Survey (Kelly and Matos, 2015), we predict the expected price of reclaimed resources over 15 years by the ARIMA model implemented with the "forecast" package of R 3.4.2. The cost parameters were obtained via determination of the actual costs for enterprises in Germany. All the parameter values used are available upon request, which include the predicted expected prices of reclaimed resources, yield, amount of collected PV modules, distance, capacity, and cost coefficients. Due to the high data volume, the data are stored in Dropbox and can be downloaded from https://www.dropbox.com/s/zauf9c6n2-qiw35w/Parameter\%20data.docx?dl ¼0.

Comparison between the OA algorithm and CPLEX
To evaluate the performance of the proposed OA algorithm, we conduct 10 instances for each combination of jNj and jLj by multiplying the values of the mean random price (l), standard deviation (C), and other parameters (p mpi , Q npt , u jmt , d nl , q nlt , g jmt , f jmt , K jm , and k mp ) by (1 þ i), where i is drawn uniformly from ½À0:1; 0:1. The instances are solved by CPLEX and the proposed OA algorithm. Both the OA algorithm and CPLEX obtain the optimal solutions within a reasonable running time. The mean and standard derivation of the CPU times are reported in Table 3 Compared with the results solved directly by CPLEX, the proposed OA algorithm is superior to CPLEX in terms of CPU time in the large-scale instances. The results with different values of risk parameters (h, and c 0 ) are reported in Section 1 of the supplementary material. 7.3. Performance analysis of the robust approach compared with the SAA approach We employ the SAA method as a benchmark approach to evaluate the performance of the resulting robust model. A scenario-based model associated with the risk-averse twostage stochastic programming model is proposed in Appendix A to implement the SAA method. The expected random prices (l) of the reclaimed resources are obtained by the ARIMA forecast model based on the historical price data. We generate samples uniformly drawn from the interval [0, 2l] to observe the impact of price volatility on the recycling network. The prices of reclaimed resources are also assumed to be mutually independent. Thus, the variance of the prices is m 2 /3, and the covariance matrix (C) is a diagonal matrix. That is, the ith diagonal entry is mi 2 /3 and the other entries are zero. We randomly selected jNj recycling centers and jLj collection centers from 10 recycling centers and 23 collection centers in the test data, respectively. The numbers of decision variables and constraints in the SAA and SOCP models are presented in Table 4 We take the pair of ðjNj; jLjÞ with (10, 5), (15, 10), and (23, 10), and three groups of tests are conducted. Due to space limitations, we only present the computational results for certain risk aversion parameters h ¼ 1, ¼ 0:01, and c 0 ¼ 0:05 in Tables 5 and 7. Moreover, the results associated with other values of h, and c 0 are reported in Section 2 of the supplementary material, which demonstrate the consistent conclusions. Table 5 reports the performance (costs, gaps, and CPU time) of the SAA and SOCP approach with different problem sizes. Although the CPU times of the approaches increase when the problem size increases, the CPU time of the SAA is higher than that of the SOCP method. Moreover, when the problem size becomes large (ðjNj; jLjÞ ¼ ð15; 10Þ; ð23; 10Þ), the transportation and processing costs obtained by the SOCP are smaller than those calculated by the SAA, which is consistent with the result in Table 6 Rows 3-6 in Table 6 present the amount of PV modules processed by the corresponding technology and row 7 presents the total transport volume (multiplication of distance and flow) of PV modules. As shown in Table 6, the amount of PV modules processed by Technology 1 in the SOCP exceeds the counterpart in the SAA, and the total transport volume in the SOCP is less than its counterpart in the SAA, which implies that the solution obtained by the SOCP is more conservative and that the investor tends to utilize technology with a cheaper unit processing cost to process PV modules and chooses a closer recycling center to be served. Consequently, the transportation and processing costs obtained by the SOCP approach are reduced. In contrast, the investment and operational costs are the same because the marginal cost of adjusting investment decisions  Table 4. The numbers of decision variables and constraints in the SAA and SOCP approach.

Continuous Binary variables
variables is significant. Thus, the investment decisions are not impacted by the robust approach. As the revenue, profit, and CVaRr are affected by the samples, we utilize the Monte Carlo simulation method to evaluate the three performance indicators of the SAA and SOCP approaches. A set of 100 000 samples of the price are uniformly drawn from [0, 2l] and the mean and variance of the revenue and profit of both solutions obtained by the SAA and SOCP approaches are then calculated. CVaRr is computed by averaging the profit of the worst of 100 000 samples. The comparison results are reported in Table 7, and the relative differences between the expectation and variance of the profits obtained by the SAA and SOCP methods are respectively defined as: As shown in Table 7, the variances of the profit and revenue of the SAA are both larger than those of the SOCP approach, whereas the averages of the profit and revenue are opposite when the problem size is (15, 10) or (23, 10). Moreover, according to the values reported in the columns labeled by Gap EP and Gap V P , we conclude that the SOCP approach can considerably reduce the profit variation at the expense of a relatively small decrease in expected profit (Gap EP Gap V P ). The robustness of the solutions obtained by the SOCP approach is validated. CVaRr obtained by the SOCP method is larger than that of the SAA, which demonstrates that the recycling network obtained by the SOCP approach has higher profit under the worst cases than that of the SAA. In addition, the profit is negative due to the mandatory recycling of all PV modules, which leads to high construction and operational costs of the recycling network and low revenue from selling reclaimed resources.

Impacts of risk preference
In this section, we conduct numerical experiments to explore managerial insights into the risk preference of investors. Three risk-preference parameters are presented in model (40). denotes the confidence level of CVaRr , which is defined in Definition 1. CVaRr quantifies the mean value of the worst of the total profits. Therefore, a smaller represents more averse behavior of investors. h is the weight of CVaRr in the two-stage stochastic model. A higher h means that the investors have more risk-averse preference. Then, we conduct numerical experiments with different risk-averse parameters (, h, and c 0 ) to examine their impacts on network performance. Figure 4 reports the impacts of risk preference on the risk measure (CVaRr ). We set three cases for ¼ (0.01, 0.05, 0.1) and h ¼ ð0:1; 0:3; :::; 2Þ, and c 0 is set to 0.05 (Due to space limitations, we only present the computational results for c 0 ¼ 0:05 in Figures 4 and 5. The results associated with other values of c 0 (0.01, 0.1 and 0.95) are reported in Figures  1-6 in Section 3 of the supplementary material, which demonstrate the consistent conclusions). According to Figures  4(a), 4(b) and 4(c), the following observations are obtained: 1. As h increases, CVaRr has an upward trend in Figures  4(a), 4(b), and 4(c), which is consistent with the meanrisk model (22). This finding indicates that investors attempt their best to enhance CVaRr to reduce risk when the risk-averse level rises (high h). 2. CVaRr increases when rises by comparing Figures 4(a), 4(b), and 4(c). In addition, the value of CVaRr slowly grows and tends to stabilize after h increases to a certain extent, which means that CVaRr is subjected to intrinsic price volatility and has an upper bound to be maximized. Marginal risk reduction decreases as the risk-averse level increases. Thus, risk-averse investors cannot limitlessly reduce risk in practice. Figure 5 demonstrates the revenue, total cost, and profit of the recycling network using various risk preference parameters. The observations are as follows: 1. The revenue, total cost, and profit decrease as h grows.
This finding implies that risk-averse investors (large h) are inclined to make conservative decisions, such as scaling back on investing in technology that can process the resource with volatile price to avoid extreme loss. Thus, the total cost decreases and the revenue and profit reduce. 2. Unlike h, the revenue, total cost, and profit increase when the risk parameter increases. This finding indicates that additional risk-seeking investors (large ) are willing to expand investments for high profit.  Figure 6. Impacts of c 0 on CVaRr and profit on (a) Initial yield (b) Initial yield and two yields with growth rate Next, we discuss the impacts of risk parameter c 0 on the recycling network. In Figure 6, and h are set to be 0.05 and 1, respectively. CVaRr and profit under different c 0 are illustrated in Figure 6. The following results are obtained: 1. From Figure 6(a), the profit has a decreasing trend when c 0 increases. As c 0 increases, risk-averse investors tend to adjust their strategy for investing in recycling technology and obtain additional reclaimed resources with low price volatility. Thus, investors can reduce volatility by sacrificing profit.
2. From Figure 6(b), CVaRr rises when c 0 increases. Indeed, additional risk-averse investors attempt to decrease revenue variance, and the recycling network becomes resistant to risk. Thus, CVaRr rises with low revenue variance.

Impacts of technological advancement
Technological advancement has a significant impact on the recycling network. In the wake of technological advancement, the yield of PV module recovery technology increases or new technology with a high yield emerges. Therefore, we further examine the impact of the yield of the technology on the recycling network. Figure 7 illustrates the trade-off curves of the profit and CVaRr under different yields by using the weighting method (Cohon, 1978). The observations are summarized as follows: 1. A trade-off between the CVaRr and profit illustrated in Figure 7(a) signifies that if investors want considerable investment profits, then they must take corresponding risks. According to the trade-off curves in Figure 7(a), investors can select an optimal investment policy that best represents their risk aversion level. 2. As shown in Figure 7(b), the trade-off curves with high yield are in the upper right of the figure whereas those with low yield are in the bottom left. This result demonstrates that the recycling network with high yield obtains a large profit at the same risk level. In other words, on the premise of having the same profit, the recycling network with high yield has large risk resistance (high CVaRr ). The managerial insight from this finding implies that intensifying efforts on technological advancement can effectively reduce risk and increase profit. Moreover, we discuss the impacts of h on the CVaRr and profit under different yields, and the results are reported in Figure 8. First, at the same level of risk aversion (same h), CVaRr and profit are larger when the recycling network has higher yield, which also confirms the conclusions obtained from Figure 7. Then, we find that CVaRr (profit) increases (decreases) sharply in each subfigure when the risk aversion parameter (h) reaches a certain level (h ¼ 0:5, h ¼ 1, and h ¼ 1:7 in Figures  8(a), 8(b) and 8(c), respectively). These sharp increases (decreases) are caused by the adjustment of technology investment, and subsequent slight changes of CVaRr and profit are caused by the adjustment of transportation and processing decisions. Investors become more risk-averse with increasing h, and they reduce their technology investment when the risk aversion is beyond a critical level that they are willing to withstand. These critical levels are different among the three subfigures with different yields, which suggests that recycling networks with higher yield are more risk-resistant and there is no need to adjust investment decision even if the investors' risk aversion slightly increases. Thus, pursuing technological advancement is an effective approach for enhancing the risk resistance capability of recycling networks regardless of the change in investors' risk preference.

Analysis of topology structure
We compare topology networks determined by the solutions of the robust models for different risk preferences. The geographic data are adapted from real data in NRW. We consider different risk-averse levels by varying the values of h and , while risk parameters c 0 ¼ 0:05 and other parameters' values are consistent with those in Sections 7.1 and 7.3 (problem size ðjNj; jLjÞ ¼ ð23; 10Þ, price data, and cost parameters, among others). Figure 9 illustrates two representative topology structures of the recycling network under different risk preferences (h and ). In general, Figure 9(a) presents the result obtained by riskneutral investors, and Figure 9(b) is the result obtained by risk-averse investors. In the figure, black dots denote the locations of the collection centers; red or blue circles indicate the recycling centers, where the color of the circles represents the technologies invested in that recovery center and the size of the circles represents the total amount of PV modules processed in that recycling center by the corresponding technology; the connecting lines illustrate the amount of PV module proportionally delivered to the recycling centers from collection centers (the black dots are not linked, which indicates that the corresponding collection centers have no PV modules to be recovered); and the purple numbers represent the investment period of the recycling technology. Reviewing Figure 9, we obtain the following observations: 1. Risk-averse investors associated with Figure 9(b) invest fewer recycling technologies than risk-neutral investors associated with Figure 9(a). Risk-averse investors generally tend to invest fewer technologies to mitigate the investment risks because technology investment is a long-term strategic decision with large fixed capital investments. 2. Several collection centers are served by a single recycling center in Figure 9(a), while other collection centers are served by two recycling centers in Figure 9(b). This finding implies that investors can mitigate the risk of the recycling network due to price uncertainty by assigning multi-recycling centers for each collection center. 3. Technology 4 is invested in the first year, and other technologies are invested after 12 years. As previously mentioned, Technology 4 can process four types of PV modules, whereas the other technologies only process one or two PV modules. The flexible technology (Technology 4) is selected in the initial period, which can satisfy the requirements of all types of PV modules. Then, as a multi-period decision, investors generally invest as late as possible when capacity is satisfied, and postponed strategies can reduce the total operational costs, such as equipment maintenance costs. 4. In addition, a certain technology (Technology 1 shown in Figure 9(b)) rather than the flexible technology (Technology 4 shown in Figure 9(a) is invested in the 12th period. Compared with Technology 1, Technology 4 is more flexible and has a higher yield for reclaimed resources with high price volatility. Risk-averse investors are inclined to invest in technology that recycles additional reclaimed resources with low price volatility when the capacity is sufficient for different technologies, thereby effectively resisting the risk due to price fluctuations.
We do not report other topology analysis experiments with different values of and h because the same observations are obtained.

Conclusions and future research
This article investigates a multi-period PV module recycling network planning problem under resource price uncertainty. The problem is first formulated as a profit maximization riskneutral model, and a risk-averse model incorporating risk measure (CVaRr ) is then proposed. Due to high degrees of resource price uncertainty in the market, which implies that precisely estimating the probability distribution of the price is difficult, robust risk-neutral/risk-averse models, which are computationally tractable, are formulated. An OA algorithm is developed to solve the robust model. Based on real data from Germany, extensive numerical experiments are conducted to demonstrate the benefits of the robust model.
Moreover, managerial insights are explored, and several helpful insights are summarized as follows: 1. As the risk-averse level of investors increases, investors are inclined to take conservative decisions and scale back on the investment in technology that recycles resources with volatile prices to avoid extreme loss. Thus, the total cost decreases and return (revenue, profit) reduces ( Figure 5). Consequently, the risk of investments is mitigated at the expense of a low profit. 2. The recycling network with a high recycling technology yield has large risk resistance and profit (Figure 7). Investors should continuously intensify efforts on technological advancement to improve their competitiveness. 3. Investors can reduce the risk of the recycling network due to price uncertainty by assigning multiple recycling centers to provide service to each collection center. In addition, as a multi-period decision, the preferred strategies are as follows: investing in a flexible technology in the early period to satisfy the availability of all types of PV module processing capabilities, and investing in certain technologies rather than a flexible one as late as possible on the premise of sufficient capacity, which can cut down total operational costs and reduce financial risk ( Figure 9).
Finally, several future research directions are outlined. First, additional sophisticated prediction technology will be selected to predict the parameter values, such as resource price and collected amount of PV modules. Second, the robust PV module recycling network planning problem considering uncertain parameters such as demands and unit costs is worth studying. Third, incorporating additional decisions, such as location and assignment decisions associated with collection centers and electricity stations, into the models is interesting. Fourth, it is worth comparing the proposed OA algorithm with the generalized Benders decomposition, which is also a well-known solution approach that solves a convex program by extending the Benders decomposition algorithms that can only solve a linear program (Geoffrion, 1972).