Computational simulations of multiphase flow are challenging because many practical applications require adequate resolution of not only interfacial physics associated with moving boundaries with possible topological changes, but also around three-dimensional, irregular solid geometries. In this paper, we highlight recent efforts made in simulating multiphase fluid dynamics around complex geometries, based on an Eulerian-Lagrangian framework. The approach uses two independent but related grid layouts to track the interfacial and solid boundary conditions, and is capable of capturing interfacial as well as multiphase dynamics. In particular, the stationary Cartesian grid with time dependent, local adaptive refinement is utilized to handle the computation of the transport equations, while the interface shape and movement are treated by marker-based triangulated surface meshes which freely move and interact with the Cartesian grid. The markers are also used to identify the location of solid boundaries and enforce the no-slip condition there. Issues related to the contact line treatment, topological changes of multiphase fronts during merger or breakup of objects, and necessary data structures and solution techniques are also highlighted. Selected test cases including spacecraft fuel tank flow management and liquid plug flow dynamics are presented.

Radial basis functions (RBFs) can be used to approximate derivatives and solve differential equations in several ways. Here, we compare one important scheme to ordinary finite differences by a mixture of numerical experiments and theoretical Fourier analysis, that is, by deriving and discussing analytical formulas for the error in differentiating exp(ikx) for arbitrary k. "Truncated RBF differences" are derived from the same strategy as Fourier and Chebyshev pseudospectral methods: Differentiation of the Fourier, Chebyshev or RBF interpolant generates a differentiation matrix that maps the grid point values or samples of a function u(x) into the values of its derivative on the grid. For Fourier and Chebyshev interpolants, the action of the differentiation matrix can be computed indirectly but efficiently by the Fast Fourier Transform (FFT). For RBF functions, alas, the FFT is inapplicable and direct use of the dense differentiation matrix on a grid of N points is prohibitively expensive (O(N²)) unless N is tiny. However, for Gaussian RBFs, which are exponentially localized, there is another option, which is to truncate the dense matrix to a banded matrix, yielding “truncated RBF differences”. The resulting formulas are identical in form to finite differences except for the difference weights. On a grid of spacing h with the RBF as φ(x)=exp(−α²(x/h)²),
 

where without approximation w_m=(−1)^m+12α²/sinh(mα²). We derive explicit formula for the differentiation of the linear function, f(X)≡X, and the errors therein. We show that Gaussian radial basis functions (GARBF), when truncated to give differentiation formulas of stencil width (2M+1), are significantly less accurate than (2M)-th order finite differences of the same stencil width. The error of the infinite series (M=∞) decreases exponentially as α→0. However, truncated GARBF series have a second error (truncation error) that grows exponentially as α → 0. Even for α ∼ O(1) where the sum of these two errors is minimized, it is shown that the finite difference formulas are always superior. We explain, less rigorously, why these arguments extend to more general species of RBFs and to an irregular grid. There are, however, a variety of alternative differentiation strategies which will be analyzed in future work, so it is far too soon to dismiss RBFs as a tool for solving differential equations.

We perform numerical simulations of the sedimentation of rigid fibers suspended in a viscous incompressible fluid at nonzero Reynolds numbers. The fiber sedimentation system is modeled as a two-dimensional immersed boundary problem, which naturally accommodates the fluid-particle interactions and which allows the simulation of a large number of suspending fibers. We study the dynamics of sedimenting fibers under a variety of conditions, including differing fiber densities, Reynolds numbers, domain boundary conditions, etc. Simulation results are compared to experimental measurements and numerical results obtained in previous studies.

A stable and explicit second order accurate finite difference method for the elastic wave equation in curvilinear coordinates is presented. The discretization of the spatial operators in the method is shown to be self-adjoint for free-surface, Dirichlet and periodic boundary conditions. The fully discrete version of the method conserves a discrete energy to machine precision. 

