Abstract

Silicon is one of the commonly used semiconductors for various industrial applications. Traditional silicon synthesis methods are often expensive and cannot meet the continuously growing demands for high-purity Si; electrodeposition is a promising and simple alternative. However, the electrodeposited products often possess nonuniform thicknesses due to various sources of uncertainty inherited from the fabrication process; to improve the quality of the coating products, it is crucial to better understand the influences of the sources of uncertainty. In this paper, uncertainty quantification (UQ) analysis is performed on the silicon electrodeposition process to evaluate the impacts of various experimental operation parameters on the thickness variation of the coated silicon layer and to find the optimal experimental conditions. To mitigate the high experimental and computational cost issues, a Gaussian process (GP) based surrogate model is constructed to conduct the UQ study with finite element (FE) simulation results as training data. It is found that the GP surrogate model can efficiently and accurately estimate the performance of the electrodeposition given certain experimental operation parameters. The results show that the electrodeposition process is sensitive to the geometric settings of the experiments, i.e., distance and area ratio between the counter and working electrodes; whereas other conditions, such as the potential of the counter electrode, temperature, and ion concentration in the electrolyte bath are less important. Furthermore, the optimal operating condition to deposit silicon is proposed to minimize the thickness variation of the coated silicon layer and to enhance the reliability of the electrodeposition experiment.

Introduction

Silicon is one of the most important and commonly used semiconductors for the applications of computer chips, photovoltaic cells, and optoelectronic devices, due to its chemical stability and unique optical and electronic properties [13]. Various methods have been used to synthesize the Si-based nanostructures, such as sputtering, plasma deposition from silane, laser ablation, and electron-beam evaporation. [47]. However, these methods are usually expansive, complicated, and time-consuming procedures, which are based on the use of high-temperature vapors and melts, and cannot uniformly deposit silicon onto complex three-dimensional substrates [8]. With the continuously growing demand for high-purity Si, the development of advanced deposition techniques to produce low-cost Si devices has become a research interest.

Electrochemistry is a valuable tool for the development of perfectly covering and conformal coatings because it allows accurate controls over the size of deposits and the deposition rate by manipulating the operating parameters, e.g., overvoltage, current density, and bath composition [9,10]. Besides, it also provides some other merits, such as the good compatibility between deposit layers and most conducting substrates, ease of scale-up, and applicability to a wide variety of different shapes and surface geometries [11,12]. Previous studies reported the potential to electrochemically synthesize silicon and other semiconductors in organic solvents and molten salts [1315]. In this study, we adopt the electrodeposition process in the organic solvent with SiCl4 as solute [8]. During the reactions, SiCl4 is electrochemically reduced to form silicon through the following equations:
$SiCl4 + 4e− → Si + 4Cl− (working electrode)$
(1)
$4Cl−→ 2Cl2+ 4e−(counter electrode)$
(2)

Nevertheless, even though electrodeposition is a relatively simple method, various sources of uncertainty may inevitably be introduced during the operating process, such as the variation of the experimental setting, deviation of material properties from the nominal values, and the edge effect, which could result in electrochemical heterogeneity and uneven deposits [16]. Therefore, to improve the quality of the electroplated products, i.e., the uniformity of deposits, it is crucial to better understand the influences of different sources of uncertainty on the uniformity of the resulting products through uncertainty quantification (UQ) analysis and further obtain the optimal experimental conditions. Previously, different experimental and numerical methods have been applied to investigate the influences of various experimental operation parameters on the nonuniformity of the electrodeposition process [1719]. However, these methods are often expensive and time-consuming, which are not suitable to be used for UQ study. In this study, five experimental operation parameters are selected as the primary sources of uncertainty [20], including the distance between the anode (counter electrode) and cathode (working electrode) d, the area ratio of cathode over anode r, the potential of the anode V, the temperature of the electrolyte bath T, and initial concentration of Si4+ in the electrolyte c. The UQ analysis is then performed on the silicon electrodeposition process, which includes two steps. Firstly, a finite element (FE) model is constructed to simulate the silicon electrodeposition process under various experimental conditions and to predict the thickness and uniformity of the deposited layer; second, a Gaussian process (GP) based surrogate model is developed using the simulation results as training data to reproduce the outputs of the FE model efficiently and accurately and to conduct the UQ on the silicon electrodeposition process.

Numerical simulation models like FE analysis have been widely adopted in the electrochemistry field to avoid conducting a large number of experiments that are time-consuming and expensive [2125]. These models can help understand various complex coupled phenomena and the underlying mechanisms, such as electrochemical metal-reduction reactions on moving deposition surfaces, and transport of charged and neutral species via diffusion, electromigration, and convection in the electrolyte, etc. More importantly, through precise prediction of the electrodeposition process, the size and microstructure of deposits can be directly mapped to the operation parameters, e.g., geometries of the electrodes, applied potential, bath composition, pH, and temperature, etc. Chen et al. [26] reported the development of FE models to simulate the two-dimensional (2D) processes of nickel electrodeposition in the lithography, electroforming and molding microfabrication. Species concentrations, electrolyte potential, flow field, and positions of the moving deposition surfaces were computed by coupling the arbitrary Lagrangian–Eulerian approach to track the entire transient deposition processes and to calculate local current densities that influenced the microstructure and functional/mechanical properties of the deposit. Lupo et al. [27] investigated the dendrite formation phenomena during the electrodeposition process using FE methods and validated the model with experimental results. It was found that the FE models could precisely estimate the growth of dendrites and their detailed size distribution over the whole electrode surface. Wheeler et al. [28] presented the modeling of superconformal electrodeposition/superfill using the level set method in the Eulerian framework to simulate the additive accumulation and conservation on a deforming copper/electrolyte interface. However, one critical disadvantage of the FE models is that they are often composed of multiple components, partial differential equations, and physical modules, with the degrees-of-freedom of million level; as a result, the computational cost is rather high, and it is not feasible to implement the UQ analysis via this technique.

