Wax Deposition Modeling with Considerations of Non-Newtonian

Mar 27, 2017 - ... of Non-Newtonian Characteristics: Application on Field-Scale Pipeline ..... accounting for the effect of shear stress at wall on th...
5 downloads 0 Views 2MB Size
Article pubs.acs.org/EF

Wax Deposition Modeling with Considerations of Non-Newtonian Characteristics: Application on Field-Scale Pipeline Sheng Zheng,† Mohamed Saidoun,‡ Thierry Palermo,§ Khalid Mateen,‡ and H. Scott Fogler*,† †

Department of Chemical Engineering, University of Michigan, 2300 Hayward Street, Ann Arbor, Michigan 48109, United States TOTAL E&P, 1201 Louisiana St. Suite 1800, Houston, Texas 77210, United States § TOTAL E&P, Avenue Larribau, 64018 Pau Cedex, France ‡

S Supporting Information *

ABSTRACT: Wax removal by pigging is costly in sub-sea oil production. Cost-effective scheduling of pigging can be achieved based on the deposition rate predicted by wax deposition models. Conventional wax deposition models predict wax deposition rates on the basis of Newtonian fluid mechanics. Such an approach can become invalid for highly waxy crude oils with nonNewtonian rheology. In this investigation, different simulation techniques, including large eddy simulation, Reynolds-averaged Naiver−Stokes equations, and the law of the wall, were applied to model non-Newtonian pipe flow. It was discovered that the law of the wall method is the best method to calculate the velocity profile, shear stress and the turbulent momentum diffusivity in turbulent non-Newtonian pipe flow of waxy oil. An enhanced wax deposition model considering the non-Newtonian characteristics of waxy oil using the law of the wall method was developed and applied to predict wax deposition rates in a fieldscale pipeline.

1. INTRODUCTION 1.1. Wax Deposition in Sub-sea Pipelines and Wax Deposition Modeling. During sub-sea oil transportation, the crude oil is cooled by the seawater. Dissolved wax molecules contained by the crude oil will precipitate when the temperature of the fluid is below the wax appearance temperature (WAT). Precipitation of wax is expected to be the most profound at the pipe wall where the lowest temperature is encountered among all radial locations, generating a radial concentration gradient of the dissolved wax and a net flux of dissolved wax toward the wall.1 Radial diffusion of wax molecules causes wax to deposit on the inner pipe wall, reducing the effective pipe diameter and posing severe risks to oil production. Wax remediation techniques including mechanical removal, thermal insulation, chemical treatment, etc. need to be implemented to address the wax deposition issue during production.2 The design of these wax remediation techniques, such as the pigging operation, is a nontrivial task. It should be noted that frequent pigging will generate significant operational costs while insufficient pigging will lead to thick and hard deposits that are impossible to remove.3 Therefore, determination of a proper pigging frequency is crucial to sub-sea oil production. Determination of a proper pigging frequency relies on estimated wax deposition and aging rates during the design phase of oil field development. Wax deposition is a process depending on various parameters, including the oil flow rate, oil and ambient temperatures, shear stress, oil composition, and so on.4 In order to investigate the effects of the interplay between these parameters on the wax deposition rate, wax deposition models need to be constructed. Academic wax deposition models have been developed on the basis of the mechanism of molecular diffusion.1 Virtually all exiting wax deposition models simulate the heat loss of the oil flow to the ambient and the radial diffusion of wax in order to predict the wax deposition rate in pipelines.5−8 A majority of the © 2017 American Chemical Society

wax deposition modeling studies have been carried out based on oils with relatively low wax contents, e.g., 15 wt%, and it is expected that such waxy oils can become highly non-Newtonian due to the precipitation of wax at temperatures below the WAT.9 The effect of the non-Newtonian fluid characteristics on the modeling of flow, heat transfer, and mass transfer is not well-understood. Applications of wax deposition models to real-world oil transportation pipelines are extremely rare in published investigations. Moreover, a large number of adjustable parameters are usually required to achieve reasonable predictions at field scale due to the large uncertainties associated with field operations.9,10 Continuing efforts have been dedicated to improving wax deposition modeling capabilities, especially upscaling from laboratory-scale predictions to field applications. 1.2. Non-Newtonian Characteristics of Waxy Crude Oil. A dynamic yield stress and shear thinning characteristics can be observed for waxy crude oils when the temperature is below the WAT. Such unique non-Newtonian characteristics originate from the precipitation of wax at temperatures below the WAT and the resulting wax-in-oil suspension. The complex rheology of waxy crude oil below the WAT has been extensively studied.11−15 These investigations show that various constitutive equations, including the Cross model,16 the Casson model,17 and the Herschel−Bulkley model,18 can usually be used to fit the nonNewtonian flow curves of waxy crude oils. Furthermore, it has been found that the non-Newtonian rheology of waxy crude oil has a significant impact in the modeling of pipeline restart.19,20 Received: February 18, 2017 Revised: March 24, 2017 Published: March 27, 2017 5011

DOI: 10.1021/acs.energyfuels.7b00504 Energy Fuels 2017, 31, 5011−5023

Article

Energy & Fuels Unfortunately, the role of non-Newtonian fluid characteristics of waxy crude oil in wax deposition modeling is not wellunderstood. The next subsection discusses the steps involved in wax deposition modeling and the potential roles of nonNewtonian fluid characteristics in these steps. 1.3. Wax Deposition Modeling Considering NonNewtonian Fluid Mechanics. A wax deposition model consists of three steps: (1) hydrodynamic calculations and (2) heat and mass transfer calculations, which then allow for (3) calculation of the deposit thickness. The hydrodynamic calculation predicts the pressure drop along the pipeline. Pressure drop predictions can be used to compare with the measured pressure drop in experiments/field operations to benchmark the wax deposition model. In addition, hydrodynamic calculations predict the radial velocity profile and the eddy momentum diffusivity to be used in the subsequent heat and mass transfer calculations. Heat and mass transfer calculations simulate the heat loss from oil to the surroundings and the radial molecular diffusion of wax, followed by calculations of deposit growth and aging. These steps in the wax deposition modeling algorithm need to be enhanced to capture the non-Newtonian characteristics of waxy crude oil. Benallal et al. and Zheng et al. first attempted to model wax deposition with non-Newtonian fluid mechanics in laminar flow regime.21,22 In this investigation, we analyzed non-Newtonian turbulent characteristics to enable turbulent applications such as industrial-scale pipe flow. The following section, i.e., section 2, discusses two existing numerical techniques, including large eddy simulation (LES) and Reynolds-averaged Navier−Stokes (RANS) simulation, for: • hydrodynamic modeling to calculate the shear stress at the wall/deposit-fluid interface as well as the radial velocity profile and eddy diffusivities • heat and mass transfer modeling Predictions from the theoretically advanced LES are used as benchmarks to assess the more commonly used RANS. It was found that RANS cannot be used to model non-Newtonian turbulent flow. In section 3, we will develop a method to model hydrodynamics, heat transfer, and mass transfer by adapting the law of the wall method. The predictions from the modified law of the wall agree well with those from LES. Section 4 discusses the modeling of wax deposit formation based on first principles of the rheology of waxy oil. It should be noted that this workflow developed for non-Newtonian flow, heat and mass transfer modeling is based on single-phase oil flow. Extension of the methodology to multiphase flow regimes is beyond the scope of this investigation. Finally, the wax deposition model developed by combining findings from sections 2−4 will be applied to simulate deposit build-up in a real-world field pipeline, discussed in section 5.

