An entropy stable fully discrete shock capturing space-time Discontinuous Galerkin (DG) method was proposed in a recent paper to approximate hyperbolic systems of conservation laws. This numerical scheme involves the solution of a very large nonlinear system of algebraic equations, by a Newton-Krylov method, at every time step. In this paper, we design efficient preconditioners for the large, non-symmetric linear system, that needs to be solved at every Newton step. Two sets of preconditioners, one of the block Jacobi and another of the block Gauss-Seidel type are designed. Fourier analysis of the preconditioners reveals their robustness and a large number of numerical experiments are presented to illustrate the gain in efficiency that results from preconditioning. The resulting method is employed to compute approximate solutions of the compressible Euler equations, even for very high CFL numbers.

Removing geometric details from the computational domain can significantly reduce the complexity of downstream task of meshing and simulation computation, and increase their stability. Proper estimation of the sensitivity analysis error induced by removing such domain details, called defeaturing errors, can ensure that the sensitivity analysis fidelity can still be met after simplification. In this paper, estimation of impacts of removing arbitrarily constrained domain details to the analysis of incompressible fluid flows is studied with applications to fast analysis of incompressible fluid flows in complex environments. The derived error estimator is applicable to geometric details constrained by either Dirichlet or Neumann boundary conditions, and has no special requirements on the outer boundary conditions. Extensive numerical examples were presented to demonstrate the effectiveness and efficiency of the proposed error estimator.

Using the incompressible isotropic turbulent fields obtained from direct numerical simulation and large-eddy simulation, we studied the statistics of oscillation structures based on local zero-crossings and their relation with inertial-range intermittency for transverse velocity and passive scalar. Our results show that for both the velocity and passive scalar, the local oscillation structures are statistically scale-invariant at high Reynolds number, and the inertial-range intermittency of the overall flow region is determined by the most intermittent structures characterized by one local zero-crossing. Local flow patterns conditioned on the oscillation structures are characterized by the joint probability density function of the invariants of the filtered velocity gradient tensor at inertial range. We demonstrate that the most intermittent regions for longitudinal velocity tend to lay at the saddle area, while those for the transverse velocity tend to locate at the vortex-dominated area. The connection between the ramp-cliff structures in passive scalar field and the corresponding saddle regions in the velocity field is also verified by the approach of oscillation structure classification.

Numerical methods based on gyrocenter gauge kinetic theory are suitable for first principle simulations of high frequency waves in magnetized plasmas. The δf gyrocenter gauge PIC simulation for linear rf wave has been previously realized. In this paper we further develop a full-f nonlinear PIC algorithm appropriate for the nonlinear physics of high frequency waves in magnetized plasmas. Numerical cases of linear rf waves are calculated as a benchmark for the nonlinear GyroGauge code, meanwhile nonlinear rf-wave phenomena are studied. The technique and advantage of the reduction of the numerical noise in this full-f gyrocenter gauge PIC algorithm are also discussed.

We develop a continuum hydrodynamic model for two-phase immiscible flows that involve electroosmotic effect in an electrolyte and moving contact line at solid surfaces. The model is derived through a variational approach based on the On-sager principle of minimum energy dissipation. This approach was first presented in the derivation of a continuum hydrodynamic model for moving contact line in neutral two-phase immiscible flows (Qian, Wang, and Sheng, J. Fluid Mech. 564, 333-360 (2006)). Physically, the electroosmotic effect can be formulated by the Onsager principle as well in the linear response regime. Therefore, the same variational approach is applied here to the derivation of the continuum hydrodynamic model for charged two-phase immiscible flows where one fluid component is an electrolyte exhibiting electroosmotic effect on a charged surface. A phase field is employed to model the diffuse interface between two immiscible fluid components, one being the electrolyte and the other a nonconductive fluid, both allowed to slip at solid surfaces. Our model consists of the incompressible Navier-Stokes equation for momentum transport, the Nernst-Planck equation for ion transport, the Cahn-Hilliard phase-field equation for interface motion, and the Poisson equation for electric potential, along with all the necessary boundary conditions. In particular, all the dynamic boundary conditions at solid surfaces, including the generalized Navier boundary condition for slip, are derived together with the equations of motion in the bulk region. Numerical examples in two-dimensional space, which involve overlapped electric double layer fields, have been presented to demonstrate the validity and applicability of the model, and a few salient features of the two-phase immiscible electroosmotic flows at solid surface. The wall slip in the vicinity ofmoving contact line and the Smoluchowski slip in the electric double layer are both investigated.

