## 1 Introduction

In the paper two-fluid model of gas-solid mixture is discussed. In this model gas is considered as a carrier phase and dust grains are considered as dispersed phase. It is assumed that the solid phase has one typical size in each volume (so this phase is monodisperse). Therefore, the continuity and motion equations for gas and dispersed phase have the following form:

(1) |

(2) |

where and are volume density of gas and dust, and are velocities of gas and dust, is pressure, is gravity acceleration, are sourses and sinks for gas and dust, are forces affecting gas and dust except for pressure, gravitaion and drag, is drag force per unit volume:

(3) |

where is a velocity relaxation time. Here is particle size, is the sound speed in gas. In this work we consider particular case when does not depend on that corresponds to Epstein and Stokes regimes (see details, e.g. in [1]).

Computing of pure gas dymanics (solution of (1) with ) using explicit schemes requires time step that satisfies Courant condition:

(4) |

where is the Courant limeter. Additional necessary condition arises during computing of dusty gas dynamisc using explicit schemes

(5) |

Violation of this condition leads to numerical instability. Condition (5) is extremally prohibitive for intense interphase interaction (for small ). Intense interphase interaction arises in many applications of gas-particle mixtures (e.g. in modelling reactors with finely-dispersed catalyst, in planet formation from gas-dust circumstellar disks etc.) and is characterized by the fact that time of momentum transfer is much less than process time.

In the case of intense interphase interaction, is called stiff relaxation term [2]. Effective schemes for problems with stiff relaxation terms are designed in way that they could preserve asymptotical solution even with . Main ideas of such design are described in [3, 4] and assosiated with using of implicit approximation of stiff relaxation term along with explisit approximation of other terms. This approaches are developed for euler methods for solving fluid dynamics equations and are used in different applications e.g. [5, 6]. Transfering this ideas to lagrangian methods (for example, to TFSPH in which gas and dust are modelled by different sets of particles) encounter new difficulties caused by the fact that carrier and disperse phases quantities are known in different points of space. For TFSPH authors of [7, 8, 9]

shown that the way of interpolation this parametres for computing

influences the method’s property to preserve asymptotic. In particular, [7] demonstrated that classical method for computing drag force in TFSPH [10] captures the asymptotic of the solution for small only with(6) |

Getting over this restriction is crucial for modelling mixtures with intense interphase interaction. For this reason authors of [8], [9] proposed other methods to computing drag forces. In this paper we present quantitative comparison of the classical [10] and new [8], [9] approaches focusing on their ability to preserve asymptotic properties of the solution. The methods are described in detail in section 2, the test problem results are given in section 3 and the summary is provided in section 4.

##
2 Methods for Computing Drag Terms in Two-Fluid

Smoothed Particle Hydrodynamics

Further we will give the schemes for solving the equations of gas and dust motion (7) in the standard SPH notation. We will consider only the schemes in which gas and dust are simulated by different sets of particles, i.e. by the two-fluid approach for smoothed particle hydrodynamics (TFSPH). Let be the number of the time step. Following the notation introduced in [10], we will use as the indices for gas particles, and as the indices for dust particles.

### 2.1 The MK Monaghan–Kocharyan explicit scheme.

A method for computing the drag force, which was proposed in [10] (hereinafter referred to as MK (Monaghan-Kocharyan Drag)), is classical for smoothed particle hydrodynamics. This method is based on the computing of the relative velocity between each pair of gas-dust particles and is employed in astrophysical and engineering applications of two-phase medium mechanics [11, 12, 13, 14] and others.

We implemented this scheme so that the summand accounting for drag uses the velocities from the previous time step:

(8) |

(9) |

(10) |

where and are the masses of gas and dust particles, respectively, , is a clipping constant, , is the radius of a spherical dust particle with index , is the constant determined by dimensionality of the problem (for one-dimensional problems, ), and is the smoothing kernel.

### 2.2 The semi-implicit ISPH scheme with interpolation of the first order approximation with respect to time.

The second method for computing the drag consists in the calculation of gas characteristics at the points where dust particles are located (and vice versa) using the SPH interpolation formulas:

(11) |

where is the dust velocity at a spatial point where the gas particle with index is located, and is the gas velocity at a spatial point where the dust particle with index is located.

