Modelling gravity anomalies from topography using the finite element method with adaptive mesh refinement

Abstract

Listen

Translate article

The first-order finite element method with unstructured hexahedral mesh is applied to gravity forward modelling problems, to calculate gravitational effect of a given density distribution over a domain. Mesh adaptivity is used to facilitate the finite element computation, to characterize small-scale features involved in the domain, and locally adapt mesh where dominant error of approximate solution exists based on an a-posteriori error estimation. Algebraic constraints are used to handle non-conformity of the finite element method at hanging nodes generated from local refinement, and the element-based conjugate gradient method is applied to obtain a solution.Gravity forward modelling is a cornerstone of gravity inversion that has a strong practical application background. In practical applications of gravity forward modelling with large-scale computational domains, complex topography and irregular internal interfaces are inevitable, and these small-scale features require particular attention in domain discretization. In order to reduce discretization error, fine elements are needed in the vicinity of small-scale features to approximate their subtle structures. On the other hand, a uniform fine mesh over the domain is an adverse factor to large-scale computations. Thus to keep a balance between discretization error and total number of elements, adaptive mesh refinement is applied in domain discretization, and fine elements are just used to approximate topography and internal interfaces. Due to the fixed location of topography and internal interfaces, the computational domain is first discretized into a coarse uniform mesh, then elements contacting with topography and internal interfaces are recursively refined till a desirable resolution is achieved. As a result a non-uniform hexahedral mesh with hanging nodes is generated as initial mesh.Domain splitting is used to separate anomalies from a background model, which has a much stronger gravitational effect than anomalies, so that gravity field generated from anomalies can be studied accurately. Numerical solution based on initial mesh is rarely accurate enough, moreover approximation error is not evenly distributed due to anomalies and the use of non-uniform initial mesh, hence global refinement is not an optimal choice to improve accuracy of approximate solution, and a better option is mesh adaptation. According to an error estimation, only the elements where major errors exist are refined in adaptive mesh refinement. The dynamic use of adaptive mesh refinement helps to gain a better balance between solution accuracy and computing cost.Memory-efficiency is an important factor of a numerical method, especially when a multi-scale problem is solved. A strategy that can be used to reduce the requirement on memory is application of the element-based conjugate gradient method. In element-based conjugate gradient method, matrix-vector multiplication is performed in an element-to-element way without explicitly assembling the global finite element matrices, hence memory is saved from accessing and manipulating large matrices. This is particularly helpful if finite element computation is carried out in an iterative process, such as the use of adaptive mesh refinement or solving nonlinear problems. Numerical experiments confirm that the computational domain with complex small-scale features can be accurately approximated by a non-uniform and unstructured hexahedral mesh, and Kelly error estimator is an efficient and accurate error estimator. Comparing with computation using global refinement, finite element computation with mesh adaptation is able to deliver more accurate solutions but with fewer number of unknowns. Using the finite element method with mesh adaptivity based on Kelly error estimator, accurate approximations of gravity field of anomalies are obtained.