# Asphaltene Formation Damage Stimulation by Ultrasound: An Analytical Approach Using Bundle of Tubes Modeling

Authors: Arash Rabbani; Mohammad Hossein Ghazanfari; Mahmood Amani
Journal: Journal of Petroleum Engineering
 (2015)

---

## Abstract

This study presents a novel approach for bundle of tubes modeling of permeability impairment due to asphaltene-induced formation damage attenuated by ultrasound which has been rarely attended in the available literature. Model uses the changes of asphaltene particle size distribution (APSD) as a function of time due to ultrasound radiation, while considering surface deposition and pore throat plugging mechanisms. The proposed model predicts the experimental data of permeability reduction during coinjection of solvent and asphaltenic oil into core with reasonable agreement. Viscosity variation due to sonication of crude oil is used to determine the fluid mobility applied in the model. The results of modeling indicate that the fluid samples exposed to ultrasound may cause much less asphaltene-induced damage inside the porous medium. Sensitivity analysis of the model parameters showed that there is an optimum time period during which the best stimulation efficiency is observed.  The results of this work can be helpful to better understand the role of ultrasound prohibition in dynamic behavior of asphaltene deposition in porous media. Furthermore, the present model could be potentially utilized for modeling of other time-dependent particle induced damages.


---

## Body

## 1. Introduction

Ultrasonic waves have a wide variety of applications in the oil field industries. Special characteristics of these waves such as high frequency and performing cavitation have brought them forward as a unique gadget in many stimulation operations [1]. Ultrasound radiation can change the molecular weight of asphaltenes. Lin et al. [2] applied ultrasonic waves to crack refinery oil residue which contains high amount of asphaltenic components in the presence of abundant hydrogen. Because of the local microscale high pressure and high temperature created by ultrasonic cavitation, asphaltene cracking can be even done in the atmospheric pressure. Dunn and Yen [3] successfully converted about 35% of asphaltene into gas oil and resin fractions under very mild condition using the ultrasound power and its cavitation feature. Kang et al. [4] published experimental results discussing the effects of ultrasound frequency on changing the average molecular weight of asphaltene content and showed that there is no direct relationship between frequency and capability of breaking down the asphaltene molecules.The viscosity of asphaltenic oil can be changed by radiating ultrasonic waves. Islam et al. [5] ran experiments to ascertain the relationship between asphaltene content molecular weight and bitumen viscosity. They showed that oil viscosity significantly increases if containing more asphaltenic heavy components. Consequently, by accepting that ultrasound can crack the asphaltene molecular chains, it can also reduce the viscosity of oil. Chakma and Berruti [6] reported 15% viscosity reduction of heavy asphaltenic oil samples after radiating the ultrasonic waves.The viscosity of colloids is related to the size distribution and stability of their dispersed particles [7]. Zhang et al. [8] studied the effect of ultrasonic waves on the stability of asphaltene emulsion and, by comparing the data before and after ultrasonic treatment, they found that the colloidal stability of asphaltene improves and it is mainly a physical process.Mousavi et al. [9, 10] discussed the applicability of ultrasonic waves on the reduction of asphaltene flocculation. Najafi et al. [11] observed the results of ultrasonic irradiation on asphaltene particle size distribution (APSD) and found that the average sizes of asphaltene particles can be mitigated if they were exposed to the ultrasonic waves for specific time. Rad et al. [12] developed a model describing the behavior of asphaltene flocculation kinetics to show that ultrasonic exposure can reduce the asphaltene flock size after passing an optimum time of sonication.Although it has been shown that ultrasonic waves positively affect the size of asphaltene particles, fundamental understanding of how formation damage induced by asphaltene deposition can be inhibited under the influence of ultrasonic waves is not well understood. Here, a bundle of tubes-based model is developed and used for predicting of formation damage inhibition under the influence of ultrasound. Most of previous dynamic models of permeability impairment cannot take into account the role of ultrasound prohibition on asphaltene precipitation; especially they suffer from the weakness in considering the role of APSD in the model behavior [13–15]. Here, the experimental data of APSD available in the literature measured by confocal microscopy at different ultrasonic exposure times is used for the modeling of permeability impairment in asphaltene-induced damage formations.


## 2. Methodology

### 2.1. Damage Prediction

In order to predict the formation damage, an appropriate permeability modeling should be done as well as considering the particle migration and deposition mechanisms. Stein [17] is almost the pioneer in combination of bundle of tube model and filtration theories to describe particle entrapment mechanisms inside the porous media. Modeling of the asphaltene-induced damage is mainly based on the theories of fine particle straining and clogging inside the porous medium. Here, to predict formation damage, bundle of tubes model is used.


#### 2.1.1. Bundle of Tubes Modeling

One of the prevailing methods in order to model the fluid flow through porous media is bundle of tubes. Considering this model, the permeability of porous media in 1D flow can be expressed as [18]
(1)
k
=
n
p
·
π
·
r
p
4
8
·
τ
·
A
,
where k is permeability, n
p is the total number of tubes, r
p is average radius of pores, A is cross-sectional area, and τ is the tortuosity of the medium. Lanfrey et al. [19] have presented a tortuosity model for randomly packed fixed beds as follows:
(2)
τ
=
0.84
3
π
4
4
/
3
·
1
-
φ
4
/
3
φ
S
2
,
where S is the mean sphericity of grain particles and φ is the porosity. Sphericity is a measure of how spherical (round) an object is. It depends on porosity and specific surface of porous media which could be experimentally obtained [20]. Considering the bundle of tubes model, medium porosity with length L is obtained as follows:
(3)
φ
=
n
p
·
π
·
r
p
2
·
τ
·
L
A
·
L
.In order to find the number and radius of tubes, combining (1) and (3) results, we get
(4)
r
p
=
2
τ
2
k
φ
,
(5)
n
p
=
8
k
·
τ
·
A
π
r
p
2
.


#### 2.1.2. Asphaltene Particle Size Distribution

It seems that the most of previous dynamic models of permeability impairment cannot take into account the role of ultrasound prohibition in asphaltene precipitation; meanwhile they rarely consider the effect of APSD in the model behavior [13–15, 21–25]. It has been expected that APSD majorly affects permeability impairment. It has been shown that the relatively large particles have a noticeable chance to plug the pore throats of the rock [24]. Plugging the throats will decrease the permeability, so considering an appropriate size distribution for asphaltene particles seems indispensable. The size distribution can be estimated by image processing methods applied on the microscopic images of the precipitated samples of asphaltenic oil [11]. Also, there are computational models which can estimate the APSD such as the thermodynamic colloidal model which is a phase behavior model capable of simulating the APSD as a function of the thermodynamics of the system [26–28]. Although many distribution functions can be fitted on asphaltene particle size data, we have chosen Gamma distribution because it is common in the literature [29, 30]. This is partly due to the fact that the gamma distribution provides high flexibility in representing such data [29]. The gamma distribution function is given as
(6)
f
x
=
x
α
-
1
e
-
x
/
β
β
α
Γ
α
,
where α and β are parameters related to the mean and standard deviation, always greater than 0. When σ is the standard deviation of data and μ is the arithmetic mean, σ
2
=
α
β
2 and μ
=
α
β can be expressed. Γ
(
α
) is the gamma function given by the formula
(7)
Γ
α
=
∫
0
∞
t
α
-
1
·
e
-
t
d
t
.


#### 2.1.3. Filtration Theories

