Transport of non-spherical particles in pipe flow with suction

(1)

Transport of non-spherical particles in pipe

flow with suction

Emil Wångby

Space Engineering, master's level 2020

Luleå University of Technology

(2)

Transport of non-spherical

particles in pipe flow with

suction

Emil W˚

angby

Division of Fluid and Experimental Mechanics Department of Engineering Sciences and Mathematics

Lule˚a University of Technology

Lule˚a, Sweden

June 2020

Supervisor: Hans ˚Akerstedt

(3)

Abstract

The interest of how small non-spherical particles transport behaviour when transported in pipe-flow is of large interest in a variety applications. This kind of theory have been used when studying composite manufacturing and how particles behaves in the human lungs. The main focus is to study the statistical deposition rate in a flow-field with and without capillary action and gravity.

Two kind of particle shapes are of main interest which are prolate and oblate spheroids. In this study the method of vector projection is used to track particle orientation instead of the more common methods of Euler-angles or quaternions. The method of tracking the particle motion used is Lagrangian tracking method which solves the equations of motion for the particles indi-vidually.

When studying particles of nano-scale the importance of the phenomenon called Brownian motion arises. The inclusion if the Brownian motion gives rise to the solving of stochastic differential equations for the particle trans-port. To solve the resulting equations of transport a MATLAB program was developed to using the numerical Euler-Maruyama scheme. Simulations is done with a large amount of particles with a varying particle size and aspect ratio. The deposition results are compared between the different particles shape and sizes.

(4)

CONTENTS

1 Introduction

The interest of micro- and nano-particle transportation behaviour today are wide and used in different kinds of applications. In composite manufac-turing nano-particles are used to reinforce fibers and to give the composite desired functionality or additional properties [7]. With the increase of usage of nano-particles in composite manufacturing the humans are also exposed more to particles of this scale as well. Increasing the human exposure to nano-particles has increased the need to understand non-spherical nano-nano-particles transports in human airways. The understanding how nano-particle trans-portation in the human airways is not only of interest for its health effect of inhaling particles in manufacturing but also medical applications [6].

The transport behaviour particles of these scales are heavily dependent of the fluid flow as well as particles size and shape. Especially Brownian mo-tion and gravitamo-tional force are directly dependent on the shape and size of the particle. In this thesis mainly non-spherical particles, prolate and oblate spheroids, will be covered to study how different shapes effect transport be-haviour. The theory behind the forces acting on these spheroids will be outlined to give an understanding of what effects the particle transportation. The forces that will act on the particle in this theory are: hydrodynamic, Brownian and gravitational. The fluid flow will act with a hydrodynamic force is on the particle which is dependent on particle shape and orienta-tion. The Brownian force is originating from the random collisions with other molecules of the fluid and gives rise to a stochastic force acting on the particle. The gravitational force will act with a constant force in one direc-tion on the particle depending the particle shape and size. In this thesis the Lagrangian tracking method will be used which the motion of each particle is tracked individually. So for the theory in this thesis each force and torque will be solved for each particle individually. Combining all forces acting on the particle will derive the equations of motion for the particle which for this theory is a set of stochastic differential equations. The numerical solving pro-cedure is done using MATLAB through a iterative numerical method called Euler-Maruyama method [9] and is described in section 3.1.

(7)

1.1 Previous work

There is previous work made on the subject of studying transport of micro-and nano-particles for different cases micro-and applications. The Division of Fluid Mechanics at the Department of Applied Physics and Mechanical Engineering

at Lule˚a university of technology has produced two doctoral theses [5, 6]

and several articles covering this subject. These works covers both tracking methods, the Lagrangian method and the Euler method. There has also been different methods used to describe the orientation of the particles, for example Euler angles and quaternions.

1.2 Aim of the thesis

(8)

2 Theory

In this section all theory and equations are being considered for one particle with spheroidal shape. All equations can also be rewritten to define a particle with spherical shape. The calculation domain for the theory is of a pipe with symmetrical cylindrical shape which can be seen in Figure 1. With x-axis along the pipe being symmetry axis and yz-axis being the radial direction.

Figure 1: Illustration of the the pipe simulated in this thesis. (Not to scale)

2.1 Flow Field

The flow field used in this thesis is parabolic, fully developed laminar flow. It is derived from Navier-Stokes equations when being in the limit of small Reynolds number and having axis symmetry[7] and given by

u = U0 1 − δx R 1 −r R 2 ˆ ex+ δU0 1 2 r R − 1 2 r R 2 ˆ er. (2.1)

The first term is the velocity component in the direction of the pipe, which is oriented in the x-direction. The second part term is the radial component

which is driven by the capillary action. U0 is the maximum velocity in the

(9)

2.2 Particle geometry

by dividing by the radial length R of the pipe. The suction velocity at the wall is then described by

V0 = δU0

1

4 δ 1. (2.2)

To describe the flow in Cartesian coordinates the radial direction of the flow is split into its y- and z-components using

ˆ

er = cos (α)ˆey+ sin (α)ˆez α = arctan

z y

. (2.3)

2.2 Particle geometry

Two kinds of particles will be considered in this thesis, prolate and oblate spheroids which can be seen in Figure 2a and 2b. These particles are defined

by an aspect ratio, β, and geometric diameter, df. β is here defined as the

ratio between the length along the symmetry axis to the length perpendicular to the symmetry axis which gives

βprolate =

a

b βoblate =

b

a. (2.4)

The geometric diameter df is the length perpendicular to symmetry axis

df,prolate = 2b and df,oblate = 2a. For both of the particles the polar-axis is

along the z0-axis, that is the major axis for prolate particles and minor axis

for oblate particles.

(a) Prolate (b) Oblate

Figure 2: Illustration of the particles used in the simulations. a) the prolate particle with semi-major axis a semi-minor axis b. b) the oblate particle with semi-major axis a and semi-minor axis b. The polar-axis direction is

(10)

2.3 Coordinate systems

To correctly describe the position and orientation of the particle, two coordi-nate systems has to be defined. One being a global coordicoordi-nate system fixed in space and defined as [x, y, z] and the other one being a local system fixed

to the particle and defined as [x0, y0, z0].

2.3.1 Local to global projection

The particle motion will be described in the global system so any local coor-dinates needs to be transformed from local to global. Transforming a vector from local to global coordinates is done trough splitting the vector into two components. The two components is parallel respectively perpendicular to

the particles orientational vector ˆn

a = a k ˆn + a ⊥ ˆn = ak+ a⊥. (2.5)

Here a is an arbitrary vector and ˆn is the particles orientational vector. To

calculate these two components a vector projection and orthogonal projection

is done onto the particles orientational vector ˆn. Since the orientational

