# A Stable and High-Precision Downward Continuation Method of Magnetic Data
## 2 Methodology of the Proposed Method
For this reason, the logarithmic barrier method is introduced to solve Equation (9), and the objective function is given by [23]  
\[\mathbf{\Phi}=\|\mathbf{G}\mathbf{k}-\mathbf{d}_{obs}\|_{2}^{2}+\lambda^{2} \|\mathbf{W}_{k}\mathbf{Z}\mathbf{k}\|_{2}^{2}-2\mu\sum_{j=1}^{M}\ln(k_{j}) \tag{13}\]  
The \(\mu\) in Equation (13) is the barrier parameter, and \(-2\mu\sum_{j=1}^{M}\ln(k_{j})\) is the barrier term which is applied to integrate physical bound constraints into the objective function. After that, the Newton method is used to minimize Equation (13) with the regularization parameter fixed during this process. The \(\mu\) is decreased with iterations, and the closer it gets to zero, the more Equation (9) approaches the final result. The nth iteration result is given as Equation (14). \[(\mathbf{G}^{T}\mathbf{G+}\mathbf{\lambda W}_{k}^{T}\mathbf{W}_{k}+\mu^{(n)} \mathbf{X}^{-2})\Delta\mathbf{k}\mathbf{=-\mathbf{G}^{T}\delta\mathbf{d}- \lambda\mathbf{W}_{k}^{T}\mathbf{W}_{k}\delta\mathbf{k}+\mu^{(n)}\mathbf{X}^ {-1}\mathbf{e} \tag{14}\]  
where \(\mathbf{X}=diag\{k_{1},\cdots,k_{M}\}\), \(\mathbf{e}=\{1,\cdots,1\}\), \(\delta\mathbf{d}=\mathbf{G}\mathbf{k}^{(n-1)}-\mathbf{d}_{obs}\) and \(\delta\mathbf{k}=\mathbf{k}^{(n-1)}-\mathbf{k}^{(0)}\). After that, for preventing the outcome of each iteration fluctuating around the reasonable result instead approaching it, the reduced step length is applied to update the new model, as Equation (15) shows. Equation (16) is the maximum permissible step length of the nth iteration. \[\mathbf{k}^{(n)}\mathbf{=}\mathbf{k}^{(n-1)}\mathbf{+}\gamma\beta^{(n)} \Delta\mathbf{k} \tag{15}\]  
\[\beta^{(n)}=\left\{ \begin{array}{ll}1,&if\quad\Delta\mathbf{k}>0\\ \min_{\Delta\mathbf{k}_{j}<0}\frac{\kappa_{j}^{(n-1)}}{|\Delta\mathbf{k}_{j}|},&\quad\text{otherwise}\end{array} \right. \tag{16}\]  
Then the barrier parameter can be updated by Equation (18) and \(\gamma=0.925\) in this paper. \[\mu^{(n+1)}\mathbf{=}[1-\min(\gamma\beta^{(n)})]\mu^{(n)} \tag{17}\]  
Using the conjugate gradient method (CG) to solve Equation (13) by iteration, when \(\mu\) satisfies the preset conditions or reaches the upper limit of the iteration, the iteration ends. Finally, \(\mathbf{K}^{(n)}\) is the final model parameter, and based on that, we can obtain the magnetic data of the target depth by forward calculating. It is assumed that the magnetic data of the downward continuation surface in Figure 1 is \(\mathbf{d}_{down}\). The distance of downward continuation is \(h\); since the target interface of the downward continuation can be a plane or an undulating surface, naturally, the \(h\) can be constant or a different distance at each point. The forward operator, Equation (4), transforms to Equation (18). \[ \begin{split}& g_{ij}=\frac{H\cdot K}{4\pi}\{d_{1}\cdot\ln[r+(x_{mj }-x_{i})]+d_{2}\cdot\ln[r+(y_{mj}-y_{i})]+d_{3}\cdot\ln[r+(z_{mj}-h_{i})]\\ &\quad+d_{4}\cdot\arctan\frac{(x_{mj}-x_{i})(y_{mj}-y_{i})}{(x_{ mi}-x_{i})^{2}+r(z_{mj}-h_{i})+(z_{mj}-h_{i})^{2}}\\ &\quad+d_{5}\cdot\arctan\frac{(x_{mj}-x_{i})(y_{mj}-y_{i})}{(y_{ mi}-y_{i})^{2}+r(z_{mj}-h_{i})+(z_{mj}-h_{i})^{2}}\\ &\quad+d_{6}\cdot\arctan\frac{(x_{mj}-x_{i})(y_{mj}-y_{i})}{r(z _{mj}-h_{i})}\big{|}c^{(x_{mj}-x_{i})(y_{mj}-y_{i})}_{x_{jm}}|c^{z_{jm}}_{y_{ jmin}}|c^{z_{jm}}_{z_{jmax}}\\ \end{split}  \tag{18}\]  
where \(r=[(x_{mj}-x_{i})^{2}+(y_{mj}-y_{i})^{2}+(z_{mj}-h_{i})^{2}]^{1/2}\), \(j=1\cdots M\), and \(i=1\cdots N\). The final result of the downward continuation can be obtained by Equation (19). The specific implementation of the proposed method is shown in Algorithm 1. Based on this algorithm, the magnetic data of arbitrary surface can be obtained. \[\mathbf{d}_{down}=\mathbf{G}\mathbf{k}^{(n)} \tag{19}\]  
```
Preparation : \(\mathbf{d}_{obs},\mathbf{k}_{min},\mathbf{k}_{max}\), maximum iterations \(N_{max}\). 1. calculate \(\mathbf{G}\), \(\mathbf{W}_{k}\), and \(\mathbf{Z}\) after adding constraint data. 2. Initialize \(\mathbf{k}^{(0)}=0.01\) and \(n=\)0. 3. Estimate the regularization parameter \(\lambda\). 4.