## 2 Materials and Methods
### SSIP Equivalent Inversion Method
where \(J_{\rho}=\left[\partial d_{pci}/\partial m_{j}\right]\) represents the partial derivative matrix of simulated data with respect to model parameters. The primary computation methods include the reciprocity principle method (computationally intensive and memory-demanding) and the quasi-Newton method (computationally efficient). In practice, these two approaches are often combined, with the choice between them depending on the scale of the inversion problem. \(\Delta\mathbf{d}_{\rho}=\mathbf{d}_{\rho a}-\mathbf{d}_{pc}\) denotes the data residual vector. The conjugate gradient method is typically employed to solve Equation (4), yielding the model parameter update \(\Delta m_{\rho}\). Through multiple inversion iterations, the subsurface resistivity model \(\mathbf{m}_{\rho}\) can be obtained and substituted into the following equation:  
\[m_{\rho 1}=m_{\rho 0}+u*\Delta m_{\rho} \tag{5}\]  
In Equation (5), \(u\) represents the correction factor for the iteration step length. The 0.618 golden section search method is employed to determine a suboptimal correction factor \(u\). Subsequently, the model parameters for the next inversion iteration can be computed using Equation (5), followed by model updating. This iterative process continues until the mean squared error converges to meet the predefined criteria, ultimately yielding the subsurface resistivity model \(m_{\rho}\). Building upon this foundation, a linear inversion method for chargeability (applicable when chargeability values are relatively small) was employed by solving a single linear inversion equation, and the subsurface chargeability model can be obtained as follows:  
\[({{J_{\rho}}^{T}}{{W_{\eta}}^{T}}{{W_{\eta}}}{{I_{\rho}}}+\lambda_{\eta s}{{S^{ T}}}{{S}}+\lambda_{\eta b}{{B^{T}}}{{B}})\mathbf{\eta}={{J_{\rho}}^{T}}{{W_{\eta}}^{T}} {{\mathbf{\eta}}_{a}}+\lambda_{\eta b}{{B^{T}}}{{B}}{{\mathbf{\eta}}_{b}} \tag{6}\]  
where \({{W_{\eta}}}\) is the chargeability weighting matrix, \(\mathbf{\eta}\) represents the chargeability model parameters, \(\mathbf{\eta}_{a}\) denotes the measured chargeability data, \(\eta_{b}\) stands for the background chargeability distribution, \(\lambda_{\eta s}\) and \(\lambda_{\eta b}\) are respectively the smoothing and background constraint factors for chargeability inversion, and \({{J_{\rho}}}\) is the partial derivative matrix that has been previously computed during resistivity inversion, resulting in significantly reduced computational cost for the linear chargeability inversion. When subsurface polarizability parameters are sufficiently large such that the measured apparent polarizability and polarizability model exhibit nonlinear relationships, the objective function for polarizability inversion is similarly constructed under the least squares criterion as follows:  
\[\Phi_{\eta}(\mathbf{m}_{\eta})=\left\|{{W_{\eta}}\left[{{\mathbf{d}_{\eta a}}-{\mathbf{d }_{\eta c}}(\mathbf{m}_{\eta})}\right]}\right\|_{2}^{2}+\lambda_{\eta s}\left\|{ {\mathbf{S}}{\mathbf{m}_{\eta}}}\right\|_{2}^{2}+\lambda_{\eta b}\left\|{{\mathbf{B}} \left[{{\mathbf{m}_{\eta}}-{\mathbf{m}_{\eta b}}}\right]}\right\|_{2}^{2} \tag{7}\]  
The parameters are analogous to those in Equation (3). By differentiating both sides of Equation (7) with respect to the model parameters and setting the derivatives to zero, the generalized linear inversion equation for polarizability can be derived thusly:  
\[({{J_{\eta}}^{T}}{{W_{\eta}}^{T}}{{W_{\eta}}}{{I_{\eta}}}+\lambda_{\eta s}{{S^ {T}}}{{S}}+\lambda_{\eta b}{{B^{T}}}{{B}})\Delta\mathbf{m}_{\eta}={{J^{T}}}{{W_{ \eta}}^{T}}\Delta{\mathbf{d}_{\eta}}-\lambda_{\eta s}{{S^{T}}}{{S}}{{\mathbf{m}_{\eta} }}+\lambda_{\eta b}{{B^{T}}}{{B}}{{\mathbf{B}}}(\mathbf{m}_{\eta b}-\mathbf{m}_{\eta}) \tag{8}\]  
By solving Equation (8), the polarizability model parameter update \(\Delta\mathbf{m}_{\eta}\) can be obtained and substituted into the following equation:  
\[m_{\eta 1}=m_{\eta 0}+u*\Delta m_{\eta} \tag{9}\]  
The 0.618 golden section search method is implemented to determine a suboptimal correction factor \(u\).