Fine migration in porous media occurs in three main patterns. For the particles larger than the average size of the pores, migration cannot be observed. As a result severe particle deposition occurs and permeability reduces significantly (Figure1, Case I) [31].Figure 1
Different mechanisms of filtration. Case I: filter cake, Case II: straining, and Case III: penetration and escaping.The smaller particles have higher chance in penetrating into the porous media. For the penetrated particles, deposition probability should be determined in order to evaluate the plugging phenomena in the porous media caused by the precipitated asphaltene particles. In such systems, solids-to-tube diameter ratios are typically in the range of 0.1–0.3 [28, 32]. In this range, according to the filtration theories, straining filtration is applied. Particle straining can decrease the permeability (Figure 1, Case II) [31]. In the other case when the particles are very small compared to the porous media pore sizes, there is no or very little retention inside the porous medium (Figure 1, Case II). Rege and Fogler [24] have modified and rearranged a probability function for depositing a particle in a pore tube which was firstly developed by Stein [17]:
(8)
P
x
=
4
·
θ
·
x
r
p
2
-
θ
·
x
r
p
3
+
θ
·
x
r
p
4
,
where x is the radius of asphaltene particle [24] and θ is a lumped parameter that takes into account the effect on deposition of several interparticle forces such as gravitational, inertial, hydrodynamic, electric double layer, and Van der Waals forces [33]. If for certain values of θ, r
p, and x, the value of (θ
·
x) is greater than r
p, the probability of capture is set equal to 1.0 and the particle is captured. Higher θ value represents higher probability of particle deposition. Its relationship with some hydrodynamic characteristics of fluid flow inside a pore tube could be expressed as
(9)
θ
=
θ
0
e
(
-
v
/
v
*
)
,
where e is the exponential function, θ
0 is a constant dependent on ionic conditions [24], v is the flow mean velocity inside a tube known as interstitial velocity which is the superficial velocity divided by porosity [32], and v
* is the minimum fluidization velocity that can be calculated by knowing the particle and fluid density [34]. It is obvious from (9) that higher flow velocities reduce the value of θ. This leads to lower deposition probability for all sizes of particles.In lower flow velocities, fluid loses the capability of carrying the asphaltene particles. There are theoretical as well as experimental relationships to estimate the minimum fluidization velocity. An equation for minimum fluidization velocity can be obtained by setting the pressure drop across the deposited bed equal to the weight of the deposited bed per unit area of cross-section, allowing for the buoyant force of the displaced fluid [34]:
(10)
Δ
p
L
=
g
1
-
ε
a
ρ
a
-
ρ
f
,
where Δ
p
/
L is the pressure drop per length of medium, g is gravitational acceleration related to the buoyant force, ρ
a is particles density which in this problem is assumed as the asphaltene density, ρ
f is fluid density, and ε
a is the void fraction of the deposited bed of asphaltene particles when lying above each other. As a simplifying assumption it is reasonable to consider the packing porosity of asphaltene particles (ε
a) as about 0.4 [34]. So, on the one hand, (10) gives the pressure drop estimation for fluid in its minimum fluidization velocity. On the other hand, the Ergun equation for pressure drop in packed beds [35] by considering asphaltene particles as spheres can be rearranged to
(11)
Δ
p
L
=
150
μ
v
x
-
2
1
-
ε
a
2
ε
a
3
+
1.75
ρ
f
1
-
ε
a
x
-
ε
a
3
,
where μ is fluid viscosity and x
- is the average size of particles (asphaltene particles). By setting (10) equal to the pressure drop term in (11) and solving the produced equation for v (superficial velocity), the minimum fluidization velocity is obtained. For very small particles, only the laminar-flow term of the Ergun equation is significant [34]. So, the minimum fluidization velocity can be written as
(12)
v
*
=
g
ρ
a
-
ρ
f
150
μ
ε
a
3
1
-
ε
a
x
-
2
.Here, if low interparticle forces between the migrating particles are assumed,θ
0 can be set equal to 1 [24]. Therefore, the initial capture probability can be predicted for each particle size directly from (8). By deposition of the fraction of asphaltene particles on the pore walls, the initial particle size distribution (PSD) will be changed into new one which contains smaller particles that have higher chance to escape from next pore tubes without being deposited. Thus, the deposited PSD (
f
d
x
) can be calculated by product of the initial precipitated distribution (
f
x
) and probability of adsorption (
P
(
x
)
) which is obtained from (8). PSD of the deposited particles can be expressed as
(13)
f
d
x
=
f
x
·
P
x
.


#### 2.1.4. Formation Damage Dynamic Modeling

To present a proper model for asphaltene deposition, porous media can be divided into some hypothetical sections and deposited PSD should be calculated for each section. When fluid moving through the porous media, a portion of particles are precipitated but not deposited because of the fluidization capability of fluid. Also, portion of the precipitated particles are deposited because of their larger size. Simultaneously, all of the particles whether deposited or precipitated are attaching to their neighbor particles and form larger flocks. Consequently, there is spatial and temporal change in PSD while flowing inside the porous media. We write all PSD functions depending on time(
t
) and section number (
i
). Thus, mass balance equation of particles can be written for the ith section as (Figure 2)
(14)
f
in
x
,
i
,
t
n
in
i
,
t
-
f
d
x
,
i
,
t
n
d
i
,
t
+
f
pr
x
,
t
n
pr
i
,
t
=
f
out
x
,
i
,
t
n
out
i
,
t
,
where x is the size of particle, i is section number, t is the time after starting the injection, f
in is the section input PSD, n
in is the total number of input particles entering the section, n
d the total number of deposited particles in the section, f
pr is the precipitated PSD in the section, n
pr is the total number of precipitated particles, f
out is the section output PSD, and n
out is the total number of section output particles. When we are multiplying f
in
x
,
i
,
t and n
in
(
i
,
t
), we are calculating the number of particles with size x which enter the ith section in time t. Here, f
pr does not depend on the hydrodynamic characteristics of the section. On the other hand, f
out, f
in, and f
d are functions of section number that indicates the dependence of these parameters on the hydraulic specifications of the section. By assuming one-dimensional flow in porous media, the input distribution of each section is equal to the output of the previous section, meaning that
(15)
f
out
x
,
i
,
t
=
f
in
x
,
i
+
1
,
t
+
Δ
t
,
where Δ
t shows the time required for fluid to go through one section of the medium and is only dependent on the rate of injection. Figure 2 illustrates the sections and diagram of asphaltene particles mass balance for one section. We have assumed equal density for all of the asphaltene particles. Since flow is linear and we have considered no backflow, the mass balance equation can be solved for finding the outlet PSD (f
out), section by section. For the first section, the inlet PSD (f
in) is known from microscopic imaging and the deposited PSD (
f
d
) can be found from (13). The total number of inlet particles is obtained by knowing the mass content of asphaltene from IP 143 test and inlet PSD (f
in). The total number of the deposited particles in the section (
n
d
) is obtained by knowing the mass content of the deposited particles and deposited PSD (
f
d
) (equation (18)). Consequently, the number of outlet particles with size x is determined. By repeating the whole procedure for other sizes of particles, we can found the size and number of particles which have not deposited in the section and are available to enter the next section. The outlet distribution of each section is the inlet of the next one and this is repeated to simulate all of the sections.Figure 2
Schematic of damage simulation sections and particle mass balance.Figure3 illustrates 10 sections of a hypothetical porous media simulated for one pore volume injection. The black strips on each section indicate blocked area by asphaltene deposition. It is expected that the first sections experience more asphaltene deposition than the last sections. The first section experiences nine time steps of asphaltene deposition and this leads to massive permeability reduction in this section. The next sections, respectively, face lighter deposition because of the less remained particles in fluid.Figure 3
Corresponding time and section numbers for one pore volume of injected fluid.It is predictable that particle deposition is more severe in the initial sections because when the large particles arrive initial sections they got deposited and do not infect the next sections. In the present study, four sections have been considered for the medium. During the injection, average particles size enhanced because of the asphaltene flocculation, so it is important to find the relation between size increment and time to determine the size of particles involved in deposition process at each moment. Figure4 presents the results of simulation on the specific rock core discussed in the paper and as it can be seen by dividing the medium into four sections, one could conclude that primitive sections contain higher number of the deposited asphaltene particles. Also, Figure 5 shows the mass content of asphaltene particles deposited in each section versus injected pore volume of solvent and crude oil. It seems that the amounts of asphaltene content in each section are approaching a maximum value which is related to the maximum possible amount of asphaltene precipitation for the specific sample of crude oil.Figure 4
Number of asphaltene flocks (particles) in different sections versus injected pore volume of solvent and crude oil predicted by the model.Figure 5
Mass of the deposited asphaltene content in different sections versus injected pore volume of solvent and crude oil predicted by the model.


#### 2.1.5. Hydraulic Radius Approach

The hydraulic radius is a measure of channel flow efficiency. Based on the constant shear stress at the boundary assumption, the hydraulic radius is defined as the ratio of the channel’s void volume to its wetted surface area [36].Hydraulic radius can be expressed asr
H
=
V
p
/
A
p, where V
p is volume of a single pore tube and A
p is its wetted surface area. The deposition effects of asphaltene particles in porous media are simulated by assuming the reduction of hydraulic radius in the medium. As asphaltene particles got deposited inside the pore channels, they reduce the channel void volume and increase the pore surface area. According to the definition of hydraulic radius, both of those events can decrease the value of the hydraulic radius by enhancing the fraction denominator and reducing its nominator. By knowing the deposited asphaltene mass content in each section and asphaltene density, the volume of deposited asphaltene can be simply estimated. Assuming asphaltene particles as spheres, the new surface area of pores as well as new (damaged) hydraulic radius can be calculated. The formulation for calculating total surface added by particle deposition and the volume of deposited asphaltene per section is provided below:
(16)
A
a
=
n
a
∫
0
∞
f
t
x
·
4
π
x
2
·
γ
·
d
x
,
(17)
V
a
=
m
a
ρ
a
,
where A
a is the total asphaltene particles’ surface area (for spherical particles), n
a the total number of asphaltene particles, V
a the total volume of asphaltene particles, m
a the total mass of asphaltene particles inside the porous media, ρ
a is asphaltene density, and γ is a defined coefficient to compensate the difference between the surface area of spheres and real asphaltene particles. Experiments indicate that γ is a number usually between 0.6 and 1.1 [28]. Although asphaltene particles are not perfect sphere, we have assumed that γ is equal to 1 in our simulation.Based on (16), in order to find the total surface of asphaltene particles, total number of particles should be specified. In this regard, by considering that the volumetric summation of all particles represents the asphaltene volumetric content in the section, the total number of asphaltene particles can be calculated from following equation:
(18)
n
a
=
V
a
∫
0
∞
f
t
x
·
4
/
3
π
x
3
d
x
.Thus, the new (damaged) hydraulic radius is obtained from (19). Also, according to equation (20), the value of the equivalent pore tube radius is two times greater than the defined hydraulic radius:
(19)
r
H
~
=
n
p
·
V
p
-
V
a
n
p
·
A
p
-
A
a
,
(20)
r
p
~
=
2
·
r
H
~
,
where r
p
~ is the radius of pore tubes after asphaltene deposition and r
H
~ is the hydraulic radius after reducing the conductivity due to deposition of asphaltene.


