Particle flows around an intruder

TL;DR

This study numerically investigates particle flows around an intruder, revealing a transition in drag force laws from Epstein's to Newton's and Stokes' laws, and identifies turbulent-like behavior at high speeds and large intruder sizes.

Contribution

It demonstrates the crossover of drag laws depending on flow speed and intruder size, and analyzes turbulent-like flow behavior behind the intruder.

Findings

Crossover from Epstein's to Newton's law based on speed ratio.

Transition from Epstein's to Stokes' law in low-speed regime.

Detection of turbulent-like flow behind the intruder at high speeds.

Abstract

Particle flows injected as beams and scattered by an intruder are numerically studied. We find a crossover of the drag force from Epstein's law to Newton's law, depending on the ratio of the speed to the thermal speed. These laws can be reproduced by a simple analysis of a collision model between the intruder and particle flows. The crossover from Epstein's law to Stokes' law is also found for the low-speed regime as the time evolution of the drag force caused by beam particles. We also show the existence of turbulent-like behavior of the particle flows behind the intruder with the aid of the second invariant of the velocity gradient tensor and the relative mean square displacement for the high-speed regime and a large intruder.

Tables

Table 1. Table 1: The set of used parameters.

$D_c/d$	$N_s$	$N$	$D_t/d$
$5$	$144$	$30,000$	$20$
$13$	$744$	$101,250$	$30$
$28$	$3,364$	$810,000$	$60$
$98$	$39,204$	$9,600,000$	$200$

Table 2. Table 2: The set of used parameters in Appendix B .

$D_c/d$	$N_s$	$N$	$D_t/d$
$5$	$144$	$30,000$	$20$
$9$	$400$	$30,000$	$40$
$13$	$784$	$30,000$	$40$
$16$	$1,156$	$30,000$	$40$
$19$	$1,600$	$30,000$	$40$

Equations

$$F={\int }_0^{\pi }\sin \theta d\theta \Delta p_{\parallel }{\Omega }_{\mathrm{c}}=\frac{\pi }{3}nm(D+d{)}^2{V}^2,$$

$$F=\frac{\pi }{3}nm(D+d{)}^2{v}_{\mathrm{T}}V,$$

$$F_{\mathrm{St}}=3\pi \eta DV,$$

$$\eta (\varphi )$$

$$\times \left[\frac{1}{g_0(\varphi )}{\left(1+\frac{8}{5}\varphi g_0(\varphi )\right)}^2+\frac{768}{25\pi }{\varphi }^2{g}_0(\varphi )\right].$$

$$F=2\pi \eta {\left(\frac{D}{2}\right)}^2{\int }_0^{\pi }d\theta {\dot{\gamma }}_{\mathrm{loc}}\sin \theta \sim \eta DV,$$

$${\omega }_{xz}=\frac{\partial v_z}{\partial x}-\frac{\partial v_x}{\partial z}.$$

$${\Delta }^2(t)\equiv ⟨|\delta {\mathbf{r}}_i(t)-\delta {\mathbf{r}}_j(t){|}^2⟩/d^2,$$

$${\theta }_0={\cos }^{-1}\left[1-{\left(\frac{D}{D_{\mathrm{t}}}\right)}^2\right],$$

$$\frac{p}{nT}=Z(\varphi )=\frac{1+\varphi +{\varphi }^2-{\varphi }^3}{(1-\varphi {)}^3},$$

$$C_{V}dT+T{\left(\frac{\partial p}{\partial T}\right)}_V\frac{d{L}^3}{N}=0,$$

$$\log \frac{p}{\varphi Z(\varphi )}-\frac{2}{3}{\int }^{\varphi }\frac{Z({\varphi }^{\prime })}{{\varphi }^{\prime }}d{\varphi }^{\prime }=const.$$

$$c(\varphi )=\sqrt{{\left(\frac{\partial p}{\partial \rho }\right)}_S}=f(\varphi )\sqrt{\frac{T}{m}},$$

$$f(\varphi )=\sqrt{\frac{5+10\varphi -3{\varphi }^2-24{\varphi }^3+37{\varphi }^4-22{\varphi }^5+5{\varphi }^6}{3(1-\varphi {)}^6}}.$$

$$F=\frac{\pi }{3}nm(D+d{)}^{2}cV,$$

