Practical Numerical Method Research Articles

Multi-field coupling analysis of mechanical responses in methane hydrate exploitation with a practical numerical approach combining T + H with DEM

The methane hydrate (MH) is widely distributed in deep-sea sedimentary layers, considered as the most promising clean energy to replace conventional oil and gas in the 21st century. Large-scale exploitation of MH may lead to a series of geotechnical engineering problems. To numerically estimate the seabed stability during MH exploitation which should involve multi-field coupling simulation of complex mechanical behaviors of methane hydrate bearing sediment (MHBS), a new practical hybrid numerical method is proposed, which combines the TOUGH + HYDRATE (T + H) code, targeting multiphase flow analysis of hydrate-bearing geologic systems, with the improved Distinct Element Method (DEM). Porosity, temperature, pressure, salinity fields are exchanged between DEM and T + H to realize the thermal-hydro-mechanical-chemical (THMC) coupling and to reproduce the progressive discontinuous failure process. The results of verification show obvious uniformity between numerical and analytical results for one-dimensional (1D) consolidation and 1D heat conduction tests on soils. The paper tries to analyze evolutions of THMC fields, gas production, hydrate bond breakage and other mechanical behaviors of actual depressurization production by using this method, focusing on geotechnical risks of exploitation. In addition, study observes significant changes within certain range from the wellbore for gas production and mechanical responses in the early stage.

A kernel-based method for fractional integro-differential equations with a weakly singular kernel in multi-dimensional complex domains

We present an efficient and practical numerical method for solving fractional integro-differential equations with a weakly singular kernel in complex two and three-dimensional domains. The reproducing kernel pseudo-spectral method is used to discretize the multidimensional regular and irregular spatial domains, while time discretization is done using a third-order backward differentiation formula. We have investigated the unconditional stability of the utilized temporal discretization. Finally, we have used the proposed method to solve the examples presented in the articles and compared the results. The obtained results show the accuracy and ability of the presented method.

Numerical analysis method of equivalent deformation parameters for grouted arch dam abutment fractured rock masses considering scale effects and anisotropy

Research on the equivalent deformation parameters of grouted arch dam abutment rock masses has recently emerged as a key in the design, construction and stability assessment of arch dams. With the aim of synthetically considering scale effects and anisotropy of the rock mass equivalent deformation parameters after grouting, in this paper, a practical numerical calculation method is developed. First, the Latin hypercube sampling method is introduced to generate a network of fractures of the arch dam abutment rock mass. Subsequently, the finite element models of the rock masses with different sizes in various directions are constructed. On this basis, equivalent continuum theory and cohesive element model are applied for numerical calculation of grouting fractured rock masses. Then, the scale effect and anisotropy of the deformation modulus of rock mass under different grouting conditions are investigated. Finally, the Jinping super-high arch dam is taken as an example to validate the theoretical significance and application effect of the proposed method. The analysis results indicate that the process of fracture and fragmentation of grouted rock masses is characterized precisely, and the deformation of material units can be described as fully as possible to calculate more realistic equivalent deformation parameters of the fractured rock masses.

The Strehl ratio as a phase histogram.

It is shown that the Strehl ratio can always be written as an integral over an apodization-weighted phase histogram. The corresponding mathematical formalism, based on Federer's co-area formula, is enumerated, and a practical numerical method to quickly and accurately calculate apodization-weighted phase histograms is detailed and compared with similar methods. Conditions for expressing the Strehl ratio as a product S=S 1 S 2 are investigated.

Numerical Approach and Verification Method for Improving the Sensitivity of Ferrous Particle Sensors with a Permanent Magnet

This study aimed to improve the sensitivity of ferrous particle sensors used in various mechanical systems such as engines to detect abnormalities by measuring the number of ferrous wear particles generated by metal-to-metal contact. Existing sensors collect ferrous particles using a permanent magnet. However, their ability to detect abnormalities is limited because they only measure the number of ferrous particles collected on the top of the sensor. This study provides a design strategy to boost the sensitivity of an existing sensor using a multi-physics analysis method, and a practical numerical method was recommended to assess the sensitivity of the enhanced sensor. The sensor’s maximum magnetic flux density was increased by around 210% compared to the original sensor by changing the core’s form. In addition, in the numerical evaluation of the sensitivity of the sensor, the suggested sensor model has improved sensitivity. This study is important because it offers a numerical model and verification technique that may be used to enhance the functionality of a ferrous particle sensor that uses a permanent magnet.

A Second-Order Accurate Numerical Approximation for a Two-Sided Space-Fractional Diffusion Equation