#### 2.1.6. Permeability Impairment

It is predictable that asphaltene deposition can reduce the number and radius of pore tubes which conduct the fluid flow. By considering ther
p
~ value for the radius of hypothetical pore tubes after asphaltene deposition and the total void volume of porous media, the new number of total pore tubes can be calculated from (21). Knowing r
p
~ and n
p
~ leads to find the permeability (
k
~
) in the new condition of porous media which asphaltene particles have deposited inside. As presented in (1), attending the Poiseuille equation, new permeability is estimated from (22) as follows:
(21)
n
p
~
=
A
·
L
·
ε
π
r
p
~
2
·
L
·
τ
,
(22)
k
~
=
n
p
~
·
π
·
r
p
~
4
8
·
τ
·
A
.To sum up, the procedure of damage modeling during asphaltene deposition is clarified in Figure6 within the flowchart. The simulation continues until mass balance equation is solved for all sections and all time steps.Figure 6
Flowchart of model calculations.


### 2.2. Ultrasonic Stimulation

#### 2.2.1. Alteration of APSD after Sonication

The size of the asphaltene flocks is an important factor in modeling the deposition because larger flocks are more prone to settle and form deposits [37]. Asphaltene exposure to ultrasonic wave can change its molecular and structural properties [6]. Besides, according to previous studies, sonication can change the kinetic of asphaltene flocculation from irreversible to reversible [12]. Figure 7 presents the PSD of precipitated asphaltene in the mixture asphaltenic oil and n-pentane after 5 min of mixing. In this case, by increasing the sonication time, the average size of particles is reduced.Figure 7
APSD of oil at different times of sonication and after 5 minutes of flocculation.Najafi et al. [11] showed experimental results of asphaltene flocculation and the kinetics of aggregates creation after ultrasonic exposure by adding specific ratio of alkane as the solvent. Ultrasonic radiation was applied using a wave generator with 45 kHz frequency and 75 W output power. Their experiments have shown that sonication firstly causes an increase in the size of aggregates. But after about 90 minutes after the beginning of flocculation, a reduction in aggregate size of sonicated samples is observed. Flocculation time starts at the moment that oil and n-pentane are mixed. Figure 8 illustrates the change of particle size versus the flocculation time in different sonication periods. Asphaltene average particle size increases uniformly for the case of nonsonicated oil by passing the flocculation time. In the sonicated samples, particles’ average size has begun to decrease after about 90 minutes of flocculation time. In a different aspect from Figure 8, Figure 9 presents the effect of sonication time on a single sample after passing different time steps of flocculation. Despite of a slight increment in asphaltene particle size by passing the sonication time, in an overall view, average particles size is decreased significantly after 40 min of sonication.Figure 8
Asphaltene average particle size during the flocculation process at different ultrasonic radiation time.Figure 9
Asphaltene particle size versus time of ultrasonic radiation.


#### 2.2.2. Viscosity Behavior of Sonicated Asphaltenic Crude Oil

Najafi et al. [11] have published experimental results showing the viscosity variation of asphaltenic crude oil under different time intervals of sonication. Figure 10 indicates that sonication reduces the initial viscosity of oil to a minimum value which appears around 30 minutes. The decrease in oil viscosity before ultrasonic radiation time of 20 minutes could be attributed to the produced ultrasonic energy, which may cause partial degradation of resin intramolecular bond sand followed by its partial separation from asphaltene particles. Also, the increase of oil viscosity after 20 minutes of sonication may be related to the breakdown of asphaltene colloids to the finer particles that consequently enhances the solubility and suspension of asphaltene in the crude oil [38]. During this test, the temperature is fixed on 26°C.Figure 10
Changes of dynamic and kinematic viscosity of Sarvak oil versus ultrasonic exposure time [11].


## 2.1. Damage Prediction

In order to predict the formation damage, an appropriate permeability modeling should be done as well as considering the particle migration and deposition mechanisms. Stein [17] is almost the pioneer in combination of bundle of tube model and filtration theories to describe particle entrapment mechanisms inside the porous media. Modeling of the asphaltene-induced damage is mainly based on the theories of fine particle straining and clogging inside the porous medium. Here, to predict formation damage, bundle of tubes model is used.


### 2.1.1. Bundle of Tubes Modeling

One of the prevailing methods in order to model the fluid flow through porous media is bundle of tubes. Considering this model, the permeability of porous media in 1D flow can be expressed as [18]
(1)
k
=
n
p
·
π
·
r
p
4
8
·
τ
·
A
,
where k is permeability, n
p is the total number of tubes, r
p is average radius of pores, A is cross-sectional area, and τ is the tortuosity of the medium. Lanfrey et al. [19] have presented a tortuosity model for randomly packed fixed beds as follows:
(2)
τ
=
0.84
3
π
4
4
/
3
·
1
-
φ
4
/
3
φ
S
2
,
where S is the mean sphericity of grain particles and φ is the porosity. Sphericity is a measure of how spherical (round) an object is. It depends on porosity and specific surface of porous media which could be experimentally obtained [20]. Considering the bundle of tubes model, medium porosity with length L is obtained as follows:
(3)
φ
=
n
p
·
π
·
r
p
2
·
τ
·
L
A
·
L
.In order to find the number and radius of tubes, combining (1) and (3) results, we get
(4)
r
p
=
2
τ
2
k
φ
,
(5)
n
p
=
8
k
·
τ
·
A
π
r
p
2
.


### 2.1.2. Asphaltene Particle Size Distribution

It seems that the most of previous dynamic models of permeability impairment cannot take into account the role of ultrasound prohibition in asphaltene precipitation; meanwhile they rarely consider the effect of APSD in the model behavior [13–15, 21–25]. It has been expected that APSD majorly affects permeability impairment. It has been shown that the relatively large particles have a noticeable chance to plug the pore throats of the rock [24]. Plugging the throats will decrease the permeability, so considering an appropriate size distribution for asphaltene particles seems indispensable. The size distribution can be estimated by image processing methods applied on the microscopic images of the precipitated samples of asphaltenic oil [11]. Also, there are computational models which can estimate the APSD such as the thermodynamic colloidal model which is a phase behavior model capable of simulating the APSD as a function of the thermodynamics of the system [26–28]. Although many distribution functions can be fitted on asphaltene particle size data, we have chosen Gamma distribution because it is common in the literature [29, 30]. This is partly due to the fact that the gamma distribution provides high flexibility in representing such data [29]. The gamma distribution function is given as
(6)
f
x
=
x
α
-
1
e
-
x
/
β
β
α
Γ
α
,
where α and β are parameters related to the mean and standard deviation, always greater than 0. When σ is the standard deviation of data and μ is the arithmetic mean, σ
2
=
α
β
2 and μ
=
α
β can be expressed. Γ
(
α
) is the gamma function given by the formula
(7)
Γ
α
=
∫
0
∞
t
α
-
1
·
e
-
t
d
t
.


### 2.1.3. Filtration Theories

