Modeling the Thixotropic Behavior of Waxy Crude - Industrial

May 16, 2013 - Waxy crude shows viscoelastic, yielding, and thixotropic behaviors below the gelation point. In this paper, to depict those complex rhe...
1 downloads 0 Views 1MB Size
Article pubs.acs.org/IECR

Modeling the Thixotropic Behavior of Waxy Crude Houxing Teng and Jinjun Zhang* National Engineering Laboratory for Pipeline Safety/Beijing Key Laboratory of Urban Oil & Gas Distribution Technology, China University of Petroleum, Beijing 102249, China S Supporting Information *

ABSTRACT: Waxy crude shows viscoelastic, yielding, and thixotropic behaviors below the gelation point. In this paper, to depict those complex rheological behaviors a new thixotropic model is proposed. In the model, the total shear stress is considered to be a combination of an elastic stress and a viscous stress. The elastic stress, assumed to be structure- and strain-dependent, is described by a nonlinear damping function together with a structural parameter. The evolution of structural parameter is described by a new kinetic equation, whose breakdown term is assumed to be dependent on energy dissipation rate rather than on shear rate. The proposed model is validated respectively by the test of stepwise increases in shear rate and the hysteresis loop test for waxy crude and further validated by predicting the transient flow curve of hysteresis loop test with the model parameters obtained from the test of stepwise increases in shear rate. The results showed that the model’s capability of fitting experimental data and prediction is satisfactory, and the model nicely predicts the transition from an “elastically” dominated region to a “viscous” dominated region without a discontinuity in the stress−strain curve.

1. INTRODUCTION Crude oil is a primary and essential energy source throughout the world. Waxy crude oils represent about 20% in the world petroleum reserves produced and pipelined among various nonconventional oils nowadays.1 In China, more than 80% of crude oil produced is waxy crude oil. Waxy crude oil shows complex temperature-dependent rheological behaviors. Below the wax appearance temperature (WAT), molecules of wax begin to precipitate and form solid wax crystals in the crude oil. As the temperature further decreases, the amount of precipitated wax crystals increases, and the precipitated wax crystals interlock with each other and form a three-dimensional spongelike network structure below the gelation temperature, completing the transition of crude oil from sol to colloidal gel.2 Consequently the gelled waxy crude oil shows complex rheological behaviors, such as viscoelasticity, yield stress, and thixotropy.2,3 The experimental results show that the mechanical response of a gelled waxy crude is linear viscoelastic at the beginning of deformation.4 As the shear strain increases, the three-dimensional network structure composed of loose clusters is disrupted into a small number of relatively large flocs or aggregates, losing the connectivity of network structure.2,4 In this process, the viscoelasticity decays, and the mechanical response gradually transfers from being dominated by the elastic effect to being dominated by the viscous effect. After shearing under some high shear rate for a period of time, the internal network structure is nearly totally destroyed, thus the viscoelastic behavior is very weak and even can be ignored, and then the rheological response is mainly thixotropic.4 In summary, the mechanical response of a gelled waxy crude can be divided into the three regions: a preyield (viscoelastic) region, a yielding (viscoelastic−thixotropic) region, and a postyield (thixotropic and/or viscous) region.4−6 The thixotropic behavior of gelled waxy crude has been extensively studied. However, most of the existing thixotropic © XXXX American Chemical Society

models for waxy crude are viscoplastic model, utilizing a yield stress to characterize the yielding behavior.7,8 The viscoplastic models assume that the material begins to flow only if the shear stress is greater than yield stress. For this reason, they can only depict the thixotropic behavior after yielding but cannot account for the initial viscoelastic behavior before the yield point.9 Some of the recent studies attempted to describe both the viscoelastic and thixotropic effects of thixotropic materials with yield stress.7,8,10−12 For the structural kinetic thixotropic model, there are generally two types of models, Type I and Type II, to take into account both the viscoelasticity and the thixotropy of materials.11The Type I model assumes that the total shear stress is comprised of an elastic and a viscous stress, using the elastic stress to represent the elastic behaviors;8,11,13,14 while the Type II model, also named the Maxwell-type model, utilizes the Maxwell9,10 or Jefferys12 mechanical analog to account for the viscoelastic property. The Type I model uses the elastic stress to account for the viscoelastic behavior. The stress equation, which is usually algebraic, does not establish on a mechanical analog but originates from the Bingham model; and to describe the evolution of elastic strain pertaining to the microstructure, an additional kinetic equation is needed. Doraiswamy et al.13 proposed a model for concentrated suspensions, which assumes that the material has an initial elastic behavior before yielding and a purely viscous behavior afterward. Later, Mujumdar et al.8 extended this model to describe the smooth transition of thixotropic materials with shear stress from an elastic dominated phase to a viscous dominated phase. Dullaert and Mewis11 developed a general structural kinetics thixotropic model for thixotropic systems based on inelastic suspending Received: March 27, 2013 Revised: May 4, 2013 Accepted: May 16, 2013