the wall shear stress of both Newtonian and non-Newtonian fluids can be analytically derived as a function of the pipe diameter, the flow rate and the rheological parameters of the fluid of interest. However, transportation of crude oil in sub-sea pipelines usually occurs in turbulent flow regime. Therefore, it is imperative to establish reliable methods to calculate the wall shear stress for industrial-scale field pipeline transporting waxy crude oils. At temperatures above the WAT, a waxy crude oil behaves as a Newtonian fluid23 and the calculation of wall shear stress for Newtonian turbulent pipe flow can be readily achieved using a friction factor correlation such as the one by Blasius,24 shown in eq 1. fD = 0.316Re−1/4

(1)

where fD = Darcy friction factor Re = Reynolds number

The wall shear stress can then be calculated using the friction factor defined in eq 1: τwall/interface = fD

ρoil U 2 (2)

8

where τwall/interface = shear stress imposed by the fluid at wall/deposit‐fluid interface (Pa) ρoil = density of the oil (kg/m 3) U = superficial velocity of the oil (m/s)

Unfortunately, methods to calculate the wall shear stress of turbulent non-Newtonian pipe flows are not well-established.25 The most common approach is to calculate the shear stress at the wall on the basis of friction factors correlations. However, friction factor correlations are available only for non-Newtonian fluids whose flow curves follow simple rheological models such as the power law model26,27 and the Herschel−Bulkley model.28,29 The complexity of the rheological behavior of a waxy crude oil can be beyond these classical rheological models.11,30 Friction factor correlations for non-Newtonian fluids with complex rheological behaviors do not exist at this time. In addition, different friction factor correlations can give drastically different predictions for the wall shear stress, even for the same non-Newtonian fluid. Under extreme conditions, the prediction of the wall shear stress obtained from one correlation can be twice as high as the prediction from another correlation.25 Therefore, it is imperative to develop reliable methods to model the hydrodynamics of nonNewtonian turbulent pipe flow. The turbulent characteristics of non-Newtonian pipe flow can be modeled from first principles with direct numerical simulation (DNS). This simulation can be achieved by numerically solving the Navier−Stokes equations, shown in eq 3.

2. HYDRODYNAMIC MODELING WITH NON-NEWTONIAN CHARACTERISTICS The hydrodynamic calculations will predict the pressure drop and the radial velocity profile as well as the eddy momentum diffusivity. The pipeline pressure drop, or equivalently, the shear stress at the pipe wall/deposit-fluid interface, is a critical parameter for pipeline design. Moreover, as will be discussed in section 4, the prediction of the deposit growth rate is achieved by comparing the shear stress imposed by the fluid and the dynamic yield stress of the deposit. In the laminar flow regime,

ρoil

∂U + U ·∇U = −∇p + ∇·τ ∂t = −∇p + ∇·μ[∇U + (∇U)T ]

5012

(3)

DOI: 10.1021/acs.energyfuels.7b00504 Energy Fuels 2017, 31, 5011−5023

Article

Energy & Fuels

can be used to account for the non-Newtonian fluid mechanics in wax deposition modeling. The governing equations of LES and RANS are included in Appendix A. 2.2.1. Input Parameters for LES. Sub-sea pipelines are usually with the length scale of tens of kilometers. Along the length of a sub-sea pipeline, the temperature of the fluid can change by as much as 50 °C, leading to a drastic change in the volume fraction of solid wax of the fluid along the axial direction. Consequently, the rheological parameters of the fluid can vary significantly along the pipe from upstream to downstream locations. In order to demonstrate this unique development of the hydrodynamics along the sub-sea pipeline associated with the change in the local rheological parameters, large eddy simulations were performed at five difference axial locations, i.e., locations 1−5, with location 1 representing the inlet of the pipeline. The Newtonian fluid mechanics approach is used in this preliminary simulation in order to estimate the solid volume fraction in the bulk and the rheological parameters associated with the local solid wax volume fraction. In the complete wax deposition model presented later in the manuscript, the heat/ mass transfer simulations are integrated with non-Newtonian hydrodynamic simulations. We begin the analysis by using the Herschel−Bulkley model for the LES, as this simple rheological model accurately fits the measured shear stress-shear rate curves of the oil of interest. The rheological parameters in the Herschel−Bulkley equation, namely the dynamic yield stress, τy, the consistency, K and the flow index, n were obtained on the basis of the measured viscosity−temperature curves at various shear rates. The Herschel−Bulkley model parameters at a particular temperature were obtained by fitting the viscosity measurements at this temperature under various shear rates using the Herschel− Bulkley equation. All input parameters discussed in this section were included in Appendix B. 2.2.2. Computational Specifications. The LES was implemented with a commercial computational fluid dynamic software, FLUENT 16. The computational mesh represents a 1.5-m-long section of pipeline with an inner diameter of 0.305 m. In order to reduce the required length of the computational domain to achieve fully developed turbulent flow, a pair of periodic boundary conditions is used at the inlet and outlet of the computational domain. The entire computational mesh contains 7 000 000 computational cells. Fine computational mesh grids were used in the vicinity of the pipe wall. A transient LES with a time step of 4 ms is performed to model the evolution of the fluctuating velocity field. With this time step, it will take 2.5 h to simulate one residence time (∼2 s) using 48 CPUs in parallel. It usually takes a simulation with five residence times to reach a statistical steady state, after which five additional residence times are simulated to collect statistics. Therefore, each simulation to obtain time-averaged velocity profile and pressure drop in a 1.5-m pipe section will require 25 h. This computational time required to model a 1.5-m pipe section prohibits the direct scale-up of LES to industrial-scale wax deposition modeling where the pipeline is usually several kilometers in length. Detailed procedures to perform LES were included in Appendix C. It should be noted that the viscosity defined by the Herschel− Bulkley model approaches infinity as the shear rate approaches zero. As a result, this mathematical singularity can cause divergence in the CFD solver. In order to overcome this numerical challenge, the Herschel−Bulkley−Papanastasiou model32 shown in eq 5, is used instead of the conventional

where t = time (s) U = velocity vector (m/s) p = pressure (Pa) τ = stress tensor (Pa) μ = viscosity of the oil (Pa·s)

The viscosity of the fluid can be described by a non-Newtonian constitutive equation, such as the Herschel−Bulkley equation shown in eq 4. μ=

τy + K |γ |̇ n |γ |̇

(4)

where τy = dynamic yield stress (Pa) K = consistency (Pa·sn) |γ |̇ = strain rate magnitude (s−1) n = flow index

