Third-order less oscillatory and less diffusive compact stencil-based upwind schemes, and their applications to incompressible flows and free surface flows

Yokoi, Kensuke

doi:10.1063/5.0112953

We propose novel third-order less oscillatory and less diffusive compact stencil-based upwind schemes for the approximation of the continuity equation. The proposed schemes are based on the constrained interpolation profile-conservative semi-Lagrangian schemes. An important feature of the proposed schemes is that the interpolation functions are constructed using only variables within one upwind cell (a cell average and two boundary values). Furthermore, the proposed schemes have third-order accuracy and are also less oscillatory, less diffusive, and fully conservative. The proposed schemes are validated through various benchmark problems and comparisons with experiments of two droplets collision/separation and droplet splashing. The numerical results have shown that the proposed schemes have third-order accuracy for smooth solution, and capture discontinuities and smooth solutions simultaneously without numerical oscillations. The proposed schemes can capture the secondary vorticity of lid-driven cavity flow of Re = 7500 with a Cartesian grid of 64 × 64. The numerical results of two droplets collision/separation of We = 40 show that the proposed schemes can reproduce droplets collision/separation with quite coarse grids. These numerical results of droplet splashing have demonstrated that proposed schemes can reduce numerical diffusions well against existing schemes and robust.

I. INTRODUCTION

In this paper, we propose novel third-order less oscillatory and less diffusive compact stencil-based upwind schemes for the approximation of the continuity equation

\frac{\partial ϕ}{\partial t} + \frac{\partial (u ϕ)}{\partial x} = 0,

(1)

where

ϕ

is the scalar and u is the velocity. The proposed method is based on the constrained interpolation profile-conservative semi-Lagrangian (CIP-CSL) schemes.^{15,16,18,27,34–36}

The CIP-CSL schemes are solvers of the conservation equation based on a multi-moment concept, which uses both cell average and boundary value as variables. In CIP-CSL schemes, the order of accuracy can be improved by using both cell average and boundary value without using wider stencils (e.g., the CIP-CSL schemes can achieve third-order accuracy only with variables within a cell). On the other hand, most of other schemes such as the essentially non-oscillatory schemes (ENO),^9,22,23 weighted essentially non-oscillatory schemes (WENO),^13,17 and the piecewise parabolic method (PPM)⁶ use only single-moment (a type of variables, for instance, cell average if finite volume methods and point value if finite difference methods). In these single-moment-based schemes, a wider stencil must be used to increase the order of accuracy. Thus, the CIP-CSL schemes are a different type of schemes from finite volume or finite different schemes.

The CIP-CSL schemes have been applied to various fluid problems^1,30,31,33 such as droplets collision/separation,⁴⁴ droplet impact on dry surface,⁴⁶ droplet splashing on dry surface,^41–44 interactions between particles and free surface flow,⁴⁴ and oceanic flows.³⁷ The fluid solvers^30,31,33,44 based on the CIP-CSL schemes are highly efficient especially for incompressible flow simulations. For instance, the fluid solvers can simulate droplet splashing on a hydrophobic substrate at least qualitatively with a coarse grid (192 × 192 × 48) and the numerical simulation was completed within 2 hours on a standard desktop computer (Intel Core i7–3820 3.6 GHz, 8 GB memory).⁴⁴ As drawbacks of the CIP-CSL schemes against single-moment-based schemes, the schemes require additional memory and calculation time for these additional degrees of freedom (DoF) because CIP-CSL schemes use both cell average and boundary value as variables. Although the drawback on additional memory cannot be resolved, these additional DoF surely improve fluid calculations (i.e., numerical resolution is increased by the additional DoF) and memory is not so expensive recently. If the calculation time is also proportionally increased due to the additional DoF, the use of the CIP-CSL schemes would not have many advantages against single-moment-based schemes because it is almost like just increased numerical resolution of single-moment-based schemes. However, if the CIP-CSL schemes are used within semi-implicit formulations (the pressure is solved implicitly and velocity explicitly) such as VSIAM3^30,31,33 or VSIAM3-TM,⁴⁴ the computational time is not increased proportionally. Such semi-implicit formulations are mainly used for incompressible flow simulations. In VSIAM3 and VSIAM3-TM, these additional DoF are used only for velocity and not used for pressure. Then, although the calculation cost of velocity is surely increased, the cost of pressure calculation is identical to that of single-moment-based schemes. As it is well known, in incompressible flow simulations based on the semi-implicit formulations, the pressure calculation dominates the total calculation time (typically >80%). Therefore, additional calculation time by the additional DoF is small for the total calculation time. However, the accuracy of velocity is meaningfully improved due to the additional DoF. If the CIP-CSL schemes are used for full explicit formulations (e.g., compressible flow simulations), the calculation cost is proportionally increased as DoF and the advantages are not so clear compared to incompressible flow cases. Therefore, in this paper, we focus only on incompressible flows as applications of the proposed schemes.

There are several compact stencil-based CIP-CSL schemes such as CIP-CSL2 (CIP-CSL with second-order polynomial function)³⁶ and CIP-CSLR (CIP-CSL with rational function).³⁵ In CSL2, which is based on a second-order polynomial interpolation function, three moments within the upwind cell (i.e., a cell average and two boundary values) are used to construct the interpolation function (the details are explained in Subsection II A). CSL2 has third-order accuracy even though the interpolation function is constructed based on moments within one upwind cell. However, as a disadvantage, CSL2 is oscillatory. In CSLR, the same moments with CSL2 are used to construct the interpolation function and oscillations have successfully been minimized by using the rational function. However, CSLR cannot maintain third-order accuracy.

In this paper, we propose the CSL2T (CSL second-order polynomial and hyperbolic tangent functions) and CSL2R (CSL second-order polynomial and rational functions) schemes, which can maintain third-order accuracy of CSL2 but can also minimize oscillations like CSLR. In this paper, we focus on solvers, which can construct interpolation functions using only variables within one upwind cell because such solvers can avoid loss of accuracy near wall boundaries (the proposed schemes can maintain third-order accuracy even near wall boundaries) and practically important. In this paper, we also propose non-oscillatory selector (NOS). NOS is used to switch between an oscillatory high-order schemes such as CSL2 and a non-oscillatory scheme (e.g., CSLT or CSLR). There are some similar concepts with NOS such as limiters,^10,21 boundary variation diminishing (BVD),^7,24 and multi-dimensional optimal order detection (MOOD).⁵ As an advantage of NOS against limiters, NOS can maintain third-order accuracy. Although we also tried some limiters, these limiters did not work well for CSL schemes. This will be because CSL schemes use not only cell averages like finite volume methods but also boundary values. In this paper, we have employed a hyperbolic tangent function^12,32,39 as a non-oscillatory solution like BVD schemes. However, the selection strategies of high-order interpolation function and the hyperbolic tangent function are totally different. Although NOS uses discrete maximum principle to select a suitable interpolation function like MOOD, there are many differences with MOOD. For instance, although MOOD is for high-order finite volume schemes, NOS is for mixed formulations, which uses both finite volume and finite difference methods. MOOD has not employed the hyperbolic tangent function as a candidate of interpolation functions.

The details of proposed CIP-CSL schemes and some existing CIP-CSL schemes are given in Sec. II. In Sec. III, numerical results of sine wave propagation test, Jiang-Shu test, Burgers' equation, Zalesak problem, lid-driven cavity flows, droplet collision and separation, and droplet splashing are given. The summary is given in Sec. IV.

II. NUMERICAL METHODS

In this section, we review CSL2 and CSLR schemes, and propose third-order less oscillatory and less diffusive compact stencil-based upwind schemes (CSL2T and CSL2R with NOS). In Secs. II A and II B, CSL2³⁶ and CSLR³⁵ are explained. In Sec. II C, we propose the CSLT scheme based on the hyperbolic tangent function. In Sec. II D, we propose NOS and the CSL2T scheme based on CSL2 and CSLT schemes. In Sec. II E, we also propose the CSL2R scheme based on CSL2, CSLR, and NOS.

A. CIP-CSL2 scheme