As a result, all features of the gas-dust medium become known for each model particle. This method and its modifications are applied in refs. [8, 15, 16]. We apply this idea to construct a semi-implicit scheme that would not require the fulfillment of condition (5) for obtaining stable solutions. In particular, parsimonious is the following ISPH scheme with the first order approximation with respect to time (the quantities calculated using interpolation formulas are marked in blue, while the quantities derived from those calculated by interpolation formulas are marked in red):

(12) |

(13) |

(14) |

(15) |

(16) |

### 2.3 A new SPH-IDIC scheme – the implicit ‘drag in cell”.

In addition, computing of the drag force can be based on the idea of the particle-in-cell method for simulation of gas-dust flows [17]. The parsimonious semi-implicit SPH-IDIC approach based on this idea was suggested and tested in our earlier paper [9]. A detailed description of this approach is presented below.

At each time instant, we will decompose the entire calculation region into disjoint volumes so that the merging of these volumes will coincide with the entire region. Suppose a separate volume contains gas particles of a similar mass and dust particles of a similar mass , with , . Introduce the volume-averaged values of and (anywise) and assume that

(17) |

thus determining

(18) |

Let us assume that in computing the drag force that acts from gas on dust, the gas velocity is constant over the entire volume and equal to , whereas dust particles have different velocities (and vice versa). In addition, we will calculate the drag factor and density using values of the quantities from the preceding time layer, and relative velocity – from the subsequent layer. The resulting scheme will have the form

(19) |

(20) |

(21) |

If the time derivative in (19)-(20) is approximated to the first order, there exists a parsimonious method to calculate . In [9], it was shown that the semi-implicit scheme (19)-(21) with the first order approximation with respect to time satisfies the momentum conservation law for each cell, i.e. the momentum lost by gas due to drag on dust completely coincides with the momentum acquired by dust due to drag on gas.

## 3 Results

To measure the ability of the schemes to preserve asymptotical properties of the solution in case of intense interphase interaction two well-known test problems with available reference solution are used. The first is a problem of sound-wave propagation in isothermal gas-dust mixture. This problem has a smooth solution and suits perfectly to study the way of drag computing in SPH. We will refer to this problem as Dustywave. The second is a problem of shock wave propagation in coupled gas-dust mixture with initial conditions known as Sod shock tube. Due to discontinuity of the solution and complexity of wave structure this problem is common and challenging test for computational gas dynamics. For this problem Dustyshock name is reserved. Both Dustywave and Dustyshock problems are described in detail [9]. Moreover, for all numerical experiments in the paper we take the same physical and numerical parameters as in [9]. In particular, for both problems we take high drag coefficient and high concentration of dust in gas

(22) |

which is guarantee intense interphase interaction and challenge for simulation.

Fig. 1

shows the solution of Dustywave problem at the time moment

with MK, ISPH, IDIC methods. In this case (22) leads to . The figure displays the dust and gas velocity obtained with different smoothing lengths . The left panels show the results of computing for explicit schemes MK. At , the time step and the number of particles are used; while at , the step and . The middle and right panels present the results of computing for semi-implicit schemes ISPH and IDIC with ; in this case, the number of particles was the same as for the left panels.Fig. 2 shows gas and dust velocity as the solution of Dustyshock problem at the time moment . The same methods MK, ISPH, IDIC as in Fig. 1 are used. For we take , , for — ,

One can see that the numerical solutions obtained by MK and ISPH methods with the smoothing length increased from (the condition (6) is satisfied) to ((6) is violated) acquire a pronounced dissipation. The observed tendency to solution dissipation was described in [8, 7]. We can see from Figs. 1, 2 and Tables 1, 2 that the maximum level of dissipation is obtained in the case of explicit schemes MK. As follows from the central panel in Fig. 1, the semi-implicit ISPH scheme gives a smaller dissipation at the same smoothing length as compared to MK. Moreover, thanks to semi-implicit approximation of the drag force, the ISPH scheme has no restrictions on the time step (5).

MK and ISPH are the fully Lagrangian methods, which means that all forces are calculated without the introduction of a spatial grid. The IDIC method is a combination of Lagrangian and Euler approaches because the drag force is computed using the decomposition of particles into Euler volumes. One can see on the right upper panel of Fig.1 that these numerical solutions obtained by IDIC are free of dissipation, and at the wave amplitude is reproduced without visible error, in distinction to MK and ISPH.

MK | ISPH | IDIC | ||
---|---|---|---|---|

0.0417 | 0.0332 | 0.0003 | ||

0.1486 | 0.079 | 0.0012 |

MK | ISPH | IDIC | ||
---|---|---|---|---|

0.1177 | 0.2101 | 0.0457 | ||