Machine learning-based surrogate models have drawn great attention in the material science field recently [29], which are constructed through data-driven, bottom-up techniques. Different from FE models, Machine learning surrogate models can directly capture the general underlying trends of the performances of the systems and efficiently approximate the outputs of large-scale complex systems given any desirable input variables through the realization of a regression model and a random process [30,31]. The hyperparameters and random process in the regression model are trained using the FE simulation results and experimental data. Surrogate models with different predictors, e.g., the response surface method [32], neural networks [33], and GP (as known as Kriging) [3436], have been applied for the structural dynamic analysis of various engineering systems. In the previous study, we adopted GP-based surrogate models to help understand the lithiation-induced stresses [22] and crack initiation phenomenon [25,37] in the Si-based lithium-ion battery anodes over the design space. It was found that the GP model could effectively estimate the lithiation-induced stress concentration and crack formation in the Si layers without losing fidelity.

In this study, the UQ analysis is performed on the silicon electrodeposition process via numerical simulation methods. The overall modeling procedures contain two steps. In the first step, a 2D multiphysics-based FE model is developed to simulate the electrodeposition process; the thickness variations of electroplated silicon layers are calculated with different combinations of operation parameters as input variables, including the distance d and area ratio r between the counter and working electrodes, the potential of the counter electrode V, the temperature of the electrolyte bath T, and initial concentration of Si4+c; the design sensitivity analyses of these parameters are conducted. In the second step, the GP-based surrogate model is trained with the FE simulation results; the propagation of uncertainty from inputs over the output space is evaluated; and the improvement strategy of the experimental setup for silicon electrodeposition is provided to minimize the thickness variations of the coated silicon layer.

Finite Element Model Development

In this section, the multiphysics-based FE model is established in the commercially available FE simulation tool (comsolmultiphysics) to simulate the silicon electrodeposition process. The electrochemical kinetic on the electrode surfaces, the deformation of the electrode boundary, the transportation of Si4+ and Cl and potential change in the electrolyte are calculated. In the following, governing equations that control the reactions and deformation are discussed; then, the model implementation in the fe software is introduced.

Governing Equations.

The electrodeposition process is modeled with tertiary current distribution [38], accounting for the effect of electrolyte composition and ionic strength on the electrochemical process, as well as solution resistance and electrode kinetics, to solve for the potential and current density in the electrodes and the electrolyte, respectively. The fluxes of Si4+ and Cl ions in the electrolyte Ni are described by the Nernst-Planck Eq. (3), including the diffusion term in concentration gradient and the migration term in an electric field. The convection of electrolyte in a flow field is neglected
$Ni=−Di∇ci−ziuiFci∇ϕl$
(3)
where Ni is the flux of ion i in mol/(m2·s), Di is the diffusion coefficient in m2/s, ci is the concentration of the ion i, in mol/m3, zi is the valence, ui is the ith ion mobility in s·mol/kg, ϕl is the electrolyte potential in V, and F = 96487 C/mol is the Faraday's constant. Furthermore, the assumption of electroneutrality is held for species concentrations [39]. The current density il and current balance in the electrolyte are described in the following equation:
$il=F∑i=12ziNi and ∇il=0$
(4)
The current density of charge transfer reactions at the electrode–electrolyte interfaces iloc is expressed as a function of the overpotential η using the Butler–Volmer equation
$iloc=i0[exp(αazFη/RT)−exp(αczFη/RT)]$
(5)

where i0 is exchanged current density in A/m2, T is the temperature in K, R = 8.314 J/(K·mol) is the universal gas constant, αc and αa are the cathodic and anodic charge transfer coefficient, respectively.

The amount of deposition is calculated using Faraday's laws of electrolysis [40], which states that the amount of silicon deposited on the working electrode is proportional to the amount of electricity used
$dloc=iloc⋅t⋅M/(z⋅F⋅ρ)$
(6)

where dloc is the coating thickness, t is the plating time in s, M =28 g/mol, and ρ = 2300 kg/m3 are the atomic weight and density of silicon, respectively.

Thickness variation of the deposited layer k is defined as
$k=[(davg−dmin)/davg]×100%$
(7)

where dmin and davg are the minimum and average coating thickness, respectively.

Model Implementation.

The silicon electrodeposition process is solved in 2D continuum space using a tertiary current distribution module. Figure 1 illustrates the schematic of the model geometry, containing the finite element mesh. The deposition at the working electrode (cathode) and the electricity supply at the counter electrode (anode) are set to take place at 100% current yield. The process is assumed to be time dependent, where the silicon electrodeposition is simulated as the moving boundary of the cathode. The model is governed by the mass conservation for Si4+ and Cl ions and the electroneutrality condition. The other boundaries are assumed to be insulated.

Fig. 1
Fig. 1
Close modal
As the baseline scenario, the distance between the two electrodes d is set to be 5 mm. The area of the anode is fixed as 12 mm; and the area ratio of cathode over anode r is 1 as the baseline condition. The potential of the anode V is 0 V, whereas the potential of the cathode is fixed at −1.7 V. The temperature T is 293 K. The initial concentration of Si4+ in the electrolyte c is 200 mol/m3. The diffusion coefficients Di of Si4+ and Cl ions are assumed to obey the Arrhenius type relationship (Eq. (8)) [41]. Table 1 summarizes the baseline condition, as well as the design space of five operation parameters as the primary sources of uncertainty. Table 2 lists the other electrochemical parameters used in the model
$Di=D0⋅ exp(−ED/RT)$
(8)

where D0 = 4.13 × 10−10 is a pre-exponential factor, and ED = 19,200 J/mol is the activation energy for diffusion.

Table 1

Baseline condition and design space of the operation parameters examined as the sources of uncertainty

Operation parametersBaselineRange
Distance between the electrodes d (mm)5.05.0–15.0
Area ratio of cathode over anode r1.00.2–1.8
Potential of the anode V (V)0.00.0–2.0
Temperature T (K)293.0293.0–333.0
Initial concentration of Si4+c (mol/m3)200.0120.0–240.0
Operation parametersBaselineRange
Distance between the electrodes d (mm)5.05.0–15.0
Area ratio of cathode over anode r1.00.2–1.8
Potential of the anode V (V)0.00.0–2.0
Temperature T (K)293.0293.0–333.0
Initial concentration of Si4+c (mol/m3)200.0120.0–240.0
Table 2

Electrochemical parameters used in the model