The CIP-CSL2 scheme³⁶ would be the simplest CIP-CSL formulation. The CIP-CSL2 scheme uses three moments in the upwind cell (i.e., one cell average

{\bar{ϕ}}_{i}

and two boundary values

ϕ_{i - 1 / 2}

and

ϕ_{i + 1 / 2}

⁠) to interpolate within the upwind cell (i.e., between

x_{i - 1 / 2}

and

x_{i + 1 / 2}

⁠) as shown in Fig. 1. Then, we can interpolate between

x_{i - 1 / 2}

and

x_{i + 1 / 2}

by a quadratic function

Φ_{i}^{CSL 2} (x)

⁠,

Φ_{i}^{CSL 2} (x) = C_{2, i}^{CSL 2} {(x - x_{i - 1 / 2})}^{2} + C_{1, i}^{CSL 2} (x - x_{i - 1 / 2}) + C_{0, i}^{CSL 2} .

(2)

FIG. 1.

Schematic figure of the CIP-CSL2 scheme. u i − 1 / 2 < 0 is assumed. The moments, which are indicated by gray color ( ϕ i − 1 / 2 , ϕ ¯ i, and ϕ i + 1 / 2), are used to construct the quadratic interpolation function Φ i CSL 2 ( x ).

View large Download slide

Schematic figure of the CIP-CSL2 scheme. $u_{i - 1 / 2} < 0$ is assumed. The moments, which are indicated by gray color (⁠ $ϕ_{i - 1 / 2}, {\bar{ϕ}}_{i}$ ⁠, and $ϕ_{i + 1 / 2}$ ⁠), are used to construct the quadratic interpolation function $Φ_{i}^{CSL 2} (x)$ ⁠.

By using the following constraints,

ϕ_{i - 1 / 2} = Φ_{i}^{CSL 2} (x_{i - 1 / 2}),

(3)

ϕ_{i + 1 / 2} = Φ_{i}^{CSL 2} (x_{i + 1 / 2}),

(4)

{\bar{ϕ}}_{i} = \frac{1}{Δ x} \int_{x_{i - 1 / 2}}^{x_{i + 1 / 2}} Φ_{i}^{CSL 2} (x) d x,

(5)

these coefficients

C_{2, i}^{CSL 2}, C_{1, i}^{CSL 2},

and

C_{0, i}^{CSL 2}

can be calculated as follows:

C_{2, i}^{CSL 2} = \frac{1}{Δ x^{2}} (- 6 {\bar{ϕ}}_{i} + 3 ϕ_{i - 1 / 2} + 3 ϕ_{i + 1 / 2}),

(6)

C_{1, i}^{CSL 2} = \frac{1}{Δ x} (6 {\bar{ϕ}}_{i} - 4 ϕ_{i - 1 / 2} - 2 ϕ_{i + 1 / 2}),

(7)

C_{0, i}^{CSL 2} = ϕ_{i - 1 / 2} .

(8)

By using the interpolation function

Φ_{i}^{CSL 2} (x)

⁠, the boundary value

ϕ_{i - 1 / 2}

is updated by the conservation equation of a differential form

\frac{\partial ϕ}{\partial t} + u \frac{\partial ϕ}{\partial x} = - ϕ \frac{\partial u}{\partial x} .

(9)

A semi-Lagrangian approach is used for Eq.

