We review, implement, and compare numerical integration schemes for spatially bounded diffusions stopped at the boundary which possess a convergence rate of the discretization error with respect to the timestep h higher than . We address specific implementation issues of the most general-purpose of such schemes. They have been coded into a single Matlab program and compared, according to their accuracy and computational cost, on a wide range of problems in up to ℝ48. The paper is self-contained and the code will be made freely downloadable.

We construct a new first-order central-upwind numerical method for solving systems of hyperbolic equations in conservative form. It applies in multidimensional structured and unstructured meshes. The proposed method is an extension of the UFORCE method developed by Stecca, Siviglia and Toro, in which the upwind bias for the modification of the staggered mesh is evaluated taking into account the smallest and largest wave of the entire Riemann fan. The proposed first-order method is shown to be identical to the Godunov upwind method in applications to a 2 x 2 linear hyperbolic system. The method is then extended to non-linear systems and its performance is assessed by solving the two-dimensional inviscid shallow water equations. Extension to second-order accuracy is carried out using an ADER-WENO approach in the finite volume framework on unstructured meshes. Finally, numerical comparison with current competing numerical methods enables us to identify the salient features of the proposed method.

Many physical processes are described by elliptic or parabolic partial differential equations. For linear stability problems associated with such equations, the inverse Laplacian provides a very effective preconditioner. In addition, it is also readily available in most scientific calculations in the form of a Poisson solver or an implicit diffusive timestep. We incorporate Laplacian preconditioning into the inverse Arnoldi method, using BiCGSTAB to solve the large linear systems. Two successful implementations are described: spherical Couette flow described by the Navier-Stokes equations and Bose-Einstein condensation described by the nonlinear Schrödinger equation.

This paper presents a new approach to verify the accuracy of computational simulations. We develop mathematical theorems which can serve as robust a posteriori error estimation techniques to identify numerical pollution, check the performance of adaptive meshes, and verify numerical solutions. We demonstrate performance of this methodology on problems from flow thorough porous media. However, one can extend it to other models. We construct mathematical properties such that the solutions to Darcy and Darcy-Brinkman equations satisfy them. The mathematical properties include the total minimum mechanical power, minimum dissipation theorem, reciprocal relation, and maximum principle for the vorticity. All the developed theorems have firm mechanical bases and are independent of numerical methods. So, these can be utilized for solution verification of finite element, finite volume, finite difference, lattice Boltzmann methods and so forth. In particular, we show that, for a given set of boundary conditions, Darcy velocity has the minimum total mechanical power of all the kinematically admissible vector fields. We also show that a similar result holds for Darcy-Brinkman velocity. We then show for a conservative body force, the Darcy and Darcy-Brinkman velocities have the minimum total dissipation among their respective kinematically admissible vector fields. Using numerical examples, we show that the minimum dissipation and total mechanical power theorems can be utilized to identify pollution errors in numerical solutions. The solutions to Darcy and Darcy-Brinkman equations are shown to satisfy a reciprocal relation, which has the potential to identify errors in the numerical implementation of boundary conditions. It is also shown that the vorticity under both steady and transient Darcy-Brinkman equations satisfy maximum principles if the body force is conservative and the permeability is homogeneous and isotropic. A discussion on the nature of vorticity under steady and transient Darcy equations is also presented. Using several numerical examples, we will demonstrate the predictive capabilities of the proposed a posteriori techniques in assessing the accuracy of numerical solutions for a general class of problems, which could involve complex domains and general computational grids.

We propose a multi-layer approach to simulate hyperpycnal and hypopycnal plumes in flows with free surface. The model allows to compute the vertical profile of the horizontal and the vertical components of the velocity of the fluid flow. The model can describe as well the vertical profile of the sediment concentration and the velocity components of each one of the sediment species that form the turbidity current. To do so, it takes into account the settling velocity of the particles and their interaction with the fluid. This allows to better describe the phenomena than a single layer approach. It is in better agreement with the physics of the problem and gives promising results. The numerical simulation is carried out by rewriting the multilayer approach in a compact formulation, which corresponds to a system with nonconservative products, and using path-conservative numerical scheme. Numerical results are presented in order to show the potential of the model.