Full text

I Introduction

To know fluid flows around an intruder is a fundamental problem Lamb ; Batchelor . When the Reynolds number is low, the drag force acting on a spherical intruder obeys Stokes’ law in which the drag is proportional to the relative speed between the intruder and the fluid, the viscosity of the fluid, and the radius of the intruder. On the other hand, the drag force satisfies Newton’s law for the high Reynolds number, in which the drag is proportional to the square of the relative speed and the cross section of the intruder.

To understand particle flows is important in various fields Hirano06 ; Adil06 ; Das10 ; Armitage ; Vallado14 ; Ramsey ; Patterson09 ; Hutzler12 ; Cho75 ; Lee87 . The perfect fluidity is a key concept to understand quark-gluon matter Hirano06 ; Adil06 . The drag coefficient of quarks through the quark-gluon plasma is evaluated theoretically Das10 . The drag force acting on dust in the protoplanetary disks is known to show the crossover from Epstein’s law Epstein24 , in which the drag force is proportional to the cross section and the relative moving speed, to Stokes’ law as the average size of dust increases Armitage . This problem is also related to the designs of artificial satellites Vallado14 . Atomic and molecular beams Ramsey ; Patterson09 ; Hutzler12 are widely used for nanotechnologies such as a beam epitaxy on a thin film Cho75 . We believe that the drag force still satisfies Stokes’ law for molecular flows of low Reynolds number Vergeles95 ; Chen06 ; Li08 ; Itami15 ; Asano18 . In atomic and molecular flows, Stokes’ law can be used only for systems in the zero Knudsen number limit, i.e., if the mean free path of molecules is much smaller than the size of the intruder Cunningham10 ; Knudsen11 ; Sone ; Takata93 ; Taguchi15 ; Taguchi17 . The correction to Stokes’ drag for rarefied gases in the low Knudsen number is theoretically confirmed by the kinetic theory of rarefied gases Cunningham10 ; Knudsen11 ; Sone ; Takata93 ; Taguchi15 ; Taguchi17 . The drag force, however, acting on a slowly moving intruder in a crowd of molecules for the high Knudsen number satisfies Epstein’s law Epstein24 ; Li03 .

There are various experimental and numerical studies on the drag force acting on an intruder in granular flows Albert99 ; Chehata03 ; Wassgren03 ; Geng04 ; Bharadwaj06 ; Potiguar13 ; Takehara10 ; Reddy11 ; Hilton13 ; Guillard13 ; Takada17 ; Candelier10 ; Gnann ; Wang14 . Note that thermal fluctuations do not play any roles for granular particles. A variety of velocity dependences of the drag force is reported depending on the protocols of experiments and simulations. Granular jet experiments Cheng07 ; Ellowitz13 and simulations Ellowitz13 ; Sano12 ; Muller14 suggest that the granular jet flows can be approximated by a perfect fluid model Ellowitz13 ; Muller14 .

It is natural to expect that a molecular flow can create turbulence if the corresponding Reynolds number of the flow is sufficiently high. To verify this conjecture, the Kármán vortices behind an intruder have already been observed in molecular dynamics (MD) simulations Komatsu14 ; Asano18 ; Rapaport86 . Nevertheless, there are no published papers, as long as we know, to reproduce a fully developed turbulent flow by MD.

In this paper, we numerically study the drag force acting on a stationary spherical intruder in particle flows by controlling the ratio of the injected speed of the particles to the thermal speed, which is proportional to the sound speed, in terms of the MD. We also try to reproduce a turbulent-like behavior of particle flows behind the intruder with the aid of the second invariant of the velocity gradient tensor and the relative mean square displacement of two moving particles if the injected speed is much larger than the thermal speed and the size of the intruder is much larger than the molecule size.