It should be noted that in order to resolve the motions of turbulent eddies, the size of the computational grid needs to be smaller than the size of the turbulent eddies.31 For an industrial pipeline, the Reynolds number can be on the order of 106. The size of the smallest eddies associated with this Reynolds number will be on the order of 10−3 m. Consequently, resolution of all the turbulent eddies in a 10-km industrial-scale pipeline with an inner diameter of 12 in. will require a computational mesh with as many as 1012 computational cells. Simulations with such computational intensity are not feasible. As a result, turbulence modeling approaches, such as DNS, are usually prohibitive in engineering applications, e.g., wax deposition modeling in subsea pipelines, due to the computational intensity. In comparison with DNS, large eddy simulation (LES) is more promising to be adapted to model the turbulence characteristics in industrialscale pipelines owing to its relatively more tolerable computational intensity. In LES, only the turbulent eddy motions larger than the computational cell size are resolved, thereby significantly reducing the computational intensity and cost of LES. With the existing computational power, LES is the best available approach to investigate the turbulence characteristics for non-Newtonian pipe flow at field scale. However, the computational intensity associated with LES still does not allow it to be applied on the entire pipeline with a length scale of kilometers. Therefore, we will first perform LES for sections of a field-scale pipeline (not the entire pipeline) in order to gain insights on the turbulent characteristics of non-Newtonian pipe flow. We will then use LES to evaluate the turbulent characteristics predicted by more computationally efficient approaches for hydrodynamic modeling including Reynoldsaverage Navier−Stokes (RANS) method and the law of the wall to see if they can be used. 2.1. Large Eddy Simulation. In this section, we will perform large eddy simulation (LES) to uncover the turbulent characteristics of the non-Newtonian oil flow in an industrial-scale pipeline. The findings from LES will be compared with the predictions from the more commonly used Reynolds-averaged Navier−Stokes (RANS) method to assess if the RANS method 5013

DOI: 10.1021/acs.energyfuels.7b00504 Energy Fuels 2017, 31, 5011−5023

Article

Energy & Fuels μ=

τy γ̇

[1 − exp(−mγ )] ̇ + Kγ ṅ − 1

the comparison between the viscosity calculated with LES and RANS simulation.

(5)

Herschel−Bulkley model. This model preserves the nonNewtonian characteristics by the conventional Herschel− Bulkley model but continuously approaches an upper bound mτy instead of infinity as the shear rate approaches zero. In this simulation, the constant m is set to be 104 s.32 2.2. Numerical Results from LES at Axial Locations 1−3: Turbulent Flow. It should be noted that the crude oil is Newtonian at the inlet, i.e., location 1, as the inlet temperature is higher than the WAT. Flow characteristics, such as velocity profile and pressure drop for Newtonian pipe flow have been extensively investigated and are readily available.33 Using the parameters at the inlet, we compared the pressure drop as well as the velocity profile predicted by LES and RANS and confirmed that the LES was correctly implemented. These comparisons are included in Appendix D. At locations 2 and 3, the flow is expected to be nonNewtonian, as the temperature is below the WAT at these two locations and thus wax particles will form, causing the fluid to become non-Newtonian. The LES results from locations 2 and 3 are qualitatively identical. We will now analyze the velocity profile predicted with rheological properties at location 3 as an example. Figure 1 shows an instantaneous velocity profile plotted at a cut-plane along the axial direction. It can be seen that the velocity field fluctuates spatially, indicating that a non-Newtonian turbulent flow regime is encountered at this location. As one approaches the centerline, the shear stress decreases and can potentially be lower than the local dynamic yield stress of the fluid.34 It was previously speculated that a plug could form around the centerline. Based on Figure 1, it can be seen that the instantaneous fluctuation constantly breaks up the fluid microstructure near the centerline and therefore preventing this central plug from forming. In comparison with LES, RANS is more computationally economical and commonly used for turbulent modeling.35−37 Therefore, it is of interest to evaluate whether RANS can be used to calculate the velocity and turbulent diffusivity profiles of nonNewtonian turbulent pipe flow. It should be noted that the timeaveraged strain rate is used to calculate the strain rate-dependent viscosity in the RANS governing equations because RANS only tracks time-averaged quantities. Therefore, the calculation for viscosity is different in LES and RANS simulation. Figure 2 shows

Figure 2. Comparison between the instantaneous viscosity predicted by LES and the time-averaged viscosity predicted by RANS.

As can be seen from Figure 2, the RANS simulation predicts a plug with an extremely high viscosity in the central region of the pipe. It should be noted that the prediction of a plug is a numerical artifact because the RANS approach cannot resolve the instantaneous fluctuating strain rate in the central region. Consequently, RANS can significantly under-estimate the intensity of turbulent mixing during heat and mass transfer simulation, leading to unreliable predictions for the temperature and concentration profiles. In order to demonstrate this shortcoming of RANS, RANS is used to solve the heat/mass transfer governing equation shown in eq 6, and the predicted radial temperature/concentration profile is compared with prediction from LES. Detailed implementations of RANS and LES to numerically solve eq 6 are included in Appendix E. ∂ϕ ∂ ⎛⎜ ∂ϕ ⎞⎟ ∂ + (Ujϕ) − Γ =0 ∂t ∂xj ∂xj ⎜⎝ ∂xj ⎟⎠

(6)

where ϕ = temperature or wax concentration Γ = thermal or mass diffusivity (m 2/s)

Figure 3 shows the predicted radial temperature/concentration profile by RANS and LES. It can be seen that the prediction by RANS differs significantly from that by LES. It can

Figure 1. Snapshot of the instantaneous velocity magnitude generated from turbulent non-Newtonian pipe flow simulation performed with rheological parameters at location 3. 5014

DOI: 10.1021/acs.energyfuels.7b00504 Energy Fuels 2017, 31, 5011−5023

Article

Energy & Fuels be observed that RANS predicts a large gradient near the central region of the pipe flow due to the under-estimation of the effective diffusivity in this region.

Figure 4. Comparison between the analytical and numerical velocity profiles of Herschel−Bulkley fluid with a laminar velocity profile expected for Newtonian pipe flow.

Moreover, the nonfluctuating velocity profile predicted by LES corresponds well with the analytical solution for the radial velocity profile of laminar Herschel−Bulkley fluid pipe flow. It should be noted that a plug is predicted in the center region of the pipe. Different from the artificial “plug” predicted due to numerical artifacts associated with RANS, this plug in laminar flow can be observed in laminar pipe flow experiments.38 It should be noted that being able to capture the transition from turbulent to laminar flow along the axial direction is critical in the modeling of heat and mass transfer as turbulent flow has significant lower heat/mass transfer resistance than laminar flow. In order to understand the impact of flow regime prediction on heat and mass transfer modeling, two heat transfer simulations were performed with (Case 1) a Newtonian viscosity model and (Case 2) a non-Newtonian viscosity model. Table 1 shows the parameters used in these two simulations.

Figure 3. Predicted radial profile of temperature/concentration with RANS and LES.

A dimensionless number, Δ, characterizing the relative rate of convective and diffusive transfer, defined in eq 7, is also calculated on the basis of the temperature/concentration profiles predicted by RANS and LES, respectively. Note that Δ is identical to the Nusselt number in the context of heat transfer and Sherwood number for mass transfer. ∂ϕ

− ∂r |wall D h D Δ = int = ϕbulk − ϕwall Γ

(7)

where Δ ≡ Nu, Nusselt number for heat transfer Δ ≡ Sh, Sherwood number for mass transfer

Table 1. Input Parameters for Heat Transfer Simulation with Newtonian and Non-Newtonian Approaches

h int = internal convective transfer coefficient (m/s) D = pipe diameter (m)

Case 1 viscosity model diameter, D (m) length, L (m) velocity (m/s) inlet temperature (°C) ambient temperature (°C) external heat transfer coefficient (W/m2/K) predicted internal heat transfer coefficient (W/m2/K)

r = radial coordinate (m) ϕbulk = temperature/concentration in the bulk ϕwall = temperature/concentration at the wall