A nonlocal continuum electrostatic model, defined as integro-differential equations, can significantly improve the classic Poisson dielectric model, but is too costly to be applied to large protein simulations. To sharply reduce the model’s complexity, a modified nonlocal continuum electrostatic model is presented in this paper for a protein immersed in water solvent, and then transformed equivalently as a system of partial differential equations. By using this new differential equation system, analytical solutions are derived for three different nonlocal ionic Born models, where a monoatomic ion is treated as a dielectric continuum ball with point charge either in the center or uniformly distributed on the surface of the ball. These solutions are analytically verified to satisfy the original integro-differential equations, thereby, validating the new differential equation system.

When gas flows through corrugated pipes, pressure waves interacting with vortex shedding can produce distinct tonal noise and structural vibration. Based on established observations, a model is proposed which couples an acoustic pipe and self-excited oscillations with vortex shedding over the corrugation cavities. In the model, the acoustic response of the corrugated pipe is simulated by connecting the lossless medium moving with a constant velocity with a source based on a discrete distribution of van der Pol oscillators arranged along the pipe. Our time accurate solutions exhibit dynamic behavior consistent with that experimentally observed, including the lock-in frequency of vortex shedding, standing waves and the onset fluid velocity capable of generating the lock-in.

A fully discrete discontinuous Galerkin method is introduced for solving time-dependent Maxwell’s equations. Distinguished from the Runge-Kutta discontinuous Galerkin method (RKDG) and the finite element time domain method (FETD), in our scheme, discontinuous Galerkin methods are used to discretize not only the spatial domain but also the temporal domain. The proposed numerical scheme is proved to be unconditionally stable, and a convergent rate is established under the L2-norm when polynomials of degree atmost r and k are used for temporal and spatial approximation, respectively. Numerical results in both 2-D and 3-D are provided to validate the theoretical prediction. An ultra-convergence of order in time step is observed numerically for the numerical fluxes w.r.t. temporal variable at the grid points.

We present DASHMM, a general library implementing multipole methods (including both Barnes-Hut and the Fast Multipole Method). DASHMM relies on dynamic adaptive runtime techniques provided by the HPX-5 system to parallelize the resulting multipole moment computation. The result is a library that is easy-to-use, extensible, scalable, efficient, and portable. We present both the abstractions defined by DASHMM as well as the specific features of HPX-5 that allow the library to execute scalably and efficiently.

This work proposes an extension to Boltzmann BGK equation for dense gases. The present model has an H-theorem and it allows choice of the Prandtl number as an independent parameter. I show that similar to Enskog equation this equation can reproduce dynamics of dense gases.

We are concerned with a model of ionic polymer-metal composite (IPMC) materials that consists of a coupled system of the Poisson and Nernst-Planck equations, discretized by means of the finite element method (FEM). We show that due to the transient character of the problem it is efficient to use adaptive algorithms that are capable of changing the mesh dynamically in time. We also show that due to large qualitative and quantitative differences between the two solution components, it is efficient to approximate them on different meshes using a novel adaptive multimesh hp-FEM. The study is accompanied with numerous computations and comparisons of the adaptive multimesh hp-FEM with several other adaptive FEM algorithms.

In this paper, we present an adaptive, analysis of variance (ANOVA)-based data-driven stochastic method (ANOVA-DSM) to study the stochastic partial differential equations (SPDEs) in the multi-query setting. Our new method integrates the advantages of both the adaptive ANOVA decomposition technique and the data-driven stochastic method. To handle high-dimensional stochastic problems, we investigate the use of adaptive ANOVA decomposition in the stochastic space as an effective dimension-reduction technique. To improve the slow convergence of the generalized polynomial chaos (gPC) method or stochastic collocation (SC) method, we adopt the data-driven stochastic method (DSM) for speed up. An essential ingredient of the DSM is to construct a set of stochastic basis under which the stochastic solutions enjoy a compact representation for a broad range of forcing functions and/or boundary conditions.