ϕ_{i - 1 / 2}^{*} = {\begin{matrix} Φ_{i - 1}^{CSL 2} (x_{i - 1 / 2} - u_{i - 1 / 2} Δ t) & if & u_{i - 1 / 2} \geq 0, \\ Φ_{i}^{CSL 2} (x_{i - 1 / 2} - u_{i - 1 / 2} Δ t) & if & u_{i - 1 / 2} < 0, \end{matrix}

(10)

\frac{\partial ϕ}{\partial t} = - ϕ^{*} \frac{\partial u}{\partial x} .

(11)

The cell average

{\bar{ϕ}}_{i}

is update by a finite volume formulation

ϕ_{i}^{n + 1} = ϕ_{i}^{n} - \frac{F_{i + 1 / 2}^{CSL 2} - F_{i - 1 / 2}^{CSL 2}}{Δ x},

(12)

where

F_{i - 1 / 2}^{CSL 2}

is the flux at the cell boundary

x_{i - 1 / 2}

⁠,

F_{i - 1 / 2}^{CSL 2} = {\begin{matrix} - \int_{x_{i - 1 / 2}}^{x_{i - 1 / 2} - u_{i - 1 / 2} Δ t} Φ_{i - 1}^{CSL 2} (x) d x & if & u_{i - 1 / 2} \geq 0, \\ - \int_{x_{i - 1 / 2}}^{x_{i - 1 / 2} - u_{i - 1 / 2} Δ t} Φ_{i}^{CSL 2} (x) d x & if & u_{i - 1 / 2} < 0. \end{matrix}

(13)

B. CIP-CSLR scheme

The CIP-CSLR scheme is briefly explained here. In the CIP-CSLR scheme,³⁵ the following function

Φ_{i}^{CSLR} (x)

⁠,

Φ_{i}^{CSLR} (x) = \frac{C_{2, i}^{CSLR} C_{1, i}^{CSLR} {(x - x_{i - 1 / 2})}^{2} + 2 C_{2, i}^{CSLR} (x - x_{i - 1 / 2}) + ϕ_{i - 1 / 2}}{{(1 + C_{1, i}^{CSLR} (x - x_{i - 1 / 2}))}^{2}},

(14)

is used to interpolate between

x_{i - 1 / 2}

and

x_{i + 1 / 2}

⁠. The coefficients,

C_{2, i}^{CSLR}, C_{1, i}^{CSLR}

⁠, and

C_{0, i}^{CSLR}

⁠, are determined as follows:

C_{2, i}^{CSLR} = C_{1, i}^{CSLR} ϕ_{i} + (ϕ_{i} - ϕ_{i - 1 / 2}) / Δ x,

(15)

C_{1, i}^{CSLR} = \frac{1}{Δ x} (\frac{| (ϕ_{i - 1 / 2} - ϕ_{i}) | + ε}{| (ϕ_{i} - ϕ_{i + 1 / 2}) | + ε} + 1),

(16)

C_{0, i}^{CSLR} = ϕ_{i - 1 / 2},

(17)

by using the same constraints with these of CSL2 []. Here, ε is a small number to avoid zero division. We used

ε = 10^{- 15}

for all results in this paper. The rest of the procedure is the same with CSL2.

C. CIP-CSLT scheme

We propose the CSLT (CSL with tangent hyperbolic function) scheme, which uses a hyperbolic tangent function as the interpolation function (Fig. 2). A similar scheme has been proposed in Ref. 16. In Ref. 16, only cell averages were used to construct the hyperbolic tangent function (like standard finite volume methods) and a wider stencil (three cells) was used. In the proposed scheme, not only cell averages but also boundary values are used (like CIP-CSL schemes) and a shorter stencil (one cell average and two boundary values) is used. The interpolation function is defined as

Φ_{i}^{CSLT} (x) = ϕ_{min} + \frac{ϕ_{max} - ϕ_{min}}{2} (1 + γ tanh (β (\frac{x - x_{i - 1 / 2}}{Δ x} - \tilde{x_{i}}))),

(18)

where

ϕ_{min} = min (ϕ_{i - 1 / 2}, ϕ_{i + 1 / 2}), ϕ_{max} = max (ϕ_{i - 1 / 2}, ϕ_{i + 1 / 2})

⁠, and

γ =

sgn (ϕ_{i + 1 / 2} - ϕ_{i + 1 / 2})

⁠.

\tilde{x_{i}}

is the jump location and can be calculated by using . β is the parameter to control jump thickness. The rest of the procedure is the same with CSL2. The CIP-CSLT scheme is a non-oscillatory scheme because the interpolation function (the hyperbolic tangent function) is a monotonic function with a minimum of

ϕ_{min}

and a maximum of

ϕ_{max}

⁠.

FIG. 2.

Schematic figure of the CIP-CSLT scheme. u i − 1 / 2 < 0 is assumed. The moments, which are indicated by gray color ( ϕ i − 1 / 2 , ϕ ¯ i, and ϕ i + 1 / 2), are used to construct the step function Φ i CSLT ( x ).

View large Download slide

Schematic figure of the CIP-CSLT scheme. $u_{i - 1 / 2} < 0$ is assumed. The moments, which are indicated by gray color (⁠ $ϕ_{i - 1 / 2}, {\bar{ϕ}}_{i}$ ⁠, and $ϕ_{i + 1 / 2}$ ⁠), are used to construct the step function $Φ_{i}^{CSLT} (x)$ ⁠.

D. Non-oscillatory selector (NOS) and CIP-CSL2T scheme

We propose non-oscillatory selector (NOS), which can prevent numerical oscillations of a high-order scheme such as CSL2 by combining with a non-oscillatory scheme such as CSLT and CSLR. The basic idea of the proposed selector is that first we solve the equation using a high-order scheme like CSL2. If oscillations are detected for some cell averages and/or boundary values, we recalculate these oscillatory cell averages and/or boundary values using a non-oscillatory scheme like CSLT. If numerical oscillations are not detected, we simply use the solution of the high-order scheme. Some schematic examples are shown in Fig. 3. For instance, if CSL2 is used as the high-order scheme, CSL2 possibly causes numerical oscillations as shown in Figs. 3(a) and 3(b) because CSL2 employs the quadratic function as the interpolation function. If the interpolation functions, which are indicated by the dotted ellipses in Figs. 3(a) and 3(b), are used to update

ϕ_{i + 1 / 2}, ϕ_{i + 1 / 2}

possibly causes some numerical oscillations. For such cases (when an oscillation was detected), we recalculate

ϕ_{i + 1 / 2}

using CSLT as shown in Figs. 3(a) and 3(b). CSLT has employed the hyperbolic tangent function (a monotonic function) as the interpolation function, so we can prevent the oscillation. If numerical oscillations are not detected like Fig. 3(c), we simply use the solution by CSL2. We also consider a successive gradient in this selector. The procedure is as follows (here we consider that CSL2 is the high-order scheme and CSLT is the non-oscillatory scheme):

FIG. 3.

View large Download slide

Schematic figure of the non-oscillatory selector (NOS) for three different cases. (a) and (b) show typical cases when CSL2 is oscillatory and (c) a typical non-oscillatory case of CSL2.

Step 1: Update all cell averages (⁠ ${\bar{ϕ}}_{i}^{CSL 2, n + 1}$ ⁠) and boundary values (⁠ $ϕ_{i - 1 / 2}^{CSL 2, n + 1}$ ⁠) using CSL2.
Step 2: If the cell average ${\bar{ϕ}}_{i}^{CSL 2, n + 1}$ is oscillatory, the fluxes $F_{i - 1 / 2}^{CSL 2}$ and $F_{i + 1 / 2}^{CSL 2}$ are replaced with $F_{i - 1 / 2}^{CSLT}$ and $F_{i + 1 / 2}^{CSLT}$ ⁠, and ${\bar{ϕ}}_{i}$ is re-updated using CSLT (neighbor cell averages associated with $F_{i - 1 / 2}^{CSLT}$ or $F_{i + 1 / 2}^{CSLT}$ are also re-updated). Then, we also re-update the boundary values, $ϕ_{i - 1 / 2}$ and $ϕ_{i + 1 / 2}$ using CSLT. We assume that ${\bar{ϕ}}_{i}^{CSL 2, n + 1}$ is oscillatory when (⁠ ${\bar{ϕ}}_{i}^{CSL 2, n + 1} > {\bar{ϕ}}_{i}^{max, n} - ε$ ⁠) or (⁠ ${\bar{ϕ}}_{i}^{CSL 2, n + 1} < {\bar{ϕ}}_{i}^{min, n} + ε$ ⁠). Here, ${\bar{ϕ}}_{i}^{max, n}$ is the maximum value of the variables (at time step n), which are used to calculate ${\bar{ϕ}}_{i}^{CSL 2, n + 1}$ ⁠, and ${\bar{ϕ}}_{i}^{min, n}$ is the minimum value of these variables.
Step 3: If the boundary value $ϕ_{i - 1 / 2}^{CSL 2, n + 1}$ is oscillatory, $ϕ_{i - 1 / 2}$ is re-updated using CSLT. Then, $F_{i - 1 / 2}^{CSL 2}$ is replaced with $F_{i - 1 / 2}^{CSLT}$ (if it was not replaced in step 2) and also re-updates cell averages, which are associated with $F_{i - 1 / 2}^{CSLT}$ ⁠. We assume that $ϕ_{i - 1 / 2}^{CSL 2, n + 1}$ is oscillatory when (⁠ $ϕ_{i - 1 / 2}^{CSL 2, n + 1} > ϕ_{i - 1 / 2}^{max, n} - ε$ ⁠) or (⁠ $ϕ_{i - 1 / 2}^{CSL 2, n + 1} < ϕ_{i - 1 / 2}^{min, n} + ε$ ⁠). Here, $ϕ_{i - 1 / 2}^{max, n}$ is the maximum value of the variables (at time step n), which are used to calculate $ϕ_{i - 1 / 2}^{CSL 2, n + 1}$ ⁠, and $ϕ_{i - 1 / 2}^{min, n}$ is the minimum value of these variables.
Exception CSLT cannot handle preexisting cusps (when ${\bar{ϕ}}_{i} > max (ϕ_{i - 1 / 2}, ϕ_{i + 1 / 2}) - ε$ or ${\bar{ϕ}}_{i} < min (ϕ_{i - 1 / 2}, ϕ_{i + 1 / 2}) + ε$ ⁠) because then $\tilde{x}$ is going to be infinity. For preexisting cusps, alternative solvers are required. In this paper, CSL2 is used for cusps when the successive gradient (r_i) is less than α; otherwise, the first-order upwind scheme is used. The successive gradient is defined as

$r_{i} = max (\frac{| {\bar{ϕ}}_{i} - ϕ_{i - 1 / 2} |}{| ϕ_{i + 1 / 2} - \bar{ϕ_{i}} |}, \frac{| ϕ_{i + 1 / 2} - {\bar{ϕ}}_{i} |}{| {\bar{ϕ}}_{i} - ϕ_{i - 1 / 2} |}) .$

(19)

In this paper, $α = 4.01$ is mainly used.

NOS is intuitive and easy to implement. We call the scheme, which combined CSL2 and CSLT as CSL2T (CSL with second-order polynomial and hyperbolic tangent functions).

E. CIP-CSL2R scheme

By using NOS, we can minimize numerical oscillations of CSL2 with any less-oscillatory scheme. In this paper, we also proposed CSL2R (CSL with second-order polynomial and rational functions) scheme, which is based on CSL2 and CSLR. The procedure is almost same with CSL2T. We can simply replace CSLT with CSLR in the CSL2T scheme.

III. NUMERICAL RESULTS

A. Sine wave propagation

The conservation equation is solved with the initial condition

ϕ (x, 0) = sin (2 π x)

⁠. The domain

[0, 1]

⁠, u(x) = 1 and periodic boundary conditions are used. Four different grid sizes (N = 40, 80, 160, and 320) are used with

Δ t = 0.4 Δ x

and

Δ x = 1 / N

⁠. Errors are defined as follows:

L_{1} = \frac{1}{N} \sum_{i = 1}^{N} | ϕ_{i} - ϕ_{exact, i} |,

(20)

L_{\infty} = max (| ϕ_{i} - ϕ_{exact, i} |) .

(21)

Table I shows the numerical results. CSL2 has third-order accuracy. CSLR cannot maintain third-order accuracy. CSL2T and CSL2R have third-order accuracy when

α ≳ 4.001

⁠. This means that NOS does not disturb the accuracy of the CSL2 scheme when

α ≳ 4.001

⁠. If

α ≲ 4.00

(for instance,

α = 3.999

⁠) is used, the order of accuracy cannot be maintained in both CSL2T and CSL2R. In this paper,

α = 4.01

is used for the rest of results.

TABLE I.

L₁ and $L_{\infty}$ errors in sine wave propagation at t = 1.

Method	N	L₁ error	L₁ order	$L_{\infty}$ error	$L_{\infty}$ order
CSL2	40	$9.79 \times 10^{- 5}$	⋯	$1.53 \times 10^{- 4}$	⋯
	80	$1.23 \times 10^{- 5}$	2.99	$1.93 \times 10^{- 5}$	2.99
	160	$1.53 \times 10^{- 6}$	3.01	$2.41 \times 10^{- 6}$	3.00
	320	$1.92 \times 10^{- 7}$	2.99	$3.01 \times 10^{- 7}$	3.00
CSLR	40	$2.29 \times 10^{- 3}$	⋯	$1.19 \times 10^{- 2}$	⋯
	80	$5.53 \times 10^{- 4}$	2.05	$4.75 \times 10^{- 3}$	1.33
	160	$1.29 \times 10^{- 4}$	2.10	$1.83 \times 10^{- 3}$	1.38
	320	$2.75 \times 10^{- 5}$	2.23	$6.91 \times 10^{- 4}$	1.41
CSL2T (any β, $α ≳ 4.001$ ⁠)	40	$9.79 \times 10^{- 5}$	⋯	$1.53 \times 10^{- 4}$	⋯
	80	$1.23 \times 10^{- 5}$	2.99	$1.93 \times 10^{- 5}$	2.99
	160	$1.53 \times 10^{- 6}$	3.01	$2.41 \times 10^{- 6}$	3.00
	320	$1.92 \times 10^{- 7}$	2.99	$3.01 \times 10^{- 7}$	3.00
CSL2T (⁠ $β = 4.0$ ⁠, $α = 3.999$ ⁠)	40	$8.50 \times 10^{- 4}$	⋯	$4.12 \times 10^{- 3}$	⋯
	80	$1.75 \times 10^{- 4}$	2.29	$1.61 \times 10^{- 3}$	1.36
	160	$3.22 \times 10^{- 5}$	2.44	$5.86 \times 10^{- 4}$	1.46
	320	$5.58 \times 10^{- 6}$	2.53	$2.00 \times 10^{- 4}$	1.55
CSL2R (⁠ $α ≳ 4.001$ ⁠)	40	$9.79 \times 10^{- 5}$	⋯	$1.53 \times 10^{- 4}$	⋯
	80	$1.23 \times 10^{- 5}$	2.99	$1.93 \times 10^{- 5}$	2.99
	160	$1.53 \times 10^{- 6}$	3.01	$2.41 \times 10^{- 6}$	3.00
	320	$1.92 \times 10^{- 7}$	2.99	$3.01 \times 10^{- 7}$	3.00
CSL2R (⁠ $α = 3.999$ ⁠)	40	$8.99 \times 10^{- 4}$	⋯	$4.39 \times 10^{- 3}$	⋯
	80	$2.00 \times 10^{- 4}$	2.17	$1.86 \times 10^{- 3}$	1.24
	160	$3.68 \times 10^{- 5}$	2.44	$6.73 \times 10^{- 4}$	1.46
	320	$7.23 \times 10^{- 6}$	2.35	$2.59 \times 10^{- 4}$	1.38

Method	N	L₁ error	L₁ order	$L_{\infty}$ error	$L_{\infty}$ order
CSL2	40	$9.79 \times 10^{- 5}$	⋯	$1.53 \times 10^{- 4}$	⋯
	80	$1.23 \times 10^{- 5}$	2.99	$1.93 \times 10^{- 5}$	2.99
	160	$1.53 \times 10^{- 6}$	3.01	$2.41 \times 10^{- 6}$	3.00
	320	$1.92 \times 10^{- 7}$	2.99	$3.01 \times 10^{- 7}$	3.00
CSLR	40	$2.29 \times 10^{- 3}$	⋯	$1.19 \times 10^{- 2}$	⋯
	80	$5.53 \times 10^{- 4}$	2.05	$4.75 \times 10^{- 3}$	1.33
	160	$1.29 \times 10^{- 4}$	2.10	$1.83 \times 10^{- 3}$	1.38
	320	$2.75 \times 10^{- 5}$	2.23	$6.91 \times 10^{- 4}$	1.41
CSL2T (any β, $α ≳ 4.001$ ⁠)	40	$9.79 \times 10^{- 5}$	⋯	$1.53 \times 10^{- 4}$	⋯
	80	$1.23 \times 10^{- 5}$	2.99	$1.93 \times 10^{- 5}$	2.99
	160	$1.53 \times 10^{- 6}$	3.01	$2.41 \times 10^{- 6}$	3.00
	320	$1.92 \times 10^{- 7}$	2.99	$3.01 \times 10^{- 7}$	3.00
CSL2T (⁠ $β = 4.0$ ⁠, $α = 3.999$ ⁠)	40	$8.50 \times 10^{- 4}$	⋯	$4.12 \times 10^{- 3}$	⋯
	80	$1.75 \times 10^{- 4}$	2.29	$1.61 \times 10^{- 3}$	1.36
	160	$3.22 \times 10^{- 5}$	2.44	$5.86 \times 10^{- 4}$	1.46
	320	$5.58 \times 10^{- 6}$	2.53	$2.00 \times 10^{- 4}$	1.55
CSL2R (⁠ $α ≳ 4.001$ ⁠)	40	$9.79 \times 10^{- 5}$	⋯	$1.53 \times 10^{- 4}$	⋯
	80	$1.23 \times 10^{- 5}$	2.99	$1.93 \times 10^{- 5}$	2.99
	160	$1.53 \times 10^{- 6}$	3.01	$2.41 \times 10^{- 6}$	3.00
	320	$1.92 \times 10^{- 7}$	2.99	$3.01 \times 10^{- 7}$	3.00
CSL2R (⁠ $α = 3.999$ ⁠)	40	$8.99 \times 10^{- 4}$	⋯	$4.39 \times 10^{- 3}$	⋯
	80	$2.00 \times 10^{- 4}$	2.17	$1.86 \times 10^{- 3}$	1.24
	160	$3.68 \times 10^{- 5}$	2.44	$6.73 \times 10^{- 4}$	1.46
	320	$7.23 \times 10^{- 6}$	2.35	$2.59 \times 10^{- 4}$	1.38

B. Jiang-Shu test

We validate proposed CSL schemes through Jiang-Shu problem¹³ u(x) = 1, the domain

[- 1, 1]

⁠, N = 200,

Δ t = 0.4 Δ x, Δ x = 2 / N

⁠, and periodic boundary conditions are used in this test. The initial condition is given as

ϕ (x, 0) = {\begin{array}{l} \frac{1}{6} (G (x, β^{J S}, z - δ) + G (x, β^{J S}, z + δ) + 4 G (x, β^{J S}, z)) & if - 0.8 \leq x < - 0.6, \\ 1 & if - 0.4 \leq x < - 0.2, \\ 1 - | 10 (x - 0.1) | & if 0.0 \leq x < 0.2, \\ \frac{1}{6} (F (x, α^{J S}, a - δ) + F (x, α^{J S}, a + δ) + 4 F (x, α^{J S}, a)) & if 0.4 \leq x < 0.6, \\ 0, & otherwise, \end{array}

(22)

where

G (x, β^{J S}, z) = e^{- β^{J S} {(x - z)}^{2}},

(23)

F (x, α^{J S}, a) = \sqrt{max (1 - {(α^{J S})}^{2} {(x - a)}^{2}, 0)},

(24)

where a = 0.5,

z = - 0.7, δ = 0.005, α^{J S} = 10

⁠, and

β^{J S} = log (2) / (36 δ^{2})

⁠.

Figure 4 shows the results at t = 2 (1 period). CSL2 is oscillatory around discontinuities as shown in Fig. 4(a). Although CSLR can prevent numerical oscillations, the result is diffusive as shown in Fig. 4(b). CSL2T (β = 4) and CSL2R can minimize numerical oscillations of CSL2 without increasing diffusions much [Figs. 4(c) and 4(d)].

FIG. 4.

View large Download slide

Numerical results of Jiang-Shu test at t = 2 (1 period). (a)–(d) are numerical results by CSL2, CSLR, CSL2T (β = 4), and CSL2R, respectively. N = 200 and CFL = 0.4 are used.

Figure 5 shows the results when the parameter β in CSLT was altered. If a smaller β [e.g., Fig. 5(a)] is used, the solution becomes slightly diffusive and errors are increased as shown in Table II. If a larger β [e.g., Figs. 5(b) and 5(c)] is used, the solution becomes less diffusive and errors are decreased as shown in Table II. If a large β (e.g., β = 10) is used, CSL2T can delete numerical oscillations almost without increasing diffusion as shown in Fig. 5(c). Although it is better to use larger β for the 1D test problem, there are some side effects of sharp sigmoid function in multidimensional fluid simulations (an example is shown in lid-driven cavity flows later). Therefore, in this paper, β = 4 is mainly used. CSL2T is superior to CSL2R when $β ≳ 3$ is used in this test.

FIG. 5.

View large Download slide

Numerical results of Jiang-Shu test at t = 2 (1 period). (a)–(c) are numerical results by CSL2T with β = 3, β = 5, and β = 10, respectively. N = 200 and CFL = 0.4 are used.

TABLE II.

L₁ and $L_{\infty}$ errors in the complex wave propagation test at t = 1 (1 period). N = 200 and CFL = 0.4 are used.

	L₁ error	$L_{\infty}$ error
CSL2	$1.49 \times 10^{- 2}$	$2.46 \times 10^{- 1}$
CSLR	$1.90 \times 10^{- 2}$	$2.60 \times 10^{- 1}$
CSL2T (β = 3)	$1.27 \times 10^{- 2}$	$2.33 \times 10^{- 1}$
CSL2T (β = 4)	$1.22 \times 10^{- 2}$	$2.33 \times 10^{- 1}$
CSL2T (β = 5)	$1.19 \times 10^{- 2}$	$2.32 \times 10^{- 1}$
CSL2T (β = 10)	$1.12 \times 10^{- 2}$	$2.31 \times 10^{- 1}$
CSL2R	$1.28 \times 10^{- 2}$	$2.33 \times 10^{- 1}$

C. Burgers' equation

In this test, we solve the inviscid Burgers' equation in its conservative formulation

\frac{\partial u}{\partial t} + \frac{\partial (u^{2} / 2)}{\partial x} = 0,

(25)

with the initial condition

u (x, 0) = 0.5 + 0.4 cos (2 π x)

⁠. N = 100, CFL = 0.2, and periodic boundary conditions are used. Figure 6 shows the results at t = 1. The reference solution is created by using CSL3CW with N = 10 000. CSL2 has numerical oscillation around the discontinuity. CSLR, CSL2T, and CSL2R can capture the discontinuity well without numerical oscillation. CSL2T and CSL2R are slightly less diffusive than CSLR around the discontinuity.

FIG. 6.

View large Download slide

Numerical results of Burgers' equation at t = 1. N = 100 and CFL = 0.2 are used.

D. Zalesak problem

Zalesak's test problem⁴⁷ in which a notched circle is rotated is widely used as a test of scalar advection schemes. The initial condition is given by

ϕ = {\begin{matrix} 1 & if & \sqrt{{(x - 0.5)}^{2} + {(y - 0.75)}^{2}} < 0.17 and \\ (y > 0.85 or | x - 0.5 | > 0.03), \\ 0 & if & others, \end{matrix}

(26)

u (x, y) = (y - 0.5, 0.5 - x) .

(27)

The one revolution is completed with 628 time steps.

Figure 7 shows numerical results (top view of 0.5-contour) by CSL2, CSLR, CSL2T (β = 4), and CSL2R after one revolution. Both results show good agreement with the exact solution. Figure 8 shows the perspective views of the numerical results by CSL2, CSLR, CSL2T (β = 4), and CSL2R. As shown in Fig. 8, the numerical result by CSL2 is oscillatory (a), CSLR can minimize numerical oscillations but the result is diffusive panel (b). CSL2T and CSL2R can minimize oscillations without increasing numerical diffusions much as shown in Figs. 8(c) and 8(d), respectively. Table III shows L₁ errors after one revolution and four revolutions. The result indicates that CSL2T is superior to CSL2R in this test and the trend is unchanged after four revolutions.

FIG. 7.

View large Download slide

Top views of numerical results of Zalesak problem (after one revolution). (a)–(d) are numerical results by CSL2, CSLR, CSL2T (β = 4), and CSL2R, respectively. The dotted and solid lines represent 0.5 contour lines of numerical results and the exact solution, respectively. A Cartesian grid of 100 × 100 was used.

FIG. 8.

View large Download slide

Numerical results of Zalesak problem (after one revolution). (a)–(d) are numerical results by CSL2, CSLR, CSL2T (β = 4), and CSL2R, respectively. A Cartesian grid of 100 × 100 was used.

TABLE III.

L₁ errors in Zalesak problem (after one revolution and four revolutions). A Cartesian grid of 100 × 100 was used.

Schemes	L₁ Error after one revolution	L₁ error after four revolutions
CSL2	$1.55 \times 10^{- 2}$	$2.21 \times 10^{- 2}$
CSLR	$1.73 \times 10^{- 2}$	$2.75 \times 10^{- 2}$
CSL2T (β = 4)	$1.28 \times 10^{- 2}$	$1.97 \times 10^{- 2}$
CSL2R	$1.30 \times 10^{- 2}$	$1.99 \times 10^{- 2}$

Schemes	L₁ Error after one revolution	L₁ error after four revolutions
CSL2	$1.55 \times 10^{- 2}$	$2.21 \times 10^{- 2}$
CSLR	$1.73 \times 10^{- 2}$	$2.75 \times 10^{- 2}$
CSL2T (β = 4)	$1.28 \times 10^{- 2}$	$1.97 \times 10^{- 2}$
CSL2R	$1.30 \times 10^{- 2}$	$1.99 \times 10^{- 2}$

E. Lid-driven cavity flows

The proposed CIP-CSL schemes were applied to lid-driven cavity flow problems, and these numerical results are compared with numerical results by Ghia et al.⁸ Figure 9 shows numerical results of lid-driven cavity flows of Re = 7500 by CSLR, CSL2T (β = 4), and CSL2R (CSL2 is not considered hereafter because such oscillatory schemes cause some issues such as unnatural oscillations in solutions, and instability of calculations). A Cartesian grid of 128 × 128 was used. All numerical results have shown good agreements with the reference solution by Ghia (a grid of $256 \times 256$ was used in the reference solution). However, as shown in Fig. 9, we cannot observe differences of these results clearly with the resolution 128 × 128. Therefore, we compare these results with a lower numerical resolution. Figure 10 shows the numerical results when a Cartesian grid of 64 × 64 was used. In this resolution, it is clear that CSL2T and CSL2R are superior to CSLR. To observe the differences among CSL2T, CSL2R, and CSLR more clearly, we plotted numerical results of CSLR, CSL2T, and CSL2R together, and magnified the area $- 0.5 \leq x \leq 0$ and $- 1 \leq y \leq - 0.5$ ⁠. As shown in Fig. 11, CSL2T is slightly superior to CSL2R. We also study the effect of β of CSL2T in this test. As shown in Fig. 12(a), if β = 3 is used, CSL2T is even less accurate than CSLR. If β is increased from β = 3 up to around $β = 4.5$ ⁠, the results are improved and are superior to results by CSL2R as well as CSLR. However, as shown in Fig. 12(b), if β is increased further around from $β = 4.5$ ⁠, the results are getting worse. Based on these results in Fig. 12, optimal β would be around $4 \leq β \leq 5$ ⁠. If around $4 \leq β \leq 5$ is used, CSL2T is superior to CSL2R. In this paper, β = 4 is used hereafter.

FIG. 9.

View large Download slide

Numerical results of lid-driven cavity flows of Re = 7500. (a)–(c) are numerical results by CSLR, CSL2T (β = 4), and CSL2R, respectively. A Cartesian grid of 128 × 128 was used. The solid and dot lines represent x- and y-components of the velocity fields on the lines x = 0 and y = 0, respectively. Dots represent numerical results by Ghia et al.⁸

FIG. 10.

View large Download slide

Numerical results of lid-driven cavity flows of Re = 7500. (a)–(c) are numerical results by CSLR, CSL2T (β = 4), and CSL2R, respectively. A Cartesian grid of 64 × 64 was used.

FIG. 11.

View large Download slide

Numerical results of lid-driven cavity flows of Re = 7500. These lines are numerical results by CSLR, CSL2T (β = 4), and CSL2R, respectively. A Cartesian grid of 64 × 64 was used.

FIG. 12.

View large Download slide

Numerical results of lid-driven cavity flows of Re = 7500. (a) Represents numerical results by CSL2T with β = 3, $β = 3.5$ ⁠, β = 4, and $β = 4.5$ ⁠. (b) Represents numerical results by CSL2T with β = 4, β = 5, β = 6, and β = 10. A Cartesian grid of 64 × 64 was used.

Figures 13(a) and 13(b) show numerical results of Re = 1000 and 5000 by CSL2T (β = 4), respectively. CSL2T can also simulate cavity flows of other Reynolds numbers well. We have also confirmed that CSLR and CSL2R schemes can also simulate cavity flows of these Reynolds number well. Figure 14 shows the velocity field of Re = 7500 by CSL2T when a Cartesian grid of 64 × 64 was used. Secondary vorticity at bottom right corner has been captured with the resolution of 64 × 64.

FIG. 13.

View large Download slide

Numerical results of lid-driven cavity flows of Re = 1000 (a) and 5000 (b) when a Cartesian grid of 128 × 128 was used.

FIG. 14.

View large Download slide

Numerical results of lid-driven cavity flows of Re = 7500 by CSL2T (β = 4). A Cartesian grid of 64 × 64 was used.

F. Two droplets collisions and separations

We performed three-dimensional numerical simulations of two droplets collision and separation, which involves topology change of liquid interfaces.³ The numerical method to simulate free surface flows is based on VSIAM3-TM,⁴⁴ the CLSVOF (coupled level set^19,26 and volume-of-fluid^11,14) method,^25,38 the THINC/WLIC (tangent of hyperbola for interface capturing/weighted line interface calculation) scheme,^12,32,39 and the density scaled balanced continuum surface force model^43,45 with level set curvature correction.⁴² For the full details of the implementation of the free surface flow solver, see Refs. 40 and 42–44.

Figures 15(a)–15(c) show snapshots of the numerical results of two droplets collision and separation (We = 40) by CSLR, CSL2R, and CSL2R, respectively. The corresponding experimental result can be found in Ref. 3. In these numerical simulations, quantitative parameters were used. The density ratio is 1.25:1000 (air:liquid). The mesh size $Δ = D / 7.4$ was used, where D is the diameter of initial droplets. In this numerical resolution (⁠ $Δ = D / 7.4$ ⁠), only CSL2T could reproduce the separation, and CSL2R and CSLR failed to do it. However, if the resolution was increased, CSLR and CSL2R could also reproduce the separation with $Δ = D / 7.5$ and $Δ = D / 8.0$ ⁠, respectively. The numerical results have demonstrated that CSL2T and CSL2R can reduce numerical diffusion in the free surface flow application compared to CSLR, and CSL2T is slightly better than CSL2R in this test. Although there are some previous numerical work of droplet collision and separations,^2,4,20,28 compared to these previous work, our simulation has captured the phenomena with a quite course grid.

FIG. 15.

Numerical results of two droplet collisions and separation by CSLR (a), CSL2T (b), and CSL2R (c) with the corresponding experimental result (d).3 The time evolution is from right to left. The mesh size is Δ = D / 7.4. Reproduced with permission from N. Ashgriz and J. Y. Poo, “Coalescence and separation in binary collisions of liquid drops,” J. Fluid Mech. 221, 183 (1990). Copyright 1990 Cambridge University Press.

View large Download slide

Numerical results of two droplet collisions and separation by CSLR (a), CSL2T (b), and CSL2R (c) with the corresponding experimental result (d).³ The time evolution is from right to left. The mesh size is $Δ = D / 7.4$ ⁠. Reproduced with permission from N. Ashgriz and J. Y. Poo, “Coalescence and separation in binary collisions of liquid drops,” J. Fluid Mech. 221, 183 (1990). Copyright 1990 Cambridge University Press.

G. Droplet splashing

We performed numerical simulations of prompt splashing using CSLR, CSL2T, and CSL2R with VSIAM3-TM.⁴⁴ The corresponding experimental result can be found in Ref. 29. In the numerical simulations, quantitative parameters, the densities $ρ_{liquid} = 1000$ ⁠, $ρ_{air} = 1.25$ kg/m³, viscosities $μ_{liquid} = 1.0 \times 10^{- 3}$ ⁠, $μ_{air} = 1.82 \times 10^{- 5}$ Pa s, surface tension $σ = 7.2 \times 10^{- 2}$ N/m, gravity 9.8 m/s², initial droplet diameter D = 1.86 mm, impact speed 2.98 m/s, and the equilibrium contact angle 163° are used. A Cartesian grid of $224 \times 224 \times 48$ is used. The mesh size $Δ = D / 45$ was used. Figure 16 shows the result. CSLR, CSL2T, and CSL2R can capture the physics of droplet splashing including satellite droplets and spikes. It seems that CSL2T and CSL2R capture satellite droplets and spikes more than CSLR and that the results by CSLR are more diffusive than these by CSL2T and CSL2R. The result has also shown the CSL2T and CSL2R as well as CSLR are robust.

FIG. 16.

View large Download slide

A comparison among the numerical results by CSLR (a), CSL2T (b), and CSL2R (c). A distilled water droplet of 1.86 (mm) impacts onto a superhydrophobic substrate (the equilibrium angle is 163°). The droplet impact speed is 2.980 (m/s). A Cartesian grid of $224 \times 224 \times 48$ and $α = 1.5 Δ x$ is used. The mesh size is $Δ = D / 45$ ⁠. The corresponding experimental result can be found in Ref. 29.

IV. SUMMARY

We proposed novel third-order less oscillatory and less diffusive compact stencil-based upwind schemes, which are called the CIP-CSL2T and CIP-CSL2R schemes for the approximation of hyperbolic conservation laws. We also proposed NOS, which can select a suitable CSL scheme between CSL2 and CSLT in CSL2T, and CSL2 and CSLR in CSL2R. Both CSL2T and CSL2R have third-order accuracy for the sine wave test and capture discontinuities and smooth solutions without numerical oscillations. The proposed schemes could capture the secondary vorticity at bottom right corner of lid-driven cavity flow of Re = 7500 with a Cartesian grid of 64 × 64. The numerical results of two droplets collision/separation of We = 40 show that the proposed schemes can reproduce droplets collision/separation with low numerical resolutions of 7.4–7.5 meshes for the diameter of initial droplets. The numerical results of droplet splashing have demonstrated that proposed schemes can reduce numerical diffusions well against CSLR and robust.

ACKNOWLEDGMENTS

This research was supported in part by FLEXIS and IROHMS, which are part-funded by the European Regional Development Fund (ERDF), through the Welsh Government. The numerical simulations were partially conducted on computers at Yukawa Institute of Theoretical Physics in Kyoto University and at HPC Wales. We acknowledge the support of the Supercomputing Wales project, which is part-funded by the European Regional Development Fund (ERDF) via Welsh Government.

AUTHOR DECLARATIONS

Conflict of Interest

The author declares that there are no conflicts to disclose.

Author Contributions

Kensuke Yokoi: Conceptualization (lead); Data curation (lead); Formal analysis (lead); Funding acquisition (lead); Investigation (lead); Methodology (lead); Project administration (lead); Resources (lead); Software (lead); Validation (lead); Visualization (lead); Writing – original draft (lead); Writing – review & editing (lead).

DATA AVAILABILITY

The data that support the findings of this study are available within the article.

REFERENCES

1.

M.

Al-Mosallam

and

K.

Yokoi

, “

Efficient implementation of volume/surface integrated average-based multi-moment method

,”

Int. J. Comput. Methods

14

(

02

),

1750010

(

2017

).

https://doi.org/10.1142/S0219876217500104

Google Scholar

Crossref

2.

A.

Amani

,

N.

Balcázar

,

E.

Gutiérrez

, and

A.

Oliva

, “

Numerical study of binary droplets collision in the main collision regimes

,”

Chem. Eng. J.

370

,

477

–

498

(

2019

).

https://doi.org/10.1016/j.cej.2019.03.188

Google Scholar

Crossref

3.

N.

Ashgriz

and

J. Y.

Poo

, “

Coalescence and separation in binary collisions of liquid drops

,”

J. Fluid Mech.

221

,

183

–

204

(

1990

).

https://doi.org/10.1017/S0022112090003536

Google Scholar

Crossref

4.

X.

Chen

and

V.

Yang

, “

Direct numerical simulation of multiscale flow physics of binary droplet collision

,”

Phys. Fluids

32

(

6

),

062103

(

2020

).

https://doi.org/10.1063/5.0006695

Google Scholar

Crossref

5.

S.

Clain

,

S.

Diot

, and

R.

Loubère

, “

A high-order finite volume method for systems of conservation laws–multi-dimensional optimal order detection mood

,”

J. Comput. Phys.

230

(

10

),

4028

–

4050

(

2011

).

https://doi.org/10.1016/j.jcp.2011.02.026

Google Scholar

Crossref

6.

P.

Colella

and

P. R.

Woodward

, “

The piecewise parabolic method (PPM) for gas-dynamical simulations

,”

J. Comput. Phys.

54

(

1

),

174

–

201

(

1984

).

https://doi.org/10.1016/0021-9991(84)90143-8

Google Scholar

Crossref

7.

X.

Deng

,

Y.

Shimizu

,

B.

Xie

, and

F.

Xiao

, “

Constructing higher order discontinuity-capturing schemes with upwind-biased interpolations and boundary variation diminishing algorithm

,”

Comput. Fluids

200

,

104433

(

2020

).

https://doi.org/10.1016/j.compfluid.2020.104433

Google Scholar

Crossref

8.

U.

Ghia

,

K.

Ghia

, and

C.

Shin

, “

High-Re solutions for incompressible flow using the Navier–Stokes equations and a multigrid method

,”

J. Comput. Phys.

48

(

3

),

387

–

411

(

1982

).

https://doi.org/10.1016/0021-9991(82)90058-4

Google Scholar

Crossref

9.

A.

Harten

,

B.

Engquist

,

S.

Osher

, and

S. R.

Chakravarthy

, “

Uniformly high order accurate essentially non-oscillatory schemes, III

,”

J. Comput. Phys.

71

(

2

),

231

–

303

(

1987

).

https://doi.org/10.1016/0021-9991(87)90031-3

Google Scholar

Crossref

10.

A.

Harten

and

S.

Osher

, “

Uniformly high-order accurate nonoscillatory schemes. I

,”

SIAM J. Numer. Anal.

24

(

2

),

279

–

309

(

1987

).

https://doi.org/10.1137/0724022

Google Scholar

Crossref

11.

C.

Hirt

and

B.

Nichols

, “

Volume of fluid (VOF) method for the dynamics of free boundaries

,”

J. Comput. Phys.

39

(

1

),

201

–

225

(

1981

).

https://doi.org/10.1016/0021-9991(81)90145-5

Google Scholar

Crossref

12.

S.

Ii

,

K.

Sugiyama

,

S.

Takeuchi

,

S.

Takagi

,

Y.

Matsumoto

, and

F.

Xiao

, “

An interface capturing method with a continuous function: The THINC method with multi-dimensional reconstruction

,”

J. Comput. Phys.

231

(

5

),

2328

–

2358

(

2012

).

https://doi.org/10.1016/j.jcp.2011.11.038

Google Scholar

Crossref

13.

G.-S.

Jiang

and

C.-W.

Shu

, “

Efficient implementation of weighted ENO schemes

,”

J. Comput. Phys.

126

(

1

),

202

–

228

(

1996

).

https://doi.org/10.1006/jcph.1996.0130

Google Scholar

Crossref

14.

J.

Li

, “

Calcul d'interface affine par morceaux

,”

C. R. Acad. Sci., Ser. II: Méc., Phys., Chim., Astron.

320

(

8

),

391

–

396

(

1995

).

Google Scholar

15.

Q.

Li

,

S.

Omar

,

X.

Deng

, and

K.

Yokoi

, “

Constrained interpolation profile conservative semi-Lagrangian scheme based on third-order polynomial functions and essentially non-oscillatory (CIP-CSL3ENO) scheme

,”

Commun. Comput. Phys.

22

(

3

),

765

–

788

(

2017

).

https://doi.org/10.4208/cicp.OA-2016-0117

Google Scholar

Crossref

16.

Q.

Li

,

J.

Xia

,

K.

Yokoi

, and

S.

Omar

, “

Boundary variation diminished conservative semi-Lagrangian method for both compressible and incompressible flows

,”

Phys. Fluids

33

(

11

),

117114

(

2021

).

https://doi.org/10.1063/5.0067226

Google Scholar

Crossref

17.

X.-D.

Liu

,

S.

Osher

, and

T.

Chan

, “

Weighted essentially non-oscillatory schemes

,”

J. Comput. Phys.

115

(

1

),

200

–

212

(

1994

).

https://doi.org/10.1006/jcph.1994.1187

Google Scholar

Crossref

18.

N.

Onodera

,

T.

Aoki

, and

K.

Yokoi

, “

A fully conservative high-order upwind multi-moment method using moments in both upwind and downwind cells

,”

Int. J. Numer. Methods Fluids

82

(

8

),

493

–

511

(

2016

).

https://doi.org/10.1002/fld.4228

Google Scholar

Crossref

19.

S.

Osher

and

J. A.

Sethian

, “

Fronts propagating with curvature-dependent speed: Algorithms based on Hamilton-Jacobi formulations

,”

J. Comput. Phys.

79

(

1

),

12

–

49

(

1988

).

https://doi.org/10.1016/0021-9991(88)90002-2

Google Scholar

Crossref

20.

Y.

Pan

and

K.

Suga

, “

Numerical simulation of binary liquid droplet collision

,”

Phys. Fluids

17

(

8

),

082105

(

2005

).

https://doi.org/10.1063/1.2009527

Google Scholar

Crossref

21.

P. L.

Roe

, “

Characteristic-based schemes for the Euler equations

,”

Annu. Rev. Fluid Mech.

18

(

1

),

337

–

365

(

1986

).

https://doi.org/10.1146/annurev.fl.18.010186.002005

Google Scholar

Crossref

22.

C.-W.

Shu

and

S.

Osher

, “

Efficient implementation of essentially non-oscillatory shock-capturing schemes

,”

J. Comput. Phys.

77

(

2

),

439

–

471

(

1988

).

https://doi.org/10.1016/0021-9991(88)90177-5

Google Scholar

Crossref

23.

C.-W.

Shu

and

S.

Osher

, “

Efficient implementation of essentially non-oscillatory shock-capturing schemes, II

,”

J. Comput. Phys.

83

(

1

),

32

–

78

(

1989

).

https://doi.org/10.1016/0021-9991(89)90222-2

Google Scholar

Crossref

24.

Z.

Sun

,

S.

Inaba

, and

F.

Xiao

, “

Boundary variation diminishing (BVD) reconstruction: A new approach to improve Godunov schemes

,”

J. Comput. Phys.

322

,

309

–

325

(

2016

).

https://doi.org/10.1016/j.jcp.2016.06.051

Google Scholar

Crossref

25.

M.

Sussman

and

E. G.

Puckett

, “

A coupled level set and volume-of-fluid method for computing 3D and axisymmetric incompressible two-phase flows

,”

J. Comput. Phys.

162

(

2

),

301

–

337

(

2000

).

https://doi.org/10.1006/jcph.2000.6537

Google Scholar

Crossref

26.

M.

Sussman

,

P.

Smereka

, and

S.

Osher

, “

A level set approach for computing solutions to incompressible two-phase flow

,”

J. Comput. Phys.

114

(

1

),

146

–

159

(

1994

).

https://doi.org/10.1006/jcph.1994.1155

Google Scholar

Crossref

27.

R.

Tanaka

,

T.

Nakamura

, and

T.

Yabe

, “

Constructing exactly conservative scheme in a non-conservative form

,”

Comput. Phys. Commun.

126

(

3

),

232

–

243

(

2000

).

https://doi.org/10.1016/S0010-4655(99)00473-7

Google Scholar

Crossref

28.

S.

Tanguy

and

A.

Berlemont

, “

Application of a level set method for simulation of droplet collisions

,”

Int. J. Multiphase Flow

31

(

9

),

1015

–

1035

(

2005

).

https://doi.org/10.1016/j.ijmultiphaseflow.2005.05.010

Google Scholar

Crossref

29.

P.

Tsai

,

S.

Pacheco

,

C.

Pirat

,

L.

Lefferts

, and

D.

Lohse

, “

Drop impact upon micro- and nanostructured superhydrophobic surfaces

,”

Langmuir

25

(

20

),

12293

–

12298

(

2009

).

https://doi.org/10.1021/la900330q

Google Scholar

Crossref

PubMed

30.

F.

Xiao

, “

Unified formulation for compressible and incompressible flows by using multi-integrated moments. I. One-dimensional inviscid compressible flow

,”

J. Comput. Phys.

195

(

2

),

629

–

654

(

2004

).

https://doi.org/10.1016/j.jcp.2003.10.014

Google Scholar

Crossref

31.

F.

Xiao

,

R.

Akoh

, and

S.

Ii

, “

Unified formulation for compressible and incompressible flows by using multi-integrated moments. II. Multi-dimensional version for compressible and incompressible flows

,”

J. Comput. Phys.

213

(

1

),

31

–

56

(

2006

).

https://doi.org/10.1016/j.jcp.2005.08.002

Google Scholar

Crossref

32.

F.

Xiao

,

Y.

Honma

, and

T.

Kono

, “

A simple algebraic interface capturing scheme using hyperbolic tangent function

,”

Int. J. Numer. Methods Fluids

48

(

9

),

1023

–

1040

(

2005

).

https://doi.org/10.1002/fld.975

Google Scholar

Crossref

33.

F.

Xiao

,

A.

Ikebata

, and

T.

Hasegawa

, “

Numerical simulations of free-interface fluids by a multi-integrated moment method

,”

Comput. Struct.

83

(

6–7

),

409

–

423

(

2005

).

https://doi.org/10.1016/j.compstruc.2004.06.005

Google Scholar

Crossref

34.

F.

Xiao

and

T.

Yabe

, “

Completely conservative and oscillationless semi-Lagrangian schemes for advection transportation

,”

J. Comput. Phys.

170

(

2

),

498

–

522

(

2001

).

https://doi.org/10.1006/jcph.2001.6746

Google Scholar

Crossref

35.

F.

Xiao

,

T.

Yabe

,

X.

Peng

, and

H.

Kobayashi

, “

Conservative and oscillation-less atmospheric transport schemes based on rational functions

,”

J. Geophys. Res.: Atmos.

107

(

D22

),

ACL 2–1

, https://doi.org/10.1029/2001JD001532 (

2002

).

Google Scholar

36.

T.

Yabe

,

R.

Tanaka

,

T.

Nakamura

, and

F.

Xiao

, “

An exactly conservative semi-Lagrangian scheme (CIP-CSL) in one dimension

,”

Mon. Weather Rev.

129

(

2

),

332

–

344

(

2001

).

https://doi.org/10.1175/1520-0493(2001)129<0332:AECSLS>2.0.CO;2

Google Scholar

Crossref

37.

S.

Yamashita

,

C.

Chen

,

K.

Takahashi

, and

F.

Xiao

, “

Large scale numerical simulations for multi-phase fluid dynamics with moving interfaces

,”

Int. J. Comput. Fluid Dyn.

22

(

6

),

405

–

410

(

2008

).

https://doi.org/10.1080/10618560802199691

Google Scholar

Crossref

38.

K.

Yokoi

, “

Numerical method for complex moving boundary problems in a Cartesian fixed grid

,”

Phys. Rev. E

65

,

055701

(

May 2002

).

https://doi.org/10.1103/PhysRevE.65.055701

Google Scholar

Crossref

39.

K.

Yokoi

, “

Efficient implementation of THINC scheme: A simple and practical smoothed VOF algorithm

,”

J. Comput. Phys.

226

(

2

),

1985

–

2002

(

2007

).

https://doi.org/10.1016/j.jcp.2007.06.020

Google Scholar

Crossref

40.

K.

Yokoi

, “

A numerical method for free-surface flows and its application to droplet impact on a thin liquid layer

,”

J. Sci. Comput.

35

(

2

),

372

–

396

(

2008

).

https://doi.org/10.1007/s10915-008-9202-z

Google Scholar

Crossref

41.

K.

Yokoi

, “

Numerical studies of droplet splashing on a dry surface: Triggering a splash with the dynamic contact angle

,”

Soft Matter

7

,

5120

–

5123

(

2011

).

https://doi.org/10.1039/c1sm05336a

Google Scholar

Crossref

42.

K.

Yokoi

, “

A practical numerical framework for free surface flows based on CLSVOF method, multi-moment methods and density-scaled CSF model: Numerical simulations of droplet splashing

,”

J. Comput. Phys.

232

(

1

),

252

–

271

(

2013

).

https://doi.org/10.1016/j.jcp.2012.08.034

Google Scholar

Crossref

43.

K.

Yokoi

, “

A density-scaled continuum surface force model within a balanced force formulation

,”

J. Comput. Phys.

278

,

221

–

228

(

2014

).

https://doi.org/10.1016/j.jcp.2014.08.034

Google Scholar

Crossref

44.

K.

Yokoi

,

M.

Furuichi

, and

M.

Sakai

, “

An efficient multi-dimensional implementation of VSIAM3 and its applications to free surface flows

,”

Phys. Fluids

29

(

12

),

121611

(

2017

).

https://doi.org/10.1063/1.4996183

Google Scholar

Crossref

45.

K.

Yokoi

,

R.

Onishi

,

X.-L.

Deng

, and

M.

Sussman

, “

Density-scaled balanced continuum surface force model with a level set based curvature interpolation technique

,”

Int. J. Comput. Methods

13

(

04

),

1641004

(

2016

).

https://doi.org/10.1142/S0219876216410048

Google Scholar

Crossref

46.

K.

Yokoi

,

D.

Vadillo

,

J.

Hinch

, and

I.

Hutchings

, “

Numerical studies of the influence of the dynamic contact angle on a droplet impacting on a dry surface

,”

Phys. Fluids

21

(

7

),

072102

(

2009

).

https://doi.org/10.1063/1.3158468

Google Scholar

Crossref

47.

S. T.

Zalesak

, “

Fully multidimensional flux-corrected transport algorithms for fluids

,”

J. Comput. Phys.

31

(

3

),

335

–

362

(

1979

).

https://doi.org/10.1016/0021-9991(79)90051-2

Google Scholar

Crossref

Published open access through an agreement with Cardiff University School of Engineering

Third-order less oscillatory and less diffusive compact stencil-based upwind schemes, and their applications to incompressible flows and free surface flows

I. INTRODUCTION

II. NUMERICAL METHODS

A. CIP-CSL2 scheme

B. CIP-CSLR scheme

C. CIP-CSLT scheme

D. Non-oscillatory selector (NOS) and CIP-CSL2T scheme

E. CIP-CSL2R scheme

III. NUMERICAL RESULTS

A. Sine wave propagation

B. Jiang-Shu test

C. Burgers' equation

D. Zalesak problem

E. Lid-driven cavity flows

F. Two droplets collisions and separations

G. Droplet splashing

IV. SUMMARY

ACKNOWLEDGMENTS

AUTHOR DECLARATIONS

Conflict of Interest

Author Contributions

DATA AVAILABILITY

REFERENCES

Citing articles via

Submit your article

Sign up for alerts

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

Third-order less oscillatory and less diffusive compact stencil-based upwind schemes, and their applications to incompressible flows and free surface flows Open Access

I. INTRODUCTION

II. NUMERICAL METHODS

A. CIP-CSL2 scheme

B. CIP-CSLR scheme

C. CIP-CSLT scheme

D. Non-oscillatory selector (NOS) and CIP-CSL2T scheme

E. CIP-CSL2R scheme

III. NUMERICAL RESULTS

A. Sine wave propagation

B. Jiang-Shu test

C. Burgers' equation

D. Zalesak problem

E. Lid-driven cavity flows

F. Two droplets collisions and separations

G. Droplet splashing

IV. SUMMARY

ACKNOWLEDGMENTS

AUTHOR DECLARATIONS

Conflict of Interest

Author Contributions

DATA AVAILABILITY

REFERENCES

Citing articles via

Submit your article

Sign up for alerts

Related Content

Resources

Explore

pubs.aip.org

Connect with AIP Publishing

This Feature Is Available To Subscribers Only

Third-order less oscillatory and less diffusive compact stencil-based upwind schemes, and their applications to incompressible flows and free surface flows