## 3 Main Methodologies
### GMM-Based APF Design  
When performing the proximity operation to the target, the complex shape of the target, as shown in Figure 1, must be considered for the requirement of collision avoidance. Therefore, GMM is applied, which can give an analytical expression and more accurate depiction of the target outline. The GMM can be regarded as a model composed of K single Gaussian probability density functions [28], as shown in Equation (10). These K sub-models are hidden variables of the mixed model. The GMM is used here because the Gaussian distribution has good mathematical properties and good computing performance. \[p(\mathbf{\mu}|\Theta)=\sum\limits_{j=1}^{K}\pi_{j}\Phi\Big{(}\mathbf{\mu}| \mu_{j},\Sigma_{j}\Big{)}=\sum\limits_{j=1}^{K}\frac{\pi_{j}}{(2\pi)^{\frac{ \mathbb{D}}{2}}\big{|}\Sigma_{j}\big{|}^{2}}e^{-\frac{1}{2}(\mathbf{\tau}-\mathbf{ \mu}_{j})^{\mathbb{T}}\Sigma_{j}^{-1}(\mathbf{\tau}-\mathbf{\mu}_{j})}, \tag{10}\]  
where, \(\Phi\Big{(}\mathbf{\mu}|\mu_{j}\mathbf{\Sigma}_{j}\Big{)}\) means the GMM's \(j\)th gaussian component. \(\mathbf{\mu}_{j}\in R^{3}\), \(\Sigma_{j}\in R^{3\times 3}\), and \(\pi_{j}\in R\) are the expectation, covariance matrix, and weight of \(j\)th gaussian component, respectively. Meanwhile, \(\pi_{j}\) satisfies the constrain \(\sum\limits_{j=1}^{K}\pi_{j}=1\). \(D\) is the dimension of position vector \(\mathbf{r}\). \(\Theta=\Big{\{}\theta_{j}=\Big{(}\pi_{j},\mathbf{\mu}_{j},\Sigma_{j}\Big{)} \Big{\}}\) contains three key parameters, which determine the GMM. To describe the outline by GMM, the three key parameters should be solved first. Utilizing servicer's camera or other sensing equipment, a series of feature points of the target's surface can be obtained in space, such as point-cloud model [26], which is recorded as \(R=\Big{\{}\mathbf{r}_{i}=[x_{i},y_{i},z_{i}]^{T},i=1,2,\cdots,N\Big{\}}\). Then, using the EM algorithm, the outline, described by the feature points, can be fitted by GMM. The steps of solving GMM by EM algorithm are shown below [43]. 1. Set the initial parameters \(\Theta\). 2. E-step: Calculate the probability \(\gamma_{ij}\) that each data \(r_{i}\) comes from sub-model \(j\), based on above parameters. \[\gamma_{ij}=\frac{\pi_{j}p\Big{(}\mathbf{\tau}_{i}\big{|}\theta_{j}\Big{)}}{p (\mathbf{\tau}_{i}|\Theta)}.\] (11)
3. M-step: Calculate the next generation parameters as follows, \[\mathbf{\mu}_{j}=\frac{\sum_{i=1}^{N}\gamma_{ij}\mathbf{\tau}_{i}}{\sum_{i=1} ^{N}\gamma_{ij}},j=1,2,\ldots,K,\] (12)  
Figure 2: Flowchart of the two-spacecraft proximity operation process. \[\Sigma_{j}=\frac{\sum_{i=1}^{N}\gamma_{ij}\Big{(}r_{i}-\mu_{j}\Big{)}\Big{(}r_{i}- \mu_{j}\Big{)}^{T}}{\sum_{i=1}^{N}\gamma_{ij}},j=1,2,\ldots,K, \tag{13}\]  
\[\pi_{j}=\frac{\sum_{i=1}^{N}\gamma_{ij}}{N},j=1,2,\ldots,K. \tag{14}\]  
Iterate the E-step and M-step until satisfying the condition in Equation (15). \(\varepsilon\) is small positive constant, and the superscript \(k\) means the \(k\)th calculation. \[\|\Theta^{(k+1)}-\Theta^{(k)}\|<\varepsilon \tag{15}\]  
In actual calculation process, the initial parameters greatly affect the calculation results, and the E-step often takes a long time. To improve the computational efficiency, in other words, to optimize initial parameters \(\mu_{j}\), the K-means clustering algorithm is employed. The main effect of K-means is to divide a set of data into K clusters according to Euclidean distance, the steps of which are as follows [44]. 1. Set the initial K means. 2. Calculate the Euclidean distance between \(r_{i}\) and each mean \(\mu_{j}\) (j = 1, 2,..., K), then, dividing \(r_{i}\) to the nearest cluster center, i.e., mean and recording as \(r_{ij}\). 3. Recalculate each cluster center by the data in same cluster. 4.