In this paper, we investigate a practical numerical method for solving a one-dimensional two-sided space-fractional diffusion equation with variable coefficients in a finite domain, which is based on the classical Crank-Nicolson (CN) method combined with Richardson extrapolation. Second-order exact numerical estimates in time and space are obtained. The unconditional stability and convergence of the method are tested. Two numerical examples are also presented and compared with the exact solution.

Further Study on One of the Numerical Methods for Pure Loss of Stability in Stern Quartering Waves

The International Maritime Organization (IMO) finalized the second-generation intact stability criteria in 2022. However, an accurate and practical numerical method for stability loss has yet to be established. Therefore, a 6 DOF numerical model is further improved based on the previous study. Firstly, the rolling motion is simulated using a seakeeping model instead of the previous maneuvering mathematical model. Secondly, the roll-restoring variation is calculated directly considering the instantaneous wet hull instead of the previous pre-calculated method. Thirdly, transferring frequency to time is used to obtain heave and pitch motions, further considering yaw angle and sway velocity. Fourthly, the dynamic forces for sway, roll, and yaw motions are calculated, further considering the effect of the speed variation. Fifthly, the 6 DOF motions are used to determine the instantaneous wet hull, and the FK force and the hydrostatic force are calculated by the body’s exact method. Finally, a new conclusion is obtained that the sway and yaw motions’ effect on the ship speed loss, the relative longitudinal wave profile by the speed loss, the rudder angles, and the accompanying rudder forces in the rolling direction are significant, and much more than their centrifugal force or coupled force in the rolling direction.

Influence of high–low hybrid towers on the wind flow and energy of complex terrain

In this paper, an Actuator Disk Model (ADM) was introduced into a commercial computational fluid dynamic software and used to simulate the wake effect of wind turbines in a wind field with complex terrain. Then, the influence of towers with difference height on wake flow and power generation can be studied, which is particularly important for the wind field design in complex terrain. The numerical results are verified by means of comparison with the measured data of Supervisory Control and Data Acquisition (SCADA). More detailed numerical results show that the optimization of the design of high and low hybrid towers with difference size wind turbines can not only make full use of the advantages of high wind speed at high altitudes, but also reasonably improve the wake flow generated by the front wind turbine and increase the power generation of the wind field. It is found that the power increment generated by the high–low hybrid towers arrangement is roughly inversely proportional to U00.877, when only evaluating the effect of tower height; moreover, if the rotor diameter is further increased, the power gain of the hybrid arrangement scheme is roughly inversely proportional to U01.167. The research provides a robust and practical numerical method for the design of wind field layout and the optimization of high and low hybrid towers.

Efficient and practical phase-field method for the incompressible multi-component fluids on 3D surfaces with arbitrary shapes

The fluid dynamics on a curved surface can be used to approximately describe many physical phenomena in the natural world, such as the atmospherical circulation on a planet. The present work develops a simple, practical, and efficient numerical method for the multi-phase fluid simulations on arbitrarily curved surfaces in 3D space. The multi-component Cahn–Hilliard type phase-field fluid model is adopted to capture the motion of multiple interfaces. The velocity field and pressure are updated by solving the incompressible Navier–Stokes equations. To achieve the decoupled computations of pressure and velocity, the Chorin's projection method with pressure correction is considered. A linear semi-implicit time-marching scheme is designed for the phase-field equations. By adopting the embedding method and appropriate pseudo-Neumann boundary conditions, the direct calculations on the surface can be extended into a 3D narrow band domain including the surface. The standard Laplacian operator can be used to replace the surface Laplace–Beltrami operators and the finite difference method is used to perform the discretization in space. In each time step, all variables can be updated in a step-by-step manner because the proposed method is totally decoupled. Therefore, the whole algorithm is highly efficient and easy to implement. Extensive numerical experiments indicate that the proposed method has good performance for the multi-component simulations on curved surfaces with complex shapes.

On the transition of the contact eccentricity of an elastic layer resting on a rigid frictionless substrate

Eccentric distributions of contact stress and area would be induced when friction is taken into account, and the problem becomes more complex in a layer-substrate system. To obtain a detailed understanding of the contact eccentricity, a comprehensive analysis on the coupling effects among Poisson's ratio, layer thickness and the coefficient of friction is conducted in this paper. By applying the Fourier transform techniques, the problem is converted to the solution of a singular integral equation of the second kind, and a practical numerical method based on an iterative process is proposed. The results obtained are compared with previous studies. It is seen that eccentric shifting of the pressure distribution and contact zone would be more apparent with increasing coefficient of friction, and the moving direction would be reversed when Poisson's ratio exceeds a critical value. Results are presented for the intriguing Poisson's ratio which could significantly reduce the eccentricity.

Simple and practical method for the simulations of two-component PFC models for binary colloidal crystals on curved surfaces

