The Moon is directly exposed to solar radiation and solar wind plasma (drifting protons and electrons) lacking an atmosphere and a global magnetic field. Consequently, the lunar surface is electrically charged by the bombardment of solar wind plasma and emission/collection of photoelectrons. Near the illuminated lunar surface, the plasma sheath is dominated by photoelectrons, thus usually referred to as "photoelectron sheath". Additionally, dust grains on the lunar surface may get charged and levitated from the surface under the influence of the electric field within the plasma sheath as well as gravity. This work is motivated by the high computational cost associated with uncertainty quantification (UQ) analysis of plasma simulations using high-fidelity kinetic models such as particle-in-cell (PIC). The main quantities of interest (QoI) of this study is the vertical structure of the photoelectron sheath and its effects on levitation of dust grains with different sizes and electric charges.
Both the electric potential (𝜙𝜙) and the electric field (𝐸𝐸) on lunar surface are determined by many parameters, such as solar wind drifting velocity (𝑣𝑣 d ), electron temperature (𝑇𝑇 e ), photoelectron temperature (𝑇𝑇 p ), density of ions at infinity (𝑛𝑛 i,∞ ), and density of photoelectrons (𝑛𝑛 p ), etc. Due to uncertain factors in lunar environment, the electric potential, electric field, and the dust levitation height, etc., are also uncertain. While many sources uncertainty may exist, they are generally categorized as either aleatory or epistemic. Uncertainties are characterized as epistemic, if the modeler sees a possibility to reduce them by gathering more data or by refining models. Uncertainties are categorized as aleatory if the modeler does not foresee the possibility of reducing them [1]. An example of the aleatory uncertainty in lunar environment is the solar wind parameters, and an example of the epistemic uncertainty is the photoelectron temperature which is obtained by limited measurement data from Apollo missions. For lunar landing missions, one needs to take into consideration the uncertainties of the electrostatic and dust environment near the lunar surface. For example, the upper and lower bounds of the electric field and dust grain levitation heights in the photoelectron sheath should be considered when determining whether it is safe for a certain area to land a spacecraft.
Determining the bounds of the electric potential, electric field, and dust levitation height, however, is computationally expensive, because the particle-based kinetic models such as particle-in-cell simulations are time-consuming to evaluate. To address this issue, we develop an adaptive Kriging method (AKM) which can determine those bounds with a small number of calculations of the model. It is straightforward to train and obtain an accurate Kriging model [2] to replace the actual model and then calculate the bounds with the model. However, it is not necessary for the Kriging model to be accurate everywhere in its input space, because it will need more training samples and hence decrease the efficiency. Since the objective is to determine those bounds, we only need the Kriging model to be partially accurate near the regions of interest, as long as it can help find those bounds accurately. This way, we can save more computational efforts. To this end, we develop a task-oriented learning function and a stopping criterion to adaptively train the Kriging model. We start with an analytic model for the 1-D photoelectron sheath near the lunar surface [3,4]. This model is computationally cheap and hence the accurate results can be obtained by brute force. With the accurate results, we can test the accuracy of the proposed method. It is noted here that the ultimate application of this method is not the simple 1-D problem presented in this work, but more complicated or computationally expensive models such as 3-D fully kinetic particle-in-cell plasma simulations.
The rest of this paper is organized as follows. Section 2 presents the 1-D photoelectron sheath and dust levitation problem on lunar surface, as well as the 1-D analytic model. Section 3 briefly introduces the Kriging method and general AKM. Section 4 presents the proposed AKM. Section 5 presents the results. Conclusions are given in Section 6.
We employ the recently derived 1-D photoelectron sheath model for the lunar surface [3,4]. As given in detail in [3,4], there are three types of electric potential profiles [3][4][5][6] in the photoelectron sheath: Type A, Type B, and Type C, as shown in Fig. 1, where 𝜙𝜙 is the electric potential and 𝑍𝑍 is the vertical coordinate. In this study, we focus on Type C sheath profile as it is expected at the polar regions of the Moon, where the next lunar landing mission will likely occur.
Both the electrical potential 𝜙𝜙 and corresponding electric field 𝐸𝐸 are functions of 𝑍𝑍 with a series of parameters 𝐏𝐏 = �𝑣𝑣 d , 𝑇𝑇 e , 𝑇𝑇 p , 𝑛𝑛 i,∞ , 𝑛𝑛 p �. To obtain 𝜙𝜙(𝑍𝑍; 𝐏𝐏) and 𝐸𝐸(𝑍𝑍; 𝐏𝐏), we need to solve an ordinary differential equation (ODE) [3]. Once the potential profile 𝜙𝜙 is obtained, it is straightforward to calculate electric field 𝐸𝐸 by
A typical Type C sample curve of 𝐸𝐸(𝑍𝑍; 𝐏𝐏) is shown in Fig. 2. Note that both 𝜙𝜙 and 𝐸𝐸 converge to zero at large values of 𝑍𝑍 where it is used as the electric potential reference (zero potential and zero field).
Subjected to the electric field force, a charged dust on lunar surface may be levitated [7,8]. Above the lunar surface, there is a position where the upward electric field force balances the downward gravity [4]. This position is referred to as equilibrium levitation height, denoted as 𝑍𝑍 * . 𝑍𝑍 * can be solved through the following equation of static equilibrium of a charged dust in an electric field:
where 𝑞𝑞 is the dust charge, 𝑚𝑚 is the mass of the dust, and 𝑚𝑚 = 1.62 m/s 2 is the gravity acceleration on lunar surface [9]. With the assumption of spherical dust grains, 𝑚𝑚 is given by
where 𝑟𝑟 is the radius of the lunar dust grain, and 𝜌𝜌 = 1.8 g/cm 3 is the mass density of dust grains [10].
For simplicity, Eq. ( 2) is rewritten as
where 𝑤𝑤 = 𝑚𝑚𝑚𝑚/𝑞𝑞. Once both 𝐸𝐸(𝑍𝑍; 𝐏𝐏) and 𝑤𝑤 have been given or determined, a root-finding scheme is employed to solve Eq. ( 4) for 𝑍𝑍 * . Fig. 3 shows an example of how to obtain 𝑍𝑍 * graphically.
Due to the lack of information, it is almost impossible to obtain the distribution functions of 𝐏𝐏.
However, this method is too expensive or even unaffordable. One reason is that solving the ODE a large number 𝑁𝑁 MCS of times is time-consuming, even when the analytic solution to the ODE is available for the 1-D problem. Another reason is that there is no analytic solution to complex 2-D or 3-D problems where kinetic particle-in-cell simulations are usually employed to solve the electrostatic field through Poisson's equation.
The objective of this study is to develop a method to determine 𝜙𝜙(𝑍𝑍), 𝜙𝜙(𝑍𝑍), 𝐸𝐸(𝑍𝑍) and 𝐸𝐸(𝑍𝑍) accurately and then calculate 𝑍𝑍 * of dust grains. It is noted here that the ultimate application of this method is not the relatively simple 1-D problem presented in this work, but more complicated or computationally expensive models such as 3-D fully kinetic particle-in-cell plasma simulations. For computationally expensive models, evaluating the model consumes the majority of computational resource, so we will use the number 𝑁𝑁 ODE of ODEs that we need to solve as a measure of the computational cost.
Before presenting the proposed method, we briefly introduce the Kriging model [12,13] and AKM [13][14][15][16][17][18][19][20][21][22][23][24][25][26][27][28], on which the proposed method is based.
Kriging model makes regression to a black-box function (BBF) using a training sample set, or a design of experiment (DoE). The main idea of Kriging is to treat the BBF as a realization of a Gaussian random field indexed by the input variables of the BBF. The theoretical foundation of Kriging model is exactly the Bayesian inference [29] . From the perspective of Bayesian interface, a prior Gaussian random field is trained by the DoE and hence a posterior Gaussian random field is generated. Then the mean value function, also indexed by the input variables of the BBF, of the posterior random field is the Kriging prediction to the BBF. Besides, the variance function, also indexed by the input variables of the BBF, of the posterior random field quantifies the local prediction uncertainty or prediction error.
The randomness, or uncertainty, of the posterior random field mainly comes from the fact that only a limit number of samples of the BBF are used to train the prior random field. In other words, only part of the information of the BBF is available, and the missing part of information leads to the epistemic uncertainty in the random field. Generally, the more training samples we use, the less epistemic uncertainty will result and with stronger confidence will we predict the BBF.
A simple yet widely used prior random field is the stationary Gaussian random field given by
where 𝜇𝜇 is an unknown parameter representing the mean value of the random field 𝐾𝐾(𝐗𝐗) and 𝜂𝜂(𝐗𝐗; 𝜉𝜉 2 , 𝜃𝜃) is a zero-mean stationary Gaussian random field indexed by 𝐗𝐗, the input variables of a BBF 𝑘𝑘(𝐗𝐗). Both the variance parameter 𝜉𝜉 2 and correlation parameters 𝛉𝛉 of 𝜂𝜂(𝐗𝐗; 𝜉𝜉 2 , 𝛉𝛉) are unknown. The parameters 𝜇𝜇, 𝜉𝜉 2 and 𝛉𝛉 fully define the prior random field 𝐾𝐾(𝐗𝐗). A DoE, or a training sample set, of 𝑘𝑘(𝐗𝐗) is used to train 𝐾𝐾(𝐗𝐗) and then determine those parameters.
The correlation function 𝐶𝐶�𝐱𝐱 (𝑖𝑖) , 𝐱𝐱 (𝑗𝑗) � of 𝜂𝜂(𝐗𝐗; 𝜉𝜉 2 , 𝛉𝛉) is given by
where 𝑅𝑅�𝐱𝐱 (𝑖𝑖) , 𝐱𝐱 (𝑗𝑗) ; 𝛉𝛉� is the correlation coefficient function of 𝜂𝜂(𝐗𝐗; 𝜉𝜉 2 , 𝛉𝛉) at two points 𝐱𝐱 (𝑖𝑖) and 𝐱𝐱 (𝑗𝑗) of 𝐗𝐗. There are many models for 𝑅𝑅�𝐱𝐱 (𝑖𝑖) , 𝐱𝐱 (𝑗𝑗) ; 𝛉𝛉�. A widely used model is known as the Gaussian model, or squared exponential model, given by
where 𝐷𝐷 is the dimension of 𝐗𝐗, 𝑥𝑥 𝑑𝑑 (𝑖𝑖) is the 𝑑𝑑 th component of 𝐱𝐱 (𝑖𝑖) , and 𝜃𝜃 𝑑𝑑 is the 𝑑𝑑 th component of 𝛉𝛉.
For a BBF 𝑘𝑘(𝐗𝐗), the Kriging model predicts 𝑘𝑘(𝐱𝐱) as 𝑘𝑘 � (𝐱𝐱), which is a normal variable whose mean value and variance are 𝑘𝑘 � (𝐱𝐱) and 𝜎𝜎 2 (𝐱𝐱), respectively. Note that 𝜎𝜎 2 (𝐱𝐱) is also termed as mean squared error (MSE). Generally, 𝑘𝑘 � (𝐱𝐱) is regarded as the deterministic prediction to 𝑘𝑘(𝐱𝐱), since a deterministic prediction is usually needed. 𝜎𝜎 2 (𝐱𝐱) measures the prediction uncertainty, or prediction error, and therefore we can validate a Kriging model simply using 𝑘𝑘 � (𝐱𝐱) and 𝜎𝜎 2 (𝐱𝐱) without employing traditional validation methods, such as the cross validation [30]. Because of this advantage, many algorithms have been proposed to adaptively train a Kriging model for expensive BBFs [14][15][16][17][18][19][20][21][22][23][24][25][26][27][31][32][33][34][35][36]. When sufficient training samples have been used for training, 𝜎𝜎 2 (𝐱𝐱) converges to 0 and the normal variable 𝑘𝑘 � (𝐱𝐱)
degenerates to a deterministic value, i.e., the exact value of 𝑘𝑘(𝐱𝐱).
The main idea of AKM is to adaptively add training samples to update the Kriging model iteratively until an expected accuracy is achieved. Fig. 6 Given a specific engineering problem, the key of employing an AKM is to make good use of all available information, such as the features of the BBF and QoI, and then design a customized or task-oriented error metric, stopping criterion and learning function.
In the UQ community, a great number of AKMs have been developed to solve varies kinds of problems, such as reliability analysis [15, 17-24, 26, 31-33, 36], robustness analysis [14],
sensitivity analysis [34], robust design [25,35], and reliability-based design [16,27], etc.
In this section, we present detailed procedures of calculating 𝜙𝜙(𝑍𝑍) and 𝜙𝜙(𝑍𝑍). Similar procedures can also apply to calculate 𝐸𝐸(𝑍𝑍) and 𝐸𝐸(𝑍𝑍).
The main idea of the proposed method is to employ the framework of AKM and customize it to calculate 𝜙𝜙(𝑍𝑍) and 𝜙𝜙(𝑍𝑍) (as well as 𝐸𝐸(𝑍𝑍) and 𝐸𝐸(𝑍𝑍)). Fig. 7 shows the brief flowchart of the proposed method. In Step 1, we evenly generate 𝑁𝑁 in initial samples of 𝐏𝐏.
Generally, 𝑁𝑁 in is much smaller than 𝑁𝑁 MCS . Details of this step will be given in Subsection 4.2. In
Step
)
In Step 5, an error metric is developed to measure the error of 𝜙𝜙(𝑍𝑍) and 𝜙𝜙(𝑍𝑍) estimated by Eq. ( 16) and Eq. ( 17).
Step 6 is about a stopping criterion. Details about Steps 5 and 6 will be given in Subsection 4.4. The learning function involved in Step 7 will be given in Subsection 4.3. The implementation of the proposed method will be given in Subsection 4.5.
There are two significant differences between most existing AKMs and the proposed method. First, the former aims at estimating a constant value, such as the structural reliability and robustness, while the latter aims at estimating two functions, i.e., 𝜙𝜙(𝑍𝑍) and 𝜙𝜙(𝑍𝑍). Second, when given a specific value of input, the output of the BBFs involved in the former methods is a single value. However, in this work, with a given realization 𝐩𝐩 of 𝐏𝐏, the output of solving the ODE is a function 𝜙𝜙(𝑍𝑍; 𝐩𝐩). With those differences, we cannot use the existing error metrics, stopping criteria or learning functions. Instead, we take into consideration those differences and design a new error metric, stopping criterion and learning function to fit the problem. This is the main contribution of the proposed algorithm.
For numerical computation, we need to evenly discretize Ω into a few points. Suppose 𝑃𝑃 𝑖𝑖 , the 𝑖𝑖 th component of 𝐏𝐏, is discretized into 𝑁𝑁 𝑖𝑖 points, then Ω will be discretized into in total
points. For convenience, we denote the set of those points as 𝐩𝐩 MCS The 𝑁𝑁 in initial training samples 𝐩𝐩 in of 𝐏𝐏 are selected such that they are distributed in Ω as even as possible. Commonly used sampling methods include random sampling, Latin hypercube sampling and Hammersley sampling [37]. In this study, we employ the Hammersley sampling method because it has better uniformity properties over a multidimensional space [38]. The Hammersley sampling method firstly generates initial training samples in a 5dimensional hypercube [0,1] 5 and then they are mapped into Ω to get the initial training samples of 𝐏𝐏. Note that the five dimensions of the hypercube are assumed to be independent, with the assumption that all variables in 𝐏𝐏 are independent. Those initial training samples, however, are not necessarily among 𝐩𝐩 MCS , so we need to round them to the nearest ones in 𝐩𝐩 MCS . Since the components of 𝐏𝐏 do not necessarily share the same dimension unit, the distances which we use to find the nearest samples should be normalized. For example, the distance 𝑑𝑑 between a sample 𝐩𝐩 (h) generated by Hammersley and a candidate sample 𝐩𝐩 (c) in 𝐩𝐩 MCS is given by
𝑃𝑃 𝑖𝑖,max -𝑃𝑃 𝑖𝑖,min � 2 5 𝑖𝑖=1 (18) where 𝑝𝑝 𝑖𝑖 (ℎ) is the 𝑖𝑖 th component of 𝐩𝐩 (h) , 𝑝𝑝 𝑖𝑖 (𝑐𝑐) is the 𝑖𝑖 th component of 𝐩𝐩 (c) Since for any 𝐩𝐩 ∈ 𝐩𝐩 MCS , it is known that 𝜙𝜙(𝑍𝑍 max ; 𝐩𝐩) ≡ 0 (Fig. 1), theoretically we also need to add all the 𝑁𝑁 𝑃𝑃 points 𝑍𝑍 max × 𝐩𝐩 MCS as input training samples so that we make good use of all known information. However, it is not practical to do so. For example, if 𝑁𝑁 𝑖𝑖 = 10, 𝑖𝑖 = 1,2, … ,5, we need to add 𝑁𝑁 𝑃𝑃 = 10 5 points as input training samples. So many training samples will make 𝜙𝜙 � (𝑍𝑍; 𝐏𝐏) complex, expensive and not accurate, losing its expected properties. To balance the need to add them and the drawback of adding all of them, we add part of them.
Specifically, we evenly generate 𝑁𝑁 𝑃𝑃 ′ samples 𝐩𝐩 ′ of 𝐏𝐏 using procedures similar to what we used to generate 𝐩𝐩 in . Then 𝐱𝐱 inp2 is given by
The input training sample set 𝐱𝐱 inp = 𝐱𝐱 inp1 ⋃ 𝐱𝐱 inp2 . Denote the corresponding electric potential 𝜙𝜙 at 𝐱𝐱 inp as 𝛟𝛟 out . The input-output training sample pairs �𝐱𝐱 inp , 𝛟𝛟 out � are used to build the initial 𝜙𝜙 � (𝑍𝑍; 𝐏𝐏). More training samples will be added to update 𝜙𝜙 � (𝑍𝑍; 𝐏𝐏) later.
Generally, the initial Kriging model is not accurate enough to get 𝜙𝜙(𝑍𝑍) or 𝜙𝜙(𝑍𝑍) accurately through Eq. ( 5) and Eq. ( 6). To improve the accuracy of 𝜙𝜙 � (𝑍𝑍; 𝐏𝐏) and hence of 𝜙𝜙(𝑍𝑍)
and 𝜙𝜙(𝑍𝑍), we need to add training samples of 𝜙𝜙(𝑍𝑍; 𝐏𝐏) to refine 𝜙𝜙 � (𝑍𝑍; 𝐏𝐏). A learning function is used to determine which sample of 𝐏𝐏, and hence of 𝜙𝜙(𝑍𝑍; 𝐏𝐏), should be added.
In our previous work [3], we used the learning function given by
where 𝐩𝐩 (next) is the next to-be-added sample of
To estimate 𝜙𝜙(𝑍𝑍), we also propose a learning function given by
In order to estimate both 𝜙𝜙(𝑍𝑍) and 𝜙𝜙(𝑍𝑍) simultaneously, we combine Eq. ( 23) and Eq. ( 25) to propose a learning function given by
Once 𝐩𝐩 (next) has been determined, we solve the ODE to numerically get a function 𝜙𝜙�𝑍𝑍; 𝐩𝐩 (next) �, from which we get (𝑁𝑁 𝑍𝑍 -1) points (the remaining one at 𝑍𝑍 max , where 𝜙𝜙 ≡ 0, is excluded) and add them into �𝐱𝐱 inp , 𝛟𝛟 out � to enrich the training samples.
Since Eq. ( 16) and Eq. ( 17) cannot obtain absolutely accurate 𝜙𝜙(𝑍𝑍) and 𝜙𝜙(𝑍𝑍) due to the prediction error of 𝜙𝜙 � (𝑍𝑍; 𝐏𝐏), we need an error metric to measure the error of currently estimated 𝜙𝜙(𝑍𝑍) and 𝜙𝜙(𝑍𝑍). Since � 𝜉𝜉(𝑧𝑧,𝐩𝐩)
� measures the absolute expected improvement rate of
� is small for any 𝑧𝑧 ∈ 𝐳𝐳 MCS and 𝐩𝐩 ∈ 𝐩𝐩 MCS , 𝜙𝜙(𝑍𝑍) is expected to sufficiently accurate. Therefore, we propose to use max 𝑧𝑧∈𝐳𝐳 MCS ,𝐩𝐩∈𝐩𝐩 MCS � 𝜉𝜉(𝑧𝑧,𝐩𝐩) 𝜙𝜙(𝑧𝑧) � to quantify the error of 𝜙𝜙(𝑍𝑍).
� is used to quantify the error of 𝜙𝜙(𝑍𝑍). Combining both, we have the error metric 𝛤𝛤, which measures the error of both 𝜙𝜙(𝑍𝑍) and 𝜙𝜙(𝑍𝑍), given by
Once 𝛤𝛤 is small enough, the estimated 𝜙𝜙(𝑍𝑍) and 𝜙𝜙(𝑍𝑍) are expected to be sufficiently accurate. Therefore, the stopping criterion shown in Fig. 7 is defined as
where 𝛾𝛾 is a threshold that controls the efficiency and accuracy of the proposed AKM. Generally speaking, a smaller 𝛾𝛾 will lead to higher accuracy but lower efficiency.
As shown in Fig. 1, 𝜙𝜙(𝑍𝑍; 𝐏𝐏) approaches zero when 𝑍𝑍 takes large value. As a result, 𝜙𝜙(𝑧𝑧)
and 𝜙𝜙(𝑧𝑧) in Eq. ( 26) and Eq. ( 27) are likely to take very small values close to zero. It leads to the singularity of the calculation of Eq. ( 26) and Eq. ( 27), doing harm to the robustness of the proposed algorithm. To solve this issue, we translate all training samples of 𝜙𝜙(𝑍𝑍; 𝐏𝐏) simply by adding a negative constant ϵ. This way, the translated 𝜙𝜙(𝑍𝑍; 𝐏𝐏) will never approach zero and the singularity issue is evitable. Trained by the translated samples of 𝜙𝜙(𝑍𝑍; 𝐏𝐏), the Kriging model 𝜙𝜙 � (𝑍𝑍; 𝐏𝐏) will also lead to the translation of 𝜙𝜙(𝑍𝑍) and 𝜙𝜙(𝑍𝑍). We can translate 𝜙𝜙(𝑍𝑍) and 𝜙𝜙(𝑍𝑍) back simply by subtracting ϵ from them. Note that there is no rigorous theory to quantify how ϵ affect the properties of the proposed AKM. We suggest determining ϵ using
where mean(•) represents mean value.
Based on all the procedures given above, we generate the pseudo codes of the proposed AKM given in Algorithm 1.
Theoretically, it is vital to validate the Kriging model to make sure that it has been trained accurately. An explicit validation, however, is not involved in the proposed AKM. There are two main reasons. First, the adaptive training focuses on the accuracy of QoI instead of the accuracy of the Kriging model. Once there is an indication that the accuracy of QoI in current training iteration is sufficient, i.e., the stopping criterion in Eq. ( 27) is satisfied, the algorithm stops adding more training samples, no matter the Kriging model is globally accurate or not. As a result, when the algorithm has converged, it is very likely that the Kriging model is accurate only on some subdomains but not accurate on other domains. Therefore, it is not suitable to do explicit cross validation. Second, the error metric 𝛤𝛤 can measure the accuracy of QoI, and therefore we in fact do validation implicitly. As long as the accuracy of QoI is sufficient, it does not matter if the Kriging model is or not accurate on the entire domain.
In this section, we illustrate the proposed AKM. MCS is used to solve the same problems with brute force. Results by MCS are treated as standard to verify the proposed AKM. We build the Kriging model and calculate the Kriging predictions using the DACE toolbox [39]. The anisotropic Gaussian kernel is used.
We consider the 1-D photoelectron sheath problem discussed in Section 2. The sun elevation angle is given as 9 degrees. The maximal and minimal values of 𝐏𝐏 = �𝑣𝑣 d , 𝑇𝑇 e , 𝑇𝑇 p , 𝑛𝑛 i,∞ , 𝑛𝑛 p � are given in Table 1. We use both MCS and the proposed AKM to estimate 𝜙𝜙(𝑍𝑍) and 𝜙𝜙(𝑍𝑍). The values of all involved parameters are given in Table 2.
The domain Ω of 𝐏𝐏 is discretized into 𝑁𝑁 𝑃𝑃 = 5 5 points, which are assembled into 𝐩𝐩 MCS .
The 𝑁𝑁 in = 5 samples in hypercube space [0,1] 5 , generated by Hammersley sampling method, are given in Table 3. Then the 5 samples are mapped into Ω, as given in Table 4. Rounding the 5 samples in Ω to the nearest ones in 𝐩𝐩 MCS , we get the initial samples 𝐩𝐩 in of 𝐏𝐏, as given in Table 5. Solving the ODE five times, each with a sample in 𝐩𝐩 in , we get five samples of 𝜙𝜙(𝑍𝑍; 𝐏𝐏) as shown in Fig. 8. Note that all points in �𝐱𝐱 inp2 , 𝛟𝛟 out2 � have 𝑍𝑍 = 𝑍𝑍 max and 𝜙𝜙 = 0. Combining � 𝐱𝐱 inp1 , 𝛟𝛟 out1 � and �𝐱𝐱 inp2 , 𝛟𝛟 out2 �, we have 345 training points in �𝐱𝐱 inp , 𝛟𝛟 out �. To do the translation mentioned in Subsection 4.5, we update 𝛟𝛟 out simply by 𝛟𝛟 out = 𝛟𝛟 out + ϵ , where ϵ = -6.97 V is obtained with Eq. ( 29). With the updated �𝐱𝐱 inp , 𝛟𝛟 out �, we build an initial Kriging model and then estimate 𝜙𝜙(𝑍𝑍) and 𝜙𝜙(𝑍𝑍) through Eq. ( 16) and Eq. ( 17). Finally, we translate 𝜙𝜙(𝑍𝑍) and 𝜙𝜙(𝑍𝑍)
back by 𝜙𝜙(𝑍𝑍) = 𝜙𝜙(𝑍𝑍) -ϵ and 𝜙𝜙(𝑍𝑍) = 𝜙𝜙(𝑍𝑍) -ϵ. Fig. 9
In this example, we still consider the same 1-D photoelectron sheath problem in Subsection 5.1, but the objective is to estimate 𝐸𝐸(𝑍𝑍) and 𝐸𝐸(𝑍𝑍) and then calculate the dust levitation height. The values of all involved parameters are given in Table 6.
The procedures used to estimate 𝐸𝐸(𝑍𝑍) and 𝐸𝐸(𝑍𝑍) are almost the same as that used to estimate 𝜙𝜙(𝑍𝑍) and 𝜙𝜙(𝑍𝑍). The only difference is that the samples of 𝐸𝐸(𝑍𝑍; 𝐏𝐏) instead of 𝜙𝜙(𝑍𝑍; 𝐏𝐏) are used. The final estimation of 𝐸𝐸(𝑍𝑍) and 𝐸𝐸(𝑍𝑍) is shown in Fig. 11. It shows that the proposed AKM method is very accurate. As for the efficiency, the proposed method needs only 𝑁𝑁 ODE = 𝑁𝑁 in + 18 = 23 ODE solutions. Compared to 𝑁𝑁 𝑃𝑃 = 3,125 ODE solutions needed in MCS, the proposed method is very efficient.
When the upper and lower bounds of the electric field have been determined, we can use the them to determine the levitation heights of the dust grains. Assuming there are two types of dust grains, A and B, in the electric field. The relevant parameters of the grains are given in Table 7, where e = 1.062 × 10 -19 C is the electric charge of an electron. The dust levitation heights are shown in Fig. 12 and given in Table 8. Given any dust grain with known 𝑤𝑤 value, we can easily determine its levitation height interval using the method shown in Fig. 12. This will help to evaluate the risk or damage caused by the levitated dust grains for lunar exploration missions.
We presented an adaptive Kriging based method to perform UQ analysis of the 1-D photoelectron sheath and dust levitation on the lunar surface. A recently derived 1-D photoelectron sheath model was used as the high-fidelity physics-based model and the blackbox function. Adaptive Kriging method, with a task-oriented learning function and stopping criterion, was utilized to improve the efficiency in calculating the upper and lower bounds of electric potential as well as dust levitation height, given the intervals of model parameters.
Experiment analysis shows that the proposed AKM method is both accurate and efficient.
Current and ongoing efforts are focused on building adaptive Kriging model for 2-D and 3-D kinetic particle simulations of lunar plasma/dust environment and perform UQ analysis.