We consider two existing FFT-based fast-convolution iterative solution techniques for the scalar T-matrix multiple-scattering equation [1]. The use of the FFT operation requires field values be expressed on a regular Cartesian grid and the two techniques differ in how to go about achieving this. The first technique [6, 7] uses the non-diagonal translation operator [1, 9] of the spherical multipole field, while the second method [11] uses the diagonal translation operator of Rokhlin [10]. Because of its use of the non-diagonal translator, the first technique has been thought to require a greater number of spatial convolutions than the second technique. We establish that the first method requires only half as many convolution operations as the second method for a comparable numerical accuracy and demonstrate, based on an actual CPU time comparison, that it can therefore perform iterations faster than the second method. We then consider the respective symmetry relations of the non-diagonal and diagonal translators and discuss a memory-reduction procedure for both FFT-based methods. In this procedure, we need to store only the minimum sets of near-field and far-field translation operators and generate missing elements on the fly using the symmetry relations. We show that the relative cost of generating the missing elements becomes smaller as the number of scatterers increases.

A hybrid finite element-Laplace transform method is implemented to analyze the time domain electromagnetic scattering induced by a 2-D overfilled cavity embedded in the infinite ground plane. The algorithm divides the whole scattering domain into two, interior and exterior, sub-domains. In the interior sub-domain which covers the cavity, the problem is solved via the finite element method. The problem is solved analytically in the exterior sub-domain which slightly overlaps the interior subdomain and extends to the rest of the upper half plane. The use of the Laplace transform leads to an analytical link condition between the overlapping sub-domains. The analytical link guides the selection of the overlapping zone and eliminates the need to use the conventional Schwartz iteration. This dramatically improves the efficiency for solving transient scattering problems. Numerical solutions are tested favorably against analytical ones for a canonical geometry. The perfect link over the artificial boundary between the finite element approximation in the interior and analytical solution in the exterior further indicates the reliability of the method. An error analysis is also performed.

In this study, a three-dimensional artificial compressibility solver based on the average-state Harten-Lax-van Leer-Contact (HLLC) [13] type Riemann solution is first proposed and developed to solve the time-dependent incompressible flow equations. To implement unsteady flow calculations, a dual time stepping strategy including the LU decomposition method is used in the pseudo-time iteration and the second-order accurate backward difference is adopted to discretize the unsteady flow term. Also a third-order accurate HLLC numerical flux is derived for approximating the inviscid terms. To verify numerical accuracy, flows over a lid-driven cavity and an oscillating flat plate are chosen as the benchmark tests. In addition, the current solver is extended to solve blood flows in a realistic human aorta measured from MRI (Magnetic Resonance Imaging). The simulation geometry was derived from a three-dimensional reconstruction of a series of two-dimensional slices obtained in vivo. Numerical results demonstrate wall stresses were highly dynamic, but were generally high along the outer wall in the vicinity of the branches and low along the inner wall, particularly in the descending aorta. The maximum wall stress distribution is presented on the aortic arch in the systole. In addition, extensive counter-clockwise secondary flows and three-dimensional helical vortex influenced considerably by the presence of vessel contraction, torsion and the branches were shown in the descending aorta in the late systole and early diastolic cycles.

In this paper we discuss the inverse scattering algorithm for predicting internal multiple reflections (reverberation artefacts), focusing our attention on the construction mechanisms. Roughly speaking, the algorithm combines amplitude and phase information of three different arrivals (sub-events) in the data set to predict one multiple reflection. The three events are conditioned by a certain relation which requires that their pseudo-depths, defined as the depths of their turning points relative to the constant background velocity, satisfy a lower-higher-lower relationship. This implicitly assumes a pseudo-depth monotonicity condition, i.e., the relation between the actual depths and the pseudo-depths of any two sub-events is the same. We study this relation in pseudo-depth and show that it is directly connected with a similar relation between the vertical or intercept times of the sub-events. The paper also provides the first multidimensional analysis of the algorithm (for a vertically varying acoustic model) with analytical data. We show that the construction of internal multiples is performed in the plane waves domain and, as a consequence, the internal multiples with headwaves sub-events are also predicted by the algorithm. Furthermore we analyze the differences between the time monotonicity condition in vertical or intercept time and total travel time and show a 2D example which satisfies the former but not the latter. Finally we discuss one case in which the monotonicity condition is not satisfied by the sub-events of an internal multiple and discuss ways of lowering these restrictions and of expanding the algorithm to address these types of multiples.

The full torus electromagnetic gyrokinetic particle simulations using the hybrid model with kinetic electrons in the presence of magnetic shear is presented. The fluid-kinetic electron hybrid model employed in this paper improves numerical properties by removing the tearing mode, meanwhile, preserves both linear and nonlinear wave-particle resonances of electrons with Alfven wave and ion acoustic wave.