Herein, we propose a practical numerical method for solving two-component phase-field crystal model with vacancy effect on curved surfaces and perform various numerical simulations. In the present work, the curved surfaces are defined by the zero level-set of signed distance functions. By adopting the pseudo-Neumann boundary condition and closest-point type method, we can extend the computation into a three-dimensional narrow band domain which contains the curved surfaces. The spatial discretization is performed by using the standard finite difference method and the seven-point Laplacian operator can be used to replace the surface Laplace–Beltrami operator. Therefore, the numerical implementation is simple and efficient. The semi-implicit discretization with stabilization term is considered to update the solutions in time. By simulating the pattern formations on curved surfaces with various shapes, we validate the good performance of the proposed method.

Simulating a pulsed-power-driven plasma with ideal MHD

We describe a simple practical numerical method for simulating plasma driven within a vacuum chamber by a pulsed power generator. Typically, in this type of simulation, the vacuum region adjacent to the plasma is approximated as a highly resistive, light fluid; this involves computationally expensive solvers describing the diffusion of the magnetic field through this fluid. Instead, we provide a recipe for coupling pulsed power generators to the magnetohydrodynamics (MHD) domain by approximating the perfectly insulating vacuum as a light, perfectly conducting, inviscid MHD fluid and discuss the applicability of this counter-intuitive technique. This much more affordable ideal MHD representation is particularly useful in situations where a plasma exhibits interesting three-dimensional phenomena, either due to the design of the experiment or due to developing instabilities. We verified that this coupling recipe works by modeling an exactly solvable flux compression generator as well as a self-similar Noh-like solution and demonstrated convergence to the theoretical solution. We also showed examples of simulating complex three-dimensional pulsed power devices with this technique. We release our code implementation to the public (see https://github.com/beresnyak/verif_coupling).

A new viscoelastic artificial boundary with improved numerical stability in explicit calculation of wave propagation problems in infinite domains

The viscoelastic artificial boundary is currently a practical numerical method for addressing the wave propagation problem in infinite domain. When applying this method in explicit dynamic calculations, however, the stability conditions of the boundary area are more stringent than those of the internal domain, leading to smaller critical integral time increments and heavy computational workloads. In this study, the stability conditions of the viscoelastic artificial boundary and their influencing factors in an explicit integral algorithm are analyzed. On that basis, a stability improvement method by introducing added mass to the artificial boundary is proposed. Then through a parametric analysis, the recommended values of the introduced added mass are presented. After adopting these recommended values, the stability conditions of the artificial boundaries get loosed, and the stability of the overall numerical system is controlled by the internal area. Therefore the stability condition of the internal domain can be used to estimate the critical time increment in the explicit calculation. Numerical examples show that the improved viscoelastic artificial boundary proposed in this study can maintain the accuracy of the original method while significantly improving the calculation efficiency, thus offering strong practicability.

Sequential convex programming for non-linear stochastic optimal control

This work introduces a sequential convex programming framework for non-linear, finitedimensional stochastic optimal control, where uncertainties are modeled by a multidimensional Wiener process. We prove that any accumulation point of the sequence of iterates generated by sequential convex programming is a candidate locally-optimal solution for the original problem in the sense of the stochastic Pontryagin Maximum Principle. Moreover, we provide sufficient conditions for the existence of at least one such accumulation point. We then leverage these properties to design a practical numerical method for solving non-linear stochastic optimal control problems based on a deterministic transcription of stochastic sequential convex programming.

Static analysis of axially loaded piles in multilayered soils using a BEM/FEM formulation

In this paper a practical numerical method for calculating the settlement response of piles embedded in multilayered soils is presented. Finite beam elements are used for the discretization of piles. The soil mass, which is considered as a semi-infinite or finite medium with each layer being elastic, homogeneous and isotropic, is modeled with the boundary element method and the Steinbrenner's approximation. In addition, the critical case in which the elastic modulus is decreasing with depth is appropriately treated by considering the soil as a half space inhomogeneous medium. The present formulation, besides yielding fully adequate results for design purposes when compared with other approaches, presents a more practical problem modeling and data preparation. Results of numerical examples considering single piles and pile groups are compared with those of other formulations in order to demonstrate the efficiency and accuracy of the proposed formulations.

Numerical method for slope stabilization design with piles

In this paper, a detailed numerical approach is introduced for the slope stabilization design with piles considering the pile capability as well as the slope stability. Two simulation methods for stabilizing piles, i.e., the solid element simulation and the composite element simulation, are introduced and discussed about their applicability. The solid element simulation is found to be potential due to good accuracy and computation efficiency. Interface is employed for pile-soil interaction. Together with strength reduction method, the numerical techniques are applied in the reinforcement for a three-dimensional slope with one row of piles. By varying pile length, pile spacing and pile location, the optimal arrangement of pile reinforcement for the slope is discussed considering both slope stability and pile mechanical behavior. Therefore, an economical as well as effective pile stabilization design is able to be achieved through the suggested practical numerical method in this paper.

Practical numerical method for erosion risk prediction on ship propellers

It is important to make predictions of cavitation-induced erosion risk on ship propellers in the design phase. Since a cavitation tunnel test on a propeller model coated by soft paint, that is, a standard experimental method for evaluating erosion risk, is costly and time-consuming, numerical methods are necessary for erosion risk predictions. DES is made for cavitating flows around the propeller with a numerically modelled hull wake at the inflow. After achieving a converged solution, an erosion risk index is computed in each cell connecting to the blade surface and accumulated over a propeller rotation. Cavitation simulations are made for two propellers designed for a single-screw ship, of which one showed an erosion indication and the other showed no indication after cavitation tunnel tests with soft paint coating. Three index formulations are compared with the experiment result. The high value region of Index 1 based on the potential energy density of collapsing bubbles corresponds better with the eroded spot indicated by partial and complete paint removals in the experiment than those of the other indices. The maximum value of Index 1 for the non-eroded propeller is lower by more than an order of magnitude than that for the eroded one, whereas the maximum values of the other indices are of the same order of magnitude for both propellers. The validation of Index 1 is in agreement with the criterion that the maximum index needs to be below 1,000 J/m3 for erosion-free propeller designs. The design evolution based on the erosion risk index and propulsive efficiency from CFD shows that it can be a practical tool for a quantitative evaluation of blade surface erosion risk in the propeller design phase.

Numerical Solution of Two Dimensional Time-Space Fractional Fokker Planck Equation With Variable Coefficients

This paper presents a practical numerical method, an implicit finite-difference scheme for solving a two-dimensional time-space fractional Fokker–Planck equation with space–time depending on variable coefficients and source term, which represents a model of a Brownian particle in a periodic potential. The Caputo derivative and the Riemann–Liouville derivative are considered in the temporal and spatial directions, respectively. The Riemann–Liouville derivative is approximated by the standard Grünwald approximation and the shifted Grünwald approximation. The stability and convergence of the numerical scheme are discussed. Finally, we provide a numerical example to test the theoretical analysis.

A sharp interface Lagrangian-Eulerian method for rigid-body fluid-structure interaction

This paper introduces a sharp interface method to simulate fluid-structure interaction (FSI) involving rigid bodies immersed in viscous incompressible fluids. The capabilities of this methodology are benchmarked using a range of test cases and demonstrated using large-scale models of biomedical FSI. The numerical approach developed herein, which we refer to as an immersed Lagrangian-Eulerian (ILE) method, integrates aspects of partitioned and immersed FSI formulations by solving separate momentum equations for the fluid and solid subdomains, as in a partitioned formulation, while also using non-conforming discretizations of the dynamic fluid and structure regions, as in an immersed formulation. A simple Dirichlet-Neumann coupling scheme is used, in which the motion of the immersed solid is driven by fluid traction forces evaluated along the fluid-structure interface, and the motion of the fluid along that interface is constrained to match the solid velocity and thereby satisfy the no-slip condition. To develop a practical numerical method, we adopt a penalty approach that approximately imposes the no-slip condition along the fluid-structure interface. In the coupling strategy, a separate discretization of the fluid-structure interface is tethered to the volumetric solid mesh via stiff spring-like penalty forces. Our fluid-structure coupling scheme relies on an immersed interface method (IIM) for discrete geometries, which enables the accurate determination of both velocities and stresses along complex internal interfaces. Numerical methods for FSI can suffer from instabilities related to the added mass effect, but computational tests indicate that the methodology introduced here remains stable for selected test cases across a broad range of solid-fluid mass density ratios, including extremely small, nearly equal, equal, and large density ratios. Biomedical FSI demonstration cases include results obtained using this method to simulate the dynamics of a bileaflet mechanical heart valve in a pulse duplicator, and to model transport of blood clots in a patient-averaged anatomical model of the inferior vena cava.

Development of a Numerical Method to Detect Multiple Positions of CO2 Seepage on the Seafloor

When Carbon dioxide Capture and Storage (CCS), which is considered as a promising way of mitigating global warming, is conducted in the sea, it is obliged to monitor CO2 concentration in the sea, considering the risk of CO2 seepage from the reservoir into the sea and consequent impacts on the oceanic environment. If abnormally high concentration is observed, the seepage positions must be found at a high cost and a lot of time. Therefore, to narrow down the search range of the seepage positions, a practical numerical method is required. In this study, a numerical method of estimating multiple CO2 source positions was used for detecting two numerically set seepage positions in a three-dimensional domain with unsteady flow and it is showed that the method is effective for the objectives. We also proposed a procedure to improve the estimation accuracy of the method.