Fine migration in porous media occurs in three main patterns. For the particles larger than the average size of the pores, migration cannot be observed. As a result severe particle deposition occurs and permeability reduces significantly (Figure1, Case I) [31].Figure 1
Different mechanisms of filtration. Case I: filter cake, Case II: straining, and Case III: penetration and escaping.The smaller particles have higher chance in penetrating into the porous media. For the penetrated particles, deposition probability should be determined in order to evaluate the plugging phenomena in the porous media caused by the precipitated asphaltene particles. In such systems, solids-to-tube diameter ratios are typically in the range of 0.1–0.3 [28, 32]. In this range, according to the filtration theories, straining filtration is applied. Particle straining can decrease the permeability (Figure 1, Case II) [31]. In the other case when the particles are very small compared to the porous media pore sizes, there is no or very little retention inside the porous medium (Figure 1, Case II). Rege and Fogler [24] have modified and rearranged a probability function for depositing a particle in a pore tube which was firstly developed by Stein [17]:
(8)
P
x
=
4
·
θ
·
x
r
p
2
-
θ
·
x
r
p
3
+
θ
·
x
r
p
4
,
where x is the radius of asphaltene particle [24] and θ is a lumped parameter that takes into account the effect on deposition of several interparticle forces such as gravitational, inertial, hydrodynamic, electric double layer, and Van der Waals forces [33]. If for certain values of θ, r
p, and x, the value of (θ
·
x) is greater than r
p, the probability of capture is set equal to 1.0 and the particle is captured. Higher θ value represents higher probability of particle deposition. Its relationship with some hydrodynamic characteristics of fluid flow inside a pore tube could be expressed as
(9)
θ
=
θ
0
e
(
-
v
/
v
*
)
,
where e is the exponential function, θ
0 is a constant dependent on ionic conditions [24], v is the flow mean velocity inside a tube known as interstitial velocity which is the superficial velocity divided by porosity [32], and v
* is the minimum fluidization velocity that can be calculated by knowing the particle and fluid density [34]. It is obvious from (9) that higher flow velocities reduce the value of θ. This leads to lower deposition probability for all sizes of particles.In lower flow velocities, fluid loses the capability of carrying the asphaltene particles. There are theoretical as well as experimental relationships to estimate the minimum fluidization velocity. An equation for minimum fluidization velocity can be obtained by setting the pressure drop across the deposited bed equal to the weight of the deposited bed per unit area of cross-section, allowing for the buoyant force of the displaced fluid [34]:
(10)
Δ
p
L
=
g
1
-
ε
a
ρ
a
-
ρ
f
,
where Δ
p
/
L is the pressure drop per length of medium, g is gravitational acceleration related to the buoyant force, ρ
a is particles density which in this problem is assumed as the asphaltene density, ρ
f is fluid density, and ε
a is the void fraction of the deposited bed of asphaltene particles when lying above each other. As a simplifying assumption it is reasonable to consider the packing porosity of asphaltene particles (ε
a) as about 0.4 [34]. So, on the one hand, (10) gives the pressure drop estimation for fluid in its minimum fluidization velocity. On the other hand, the Ergun equation for pressure drop in packed beds [35] by considering asphaltene particles as spheres can be rearranged to
(11)
Δ
p
L
=
150
μ
v
x
-
2
1
-
ε
a
2
ε
a
3
+
1.75
ρ
f
1
-
ε
a
x
-
ε
a
3
,
where μ is fluid viscosity and x
- is the average size of particles (asphaltene particles). By setting (10) equal to the pressure drop term in (11) and solving the produced equation for v (superficial velocity), the minimum fluidization velocity is obtained. For very small particles, only the laminar-flow term of the Ergun equation is significant [34]. So, the minimum fluidization velocity can be written as
(12)
v
*
=
g
ρ
a
-
ρ
f
150
μ
ε
a
3
1
-
ε
a
x
-
2
.Here, if low interparticle forces between the migrating particles are assumed,θ
0 can be set equal to 1 [24]. Therefore, the initial capture probability can be predicted for each particle size directly from (8). By deposition of the fraction of asphaltene particles on the pore walls, the initial particle size distribution (PSD) will be changed into new one which contains smaller particles that have higher chance to escape from next pore tubes without being deposited. Thus, the deposited PSD (
f
d
x
) can be calculated by product of the initial precipitated distribution (
f
x
) and probability of adsorption (
P
(
x
)
) which is obtained from (8). PSD of the deposited particles can be expressed as
(13)
f
d
x
=
f
x
·
P
x
.


### 2.1.4. Formation Damage Dynamic Modeling

To present a proper model for asphaltene deposition, porous media can be divided into some hypothetical sections and deposited PSD should be calculated for each section. When fluid moving through the porous media, a portion of particles are precipitated but not deposited because of the fluidization capability of fluid. Also, portion of the precipitated particles are deposited because of their larger size. Simultaneously, all of the particles whether deposited or precipitated are attaching to their neighbor particles and form larger flocks. Consequently, there is spatial and temporal change in PSD while flowing inside the porous media. We write all PSD functions depending on time(
t
) and section number (
i
). Thus, mass balance equation of particles can be written for the ith section as (Figure 2)
(14)
f
in
x
,
i
,
t
n
in
i
,
t
-
f
d
x
,
i
,
t
n
d
i
,
t
+
f
pr
x
,
t
n
pr
i
,
t
=
f
out
x
,
i
,
t
n
out
i
,
t
,
where x is the size of particle, i is section number, t is the time after starting the injection, f
in is the section input PSD, n
in is the total number of input particles entering the section, n
d the total number of deposited particles in the section, f
pr is the precipitated PSD in the section, n
pr is the total number of precipitated particles, f
out is the section output PSD, and n
out is the total number of section output particles. When we are multiplying f
in
x
,
i
,
t and n
in
(
i
,
t
), we are calculating the number of particles with size x which enter the ith section in time t. Here, f
pr does not depend on the hydrodynamic characteristics of the section. On the other hand, f
out, f
in, and f
d are functions of section number that indicates the dependence of these parameters on the hydraulic specifications of the section. By assuming one-dimensional flow in porous media, the input distribution of each section is equal to the output of the previous section, meaning that
(15)
f
out
x
,
i
,
t
=
f
in
x
,
i
+
1
,
t
+
Δ
t
,
where Δ
t shows the time required for fluid to go through one section of the medium and is only dependent on the rate of injection. Figure 2 illustrates the sections and diagram of asphaltene particles mass balance for one section. We have assumed equal density for all of the asphaltene particles. Since flow is linear and we have considered no backflow, the mass balance equation can be solved for finding the outlet PSD (f
out), section by section. For the first section, the inlet PSD (f
in) is known from microscopic imaging and the deposited PSD (
f
d
) can be found from (13). The total number of inlet particles is obtained by knowing the mass content of asphaltene from IP 143 test and inlet PSD (f
in). The total number of the deposited particles in the section (
n
d
) is obtained by knowing the mass content of the deposited particles and deposited PSD (
f
d
) (equation (18)). Consequently, the number of outlet particles with size x is determined. By repeating the whole procedure for other sizes of particles, we can found the size and number of particles which have not deposited in the section and are available to enter the next section. The outlet distribution of each section is the inlet of the next one and this is repeated to simulate all of the sections.Figure 2
Schematic of damage simulation sections and particle mass balance.Figure3 illustrates 10 sections of a hypothetical porous media simulated for one pore volume injection. The black strips on each section indicate blocked area by asphaltene deposition. It is expected that the first sections experience more asphaltene deposition than the last sections. The first section experiences nine time steps of asphaltene deposition and this leads to massive permeability reduction in this section. The next sections, respectively, face lighter deposition because of the less remained particles in fluid.Figure 3
Corresponding time and section numbers for one pore volume of injected fluid.It is predictable that particle deposition is more severe in the initial sections because when the large particles arrive initial sections they got deposited and do not infect the next sections. In the present study, four sections have been considered for the medium. During the injection, average particles size enhanced because of the asphaltene flocculation, so it is important to find the relation between size increment and time to determine the size of particles involved in deposition process at each moment. Figure4 presents the results of simulation on the specific rock core discussed in the paper and as it can be seen by dividing the medium into four sections, one could conclude that primitive sections contain higher number of the deposited asphaltene particles. Also, Figure 5 shows the mass content of asphaltene particles deposited in each section versus injected pore volume of solvent and crude oil. It seems that the amounts of asphaltene content in each section are approaching a maximum value which is related to the maximum possible amount of asphaltene precipitation for the specific sample of crude oil.Figure 4
Number of asphaltene flocks (particles) in different sections versus injected pore volume of solvent and crude oil predicted by the model.Figure 5
Mass of the deposited asphaltene content in different sections versus injected pore volume of solvent and crude oil predicted by the model.


### 2.1.5. Hydraulic Radius Approach

The hydraulic radius is a measure of channel flow efficiency. Based on the constant shear stress at the boundary assumption, the hydraulic radius is defined as the ratio of the channel’s void volume to its wetted surface area [36].Hydraulic radius can be expressed asr
H
=
V
p
/
A
p, where V
p is volume of a single pore tube and A
p is its wetted surface area. The deposition effects of asphaltene particles in porous media are simulated by assuming the reduction of hydraulic radius in the medium. As asphaltene particles got deposited inside the pore channels, they reduce the channel void volume and increase the pore surface area. According to the definition of hydraulic radius, both of those events can decrease the value of the hydraulic radius by enhancing the fraction denominator and reducing its nominator. By knowing the deposited asphaltene mass content in each section and asphaltene density, the volume of deposited asphaltene can be simply estimated. Assuming asphaltene particles as spheres, the new surface area of pores as well as new (damaged) hydraulic radius can be calculated. The formulation for calculating total surface added by particle deposition and the volume of deposited asphaltene per section is provided below:
(16)
A
a
=
n
a
∫
0
∞
f
t
x
·
4
π
x
2
·
γ
·
d
x
,
(17)
V
a
=
m
a
ρ
a
,
where A
a is the total asphaltene particles’ surface area (for spherical particles), n
a the total number of asphaltene particles, V
a the total volume of asphaltene particles, m
a the total mass of asphaltene particles inside the porous media, ρ
a is asphaltene density, and γ is a defined coefficient to compensate the difference between the surface area of spheres and real asphaltene particles. Experiments indicate that γ is a number usually between 0.6 and 1.1 [28]. Although asphaltene particles are not perfect sphere, we have assumed that γ is equal to 1 in our simulation.Based on (16), in order to find the total surface of asphaltene particles, total number of particles should be specified. In this regard, by considering that the volumetric summation of all particles represents the asphaltene volumetric content in the section, the total number of asphaltene particles can be calculated from following equation:
(18)
n
a
=
V
a
∫
0
∞
f
t
x
·
4
/
3
π
x
3
d
x
.Thus, the new (damaged) hydraulic radius is obtained from (19). Also, according to equation (20), the value of the equivalent pore tube radius is two times greater than the defined hydraulic radius:
(19)
r
H
~
=
n
p
·
V
p
-
V
a
n
p
·
A
p
-
A
a
,
(20)
r
p
~
=
2
·
r
H
~
,
where r
p
~ is the radius of pore tubes after asphaltene deposition and r
H
~ is the hydraulic radius after reducing the conductivity due to deposition of asphaltene.