ParameterValue
Potential of the cathode (V)−1.7
Equilibrium potential at the silicon surface (V)0.0
Cathodic charge transfer coefficient αc0.5
Anodic charge transfer coefficient αa0.5
Valence of Si4+4.0
Valence of Cl−1.0
Electrical conductivity of silicon1000 S/m
ParameterValue
Potential of the cathode (V)−1.7
Equilibrium potential at the silicon surface (V)0.0
Cathodic charge transfer coefficient αc0.5
Anodic charge transfer coefficient αa0.5
Valence of Si4+4.0
Valence of Cl−1.0
Electrical conductivity of silicon1000 S/m

The mesh size of the cathode boundary is set to be 0.05 mm to precisely predict the growth of the deposited silicon layer. The electrolyte domain is meshed using free triangular elements. Complete geometry consists of 7041 domain elements and 373 boundary elements.

Gaussian Process-Based Surrogate Model

To further improve the efficiency to explore the operation parameter space without losing accuracy, a GP-based surrogate model is constructed. GP surrogate model is a statistical procedure that generates an estimated surface from a scattered set of points via a covariance governed Gaussian process interpolation method. It requires training samples for model construction and predicting responses of new sample points [42].

Assume the true performance function of thickness variation G(x) has a form of
$G(x)≈GGP(x)=fT(x)⋅α+S(x)$
(9)
where $fT(x)=[f1(x), …,fb(x)]$ is a basis function, $α=[α1, …, αb]$ is a regression coefficient vector , and S(x) is a Gaussian stochastic process with zero mean and certain covariance matrix. The covariance function between two inputs xi and xj is expressed as the squared-exponential form
$Cov(i,j)=σ2R(xi,xj)=σ2⋅ exp [−∑p=1kap|xip−xjp|bp]$
(10)
wher, ap and bp are parameters of the GP model, and k is the number of the input parameters x. R(xi, xj) is an n × n symmetric correlation matrix and is defined based on the distance between two sample points xi and xj. It is worth noting that, other forms of the covariation function, like the Mercer kernel, can be utilized to address more complex problems [43]. The GP model is then trained with the imported simulation results (XS, GS) to generate the initial thickness variation function G(x). For any new point $x′$, its performance can be estimated via maximizing the likelihood function [44] as a normal distributed variable with mean $Ĝ(x′)$ and variance $e⌢(x′)$
$Ĝ(x′)=μ+rTR−1(G−Aμ)$
(11)
$e⌢(x′)=σ2[1−rTR−1r+(1−ATR−1r)2/ATR−1A]$
(12)

where r is the correlation vector between $x′$ and the input parameters x, A is an n × 1 unit vector.

Influenced by the sampling situations of the data pool, the initial GP model may potentially have low fidelity in some regions of the electrodeposition operation parameters. Hence, the maximum expected information value-based adaptive sampling method is further used to improve the model fidelity [44]. Generally, three consecutive steps are involved in the approach: (1) an initial set of sample points (XS, GS) are generated by the FE model to develop the initial GP model, where XS is the model input parameters, i.e., the five operation parameters in the study, and GS is the corresponding output performance, i.e., the thickness variation. In this study, a uniform grid sampling strategy is used to select the input XS. (2) In the next step, the simulation results (XS, GS) are stored in a sample data pool and used to construct the predictive models M employing the GP-based surrogate modeling technique [4446]. (3) The most important sample points x* are then selected using the maximum expected information value-based adaptive sampling technique. The performances of these points G* are evaluated with the FE model, and the new sample points (x*, G*) are imported into the sample data pool to upgrade the surrogate model. (4) Cumulative confidence level (CCL) is used to qualify the accuracy of reliability of the surrogate model via Eqs. (13) and (14) [44]
$CCL(M,Xm)=∑i=1NCL(xm,i)N$
(13)
$CL(xm,i)=Φ(|Ĝ(xm,i)|e⌢(xm,i))$
(14)

where CL(xm,i) represents the confidence level of the classification for sample xm,i, Φ is a standard normal cumulative distribution function, and |.| is the absolute operator. Note that CCL(M, Xm) is always a positive value within [0.5, 1]. Eventually, the updated GP model acquires enhanced fidelity, which can be used for the thickness variation prediction and optimization design.

Results and Discussion

In this section, the simulation results from the FE model are first discussed. Then, the design sensitivity analysis is applied to the five operation parameters to evaluate their influences. Afterward, the GP-based surrogate model results are introduced; and the uncertainty quantification study results of the silicon electrodeposition process are discussed. Finally, a design improvement strategy for the experimental condition is provided based upon the surrogate model.

Finite Element Model Results.

Figure 2 illustrates the thickness change of the electrodeposited silicon layer along the cathode surface from the developed FE model when the simulation runs from 0 to 600 s. As can be seen, a layer of silicon is steadily growing on the cathode substrate, whereas a rather uneven coating is shown: the deposits at the left (<2 mm) and right (>18 mm) edges of the cathode substrate are much thicker than those in the center. For instance, after 600 s deposition, the center area is about 12.7 mm thick, and remains relatively consistent and uniform; however, the two edge regions can reach up to 56.6 mm. Consequently, the thickness variation of the deposition process is calculated to be 11.1% via Eq. (7). The nonuniformity phenomenon can be explained by the edge effect of electrode current density [47]. As illustrated in the inset of Fig. 2, there are concentrated electrode and electrolyte current densities near the edges, leading to high electrochemical reaction rates and thick coatings; on the contrary, the electric field and current density flow in the center area keep consistent, resulting in the uniform coating. It is worth noting that, the nonuniformity is not getting exaggerated with the increase of simulation time. This may be because the electrical conductivity of silicon is low compared to the metal materials; as the silicon layer becomes thicker near the edges, it becomes more and more difficult to transfer charges from the silicon–electrolyte interface to the cathode substrate, which gradually reduces the electrode current density and limits the electrodeposition reaction rates. In the following discussions, the thickness variation of the coating layer is treated as the target property to investigate the influences of different operation parameters and to conduct the uncertainty quantification analysis on the silicon deposition process.

Fig. 2
Fig. 2
Close modal

Sensitivity Analysis.

Before the construction of the surrogate model, a sensitivity analysis is conducted on the silicon deposition process to evaluate the influences of different operation parameters on the thickness variation of coated Si layer, to identify the critical parameters with high sensitivity, and to minimize the discrepancy between simulation results and surrogate model.

