## 3 Methodologies  
### Variational Modal Decomposition Algorithm (VMD)  
Variational mode decomposition (VMD) is a new signal adaptive-processing technique used to decompose a non-smooth, nonlinear signal into several discrete intrinsic mode functions (IMFs) [32]. Compared with traditional modal decomposition methods such as EMD, EEMD, etc., the greatest superiority of VMD lies in its adaptive tuning, i.e., a parameter is used to control the number of decomposed modes, and this adaptivity makes VMD more flexible in pre-processing signals of a different nature and complexity [33]. In the case of a preprocessed signal consisting of K bandwidth limited modal components \(u_{k}(t)\), the center frequency of each IMF is denoted by \(\omega(t)\), and the modal sum pairs are approximated to be equal to the original signal, subject to certain error allowances. The next step in the construction of the VMD is as follows:  
By applying the Hilbert transform, the analytical signal of \(u_{k}(t)\) is obtained, its one-sided spectrum is calculated, and an exponential term is added to correct the estimated frequency of the center of each modal function and multiply it with the exponential signal of the one-sided interspectrum to modulate the center band to the corresponding baseband [34]. The square of the gradient of the demodulated signal is computed to determine the estimated fundamental frequency range of each modal function. This process establishes the constrained variational problem, leading to the following mathematical model:  
\[\begin{split}\text{G}=\min\limits_{\{u_{k}\},\{\omega_{k}\}}\left\{ \sum\limits_{k=1}^{K}\partial_{t}\left[\left(\delta(t)+\frac{\text{j}}{\pi t }\right)*u_{k}(t)\right]\text{e}^{-\text{j}\omega_{k}t}\right\|_{2}^{2}\right\} \\ \text{s.t.}\sum\limits_{k=1}^{K}u_{k}=h(t)\end{split} \tag{1}\]  
where G is the minimum of the sum of the \(u_{k}(t)\) estimated bandwidths, K is the number of modal decompositions, \(\partial_{t}\) is the bias of the sought function, \(\delta(t)\) is the Dirac distribution function, \"\(\omega\)\" is the convolution operator, and \(h(t)\) is the original load signal [35].  
Introducing the Lagrange multiplier \(\varepsilon(t)\) and the quadratic penalty factor \(\alpha\), the constrained problem is transformed into an unconstrained problem, and the expression of the augmented Lagrange function \(L\) is obtained as shown in Equation \(\{u_{k}(t)\}\), \(\{\omega_{k}\}\) is the set of the decomposed modal components and the center frequency:  
\[L(\{u_{k}\},\{\omega_{k}\},\varepsilon)=\alpha\sum\limits_{k=1}^{K}\left\| \partial_{t}\left[\left(\delta(t)+\frac{j}{\pi t}\right)*u_{k}(t)\right]e^{-j \omega_{k}t}\right\|_{2}^{2}+\left\|h(t)-\sum\limits_{k=1}^{K}u_{k}(t)\right\| _{2}^{2}+\left\langle\varepsilon(t),h(t)-\sum\limits_{k=1}^{K}u_{k}(t)\right\rangle \tag{2}\]  
The modal components and center frequencies are solved iteratively by the alternating direction multiplier method (ADMM) to find the \"saddle point\" of the improved Lagrangian. Then, the updated \(u_{k}\), \(\omega_{k}\) and \(\lambda\) are:  
\[u_{k}^{n+1}(\omega)=\frac{f(\omega)-\sum\limits_{i\
eq k}u_{i}(\omega)+\frac{ \lambda(\omega)}{2}}{1+2\alpha(\omega-\omega_{k})^{2}} \tag{3}\]\[\omega_{k}^{n+1}=\frac{\int_{0}^{\infty}\omega|\hat{u}_{k}^{n+1}(\omega)|^{2} \mathrm{d}\omega}{\int_{0}^{\infty}|\hat{u}_{k}^{n+1}(\omega)|^{2}\mathrm{d}\omega} \tag{4}\]  
\[\tilde{\lambda}^{n+1}(\omega)=\tilde{\lambda}^{n}(\omega)+\gamma[\hat{f}( \omega)-\sum_{k=1}^{n}\hat{u}_{k}^{n+1}(\omega)] \tag{5}\]  
The judgement accuracy is \(\varepsilon>0\), where \(\varepsilon\) is the iteration error iteration stop condition, reached that stops:  
\[\sum_{k}\left\|\hat{u}_{k}^{n+1}-\hat{u}_{k}^{n}\right\|_{2}^{2}/\left\|\hat{u} _{k}^{n}\right\|_{2}^{2}<\varepsilon \tag{6}\]