### 2.1.6. Permeability Impairment

It is predictable that asphaltene deposition can reduce the number and radius of pore tubes which conduct the fluid flow. By considering ther
p
~ value for the radius of hypothetical pore tubes after asphaltene deposition and the total void volume of porous media, the new number of total pore tubes can be calculated from (21). Knowing r
p
~ and n
p
~ leads to find the permeability (
k
~
) in the new condition of porous media which asphaltene particles have deposited inside. As presented in (1), attending the Poiseuille equation, new permeability is estimated from (22) as follows:
(21)
n
p
~
=
A
·
L
·
ε
π
r
p
~
2
·
L
·
τ
,
(22)
k
~
=
n
p
~
·
π
·
r
p
~
4
8
·
τ
·
A
.To sum up, the procedure of damage modeling during asphaltene deposition is clarified in Figure6 within the flowchart. The simulation continues until mass balance equation is solved for all sections and all time steps.Figure 6
Flowchart of model calculations.


## 2.1.1. Bundle of Tubes Modeling

One of the prevailing methods in order to model the fluid flow through porous media is bundle of tubes. Considering this model, the permeability of porous media in 1D flow can be expressed as [18]
(1)
k
=
n
p
·
π
·
r
p
4
8
·
τ
·
A
,
where k is permeability, n
p is the total number of tubes, r
p is average radius of pores, A is cross-sectional area, and τ is the tortuosity of the medium. Lanfrey et al. [19] have presented a tortuosity model for randomly packed fixed beds as follows:
(2)
τ
=
0.84
3
π
4
4
/
3
·
1
-
φ
4
/
3
φ
S
2
,
where S is the mean sphericity of grain particles and φ is the porosity. Sphericity is a measure of how spherical (round) an object is. It depends on porosity and specific surface of porous media which could be experimentally obtained [20]. Considering the bundle of tubes model, medium porosity with length L is obtained as follows:
(3)
φ
=
n
p
·
π
·
r
p
2
·
τ
·
L
A
·
L
.In order to find the number and radius of tubes, combining (1) and (3) results, we get
(4)
r
p
=
2
τ
2
k
φ
,
(5)
n
p
=
8
k
·
τ
·
A
π
r
p
2
.


## 2.1.2. Asphaltene Particle Size Distribution

It seems that the most of previous dynamic models of permeability impairment cannot take into account the role of ultrasound prohibition in asphaltene precipitation; meanwhile they rarely consider the effect of APSD in the model behavior [13–15, 21–25]. It has been expected that APSD majorly affects permeability impairment. It has been shown that the relatively large particles have a noticeable chance to plug the pore throats of the rock [24]. Plugging the throats will decrease the permeability, so considering an appropriate size distribution for asphaltene particles seems indispensable. The size distribution can be estimated by image processing methods applied on the microscopic images of the precipitated samples of asphaltenic oil [11]. Also, there are computational models which can estimate the APSD such as the thermodynamic colloidal model which is a phase behavior model capable of simulating the APSD as a function of the thermodynamics of the system [26–28]. Although many distribution functions can be fitted on asphaltene particle size data, we have chosen Gamma distribution because it is common in the literature [29, 30]. This is partly due to the fact that the gamma distribution provides high flexibility in representing such data [29]. The gamma distribution function is given as
(6)
f
x
=
x
α
-
1
e
-
x
/
β
β
α
Γ
α
,
where α and β are parameters related to the mean and standard deviation, always greater than 0. When σ is the standard deviation of data and μ is the arithmetic mean, σ
2
=
α
β
2 and μ
=
α
β can be expressed. Γ
(
α
) is the gamma function given by the formula
(7)
Γ
α
=
∫
0
∞
t
α
-
1
·
e
-
t
d
t
.


## 2.1.3. Filtration Theories

Fine migration in porous media occurs in three main patterns. For the particles larger than the average size of the pores, migration cannot be observed. As a result severe particle deposition occurs and permeability reduces significantly (Figure1, Case I) [31].Figure 1
Different mechanisms of filtration. Case I: filter cake, Case II: straining, and Case III: penetration and escaping.The smaller particles have higher chance in penetrating into the porous media. For the penetrated particles, deposition probability should be determined in order to evaluate the plugging phenomena in the porous media caused by the precipitated asphaltene particles. In such systems, solids-to-tube diameter ratios are typically in the range of 0.1–0.3 [28, 32]. In this range, according to the filtration theories, straining filtration is applied. Particle straining can decrease the permeability (Figure 1, Case II) [31]. In the other case when the particles are very small compared to the porous media pore sizes, there is no or very little retention inside the porous medium (Figure 1, Case II). Rege and Fogler [24] have modified and rearranged a probability function for depositing a particle in a pore tube which was firstly developed by Stein [17]:
(8)
P
x
=
4
·
θ
·
x
r
p
2
-
θ
·
x
r
p
3
+
θ
·
x
r
p
4
,
where x is the radius of asphaltene particle [24] and θ is a lumped parameter that takes into account the effect on deposition of several interparticle forces such as gravitational, inertial, hydrodynamic, electric double layer, and Van der Waals forces [33]. If for certain values of θ, r
p, and x, the value of (θ
·
x) is greater than r
p, the probability of capture is set equal to 1.0 and the particle is captured. Higher θ value represents higher probability of particle deposition. Its relationship with some hydrodynamic characteristics of fluid flow inside a pore tube could be expressed as
(9)
θ
=
θ
0
e
(
-
v
/
v
*
)
,
where e is the exponential function, θ
0 is a constant dependent on ionic conditions [24], v is the flow mean velocity inside a tube known as interstitial velocity which is the superficial velocity divided by porosity [32], and v
* is the minimum fluidization velocity that can be calculated by knowing the particle and fluid density [34]. It is obvious from (9) that higher flow velocities reduce the value of θ. This leads to lower deposition probability for all sizes of particles.In lower flow velocities, fluid loses the capability of carrying the asphaltene particles. There are theoretical as well as experimental relationships to estimate the minimum fluidization velocity. An equation for minimum fluidization velocity can be obtained by setting the pressure drop across the deposited bed equal to the weight of the deposited bed per unit area of cross-section, allowing for the buoyant force of the displaced fluid [34]:
(10)
Δ
p
L
=
g
1
-
ε
a
ρ
a
-
ρ
f
,
where Δ
p
/
L is the pressure drop per length of medium, g is gravitational acceleration related to the buoyant force, ρ
a is particles density which in this problem is assumed as the asphaltene density, ρ
f is fluid density, and ε
a is the void fraction of the deposited bed of asphaltene particles when lying above each other. As a simplifying assumption it is reasonable to consider the packing porosity of asphaltene particles (ε
a) as about 0.4 [34]. So, on the one hand, (10) gives the pressure drop estimation for fluid in its minimum fluidization velocity. On the other hand, the Ergun equation for pressure drop in packed beds [35] by considering asphaltene particles as spheres can be rearranged to
(11)
Δ
p
L
=
150
μ
v
x
-
2
1
-
ε
a
2
ε
a
3
+
1.75
ρ
f
1
-
ε
a
x
-
ε
a
3
,
where μ is fluid viscosity and x
- is the average size of particles (asphaltene particles). By setting (10) equal to the pressure drop term in (11) and solving the produced equation for v (superficial velocity), the minimum fluidization velocity is obtained. For very small particles, only the laminar-flow term of the Ergun equation is significant [34]. So, the minimum fluidization velocity can be written as
(12)
v
*
=
g
ρ
a
-
ρ
f
150
μ
ε
a
3
1
-
ε
a
x
-
2
.Here, if low interparticle forces between the migrating particles are assumed,θ
0 can be set equal to 1 [24]. Therefore, the initial capture probability can be predicted for each particle size directly from (8). By deposition of the fraction of asphaltene particles on the pore walls, the initial particle size distribution (PSD) will be changed into new one which contains smaller particles that have higher chance to escape from next pore tubes without being deposited. Thus, the deposited PSD (
f
d
x
) can be calculated by product of the initial precipitated distribution (
f
x
) and probability of adsorption (
P
(
x
)
) which is obtained from (8). PSD of the deposited particles can be expressed as
(13)
f
d
x
=
f
x
·
P
x
.


## 2.1.4. Formation Damage Dynamic Modeling