The value of Δ predicted by LES is 823, whereas that predicted by RANS is 80. It can be seen that, for the non-Newtonian case with the RANS method, the predicted Δ value is significantly lower, indicating under-estimation of the heat/mass transfer in the central region due to the artificial plug. As a result, RANS cannot be used for heat/mass transfer modeling. 2.3. Numerical Results from LES Simulation at Axial Locations 4 and 5: Transition from Turbulent to Laminar Flow. Figure 4 shows the instantaneous steady-state radial velocity profile obtained with rheological properties at location 5. It can be observed that the velocity profile at axial location 5 is nonfluctuating, indicating that the flow field has become laminar at this downstream location. The significant amount of solid wax in the bulk at this location leads to a high viscosity in the bulk, dampening the turbulent eddies and causing transition from turbulent to laminar flow. In order to support this claim, the dampening of the turbulent kinetic energy of a fluid packet traveling along the centerline is included in Appendix F.

Case 2

Arrhenius

non-Newtonian 0.3048 1.524 (L = 5D) 0.75 30 5 20 460 60

The conventional wax deposition model predicts a turbulent flow regime with the Newtonian viscosity−temperature dependence calculated using a form of the Arrhenius equation. However, the actual viscosity of the wax-in-oil suspension is significantly higher than the value calculated by the Arrhenius temperature dependence, as the Arrhenius temperature dependence does not account for the effect of suspended solid on viscosity. Due to the under-estimation of the viscosity with the Arrhenius equation, conventional wax deposition models fail to predict the correct (i.e., laminar) flow regime. The internal heat transfer coefficient associated with the turbulent flow predicted from Case 1 using Newtonian fluid mechanics is 460 W/m2/K while that from a laminar flow in Case 2 is 60 W/m2/K. We note that an 5015

DOI: 10.1021/acs.energyfuels.7b00504 Energy Fuels 2017, 31, 5011−5023

Article

Energy & Fuels

This friction factor correlation is first utilized to predict the flow regime. The transition point is defined as the intersection of the wall shear stress−mean velocity relationship generated by assuming laminar and turbulent flows, respectively. The rheological parameters taken at location 5 of the pipeline were used to generate an example for the flow pattern prediction, shown in Figure 5.

approximately 700% over-estimation in the internal heat transfer coefficient is observed by neglecting the transition from turbulent to laminar discussed above. This over-estimation can translate into an over-estimation in the rate of heat transfer. As a result, it is essential to incorporate the modeling of nonNewtonian fluid mechanics in wax deposition modeling.

3. A COMPUTATIONALLY EFFICIENT METHOD FOR HYDRODYNAMIC MODELING TO BE USED IN WAX DEPOSITION MODELING As was discussed in the previous section, neither LES nor RANS can be applied to model wax deposition in a pipeline. LES is unacceptably time-consuming, and RANS can generate numerical artifacts. Therefore, it is imperative to develop a reliable and computationally efficient method to predict the wall shear stress, velocity, temperature and concentration profiles to be used to predict wax deposition. This desired method should be able to capture the flow characteristics uncovered by LES simulations, recapitulated below: • a flow pattern that can transit from turbulent to laminar flow due to the increase in viscosity • a plug in the near center line region that can be expected in the laminar flow regime but is likely to be broken up by turbulent eddies in the turbulent flow regime • a laminar boundary layer that can be observed in the vicinity of the wall in the turbulent flow regime, while the local velocity fluctuates with time outside this boundary layer We now propose a reliable and computationally efficient method to perform hydrodynamic, heat transfer, and mass transfer modeling of non-Newtonian pipe flow. This method consists of three steps: • prediction of flow regime • calculation of velocity profile • calculation of turbulent diffusivity for heat and mass transfer modeling These three steps will be discussed individually in the three subsections to follow. 3.1. Prediction of Flow Regime. Chilton and Stainsby developed a friction factor correlation that can predict the wall shear stress based on the fluid rheological parameters, the pipeline diameter and the flow rate,29 shown in eq 8. ⎤−0.25 ⎡ R fHB = 0.079⎢ 2 HB 4 ⎥ ⎣ n (1 − X ) ⎦

Figure 5. Prediction of wall shear stress based on the laminar and turbulent Chilton−Stainsby friction factor correlations and the definition of the onset for turbulent−laminar transition.

A sensitivity analysis was also performed to understand the effect of rheological parameters, τy, n, and K, on the prediction of laminar−turbulent transition. It was found that the transition occurs at a higher mean velocity with larger values of τy and K and with smaller values of n. Details with respect to this sensitivity analysis are included in Appendix G. Table 2 shows the flow pattern predictions (laminar or turbulent) based on the Chilton−Stainsby friction factor correlation as well as predictions by LES. Table 2. Predictions of Flow Regime with LES and Chilton− Stainsby (C-S) Friction Factor Correlation flow pattern

(8)

location

LES

C−S

axial location (% total length)

1 2 3 4 5

turbulent turbulent turbulent laminar laminar

turbulent turbulent turbulent laminar laminar

0 16 21 38 55

where RHB =

X=

μwall

τy τwall

c=

3.2. Calculation of Velocity Profile. For the laminar flow regimes, an analytical solution for the velocity profile and pressure drop can be obtained. Therefore, we will now focus on developing a computationally efficient approach for the hydrodynamic modeling of turbulent flow. In this subsection, we will explore modifications of the conventional law of the wall that are necessary in order for it to be applicable to non-Newtonian pipe flow. The conventional law of the wall can be expressed with eqs 9 and 10.

ρoil UD

,

( 3n4+n 1 )( 1 − aX −1bX − cX ) 2

a=

1 , 2n + 1

3

b=

2n , (n + 1)(2n + 1)

2n2 (n + 1)(2n + 1)

⎧ y+ y+ ≤ 5 ⎪ ⎪ U z+ = ⎨ 5 ln y+ − 3.05 5 < y+ ≤ 30 ⎪ ⎪ 2.5 ln y+ + 5.5 y+ ≥ 30 ⎩

fHB = friction factor for a Herschel−Bulkley fluid RHB = Reynolds number for a Herschel−Bulkley fluid μwall = viscosity of oil at wall (Pa·s) 5016

(9)

DOI: 10.1021/acs.energyfuels.7b00504 Energy Fuels 2017, 31, 5011−5023

Article

Energy & Fuels U z+ =

Uz τwall /ρoil

,

y+ =

y ν

τwall ρoil

shear stress or equivalently the pressure gradient needs to be predicted. As is shown in Appendix H, it was verified that reliable wall shear stress/pressure gradient predictions for non-Newtonian pipe flow can be achieved with the Chilton−Stainsby friction factor correlation. Figure 7 shows the dimensionless radial velocity profile predicted by the modified law of the wall and LES. As can be seen, the dimensionless radial velocity profile predicted by LES corresponds well with the law of the wall in the laminar region as well as the buffer region but deviates in the outer region. This small deviation is potentially due to the fact that viscosity increases as the distance to the wall distance increases.39

(10)

where Uz = axial velocity magnitude (m/s) y = distance from the wall (m) ν = kinematic viscosity of oil (m 2/s)