The design optimization and analysis of charged particle beam systems employing intense beams requires a robust and accurate Poisson solver. This paper presents a new type of Poisson solver which allows the effects of space charge to be elegantly included into the system dynamics. This is done by casting the charge distribution function into a series of basis functions, which are then integrated with an appropriate Green's function to find a Taylor series of the potential at a given point within the desired distribution region. In order to avoid singularities, a Duffy transformation is applied, which allows singularity-free integration and maximized convergence region when performed with the help of Differential Algebraic methods. The method is shown to perform well on the examples studied. Practical implementation choices and some of their limitations are also explored.

Applicability of Feynman path integral approach to numerical simulations of quantum dynamics of an electron in real time domain is examined. Coherent quantum dynamics is demonstrated with one dimensional test cases (quantum dot models) and performance of the Trotter kernel as compared with the exact kernels is tested. Also, a novel approach for finding the ground state and other stationary sates is presented. This is based on the incoherent propagation in real time. For both approaches the Monte Carlo grid and sampling are tested and compared with regular grids and sampling. We asses the numerical prerequisites for all of the above.

In this paper, we will develop a fast iterative solver for the system of linear equations arising from the local discontinuous Galerkin (LDG) spatial discretization and additive Runge-Kutta (ARK) time marching method for the KdV type equations. Being implicit in time, the severe time step , with the k-th order of the partial differential equations (PDEs)) restriction for explicit methods will be removed. The equations at the implicit time level are linear and we demonstrate an efficient, practical multigrid (MG) method for solving the equations. In particular, we numerically show the optimal or sub-optimal complexity of the MG solver and a two-level local mode analysis is used to analyze the convergence behavior of the MG method. Numerical results for one-dimensional, two-dimensional and three-dimensional cases are given to illustrate the efficiency and capability of the LDG method coupled with the multigrid method for solving the KdV type equations.

Acupuncture has been in clinical practice in China for thousands of years and its analgesia effect is worldwide accepted. However, the mechanism of acupuncture effect is not well understood. The study focus on signaling pathways induced by acupuncture, analyzes the cooperative action of the acupoints’ structure and the associated chemical mediators during acupuncture, establishes a mathematical model clarifying the roadmap of electroneurographic signal startup and transmission mechanism induced by acupuncture, quantitatively analyzing the response in acupoints to acupuncture. These work contribute to reveal the activation and transmission mechanism of neural signals induced by acupuncture from systems biology perspective, lay the foundation for the integration of acupuncture theory and modern science and further guide the clinical treatment and experimental research of acupuncture.

Defects arise when nematic liquid crystals are under topological constraints at the boundary. Recently the study of defects has drawn a lot of attention because of the growing theoretical and practical significance. In this paper, we investigate the relationship between two-dimensional defects and three-dimensional defects within nematic liquid crystals confined in a shell. A highly accurate spectral method is used to solve the Landau-de Gennes model to get the detailed static structures of defects. Interestingly, the solution is radial-invariant when the thickness of the shell is sufficiently small. As the shell thickness increases, the solution undergoes symmetry break to reconfigure the disclination lines. We study this three-dimensional reconfiguration of disclination lines in detail under different boundary conditions. In particular, we find that the temperature plays an important role in deciding whether the transition between two-dimensional defects and three-dimensional defects is continuous or discontinuous for the shell with planar anchoring condition on both inner and outer surfaces. We also discuss the characterization of defects in two- and three-dimensional spaces within the tensor model.

