# A Fast Method for Calculation of Marine Gravity Anomaly
## 2 Method
### The Proposed Fast Method  
In the proposed fast method, the calculation process is improved to avoid judgment statements and GPU parallel computing technology is used to implement the calculation process. Generally, Equation (3) uses gridded data of global geoid height undulation to compute the marine gravity anomaly of one point. Therefore, a calculation exists of a double loop along the latitude and longitude directions. The marine gravity anomalies of the study area are also computed into gridded data, which leads to another double loop. To reduce the number of loops and simplify the program, the gridded data are transformed into a row vector or a column vector, which is shown in Figure 1, and their latitude and longitude coordinates are transformed into row vectors or column vectors. The classical method needs to judge whether the calculation point \(P\) and the integration point \(Q\) coincide. If the two points do not coincide, the integration point \(Q\) is incorporated into the calculation of Equation (3). However, the judgment nested within the loop increases the computing time. In the ISF method, the integral area includes the study area, and the integration point \(Q_{p}\) that coincides with calculation point \(P\) can be found easily. Therefore, in the proposed fast method, we generate an index vector that records the location (index) of each calculation point \(P\) in the integration points (or integral grids). In the calculation of \(\Delta g_{P}(\varphi_{p},\lambda_{p})\), we first calculate each component term \((N_{q}-N_{p})\cos\varphi_{q}/\sin^{3}(\psi/2)\) of the cumulative summation part in Equation (3) directly, regardless of the current integration point \(Q\) being point \(Q_{p}\) or not. Then, we use the index vector to find the component term \((N_{q}-N_{p})\cos\varphi_{q}/\sin^{3}(\psi/2)\) that corresponds to the integration point \(Q_{p}\) and it is set as 0. Finally, we use all of the component terms \((N_{q}-N_{p})\cos\varphi_{q}/\sin^{3}(\psi/2)\) to calculate the final value of Equation (3). These operations can solve the singularity problem without judgment statements, which can accelerate the computing process. Each component term \((N_{q}-N_{p})\cos\varphi_{q}/\sin^{3}(\psi/2)\) of the cumulative summation part in Equation (3) can be computed independently, and each marine gravity anomaly in the study area can be computed independently as well. Therefore, GPU parallel computing technology can be used in the proposed fast method. GPUs were initially used for computer games. Up to now, GPUs have been developed to provide computational power for scientific applications [21]. Parallel computing based on GPUs has been widely applied in geoscience [19; 20; 21; 22; 23; 24]. MATLAB language has the advantages of powerful data processing capability, easy programming, and numerous multifunctional integrated functions. It can allow the programmer to directly perform GPU computing without knowledge of OpenMP architecture and low-level GPU computing libraries [24]. Therefore, in this study, we used MATLAB language to implement the proposed fast method. The MATLAB code for the proposed fast method is shown in Figure 2. The code in Figure 2 uses the element-wise operation to compute each component term \((N_{q}-N_{p})\cos\varphi_{q}/\sin^{3}(\psi/2)\) of the cumulative summation part, which can avoid loops and facilitate parallelization. After setting the component term \((N_{q}-N_{p})\cos\varphi_{q}/\sin^{3}(\psi/2)\) corresponding to the integration point \(Q_{p}\) to 0, the sum built-in function of MATLAB is used to compute the cumulative summation part. In the loop, resizing the vector of gravity anomaly each time involves memory deallocation or allocation and value copies, which is time-consuming [25]. Therefore, we preallocated the maximum amount of space required for the vector of gravity anomaly before the loop, and the computing time was reduced further. Given the precision of the results, this code should be performed with double precision.