The conventional law of the wall is derived on the basis of the ideas that • A thin laminar sub-layer exists in the immediate vicinity of the wall within which velocity changes drastically. In the laminar layer, the shear stress is mainly contributed by viscous stresse. • Outside this laminar layer, velocity fluctuation is significant, and the shear stress is mainly contributed by Reynolds stress. In order to check these two assumptions for non-Newtonian turbulent flow, a representative instantaneous radial velocity profile is generated from LES performed with non-Newtonian turbulent flow, shown in Figure 6. It can be seen that these two qualitative turbulence characteristics are preserved in the nonNewtonian turbulent pipe flow, with the laminar-like boundary sub-layer in the vicinity of the pipe wall and a turbulent core in the central region of the flow.

Figure 7. Comparison between the dimensionless velocity profiles predicted by LES and by the law of the wall.

Figure 8 shows the predicted velocity profiles by LES and the law of the wall. It can be seen that the radial velocity profile

Figure 6. Representative instantaneous radial velocity profiles for a turbulent non-Newtonian pipe flow.

Consequently, we can analogously derive a modified law of the wall for non-Newtonian turbulent flow. It should be noted that in the laminar sub-layer, the shear stress equals viscous stress. Therefore, velocity profile of non-Newtonian flow obeys the conventional law of the wall in the laminar region if the viscosity at wall (r = R) is used in place of the Newtonian viscosity in the definition for the dimensionless distance-to-wall, shown in eq 11. y+ =

y νwall

τwall ρoil

Figure 8. Radial velocity profiles predicted by LES and the law of the wall.

predicted by the law of the wall only differs from LES predictions by ∼10%, showing that the law of the wall method is reliable to generate predictions for the velocity profile to be used in heat and mass transfer modeling. In this section, we developed a rapid method to perform hydrodynamic modeling in order to generate input parameters for the heat and mass transfer modeling. The key findings from this section are summarized below: • The flow regime of non-Newtonian pipe flow can be readily predicted by the Chilton−Stainsby correlation.

(11)

where νwall = kinematic viscosity at wall (m 2/s)

In the conventional law of the wall with this modified definition of the dimensionless distance-to-wall, y+ will be implemented to predict the velocity profile of non-Newtonian turbulent pipe flow. In order to initiate this prediction, the wall 5017

DOI: 10.1021/acs.energyfuels.7b00504 Energy Fuels 2017, 31, 5011−5023

Article

Energy & Fuels • The modified law of the wall can generate reliable predictions for the radial velocity profile. 3.3. Calculation of Turbulent Diffusivities for Heat and Mass Transfer Modeling. Solving the heat/mass transfer equation is a critical step in wax deposition modeling. In order to generate reliable predictions for the temperature and concentration profiles, the turbulent heat/mass diffusivity needs to be calculated accurately. As was shown in section 2.2, RANS fails to predict the correct temperature and concentration profiles as it significantly under-predicts the turbulent diffusivities. The turbulent heat/mass diffusivity depends on the turbulent momentum diffusivity, shown in eqs 12 and 13, as the mixing of heat and mass is caused by turbulent eddy motions. εM Sc νT = Dwo Sc T νwall (12)

εH Pr νT = α PrT νwall

Figure 9. Predicted radial profile of temperature/concentration with the law of the wall and LES. ∂ϕ

(13)

− ∂r D h D wall Δ = int = ϕbulk − ϕwall Γ

where εM = eddy mass diffusivity (m 2/s)

For heat transfer, Δ is the Nusselt number, and for mass transfer, Δ is the Sherwood number. The values of Δ calculated on the basis of LES and the law of the wall, respectively 823 and 1116, agree well with each other. In summary, the modified law of the wall discussed in this section is appropriate for hydrodynamic, heat transfer, and mass transfer modeling for non-Newtonian turbulent pipe flow. The algorithm based on the modified law of the wall to be incorporated in the wax deposition model is summarized in Figure 10. It should be noted that Figure 10 includes another advancement in the wax deposition modeling, i.e., accounting for the effect of shear stress at wall on the calculation of deposition rate, which will be discussed in detail in the following section.

Dwo = molecular diffusivity of wax in oil (m 2/s) Sc = Schmidt number Sc T = turbulent Schmidt number εH = eddy thermal diffusivity (m 2/s) α = material thermal diffusivity of oil (m 2/s) Pr = Prandtl number PrT = turbulent Prandtl number νT = turbulent diffusivity (m 2/s)

The modified law of the wall is again used to calculate the turbulent momentum diffusivity, shown in eqs 14

νT νwall

⎧ + ⎤2 ⎡ ⎪(κy+ )2 ⎢1 − exp⎛⎜ − y ⎞⎟⎥ y+ ≤ 5 ⎪ ⎥ ⎢⎣ A ⎠⎦ ⎝ ⎪ 2 ⎪ ⎡ ⎛ y + ⎞⎤ 5 ⎪ = ⎨(κy+ )2 ⎢1 − exp⎜ − ⎟⎥ + 5 < y+ ≤ 30 ⎢⎣ ⎪ ⎝ A ⎠⎥⎦ y ⎪ 2 ⎪ ⎡ ⎛ y+ ⎞⎤ 2.5 + + 2 ⎪(κy ) ⎢1 − exp⎜ − ⎟⎥ + y ≥ 30 ⎪ ⎢⎣ ⎝ A ⎠⎥⎦ y ⎩ (14)

with κ = 0.4, A = 26. Figure 9 shows the predict temperature/concentration profile using the modified law of the wall. As can be seen, the temperature/concentration profile predicted by law of the wall is virtually identically with that predicted by LES, indicating that the law of the wall is appropriate for the calculation of turbulent diffusivity to be used in heat/mass transfer modeling. In order to compare results from the law of the wall and LES to further validate the law of the wall method, the important dimensionless number, Δ, characterizing the relative rate of convective and diffusive transfer, defined in eq 7 and repeated below, is again calculated on the basis of the temperature/ concentration profiles predicted by the law of the wall and LES.

Figure 10. Summary of the hydrodynamic, heat transfer, and mass transfer modeling algorithm based on the law of the wall. 5018

DOI: 10.1021/acs.energyfuels.7b00504 Energy Fuels 2017, 31, 5011−5023

Article

Energy & Fuels

4. MODELING OF GELATION IN THE IMMEDIATE VICINITY OF THE WALL/DEPOSIT-FLUID INTERFACE Conventional wax deposition models simulate deposit growth solely on the basis of the molecular diffusion mechanism. Wax deposition is initiated by the precipitation of wax molecules at the wall (t = 0+)/deposit-fluid interface (t > 0) when the surface temperature is below the wax appearance temperature (WAT). Precipitation of wax generates a radial concentration gradient and a net diffusive flux of dissolved wax molecules to the cold solid surface. A fraction of the dissolved wax molecules reaching the deposit-fluid interface precipitate and form an interlocking network with entrapped oil, leading to deposit growth.1 The rest of the wax molecules continue to diffuse into the existing deposit layer, causing an increase in the wax content of the deposit.1 Modeling efforts have been mostly concentrated on the calculation of radial diffusion of wax molecules toward the pipe wall (t = 0+)/deposit-fluid interface (t > 0) and the internal diffusion of wax molecules into the deposit, as these variables are essential for the calculation of the growth rate of wax deposit, shown in eq 15. ρwax Ω wax

dδ = Jwax,to interface − Jwax,into the deposit dt

Figure 11. Comparison between the deposit solid wax contents from two runs with similar concentration driving forces but drastically different wall shear stresses.