The organization of this paper is as follows: In the next section, we explain our model and setup of our study. In Sec. III, we present the results of our simulation for the drag force acting on the intruder based on not too large systems. In the first part, we show the crossover from Epstein’s drag to Newtonian drag depending on the ratio of the colliding speed to the thermal speed. In the second part, we discuss how results depend on the boundary condition on the surface of the intruder. In Sec. IV, we illustrate the existence of a crossover from Epstein’s drag to Stokes’ drag as time goes on. This section consists of two parts. In the first subsection, we show the numerical results to exhibit the crossover. In the second subsection, we explain the mechanism of the crossover. In Sec. V, we demonstrate the existence of a turbulent-like flow behind the intruder in which the relative motion of two colliding particles is super-ballistic if the size of the intruder is much larger than the molecule size. In the first subsection, we illustrate the conditions to observe turbulent-like flows produced by molecule flows. We also discuss the angle distribution of scattered particles by the intruder. In Sec. VI, we conclude our results. In Appendix A, we examine whether the sound speed is applicable to characterize the drag. In Appendix B, we briefly summarize our numerical results when the intruder is initially put inside the beam particles under the periodic boundary condition in the flow direction.

II Model

In this section, let us explain our model and the setup of our simulation. The system consists of two parts: one is a fixed intruder, and another is a collection of the mobile particles, as shown in Fig. 1. The intruder is made by one core particle whose diameter is $D_{\rm c}$. In most cases, the intruder is covered by $N_{\rm s}$ identical particles whose diameters are $d$ on the surface of the core particle. We examine four sizes of the intruders as $(D_{\rm c}/d,N_{\rm s})=(5,144)$, $(13,744)$, $(28,3{,}364)$, and $(98,39{,}204)$. The diameter of the intruder is given by $D=D_{\rm c}+2d$ in this case. We also examine the case of $N_{\rm s}=0$, i.e., $D=D_{\rm c}$ in Sec. III.2. Then, we find that the results are almost independent of the boundary condition on the surface of the intruder. The intruder is fixed at the origin, which has an infinitely large mass. We simulate systems containing $N=30{,}000$, $101{,}250$, $810{,}000$, and $9{,}600{,}000$ monodisperse mobile particles depending on the size of the intruder. The mass and diameter for each molecule are $m$ and $d$, respectively. Throughout this paper, collisions between particles are assumed to be elastic. We examine various initial volume fractions of beam particles ranged from $\phi=0.20$ to $0.55$. Before starting our simulation, the mobile particles are thermalized with the temperature $T$, and are moved with the translational speed $V$. The equation of motion of $i$-th particle at the position $\bm{r}_{i}$ is given by $md^{2}\bm{r}_{i}/dt^{2}=\sum_{j}\bm{F}_{ij}$ with the interparticle force $\bm{F}_{ij}=\Theta(d-r_{ij})k(d-r_{ij})\hat{\bm{r}}_{ij}$, where we have introduced $r_{ij}=|\bm{r}_{ij}|$, $\hat{\bm{r}}_{ij}=\bm{r}_{ij}/r_{ij}$, the spring constant $k$, and the step function $\Theta(x)$, i.e., $\Theta(x)=1$ for $x\geq 0$ and $\Theta(x)=0$ otherwise. Initially, we confine the mobile particles in tubes whose diameters are $D_{\rm t}=20d$, $30d$, $60d$, and $200d$ for $D_{\rm c}=5d$, $13d$, $28d$, and $98d$, respectively. Since the drag force is known to depend on the distance from the boundary of a container through a simulation Vergeles95 and a theory by fluid mechanics Brenner61 , we keep the ratio of $D_{\rm t}/D\simeq 2$ in our simulation. We assume that the interaction force between the wall of the tube and the particles is identical to the interparticle force. The used parameters are listed in Table 1.

When the mobile particles collide with the intruder, we examine two cases for the reflection: one is the random reflection of the angle with the temperature $T_{\rm w}$, and another is the simple reflection rule according to the equation of motion. For the former condition, when the colliding particle leaves the intruder, we give the velocity to the particle, whose magnitude is the thermal speed $v_{\rm T}=\sqrt{2T/m}$, and the scattering direction of the particle is random on the surface of the intruder. Here, we set $T_{\rm w}=T$ for simplicity for the former condition. We also examine two cases for the systems behind the intruder: one is a free scattering case, and another is a confined case, where the scattered particles are still confined in the tube. In the following, we use the time-averaged drag force in a certain time window, where the instantaneous force acting on the intruder is measured via $\bm{F}=\sum_{j}\bm{F}_{j}$, where $\bm{F}_{j}$ is the contact force acting on the intruder. In addition, we fix the speed $V=0.10d\sqrt{k/m}$. We have verified that the results of our simulation for $V=0.10d\sqrt{k/m}$ are consistent with those of the hard-core particles. We also note that the time increment is fixed to be $\Delta t=1.0\times 10^{-3}\sqrt{m/k}$.