A

dx.doi.org/10.1021/ie400983e | Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Industrial & Engineering Chemistry Research

Article

fluid to describe both the structure-dependent elastic and viscous contributions of particles and medium. At the very beginning of deformation, the shear stress predicted by the Dullaert model11 has elastic and viscous components; while for the Mujumdar model,8 the shear stress only has the elastic component without the viscous component. This is one distinction between them. Nevertheless, many papers show that for waxy crude2−6 and other thixotropic materials with yield stress,8,15,16 the initial mechanical response is elastic. From this aspect, the initial mechanical response predicted by the Mujumdar model seems more realistic and plausible for these materials. The Type II model usually establishes on a mechanical analog composed of a series of ideal spring and dashpot and incorporates structural kinetic theory to describe the thixotropic behavior.17 In this case, thixotropy can be considered as a special case of more general nonlinear viscoelasticity.9 For the viscoelastic information is already embedded in the stress equation; an additional evolution equation for elastic strain is not needed. The mechanical analogs of different models are usually different, and many of the published models have been reviewed in the literature.8,9 More recently, de Souza Mendes10 proposed a novel Maxwell model for elasto-viscoplastic thixotropic fluids. For the kinetic equation for the structural parameter, the steady-state flow curve is used as an input and the breakdown term is assumed to be dependent on stress rather than on shear rate. Subsequently, another similar model is developed by replacing the Maxwell mechanical analog with the Jeffreys mechanical analog.12 Of Type I and Type II methods, the thixotropy review of Mewis and Wagner18 recommends models of Type I; while the recent review of de Souza Mendes and Thompson19 argues that the Type II method is more robust. One main argument that de Souza Mendes and Thompson are against Mewis and Wagner as well as the Type I method is that “yielding in the classical sense does not depend on shear history”. Indeed, there is some controversy among rheologists regarding whether yielding is dependent on shear histories or not. Some scholars,20−22 besides Mewis and Wagner, argue that “yield stress depends on the history of the sample and there is a close link between thixotropy and yield stress, because both phenomena are believed to be caused by the same underlying physics - the internal microstructure”. More recently, Moller et al. showed that this controversy disappears when a clear distinction is made between two types of yield stress fluids: thixotropic and nonthixotropic fluids; the yield stress of thixotropic fluids depends on the shear history of the sample, while for the nonthixotropic ones the yield stress depends only on the shear rate.23 As to waxy crude, it shows yield stress and thixotropy, and the experimental data show that the yield stress does depend on the shear history of the sample.2−4,24 The viscoelasticity of waxy crude, very sensitive and fragile, is different from that of polymer-like materials (e.g., polymer solution, polymer melt), for the internal structure is destroyed very sharply nearby the yield point. The waxy crude’s viscoelasticity mainly refers to the viscoelastic behavior before (and/or nearby) the yield point; the most dominate and important rheological behavior of waxy crude is still thixotropy. The Maxwell-type model seems inappropriate to describe this viscoelasto-plastic behavior. In this paper, a structural kinetics model, belonging to Type I, is proposed to fully depict the viscoelasto-plastic behavior of waxy crude. The new model introduces a damping function to

describe the nonlinear viscoelastic behavior of a gelled waxy crude. The evolution of elastic strain with time and shear rate is described by a specially proposed kinetic equation. For the structural parameter, its evolution is described by a new kinetic equation whose breakdown term is assumed to be dependent on energy dissipation rate, to overcome the shortcoming of the kinetic equation whose breakdown term is assumed to be dependent on shear rate; and this new kinetic equation is compared with the other two kinetic equations whose breakdown term is assumed to be related with shear rate and shear stress. The model is validated by the test of stepwise increases in shear rate and the hysteresis loop test for waxy crude. The yielding process of waxy crude is also analyzed with the proposed model.

