First-order Finite Element Method Research Articles

A fast direct solver for surface-based whole-head modeling of transcranial magnetic stimulation

When modeling transcranial magnetic stimulation (TMS) in the brain, a fast and accurate electric field solver can support interactive neuronavigation tasks as well as comprehensive biophysical modeling. We formulate, test, and disseminate a direct (i.e., non-iterative) TMS solver that can accurately determine global TMS fields for any coil type everywhere in a high-resolution MRI-based surface model with ~ 200,000 or more arbitrarily selected observation points within approximately 5 s, with the solution time itself of 3 s. The solver is based on the boundary element fast multipole method (BEM-FMM), which incorporates the latest mathematical advancement in the theory of fast multipole methods—an FMM-based LU decomposition. This decomposition is specific to the head model and needs to be computed only once per subject. Moreover, the solver offers unlimited spatial numerical resolution. Despite the fast execution times, the present direct solution is numerically accurate for the default model resolution. In contrast, the widely used brain modeling software SimNIBS employs a first-order finite element method that necessitates additional mesh refinement, resulting in increased computational cost. However, excellent agreement between the two methods is observed for various practical test cases following mesh refinement, including a biophysical modeling task. The method can be readily applied to a wide range of TMS analyses involving multiple coil positions and orientations, including image-guided neuronavigation. It can even accommodate continuous variations in coil geometry, such as flexible H-type TMS coils. The FMM-LU direct solver is freely available to academic users.

A vertex scheme for two-phase flow in heterogeneous media

This paper presents the numerical solution of immiscible two-phase flows in porous media, obtained by a first-order finite element method equipped with mass-lumping and flux upwinding. The unknowns are the physical phase pressure and phase saturation. Our numerical experiments confirm that the method converges optimally for manufactured solutions. For both structured and unstructured meshes, we observe the high-accuracy wetting saturation profile that ensures minimal numerical diffusion at the front. Performing several examples of quarter-five spot problems in two and three dimensions, we show that the method can easily handle heterogeneities in the permeability field. Two distinct features that make the method appealing to reservoir simulators are: (i) maximum principle is satisfied, and (ii) local mass error is small and remains below nonlinear solver tolerance.

Planetary Normal Mode Computation: Parallel Algorithms, Performance, and Reproducibility

This article is an extension of work entitled “Computing planetary interior normal modes with a highly parallel polynomial filtering eigensolver.” by Shi <italic xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:xlink="http://www.w3.org/1999/xlink">et al.</i> , <xref ref-type="bibr" rid="ref1" xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:xlink="http://www.w3.org/1999/xlink">[1]</xref> originally presented at the SC18 conference. A highly parallel polynomial filtered eigensolver was developed and exploited to calculate the planetary normal modes. The proposed method is ideally suited for computing interior eigenpairs for large-scale eigenvalue problems as it greatly enhances memory and computational efficiency. In this article, the second-order finite element method is used to further improve the accuracy as only the first-order finite element method was deployed in the previous work. The parallel algorithm, its parallel performance up to 20k processors, and the great computational accuracy are illustrated. The reproducibility of the previous work was successfully performed on the Student Cluster Competition at the SC19 conference by several participant teams using a completely different Mars-model dataset on different clusters. Both weak and strong scaling performances of the reproducibility by the participant teams were impressive and encouraging. The analysis and reflection of their results are demonstrated and future direction is discussed.

Combining the modified discrete element method with the virtual element method for fracturing of porous media

Simulation of fracturing processes in porous rocks can be divided into two main branches: (i) modeling the rock as a continuum enhanced with special features to account for fractures or (ii) modeling the rock by a discrete (or discontinuous) approach that describes the material directly as a collection of separate blocks or particles, e.g., as in the discrete element method (DEM). In the modified discrete element (MDEM) method, the effective forces between virtual particles are modified so that they reproduce the discretization of a first-order finite element method (FEM) for linear elasticity. This provides an expression of the virtual forces in terms of general Hook’s macro-parameters. Previously, MDEM has been formulated through an analogy with linear elements for FEM. We show the connection between MDEM and the virtual element method (VEM), which is a generalization of FEM to polyhedral grids. Unlike standard FEM, which computes strain-states in a reference space, MDEM and VEM compute stress-states directly in real space. This connection leads us to a new derivation of the MDEM method. Moreover, it enables a direct coupling between (M)DEM and domains modeled by a grid made of polyhedral cells. Thus, this approach makes it possible to combine fine-scale (M)DEM behavior near the fracturing region with linear elasticity on complex reservoir grids in the far-field region without regridding. To demonstrate the simulation of hydraulic fracturing, the coupled (M)DEM-VEM method is implemented using the Matlab Reservoir Simulation Toolbox (MRST) and linked to an industry-standard reservoir simulator. Similar approaches have been presented previously using standard FEM, but due to the similarities in the approaches of VEM and MDEM, our work provides a more uniform approach and extends these previous works to general polyhedral grids for the non-fracturing domain.