III Crossover from Epstein’s law to Newtonian law

In this section, we present the results of our simulations for not large system sizes. This section consists of two parts. In Sec. III.1, we show the drag force against the translational speed as a crossover from Epstein’s law to Newtonian law, depending on the ratio of the colliding speed to the thermal speed. In Sec. III.2, we investigate whether the geometry of the intruder affects the drag law.

III.1 Crossover of the drag force from Epstein’s law to Newtonian law

We present the results of our simulation of the drag force acting on the intruder against the normalized colliding speed by the thermal speed for both the free scattering (with fixing $D/d$ and various $\phi$) and the confined cases (with fixing $\phi$ and various $D/d$) in Fig. 2. Here, the data are obtained by the time average of the instantaneous drag force in the range between $t=400\sqrt{m/k}$ and $800\sqrt{m/k}$. Because the drag force is determined from the total force acting on the intruder by all collisions of particles in front of the intruder, the drag force is unchanged whether the particle flow behind the intruder is confined or not. In this subsection, we assume that the intruder is covered by small particles whose sizes are the same as those of colliding particles and the reflection rule between colliding particles and the intruder is thermal in which reflecting particles are stochastically scattered with satisfying the Maxwell distribution.

We note that our obtained drag forces are proportional to $(D+d)^{2}$, i.e., the collision cross section Wassgren03 ; Takada19 . The drag force is proportional to the square of the speed for $V/v_{\rm T}\gg 1$, while it is proportional to $v_{\rm T}V$ for $V/v_{\rm T}\ll 1$.

The former, on the one hand, can be understood by a simple collision model. Because each momentum change in the flow direction for a collision is $\Delta p_{\parallel}=2mV\cos^{2}\theta$ where $\theta$ is the angle between the intruder and the colliding particle as shown in Fig. 3, and the collision frequency is $\Omega_{\rm c}=(\pi/4)n(D+d)^{2}V$ with the number density $n=6\phi/(\pi d^{3})$ Wassgren03 , the force acting on the intruder is given by

[TABLE]

which agrees well with the simulation results for $V\gtrsim v_{\rm T}$ (see Fig. 2). This is the simple Newtonian drag law for the high-speed regime.

On the other hand, the front of the beam of mobile particles diffuses before colliding with the intruder for $V/v_{\rm T}\ll 1$ as shown in Fig. 4(a). This expansion speed of the diffusive front is approximately given by $V_{\rm p}=\sqrt{Vv_{\rm T}}$ as shown in Fig. 4(b). When we replace $V$ in Eq. (1) by $V_{\rm p}$, we obtain

[TABLE]

which agrees well with the simulation results for $V\lesssim v_{\rm T}$ as shown in Fig. 2. Here, the sound speed is another characteristic speed that appears in this system. However, it is evident that the thermal speed is superior to the sound speed to characterize the drag force as shown in Appendix A.

III.2 Effects of the boundary conditions on the intruder

In this subsection, let us check whether the results depend on the boundary conditions between the intruder and the mobile particles. Here, we refer to the intruder introduced in the previous subsection as (i) “thermal and bumpy”, because the mobile particles reflect at random when they collide with the intruder, and the small particles are attached on the surface of the core particle. We examine the other three different types of intruders: (ii) “reflective and smooth,” (iii) “reflective and bumpy,” and (iv) “thermal and smooth.” The intruder in the condition of (ii) “reflective and smooth” consists of only one core particle, and the reflection between the intruder and the mobile particles is described by a simple elastic collision rule. The intruder in (iii) “reflective and bumpy” adopts a simple elastic collision rule, and the intruder consists of the core particle and the small particles on its surface as that used in the case (i). The intruder in (iv) “thermal and smooth” adopts the thermal scattering as used in (i), and the intruder consists of only one core intruder. Figure 5 plots the results of the drag forces under various boundary conditions for $D/d=7$ and $\phi=0.40$ in the free scattering condition. The results clearly indicate that the choice of the boundary condition on the intruder is not important for the drag force.

IV Crossover from Epstein’s law to Stokes’ law

