We present an unconditionally energy stable and uniquely solvable finite difference scheme for the Cahn-Hilliard-Brinkman (CHB) system, which is comprised of a Cahn-Hilliard-type diffusion equation and a generalized Brinkman equation modeling fluid flow. The CHB system is a generalization of the Cahn-Hilliard-Stokes model and describes two phase very viscous flows in porous media. The scheme is based on a convex splitting of the discrete CH energy and is semi-implicit. The equations at the implicit time level are nonlinear, but we prove that they represent the gradient of a strictly convex functional and are therefore uniquely solvable, regardless of time step size. Owing to energy stability, we show that the scheme is stable in the time and space discrete $ℓ^∞$(0,$T$;$H^1_h$) and $ℓ^2$(0,$T$;$H^2_h$) norms. We also present an efficient, practical nonlinear multigrid method – comprised of a standard FAS method for the Cahn-Hilliard part, and a method based on the Vanka smoothing strategy for the Brinkman part – for solving these equations. In particular, we provide evidence that the solver has nearly optimal complexity in typical situations. The solver is applied to simulate spinodal decomposition of a viscous fluid in a porous medium, as well as to the more general problems of buoyancy- and boundary-driven flows.

In event-driven algorithms for simulation of diffusing, colliding, and reacting particles, new positions and events are sampled from the cumulative distribution function (CDF) of a probability distribution. The distribution is sampled frequently and it is important for the efficiency of the algorithm that the sampling is fast. The CDF is known analytically or computed numerically. Analytical formulas are sometimes rather complicated making them difficult to evaluate. The CDF may be stored in a table for interpolation or computed directly when it is needed. Different alternatives are compared for chemically reacting molecules moving by Brownian diffusion in two and three dimensions. The best strategy depends on the dimension of the problem, the length of the time interval, the density of the particles, and the number of different reactions.

Numerical methods are developed to simulate the wave propagation in heterogeneous 2D fluid/poroelastic media. Wave propagation is described by the usual acoustics equations (in the fluid medium) and by the low-frequency Biot's equations (in the porous medium). Interface conditions are introduced to model various hydraulic contacts between the two media: open pores, sealed pores, and imperfect pores. Well-posedness of the initial-boundary value problem is proven. Cartesian grid numerical methods previously developed in porous heterogeneous media are adapted to the present context: a fourth-order ADER scheme with Strang splitting for time-marching; a space-time mesh-refinement to capture the slow compressional wave predicted by Biot's theory; and an immersed interface method to discretize the interface conditions and to introduce a subcell resolution. Numerical experiments and comparisons with exact solutions are proposed for the three types of interface conditions, demonstrating the accuracy of the approach. 

A new approach to high-order accuracy for the numerical solution of conservation laws introduced by Huynh and extended to simplexes by Wang and Gao is renamed CPR (correction procedure or collocation penalty via reconstruction). The CPR approach employs the differential form of the equation and accounts for the jumps in flux values at the cell boundaries by a correction procedure. In addition to being simple and economical, it unifies several existing methods including discontinuous Galerkin, staggered grid, spectral volume, and spectral difference. To discretize the diffusion terms, we use the BR2 (Bassi and Rebay), interior penalty, compact DG (CDG), and I-continuous approaches. The first three of these approaches, originally derived using the integral formulation, were recast here in the CPR framework, whereas the I-continuous scheme, originally derived for a quadrilateral mesh, was extended to a triangular mesh. Fourier stability and accuracy analyses for these schemes on quadrilateral and triangular meshes are carried out. Finally, results for the Navier-Stokes equations are shown to compare the various schemes as well as to demonstrate the capability of the CPR approach.

The commonly used incompressible phase field models for non-reactive, binary fluids, in which the Cahn-Hilliard equation is used for the transport of phase variables (volume fractions), conserve the total volume of each phase as well as the material volume, but do not conserve the mass of the fluid mixture when densities of two components are different. In this paper, we formulate the phase field theory for mixtures of two incompressible fluids, consistent with the quasi-compressible theory [28], to ensure conservation of mass and momentum for the fluid mixture in addition to conservation of volume for each fluid phase. In this formulation, the mass-average velocity is no longer divergence-free (solenoidal) when densities of two components in the mixture are not equal, making it a compressible model subject to an internal constraint. In one formulation of the compressible models with internal constraints (model 2), energy dissipation can be clearly established. An efficient numerical method is then devised to enforce this compressible internal constraint. Numerical simulations in confined geometries for both compressible and the incompressible models are carried out using spatially high order spectral methods to contrast the model predictions. Numerical comparisons show that (a) predictions by the two models agree qualitatively in the situation where the interfacial mixing layer is thin; and (b) predictions differ significantly in binary fluid mixtures undergoing mixing with a large mixing zone. The numerical study delineates the limitation of the commonly used incompressible phase field model using volume fractions and thereby cautions its predictive value in simulating well-mixed binary fluids.