Our ANOVA-DSM consists of offline and online stages. In the offline stage, the original high-dimensional stochastic problem is decomposed into a series of low-dimensional stochastic subproblems, according to the ANOVA decomposition technique. Then, for each subproblem, a data-driven stochastic basis is computed using the Karhunen-Loève expansion (KLE) and a two-level preconditioning optimization approach. Multiple trial functions are used to enrich the stochastic basis and improve the accuracy. In the online stage, we solve each stochastic subproblem for any given forcing function by projecting the stochastic solution into the data-driven stochastic basis constructed offline. In our ANOVA-DSM framework, solving the original highdimensional stochastic problem is reduced to solving a series of ANOVA-decomposed stochastic subproblems using the DSM. An adaptive ANOVA strategy is also provided to further reduce the number of the stochastic subproblems and speed up our method. To demonstrate the accuracy and efficiency of our method, numerical examples are presented for one- and two-dimensional elliptic PDEs with random coefficients.

Our recently developed lattice Boltzmann model is used to simulate droplet dynamical behaviour governed by thermocapillary force in microchannels. One key research challenge for developing droplet-based microfluidic systems is control of droplet motion and its dynamic behaviour. We numerically demonstrate that the thermocapillary force can be exploited for microdroplet manipulations including synchronisation, sorting, and splitting. This work indicates that the lattice Boltzmann method provides a promising design simulation tool for developing complex droplet-based microfluidic devices.

Proton transport is one of the most important and interesting phenomena in living cells. The present work proposes a multiscale/multiphysics model for the understanding of the molecular mechanism of proton transport in transmembrane proteins. We describe proton dynamics quantum mechanically via a density functional approach while implicitly model other solvent ions as a dielectric continuum to reduce the number of degrees of freedom. The densities of all other ions in the solvent are assumed to obey the Boltzmann distribution. The impact of protein molecular structure and its charge polarization on the proton transport is considered explicitly at the atomic level. We formulate a total free energy functional to put proton kinetic and potential energies as well as electrostatic energy of all ions on an equal footing. The variational principle is employed to derive nonlinear governing equations for the proton transport system. Generalized Poisson-Boltzmann equation and Kohn-Sham equation are obtained from the variational framework. Theoretical formulations for the proton density and proton conductance are constructed based on fundamental principles. The molecular surface of the channel protein is utilized to split the discrete protein domain and the continuum solvent domain, and facilitate the multiscale discrete/continuum/quantum descriptions. A number of mathematical algorithms, including the Dirichlet to Neumann mapping, matched interface and boundary method, Gummel iteration, and Krylov space techniques are utilized to implement the proposed model in a computationally efficient manner. The Gramicidin A (GA) channel is used to demonstrate the performance of the proposed proton transport model and validate the efficiency of proposed mathematical algorithms. The electrostatic characteristics of the GA channel is analyzed with a wide range of model parameters. The proton conductances are studied over a number of applied voltages and reference concentrations. A comparison with experimental data verifies the present model predictions and validates the proposed model.

In this paper, we compute a phase field (diffuse interface) model of Cahn-Hilliard type for moving contact line problems governing the motion of isothermal multiphase incompressible fluids. The generalized Navier boundary condition proposed by Qian et al. [1] is adopted here. We discretize model equations using a continuous finite element method in space and a modified midpoint scheme in time. We apply a penalty formulation to the continuity equation which may increase the stability in the pressure variable. Two kinds of immiscible fluids in a pipe and droplet displacement with a moving contact line under the effect of pressure driven shear flow are studied using a relatively coarse grid. We also derive the discrete energy law for the droplet displacement case, which is slightly different due to the boundary conditions. The accuracy and stability of the scheme are validated by examples, results and estimate order.