vector ˆn is represented in the global system, the vector is projected from the

local system onto to the global.

The vector projection of a onto ˆn is given by

ak = (a · ˆn) ˆn = ˆn ˆn · a (2.6)

and the orthogonal projection of a onto ˆn is given by

a⊥= a − ( ˆn · a) ˆn = (I − ˆn ˆn) · a. (2.7)

Combining the vector projection and orthogonal vector projection respec-tively the vector can be projected from the local system onto the global system using

a = ˆn ˆn · a + (I − ˆn ˆn) · a. (2.8)

2.4 Equations of motion

(11)

2.4 Equations of motion

2.4.1 Equation of translational motion

The equation of motion for an arbitrary particle in the global coordinate system is given by

d

dt(mv) = ΣF . (2.9)

The right side is the sum of all forces acting on the particle and will be discussed more further on. On the left side, m is the particle mass and v is particles velocity.

Equation (2.9) can be expressed in local coordinate system by:

mdvi0

dt = F

h

i0 + F_ig0 + f_iBr0 (2.10)

Forces acting on the particle are the hydrodynamic F_ih0, Brownian f_iBr0 and

gravitational F_ig0. These forces will be discussed in section 2.5.1, 2.6 and 2.7

respectively.

In the case that Stokes number(St) are small, St 1, the particle inertia may be neglected. Stokes number characterises how a particles behave when suspended in fluid flow. Particles with St 1 do follow the streamlines of the fluid, while particles with St ≈ 1 or higher do not follow the streamlines of the fluid. In this study Stokes number is small and the particle inertia may be neglected [6]. Equation (2.10) may then be simplified to

0 = F_ih0 + F_ig0 + f_iBr0 . (2.11)

2.4.2 Equation of rotational motion

The rotation of an arbitrary particle is given by d

dt(I · ω) = ΣT (2.12)

in the global coordinate system. The right side is the sum of all torques acting on the particle and will be discussed more further on. On the left side, I is the inertia and ω is angular velocity.

The expression of equation (2.12) in the local coordinate system is given by Goldstein [3] as equation (2.13). Torques acting on the particle are

hydro-dynamic, Th i0 and Brownian, T_iBr0 Ii0 dωi0 dt − (Ij0 − Ik0)ωj0ωk0 = T h i0 + T_iBr0 (2.13)

With inertia being neglected due to small Stokes number, equation (2.13) can be simplified to

(12)

2.5 Hydrodynamic force and torque

The fluid in the pipe will contribute with an hydrodynamic force as well as torque on the particle.

2.5.1 Force

The fluid that flows in the pipe will act on the particle with an hydrodynamic force

F_ih0 = πµaK_i0(u_i0 − v_i0). (2.15)

Here µ is the dynamic viscosity of the fluid, a is the geometric radius and ui0

is the fluid velocity and vi0 is the particle velocity.

Shape factor Ki0 depends on the shape, size and orientation of the particle

and is the entries in the diagonal matrix

K =   K⊥ 0 0 0 K⊥ 0 0 0 Kk  . (2.16)

With Ki0 being defined in the local coordinate system, using

K =Kkn ˆˆn + K⊥(I − ˆn ˆn)

(2.17) it is transformed into the global coordinate, thus hydrodynamic force can be expressed as

Fh = πµaKkn ˆˆn + K⊥(I − ˆn ˆn) (u − v). (2.18)

The shape factor K is dependent on the shape, size and orientation of the particle. K for prolate spheroids are given by

K⊥,prolate = 16(β2_{− 1)} 2β2₋₃ √ β2₋₁ln(β +pβ 2_{− 1) + β} (2.19) Kk,prolate = 8(β2_{− 1)} 2β2₋₁ √ β2₋₁ln(β +pβ 2_{− 1) − β} (2.20)

with the orientation perpendicular and parallel respectively. The shape factor for oblate spheroids is given by

(13)

2.5 Hydrodynamic force and torque Kk,oblate= 8(B2_{− 1)} B(B2₋₂₎ √ B2₋₁ tan −1₍√_B2_{− 1) + B} (2.22)

with the orientation perpendicular and parallel respectively. Note in equation

(2.21) and (2.22) that B is used instead of β, were is B = 1

β. This conversion

is needed since the equations for K are defined for an aspect ratio β > 1 or major-axis to minor-axis. Oblate particles are here defined with the aspect ratio β < 1 or minor-axis to major axis [6].

2.5.2 Torque

The fluid will also contribute with an hydrodynamic torque acting on the particle because of its ellipsoid shape. For the case in this thesis the equation for the hydrodynamic torque is given by the following expression

Th = 1 Br ⊥ (1 − ˆn ˆn)(Ω − ω) − 1 Br ⊥ β2_{− 1} β2_{+ 1} ( ˆn × S) · ˆn, (2.23)

which is a derivation based upon the original work by Jeffery [8]. Rotation about the symmetry axis is not included since it is not affecting transport. The matrix S is the strain rate tensor and Ω being the spin tensor and both are derived from the global flow field. The particles angular velocity is given by ω and the particles perpendicular angular mobility in the local coordinate

system is given by B_⊥r. The angular mobility for prolate spheroids is given

from [6] as B_⊥,prolater = 3 2β2₋₁ √ β2₋₁ln(β +pβ 2_{− 1) − β} 16πµa3_(β4_{− 1)} (2.24)

and for oblate spheroids it is given as

B_⊥,oblater = 3h(2−B√ 2)B B2₋₁ tan −1₍√_B2_{− 1) − B}i 16πµa3 1 B2 − B2 . (2.25)

Note in equation (2.25) that B is used instead of β, were is B = _β1. This

conversion is needed since the equation for Br

⊥are defined for an aspect ratio

(14)

2.6 Brownian motion

Brownian motion is the random fluctuating motion caused by the particle interaction with the molecules of the fluid. It will contribute with an Brow-nian force and torque which can be modelled as an independent zero-mean Gaussian white noise process given as [2]:

hfBr

i0 (t)i = 0 hf_iBr0 (t₂)f_jBr0 (t₁)i = πS_it0δ(t₂− t₁) (2.26)

hT_iBr0 (t₁)i = 0 hT_iBr0 (t₁)T_jBr0 (t₂)i = πS_ir0δ(t₁− t₂) (2.27)

With ∆t = (t2− t1) being Brownian time increment which is related to the

correlation time in the collision between the particle and the molecule. Were δ() is the Dirac delta function (Dirac delta of the Brownian increment can be

approximated as δ(t1−t2) ≈ _∆t1 ) Sit0 and S_ir0 is the translational and rotational

spectral intensity respectively given as

S_it0 = 2σT µaK_i0 Sr i0 = 2σT µa Br ⊥ . (2.28)