To present a proper model for asphaltene deposition, porous media can be divided into some hypothetical sections and deposited PSD should be calculated for each section. When fluid moving through the porous media, a portion of particles are precipitated but not deposited because of the fluidization capability of fluid. Also, portion of the precipitated particles are deposited because of their larger size. Simultaneously, all of the particles whether deposited or precipitated are attaching to their neighbor particles and form larger flocks. Consequently, there is spatial and temporal change in PSD while flowing inside the porous media. We write all PSD functions depending on time(
t
) and section number (
i
). Thus, mass balance equation of particles can be written for the ith section as (Figure 2)
(14)
f
in
x
,
i
,
t
n
in
i
,
t
-
f
d
x
,
i
,
t
n
d
i
,
t
+
f
pr
x
,
t
n
pr
i
,
t
=
f
out
x
,
i
,
t
n
out
i
,
t
,
where x is the size of particle, i is section number, t is the time after starting the injection, f
in is the section input PSD, n
in is the total number of input particles entering the section, n
d the total number of deposited particles in the section, f
pr is the precipitated PSD in the section, n
pr is the total number of precipitated particles, f
out is the section output PSD, and n
out is the total number of section output particles. When we are multiplying f
in
x
,
i
,
t and n
in
(
i
,
t
), we are calculating the number of particles with size x which enter the ith section in time t. Here, f
pr does not depend on the hydrodynamic characteristics of the section. On the other hand, f
out, f
in, and f
d are functions of section number that indicates the dependence of these parameters on the hydraulic specifications of the section. By assuming one-dimensional flow in porous media, the input distribution of each section is equal to the output of the previous section, meaning that
(15)
f
out
x
,
i
,
t
=
f
in
x
,
i
+
1
,
t
+
Δ
t
,
where Δ
t shows the time required for fluid to go through one section of the medium and is only dependent on the rate of injection. Figure 2 illustrates the sections and diagram of asphaltene particles mass balance for one section. We have assumed equal density for all of the asphaltene particles. Since flow is linear and we have considered no backflow, the mass balance equation can be solved for finding the outlet PSD (f
out), section by section. For the first section, the inlet PSD (f
in) is known from microscopic imaging and the deposited PSD (
f
d
) can be found from (13). The total number of inlet particles is obtained by knowing the mass content of asphaltene from IP 143 test and inlet PSD (f
in). The total number of the deposited particles in the section (
n
d
) is obtained by knowing the mass content of the deposited particles and deposited PSD (
f
d
) (equation (18)). Consequently, the number of outlet particles with size x is determined. By repeating the whole procedure for other sizes of particles, we can found the size and number of particles which have not deposited in the section and are available to enter the next section. The outlet distribution of each section is the inlet of the next one and this is repeated to simulate all of the sections.Figure 2
Schematic of damage simulation sections and particle mass balance.Figure3 illustrates 10 sections of a hypothetical porous media simulated for one pore volume injection. The black strips on each section indicate blocked area by asphaltene deposition. It is expected that the first sections experience more asphaltene deposition than the last sections. The first section experiences nine time steps of asphaltene deposition and this leads to massive permeability reduction in this section. The next sections, respectively, face lighter deposition because of the less remained particles in fluid.Figure 3
Corresponding time and section numbers for one pore volume of injected fluid.It is predictable that particle deposition is more severe in the initial sections because when the large particles arrive initial sections they got deposited and do not infect the next sections. In the present study, four sections have been considered for the medium. During the injection, average particles size enhanced because of the asphaltene flocculation, so it is important to find the relation between size increment and time to determine the size of particles involved in deposition process at each moment. Figure4 presents the results of simulation on the specific rock core discussed in the paper and as it can be seen by dividing the medium into four sections, one could conclude that primitive sections contain higher number of the deposited asphaltene particles. Also, Figure 5 shows the mass content of asphaltene particles deposited in each section versus injected pore volume of solvent and crude oil. It seems that the amounts of asphaltene content in each section are approaching a maximum value which is related to the maximum possible amount of asphaltene precipitation for the specific sample of crude oil.Figure 4
Number of asphaltene flocks (particles) in different sections versus injected pore volume of solvent and crude oil predicted by the model.Figure 5
Mass of the deposited asphaltene content in different sections versus injected pore volume of solvent and crude oil predicted by the model.


## 2.1.5. Hydraulic Radius Approach

The hydraulic radius is a measure of channel flow efficiency. Based on the constant shear stress at the boundary assumption, the hydraulic radius is defined as the ratio of the channel’s void volume to its wetted surface area [36].Hydraulic radius can be expressed asr
H
=
V
p
/
A
p, where V
p is volume of a single pore tube and A
p is its wetted surface area. The deposition effects of asphaltene particles in porous media are simulated by assuming the reduction of hydraulic radius in the medium. As asphaltene particles got deposited inside the pore channels, they reduce the channel void volume and increase the pore surface area. According to the definition of hydraulic radius, both of those events can decrease the value of the hydraulic radius by enhancing the fraction denominator and reducing its nominator. By knowing the deposited asphaltene mass content in each section and asphaltene density, the volume of deposited asphaltene can be simply estimated. Assuming asphaltene particles as spheres, the new surface area of pores as well as new (damaged) hydraulic radius can be calculated. The formulation for calculating total surface added by particle deposition and the volume of deposited asphaltene per section is provided below:
(16)
A
a
=
n
a
∫
0
∞
f
t
x
·
4
π
x
2
·
γ
·
d
x
,
(17)
V
a
=
m
a
ρ
a
,
where A
a is the total asphaltene particles’ surface area (for spherical particles), n
a the total number of asphaltene particles, V
a the total volume of asphaltene particles, m
a the total mass of asphaltene particles inside the porous media, ρ
a is asphaltene density, and γ is a defined coefficient to compensate the difference between the surface area of spheres and real asphaltene particles. Experiments indicate that γ is a number usually between 0.6 and 1.1 [28]. Although asphaltene particles are not perfect sphere, we have assumed that γ is equal to 1 in our simulation.Based on (16), in order to find the total surface of asphaltene particles, total number of particles should be specified. In this regard, by considering that the volumetric summation of all particles represents the asphaltene volumetric content in the section, the total number of asphaltene particles can be calculated from following equation:
(18)
n
a
=
V
a
∫
0
∞
f
t
x
·
4
/
3
π
x
3
d
x
.Thus, the new (damaged) hydraulic radius is obtained from (19). Also, according to equation (20), the value of the equivalent pore tube radius is two times greater than the defined hydraulic radius:
(19)
r
H
~
=
n
p
·
V
p
-
V
a
n
p
·
A
p
-
A
a
,
(20)
r
p
~
=
2
·
r
H
~
,
where r
p
~ is the radius of pore tubes after asphaltene deposition and r
H
~ is the hydraulic radius after reducing the conductivity due to deposition of asphaltene.


## 2.1.6. Permeability Impairment

It is predictable that asphaltene deposition can reduce the number and radius of pore tubes which conduct the fluid flow. By considering ther
p
~ value for the radius of hypothetical pore tubes after asphaltene deposition and the total void volume of porous media, the new number of total pore tubes can be calculated from (21). Knowing r
p
~ and n
p
~ leads to find the permeability (
k
~
) in the new condition of porous media which asphaltene particles have deposited inside. As presented in (1), attending the Poiseuille equation, new permeability is estimated from (22) as follows:
(21)
n
p
~
=
A
·
L
·
ε
π
r
p
~
2
·
L
·
τ
,
(22)
k
~
=
n
p
~
·
π
·
r
p
~
4
8
·
τ
·
A
.To sum up, the procedure of damage modeling during asphaltene deposition is clarified in Figure6 within the flowchart. The simulation continues until mass balance equation is solved for all sections and all time steps.Figure 6
Flowchart of model calculations.


## 2.2. Ultrasonic Stimulation

### 2.2.1. Alteration of APSD after Sonication

The size of the asphaltene flocks is an important factor in modeling the deposition because larger flocks are more prone to settle and form deposits [37]. Asphaltene exposure to ultrasonic wave can change its molecular and structural properties [6]. Besides, according to previous studies, sonication can change the kinetic of asphaltene flocculation from irreversible to reversible [12]. Figure 7 presents the PSD of precipitated asphaltene in the mixture asphaltenic oil and n-pentane after 5 min of mixing. In this case, by increasing the sonication time, the average size of particles is reduced.Figure 7
APSD of oil at different times of sonication and after 5 minutes of flocculation.Najafi et al. [11] showed experimental results of asphaltene flocculation and the kinetics of aggregates creation after ultrasonic exposure by adding specific ratio of alkane as the solvent. Ultrasonic radiation was applied using a wave generator with 45 kHz frequency and 75 W output power. Their experiments have shown that sonication firstly causes an increase in the size of aggregates. But after about 90 minutes after the beginning of flocculation, a reduction in aggregate size of sonicated samples is observed. Flocculation time starts at the moment that oil and n-pentane are mixed. Figure 8 illustrates the change of particle size versus the flocculation time in different sonication periods. Asphaltene average particle size increases uniformly for the case of nonsonicated oil by passing the flocculation time. In the sonicated samples, particles’ average size has begun to decrease after about 90 minutes of flocculation time. In a different aspect from Figure 8, Figure 9 presents the effect of sonication time on a single sample after passing different time steps of flocculation. Despite of a slight increment in asphaltene particle size by passing the sonication time, in an overall view, average particles size is decreased significantly after 40 min of sonication.Figure 8
Asphaltene average particle size during the flocculation process at different ultrasonic radiation time.Figure 9
Asphaltene particle size versus time of ultrasonic radiation.