0.1786 | 0.29 | 0.0605 |

## 4 Summary

Simulation of the dynamics of gas-aerosol particle mixtures is computationally challenging, especially in particle methods as Smoothed particle hydrodynamics. In the paper we compared ability of fully lagrangian methods MK [10] and ISPH [8] and euler-lagrangian method IDIC [9] to reproduce asymptotical properties of the solution for dynamics of gas-dust mixtures with intense interphase interaction. We found that IDIC method where drag is computed using euler cells is asymptotic preserving and allows to use timestep and smoothing length independent on drag intensity. MK and ISPH methods require fine spatial resolution for intense interphase interaction.

Acknowledgements. This work was supported by the Russian Science Foundation grant 19-71-10026.

## References

- [1] T. Saito, M. Marumoto, K. Takayama. Numerical investigations of shock waves in gas-particle mixtures. Evaluation of numerical methods for dusty-gas shock wave phenomena. Shock Waves 13: 299–322 (2003)
- [2] S. Jin, C.D. Livermore. Numerical Schemes for Hyperbolic Conservation Laws with Stiff Relaxation Terms. Journal of Computational Physics 126, 449–467 (1996)
- [3] S. Jin. Asymptotic preserving schemes for multiscale kinetic and hyperbolic equations: a review. Riv. Mat. Univ. Parma, 3, 177–216 (2012)
- [4] G. Albi, G. Dimarco, L. Pareschi. Implicit-Explicit multistep methods for hyperbolic systems with multiscale relaxation. arXiv e-prints: 1904.03865
- [5] O. P. Stoyanovskaya, V. N. Snytnikov, E. I. Vorobyov. Analysis of methods for computing the trajectories of dust particles in a gas-dust circumstellar disk. Astron.Rep., Vol. 61, No. 12:1044–1060 (2017)
- [6] P. Degond, F. Deluzet. Asymptotic-Preserving methods and multiscale models for plasma physics. Journal of Computational Physics. 336, 429–457 (2017)
- [7] G. Laibe, D. J. Price. DUSTYBOX and DUSTYWAVE: two test problems for numerical simulations of two-fluid astrophysical dust-gas mixtures. MNRAS, 418:1491–1497 (2011)
- [8] P. Lorén-Aguilar, M. R. Bate. Two-fluid dust and gas mixtures in smoothed particle hydrodynamics: a semi-implicit approach. MNRAS, 443:927–945 (2014)
- [9] O. P. Stoyanovskaya, T. A. Glushko, N. V. Snytnikov, V. N. Snytnikov, Two-Fluid Dusty Gas in Smoothed Particle Hydrodynamics: Fast and Implicit Algorithm for Stiff Linear Drag. Astronomy and Computing, 25:25–37 (2018)
- [10] J. J. Monaghan, A. Kocharyan. SPH simulation of multi-phase flow. Computer Physics Communications, 87:225–235 (1995)
- [11] S. T. Madisson, R. J. Humble, J. R. Murray. Building Planets with Dusty Gas. Bioastronomy 2002: Life Among the Stars, 213:231 (2004)
- [12] L. Barriere-Fouchet, J. F. Gonzalez, J. R. Murray, R. J. Humble, S. T. Maddison. Dust distribution in protoplanetary disks. Vertical settling and radial migration. AAp, 443:185–194 (2005)
- [13] J. F. Gonzalez, G. Laibe, S. T. Maddison. Self-induced dust traps: overcoming planet formation barriers. MNRAS, 467:1984–1996 (2017)
- [14] Xiong, Qingang and Deng, Lijuan and Wang, Wei and Ge, Wei. SPH method for two-fluid modeling of particle–fluid fluidization. Chemical Engineering Science - CHEM ENG SCI, 66:1859–1865 (2011)
- [15] R. A. Booth, D. Sijacki, C. J. Clarke. Smoothed particle hydrodynamics simulations of gas and dust mixtures. MNRAS, 452:3932–3947 (2015)
- [16] W. K. M. Rice, G. Lodato, J. E. Pringle, P. J. Armitage, I. A. Bonnell. Accelerated planetesimal growth in self-gravitating protoplanetary discs. MNRAS, 355:543–552 (2004)
- [17] M. J. Andrews, P. J. O’Rourke. The multiphase Particle-in-Cell (MP-PIC) method for dense particulate flows. Int. J. of Multiphase Flow, 22:379–402 (1996)

Comments

There are no comments yet.