Here Ki0 is the shape factor defined in section 2.5.1, B_⊥r being the

perpen-dicular angular mobility defined in section 2.5.2, σ the Boltzmann constant and T is the absolute temperature. The equation for the Brownian force in the global coordinate system is then given as

fBr = χt r St ∆t = χ t r 2σT πµa ∆t h pKkn ˆˆn + p K⊥(1 − ˆn ˆn) i . (2.29)

The Brownian torque in the global coordinate system is then given as

TBr = χr r Sr δt = χ r r 2σT πµa δt s 1 Br ⊥ (I − ˆn ˆn). (2.30)

The Brownian time increment was chosen by studying the Brownian diffusion of the particle. This was done by comparing the root mean square of the Brownian drift [6] Drrms = v u u t2 2Dt k+ D⊥t 3 ! in∆t (2.31)

were in is the number of simulation time-steps and with the Brownian

diffu-sion is given as

D_it0 =

σT

πµaKi0

. (2.32)

From comparing the Drrmsto a theoretical value [4] a suitable time increment

(15)

2.7 Gravitation

The force due to gravitation is here defined to be in the positive z-direction. Since the density of the particle is much larger than the density of the sur-rounding fluid, buoyancy is in this case neglected [7].

Gravitational force for prolate and oblate spheroids is then given as equation

F_zg = πρfgd

3

fβ

6 , βprolate > 1 βoblate < 1. (2.33)

With ρf being the density of the particle and g being the gravitational

ac-celeration. Note that β is defined as: β > 1 for prolate spheroids and β < 1 for oblate spheroids.

2.8 Complete set of equations

With all the forces and torques now defined in previous sections, the com-plete set of equations can be defined. Equation (2.11) and (2.14) given in global coordinates can be solved for the particle velocity and rate of change for the normal direction respectively. The result is two different stochastic differential equations(SDEs).

The SDE for translational motion is given by

v = dx dt = u + " 1 pKk ˆ n ˆn + √1 K⊥ (I − ˆn ˆn) # ·   0 0 k   + s 2σT πµa∆t " 1 pKk ˆ n ˆn + √1 K⊥ (I − ˆn ˆn) # · χt_. (2.34) With k being k = gρfd 2 fβ 6µ , βprolate> 1 βoblate< 1 (2.35)

for prolate and oblate spheroids. If the gravitational force is not included, then equation (2.34) is changed to the SDE given by

(16)

2.8 Complete set of equations

The angular motion is expressed as an SDE for the rate of change for the normal direction given as

d ˆn dt = A · ˆn − β2_{− 1} β2_{+ 1} ( ˆnT · S · ˆn) ˆn − S · ˆn + r (2σT ) ∆t pB r ⊥(χr× ˆn) . (2.37)

The variables, χt _{and χ}r_{, are two independent zero-mean Gaussian}

dis-tributed random numbers which is generated at every time-step. And S (referred to strain rate tensor earlier) is the symmetric part and A is the asymmetric part of the gradient of the velocity field were S is given as

S = 1 2 ∇u + (∇u)T= 1 2 ∂ui ∂xj +∂uj ∂xi (2.38) and with A given as

(17)

3 Simulation Procedure

3.1 Euler–Maruyama numerical method

The solving of equation (2.34), (2.36) and (2.37) is done numerically using MATLAB. With the normal direction being solved first and the new normal

direction is updated as well as re-normalized with every time-step, hn. The

re-normalization is done to prevent numerical errors to accumulate with every new time-step calculation. With the new normal direction the position vector is solved and updated for each time-step. The numerical method considered is the Euler-Maruyama which is a method developed for stochastic differential equations (SDE). A general SDE contains an deterministic part and a non-deterministic part and is given by

Xt

dt = f (Xt, t)dt + g(Xt, t)dWt. (3.1)

With f (Xt, t) being the drift coefficient which is the deterministic part and

g(Xt, t) is the diffusion coefficient being the non-deterministic part of the

equation. The non-deterministic part of the equation originate from the stochastic process called standard Weiner process. If the stochastic part of the equation would be excluded it would end up being an ODE. To solve the SDE given in equation (3.1) the Euler-Maruyama method is applied and the solution is on the form [9]

Yn+1 = Yn+ f (Yn)hn+ g(Yn)∆Wn. (3.2)

Were ∆Wn originates from the standard Weiner process and is a random

variable from the normal distribution, N , with mean value of µ = 0 and

variance σ2 = hn.

∆Wn = [Wt+h− Wt] ∼ N (0, hn) ∼

p

hnN (0, 1) (3.3)

3.2 Numerical Setup

In these simulation the simulated pipe is of a length 0.4 m and with a radius

of R = 0.1 mm. The flow field has a maximum velocity of U0 = 0.2 mm/s

and when capillary action is included two different suction strengths are

considered, δ = 5 · 10−5 and δ = 10 · 10−5. These two different suction

strengths result in two different radial velocities at the wall, V0 = 2.5 · 10−9

m/s and V0 = 5 · 10−9 m/s respectively. The fluid viscosity is set to µ = 0.1

Pa with a density of 1000 kg/m3_{. The particles density is ρ}

(18)

3.2 Numerical Setup

and with a mean free path of 0.2 nm. The time-step for the simulation was set to h = 10 ms.

To compare between prolate and oblate particles they are assigned a repre-sentative diameter given as

dre= dfβ

1

3, β

prolate > 1 βoblate< 1 (3.4)

which is the diameter of a sphere with the same volume. In this thesis three

different diameters is studied, dre = 10 µm, 100 nm and 10 nm.

(19)

4 Result

4.1 Jeffery Orbits

When studying the the particles angular motion a distinctive behaviour called Jeffery orbits can be observed for both prolate and oblate particles [8]. This behaviour can be seen by recording how the angle θ changes over time. The angle θ is the angle between the global x-axis along the pipe and the particles

orientational vector ˆn (which also coincides with particles polar axis z0).

The θ and φ angles are illustrated in Figure 3 in relation to the particles

orientation vector ˆn. x y z φ θ x0 y0 z0 ˆ n

Figure 3: Illustration of how the angles θ and φ relate to the particles

orien-tation vector ˆn in the global coordinate system.

(20)

4.1 Jeffery Orbits 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4

x*[m]

0 0.5 1 1.5 2 2.5 3

radians

Angle between local z -axis and global x-axis

Figure 4: Prolate particle released at z = −0.75R and starting angle θ = π/2. Only hydrodynamics is acting on the particle.