In this section, we present the crossover from Epstein’s law to Stokes’ law for the low-speed regime. This section consists of two parts. In the first subsection, we present numerical results to exhibit the crossover from Epstein’s law to Stokes’ law. In the second subsection, we explain the physical mechanism of this crossover.

IV.1 How can we get crossover from Epstein’s law to Stokes’ law?

In the previous section, we have reported the crossover from Epstein’s law to Newtonian law when we control the ratio of the colliding speed to the thermal speed $V/v_{\rm T}$. Although Newtonian law for $V/v_{\rm T}\gg 1$ is expected,111We have confirmed that the drag law in this regime remains Newtonian in the range of $200\sqrt{m/k}\lesssim t\lesssim 1500\sqrt{m/k}$. Epstein’s law for $V/v_{\rm T}\ll 1$ is a little unexpected because Stokes’ law might be expected for slow flows. To clarify the condition to emerge Stokes’ law, in this section, we illustrate the existence of a crossover from Epstein’s law to Stokes’ law with the time evolution by fixing $V/v_{\rm T}\ll 1$.

When we focus on the region in front of the intruder, the packing fraction of this local region changes as time goes on (see Fig. 6(a)). Here, we measure the local packing fraction in the region $-D/2-10d\leq z\leq-D/2$ (see Fig. 3 for the information of the geometry). After the early stage of the collision process, there is a relatively long metastable state in which the local packing fraction is almost equal to the initial packing fraction of the beam particles. Then, $\phi_{\rm loc}$ decreases with time for $t/\sqrt{m/k}\gtrsim 600$. When we substitute the local packing fraction $\phi_{\rm loc}(t)$ into Eq. (2), we can observe the crossover from Epstein’s drag to Stokes’ drag as shown in Fig. 6(b). Then, the drag force reaches Stokes’ law in the late stage, which is given by

[TABLE]

where the shear viscosity $\eta$ is estimated from the well-known result from the Enskog theory Ferziger ; Garzo99 as

[TABLE]

Here, $g_{0}(\phi)$ is the radial distribution function at contact, which is approximately given by Carnahan and Starling formula $g_{0}(\phi)=(1-\phi/2)/(1-\phi)^{3}$, which is valid for $\phi<0.49$ Carnahan69 . We also confirm that the drag force is proportional to the cross section in Epstein’s regime, while it is proportional to $D$ in Stokes’ regime as shown in Figs. 6(c) and (d). This is also another evidence of the crossover from Epstein’s law to Stokes’ law. We note that Stokes’ law is also observed when we adopt the periodic boundary condition in the flow direction, which is investigated in Appendix B.

IV.2 Mechanism of the crossover

Let us discuss the reason why the time-dependent crossover from Epstein’s to Stokes’ laws appears. The flow is stacked around the intruder because the particles cannot move due to the existence of the outer particles near the intruder. These stacked particles can be regarded as a “boundary layer” around the intruder, which prohibits mobile particles from direct collisions with the intruder. In this case, the local shear rate can be approximated as $\dot{\gamma}_{\rm loc}\simeq(V/D)\sin\theta$, where $\sin\theta$ represents the projection parallel to the tangential direction of the surface. This leads to that the shear stress acts on the intruder whose magnitude in the flow direction is $\dot{\gamma}_{\rm loc}\eta\sin\theta$ if the viscosity $\eta$ is well defined. Therefore, we can evaluate the force acting on the intruder by integrating the shear stress over the surface of the intruder, and we obtain

[TABLE]

which is the origin of Stokes’ law.

The particles, however, are not stacked near the intruder in the early stage of the simulation. Therefore, the beam particles can directly collide with the intruder, whose speed is approximately given by the thermal speed. Then, we obtain Epstein’s law as in Eq. (2).

V Turbulent-like flows for $V/v_{\rm T}\gg 1$ and $D/d\gg 1$

In Sec. III, we did not discuss what flows can be observed for $V/v_{\rm T}\gg 1$. It is natural to expect turbulent-like flows can be observed in this regime, because $V/v_{\rm T}$ might correspond to the Reynolds number of the fluid flows. We also examine the effect of $D/d$ which has not been investigated in the previous sections.