Extended hydrodynamic models for carrier transport are derived from the semiconductor Boltzmann equation with relaxation time approximation of the scattering term, by using the globally hyperbolic moment method and the moment-dependent relaxation time. Incorporating the microscopic relaxation time and the applied voltage bias, a formula is proposed to determine the relaxation time for each moment equation, which sets different relaxation rates for different moments such that higher moments damp faster. The resulting models would give more satisfactory results of macroscopic quantities of interest with a high-order convergence to those of the underlying Boltzmann equation as the involved moments increase, in comparison to the corresponding moment models using a single relaxation time. In order to simulate the steady states efficiently, a multigrid solver is developed for the derived moment models. Numerical simulations of an n +-n-n + silicon diode are carried out to demonstrate the validation of the presented moment models, and the robustness and efficiency of the designed multigrid solver.

The formation of NaH molecules in the lowest rovibrational level of the ground electronic state is investigated using a pump-dump photoassociation (PA) scheme. In short-range region, two colliding atoms Na and H are efficiently associated into the NaH molecule in the rovibrational |0,0〉 state of the ground electronic state via the intermediately rovibrational |10,1〉 state of the excited electronic state. The changes of populations with the electric field amplitudes, frequency detunings, dump pulse duration and delay time between two laser pulses are calculated and discussed. The PA probability reaches 0.623 with a high state-selectivity.

A two-dimensional numerical scheme for the compressible Euler equations is presented and applied here to the simulation of exemplary compressible vortical flows. The proposed approach allows to perform computations on unstructured moving grids with adaptation, which is required to capture complex features of the flow-field. Grid adaptation is driven by suitable error indicators based on the Mach number and by element-quality constraints as well. At the new time level, the computational grid is obtained by a suitable combination of grid smoothing, edge-swapping, grid refinement and de-refinement. The grid modifications—including topology modification due to edge-swapping or the insertion/deletion of a new grid node—are interpreted at the flow solver level as continuous (in time) deformations of suitably-defined node-centered finite volumes. The solution over the new grid is obtained without explicitly resorting to interpolation techniques, since the definition of suitable interface velocities allows one to determine the new solution by simple integration of the Arbitrary Lagrangian-Eulerian formulation of the flow equations. Numerical simulations of the steady oblique-shock problem, of the steady transonic flow and of the start-up unsteady flow around the NACA 0012 airfoil are presented to assess the scheme capabilities to describe these flows accurately.

We consider multi-physics computations where the Navier-Stokes equations of compressible fluid flow on some parts of the computational domain are coupled to the equations of elasticity on other parts of the computational domain. The different subdomains are separated by well-defined interfaces. We consider time accurate computations resolving all time scales. For such computations, explicit time stepping is very efficient. We address the issue of discrete interface conditions between the two domains of different physics that do not lead to instability, or to a significant reduction of the stable time step size. Finding such interface conditions is non-trivial.

We discretize the problem with high order centered difference approximations with summation by parts boundary closure. We derive L2 stable interface conditions for the linearized one dimensional discretized problem. Furthermore, we generalize the interface conditions to the full non-linear equations and numerically demonstrate their stable and accurate performance on a simple model problem. The energy stable interface conditions derived here through symmetrization of the equations contain the interface conditions derived through normal mode analysis by Banks and Sjögreen in [8] as a special case.

Nonlinear filter problems arise in many applications such as communications and signal processing. Commonly used numerical simulation methods include Kalman filter method, particle filter method, etc. In this paper a novel numerical algorithm is constructed based on samples of the current state obtained by solving the state equation implicitly. Numerical experiments demonstrate that our algorithm is more accurate than the Kalman filter and more stable than the particle filter.

In this article, we present two new novel finite difference approximations of order two and four, respectively, for the three dimensional non-linear triharmonic partial differential equations on a compact stencil where the values of u, ∂2u/∂n2 and ∂4u/∂n4 are prescribed on the boundary. We introduce new ideas to handle the boundary conditions and there is no need to discretize the derivative boundary conditions. We require only 7- and 19-grid points on the compact cell for the second and fourth order approximation, respectively. The Laplacian and the biharmonic of the solution are obtained as by-product of the methods. We require only system of three equations to obtain the solution. Numerical results are provided to illustrate the usefulness of the proposed methods.