2. DEVELOPMENT OF THE THIXOTROPIC MODEL This section describes the assumptions and gives the equations that compose the proposed thixotropic model for waxy crude. Just like the common practice for the structural kinetics models, the proposed model uses a non-negative scaled structural parameter λ to represent the structuring level of the microstructure of a material. It varies between the values of 0 for the completely broken-down structure and 1 for the fully developed structure. As a common practice,11 this paper only considers the one-dimensional version of the model. 2.1. Constitutive Equation. As the model proposed by Mujumdar et al.,8 the total shear stress τ is also supposed to be composed of an elastic stress τe and a viscous stress τv, that is τ = τv + τe. When the network structure is near-continuous (i.e., λ ≈ 1), as argued in the Introduction, the mechanical response of a gelled waxy crude would be similar to that of a Hookean solid, without viscous behavior. Therefore, initially the elastic stress is assumed to obey Hooke’s law, i.e. τe = G0γe, where G0 is the shear modulus of the completely structured material (λ = 1), and γe is the elastic strain of the continuous network structure. As the shear strain increases, the network structure starts to be destroyed, resulting in the decay of shear modulus G. Most of the existing models8,11,17,25 assumed that the shear modulus is dependent on the structural parameter only. Nevertheless, the model proposed by Huang and Lu14 assumes that the shear modulus at different strains is not only controlled by a structural parameter λ but also controlled by a nonlinear damping function. In this work, the shear modulus is assumed to vary proportionally to the structural parameter and a nonlinear damping function, that is G = λG0h(γe). For the viscous stress, the viscosity is characterized by a structure-dependent consistency Δk and a completely unstructured consistency k, and the viscous stress’s dependence on shear rate is described by a kinetic index n1. The relative importance of the elastic stress to the viscous stress is determined by λ.8 Finally, the constitutive equation takes the form as follows: τ = τe + τv = λG0h(γe)γe + (1 − λ) ·(λΔk + k)γ ṅ 1

(1)

The damping function h(γe) is usually adopted to describe the nonlinear viscoelasticity in the separable KBKZ equation.26 Many forms of the damping function have been put forward,27 and the most used are the three proposed by Wagner,28 Laun,29 and Soskey and Winter.30 For waxy crude, we compared these damping functions and found that the Soskey and Winter equation is better than others. This equation also has the advantage of a simple form, with only two parameters. So, here B

dx.doi.org/10.1021/ie400983e | Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Industrial & Engineering Chemistry Research

Article

mainly described by the term g2(λ̇), i.e. in this case eq 3 would reduce to dγe/dt ≈ g2(λ̇)γ̇. Thus, eq 3 shows time dependency p via the g2(λ̇) term. The equations of g1(γ) = e−(p3.γ 4) and g2(λ̇) = −sλ̇β are found to be consistent with the above assumptions, where p3 and p4 are dimensionless parameters, and s and β are empirical parameters. Finally, the evolution of γe is described by the following kinetic equation:

we adopt the Soskey and Winter equation to describe the nonlinear viscoelasticity of waxy crude, and its equation is given as follows h(γe) =

1 1 + p1 ·γe p2

(2)