In Figure 5 one oblate particle is released at z = −0.75R at the inlet and at an angle θ = 0 which is perpendicular to the flow. Here it is observed that the particle quickly orients itself to be parallel to the flow at an angle θ = π/2 and slows down its rotation. Shortly after the rotations has slowed down the intensifies again and continues to rotate to be perpendicular to the flow again at an angle θ = π. Then quickly rotate back to be parallel to the flow to then again slow down its rotation. The rotations slows down when the particle is parallel to the flow at an angle θ = π/2. This behaviour though unsteady is periodic with in increasing periodicity just as with the prolate particles, the closer the particle gets to the pipe-wall.

0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4

x*[m]

0 0.5 1 1.5 2 2.5 3

radians

Angle between local z -axis and global x-axis

(21)

4.1 Jeffery Orbits

The effect of when both hydrodynamics and Brownian motion is included for prolate and oblate particle respectively is seen in Figure 6 and Figure 7. The new behaviour of a prolate particle can be seen in Figure 6 were it is clear that Brownian motion is effecting the particle rotation. The Jeffery orbits can still be distinguished but rotations are a bit irregular. Though the frequency of the rotations still increases as the particle is getting closer to the wall. 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4

x*[m]

0 0.5 1 1.5 2 2.5 3

radians

Angle between local z -axis and global x-axis

Figure 6: Prolate particle released at z = −0.75R and starting angle θ = π/2. Here both hydrodynamics and Brownian motion is acting on the particle.

(22)

4.1 Jeffery Orbits 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 x*[m] 0 0.5 1 1.5 2 2.5 3 radians

Angle between local z -axis and global x-axis

Figure 7: Oblate particle released at z = −0.75R and starting angle θ = 0. Here both hydrodynamics and Brownian motion is acting on the particle.

0 0.002 0.004 0.006 0.008 0.01 0.012 0.014 0.016 0.018 0.02 x*[m] 0 0.5 1 1.5 2 2.5 3 radians

(a) 0 0.002 0.004 0.006 0.008 0.01 0.012 0.014 0.016 0.018 0.02 x*[m] 0 0.5 1 1.5 2 2.5 3 radians

(b)

(23)

4.2 Deposition rate

The deposition rates seen in the graphs below are determined by calculating the percentage of particles that is considered deposited along the length of the tube. The numerical setup is described in further detail section 3.2.

4.2.1 Spherical particles

In the case when the particle is in the shape of a sphere and have an as-pect ratio of 1, the statistical deposition rate can be seen in Figure 9. It is clear from the graphs that for spherical particles deposition rate increases dramatically with the decrease in size of the particle.

From the graphs it is seen that the overall deposition rate for the particles are increasing the smaller the size of the particle is. While for the micro-scale particle the overall deposition rate lies much lower compared to the two smaller particles. Including capillary action or gravity doesn’t affect the nano-scale particles significantly and there is no clear visible difference seen in the graphs. While for the micro-scale particle the deposition rate increases with a few percentages when including capillary action or gravity.

0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 x*[m] 0 10 20 30 40 50 60 70 80 90 100 Depostion rate[%]

Spherical Capillary Action

df = 1 m

df = 100 nm

df = 10 nm

(a) Capillary action

Spherical Capillary Action with Gravity

df = 1 m

df = 100 nm

df = 10 nm

(b) Capillary action with gravity

Spherical No Capillary Action

df = 1 m df = 100 nm df = 10 nm (c) No capillary action 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 x*[m] 0 10 20 30 40 50 60 70 80 90 100 Depostion rate[%]

Spherical No Capillary Action with Gravity

df = 1 m

df = 100 nm

df = 10 nm

(d) No capillary action with gravity

(24)

4.2 Deposition rate

4.2.2 Prolate particles

For prolate particles the graphs in Figure 10 shows the deposition rate for the different cases. The curves show similar pattern as with spherical particles but with a lower rate of deposition overall.

When comparing the different particle sizes and shape there is a clear dif-ference in deposition behaviour. The most noticeable difdif-ference is the much higher deposition rate when decreasing the particle size. It also seen that an lower aspect ratio has an overall higher deposition rate.

When including capillary action the 10 nm particles the curves doesn’t show a significant difference in deposition behaviour. Meanwhile the micro-scale and 100 nm particles show a difference with a few percentages when including capillary action.

The inclusion of gravity changes the overall deposition rate with a very small change percentage for micro-scale and 100 nm particles. While the smallest particles doesn’t show any significant effect from the inclusion of gravity.

Prolate Capillary Action

= 10 dre = 1 m = 10 dre = 100 nm = 10 d re = 10 nm = 100 dre = 1 m = 100 dre = 100 nm = 100 d re = 10 nm

Prolate Capillary Action with Gravity

= 10 dre = 1 m = 10 dre = 100 nm = 10 d re = 10 nm = 100 dre = 1 m = 100 dre = 100 nm = 100 d re = 10 nm

Prolate No Capillary Action

= 10 dre = 1 m = 10 dre = 100 nm = 10 dre = 10 nm = 100 dre = 1 m = 100 dre = 100 nm = 100 dre = 10 nm (c) No capillary action 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 x*[m] 0 10 20 30 40 50 60 70 80 90 100 Depostion rate[%]

Prolate No Capillary Action with Gravity

= 10 dre = 1 m = 10 dre = 100 nm = 10 dre = 10 nm = 100 dre = 1 m = 100 dre = 100 nm = 100 dre = 10 nm

(25)

4.2 Deposition rate

4.2.3 Oblate particles

For oblate particles the graphs in Figure 11 shows the deposition rate for the different cases. Similarities from the prolate graphs is seen in the curvature pattern and higher deposition rate with smaller particle size. The pattern that a lower aspect ratio results in a higher deposition rate is also repeated here.

When capillary action is included deposition behaviour changes in a similar way as for the prolate particles. The largest impact is seen for micro-scale and 100 nm particles with the inclusion of capillary action. While the 10 nm particles doesn’t show much of a change in behaviour with or without capillary action.

Difference in deposition rate due to gravity very small The micro-scale and 100 nm particles show an impact of a few percentages while the 10 nm par-ticles doesn’t show any significant impact.

Oblate Capillary Action

B = 10 d re = 1 m B = 10 dre = 100 nm B = 10 dre = 10 nm B = 100 d re = 1 m B = 100 dre = 100 nm B = 100 dre = 10 nm

Oblate Capillary Action with Gravity

B = 10 d re = 1 m B = 10 dre = 100 nm B = 10 dre = 10 nm B = 100 d re = 1 m B = 100 dre = 100 nm B = 100 dre = 10 nm

Oblate No Capillary Action

B = 10 dre = 1 m B = 10 d re = 100 nm B = 10 dre = 10 nm B = 100 dre = 1 m B = 100 d re = 100 nm B = 100 dre = 10 nm (c) No capillary action 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 x*[m] 0 10 20 30 40 50 60 70 80 90 100 Depostion rate[%]