Numerical oscillation of the total energy can be observed when the Kohn- Sham equation is solved by real-space methods to simulate the translational move of an electronic system. Effectively remove or reduce the unphysical oscillation is crucial not only for the optimization of the geometry of the electronic structure, but also for the study of molecular dynamics. In this paper, we study such unphysical oscillation based on the numerical framework in [G. Bao, G. H. Hu, and D. Liu, An h-adaptive finite element solver for the calculations of the electronic structures, Journal of Computational Physics, Volume 231, Issue 14, Pages 4967–4979, 2012], and deliver some numerical methods to constrain such unphysical effect for both pseudopotential and all-electron calculations, including a stabilized cubature strategy for Hamiltonian operator, and an a posteriori error estimator of the finite element methods for Kohn-Sham equation. The numerical results demonstrate the effectiveness of our method on restraining unphysical oscillation of the total energies.

The lobe dynamics andmass transport between separation bubble and main flow in flow over airfoil are studied in detail, using Lagrangian coherent structures (LCSs), in order to understand the nature of evolution of the separation bubble. For this problem, the transient flow over NACA0012 airfoil with low Reynolds number is simulated numerically by characteristic based split (CBS) scheme, in combination with dual time stepping. Then, LCSs and lobe dynamics are introduced and developed to investigate themass transport between separation bubble and main flow, from viewpoint of nonlinear dynamics. The results show that stable manifolds and unstable manifolds could be tangled with each other as time evolution, and the lobes are formed periodically to induce mass transport between main flow and separation bubble, with dynamic behaviors. Moreover, the evolution of the separation bubble depends essentially on the mass transport which is induced by lobes, ensuing energy and momentum transfers. As the results, it can be drawn that the dynamics of flow separation could be studied using LCSs and lobe dynamics, and could be controlled feasibly if an appropriate control is applied to the upstream boundary layer with high momentum.

In the present work adaptation in meshless framework is proposed. The grid adaptation or mesh adaptation is quite well developed area in case of conventional grid based solvers and is popularly known as Adaptive mesh refinement (AMR). In such cases the adaptation is done by subdividing the cells or elements into finer cells or elements. In case of meshless methods there are no cells or elements but only a cloud of points. In this work we propose to achieve the meshless adaptation by locally refining the point density in the regions demanding higher resolution. This results into an adaptive enriched cloud of points. We call this method as Adaptive Cloud Refinement (ACR). The meshless solvers need connectivity information, which is a set of neighboring nodes. It is crucial part of meshless solvers. Obviously because of refining point density, the connectivity of nodes in such regions gets modified and hence has to be updated. An efficient connectivity update must exploit the fact that the node distribution would be largely unaffected except the region of adaptation. Hence connectivity updating needs to be done locally, only in these regions. In this paper we also present an extremely fast algorithm to update connectivity over adapted cloud called as ACU (Automatic Connectivity Update).

Mesh generation is a bottleneck for finite element simulations of biomolecules. A robust and efficient approach, based on the immersed boundary method proposed in [8], has been developed and implemented to generate large-scale mesh body-fitted to molecular shape for general parallel finite element simulations. The molecular Gaussian surface is adopted to represent the molecular surface, and is finally approximated by piecewise planes via the tool phgSurfaceCut in PHG [43], which is improved and can reliably handle complicated molecular surfaces, through mesh refinement steps. A coarse background mesh is imported first and then is distributed into each process using a mesh partitioning algorithm such as space filling curve [5] or METIS [22]. A bisection method is used for the mesh refinements according to the molecular PDB or PQR file which describes the biomolecular region. After mesh refinements, the mesh is optionally repartitioned and redistributed for load balancing. For finite element simulations, the modification of region mark and boundary types is done in parallel. Our parallel mesh generation method has been successfully applied to a sphere cavity model, a DNA fragment, a gramicidin A channel and a huge Dengue virus system. The results of numerical experiments show good parallel efficiency. Computations of electrostatic potential and solvation energy also validate the method. Moreover, the meshing process and adaptive finite element computation can be integrated as one PHG project to avoid the mesh importing and exporting costs, and improve the convenience of application as well.