Let us characterize the particle flows behind the intruder for large $V/v_{\rm T}$. We introduce the second invariant of the velocity gradient tensor $Q=(1/2)(-S_{\alpha\beta}S_{\alpha\beta}+W_{\alpha\beta}W_{\alpha\beta})$ where $S_{\alpha\beta}=(1/2)(\partial_{\beta}u_{\alpha}+\partial_{\alpha}u_{\beta})$ and $W_{\alpha\beta}=(1/2)(\partial_{\beta}u_{\alpha}-\partial_{\alpha}u_{\beta})$ Hunt88 . Here, we adopt Einstein’s rule for $\alpha$ and $\beta$ where duplicated indices take summation over $x$, $y$, and $z$. Figure 7 shows the contours of $Q=0$ for (a) $D=7d$ and (b) $D=100d$, where the field is coarse-grained with the scale $2d$ for visibility Zhang10 . The vortex rich regions $Q>0$ emerge behind the intruder. Such domain structure becomes complicated for large $D/d$. This behavior is similar to that observed in a turbulent flow induced by an intruder Ott00 . In particular, huge number of vortex rich domains appear for $D/d\gg 1$ (Fig. 7). This suggests that the particles flow in Fig. 7(b) is turbulent-like.

We also characterize the vortex structure by plotting the vorticity induced by the scattering of the intruder. Because the twisted structure of the flow field is not observed in our simulation, we focus on the flow structure in $xz$-plane (see Fig. 3 for the information of the geometry). Let us introduce the vorticity in $xz$-plane with $|y_{i}|\leq D/2$ as

[TABLE]

Figure 8 shows a typical snapshot of the vorticity field in $xz$-plane for $D=100d$. The positive and negative vorticity regimes are generated in the vicinity of the intruder and move toward the downstream. See also the movie in the Supplemental Material Movie . These complicated structures of vorticities are also similar to those observed in turbulent flows.

Let us study how the particles are scattered after collisions by the intruder. We focus on the relative motion of the particles which collide with the intruder almost simultaneously through the mean square displacement

[TABLE]

with $\delta\bm{r}_{i}(t)\equiv\bm{r}_{i}(t+t_{\rm c})-\bm{r}_{i}(t_{\rm c})$, where we only select two particles ($i$ and $j$) within the interval $|t_{i}-t_{j}|<\Delta t_{\rm th}\equiv 10\sqrt{m/k}$ with the collision times $t_{i}$ and $t_{j}$ with the intruder for $i$-th and $j$-th particles, respectively. Note that we choose $t_{\rm c}$ in $\delta\bm{r}_{i}(t)$ or $\delta\bm{r}_{j}(t)$ in Eq. (7) to be much larger time of $t_{i}$ and $t_{j}$. This means that these two particles exist sufficiently far from the intruder at initial. We also note that we select particles whose position in the $z$-direction is located in $-135\leq z_{i}/d\leq-100$ at time $t=0$ for $D/d=100$. This is because we try to clarify how localized particles behave as time goes on. Figure 9 shows the super-ballistic behavior $\Delta^{2}(t)\sim t^{\beta}$ where the best fitted values of exponent $\beta$ are $\beta=2.2$ for $D=7d$ and $\beta=3.0$ for $D=100d$, respectively. We have checked that this result is insensitive to the choice of $\Delta t_{\rm th}$ in the range $1.0\sqrt{m/k}\lesssim\Delta t_{\rm th}\lesssim 20\sqrt{m/k}$. The behavior for $D=100d$ is analogous to the relative motion of two tracer particles in turbulent flows, which is known as Richardson’s law $\Delta^{2}(t)\sim t^{3}$ Richardson26 ; Monin ; Beffetta02 . This suggests that the flow is in a fully turbulent state for $D/d=100$.

Next, we consider the scattering angle distribution of the mobile particles for $D=7d$. The azimuthal angle is stored when the particles reach the region $|\bm{r}_{i}|=15d$. The behavior of the angular distribution $\rho(\theta)$ for $V/v_{\rm T}\gtrsim 1$ completely differs from that for $V/v_{\rm T}\lesssim 1$ as shown in Fig. 10.

For $V/v_{\rm T}\gg 1$, the angular distribution of scattered particles has a sharp peak around the opening angle, which is the half of the apex angle of the cone of the beam scattered after collisions with the intruder Cheng07 ; Sano12 as shown in Fig. 11. We also note that this opening angle can be explained by phenomenology as in Ref. Cheng07 . Because we use repulsive and elastic particles, the opening angle is expressed as