Oblate No Capillary Action with Gravity

B = 10 dre = 1 m B = 10 d re = 100 nm B = 10 dre = 10 nm B = 100 dre = 1 m B = 100 d re = 100 nm B = 100 dre = 10 nm

(26)

4.2 Deposition rate

4.2.4 Different suction strengths

In Figure 12 the influence the suction strength has on the 10 nm particles is shown. Here the diameter is kept the same in all four graphs with the aspect ratio being different in each graph while changing the suction strength for each line.

It is seen in Figure 12 that the oblate curves tend to be more the left in each graph with an higher rate of deposition. The prolate curves are more to the right with an overall slightly lower rate of deposition. From this it is seen that neither the prolate or oblate curves changes significantly when changing the suction strength.

0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 x*[m] 0 10 20 30 40 50 60 70 80 90 100 Depostion rate[%] = 10 and B = 10, dre = 10 nm prolate = 0 prolate = 5e-5 prolate = 10e-5 oblate = 0 oblate = 5e-5 oblate = 10e-5

(a) Aspect ratio 10

0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 x*[m] 0 10 20 30 40 50 60 70 80 90 100 Depostion rate[%] = 50 and B = 50, dre = 10 nm prolate = 0 prolate = 5e-5 prolate = 10e-5 oblate = 0 oblate = 5e-5 oblate = 10e-5 (b) Aspect ratio 50 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 x*[m] 0 10 20 30 40 50 60 70 80 90 100 Depostion rate[%] = 75 and B = 75, d_re = 10 nm prolate = 0 prolate = 5e-5 prolate = 10e-5 oblate = 0 oblate = 5e-5 oblate = 10e-5 (c) Aspect ratio 75 0 0.05 0.1 0.15 0.2 0.25 0.3 0.35 0.4 x*[m] 0 10 20 30 40 50 60 70 80 90 100 Depostion rate[%] = 100 and B = 100, d_re = 10 nm prolate = 0 prolate = 5e-5 prolate = 10e-5 oblate = 0 oblate = 5e-5 oblate = 10e-5 (d) Aspect ratio 100

(27)

5 Discussion

From the results presented it can be seen that the deposition rate for oblate particles are overall higher then for prolate particles. The highest deposition rate is seen for the spherical particles. It also follows that a lower aspect ratio increases the deposition rate coinciding with the higher rate for spherical particles.

The effect of the capillary action is overall very low and is greater for larger particles compared to the smaller particles. From this it seems that the Brownian motion would be dominating reason for the particle deposition and that its significance increases with smaller particle size. Comparing the two types of spheroids it shows that the Brownian force is greater for oblate particles than it is for prolate particles, since the overall higher deposition of oblate particles.

The reasoning that the Brownian motion would be greater for oblate particles is also strengthened when looking at the results of the Jeffery orbits. When including Brownian motion the resemblance of Jeffery orbits completely dis-appear. When then looking at the prolate particles with Brownian motion included there is still some resembles of Jeffery orbits left.

The effects from gravity is unsurprisingly most noticeable for the larger par-ticles and even then very small in comparison to the effect of the capillary action. Because of the volume difference between the largest and smallest particle, it it makes sense that the smaller particles isn’t impacted as much as the larger particles by the inclusion of gravity.

Surprising is that the comparisons of the the graphs with different suction strengths doesn’t show much of difference. This leads to more credibility that the Brownian motion would be the dominating cause for the deposition and especially for the smaller particles.

When comparing the current results from Figure 9, 10 and 11 given in this study to the results given in reference [7], the results are vastly different. The results in the reference shows that geometry or size doesn’t impact the deposition behaviour in a major way. This is a great contrast to the current results and it seems more reasonable that there would be a varying deposition rate depending on size and shape of the particle.

(28)

From the comparisons made it is seen that the results in this study are in contrast to both the references compared to. In the the current study the dominating reason for deposition seems to be Brownian motion and that the suction strength doesn’t seem to have a major effect. Were in reference [7] the capillary action is stated to be the dominated reason and [1] shows that suction strength has a major impact.

6 Conclusion

Numerical simulations of particle transport in pipe-flow for different geome-tries and sizes has been done and then compared. The simulations conditions has included the effects of hydrodynamics, Brownian motion and gravitation. It is seen that the shape of the particle has a great importance in the trans-port behaviour of the particle. From the simulations done it is concluded that both capillary action and gravity has an small impact on the deposition rate and that factor for deposition is the Brownian motion. When looking at smallest particles the effect of the capillary action the gravity is not visible at all.

(29)

7 Future work

(30)

A

Appendix

A.1 Sim spherical.m

1 %% Simulation with spherical particles

2 % Program to simulate the depostion rate of spheres in pipe-flow.

3 % This simulation includes hydrodynamics, Brownian motion and gravitational

4 % force. For this script the functions, "flow.m, flow_no_cap.m, and

5 % point_gen.m are needed.

6 % The class const.m is also nedded since it contains global variables used

7 % in many functions.

8

9 % This program is written by Emil W˚angby as part of Master Thesis

10

11 %% Add folder with functions

12 addpath('functions')

13

14 %% Constants

15 df = [1e-6 1e-7 1e-8]; % Geometric diameter of particle, m

16 mu = 0.1; % Dynamic viscosity, Pa*s or kg*m/s

17 T_abs = 293; % Absolute temperature, K

18 kbT = 1.380649e-23*T_abs; % Boltzman constant*Absolute temperature, J

19 I = eye(3); % Identity matrix

20 rho = 1950; % Particle density, kg/m^3

21 g = 9.82; % Gravitational accleration, m/s^2 22

23 %% Loop constants

24 p = 1e3; % Number of particles to simulate

25 h = 10e-3; % Time step

26 breakpoint = 5e8; % Breakpoint index to stop simulate current particle

27

28 %% Loop for each diameter

29 for d=1:3

30 %% Check if save folder exists

31 save_folder = save_path; % Path to save folder

32 if ~exist(save_folder,'dir')

33 mkdir(save_folder) % Create save folder if it doesn't exist

34 end

35

(31)

A.1 Sim spherical.m

37

38 % Geometric radius of particle, m

39 a = df(d)/2;

40 % Gravitational force

41 grav = [0 0 (g*rho*df(d)^2)/(18*mu)]';

42

43

44 %% Preallocation for each particle

45 for j=1:p

46

47 % Counting variables

48 hit_wall=0; % Hit the wall

49 hit_end=0; % End of tube

50 stop_index = 0; % Index at which particle hit the wall

51 stop_pos = 0; % x-position were particle hit wall