Efficient Discrete Lagrange Multipliers in three first-order finite element discretizations for the A Posteriori Error Control in an Obstacle Problem

The efficient design of a discrete Lagrange multiplier is a key ingredient in the a posteriori error analysis for the obstacle problem. Affirmative examples exist for three different first-order finite element methods (FEMs), namely, the $P_1$ conforming Courant, the $P_1$ nonconforming Crouzeix--Raviart, and the lowest-order mixed Raviart--Thomas. With those discrete Lagrange multipliers, a general reliable and efficient a posteriori error analysis for the error in the energy norm of the displacement variables applies to all those discretizations for affine obstacles under minimal assumptions.

Reliable Averaging for the Primal Variable in the Courant FEM and Hierarchical Error Estimators on Red-Refined Meshes

Abstract A hierarchical a posteriori error estimator for the first-order finite element method (FEM) on a red-refined triangular mesh is presented for the 2D Poisson model problem. Reliability and efficiency with some explicit constant is proved for triangulations with inner angles smaller than or equal to π 2 ${\frac{\pi }{2}}$ . The error estimator does not rely on any saturation assumption and is valid even in the pre-asymptotic regime on arbitrarily coarse meshes. The evaluation of the estimator is a simple post-processing of the piecewise linear FEM without any extra solve plus a higher-order approximation term. The results also allow the striking observation that arbitrary local averaging of the primal variable leads to a reliable and efficient error estimation. Several numerical experiments illustrate the performance of the proposed a posteriori error estimator for computational benchmarks.

Medius analysis and comparison results for first-order finite element methods in linear elasticity

This paper enfolds a medius analysis for first-order nonconforming finite element methods (FEMs) in linear elasticity named after Crouzeix–Raviart and Kouhia–Stenberg, which are robust with respect to the incompressible limit as the Lamé parameter |$\lambda$| tends to infinity. The new result is a best-approximation error estimate for the stress error in |$L^2$| up to data-oscillation terms. Even for very coarse shape-regular triangulations, two comparison results assert that the errors of the nonconforming FEM are equivalent to those of the conforming first-order FEM. The explicit role of the parameter |$\lambda$| in those equivalence constants leads to an advertisement of the robust and quasi-optimal Kouhia–Stenberg FEM, in particular for nonconvex polygons. The proofs are based on conforming companions, a new discrete Helmholtz decomposition and a new discrete-plus-continuous Korn inequality for Kouhia–Stenberg finite element functions. Numerical evidence strongly supports the robustness of the nonconforming FEMs with respect to incompressibility locking and with respect to singularities, and underlines that the dependence of the equivalence constants on |$\lambda$| in the comparison of conforming and nonconforming FEMs cannot be improved. This work therefore advertises the Kouhia–Stenberg FEM as a first-order robust discretization in linear elasticity in the presence of Neumann boundary conditions.

Corrector Analysis of a Heterogeneous Multi-scale Scheme for Elliptic Equations with Random Potential

This paper analyzes the random fluctuations obtained by a heterogeneous multi-scale first-order finite element method applied to solve elliptic equations with a random potential. We show that the random fluctuations of such solutions are correctly estimated by the heterogeneous multi-scale algorithm when appropriate fine-scale problems are solved on subsets that cover the whole computational domain. However, when the fine-scale problems are solved over patches that do not cover the entire domain, the random fluctuations may or may not be estimated accurately. In the case of random potentials with short-range interactions, the variance of the random fluctuations is amplified as the inverse of the fraction of the medium covered by the patches. In the case of random potentials with long-range interactions, however, such an amplification does not occur and random fluctuations are correctly captured independent of the (macroscopic) size of the patches. These results are consistent with those obtained by the authors for more general equations in the one-dimensional setting and provide indications on the loss in accuracy that results from using coarser, and hence less computationally intensive, algorithms.

Reconstruction convergence and speed enhancement in electrical impedance tomography for domains with known internal boundaries

An improved approach for electrical impedance tomography (EIT) image reconstruction, based on modifying the forward and inverse solutions, is proposed. In this approach, the EIT forward problem is solved via the finite element method (FEM) using two types of elements. The inverse problem is solved by the modified Newton–Raphson method, whereas the condition number of the Hessian matrix is being monitored. At the early stage of the reconstruction, first-order elements are used, and if the condition number exceeds the allowable limit, the algorithm restarts. Otherwise, if the reconstruction error becomes lower than a predefined threshold, second-order elements are employed in the forward solution in order to preserve the precision of the final results. The latter stage converges in very few iterations. Since the solution speed with the first-order FEM is considerably higher than the second-order FEM, the reconstruction speed improves considerably by this approach, whereas the accuracy of the results is guaranteed by the well-conditioned Hessian matrix. Numerical simulations and experiments are followed by comparisons with other reconstruction methods which demonstrate the reliability and high solution speed of this approach. According to the results, the convergence of the proposed method is significantly improved, and its speed is 2–200 times higher than the previously developed methods with the same level of precision.