compositions of the solid deposits in these two experiments are drastically different. A deposit layer forms when its dynamic yield stress exceeds the shear stress imposed by the flow. As the dynamic yield stress of the deposit increases with the solid content of the deposit, the deposit formed at the higher flow rate (high wall shear stress) contains a significantly larger amount of solid wax than the deposit formed at the lower flow rate. Both rheometric and flow loop experiments suggest that the deposit solid content can vary significantly with shear stress. Despite the overwhelming experimental evidence suggesting that the deposit solid content changes with imposed shear stress, it has remained unclear how the solid volume fraction the depositing layer can be modeled rigorously. We now introduce a rigorous method to model the formation of a gel deposit and predict the solid volume fraction of the depositing layer on the basis of first principles of rheology. The solid wax-in-oil suspension in the immediate vicinity of the wall gradually develops a yield stress (τy) as the solid volume fraction, Ωwax, increases due to molecular diffusion of wax toward the wall and into the deposit. The wax-in-oil suspension in the immediate vicinity of the wall cease to flow when the yield stress, τy of the suspension exceeds the shear stress (τs) imposed by the fluid flow. Figure 12 provides an illustration of the process of gelation described above associated with the mathematical implementations.

(15)

where ρwax = density of the wax (kg/m 3) Ω wax = volume fraction of solid wax in the gel δ = deposit thickness ( m) t = time (s) Jwax to interface = wax diffusive flux to the interface (kg/m 2/s) Jwax into the deposit = wax diffusive flux into the deposit (kg/m 2/s)

In existing wax deposition models, the solid volume fraction of the depositing layer, Ωwax, is either assumed to be a constant or determined empirically. As a result, existing models cannot capture the variation of deposit solid wax volume fraction with changing operating conditions, such as the imposed shear stress by the fluid and fail to predict the correct deposit growth rate and solid fraction at high and low oil flow rates.40 Multiple experimental evidence suggests that the deposit solid volume fraction increases with increasing shear stress imposed by the fluid. Singh et al. measured the gelation temperature of a model wax-in-oil mixture under various shear stresses and cooling rates with a controlled stress rheometer.41 In Singh et al.’s study, the gelation temperature is defined as the point when the loss and storage modulus equal each other. It was discovered that the gelation temperature of the oil of interest decreases with increasing imposed shear stress, suggesting that a larger amount of solid has precipitated at gelation under a higher shear stress. In addition to rheometric characterizations, the variation of the solid content in the deposit was also observed from flow loop wax deposition experiments. 42 Figure 11 shows the deposit composition obtained from two flow loop wax deposition experiments performed with the same crude oil and with similar diffusive fluxes of dissolved wax, but with different oil and wall temperatures. As can be seen from Figure 11, despite the similar mass flux of dissolved wax between the two experimental runs, the

Figure 12. Illustration of gelation process in the immediate vicinity of the wall and the mathematical implementation. 5019

DOI: 10.1021/acs.energyfuels.7b00504 Energy Fuels 2017, 31, 5011−5023

Article

Energy & Fuels

Figure 13. (a) Evolution of solid volume fraction in the immediate vicinity of the wall. (b) Evolution of dynamic yield stress of wax-in-oil suspension in the immediate vicinity of the wall to the point of deposit formation.

Figure 14. Comparison between the predicted pressure drop evolutions and that observed in field.

5. APPLICATION OF THE ENHANCED WAX DEPOSITION MODEL ON A FIELD-SCALE PIPELINE The enhanced wax deposition model is applied to predict wax deposition in a real-world field-scale pipeline off the shore of Indonesia.9 This pipeline transports waxy crude oil from a central processing platform to floating production storage offloading. The sub-sea pipeline has an inner diameter of 12 in. and a length of 23 km. The pressure drop across the entire pipeline is monitored over the production period. It was observed that wax deposition causes the pressure drop to increase from ∼200 to ∼300 psi over a production period of 7 days.9 Detailed description of this field as well as the monitored operating variables, such as flow rate, temperature, and pressure drop can be found in the article by Singh et al.9 The simulation parameters and configurations are summarized in Appendix I. The wax deposition model is first used to calculate the pressure drop along the pipeline without wax deposit attached to the pipe wall. The calculated pressure drop, 180 psi, corresponds well with the pressure drop recorded right after pigging, ∼200 psi, thereby validating the hydrodynamic calculation. The heat transfer simulation predicts an average outlet temperature of 30 °C, which corresponds well with the field observation of an outlet temperature varying from 27 to 29 °C during production.9 Wax deposition modeling was then performed after the hydrodynamic and heat transfer calculations were validated by comparison between predictions and field observations. Figure 14a shows the predicted pressure drop−time trajectory by assuming two

Figure 13a shows the predicted increase of the solid volume fraction in the oil in the immediate vicinity of the wall in a fieldscale simulation. In order to calculate the change of dynamic yield stress as a function of time, the solid volume fraction at each time step shown in Figure 13a is entered into the relationship between the dynamic yield stress and the solid volume fraction. The relationship between the dynamic yield stress and the solid volume fraction was obtained by fitting of the viscosity− temperature curves measured at different shear rates using the Herschel−Bulkley equation. The increase in the dynamic yield stress of the wax−oil suspension is shown in Figure 13b. This field-scale simulation will be described in detail in the following Section 5. As can be seen from Figure 13a, the solid volume fraction of the wax-in-oil suspension in the immediate vicinity of the wall continuously increases as time elapses due to accumulation of solid wax, leading to a continuous increase in the dynamic yield stress of the solid−liquid suspension, shown in Figure 13b. This layer gels to form a deposit when the dynamic yield stress reaches or exceeds the shear stress imposed by the fluid at the interface of ∼3.5 Pa. At this point, the modeling of hydrodynamics, heat transfer, and mass transfer as well as the deposit formation was enhanced with non-Newtonian fluid mechanics. In the following section, the enhanced wax deposition model will be applied on a realworld field-scale pipeline. 5020

DOI: 10.1021/acs.energyfuels.7b00504 Energy Fuels 2017, 31, 5011−5023

Article

Energy & Fuels

shown in Table 3. Figure 15 shows the predicted average wax deposit thickness over the entire length of the pipeline as well as the axial deposit thickness profile. It should be noted that predictions such as those shown in Figure 15 are virtually impossible to measure in a field-scale pipeline. Therefore, such insights on the average deposit thickness evolution and the axial deposit thickness profile are extremely valuable for pigging design. In order to highlight the difference between the Newtonian and non-Newtonian wax deposition models, the change of pressure drop as a function of time as well as the axial deposit thickness profile were generated with the Newtonian wax deposition model, shown in Figure 16. It should be noted that exact same input parameters, except the viscosity, used for wax deposition modeling with the nonNewtonian model are also used in the Newtonian model. An Arrhenius temperature dependency was used for the viscosity in the Newtonian model without the consideration for the increase in viscosity due to suspended solid. As can be seen from Figure 16a, the Newtonian model significantly under-estimates the pressure drop increase over time due to the under-estimation in the viscosity. Misinterpretation of this under-estimation of the pressure drop build-up generated from the Newtonian model can lead to overly optimistic assessment of the wax deposition risk. As can be seen from Figure 16b, the Newtonian model also significantly under-predicts the deposit thickness. This under-estimation is due to the fact that the Newtonian model assumes the solid fraction of the deposit to be always greater than total wax content of the oil, i.e., 17 wt% for this particular crude, while based on analysis from rheology, it takes only ∼6 wt% of solid wax to immobilize the wax-in-oil suspension and form a deposit.