Figure 3 shows the effects of five operation parameters, including the distance between the two electrodes d, area ratio of cathode over anode r, the potential of the anode (counter electrode) V, temperature T, and initial concentration of Si4+c, on the silicon layer thickness variation k after 600 s simulation. It can be observed that five different operation parameters show distinct impacts.

Fig. 3
Fig. 3
Close modal

First, the thickness variation k has a nonlinear relation with the distance between two electrodes d, as shown in Fig. 3(a). When d =5 mm, k is rather small, only around 11.1%; with the increase of d, k increases simultaneously; when d 20 mm, however, further extending d will no longer affect k, which becomes a constant, 23.9%. This can be explained as follows. While the distance between the working and counter electrodes is small, the electric field in the electrolyte is strong so that the charge transfer and ion transportation between two electrodes are constrained and keep uniform, specifically near the cathode substrate; therefore, the current density concentration and sharp spikes of deposited silicon near the electrode edges are limited. On the contrary, extending the distance between electrodes could lower the electric field strength, release the constraints on ion movements, and enlarge the gradient of electrolyte potential field near the electrode edges, leading to an enhanced edge effect and a high thickness variation. Nevertheless, when the electrodes are separated more than 20 mm, the effects of the reduction of the electric field strength become minimal, thus the uneven distribution of the electrode current density will no longer be exaggerated; as a result, k becomes independent of d.

Second, the area ratio between two electrodes r shows a significant impact on the thickness variation k, as illustrated in Fig. 3(b). On one hand, when the area of cathode is smaller than that of the anode, i.e., r <1, the edge effect is severe during the electrodeposition process, causing large electrolyte potential field gradients (inset of Fig. 3(b)) and thicker silicon coatings around the edges of cathode substrate. The edge effect can be greatly mitigated when the cathode has a slightly larger area (i.e., r 1.25); consequently, k is minimized, approaching 9.8%. On the other hand, when the cathode is oversized compared to the anode (r >1.5), the electrolyte potential gradient outside the anode range (circled in the inset of Fig. 3(b)) is rather small, suggesting that the driven force for the silicon electrodeposition reaction is largely reduced; whereas, the electrolyte potential gradient in the center region remains high. Eventually, different from the previous scenarios presented in Fig. 2, the coating layer has a thicker center but thinner edges; furthermore, due to the large variation of the potential gradient between the center and edge regions, the thickness variation k of the coating layer is significantly increased.

Third, the potential of the anode V and initial concentration of Si4+c have little influences on the thickness variation k, as shown in Figs. 3(c) and 3(e), respectively, where k is fluctuating within the examined ranges and shows nonlinear relationships with V and c. This is because, even though higher V or c could lead to the improved driven force for the electrochemical reaction and correspondingly higher silicon coating rate, they will barely change the electrolyte potential gradient pattern, and thus the thickness variation is only merely affected.

Lastly, as presented in Fig. 3(d), the thickness variation k increases along with the rising temperature T in a nearly linear manner. This is caused by two reasons: first, as introduced in Eq. (8), the diffusion coefficients of different ions will increase at elevated temperature, indicating easier transportation of these species in the electrolyte; in addition, higher T could enhance the electrochemical reaction rate and increase the current density at the electrode–electrolyte interface (Eq. (5)). Both these factors will exaggerate the edge effect, resulting in the nonuniformity of the electrodeposited products.

Furthermore, the finite difference method is employed to quantify the sensitivity of the thickness variation k with respect to the five operation parameters. To avoid scaling issues, the differences are normalized to the values used in the baseline design [48]. Equation (15) is used to calculate the normalized sensitivity σ of the thickness variation k
$σ=⟨[ki(qi+Δqi)−ki(qi)]/k¯iΔqi/q¯i⟩$
(15)

where $q¯i$ is the value of the ith operation parameter in the baseline design, $k¯i$ is the calculated thickness variation with respect to the ith operation parameter, and $Δqi$ is the difference of the parameter. The sensitivity results are listed in Table 3. The results reveal that the silicon electrodeposition process is sensitive to the geometric setting of the experiments, including the distance and area ratio between two electrodes. Meanwhile, the other parameters, i.e., the potential of the anode, initial Si4+ concentration, and the temperature of the electrolyte bath, show relatively low influences. In the following, these five operation parameters are used to construct the GP-based surrogate model.

Table 3

The normalized sensitivity of the thickness variation k with respect to the five operation parameters

Operation parametersBaseline valueSensitivity
d5.0 mm0.463
r1.01.611
V0.0 V0.003
T293.0 K0.015
c200.0 mol/m30.001
Operation parametersBaseline valueSensitivity
d5.0 mm0.463
r1.01.611
V0.0 V0.003
T293.0 K0.015
c200.0 mol/m30.001

Surrogate Model Prediction Results.

The GP-based surrogate model is constructed based upon the FE simulation and the sensitivity analysis results using the uniform grid sampling within the design space of the operation parameters. Particularly, from the sensitivity analysis, two critical variables are identified that have big impacts on the uncertainty of the system performance. Afterward, more densely arranged samples for them within the examined ranges are chosen to calculate thickness variation k through the FE simulations, to model the complex nonlinear relationship between the inputs and outputs. Therefore, the sensitivity analysis can support the construction of the surrogate model and reduce the error. To have better visualization and explanation of the surrogate model, two predicted performance surface figures are presented as examples. Figures 4(a) and 4(b) illustrate the performance surfaces of silicon electrodeposition thickness variation k with respect to two selected combinations of operation parameters, i.e., distance d and area ratio r between two electrodes, as well as the potential of anode V and temperature of the electrolyte bath T, respectively. The scatters represent the FE simulation results. A highly nonlinear performance surface of k in relation to the geometric setting of the experiments, i.e., d and r, is observed in Fig. 4(a). First, when r is less than 0.75, k remains almost constant at a high level, approximately 25%, and becomes independent of d. Besides, with the variation of d, r has distinct influences: when the electrodes are close to each other, the valley of the performance surface lies at r = ∼1.25, and a sharp peak is shown at the corner of d =0 mm and r =2; however, when two electrodes are far apart, the performance surface valley is found at r = ∼1.75, and k will monotonically increase with the decrease of r. In addition, as shown in Fig. 4(b), k is fluctuating within the range of 11.5%–12% in the examined parameter space and is insensitive to the change of T and V. By reducing T, k can be slightly lowered.