In this paper, we consider variational approaches to handle the multiplicative noise removal and deblurring problem. Based on rather reasonable physical blurring-noisy assumptions, we derive a new variational model for this issue. After the study of the basic properties, we propose to approximate it by a convex relaxation model which is a balance between the previous non-convex model and a convex model. The relaxed model is solved by an alternating minimization approach. Numerical examples are presented to illustrate the effectiveness and efficiency of the proposed method.

We propose a mathematically rigorous method to measure the spontaneous curvature of a bilayer membrane by molecular dynamics (MD) simulation, which provides description of the molecular mechanisms that cause the spontaneous curvature. As a main result, for the membrane setup investigated, the spontaneous curvature is proved to be a constant plus twice the mean curvature of the membrane in its tensionless ground state. The spontaneous curvature due to the built-in transbilayer asymmetry of the membrane in terms of lipid shape is studied by the proposed method. A linear dependence of the spontaneous curvature with respect to the head-bead diameter difference and the lipid mixing ratio is discovered. The consistency with the theoretical results provides evidence supporting the validity of our method.

In this paper an explicit finite-difference time-domain scheme for solving the Maxwell's equations in non-staggered grids is presented. The proposed scheme for solving the Faraday's and Ampère's equations in a theoretical manner is aimed to preserve discrete zero-divergence for the electric and magnetic fields. The inherent local conservation laws in Maxwell's equations are also preserved discretely all the time using the explicit second-order accurate symplectic partitioned Runge-Kutta scheme. The remaining spatial derivative terms in the semi-discretized Faraday's and Ampère's equations are then discretized to provide an accurate mathematical dispersion relation equation that governs the numerical angular frequency and the wavenumbers in two space dimensions. To achieve the goal of getting the best dispersive characteristics, we propose a fourth-order accurate space centered scheme which minimizes the difference between the exact and numerical dispersion relation equations. Through the computational exercises, the proposed dual-preserving solver is computationally demonstrated to be efficient for use to predict the long-term accurate Maxwell's solutions.

High-order interpolation algorithms for charge conservation in Particle-in-Cell (PIC) simulations are presented. The methods are valid for the case that a particle trajectory is a zigzag line. The second-order and third-order algorithms which can be applied to any even-order and odd-order are discussed in this paper, respectively. Several test simulations are performed to demonstrate their validity in two-dimensional PIC code. Compared with the simulation results of one-order, high-order algorithms have advantages in computation precision and enlarging the grid scales which reduces the CPU time. 

A comparative study is conducted to evaluate three types of lattice Boltzmann equation (LBE) models for fluid flows with finite-sized particles, including the lattice Bhatnagar-Gross-Krook (BGK) model, the model proposed by Ladd [Ladd AJC, J. Fluid Mech., 271, 285-310 (1994); Ladd AJC, J. Fluid Mech., 271, 311-339 (1994)], and the multiple-relaxation-time (MRT) model. The sedimentation of a circular particle in a two-dimensional infinite channel under gravity is used as the first test problem. The numerical results of the three LBE schemes are compared with the theoretical results and existing data. It is found that all of the three LBE schemes yield reasonable results in general, although the BGK scheme and Ladd's scheme give some deviations in some cases. Our results also show that the MRT scheme can achieve a better numerical stability than the other two schemes. Regarding the computational efficiency, it is found that the BGK scheme is the most superior one, while the other two schemes are nearly identical. We also observe that the MRT scheme can unequivocally reduce the viscosity dependence of the wall correction factor in the simulations, which reveals the superior robustness of the MRT scheme. The superiority of the MRT scheme over the other two schemes is also confirmed by the simulation of the sedimentation of an elliptical particle.

The polynomial chaos expansion (PCE) is an efficient numerical method for performing a reliability analysis. It relates the output of a nonlinear system with the uncertainty in its input parameters using a multidimensional polynomial approximation (the so-called PCE). Numerically, such an approximation can be obtained by using a regression method with a suitable design of experiments. The cost of this approximation depends on the size of the design of experiments. If the design of experiments is large and the system is modeled with a computationally expensive FEA (Finite Element Analysis) model, the PCE approximation becomes unfeasible. The aim of this work is to propose an algorithm that generates efficiently a design of experiments of a size defined by the user, in order to make the PCE approximation computationally feasible. It is an optimization algorithm that seeks to find the best design of experiments in the D-optimal sense for the PCE. This algorithm is a coupling between genetic algorithms and the Fedorov exchange algorithm. The efficiency of our approach in terms of accuracy and computational time reduction is compared with other existing methods in the case of analytical functions and finite element based functions.