Table 3. Summary of Non-Newtonian Model Predictions and Field Observations

pressure drop right after pigging (psi) pressure drop with Newtonian approach (psi) pipeline outlet temperature (°C) deposit wax content (wt%)

prediction

field observation

180 160 30 24

200 200 27−29 27 ± 2

limiting scenarios for the bulk precipitation kinetic constant, i.e., the Chilton−Colburn approach9 and the solubility approach.43 As can be seen from Figure 14a, the Chilton−Colburn analogy significantly over-predicts the pressure drop build-up rate while the solubility method is more representative for the bulk precipitation kinetics in this field-scale simulation. The bulk precipitation kinetic constant is adjusted to achieve agreement between the predicted pressure drop−time trajectory and that observed in the field, shown in Figure 14b. As can be seen from Figure 14b, the predicted pressure drop matches the observed pressure drop evolution with an adjusted bulk precipitation kinetic constant of 2 × 10−2 s−1. This value is consistent with the range of bulk precipitation kinetic constant reported by Lee et al., i.e., 10−2−1 s−1.44 It is noticed that the predicted starting pressure drop is slightly lower (10%) than the observed pressure drop. This discrepancy is likely due to the fact that a certain amount of wax deposit remains attached to the pipeline after the pigging, therefore leading to a higher pressure drop right after the pigging than the pressure drop associated with a completely bare pipeline. This enhanced wax deposition model is the first model that can generate reasonable field predictions from first principles of fluid mechanics, transport and rheology. The enhanced wax deposition model predicts a wax content of 24 wt%, which corresponds well with the wax content of the pig return deposit, i.e., 27 ± 2 wt%. Table 3 summarizes the comparisons between predictions and field observations. As can be seen from Table 3, excellent agreement between the model predictions and field observations of pressure drop, temperature as well as the wax content presented in this section indicate that the enhanced wax deposition model can provide wax deposit predictions with a high level of confidence. If the Newtonian approach is used to predict the pressure drop right after pigging, a more severe under-prediction can be observed, as

6. CONCLUSIONS In this investigation, we assessed the role of non-Newtonian characteristics on wax deposition modeling. Non-Newtonian hydrodynamics as well as heat and mass transfer modeling were investigated with large eddy simulations (LES). With LES, the effects of non-Newtonian turbulent characteristics on hydrodynamics as well as heat and mass transfer models were studied from first principles. The computational cost of applying LES to model wax deposition in a field-scale pipeline was calculated. It was discovered that LES is too computationally intensive for the

Figure 15. (a) Predicted average deposit thickness evolution over a production period of 7 days. (b) Predicted axial deposit thickness profile after 7 days of production. 5021

DOI: 10.1021/acs.energyfuels.7b00504 Energy Fuels 2017, 31, 5011−5023

Article

Energy & Fuels

Figure 16. (a) Comparison between the pressure drop−time trajectories predicted by non-Newtonian and Newtonian approaches. (b) Comparison between the axial deposit thickness profiles predicted by non-Newtonian and Newtonian approaches.



modeling of wax deposition in a field-scale pipeline. The option of using common turbulent modeling technique, i.e., Reynoldsaveraged Navier−Stokes (RANS), for wax deposition modeling in industrial-scale pipeline was explored. Comparison between LES and RANS predictions revealed that commonly used RANS models lead to numerical artifacts in the hydrodynamic and transport modeling results. Based on the insights provided by LES, computationally efficient and reliable methodologies, i.e., the modified law of the wall, for the calculations of velocity profile, temperature, and concentration profiles were developed to enhance hydrodynamic, heat transfer, and mass transfer calculations in existing wax deposition models. The key findings from the hydrodynamic modeling as well as heat and mass transfer modeling using the three simulation techniques, i.e., LES, RANS, and the law of the wall, are summarized as follows: Key point 1: The previously suspected plug in the central region of the flow is constantly destroyed by turbulent eddies Key point 2: Transition of turbulent flow to laminar flow can occur as the fluid flows and cools along pipeline. Key point 3: The computational cost of using LES for wax deposition modeling in an industrial-scale pipeline is intolerable. Key point 4: RANS generates numerical artifacts in both velocity profile prediction and temperature/concentration profile predictions due to its limitation that it cannot resolve the instantaneous fluctuating shear rates near the central region of the flow. Key point 5: The modified law of the wall approach can generate rapid and accurate predictions for the velocity profile, temperature and concentration profiles. Therefore, the modified law of the wall approach is used for wax deposition modeling. The enhanced wax deposition model is applied to predict wax deposition in a real-world field-scale pipeline. The wax deposition model generates reliable predictions for the increase of pressure drop due to deposit build-up over time, the temperature at the outlet of the pipe as well as the wax content of the deposit. The good agreement between the model predictions and field observations validates that the predictions made by this enhanced model are of high confidence and reliability.

ASSOCIATED CONTENT

* Supporting Information S

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.energyfuels.7b00504. Appendix A, governing equations for large eddy simulations and Reynolds-averaged navier stokes simulations; Appendix B, input parameters for LES and RANS; Appendix C, steps to perform large eddy simulations; Appendix D, validation of the implementation of LES; Appendix E, governing equations for heat and mass transfer modeling using LES and RANS; Appendix F, time evolution of the turbulent kinetic energy leading to turbulent-to-laminar transition; Appendix G, effects of rheological parameters on flow regime transition; Appendix H, validation of the pressure gradient prediction from the Chilton−Stainsby correlation; and Appendix I, input parameters for field-scale wax deposition simulation (PDF)



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Sheng Zheng: 0000-0001-5948-7325 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS The authors would like to thank Total for the permission to publish this work. The authors acknowledge financial support from the University of Michigan Industrial Affiliates Program sponsored by Assured Flow Solutions LLC, Chevron, ConocoPhillips, Multichem a Hallibruton Service, Phillips66, Statoil, and Total. The authors would also like to thank Prof. Ronald G. Larson for the helpful discussion during the preparation of this manuscript.



REFERENCES

(1) Singh, P.; Venkatesan, R.; Fogler, H. S.; Nagarajan, N. Formation and aging of incipient thin film wax-oil gels. AIChE J. 2000, 46, 1059− 1074.

5022

DOI: 10.1021/acs.energyfuels.7b00504 Energy Fuels 2017, 31, 5011−5023

Article