### 2.2.2. Viscosity Behavior of Sonicated Asphaltenic Crude Oil

Najafi et al. [11] have published experimental results showing the viscosity variation of asphaltenic crude oil under different time intervals of sonication. Figure 10 indicates that sonication reduces the initial viscosity of oil to a minimum value which appears around 30 minutes. The decrease in oil viscosity before ultrasonic radiation time of 20 minutes could be attributed to the produced ultrasonic energy, which may cause partial degradation of resin intramolecular bond sand followed by its partial separation from asphaltene particles. Also, the increase of oil viscosity after 20 minutes of sonication may be related to the breakdown of asphaltene colloids to the finer particles that consequently enhances the solubility and suspension of asphaltene in the crude oil [38]. During this test, the temperature is fixed on 26°C.Figure 10
Changes of dynamic and kinematic viscosity of Sarvak oil versus ultrasonic exposure time [11].


## 2.2.1. Alteration of APSD after Sonication

The size of the asphaltene flocks is an important factor in modeling the deposition because larger flocks are more prone to settle and form deposits [37]. Asphaltene exposure to ultrasonic wave can change its molecular and structural properties [6]. Besides, according to previous studies, sonication can change the kinetic of asphaltene flocculation from irreversible to reversible [12]. Figure 7 presents the PSD of precipitated asphaltene in the mixture asphaltenic oil and n-pentane after 5 min of mixing. In this case, by increasing the sonication time, the average size of particles is reduced.Figure 7
APSD of oil at different times of sonication and after 5 minutes of flocculation.Najafi et al. [11] showed experimental results of asphaltene flocculation and the kinetics of aggregates creation after ultrasonic exposure by adding specific ratio of alkane as the solvent. Ultrasonic radiation was applied using a wave generator with 45 kHz frequency and 75 W output power. Their experiments have shown that sonication firstly causes an increase in the size of aggregates. But after about 90 minutes after the beginning of flocculation, a reduction in aggregate size of sonicated samples is observed. Flocculation time starts at the moment that oil and n-pentane are mixed. Figure 8 illustrates the change of particle size versus the flocculation time in different sonication periods. Asphaltene average particle size increases uniformly for the case of nonsonicated oil by passing the flocculation time. In the sonicated samples, particles’ average size has begun to decrease after about 90 minutes of flocculation time. In a different aspect from Figure 8, Figure 9 presents the effect of sonication time on a single sample after passing different time steps of flocculation. Despite of a slight increment in asphaltene particle size by passing the sonication time, in an overall view, average particles size is decreased significantly after 40 min of sonication.Figure 8
Asphaltene average particle size during the flocculation process at different ultrasonic radiation time.Figure 9
Asphaltene particle size versus time of ultrasonic radiation.


## 2.2.2. Viscosity Behavior of Sonicated Asphaltenic Crude Oil

Najafi et al. [11] have published experimental results showing the viscosity variation of asphaltenic crude oil under different time intervals of sonication. Figure 10 indicates that sonication reduces the initial viscosity of oil to a minimum value which appears around 30 minutes. The decrease in oil viscosity before ultrasonic radiation time of 20 minutes could be attributed to the produced ultrasonic energy, which may cause partial degradation of resin intramolecular bond sand followed by its partial separation from asphaltene particles. Also, the increase of oil viscosity after 20 minutes of sonication may be related to the breakdown of asphaltene colloids to the finer particles that consequently enhances the solubility and suspension of asphaltene in the crude oil [38]. During this test, the temperature is fixed on 26°C.Figure 10
Changes of dynamic and kinematic viscosity of Sarvak oil versus ultrasonic exposure time [11].


## 3. Results and Discussion

### 3.1. Permeability Impairment Predicted in Nonsonication

Hematfar et al. [16] presented experimental results of permeability reduction due to the asphaltene deposition during coinjection of alkane and oil in reservoir rock samples. Fortunately, the PSD of asphaltene versus flocculation time for this oil sample is available [11]. The input data of current asphaltene deposition modeling is obtained from combination of results of Hematfar et al. [16] and Najafi et al. [11]. The properties of rock and fluid are provided in Tables 1 and 2. To simulate their experiment core rock is divided into four sections. While the interstitial velocity is low because of the little value of flow rate (about 1 mL/min), θ is approximately set equal to θ
0, that is, equal to 1. Figure 11 presents the results of modeling in comparison to the experimental results of Hematfar et al. [16]. The results of modeling are in reasonable agreement with the experiments. k
0 is the initial value of permeability in which core is not damaged by the asphaltene particles. Obviously, more injection of oil and solvent leads to more severe formation damage. However, after about 3- or 4-pore volume injection, the permeability of porous media does not decrease any more. This can be because of the erosion effect due to higher interstitial velocity.Table 1
Porous media specifications [16].
Initial permeability
66.3 mD
Initial porosity
0.16
Irreducible water saturation
0.237
Sphericity
0.8
Length
13.8 cm
Core cross area
11.46 cm2Table 2
Specifications of asphaltene [16].
Crude oil API
17.45
Asphaltene content
12.94 wt%
Asphaltene density
1.26 gr/mL
Surface increment coefficient (γ)
1Figure 11
Comparison of model and experimental results in absence of sonication.


### 3.2. Permeability Impairment Predicted by Model for Sonicated Oil

By considering the changes of APSD during the ultrasound exposure, the modeling is again capable of predicting the permeability. Generally, sonication reduces the particle size and this leads to reducing the severity of damage. It should be noticed that, in this simulation, fluid goes under ultrasonic radiation first; then it is injected to the core and these processes do not happen simultaneously. Figure12 presents the results of permeability prediction in different time steps of sonication versus time of flocculation. The model is run for the five-pore volume coinjection of asphaltenic crude oil and solvent and it has been assumed that the ratio of solvent is sufficiently high in order to have all asphaltene particles precipitated in the condition. As it can be seen in Figure 12 the nonsonicated sample shows continuous decrease in the permeability ratio of core whereas the sonicated samples show improvement in the permeability after about 90 min of flocculation. This can be because of the fact that asphaltene flocculation can be reversible when particles are exposed to ultrasound [11]. From a different point of view, Figure 13 shows the permeability value calculated by the model for each time step of flocculation by passing the ultrasonic radiation time. As it is shown in the figure, the sonication time of 20 minutes does not represent an appropriate permeability improvement and in some cases it may intensify the damage degree. There is no deterministic explanation for this behavior around 20 min of sonication; however, it can be said that according to Figure 9, after 20 min of sonication, a local maxima is visible in asphaltene particle average size. Larger particles block the flow path efficiently and reduce the permeability.Figure 12
Changes of predicted permeability ratio versus flocculation time for different sonication periods.Figure 13
Variations of the predicted permeability ratio versus sonication period of samples for different flocculation times.


### 3.3. Effects of Uncertain Parameters of Distribution and Adsorption

It was mentioned that permeability impairment is a function of particle size distribution parameters.α and β are gamma distribution parameters, where β is known as the shape factor of distribution and α is related to the relative most probable point of distribution. Figure 14 illustrates the sensitivity analysis of permeability ratio over α and β. When β increases, permeability ratio falls. It means that more uniform and identical asphaltene particles cause less formation damage and vice versa. α is a parameter which indicates the centralization point of data. Larger alpha means larger particle size which leads to pore plugging and permeability reduction.Figure 14
Sensitivity analysis overα and β, parameters of gamma distribution function for asphaltene size (5-pore volume injection).It seems reasonable to repeat this sensitivity analysis on the parameters containing uncertainty in order to evaluate their effects on the predicted permeability reduction ratio.γ and θ
0 are two uncertain parameters related to asphaltene adsorption. Figure 15 shows the sensitivity analysis of permeability ratio over γ and θ
0. θ
0 indicates the tendency of particles to be adsorbed to the pore wall surface and γ is surface area correction factor when asphaltene particles attach to each other. In a constant permeability ratio, gamma and theta act inversely. It means that, by increasing theta, gamma decreases. Greater theta value means greater adsorption probability and this reduces permeability ratio. When γ decreases it means that asphaltene particles are packed and compressed to each other or to the pore surface and they do not increase the pore surface area during deposition. Thus, hydraulic radius is not significantly affected and permeability ratio increases.Figure 15
Sensitivity analysis on the adsorption parameters,α and β are considered as 3 and 1, respectively (5-pore volume injection).


### 3.4. Mobility Changes after Sonication