Fig. 4
Fig. 4
Close modal

The results from the GP surrogate model are in consistent with FE model predictions and sensitivity analysis results, suggesting the high fidelity of the surrogate model. In the next step, the GP surrogate model is used for the uncertainty quantification analysis of the silicon electrodeposition process.

Uncertainty Quantification Analysis.

In this study, the Monte Carlo simulation technique is adopted to perform the UQ study using the GP surrogate model. This approach provides several merits [49]. First, the Monte Carlo technique could directly map the input parameter dataset to the output performance space, and acquire the statistical distribution of the output given the distribution of the inputs, which is the most reliable approach to handle the complex nonlinear relationship between input and output [50]. Second, the GP surrogate model is used as a mapping function instead of the computationally expensive FE simulations, which could greatly reduce the computational burden without losing much accuracy when dealing with a large amount of data points.

The five input operation parameters in Table 1 are assumed to follow the truncated normal distribution [51]. Figure 5 shows the histograms of two input parameters, i.e., distance d and area ratio r between two electrodes, as examples. Correspondingly, 105 data sets are formed using the random number generator, each consisting of the aforementioned five operation parameters, where the combination of the value for each parameter has no preference and is totally random. Afterward, the distribution of the output, i.e., thickness variation k, is obtained by feeding the input data sets into the GP surrogate model. Figure 6 illustrates a histogram of the estimated output from the GP surrogate model given the generated input data sets. The output distribution is not Gaussian due to the nonlinear relationship between the input parameters and output performance. The average value and standard deviation are 0.185 and 0.030, respectively.

Fig. 5
Fig. 5
Close modal
Fig. 6
Fig. 6
Close modal

The probability distribution of the silicon electrodeposition thickness variation k is helpful to understand the uncertainty and variability of the electrodeposition process and could provide a guideline for the design improvement strategies to reduce the noise and enhance the stability. Figures 7(a) and 7(b) present the histogram of the silicon electrodeposition thickness variation k with respect to the change of coefficient of variations of the distance between the electrodes d and initial concentration of Si4+c, respectively, as examples. As can be seen, these two variables show distinct influences on the uncertainty of the electrodeposition process. On one hand, with the decrease of coefficient of variation of d from 14.4% to 3.6%, the average, and deviation of thickness variation can be decreased greatly. On the other hand, by reducing the coefficient of variation of c, the uncertainty of the system is not much affected. These indicate that different input variables have various contributions to the uncertainty of the Si anode; as a result, in order to acquire a better design of the silicon electrodeposition process with more stable performances and less uncertainty, the most effective approaches are to control the variations of distance d and area ratio r between two electrodes.

Fig. 7
Fig. 7
Close modal

Design Improvement Strategy.

In addition to the uncertainty quantification study, which can quantify the uncertainty of the silicon electrodeposition process, the surrogate model could offer another critical advantage: to guide the design optimization of the process to minimize the thickness variation k and enhance the reliability performances. By scanning the surrogate model results within the examined parameter space, the optimal operation parameters to electroplate silicon can be identified, which could lead to the smallest k. Table 4 compares k between the predicted optimal operation parameters and the baseline condition. It can be observed that, compared to the baseline condition, the optimal condition is updated, where r, V, and c are hardly changed, d is slightly extended and T is modified to be 322.3 K; as a result, the k is decreased from 11.1% to 9.97%, suggesting that the nonuniformity phenomenon can be effectively mitigated so that the silicon deposition process is upgraded. Moreover, the optimal condition is imported into the FE model to validate the GP surrogate model estimation. The thickness variation k from the FE model is almost identical to the predicted value from the surrogate model, further demonstrating the accuracy of the GP-based surrogate model.

Table 4

Comparison of the input operation parameters and output thickness variation k between the optimal and baseline designs

Optimal conditionBaseline
Input variablesDistance between the electrodes d (mm)5.55.0
Area ratio of cathode over anode r1.11.0
Potential of the anode V (V)0.20.0
Temperature T (K)322.3293.0
Initial concentration of Si4+c (mol/m3)200.7200.0
Thickness variation k from surrogate model9.97%11.10%
Thickness variation k from FE simulation9.96%11.10%
Optimal conditionBaseline
Input variablesDistance between the electrodes d (mm)5.55.0
Area ratio of cathode over anode r1.11.0
Potential of the anode V (V)0.20.0
Temperature T (K)322.3293.0
Initial concentration of Si4+c (mol/m3)200.7200.0
Thickness variation k from surrogate model9.97%11.10%
Thickness variation k from FE simulation9.96%11.10%

Conclusion

In this paper, we conduct the UQ analysis on the silicon electrodeposition process to explore the influences of various sources of uncertainty on the quality of the electrodeposited products and to find the optimal experimental conditions providing minimized thickness variation of the coated silicon layer. Five experimental operation parameters are chosen as the primary sources of uncertainty and investigated, including the distance d and area ratio r between the working (cathode) and counter (anode) electrodes, the potential of the counter electrode V, temperature of the electrolyte bath T, and initial concentration of Si4+ in the electrolyte c. The overall modeling procedures contain two steps. In the first step, a 2D multiphysics-based FE model is developed to simulate the electrodeposition process with different combinations of operation parameters as input; the thickness variation of the electroplated silicon layer is calculated; the design sensitivity analyses are then conducted to identify the critical parameters with high sensitivity. In the second step, the GP-based surrogate model is constructed with the FE simulation results as training data; the propagation of uncertainty from inputs over the output space is analyzed, and the optimal experimental condition is proposed to minimize the thickness variations of the coated silicon layer.