Energy & Fuels (2) Lashkarbolooki, M.; Seyfaee, A.; Esmaeilzadeh, F.; Mowla, D. Energy Fuels 2010, 24 (2), 1234. (3) Venkatesan, R.. The Deposition and Rheology of Organic Gels. Ph.D. Dissertation, University of Michigan, 2004. (4) Seyfaee, A.; Lashkarbolooki, M.; Esmaeilzadeh, F.; Mowla, D. J. Dispersion Sci. Technol. 2011, 32 (3), 312. (5) Huang, Z.; Lee, H. S.; Senra, M.; Fogler, H. S. AIChE J. 2011, 57 (11), 2955. (6) Huang, Z.; Lu, Y.; Hoffmann, R.; Amundsen, L.; Fogler, H. S. Energy Fuels 2011, 25 (11), 5180. (7) Lu, Y.; Huang, Z.; Hoffmann, R.; Amundsen, L.; Fogler, H. S. Energy Fuels 2012, 26 (7), 4091. (8) Panacharoensawad, E. Wax Deposition under Two-Phase OilWater Flowing Conditions. Ph.D. Dissertation, The University of Tulsa, 2012. (9) Singh, A.; Lee, H.; Singh, P.; Sarica, C. Flow Assurance: Validation of Wax Deposition Models Using Field Data from a Subsea Pipeline. Proceedings of Offshore Technology Conference, Houston, TX, 2011; OTC 21641SS. (10) Noville, I.; Naveira, L. Comparison between Real Field Data and the Results of Wax Deposition Simulation. Proceedings of SPE Latin America and Caribbean Petroleum Engineering Conference, Mexico City, Mexico, April 16−18, 2012; SPE152575. (11) Dimitriou, C. J.; McKinley, G. H. A Comprehensive Constitutive Law for Waxy Crude Oil: A Thixotropic Yield Stress Fluid. Soft Matter 2014, 10 (35), 6619−6644. (12) Dimitriou, C. J.; McKinley, G. H.; Venkatesan, R. Rheo-PIV analysis of the yielding and flow of model waxy crude oils. Energy Fuels 2011, 25 (7), 3040−3052. (13) Kané, M.; Djabourov, M.; Volle, J.-L. Rheology and structure of waxy crude oils in quiescent and under shearing conditions. Fuel 2004, 83 (11−12), 1591−1605. (14) Pedersen, K. S.; Rønningsen, H. P. Effect of Precipitated Wax on Viscosities: A Model for Predicting Non-Newtonian Viscosity of Crude Oils. Energy Fuels 2000, 14 (7), 43−51. (15) Marchesini, F. H.; Alicke, A. A.; De Souza Mendes, P. R.; Ziglio, C. M. Rheological Characterization of Waxy Crude Oils: Sample Preparation. Energy Fuels 2012, 26 (5), 2566−2577. (16) Cross, M. M. Rheology of non-Newtonian fluids: A new flow equation for pseudoplastic systems. J. Colloid Sci. 1965, 20, 417. (17) Casson, N. In Rheology of disperse systems, Proceedings of a Conference Organized by the British Society of Rheology, University College of Swansea, September 1959; Mill, C. C., Ed.; Pergamon Press: New York, 1957 (18) Herschel, W. H.; Bulkley, R. Konsistenzmessungen von GummiBenzollosungen. Colloid Polym. Sci. 1926, 39, 291−300. (19) Kumar, L.; Paso, K.; Sjöblom, J. Numerical study of flow restart in the pipeline filled with weakly compressible waxy crude oil in nonisothermal condition. J. Non-Newtonian Fluid Mech. 2015, 223, 9−19. (20) Kumar, L.; Zhao, Y.; Paso, K.; Grimes, B.; Sjöblom, J.; Lawrence, C. Numerical Study of Pipeline Restart of Weakly Compressible Irreversibly thixotropic Waxy Crude Oils. AIChE J. 2015, 61 (8), 2657− 2671. (21) Benallal, A.; Maurel, P.; Agassant, J.-F.; Darbouret, M.; Avril, G.; Peuriere, E. Wax Deposition in Pipelines: Flow-Loop Experiments and Investigations on a Novel Approach. Proceedings of SPE Annual Technical Conference and Exhibition, Denver, CO, Sept 21−24, 2008. (22) Zheng, S.; Saidoun, M.; Mateen, K.; Palermo, T.; Ren, Y.; Fogler, H. S. Wax Deposition Modeling with Considerations of NonNewtonian Fluid. Proceedings of Offshore Technology Conference, Houston, TX, 2016. (23) Li, H.; Zhang, J. Fuel 2003, 82 (11), 1387. (24) Blasius, P. R. H. Das Aehnlichkeitsgesetz bei Reibungsvorgangen in Flussigkeiten. Forschungsheft 1913, 131, 1−41. (25) Heywood, N. I.; Cheng, D. C.-H. Comparison of methods for predicting head loss in turbulent pipe flow of non-Newtonian fluids. Trans. Inst. Meas. Control 1984, 6, 33−45. (26) Dodge, D. W.; Metzner, a. B. Turbulent flow of non-newtonian systems. AIChE J. 1959, 5, 189−204.

(27) Kemblowski, Z.; Kolodziejski, J. Flow resistances of nonNewtonian fluids in transitional and turbulent flow. Int. Chem. Eng. 1973, 13 (2), 265−279. (28) Hanks, R. W., Low Reynolds Number Turbulent Pipeline Flow of Pseudohomogeneous Slurries. Proceedings of Hidrotransport 5, 1978; Paper C2, pp C23−C34, (). (29) Chilton, R. A.; Stainsby, R. Pressure Loss Equation for Laminar and Turbulent non-Newtonian Pipe Flow. J. Hydraul. Eng. 1998, 124, 522−529. (30) Palermo, T.; Tournis, E. Viscosity Prediction of Waxy Oils: Suspension of Fractal Aggregates (SoFA) Model. Ind. Eng. Chem. Res. 2015, 54, 4526−4534. (31) Pope, S. B. Turbulent Flows; Cambridge University Press: Cambridge, UK, 2000. (32) Papanastasiou, T. C. Flow of materials with yield. J. Rheol. 1987, 31, 385−404. (33) Wilkes, J. O. Fluid Mechanics for Chemical Engineers, 2nd ed.; Prentice Hall: Upper Saddle River, NJ, 2005. (34) Slatter, P. T. The laminar-turbulent transition in large pipes. Problems in Fluid Mechanics and Hydrology, Proceedings of the 12th International Conference on Transport and Sedimentation of Solid Particles, Prague, June 23−26, 1999; pp 247−256. (35) Mohammadi, B.; Pironneau, O. Analysis of the k-epsilon turbulence model; Wiley: Chichester, 1993. (36) Wu, B.; Chen, S. CFD simulation of non-Newtonian fluid flow in anaerobic digesters. Biotechnol. Bioeng. 2008, 99, 700−711. (37) Wu, B. CFD simulation of gas and non-Newtonian fluid twophase flow in anaerobic digesters. Water Res. 2010, 44, 3861−3874. (38) Poumaere, A.; Moyers-González, M.; Castelain, C.; Burghelea, T. Unsteady laminar flows of a Carbopol® gel in the presence of wall slip. J. Non-Newtonian Fluid Mech. 2014, 205, 28−40. (39) Rudman, M.; Blackburn, H. M. Direct numerical simulation of turbulent non-Newtonian flow using a spectral element method 2006, 30, 1229−1248. (40) Huang, Z. Application of The Fundamentals of Heat and Mass Transfer to The Investigation of Wax Deposition in Subsea Pipelines. Ph.D. Dissertation, University of Michigan, 2011. (41) Singh, P.; Fogler, H. S.; Nagarajan, N. Prediction of the wax content of the incipient wax-oil gel in a pipeline: An application of the controlled-stress rheometer. J. Rheol. 1999, 43 (6), 1437. (42) Hoffmann, R.; Amundsen, L. Single-Phase Wax Deposition Experiments. Energy Fuels 2010, 24, 1069−1080. (43) Venkatesan, R.; Fogler, H. S. Comments on Analogies for Correlated Heat and Mass Transfer in Turbulent Flow. AIChE J. 2004, 50, 1623−1626. (44) Lee, H. S. Computational and Rheological Study of Wax Deposition. Ph.D. Dissertation, University of Michigan, 2008.

5023

DOI: 10.1021/acs.energyfuels.7b00504 Energy Fuels 2017, 31, 5011−5023