This paper presents a new Lagrangian type scheme for solving the Euler equations of compressible gas dynamics. In this new scheme the system of equations is discretized by Runge-Kutta Discontinuous Galerkin (RKDG) method, and the mesh moves with the fluid flow. The scheme is conservative for the mass, momentum and total energy and maintains second-order accuracy. The scheme avoids solving the geometrical part and has free parameters. Results of some numerical tests are presented to demonstrate the accuracy and the non-oscillatory property of the scheme.

Parallel tempering simulation is widely used in enhanced sampling of systems with complex energy surfaces. We hereby introduce generalized canonical ensemble (GCE) instead of the usual canonical ensemble into the parallel tempering to further improve abilities of the simulation technique. GCE utilizes an adapted weight function to obtain a unimodal energy distribution even in phase-coexisting region and then the parallel tempering on GCE yields the steady swap acceptance rates (SARs) instead of the fluctuated SARs in that on canonical ensemble. With the steady SARs, we can facilitate assign the parameters of the parallel tempering simulation to more efficiently reach equilibrium among different phases. We illustrate the parallel tempering simulation on GCE in the phase-coexisting region of 2-dimensional Potts model, a benchmark system for new simulation method developing. The result indicates that the new parallel tempering method is more efficient to estimate statistical quantities (i.e., to sample the conformational space) than the normal parallel tempering, specially in phase-coexisting regions of larger systems.

A higher-order compact scheme on the nine point 2-D stencil is developed for the steady stream-function vorticity form of the incompressible Navier-Stokes (N-S) equations in spherical polar coordinates, which was used earlier only for the cartesian and cylindrical geometries. The steady, incompressible, viscous and axially symmetric flow past a sphere is used as a model problem. The non-linearity in the N-S equations is handled in a comprehensive manner avoiding complications in calculations. The scheme is combined with the multigrid method to enhance the convergence rate. The solutions are obtained over a non-uniform grid generated using the transformation r=e? while maintaining a uniform grid in the computational plane. The superiority of the higher order compact scheme is clearly illustrated in comparison with upwind scheme and defect correction technique at high Reynolds numbers by taking a large domain. This is a pioneering effort, because for the first time, the fourth order accurate solutions for the problem of viscous flow past a sphere are presented here. The drag coefficient and surface pressures are calculated and compared with available experimental and theoretical results. It is observed that these values simulated over coarser grids using the present scheme are more accurate when compared to other conventional schemes. It has also been observed that the flow separation initially occurred at Re = 21.

In this paper, we study numerically quantized vortex dynamics and their interaction in the two-dimensional (2D) Ginzburg-Landau equation (GLE) with a dimensionless parameter ε > 0 on bounded domains under either Dirichlet or homogeneous Neumann boundary condition. We begin with a review of the reduced dynamical laws for time evolution of quantized vortex centers in GLE and show how to solve these nonlinear ordinary differential equations numerically. Then we present efficient and accurate numerical methods for discretizing the GLE on either a rectangular or a disk domain under either Dirichlet or homogeneous Neumann boundary condition. Based on these efficient and accurate numerical methods for GLE and the reduced dynamical laws, we simulate quantized vortex interaction of GLE with different ε and under different initial setups including single vortex, vortex pair, vortex dipole and vortex lattice, compare them with those obtained from the corresponding reduced dynamical laws, and identify the cases where the reduced dynamical laws agree qualitatively and/or quantitatively as well as fail to agree with those from GLE on vortex interaction. Finally, we also obtain numerically different patterns of the steady states for quantized vortex lattices under the GLE dynamics on bounded domains.

This paper concerns the numerical stability of an eikonal transformation based splitting method which is highly effective and efficient for the numerical solution of paraxial Helmholtz equation with a large wave number. Rigorous matrix analysis is conducted in investigations and the oscillation-free computational procedure is proven to be stable in an asymptotic sense. Simulated examples are given to illustrate the conclusion.