It is found that a rather uneven coating is obtained, where the deposits near the left and right edges of the working electrode are much thicker than those in the center. This is mainly due to the edge effect, i.e., the concentrated electrode and electrolyte current densities near the electrode edges during the reaction. Besides, the electrodeposition process is rather sensitive to the geometric settings of the experiments, i.e., distance and area ratio between the two electrodes; both extending the distance between two electrodes and using oversized or undersized working electrode compared to the counter electrode can significantly exaggerate the thickness variation of the electrodeposited product. On the contrary, other conditions, such as the potential of the counter electrode, temperature, and ion concentration in the electrolyte bath, show less impacts. In addition, the developed GP surrogate model can efficiently and accurately estimate the performance of the electrodeposition given the experimental operation parameters. Different operation parameters are found to have distinct contributions to the uncertainty of the system; by assuming that the parameters follow a normal distribution, the average thickness variation and standard deviation are 0.185 and 0.030, respectively. Furthermore, the optimal operating condition is proposed by the GP model and validated using FE simulation to minimize the thickness variation of the coated silicon layer and to enhance the reliability of the electrodeposition experiment. To narrow the distance between two electrodes, to use a slightly larger working electrode, and elevating the temperature of the electrolyte bath are key to reduce the thickness variation.

Funding Data

• Office of Naval Research (ONR) through the Defense University Research-to-Adoption (DURA) Initiative (Grant No. N00014-18-S-F004; Funder ID: 10.13039/100000006).

• National Science Foundation (NSF) Engineering Research Center for Power Optimization of Electro-Thermal Systems (POETS) with Cooperative Agreement (Grant No. EEC-1449548; Funder ID: 10.13039/100000001).

References