A second-order finite element algorithm for solving the three-dimensional EEG forward problem

A finite element algorithm has been developed to solve the electroencephalogram (EEG) forward problem. A new computationally efficient approach to calculate the stiffness matrix of second-order tetrahedral elements has been developed for second-order tetrahedral finite element models. The present algorithm has been evaluated by means of computer simulations, by comparing with analytic solutions in a multi-spheres concentric head model. The developed finite element method (FEM) algorithm has also been applied to address questions of interest in the EEG forward problem. The present simulation study indicates that the second-order FEM provides substantially enhanced numerical accuracy and computational efficiency, as compared with the first-order FEM for comparable numbers of tetrahedral elements. The anisotropic conductivity distribution of the head tissue can be taken into account in the present FEM algorithm. The effects of dipole eccentricity, size of finite elements and local mesh refinement on solution accuracy are also addressed in the present simulation study.

Isotropic and Anisotropic Electrostatic Field Computation by Means of the Cell Method

The Cell Method (CM) is an approach to the numerical solution of field physical problems that has been recently applied to the electromagnetic computation. It is gaining more and more interest as it allows easy formulations based on integral (finite) quantities, with no need to resort to local quantities related by differential equations to be suitably regularized in numerical form. The paper deals with an application of the CM to the solution of electrostatic problems. Two formulations are presented: The first is suitable for isotropic electrostatic media, and the second is specifically suited for the anisotropic ones. The main items of the method and the structure of the two algorithms are described, and some results obtained for reference cases are reported. A comparison is also shown with the results from a first-order finite element method (FEM) algorithm.

On magnetic lens computations with FEM and BEM

The aim of the paper is to compare the Boundary Element Method and the first-order Finite Element Method (FEM) for the computation of rotationally symmetric magnetic electron lenses. It is shown that the FEM is accurate for the computation of open lenses and the computation accuracy at off-axis points is sufficient to get aberrations from ray tracing. Open lenses can be satisfactorily computed with the FEM if we use meshes with a large number of points and if we define sufficiently distant boundaries, eventually if we apply both the Dirichlet and Neumann boundary condition. It is also shown that the FEM allows fast and accurate computation of saturated magnetic electron lenses.

Accuracy of electrostatic lens computations with FOFEM

The increased speed of personal computers enables fast computation of rotationally symmetric electrostatic lenses with the first-order finite element method in meshes with a large number of mesh points. In order to produce an estimate of accuracy of the computed potential, we propose a simple procedure based on doubling the number of mesh points in each coordinate. In this way, we can produce for the lower-density mesh at each point the information about error of the potential and visualize the sources of the computation errors.

Accurate computation of magnetic lenses with FOFEM

Abstract First-order finite element method (FOFEM) is used for computations of rotationally symmetric magnetic electron lenses with saturated materials. An overview of local and global errors in the axial flux density distribution in the first-order finite element method is given, and the coefficient computation in FOFEM providing the best axial fields are shown. The quality of computed results is illustrated on several typical problems. It is shown that with graded dense mesh accurate results can be obtained fast and with small memory requirements.

Second‐order finite element method and its practical application in charged particle optics

SUMMARYThe first‐order finite element method (FOFEM) has been widely used in computing electrostatic and magnetic fields of lenses and deflectors. With the increasing requirement of high‐accuracy field computation and the advances in the power of computers, the second‐order finite element method (SOFEM), which will be demonstrated in this paper to be more accurate than FOFEM, is becoming necessary and practical. In this paper, the concept of SOFEM will be described. Several applications of SOFEM in charged particle optics will also be illustrated, such as the simulation of electrostatic and magnetic lenses and deflectors, the modelling of electron and ion sources, and the evaluation of high‐order geometric and chromatic aberrations of charged particle columns.

Finite element analysis of the performance of a ferromagnetic solid-rotor motor

It has been found that use of a real distribution of a three-dimensional electromagnetic field in a solid rotor of finite dimensions will increase the accuracy of the electromagnetic field calculations. This paper presents an analysis of a diphase asynchronous motor with a solid rotor. The calculations have been performed using the first-order finite element method.

High-order Finite-element Method to Solve the Nonlinear Coupled-wave Equations for Degenerate Two-wave and Four-wave Mixing

Abstract We use the one-dimensional finite-element method (FEM) to solve the nonlinear coupled-wave equations for degenerate two-wave mixing (DTWM) and four-wave mixing (DFWM) in a reflection geometry. High-order (quadratic) elements have been used. A comparison is made between the results obtained using the high-order FEM, and those using the first-order FEM, which have been reported previously. High accuracy is achieved using the high-order FEM which employs substantially fewer elements.