52 i = 0; % Count each time step

53

54

55 % Generating random inital positions for y and z

56 [y0,z0] = point_gen();

57

58 % Initial positions mass centre coordinates

59 x = 0; % At circle disk

60 y = y0; % Inside circle disc made by radius R

61 z = z0; % Inside circle disc made by radius R

62

63 %% Calculation of each step

64 while 1

65 % Counting each step

66 i = i + 1;

67

68 % Calculation of the angle to radial vector

69 alpha = atan2(z,y);

70

71 % Flow calculation with capillary action

72 % [uX,uY,uZ] = flow(x,y,z,alpha);

73 % Flow calculation without capillary action

74 % [uX,uY,uZ] = flow_no_cap(y,z);

75

76 % Make temporary vectors for calculation purpose

(32)

A.1 Sim spherical.m

78 U = [uX,uY,uZ]';

79

80 % Calculation of direction vector without gravity

81 X_new = X + (h*U)...

82 + (((kbT/(3*pi*mu*a*h))^(1/2))*((h^(1/2))*randn(1,3)));

83 % Calculation of direction vector with gravity

84 % X_new = X + (h*(U + grav))...

85 % + (((kbT/(3*pi*mu*a*h))^(1/2))*(h^(1/2))*randn(1,3));

86

87 % Return to the original vectors

88 x = X_new(1);

89 y = X_new(2);

90 z = X_new(3);

91

92 % Check if particle touched wall

93 r = ((y^2 + z^2)^(1/2));

94 if r > const.R

95 % Store relevent data

96 hit_wall = 1;

97 stop_index = i + 1;

98 stop_pos = x;

99

100 % Save every particle run in a .mat file

101 save(save_folder + "\sphere_"...

102 + df(d) + "_nr" + j + ".mat",...

103 'stop_index','stop_pos','hit_wall','hit_end');

104 % Stop current particle loop

105 break

106 end

107

108 % Check if particle reached end of tube

109 if x > 0.4

111 hit_end = 1;

113

115 save(save_folder + "\sphere_"...

116 + df(d) + "_nr" + j + ".mat",...

(33)

A.2 Sim prolate.m

120 break

121 end

122

123 if i == breakpoint

126

128 if ~exist(save_folder + '\Not completed runs','dir')

129 mkdir(save_folder + '\Not completed runs')

130 end

131 save(save_folder + "\Not completed runs\sphere_"...

132 + d_re(d) + "_beta" + beta + "_nr" + j + ".mat",...

133 'stop_index','y0','z0','x','y','z','uX','uY','uZ'); 134

136 break 137 end 138 139 end 140 141 end 142 143 end

A.2 Sim prolate.m

1 %% Program to simulate prolate spherioids

2 % Program to simulate the depostion rate of prolate spheroids in pipe-flow.

4 % force. The particle orientation is described using vector projection. For

5 % this script the functions, "flow.m, flow_var.m and point_gen.m are needed.

6 % The class const.m is also need since it contains global variables

7 % used in many functions.

8

10

(34)

A.2 Sim prolate.m

13

14 %% Constants

15 beta_all = [10 100 50 75]; % Aspect ratio of particle, major to minor axis

16 d_re = [1e-6 1e-7 1e-8]; % Representative diameter of particle, m

17 mu = 0.1; % Dynamic viscosity, Pa*s or kg/m*s

24 %% Loop constants

28

29 %% Loop for each aspect ratio

30 for b = 1:2

31 beta = beta_all(b); % Current aspect ratio for loop

32

34 for d=1:3

36 save_folder = save_path;% Path to save folder

38 mkdir(save_folder)

39 end

40

41 %% Particle specific constants

42

43 % Representative diameter to gemoteric diameter, m

44 df = d_re(d)*beta^(-1/3);

46 a = df/2;

47

48 % Gravitational force

49 grav = [0 0 (g*rho*df^2*beta)/(6*mu)]';

50

51 % Stokes diameter perpendicular

52 k_perp = 16*(beta^2-1)/...

(35)

A.2 Sim prolate.m

54 *log(beta+(beta^2-1)^(1/2))+beta);

55 % Stokes diameter parallel

56 k_par = 8*(beta^2-1)/...

57 ((2*beta^2-1)/(beta^2-1)^(1/2)...

58 *log(beta+(beta^2-1)^(1/2))-beta);

59

60 % Just for rotation

61 B_perp = 3*(((2*beta^2-1)/(beta^2-1)^(1/2))...

62 *log(beta+(beta^2-1)^(1/2))-beta)/(16*pi*mu*a^3*(beta^4-1));

63

65

66 for j=1:p

67

74

77

82

83 % Initial position of normal vector

84 nX = 0;

85 nY = 0;

86 nZ = 1;

87

89 while 1

90

92 i = i + 1;

93

(36)

A.2 Sim prolate.m

96

97 % Flow calculation with capillary action

98 [uX,uY,uZ] = flow(x,y,z,alpha);

99

100 % Asymmtric and Symmetric part of the velocity gradient

101 [A,S,jac] = flow_var(x,y,z,alpha);

102

104 N = [nX,nY,nZ]';

105 X = [x,y,z]';

106 U = [uX,uY,uZ]';

107

108 % Calculation of normal vector