[TABLE]

for $D<D_{\rm t}$, where $D_{\rm t}$ is the radius of tube of beam particles (see Fig. 1) Cheng07 . When we substitute $D=7d$ and $D_{\rm t}=20d$ into Eq. (8), we obtain $\theta_{0}=0.50$ (rad), which roughly agrees with the simulation results (see Fig. 10). For $V/v_{\rm T}\ll 1$, on the other hand, the particles are scattered in various directions, as shown in Fig. 10.

VI Conclusion

We have numerically studied the particle flows injected as a beam and scattered by a spherical intruder. We have found the crossover from Epstein’s law to Newton’s law, depending on the ratio of the speed to the thermal speed $V/v_{\rm T}$. This crossover can be explained by a simple collision model. The crossover from Epstein’s law to Stokes’ law is also verified as the time evolution of the drag force acting on the intruder. The turbulent-like behavior has also been observed for $V/v_{\rm T}\gg 1$ and $D/d\gg 1$, where the relative displacement between two tracer particles is super-ballistic, which satisfies Richardson’s law.

Although we mainly stress that $V/v_{\rm T}$ is an important control parameter to characterize the particles flows in this paper, it is obvious that $D/d$ is also another important parameter to characterize the flows, particularly for turbulent-like flows. The systematic studies for $D/d$ dependent flows will be one of our future research.

Acknowledgments

The authors thank Hiroshi Watanabe, Kuniyasu Saitoh, Takeshi Kawasaki, and Michio Otsuki for their useful comments. The work of S.T. has been partially supported by the Grant-in-Aid of MEXT for Scientific Research (Grant No. 20K14428). The research of H.H. has been partially supported by the Grant-in-Aid of MEXT for Scientific Research (Grant No. 16H04025) and the Scholarship ISHIZUE 2020 of Kyoto University Research Development Program. The research of the authors was partially supported by the YITP activities (YITP–T–18–03 and YITP–W–18–17). Numerical computation in this work was partially carried out at the Yukawa Institute Computer Facility.

Appendix A Volume fraction dependence of the sound speed

In this Appendix, we show how the sound speed depends on the volume fraction. We also examine whether $v_{\rm T}$ can be replaced by the sound speed.

Because the volume fraction of the mobile particles is finite, the equation of state deviates from that for the ideal gas. Thus, the equation of state at finite density is fitted by the Carnahan–Starling equation Carnahan69

[TABLE]

where $p$ is the pressure. For the adiabatic process, the first law of thermodynamics becomes

[TABLE]

where $L^{3}$ is the volume, $C_{V}=3/2$ is the heat capacity at constant volume. Substituting the equation of state (9) into Eq. (10), the following quantity is conserved along the streamline:

[TABLE]

Then, the sound speed in the adiabatic process is given by

[TABLE]

with

[TABLE]

If we adopt the expansion speed $\sqrt{c(\phi)V}$, the drag force is given by

[TABLE]

which cannot capture the simulation results for the low $V/c$ as shown in Fig. 12, where $c/v_{\rm T}=2.07$ ($\phi=0.20$), $3.22$ ($\phi=0.30$), $5.25$ ($\phi=0.40$), $9.21$ ($\phi=0.50$), and $12.7$ ($\phi=0.55$), respectively.

Appendix B Diameter of the intruder dependence of the drag force for the periodic systems

In the main text, we show how the drag force depends on the diameter when the beam particles are far from the intruder. In this Appendix, however, we investigate the diameter dependence of the intruder on the drag force in the low-speed regime when the intruder is initially located inside stationary beam particles. In addition, we adopt the periodic boundary condition in the flow direction to keep the packing fraction throughout the simulations in this Appendix. The used parameters in this Appendix are listed in Table 2.

Figure 13 shows the drag force to Stokes’ drag. As the system size increases, the drag force approaches Stokes’ law (3). This supports that the drag for the fully confined and periodic system is described by Stokes’ law, which is consistent with the previous studies Vergeles95 . The results in this Appendix are suggestive. Indeed, if we adopt the periodic boundary condition, Stokes’ law is obtained without any difficulty. In other words, Epstein’s law reported in the main text only appears as a transient from the initial hit of beam particles to approaching a steady flow.