1.
Green
,
M. A.
,
2004
, “
Recent Developments in Photovoltaics
,”
Sol. Energy
,
76
(
1–3
), pp.
3
8
.10.1016/S0038-092X(03)00065-3
2.
El Abedin
,
S. Z.
,
Borissenko
,
N.
, and
Endres
,
F.
,
2004
, “
Electrodeposition of Nanoscale Silicon in a Room Temperature Ionic Liquid
,”
Electrochem. Commun.
,
6
(
5
), pp.
510
514
.10.1016/j.elecom.2004.03.013
3.
Cui
,
Y.
, and
Lieber
,
C. M.
,
2001
, “
Functional Nanoscale Electronic Devices Assembled Using Silicon Nanowire Building Blocks
,”
Science
,
291
(
5505
), pp.
851
853
.10.1126/science.291.5505.851
4.
Zhao
,
X. W.
,
Schoenfeld
,
O.
,
Kusano
,
J.
,
Aoyagi
,
Y.
, and
Sugano
,
T.
,
1994
, “
Violet and Blue-Light Emissions From Nanocrystalline Silicon Thin-Films
,”
Jpn. J. Appl. Phys. Part 2-Lett.
,
33
(
Part 2, No. 5A
), pp.
L649
L651
.10.1143/JJAP.33.L649
5.
Tong
,
S.
,
Liu
,
X. N.
, and
Bao
,
X. M.
,
1995
, “
Study of Photoluminescence in Nanocrystalline Silicon Amorphous-Silicon Multilayers
,”
Appl. Phys. Lett.
,
66
(
4
), pp.
469
471
.10.1063/1.114059
6.
Vijayalakshmi
,
S.
,
George
,
M. A.
,
Sturmann
,
J.
, and
Grebel
,
H.
,
1998
, “
Pulsed-Laser Deposition of Si Nanoclusters
,”
Appl. Surf. Sci.
,
127-129
, pp.
378
382
.10.1016/S0169-4332(97)00659-4
7.
Xie
,
X. Y.
,
Wan
,
Q.
,
Liu
,
W. L.
,
Men
,
C. L.
,
Lin
,
Q.
, and
Lin
,
C. L.
,
2003
, “
Nanoscale Silicon Prepared on Different Substrates Using Electron-Beam Evaporation and Their Field-Emission Property
,”
Appl. Surf. Sci.
,
217
(
1–4
), pp.
39
42
.10.1016/S0169-4332(03)00561-0
8.
Chen
,
X. L.
,
Gerasopoulos
,
K.
,
Guo
,
J. C.
,
Brown
,
A.
,
Wang
,
C. S.
,
Ghodssi
,
R.
, and
Culver
,
J. N.
,
2011
, “
A Patterned 3D Silicon Anode Fabricated by Electrodeposition on a Virus-Structured Current Collector
,”
,
21
(
2
), pp.
380
387
9.
Lupan
,
O.
,
Pauporte
,
T.
, and
Viana
,
B.
,
2010
, “
Low-Temperature Growth of ZnO Nanowire Arrays on p-Silicon (111) for Visible-Light-Emitting Diode Fabrication
,”
J. Phys. Chem. C
,
114
(
35
), pp.
14781
14785
.10.1021/jp104684m
10.
Pauporté
,
T.
,
Qi
,
S.
, and
Viana
,
B.
,
2018
, “
Low Temperature Electrodeposition of Silicon Layers
,”
Oxide-Based Materials and Devices IX
,
International Society for Optics and Photonics
, 10533, p. 105332S.10.1117/12.2294944
11.
Colletti
,
L. P.
,
Flowers
,
B. H.
, and
Stickney
,
J. L.
,
1998
, “
Formation of Thin Films of CdTe, CdSe, and CdS by Electrochemical Atomic Layer Epitaxy
,”
J. Electrochem. Soc.
,
145
(
5
), pp.
1442
1449
.10.1149/1.1838502
12.
Pauporte
,
T.
,
Yoshida
,
T.
,
Goux
,
A.
, and
Lincot
,
D.
,
2002
, “
One-Step Electrodeposition of ZnO/Eosin Y Hybrid Films From a Hydrogen Peroxide Oxygen Precursor
,”
J. Electroanal. Chem.
,
534
(
1
), pp.
55
64
.10.1016/S0022-0728(02)01105-1
13.
Pandey
,
R. K.
,
Sahu
,
S.
, and
Chandra
,
S.
,
1996
,
Handbook of Semiconductor Electrodeposition
, CRC Press, Boca Raton, FL.
10.1201/9781315213989
14.
Munisamy
,
T.
, and
Bard
,
A. J.
,
2010
, “
Electrodeposition of Si From Organic Solvents and Studies Related to Initial Stages of Si Growth
,”
Electrochim. Acta
,
55
(
11
), pp.
3797
3803
.10.1016/j.electacta.2010.01.097
15.
Katayama
,
Y.
,
Yokomizo
,
M.
,
Miura
,
T.
, and
Kishi
,
T.
,
2001
, “
Preparation of a Novel Fluorosilicate Salt for Electrodeposition of Silicon at Low Temperature
,”
Electrochemistry
,
69
(
11
), pp.
834
836
.10.5796/electrochemistry.69.834
16.
Wahab
,
H. A.
,
Noordin
,
M. Y.
,
Izman
,
S.
, and
Kurniawan
,
D.
,
2013
, “
Quantitative Analysis of Electroplated Nickel Coating on Hard Metal
,”
Sci. World J.
,
2013
, pp.
1
6
.10.1155/2013/631936
17.
Broughton
,
J.
, and
Brett
,
M.
,
2005
, “
Variations in MnO2 Electrodeposition for Electrochemical Capacitors
,”
Electrochim. Acta
,
50
(
24
), pp.
4814
4819
.10.1016/j.electacta.2005.03.006
18.
Illy
,
B. N.
,
Cruickshank
,
A. C.
,
Schumann
,
S.
,
Da Campo
,
R.
,
Jones
,
T. S.
,
Heutz
,
S.
,
McLachlan
,
M. A.
,
McComb
,
D. W.
,
Riley
,
D. J.
, and
Ryan
,
M. P.
,
2011
, “
Electrodeposition of ZnO Layers for Photovoltaic Applications: Controlling Film Thickness and Orientation
,”
J. Mater. Chem.
,
21
(
34
), pp.
12949
12957
.10.1039/c1jm11225b
19.
Matlosz
,
M.
,
Vallotton
,
P. H.
,
West
,
A.
, and
Landolt
,
D.
,
1992
, “
Nonuniform Current Distribution and Thickness During Electrodeposition Onto Resistive Substrates
,”
J. Electrochem. Soc.
,
139
(
3
), pp.
752
761
.10.1149/1.2069297
20.
Solovjev
,
D.
,
Solovjeva
,
I.
, and
Litovka
,
Y. V.
,
2020
, “
Synthesis of the Suboptimal Control Algorithm for Electroplating Processes Under Conditions of Uncertainty in the Range of Processed Products
,”
IOP Conf. Series Mater. Sci. Eng.
,
709
(
2
),
p. 022063
. 10.1088/1757-899X/709/2/022063
21.
Bai
,
G. X.
, and
Wang
,
P. F.
,
2016
, “
An Internal State Variable Mapping Approach for Li-Plating Diagnosis
,”
J. Power Sources
,
323
, pp.
115
124
.10.1016/j.jpowsour.2016.05.040
22.
Zheng
,
Z.
,
Chen
,
B.
,
Gurumukhi
,
Y.
,
Cook
,
J.
,
Ates
,
M. N.
,
Miljkovic
,
N.
,
Braun
,
P. V.
, and
Wang
,
P.
,
2019
, “
Surrogate Model Assisted Design of Silicon Anode Considering Lithiation Induced Stresses
,”
IEEE International Reliability Physics Symposium (IRPS)
, Monterey, CA, Mar. 31–Apr. 4, pp.
1
6
.10.1109/IRPS.2019.8720601
23.
Zheng
,
Z. Y.
,
Chen
,
B.
,
Fritz
,
N.
,
Gurumukhi
,
Y.
,
Cook
,
J.
,
Ates
,
M. N.
,
Miljkovic
,
N.
,
Braun
,
P. V.
, and
Wang
,
P. F.
,
2019
, “
Lithiation Induced Stress Concentration for 3D Metal Scaffold Structured Silicon Anodes
,”
J. Electrochem. Soc.
,
166
(
10
), pp.
A2083
A2090
.10.1149/2.1031910jes
24.
Zheng
,
Z.
,
Chen
,
B.
,
Fritz
,
N.
,
Gurumukhi
,
Y.
,
Cook
,
J.
,
Ates
,
M. N.
,
Miljkovic
,
N.
,
Braun
,
P. V.
, and
Wang
,
P.
,
2020
, “
The Impact of Non-Uniform Metal Scaffolds on the Performance of 3D Structured Silicon Anodes
,”
J. Energy Storage
,
30
, p.
101502
.10.1016/j.est.2020.101502
25.
Zheng
,
Z.
,
Chen
,
B.
,
Xu
,
Y.
,
Fritz
,
N.
,
Gurumukhi
,
Y.
,
Cook
,
J.
,
Ates
,
M. N.
,
Miljkovic
,
N.
,
Braun
,
P. V.
, and
Wang
,
P.
,
2021
, “
A Gaussian Process Based Crack Pattern Modeling Approach for Battery Anode Materials Design
,”
J. Electrochem. Energy Convers. Storage
,
18
(
1
), p.
011011
.10.1115/1.4046938
26.
Chen
,
K. S.
, and
Evans
,
G. H.
,
2004
, “
Two-Dimensional Modeling of Nickel Electrodeposition in LIGA Microfabrication
,”
Microsyst. Technol.-Micro- Nanosyst.-Inf. Storage Process. Syst.
,
10
(
6–7
), pp.
444
450
.10.1007/s00542-004-0373-8
27.
Lupo
,
C.
, and
Schlettwein
,
D.
,
2019
, “
Modeling of Dendrite Formation as a Consequence of Diffusion-Limited Electrodeposition
,”
J. Electrochem. Soc.
,
166
(
1
), pp.
D3182
D3189
.10.1149/2.0231901jes
28.
Wheeler
,
D.
,
Josell
,
D.
, and
Moffat
,
T. P.
,
2003
, “
Modeling Superconformal Electrodeposition Using the Level Set Method
,”
J. Electrochem. Soc.
,
150
(
5
), pp.
C302
C310
.10.1149/1.1562598
29.
Ng
,
S. S. Y.
,
Xing
,
Y.
, and
Tsui
,
K. L.
,
2014
, “
A Naive Bayes Model for Robust Remaining Useful Life Prediction of Lithium-Ion Battery
,”
Appl. Energy
,
118
, pp.
114
123
.10.1016/j.apenergy.2013.12.020
30.
Sacks, J., Welch, W. J., Mitchell, T. J., and Wynn, H. P., 1989, “Design and Analysis of Computer Experiments,”
Statist. Sci.
, 4(4), pp. 409–423.
10.1214/ss/1177012413
31.
Wang
,
L. P.
,
Beeson
,
D.
,
Akkaram
,
S.
, and
Wiggs
,
G.
,
2005
, “
Gaussian Process Meta-Models for Efficient Probabilistic Design in Complex Engineering Design Spaces
,”
Proceedings of the ASME 2005 International Design Engineering Technical Conferences and Computers and Information in Engineering Conference
, Long Beach, CA, Sept. 24–28, pp. 785–798.10.1115/DETC2005-85406
32.
Jensen
,
W. A.
,
2017
, “
Response Surface Methodology: Process and Product Optimization Using Designed Experiments 4th Edition
,”
J. Qual. Technol.
,
49
(
2
), pp.
186
187
.10.1080/00224065.2017.11917988
33.
Haykin
,
S.
, and
Network
,
N.
,
2004
, Neural Networks: A Comprehensive Foundation, Prentice Hall, Upper Saddle River, NJ.
34.
Lophaven
,
S. N.
,
Nielsen
,
H. B.
,
Sondergaard
,
J.
, and
Dace
,
A.
,
2002
,
A Matlab Kriging Toolbox
,
Technical University of Denmark
,
Kongens Lyngbyand
, Report No. IMMTR-2002.
35.
Forrester
,
A.
,
Sobester
,
A.
, and
Keane
,
A.
,
2008
,
Engineering Design Via Surrogate Modelling: A Practical Guide
,
Wiley
, Hoboken, NJ.10.1002/9780470770801
36.
Fan
,
X.
,
Wang
,
P.
, and
Hao
,
F.
,
2019
, “
Reliability-Based Design Optimization of Crane Bridges Using Kriging-Based Surrogate Models
,”
Struct. Multidiscip. Optim.
,
59
(
3
), pp.
993
1005
.10.1007/s00158-018-2183-0
37.
Zheng
,
Z.
,
Xu
,
Y.
,
Chen
,
B.
, and
Wang
,
P.
,
2019
, “
Gaussian Process Based Crack Initiation Modeling for Design of Battery Anode Materials
,”
ASME
Paper No. DETC2019-97547.10.1115/DETC2019-97547
38.
Bird
,
R. B.
,
2002
, “
Transport Phenomena
,”
ASME Appl. Mech. Rev
,
55
(
1
), pp.
R1
R4
.10.1115/1.1424298
39.
Newman
,
J.
, and
Thomas-Alyea
,
K. E.
,
2012
,
Electrochemical Systems
,
Wiley
, Hoboken, NJ.
40.
Kamata
,
M.
, and
Paku
,
M.
,
2007
, “
Exploring Faraday's Law of Electrolysis Using Zinc–Air Batteries With Current Regulative Diodes
,”
J. Chem. Educ.
,
84
(
4
), p.
674
.10.1021/ed084p674
41.
Moats
,
M. S.
,
Hiskey
,
J. B.
, and
Collins
,
D. W.
,
2000
, “
The Effect of Copper, Acid, and Temperature on the Diffusion Coefficient of Cupric Ions in Simulated Electrorefining Electrolytes
,”
Hydrometallurgy
,
56
(
3
), pp.
255
268
.10.1016/S0304-386X(00)00070-0
42.
Rubinstein
,
R. Y.
, and
Kroese
,
D. P.
,
2016
,
Simulation and the Monte Carlo Method
, Vol.
10
,
Wiley
, Hoboken, NJ.10.1002/9781118631980
43.
Melkumyan
,
A.
, and
Ramos
,
F.
,
2009
, “
A Sparse Covariance Function for Exact Gaussian Process Inference in Large Datasets
IJCAI
, Menlo Park, CA, pp.
1936
1942
.10.5555/1661445.1661755
44.
Wang
,
Z. Q.
, and
Wang
,
P. F.
,
2014
, “
A Maximum Confidence Enhancement Based Sequential Sampling Scheme for Simulation-Based Design
,”
ASME J. Mech. Des.
,
136
(
2
), p.
021006
.10.1115/1.4026033
45.
Wang
,
P. F.
,
Wang
,
Z. Q.
, and
Almaktoom
,
A. T.
,
2014
, “
Dynamic Reliability-Based Robust Design Optimization With Time-Variant Probabilistic Constraints
,”
Eng. Optim.
,
46
(
6
), pp.
784
809
.10.1080/0305215X.2013.795561
46.
Wang
,
Z. Q.
, and
Wang
,
P. F.
,
2015
, “
A Double-Loop Adaptive Sampling Approach for Sensitivity-Free Dynamic Reliability Analysis
,”
Reliab. Eng. Syst. Saf.
,
142
, pp.
346
356
.10.1016/j.ress.2015.05.007
47.
Quemper
,
J.-M.
,
Dufour-Gergam
,
E.
,
Frantz-Rodriguez
,
N.
,
Gilles
,
J.-P.
,
Grandchamp
,
J.-P.
, and
Bosseboeuf
,
A.
,
2000
, “
Effects of Direct and Pulse Current on Copper Electrodeposition Through Photoresist Molds
,”
J. Micromech. Microeng.
,
10
(
2
), pp.
116
119
.10.1088/0960-1317/10/2/303
48.
,
J.
,
Mares
,
C.
,
Friswell
,
M.
, and
James
,
S.
,
2000
, “
Selection and Updating of Parameters for an Aluminium Space-Frame Model
,”
Mech. Syst. Signal Process.
,
14
(
6
), pp.
923
944
.10.1006/mssp.2000.1303
49.
Nobari
,
A.
,
Ouyang
,
H. J.
, and
Bannister
,
P.
,
2015
, “
Uncertainty Quantification of Squeal Instability Via Surrogate Modelling
,”
Mech. Syst. Signal Process.
,
60–61
, pp.
887
908
.10.1016/j.ymssp.2015.01.022
50.
Schenk
,
C. A.
,
2005
,
And G.I. Schuëller, Uncertainty Assessment of Large Finite Element Systems
, Vol.
24
,
, Secaucus, NJ.
51.
Lundstedt
,
T.
,
Seifert
,
E.
,
Abramo
,
L.
,
Thelin
,
B.
,
Nyström
,
Å.
,
Pettersen
,
J.
, and
Bergman
,
R.
,
1998
, “
Experimental Design and Optimization
,”
Chemomet. Intell. Lab. Syst.
,
42
(
1–2
), pp.
3
40
.10.1016/S0169-7439(98)00065-3