109 N_new = N + (h*(A*N - (((beta^2-1)/(beta^2+1))...

110 *((N'*S*N).*N'-S*N)))) + ((2*kbT/h)^(1/2)...

111 *(B_perp)^(1/2).*cross(h^(1/2).*randn(3,1),N));

112 % Re-normalize

113 N_new = N_new/norm(N_new);

114

116 X_new = X + (h*U)...

117 + (((((2*kbT)/(pi*mu*a*h))^(1/2))...

118 *((1/((k_par)^(1/2)))*(N*N')+(1/((k_perp)^(1/2)))...

119 *(I-(N*N'))))*(h^(1/2).*randn(3,1)));

121 % X_new = X + (h*(U + ((1/k_par)*(N*N')+(1/k_perp)...

122 % *(I-(N*N')))*grav)) + (((((2*kbT)/(pi*mu*a*h))^(1/2))...

123 % *((1/((k_par)^(1/2)))*(N*N')+(1/((k_perp)^(1/2)))...

124 % *(I-(N*N'))))*(h^(1/2).*randn(3,1)));

125

127 x = X_new(1); 128 y = X_new(2); 129 z = X_new(3); 130 nX = N_new(1); 131 nY = N_new(2); 132 nZ = N_new(3); 133

(37)

A.2 Sim prolate.m

136 if r > const.R

138 hit_wall = 1;

140 stop_pos = x;

141

143 save(save_folder + "\prolate_"...

146

148 break

149 end

150

152 if x > const.l_tube

154 hit_end = 1;

156

158 save(save_folder + "\prolate_"...

161

163 break

164 end

165

166 % If loop runs for too long, stop it

170

174 end

175 save(save_folder + "\Not completed runs\prolate_"...

(38)

A.3 Sim oblate.m

180 break 181 end 182 183 end 184 185 end 186 187 end 188 189 end

A.3 Sim oblate.m

1 %% Program to simulate oblate spherioids

2 % Program to simulate the depostion rate of oblate spheroids in pipe-flow.

4 % force. The particle orientation is described using vector projection. For

5 % this script the functions, "flow.m, flow_var.m and point_gen.m are needed.

6 % The class const.m is also need since it contains global variables

7 % used in many functions.

8

10

12 addpath('functions')

13

14 %% Constants

15 Beta_all = [10 100 50 75]; % Aspect ratio of particle, major to minor axis

16 d_re = [1e-6 1e-7 1e-8]; % Representative diameter of particle, m

17 mu = 0.1; % Dynamic viscosity, Pa*s or kg/m*s

(39)

A.3 Sim oblate.m

28

29 %% Loop for each aspect ratio

30

31 for b = 1:2

32 Beta = Beta_all(b); % Current aspect ratio for loop

33 beta = 1/Beta; % Aspect ratio length to diameter

34

36 for d=1:3

38 save_folder = save_path;% Path to save folder

40 mkdir(save_folder)

41 end

42

43 %% Particle specific constants

44

45 % Representative diameter to geometric diameter, m, beta < 1

46 df = d_re(d)*beta^(-1/3);

48 a = df/2;

49

50 % Gravitational force, beta < 1

51 grav = [0 0 (g*rho*df^2*beta)/(mu*6)]';

52

53 % Stokes diameter perpendicular, beta > 1

54 k_perp = (16*(Beta^2-1)/...

55 ((Beta*(3*Beta^2-2))/...

56 (Beta^2-1)^(1/2)*atan((Beta^2-1)^(1/2))-Beta));

57

58 % Stokes diameter parallel, beta > 1

59 k_par = (8*(Beta^2-1)/...

60 (Beta*((Beta^2-2))/(Beta^2-1)^(1/2)...

61 *atan((Beta^2-1)^(1/2))+Beta));

62

63 % Angular mobility perpendicular, beta > 1

64 B_perp = (3*((((2-Beta^2)*Beta)/(Beta^2-1)^(1/2))...

(40)

A.3 Sim oblate.m

66 (16*pi*mu*a^3*((1/Beta^2)-Beta^2)));

67

69

70 for j=1:p

71

78

81

86

87 % Initial position of normal vector

88 nX = 1;

89 nY = 0;

90 nZ = 0;

91

92

94 while 1

95

97 i = i + 1;

98

99 % Calculation of the angle to radial vector

101

102 % Flow calculation

103 [uX,uY,uZ] = flow(x,y,z,alpha);

104

105 % Asymmtric and Symmetric part of the velocity gradient

(41)

A.3 Sim oblate.m

107

109 N = [nX,nY,nZ]';

110 X = [x,y,z]';

111 U = [uX,uY,uZ]';

112

113 % Calculation of normal vector

114 N_new = N + (h*(A*N - (((beta^2-1)/(beta^2+1))...

115 *((N'*S*N).*N'-S*N)))) + ((2*kbT/h)^(1/2)...

116 *(B_perp)^(1/2).*cross(h^(1/2).*randn(3,1),N));

117 % Re-normalize

118 N_new = N_new/norm(N_new);

119

121 X_new = X + (h*U)...

122 + (((((2*kbT)/(pi*mu*a*h))^(1/2))...

123 *((1/((k_par)^(1/2)))*(N*N')+(1/((k_perp)^(1/2)))...

124 *(I-(N*N'))))*(h^(1/2).*randn(3,1)));

126 % X_new = X + (h*(U + ((1/k_par)*(N*N')+(1/k_perp)...

127 % *(I-(N*N')))*grav)) + (((((2*kbT)/(pi*mu*a*h))^(1/2))...

128 % *((1/((k_par)^(1/2)))*(N*N')+(1/((k_perp)^(1/2)))...

129 % *(I-(N*N'))))*(h^(1/2).*randn(3,1)));

130

132 x = X_new(1); 133 y = X_new(2); 134 z = X_new(3); 135 nX = N_new(1); 136 nY = N_new(2); 137 nZ = N_new(3); 138

140 r = ((y^2 + z^2)^(1/2));

141 if r > const.R

143 hit_wall = 1;

145 stop_pos = x;

146

(42)

A.3 Sim oblate.m

148 save(save_folder + "\oblate_"...

149 + d_re(d) + "_beta" + Beta + "_nr" + j + ".mat",...

151

153 break

154 end

155

157 if x > const.l_tube

159 hit_end = 1;

161

163 save(save_folder + "\oblate_"...

166

168 break

169 end

170

171 % If loop runs for too long, stop it

175

179 end

180 save(save_folder + "\Not completed runs\oblate_"...

185 break

186 end

187

(43)

A.4 functions/flow.m 189 190 end 191 192 end 193 194 end

A.4 functions/flow.m

1 function [ux,uy,uz] = flow(xi,yi,zi,alpha)

2 %FLOW Caclulating flow in x,y,z-direction

3 % From given radius of tube, inlet velocity and capilary action, the

4 % flow in each point is calculated. The radial velocity is split into

5 % cartesian coordinates from cylindrical with use of the angle alpha

6

7 % Redefine variables to make code more easily readable

8 R = const.R; % Radius of tube, [m]

9 U = const.U; % Inlet flow velocity [m/s]

10 delta = const.delta; % Suction strength <<1 []

11

12 % Calculating radial position from y- and z-position

13 r = (yi^2+zi^2)^(1/2);

14

15 % Flow in the x-direction

16 ux = U*(1-delta*(xi/R))*(1-(r/R)^2);

17 % Flow in the y-direction

18 uy = (0.5*delta*U*((r/R)-0.5*(r/R)^3))*cos(alpha);

19 % Flow in the z-direction

20 uz = (0.5*delta*U*((r/R)-0.5*(r/R)^3))*sin(alpha);

21

22 end

A.5 functions/flow var.m

1 function [A,S,jac] = flow_var(xi,yi,zi,alpha)

2 %FLOW_VAR Calculation of Symmetric and Asymmetric part

3 % of the velocity gradient

(44)

A.6 functions/point gen.m

5 % the numerical Jacobi matrix for the flow field.

6 % The Jacobi matrix is calculated from the symbolic toolbox in MATLAB.

7 % From the numerical Jacobi matrix the Symmetric and Asymmetric part

8 % of the velocity gradient can be calculated as output.

9

10 % Redefine variables to make code more easily readable

11 R = const.R; % Radius of tube, [m]

12 U = const.U; % Inlet flow velocity [m/s]

13 delta = const.delta; % Suction strength <<1 []

14

15 % Calculation of the numerical Jacobi matrix with the results from the

16 % symbolic calculation of the Jacobi matrix

17 jac=[ (U*delta*((yi^2 + zi^2)/R^2 - 1))/R,(2*U*yi*((delta*xi)/R - 1))...

18 /R^2,(2*U*zi*((delta*xi)/R - 1))/R^2; 0, (U*delta*cos(alpha)...

19 *(yi/(R*(yi^2 + zi^2)^(1/2)) - (3*yi*(yi^2 + zi^2)^(1/2))/...

20 (2*R^3)))/2, (U*delta*cos(alpha)*(zi/(R*(yi^2 + zi^2)^(1/2)) -...

21 (3*zi*(yi^2 + zi^2)^(1/2))/(2*R^3)))/2; 0, (U*delta*sin(alpha)...

22 *(yi/(R*(yi^2 + zi^2)^(1/2)) - (3*yi*(yi^2 + zi^2)^(1/2))/...

23 (2*R^3)))/2, (U*delta*sin(alpha)*(zi/(R*(yi^2 + zi^2)^(1/2))...

24 - (3*zi*(yi^2 + zi^2)^(1/2))/(2*R^3)))/2];

25

26 % From the numerical Jacobi matrix the Symmetric and Asymmetric part

27 % of the velocity gradient can be calculated

28 S = (1/2)*(jac+jac'); % Symmetric part of the velocity gradient

29 A = (1/2)*(jac-jac'); % Asymmetric part of the velocity gradient

30

31 end

A.6 functions/point gen.m

1 function [y0,z0] = point_gen()

2 %POINT_GEN Generate a random point on circle disc

3 % Using random number generator rand, wich generates a random number on

4 % interval (0,1) to create a random point inside a circle with given

5 % radius. Given is a centre point at origo.

6

7

8 % Creating a random distance to point given radius

9 r = const.R*rand;

(45)

A.7 functions/const.m

11 angle = 2*pi*rand;

12

13 % Creating a random y-coordinate given radom point and angle

14 y0 = r*cos(angle) + 0;

15 % Creating a random z-coordinate given radom point and angle

16 z0 = r*sin(angle) + 0;

17

18 end

A.7 functions/const.m

1 classdef const

2 % A collection of variables that is to be used in main

3 % script and in functions made to be able to acces those

4 % variables more easily.

5

6 properties (Constant)

7 R = 0.1e-3; % Radius of tube, [m]

8 U = 0.2e-3; % Inlet flow velocity [m/s]

9 l_tube = 0.4; % Length of tube [m]

10 delta = 5e-5; % Suction strength <<1 []

11 % delta = 10e-5 % Suction strength <<1 []

12 % delta = 0; % Suction strength <<1 []

13 end

(46)

REFERENCES

References

[1] Hans O. ˚Akerstedt. “Transport and Deposition of Large Aspect Ratio

Prolate and Oblate Spheroidal Nanoparticles in Cross Flow”. In: Pro-cesses 7.12 (2019), p. 886. issn: 2227-9717. doi: 10.3390/pr7120886. url: https://browzine.com/articles/359684367.

[2] Fa-Gung Fan and Goodarz Ahmadi. “WALL DEPOSITION OF SMALL

ELLIPSOIDS FROM TURBULENT AIR FLOWS—A BROWNIAN DYNAMICS SIMULATION”. In: Journal of Aerosol Science 31.10 (2000), pp. 1205–1229. issn: 0021-8502. doi: https : / / doi . org / 10 . 1016 / S0021 - 8502(00 ) 00018 - 5. url: http : / / www . sciencedirect . com / science/article/pii/S0021850200000185.

[3] H. Goldstein, C.P. Poole, and J.L. Safko. Classical Mechanics. Addison

Wesley, 2002. isbn: 9780201657029.

[4] P. S. Grassia, E. J. Hinch, and L. C. Nitsche. “Computer simulations of

Brownian motion of complex systems”. In: Journal of Fluid Mechanics 282 (1995), pp. 373–403. doi: 10.1017/S0022112095000176.

[5] Sofie H¨ogberg. “Modeling nanofiber transport and deposition in human

airways”. PhD thesis. 2010. isbn: 978-91-7439-171-8. url: http://urn. kb.se/resolve?urn=urn:nbn:se:ltu:diva-17105.

[6] Elise Holmstedt. “Modelling Transport of Non-Spherical Particles in

Small Channel Flow”. PhD thesis. 2016. isbn: 978-91-7583-765-9. url: http://urn.kb.se/resolve?urn=urn:nbn:se:ltu:diva-60402.

[7] Elise Holmstedt, Hans O ˚Akerstedt, and T Staffan Lundstr¨om.

“Mod-elling transport and deposition of non-spherical micro- and nano-particles in composites manufacturing”. In: Journal of Reinforced Plastics and Composites 37.8 (2018), pp. 507–519. doi: 10.1177/0731684417753741. eprint: https://doi.org/10.1177/0731684417753741. url: https: //doi.org/10.1177/0731684417753741.

[8] George Barker Jeffery and Louis Napoleon George Filon. “The motion of

(47)

REFERENCES

[9] Thomas Schaffter. “Numerical Integration of SDEs: A Short Tutorial”.

Transport of non-spherical particles in pipe flow with suction

Transport of non-spherical particles in pipe

flow with suction

Emil Wångby

Transport of non-spherical

particles in pipe flow with

suction

Emil W˚

angby

Abstract

Contents

1

Introduction

1.1

Previous work

1.2

Aim of the thesis

2

Theory

2.1

Flow Field

2.2

Particle geometry

2.3

Coordinate systems

2.4

Equations of motion

2.5

Hydrodynamic force and torque

2.6

Brownian motion

2.7

Gravitation

2.8

Complete set of equations

3

Simulation Procedure

3.1

Euler–Maruyama numerical method

3.2

Numerical Setup

4

Result

4.1

Jeffery Orbits

x*[m]

radians

Angle between local z -axis and global x-axis

x*[m]

radians

Angle between local z -axis and global x-axis

x*[m]

radians

Angle between local z -axis and global x-axis

4.2

Deposition rate

5

Discussion

6

Conclusion

7

Future work

A

Appendix

A.1

Sim spherical.m

A.2

Sim prolate.m

A.3

Sim oblate.m

A.4

functions/flow.m

A.5

functions/flow var.m

A.6

functions/point gen.m

A.7

functions/const.m

References