Although permeability is important parameter when assessing the sonication positive effects, viscosity changes should be considered. In this regard, we aim to calculate mobility which is the permeability divided by viscosity(
M
=
k
/
μ
). Therefore the optimization of the sonication time should be processed by attending the fact that the goal is to maximize the mobility of flow in the formation. As it has been discussed the viscosity of the sonicated crude oil after about 20 minutes of sonication reaches its minimum value. Also, after 20 min of flocculation, permeability ratio has a minimum. So, the optimum time of sonication cannot be found straightforwardly and requires deeper investigation. Figure 16 shows the mobility of fluid inside the formation versus time of sonication. In the 20th minute of sonication, mobility reaches its local minimum and this is not favorable. It means that, in the 20th minute of sonication, although viscosity is lower than other sonication times (Figure 10), lower permeability of formation restricts the flow (Figure 13) and consequently reduces the mobility.Figure 16
Mobility changes for different times of sonication for five-pore volume injection after 120 min of flocculation.In the field scale, in order to perform ultrasonic stimulation on oil wells, powerful ultrasound generators are sent down hole and start to propagate waves through the near well bore area. The main goal of this treatment is to open the plugged pores and reduce the skin factor. In order to have a successful ultrasonic stimulation job, the probable decrement of mobility which is illustrated in Figure16 should be carefully studied. The model presented in this study can be applicable to predict the effect of ultrasonic stimulation job on asphaltenic reservoir fluid near the well bore area. In order to apply in the field, it is suggested to upgrade the geometry of model from 1D to 3D and consider the attenuation of ultrasound while passing the through porous media.


## 3.1. Permeability Impairment Predicted in Nonsonication

Hematfar et al. [16] presented experimental results of permeability reduction due to the asphaltene deposition during coinjection of alkane and oil in reservoir rock samples. Fortunately, the PSD of asphaltene versus flocculation time for this oil sample is available [11]. The input data of current asphaltene deposition modeling is obtained from combination of results of Hematfar et al. [16] and Najafi et al. [11]. The properties of rock and fluid are provided in Tables 1 and 2. To simulate their experiment core rock is divided into four sections. While the interstitial velocity is low because of the little value of flow rate (about 1 mL/min), θ is approximately set equal to θ
0, that is, equal to 1. Figure 11 presents the results of modeling in comparison to the experimental results of Hematfar et al. [16]. The results of modeling are in reasonable agreement with the experiments. k
0 is the initial value of permeability in which core is not damaged by the asphaltene particles. Obviously, more injection of oil and solvent leads to more severe formation damage. However, after about 3- or 4-pore volume injection, the permeability of porous media does not decrease any more. This can be because of the erosion effect due to higher interstitial velocity.Table 1
Porous media specifications [16].
Initial permeability
66.3 mD
Initial porosity
0.16
Irreducible water saturation
0.237
Sphericity
0.8
Length
13.8 cm
Core cross area
11.46 cm2Table 2
Specifications of asphaltene [16].
Crude oil API
17.45
Asphaltene content
12.94 wt%
Asphaltene density
1.26 gr/mL
Surface increment coefficient (γ)
1Figure 11
Comparison of model and experimental results in absence of sonication.


## 3.2. Permeability Impairment Predicted by Model for Sonicated Oil

By considering the changes of APSD during the ultrasound exposure, the modeling is again capable of predicting the permeability. Generally, sonication reduces the particle size and this leads to reducing the severity of damage. It should be noticed that, in this simulation, fluid goes under ultrasonic radiation first; then it is injected to the core and these processes do not happen simultaneously. Figure12 presents the results of permeability prediction in different time steps of sonication versus time of flocculation. The model is run for the five-pore volume coinjection of asphaltenic crude oil and solvent and it has been assumed that the ratio of solvent is sufficiently high in order to have all asphaltene particles precipitated in the condition. As it can be seen in Figure 12 the nonsonicated sample shows continuous decrease in the permeability ratio of core whereas the sonicated samples show improvement in the permeability after about 90 min of flocculation. This can be because of the fact that asphaltene flocculation can be reversible when particles are exposed to ultrasound [11]. From a different point of view, Figure 13 shows the permeability value calculated by the model for each time step of flocculation by passing the ultrasonic radiation time. As it is shown in the figure, the sonication time of 20 minutes does not represent an appropriate permeability improvement and in some cases it may intensify the damage degree. There is no deterministic explanation for this behavior around 20 min of sonication; however, it can be said that according to Figure 9, after 20 min of sonication, a local maxima is visible in asphaltene particle average size. Larger particles block the flow path efficiently and reduce the permeability.Figure 12
Changes of predicted permeability ratio versus flocculation time for different sonication periods.Figure 13
Variations of the predicted permeability ratio versus sonication period of samples for different flocculation times.


## 3.3. Effects of Uncertain Parameters of Distribution and Adsorption

It was mentioned that permeability impairment is a function of particle size distribution parameters.α and β are gamma distribution parameters, where β is known as the shape factor of distribution and α is related to the relative most probable point of distribution. Figure 14 illustrates the sensitivity analysis of permeability ratio over α and β. When β increases, permeability ratio falls. It means that more uniform and identical asphaltene particles cause less formation damage and vice versa. α is a parameter which indicates the centralization point of data. Larger alpha means larger particle size which leads to pore plugging and permeability reduction.Figure 14
Sensitivity analysis overα and β, parameters of gamma distribution function for asphaltene size (5-pore volume injection).It seems reasonable to repeat this sensitivity analysis on the parameters containing uncertainty in order to evaluate their effects on the predicted permeability reduction ratio.γ and θ
0 are two uncertain parameters related to asphaltene adsorption. Figure 15 shows the sensitivity analysis of permeability ratio over γ and θ
0. θ
0 indicates the tendency of particles to be adsorbed to the pore wall surface and γ is surface area correction factor when asphaltene particles attach to each other. In a constant permeability ratio, gamma and theta act inversely. It means that, by increasing theta, gamma decreases. Greater theta value means greater adsorption probability and this reduces permeability ratio. When γ decreases it means that asphaltene particles are packed and compressed to each other or to the pore surface and they do not increase the pore surface area during deposition. Thus, hydraulic radius is not significantly affected and permeability ratio increases.Figure 15
Sensitivity analysis on the adsorption parameters,α and β are considered as 3 and 1, respectively (5-pore volume injection).


## 3.4. Mobility Changes after Sonication

Although permeability is important parameter when assessing the sonication positive effects, viscosity changes should be considered. In this regard, we aim to calculate mobility which is the permeability divided by viscosity(
M
=
k
/
μ
). Therefore the optimization of the sonication time should be processed by attending the fact that the goal is to maximize the mobility of flow in the formation. As it has been discussed the viscosity of the sonicated crude oil after about 20 minutes of sonication reaches its minimum value. Also, after 20 min of flocculation, permeability ratio has a minimum. So, the optimum time of sonication cannot be found straightforwardly and requires deeper investigation. Figure 16 shows the mobility of fluid inside the formation versus time of sonication. In the 20th minute of sonication, mobility reaches its local minimum and this is not favorable. It means that, in the 20th minute of sonication, although viscosity is lower than other sonication times (Figure 10), lower permeability of formation restricts the flow (Figure 13) and consequently reduces the mobility.Figure 16
Mobility changes for different times of sonication for five-pore volume injection after 120 min of flocculation.In the field scale, in order to perform ultrasonic stimulation on oil wells, powerful ultrasound generators are sent down hole and start to propagate waves through the near well bore area. The main goal of this treatment is to open the plugged pores and reduce the skin factor. In order to have a successful ultrasonic stimulation job, the probable decrement of mobility which is illustrated in Figure16 should be carefully studied. The model presented in this study can be applicable to predict the effect of ultrasonic stimulation job on asphaltenic reservoir fluid near the well bore area. In order to apply in the field, it is suggested to upgrade the geometry of model from 1D to 3D and consider the attenuation of ultrasound while passing the through porous media.


## 4. Conclusions

This work presents a bundle of tubes model for permeability impairment due to asphaltene deposition attenuated by sonication. It can be helpful to better understand the role of ultrasound exposure time in stimulation efficiency and fluid mobility enhancement. Model uses the changes of asphaltene particle size distribution as a function of time due to sonication, while assuming surface deposition and pore throat plugging mechanisms. The probability of asphaltene adsorption on the pore walls has been estimated by using a probability function which depends on pore size as well as the particle size. Based on the results obtained in this work, the following conclusions can be drawn.(i)
By affecting the asphaltene particle size distribution, ultrasounds are capable of remediating the permeability impairment caused by asphaltene formation damage.
(ii)
The proposed model predicts the experimental data of permeability reduction during coinjection of solvent and asphaltenic oil into a core sample with a reasonable agreement.
(iii)
Mobility reaches a minimum value after a specific time period of sonication. Here it is found to be around 20 minutes.
(iv)
Sensitivity analysis showed that the parameters of gamma distribution,α and β, can affect the permeability ratio. Increment of both α and β will reduce the permeability of porous media.
(v)
The results of this work can be helpful to better understand the role of ultrasound prohibition in dynamic behavior of asphaltene deposition in porous media. Also, the approach presented in this work could be utilized for the modeling of similar time-dependent particle induced formation damage such as invasion of mud solid particles into formation.



---

*Source: 804267-2015-02-24.xml*