where p1 and p2 are dimensionless parameters related to the viscoelastic property. 2.2. Kinetic Equation for the Elastic Strain. As for the elastic strain γe, the model proposed by Doraiswamy et al.13 assumes that γe will increase linearly from zero to some maximum value equal to the critical strain at which the yield point locates, and then it will remain constant as long as the deformation process continues in the same direction. The model proposed by Mujumdar et al.8 assumes that γe increases as the material is strained, but it would always be less than, or equal to, the instantaneous critical strain whose evolution with shearing conditions is described by a simple phenomenological equation. While for the model proposed by Dullaert and Mewis,11 the evolution of γe is described by a specially proposed differential equation. As for waxy crude, initially the mechanical response is viscoelastic and similar to that of a Hookean solid, hence it is physical to assume that the elastic strain γe equals the total shear strain γ at the very beginning of deformation. With the increase of γ, the network structure is disrupted into smaller flocs or aggregates. Consequently, the structural parameter λ and shear modulus decreases, and then γe becomes less than γ. Afterward, it can be deduced that (1) when the λ reaches its equilibrium value, γe should also come to its equilibrium, otherwise the shear stress will not level off; (2) when the internal structure rebuilds up as a result of flocs association to the network structure due to Brownian motions, λ and the shear modulus G should increase, while γe should recover as can be inferred. Therefore, the evolution of elastic strain γe is related with λ, shear rate, and time, i.e. having the following form: dγe/dt ∝ f(λ,dλ/dt,γ̇,t). From the above analysis, it can be presumed that (i) when the network structure is near-continuous (λ ≈ 1), the mechanical response would be similar to that of a Hookean solid, thus γe is nearly equal to γ, that is dγe/dt = γ̇; (ii) as the network structure is being disrupted into smaller flocs or aggregates, γe becomes less than γ and gradually transforms from being equal to γ (that is dγe/dt = γ̇) to being related with the evolution of λ (that is dγe/dt ∝ f(dλ/dt,γ̇); (iii) from some shear strain away from the yield point, the evolution of γe is mainly governed by the variation of λ, that is, γe increases with the breakdown of network structure, decreases with the buildup of structure, and levels off when λ comes to its equilibrium value.8 Based on the above interpretation, here we propose the following kinetic equation to describe the evolution of γe: dγe dt

̇ γ̇ = {g1(γ ) + [1 − g1(γ )]g2(λ)}

dγe dt

β

= {g1(γ ) − [1 − g1(γ )]sλ ̇ }γ ̇

with g1(γ ) = e−(p3 ·γ

p4 )

γ=

∫0

t

γ(̇ t ′)dt ′

(4)

In order to introduce the lowest possible number of empirical parameters in the model, it was found that for waxy crude the proposed kinetic equation can be further reduced to the following form: dγe/dt = {g1(γ) − [1 − g1(γ)]sλ̇}γ̇, i.e. to take β = 1; and the value of the parameter s can also be set equal to one, while it is still remained here to ensure the dimensional homogeneity of the equation. Besides, we also found that if p3 and p4 are replaced with p1 and p2, the model’s predictive capacity only worsens very little, but the number of model parameter can be lessened, which also weakens the nonlinear coupling of model parameters. This may be attributed to the fact that the decay of shear modulus G and the evolution of γe are interrelated and even maybe controlled by the some mechanism. Therefore, for the sake of simplicity, parameter p3 and p4 are replaced by p1 and p2 in this paper for waxy crude. Once the shear flow is stopped, the viscous stress τv would reduce to zero instantaneously, and the total shear stress would immediately drop to τe. Then occurs strain recovery and stress relaxation. In this case, for thixotropic materials with yield stress including waxy crude, the viscoelastic effects, i.e. stress relaxation and first normal stress difference, are often marginal phenomena and therefore are often neglected.11 If one wants to consider the evolution of γe and the stress relaxation behavior after the flow is stopped, one can refer to the treatment method employed by Mujumdar et al.8 2.3. Kinetic Equation for the Structural Parameter. As for the kinetic equation describing the evolution of structural parameter λ, in analogy with chemical reaction kinetics, it is usually assumed that its time evolution is controlled by the combined result of structure buildup and breakdown rates. So far many kinetic equations have been put forward, and some of them introduce a common prefactor for the buildup and breakdown terms.18 Phenomenologically, the prefactor generates stretched exponentials for the evolution of λ. The prefactor does not change the equilibrium value of λ at each shear rates. It only changes the descending rate of λ coming to its equilibrium value. By this way, it substantially improves the model’s fitting capability. For example, the model proposed by Dullaert and Mewis11 introduces a prefactor with time as variable, that is 1/tn2. As a consequence, the value of the prefactor at the same time point is the same no matter what the shear rate of start-up flow is. Nevertheless, the descending rate of λ at a higher shear rate of start-up flow is obviously different from that at a lower shear rate. So it seems somewhat inconsistent with the actual evolution of λ. To solve this problem, here we introduce a new prefactor 1/ (1+γn2), which takes shear strain as a variable. This prefactor follows the same functionality of the Perkins-Turner correlation.31 At the same time point, as the shear rate of

(3)

Here, the term g1(γ) on the right-hand side is used to describe the evolution of γe when the mechanical response is mainly dominated by the elastic effect, and its value is required to be close to one when the network structure is nearcontinuous (λ ≈ 1), i.e. in this case eq 3 would reduce to dγe/dt ≈ γ̇. When the elastic behavior is very weak and the mechanical response is mainly dominated by the viscous effect, the g1(γ) value should be very close to zero and the evolution of γe is C

dx.doi.org/10.1021/ie400983e | Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Industrial & Engineering Chemistry Research

Article

start-up flow increases, the shear strain γ gets a larger value, hence the value of the prefactor becomes smaller. Therefore, the descending rate of λ at a higher shear rate is different from that at a lower shear rate. In this sense, the prefactor that takes γ as a variable seems more plausible than that takes time as a variable. The structural parameter λ is now assumed to obey the following kinetic equation dλ 1 = [a(1 − λ) − λf (τ , γ )] ̇ dt 1 + γ n2

the viscoplastic models’ kinetic equations assume that the breakdown term is dependent on shear rate. To take consideration of the initial elastic dominated regime and the thixotropic (and/or viscous) dominated regime (away from the yield point), here we assume that f(τ,γ̇) is a function of a combination of shear stress and shear rate, e.g. the rate of energy dissipation ϕ, which is defined as ϕ = τγ̇ in the simple shear flow. For the constant shear rate condition, the shear stress will increase with time during the elastic regime until it reaches its peak value, thus the energy dissipation rate increases with time, so are the values of f(τ,γ̇) and λf(τ,γ̇). As a result, the breakdown term λf(τ,γ̇) will arrive at its maximum value nearby the peak shear stress point. This is in line with the actual evolution of λ. In view of the above discussion, a new kinetic equation, named kinetic equation I for convenience of further discussion, is proposed as follows

(5)

where γ is the total shear strain; n2 is a positive dimensionless constant; and a is a kinetic constant for structure buildup. In this equation, the first term on the right-hand side exhibits the structure buildup, while the second term expresses the structure breakdown. For the thixotropic materials with yield stress, it is physically reasonable to expect that for the structural parameter (i) the breakdown rate is low at early times when the microstructure is nearly undeformed, and as time elapses and the shear strain increase, the breakdown rate increases; (ii) nearby the yield point the breakdown rate comes to its maximum value; (iii) after that it begins to decay. Yziquel et al.25 as well as Mewis and Wagner18 have summarized that the evolution of λ can be categorized into the following three types of possible kinetic equation. The mostly used one assumes that the structure change is due to the shear rate,8,11,32 and some suppose that the structure evolution is associated with the shear stress,10,12 while only a few kinetic equations assume that the variation of λ is related to a combination of the two, e.g. the rate of energy dissipation.25,33 As analyzed by de Souza Mendes,10,12 if the function f(τ,γ̇) is taken to be dependent on shear rate (e.g., f(τ,γ̇) ∝ γ̇m), f(τ,γ̇) is zero at γ̇ = 0, and increases monotonically with the increase of γ̇. For a completely structured material initially at rest, if it is suddenly subjected to a constant shear rate (artificially assumed here), the stress is zero at t = 0 and increases with time during the elastic regime until it reaches the yield stress. However, the function f(τ,γ̇) does not change with time. As a consequence, the breakdown term λf(τ,γ̇) is maximum at the onset of applied shear rate (because λ = 1, and f(τ,γ̇) is constant) and decreases as λ decreases. This is not consistent with the actual evolution of λ. Therefore, de Souza Mendes10,12 came to the conclusion that it is nonphysical to assume that the destruction of microstructure is a function of shear rate. One thing to note is that the nonphysical responses mentioned above lie in the viscoelastic part before the yield point. Therefore, it does not show up in the thixotropic models (e.g., visco-plastic models) without considering the viscoelastic part before the yield point. To overcome this shortcoming, the two models proposed by de Souza Mendes10,12 assume that the breakdown term of the structural kinetic equation is a function of the instantaneous stress. This is physical because for the fully structured gelled waxy crude, it is the bond between structural units that determines the gel strength, and the external force is the agent that causes structure breaking. However, for the thixotropic dominated regime (away from the yield point) when the network structure is broken down into isolated flocs suspending in the liquid, it is the surface shear stress experienced by an isolated floc that causes flocs breaking, and it has shown that the surface shear stress is proportional to the shear rate raised to a power.7 Therefore, in this regime it is physical to assume that the breakdown term of the structural kinetic equation is dependent on shear rate. This is also the reason why most of

dλ 1 = [a(1 − λ) − bλϕm] dt 1 + γ n2

(6)

where b is a kinetic constant for shear-induced breakdown, and m is a dimensionless constant. 2.4. Parameters of the Model. In summary, the following equations compose the new thixotropic model. The equation of state: τ = λG0h(γe)γe + (1 − λ) ·(λΔk + k)γ ṅ 1

h(γe) =

1 1 + p1 ·γe p 2

The kinetic equation for structural parameter λ: dλ 1 = [a(1 − λ) − bλϕm] dt 1 + γ n2

The kinetic equation for elastic strain γe: dγe dt

= [g1 − (1 − g1)sλ]̇ γ ̇

with g1 = e−(p1·γ

p

2)

γ=

∫0

t

γ(̇ t ′)dt ′

As mentioned above, the parameter s is set to one, and its existence in the equation is to ensure the dimensional homogeneity of the equation. In summary, the proposed model contains ten adjustable parameters, i.e. G0, k, Δk, n1, n2, a, b, m, p1, p2. The value of G0 may be roughly estimated from the initial slope of shear stress vs shear strain curve34 and/or may be estimated by the low frequency plateau value of storage modulus G′.11 Via substituting the equilibrium structural parameter obtained from the structural kinetic equation into the equation of state, the constants in the equation of state (k, Δk, and n1), the ratio b/a, and m can be determined through a least-squares fitting to the steady-state flow curve at shear rates high enough when the elastic effect is vanished (i.e., τe ≈ 0). The values of p1 and p2 can be obtained from a series of largestrain relaxation tests.26 In principle, the values of the remaining parameters (and/or all of the model parameters) can be obtained via fittings to data pertaining to transient flows through a computer program with the nonlinear least-squares method, where the two kinetic equations for λ and elastic strain γe are solved numerically by the fourth order Runge−Kutta discretization method. In this work, the values of all model D

dx.doi.org/10.1021/ie400983e | Ind. Eng. Chem. Res. XXXX, XXX, XXX−XXX

Industrial & Engineering Chemistry Research

Article

parameters are obtained from a least-squares fitting to a transient flow curve.

of G′ is so small that it is verified to zero at the beginning of experiments, which indicates that the oil is purely viscous. When the oil cools to temperatures below the WAT, both of G′ and G″ increase sharply, but G′ increases more quickly than G″. Below the gelation temperature (at which G′ = G″, that is 34.5 °C in this case), G′ becomes greater than G″, meaning that the oil shows more elastic behavior than viscous behavior. At the end of the isothermal process, G′ is well above G″ and the loss angle is within 25° for the four test temperatures, which indicates that the oil has already formed a gel structure before the isothermal thixotropic test. The test result of stepwise increases and decreases in shear rate for Daqing waxy crude at 33 °C is shown in Figure 2. It can

3. EXPERIMENTAL SECTION 3.1. Materials. To investigate the viscoelastic and thixotropic behaviors of waxy crude and to validate the proposed model, a waxy crude oil from China, i.e. the Daqing waxy crude, is used. This crude has a wax content of 24.1%, a wax appearance temperature (WAT) of 39.0 °C, and a specific gravity of 0.8624 at 20 °C. In consideration of the thermal- and shear-history-dependence of rheological behaviors of waxy crude, the oil specimens were pretreated to a temperature of 80 °C to resolve all wax crystals so as to remove the “memory” of the oil for better repeatability and comparability of the experimental data, and then placed at room temperature statically for 48 h before testing. 3.2. Rheological Measurements. All measurements were performed with the coaxial cylinder sensor system (Z41Ti) of HAAKE MARS III rheometer, and the temperature of samples was controlled by a programmable water bath (AC 200) with temperature control accuracy of 0.01 °C. Pretreated oil specimen was heated to 50 °C and held at that temperature for 20 min. The specimen was then loaded into the measuring cylinder preheated at 50 °C, isothermally held for 10 min, and then statically cooled to the test temperature at a cooling rate of 0.5 °C/min. Then the oil specimen was isothermally held for 45 min before testing to let the wax crystal structure fully develop. Meanwhile, oscillatory measurement at a small strain amplitude was carried out in order to obtain information on the gel strength evolution. In this work, a frequency of 1 Hz and a strain amplitude of 0.001 were adopted, which are small enough not to disturb the wax crystal structure. After that, thixotropic measurements were performed isothermally at the test temperature. In this work, the chosen test temperatures are 32, 33, 34, and 35 °C, and a fresh specimen was used for each test. The test temperatures are nearby the gelation temperature, and the moderate deposition of wax crystals on the cylinder surface might give rise to a thin roughness layer, reducing the wall slip tendency.35 So the wall slip phenomenon is not obvious and can be ignored. 3.3. Viscoelastic and Thixotropic Behaviors of Waxy Crude. Figure 1 shows the variation of the storage modulus G′ and loss modulus G″ during the cooling and isothermal process at 33 °C. It can be seen that at temperatures higher than the WAT (39−50 °C), G″ is relatively small (