Can a Cubic Equation of State Model BitumenâSolvent Phase Behavior?

Subscriber access provided by UNIVERSITY OF CONNECTICUT

Article

Can a Cubic Equation of State Model Bitumen-Solvent Phase Behavior? Kim Johnston, Marco Satyro, Shawn David Taylor, and Harvey William Yarranton Energy Fuels, Just Accepted Manuscript • DOI: 10.1021/acs.energyfuels.7b01104 • Publication Date (Web): 03 Jul 2017 Downloaded from http://pubs.acs.org on July 9, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

Energy & Fuels is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

Can a Cubic Equation of State Model Bitumen-Solvent Phase Behavior? K. A. Johnston1, M.A. Satyro2, S.D. Taylor3, H.W. Yarranton1 1. Department of Chemical and Petroleum Engineering, University of Calgary, Calgary, AB, Canada; 2. Department of Chemical and Biomolecular Engineering, University of Clarkson, Potsdam, NY, USA; 3. Schlumberger-Doll Research, 1 Hampshire St., Cambridge, Massachusetts, USA, 02139

Abstract Cubic equations of state (CEoS), such as the Advanced Peng Robinson (APR) EoS, are convenient for use in commercial simulators and have successfully fit saturation pressures and asphaltene onset points for bitumen-solvent systems using simple quadratic mixing rules. However, this approach does not accurately predict asphaltene precipitation yields. In this study, the APR EoS with several sets of asymmetric mixing rules is evaluated against saturation pressure and asphaltene yield data for n-pentane diluted bitumen. The asymmetric van der Waals, Sandoval et al., and two forms of Huron-Vidal mixing rules with an NRTL activity coefficient model are considered. The use of asymmetric mixing rules significantly improves the match to asphaltene yield data; however, the yields are still under-predicted at high solvent contents and the tuning parameters that give the best match for asphaltene yield data are not predictive or easily correlated for other solvents. The APR EoS with symmetric van der Waals mixing rules is also evaluated with compositionally dependent binary interaction parameters. The use of compositionally dependent solvent/asphaltene binary interaction parameters allows the model to fit asphaltene yield data over the entire composition range. A set of interaction parameters is recommended that fits both asphaltene yield and saturation pressure data. The merits of this methodology as a practical option for modeling heavy oil-solvent behavior are discussed.

1 ACS Paragon Plus Environment

Energy & Fuels

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 37

1. Introduction The phase behavior of mixtures of bitumens (or heavy oils) and n-alkanes is a key factor in the design and optimization of solvent-assisted in situ heavy oil recovery processes, solvent based oil sand extraction processes, and solvent deasphalting processes. These mixtures are asymmetric and can form multiple phases including vapor-liquid (VL), liquid-liquid (LL), VLL, and possibly VLLL regions

(1-5)

. Depending on the solvent and the conditions, the mixture can split into a

bitumen-rich phase and a light solvent-rich phase or a solvent-rich phase and a heavy asphaltenerich phase. A model that captures the full range of this phase behavior, including phase compositions, has proven elusive.

To date, the most successful models for bitumen-solvent phase behavior treat all of the oil, including the asphaltenes, as a solution and treat all phase behavior, including asphaltene precipitation, as a chemical phase separation. These models include the modified regular solution (MRS) model, PC-SAFT EoS, and the cubic plus association (CPA) EoS. The MRS model can fit and, to some extent, predict asphaltene yield and onset points for various solvent-diluted bitumen systems

(6-8)

. However, this version of regular solution theory is unable to accurately

model component partitioning between an asphaltene-rich phase and a solvent-rich phase. In addition, it has not been applied to other phase transitions including VLE.

The CPA EoS can capture vapor-liquid equilibrium as well as asphaltene precipitation onsets in conventional oils and bitumen over a range of temperatures, pressures, and compositions

(9,10)

.

The PC-SAFT equation of state can also model both vapor-liquid equilibrium and the onset of asphaltene precipitation from a depressurized conventional crude oil

(11,12)

. More recently, the

PC-SAFT has been used to model asphaltene yield curves from solvent diluted bitumen and has been shown to capture the effects of asphaltene polydiversity independently

(14,15)

(13)

. Both models were evaluated

for asphaltene precipitation from conventional oils and bitumens and were

confirmed to give good predictions for phase boundaries and the amount of asphaltene precipitation. The PC-SAFT predictions for asphaltene precipitation were slightly closer to the experimental data. A disadvantage of both the PC-SAFT and the CPA equation of state is that they can have more than three roots resulting in a more complex flash and longer computation time. 2 ACS Paragon Plus Environment

Page 3 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

Cubic equations of state (CEoS) are used in most commercial reservoir and process simulation software due to their simplicity of application and relatively fast computation times. CEoS models and associated correlations have been developed and tuned for light hydrocarbon systems. While CEoS have been applied to bitumen systems, they have not yet been successfully tested on a broad range of phase behavior for bitumens and solvents. Mehrotra and Svrcek

(2)

successfully modeled solubility data for mixtures of Cold Lake bitumen and nitrogen, carbon dioxide, methane, and ethane using the Peng-Robinson (PR) EoS. Jamaluddin et al.

(15)

used the

Martin EoS to fit solubility and saturated liquid density for systems of bitumen and carbon dioxide. Shaw and coworkers

(16-18)

have used a group contribution method in conjunction with

the PR EoS to capture LLV phase behavior of Athabasca vacuum residue and solvent. Castellanos-Díaz et al.

(19)

(20)

and Agrawal et al.

modeled the phase behavior of bitumen-

propane, bitumen-carbon dioxide and bitumen-pentane systems using the Advanced PengRobinson (APR) EoS. Their model correctly predicted saturation pressures and asphaltene onset points, but did not accurately predict asphaltene yields. The model predicted that the asphaltenes become soluble at high solvent dilutions while experimental data showed that the asphaltenes remain insoluble.

A possible remedy, recommended by Agrawal et al.

(20)

, is to apply asymmetric mixing rules

rather than the symmetric van der Waals mixing rules used in the APR model. The first asymmetric mixing rules proposed were empirically developed compositional dependent expressions for the binary interaction parameters (kij) used in the symmetric van der Waals mixing rules (21-24). Another class of asymmetric mixing rules are derived by equating the excess free energy (either Gibbs or Helmholtz) from an equation of state to the excess free energy from an activity coefficient model. The Huron-Vidal and Wong-Sandler mixing rules are well known examples of excess energy mixing rules

(25-26)

. The Huron-Vidal asymmetric mixing rules are

theoretically applicable to mixtures of nonpolar and polar compounds model mixtures of water and hydrocarbons

(28-30)

(27)

and have been used to

. Gregorowicz and de Loos

(31)

used the Wong-

Sandler mixing rules to model LLV behavior in asymmetric hydrocarbon mixtures. To date, a CEoS with asymmetric mixing rules has not been applied to mixtures of bitumen and solvents.


Energy & Fuels

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 37

The objective of this study is to evaluate the applicability of asymmetric mixing rules with a CEoS to model bitumen-solvent phase behavior. The Peng-Robinson CEoS is used to model saturation pressures and asphaltene precipitation yield data from Johnston et al. (32) for n-pentane diluted bitumen. First, the asymmetric van der Waals mixing rules (kij ≠ kij) are evaluated. Two forms of Huron-Vidal mixing rules rule. Then, the Sandoval et al.

(23)

(25)

are evaluated to test the excess energy type of mixing

mixing rules are evaluated to test the composition-dependent

van der Waals type of mixing rule. Finally, the form of the compositional dependence required to fit asphaltene yield data is examined. The advantages and shortcomings of each type of model are discussed.

2. CEoS Model 2.1 Advanced Peng Robinson CEoS The Peng Robinson cubic equation of state (APR EoS) is given by (33):

P=

aα (Tr , ω ) RT − V − b V (V + b) + b(V − b)

(1)

where P is absolute pressure, T is absolute temperature, V is the molar volume, R is the universal gas constant, a is the constant that corrects for attractive potential of molecules, b is the constant that corrects for volume of the molecules, α(Tr,ω) is a function specific to the equation of state, Tr is the reduced temperature, and ω is the acentric factor. The constants a and b for a pure component are related to its critical properties as follows:

ai = 0.457235 R 2Tci2 Pci

(2)

bi = 0.0777969 RTci Pci

(3)

where Tc is the critical temperature, Pc is the critical pressure, and the subscript i denotes the component.

CEoS such as the Peng-Robinson EoS do not accurately predict liquid phase densities and are typically corrected using volume translation. The Peneloux volume translation (34) was applied to the PR EoS by Jhaveri and Youngren (35) as follows:

=

,

−

(4)


Page 5 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

where c is the volume translation. The Advanced Peng Robinson equation of state (APR EoS) (36) is the Peng Robinson EoS with volume translation implemented in the VMGSimTM software used for this study.

2.2 Mixing Rules The mixing rules considered in this study are the asymmetric van der Waals (AvdW), compositionally dependent van der Waals (CDvdW), the Sandoval/Wilczek-Vera/Vera Mixing Rules (SWV), and the Huron-Vidal (HV) mixing rules. The symmetric van der Waals mixing rules have been tested in a previous study (20) and are used as basis of comparison for this study. Each is described below.

2.2.1 van der Waals Mixing Rules The classical van der Waals mixing rules are given by:

a m = ∑ ∑ x i x j a i a j (1 − k ij )

(5)

bm = ∑ xi bi

(6)

i

j

i

where am and bm are the EoS constants for the mixture, xi is the mole fraction of the ith component, and kij is a binary interaction parameter. Symmetric (SvdW): When kij = kji, the van der Waals mixing rules are symmetric and the model has only one interaction parameter per binary pair. For bitumen-solvent systems, CastellanosDíaz et al. (19) recommended the modified Gao et al. (37) correlation to estimate symmetric binary interaction parameters as follows:

= 1 −

(7)

where nij is an adjustable exponent with a default value of 0.27. Agrawal et al. (20) recommended temperature dependent binary interaction parameters of the following form: = 1+

# "

+ $%&'

(8)

where kij1 and kij2 are constants and is given by:


Energy & Fuels

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

=1−(

)

Page 6 of 37

(9)

In this case, n, kij1, and kij2 are adjusted to fit the data. Asymmetric (AvdW): If kij ≠ kji in Equation 5, the mixing rules are asymmetric. The asymmetric van der Waals mixing rules have two interaction parameters per binary pair and are more capable of matching the phase behavior of highly non-ideal mixtures.

Compositionally Dependent (CDvdW): When the kij in Equation 5 is defined as a function of composition, the mixing rules are considered to be asymmetric, even though kij can be set equal to kji at a fixed composition. In this case, the asymmetry is with composition and the mixture can be traditionally symmetric or asymmetric at any fixed composition depending on how the interaction parameters are calculated. The fitting of the asymmetric binary interaction parameters specifically for asphaltene-solvent binaries is discussed further in Section 3.2.

2.2.2 Sandoval/Wilczek-Vera/Vera Mixing Rules (SWV) Sandoval et al.

(23)

proposed several new forms of asymmetric, compositionally dependent

mixing rules to better match data for asymmetric ternary systems. The Sandoval/WilczekVera/Vera form of mixing rule (SWV) tested in this study is expressed as follows:

a m = ∑ ∑ x i x j a i a j (1 − δ ij )

(10)

* = + ,1 − - − - . + - + -

(11)

i

where δij is given by:

j

" " + =

(12)

The SWV model has two binary interaction parameters per binary pair, kij and kji. If kij = kji, the mixing rules reduce to the symmetric van der Waals mixing rules. This mixing rule can also be defined as a compositionally dependent van der Waals mixing rule but is labeled as SWV to distinguish it from the CDvdW mixing rule defined previously.

2.2.3 Huron-Vidal Mixing Rules (HV)


Page 7 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

For a CEoS, the limit of the Gibbs energy as pressure approaches infinity (gE∞) is defined as follows:

 a (T ) a (T )  Λ g ∞E = − m − ∑ xi i b b i  m 

(13)

where Λ is a function of the equation of state and, for the Peng-Robinson EoS, is given by: Λ=

2+ 2  ln 2 2  2 − 2  1

(14)

Equation 13 can be rearranged to obtain the following expression for am as a function of gE∞ (25,36)

:

 a gE a m = bm  ∑ i xi − ∞ Λ  bi

  

(15)

The NRTL activity coefficient model can be used to calculate gE∞ in Equation 15. The NRTL model is used to represent highly non-ideal mixtures where composition at a microscopic level deviates from the overall composition. Hence, it is referred to as a local composition model. The Gibbs excess energy at infinite pressure from the NRTL model is given by: n

n

g ∞E = ∑ xi i =1

∑x j =1

g ji − g ii   (g ji − g ii ) exp − α ji RT   n g − g ii   xk exp − α ki ki  ∑ RT   k =1

j

(16)

where gji is the interaction energy between component j and i, and αji is a constant characteristic of the mixture called the non-randomness parameter. αji is analogous to the inverse of the coordination number of a lattice. In the NRTL model, there are three interaction parameters per binary pair that can be tuned to match experimental data: (gji - gii), (gij - gjj), and αij.

Huron and Vidal introduced the parameter to correct for the volume of molecules, bj, into the calculation of gE∞ as follows:


Energy & Fuels

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

n

n

g ∞E = ∑ xi

g ji − g ii   (g ji − g ii ) exp − α ji RT   n g − g ii   xk bk exp − α ki ki  ∑ RT   k =1

∑x b j

j =1

i =1

Page 8 of 37

j

(17)

This addition allows the Huron-Vidal mixing rules to be reduced to the classical van der Waals mixing rules by substituting the following values for αij, gii, and gjj:

/ = 0

(18)

1 = −$%2

1 = −2

1 1 1 −

(19) (20)

In this work, the mixing rules obtained by using Equation 16 with Equation 15 are referred to as the mole-based Huron-Vidal mixing rules and those obtained by using Equation 17 with Equation 15 as the volume-based Huron-Vidal mixing rules. Both the mole-based and volumebased Huron-Vidal mixing rules were evaluated. The mole-based Huron-Vidal mixing rules were tested using the VMGSim Gibbs Excess Peng Robinson (GEPRTM) property package, which is commercially available

(36)

. The volume-based Huron-Vidal mixing rules were tested using an

unreleased beta version of VMGSim. In the VMGSim implementation of both versions of the Huron-Vidal mixing rules (mole-based and volume-based), the temperature dependence of the interaction energies is given by: 3 3

= 4 +

5

+ 6 $%&

(21)

where Aij, Bij and Cij are tuning parameters. In this work, Aij, Aji, Cij, and Cji were all set to zero, and only Bij, Bji, and αij were tuned to fit experimental data.

3.0 Modeling Methodology Figure 1 summarizes the overall modeling methodology. The fluid is characterized into pseudocomponents with assigned properties based on a distillation assay. The characterization is input into the APR EoS model, the binary interaction parameters are set, and a flash calculation is performed. The interaction parameters are adjusted iteratively to fit the experimental data. The methodology is discussed in more detail below. 8 ACS Paragon Plus Environment

Page 9 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

Figure 1. Modeling methodology algorithm.

3.1 Oil Characterization In order to model bitumen/solvent systems using an EoS, the fluid was divided into components and pseudo-components, each with an assigned mole fraction, density, critical properties, and acentric factor. Each solvent was treated as an individual component with known properties. The maltene and C5-asphaltene (n-pentane extracted) fractions of the bitumen were characterized separately following the methodology recommended by Castellanos-Díaz et al.

(19)

and Agrawal

et al. (20). The maltene fraction was divided into pseudo-components based on a distillation assay. The asphaltene characterization was based on a gamma distribution of molecular weights. This 9 ACS Paragon Plus Environment

Energy & Fuels

characterization has been shown to match vapor-liquid and liquid-liquid equilibrium data for bitumens diluted with a variety of solvents at temperatures from 20 to 180°C and pressures up to 10 MPa

(20)

. The characterization of the maltenes and C5-asphaltenes are discussed in detail

below.

3.1.1 Maltene Characterization The bitumen characterization was prepared for WC-B-B2 bitumen from Western Canada (the bitumen properties are provided in Section 4). The C5-asphaltene content of this bitumen is 19.4 wt%.

The remaining 80.6 wt% of the bitumen is defined as maltenes. The maltene

characterization was based on distillation data, in this case, a spinning band distillation assay from Agrawal et al. (20). The distillable fraction represents approximately 30% of the bitumen by volume. The remainder of the normal boiling point curve was extrapolated from the distillation data using a Gaussian distribution (linearly extrapolating on probability coordinates), Figure 2.

600 Distillation Data 550

Extrapolation

500 450

NBP, °C

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 37

400 350 300 250 200 0

20 40 60 Cumulative Mass % Distilled

80

Figure 2. Extrapolated NBP data for the maltene fraction of the WC-B-B2 bitumen.

The next step was to divide the distillation curve into pseudo-components. Castellanos-Díaz et al.

(19)

found that number of pseudo-components did not have a significant impact on the

prediction of saturation pressures as long as there were are least two pseudo-components in the 10 ACS Paragon Plus Environment

Page 11 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

light oil section and two in the heavy oil section. In this work, ten maltene pseudo-components were used to ensure that the entire range of maltene properties is accurately represented when modeling liquid-liquid phase behavior: three in the light oil section (200°C – 375°C), five in the medium oil section (375°C- 515°C), and two in the heavy oil section (515°C – 585°C). The molecular weight and liquid density of each pseudo-component were constrained so that the mixture matched the average maltene molecular weight and liquid density of 450 g/mol and 1005 kg/m³, respectively. The Lee-Kesler correlations (summarized in Appendix B) were used to estimate the critical properties and acentric factor from the average boiling point of each pseudocomponent. The pseudo-component properties are listed in Table 1 and shown in Figure 3.

Table 1. Bitumen pseudo-components and their properties. Mole MW NBP Density Mass Fraction Fraction (g/mol) (°C) (kg/m3) Malt1 0.059 0.127 246 236 901 Malt2 0.082 0.140 308 295 934 Malt3 0.098 0.136 382 355 967 Malt4 0.085 0.102 439 397 988 Malt5 0.075 0.083 479 425 1001 Malt6 0.081 0.082 520 452 1015 Malt7 0.082 0.077 563 480 1028 Malt8 0.084 0.074 607 508 1040 Malt9 0.085 0.068 658 540 1054 Malt10 0.075 0.056 709 571 1067 Asph1 0.025 0.011 1190 667 1070 Asph2 0.025 0.009 1491 706 1075 Asph3 0.031 0.010 1758 730 1080 Asph4 0.053 0.014 2084 751 1085 Asph5 0.060 0.011 2790 777 1090

Name

Pc Tc (kPa) (°C) 2655 445 2269 503 1948 561 1750 600 1633 625 1523 650 1421 675 1326 700 1227 727 1136 754 758 817 655 845 603 863 564 879 515 897


Energy & Fuels

3000 MW (g/mol) NBP (°C) Density (kg/m³) Tc (°C) Pc (kPa)

2500

Property Value

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 37

2000 1500 1000 500 0

Pseudocomponent

Figure 3. Bitumen pseudo-component properties.

3.1.2 Asphaltene Characterization The C5-asphaltenes were characterized one of two ways; 1) based on a gamma distribution of molecular weights as recommended by Agrawal et al.

(20)

, or; 2) as a single lumped pseudo-

component. The gamma distribution characterization described below was used in the AvdW and HV models for consistency with previous studies. A single pseudo-component was used to represent the asphaltene fraction in the SWV and CDvdW models to simplify the tuning process for these two models. For all the asymmetric models investigated, similar model results can be generated with either one or five asphaltene pseudo-components although the kij values do need to be slightly adjusted for the different characterizations.

Gamma Distribution Characterization: The Gamma distribution is given by: ( MW − MWm ) β −1 f ( MW ) = Γ( β )

β

   MW − MWm  β   exp β   MWavg − MWm   MWavg − MWm 

(22)

where f is molar frequency, MW is molecular weight, β is a shape factor, Γ is the gamma function, and subscripts m and avg indicate the monomer and average, respectively. To generate the molecular weight distribution, the gamma distribution parameters, MWm, MWavg and β, were


Page 13 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

set to 800 g/mol, 2000 g/mol and 2.5, respectively, as recommended by Agrawal et al. (20). Once the MW distribution was generated, the asphaltenes were divided into five pseudo-components of approximately equal weight fraction.

The true boiling point (TBP) of the asphaltene pseudo-components was calculated using the Søreide correlation Twu correlations et al.

(19)

(38)

, while critical properties and acentric factor were calculated using the

(39)

. The characterizations reported by Agrawal et al.

(20)

and Castellanos-Díaz

showed an increase in critical pressure with an increase in molecular weight for some

pseudo-components, contrary to the overall trend in the maltene properties. To correct this problem, the critical properties and acentric factor were recalculated using the aggregate molecular weight from the gamma distribution and the monomer density as inputs. The monomer density was estimated based on density data for fractionated asphaltenes from Barrera et al.

(40)

.

The values used for each asphaltene pseudo-component are summarized in Table1. The use of these values for molecular weight and density results in a continuously decreasing trend in critical pressure for the whole bitumen, Figure 3. Since the monomers in an asphaltene nanoaggregate are held together by van der Waals forces (which are very weak compared to the covalent bonds in the monomers), it is assumed that the properties of a pseudo-component representing a nano-aggregate will be closer to the monomer density than a density calculated based on the aggregate molecular weight.

Single Component Characterization: The single lumped pseudo component was characterized by defining molecular weight and density as 1800 g/mol and 1120 kg/m3 respectively. TBP was calculated using the Søreide correlation (38). The critical properties and acentric factor calculated from Twu correlations. These correlations are summarized in Appendix B. Table 2 compares the single component characterization with the range of properties in the Gamma distribution based characterization.


Energy & Fuels

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 37

Table 2. Asphaltene pseudo-components and their properties for the Gamma distribution based and single component characterizations. Property Gamma Distribution Characterization Single Component Characterization MW (g/mol) 1190 - 2790 1800 NBP (°C) 667 - 777 721 Density (kg/m3) 1070 - 1090 1120 Tc (°C) 817 - 897 906 Pc (kPa) 758 - 515 1057

3.2 Initialization, Flash Calculation, and Binary Interaction Parameters After the oil was characterized into pseudo-components, initial values were generated for the binary interaction parameters. Flash calculations were then performed using VMGSimTM Versions 8.0 to 9.5. VMGSimTM combines the Rachford-Rice algorithm to perform flash calculations with a stability test to minimize the Gibbs free energy. The binary interaction parameters were iteratively modified to obtain the best match to experimental data.

Note that the interaction parameters between the pseudo-components were set to zero to simplify the comparison between different sets of mixing rules. In other words, the only non-zero interaction parameters were between the solvent and each oil component. This assumption had little effect on the model predictions. The effect of the temperature dependence on the saturation pressures considered in this study was also negligible at temperatures up to 180°C (see appendix). Therefore, the temperature dependence was also set to zero (kij1=0 and kij2=0) in order to simplify the preliminary fitting exercises used to compare the different types of mixing rule. Temperature dependent binary interaction parameters were only implemented in the most successful model (CDvdW) in order to generate P-X diagrams ranging from 23 to 280°C. To tune the model, the asphaltene/solvent binary interaction parameters were adjusted to match asphaltene onset and yield data. If a good match for asphaltene onset and yield data was found, the maltene/solvent binary interaction parameters were adjusted to match saturation pressure data.

4. Dataset


Page 15 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

The modeling approach was tested on phase behavior data for n-pentane diluted bitumen reported by Johnston et al.

(32)

. The data were collected for two bitumen samples from the same

Western Canadian source reservoir (WC-B-B2 and WC-B-B3); selected bitumen properties are provided in Table 3. The dataset includes: saturation pressures at 90 to 280°C for feed compositions up to 65 wt% n-pentane; asphaltene yields and onsets at ambient conditions up to 250°C and 13.8 MPa, and; phase compositions for the asphaltene-rich heavy phase at 180°C and 4.8 MPa. The phase boundaries and yields are shown in Figures 4 and Figure 5, respectively. The phase compositions are provided in Table 4. In addition, asphaltene yield data for n-heptane diluted WC-B-B2 bitumen were available at ambient conditions, Table 4, and were used to evaluate the recommended model.

Table 3: Selected properties of WC-B-B2 and WC-B-B3 bitumen. Property

WC-B-B2

WC-B-B3

Specific Gravity

1.015

1.020

Viscosity at 50°C, 1 atm

2,900

3,100

Saturates, wt%

17

-

Aromatics, wt%

46.9

-

Resins, wt%

16.7

-

C5-asphaltenes, wt%

19.4

19.2


Energy & Fuels

1

6

CDvdW model

CDvdW model

0.6

L

L1L2

0.4

0.2

4

3

L

L1L2

2

1

L1V

0 0

10

L1V

0

20 30 40 50 60 Feed Pentane Content, wt%

70

0

10


70

15

6 Psat 5

Psat

c) 120°C

140°C onset

onset

CDvdW model

CDvdW model

4

L

Pressure (MPa)

Pressure (MPa)

b) 90°C

onset

5

Pressure (MPa)

Pressure (MPa)

0.8

Psat

a) 23°C

onset

L1L2

3

2

d) 180°C

10

L

L1L2

5

1

L1V

L1V 0

0

0

10


70

0

10


70

15

15 Psat

Psat

e) 230°C

250°C onset

f) 280°C

CDvdW model

CDvdW model 10

L

Pressure (MPa)

Pressure (MPa)

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 37

L1L2

5

L

10

L1L2

5

L1V

L1V

0

0

0

10


70

0

10


70

Figure 4. Pressure-Composition diagrams for n-pentane diluted bitumen at a) 23°C and b) 90°C c) 120°C d) 180°C e) 230°C f) 280°C. Symbols are experimental data from Johnston et al. (32) and lines are CDvdW model. 16 ACS Paragon Plus Environment

Page 17 of 37

20

C5-Asphaltene Yield, wt%


20

15

10 Data 5

CDvdW Model

15

10

Data 5 CDvdW Model

a) 23°C, 0.1 MPa

b) 90°C, 4.8 MPa

0

0 40


40

100

20

50 60 70 80 Feed Pentane Content, wt%

90

20



15

10

Data 5

CDvdW Model

15

10 4.8 MPa 13.8 MPa

5

CDvdW Model

d) 180°C

c) 140°C, 4.8 MPa 0

0 40


90

40


90

20


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

15

10

Data 5

CDvdW Model

e) 250°C, 10.3 MPa 0 40


90

Figure 5. Asphaltene yields for n-pentane diluted bitumen at a) 23°C b) 90°C c) 140°C d) 180°C and e) 250°C. Symbols are experimental data from Johnston et al. (32) and lines are the CDvdW model, and ‘C5-’ denotes pentane insoluble asphaltenes. 17 ACS Paragon Plus Environment

Energy & Fuels

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 37

Table 4. Heavy phase compositions of n-pentane diluted bitumen, measured using the CDvdW model at 180°C and 4.8 MPa. Solvent Content in Feed (wt%) 59.2 63.7 72.5

Measured Solvent (wt%) 31.4 19.2 20.3

Maltene (wt%) 20.4 24.0 23.7

(32)

and modeled

CDvdW Model C5-Asph. (wt%) 48.2 56.8 56.0

Solvent (wt%) 3.1 2.5 1.8

Maltene (wt%) 6.0 4.6 3.2

C5-Asph. (wt%) 90.9 92.9 95.0

Table 5. Asphaltene yields for n-heptane diluted bitumen (WC-B-B2) at ambient conditions (41). n-Heptane Content wt % 57.5 58.0 59.5 59.7 61.5 62.0 64.5 65.5 64.3 64.8 66.7 71.4 75.0 80.0 88.3 91.1

Asphaltene Yield wt % 1.2 1.4 2.7 2.1 3.8 3.5 5.7 5.2 5.1 5.0 6.5 9.3 11.5 13.1 15.5 16.0


Page 19 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

5. Results The modeling results and discussion begins with the symmetric (SvdW) and asymmetric van der Waals (AvdW) mixing rules followed by the mole-based and volume-based Huron-Vidal (HV) excess energy mixing rules. Then, the Sandoval et al. (SWV) compositional-dependent van der Waals mixing rules are evaluated and, finally, the compositional-dependence required to fit asphaltene yield with van der Waals mixing rules (CDvdW) data is investigated and an expression for kij is recommended.

In order to minimize the time required for the model evaluations, the mixing rules for most cases were first tested with a single component asphaltene and single component maltene characterization. The model with this simplified characterization gives similar results to the full characterization as shown in Appendix A. If the preliminary model was able to qualitatively fit asphaltene yield data, the model was then evaluated on the ambient condition data with the full characterization described in Section 3. The interaction parameters were tuned to fit both the asphaltene yield curve and the saturation pressures. Only the most successful model (CDvdW) was extended to the full range of temperatures and pressures in the dataset.

5.1 Symmetric van der Waals Mixing Rules Agrawal et al.

(20)

and Johnston et al.

(32)

demonstrated that the APR with symmetric van der

Waals mixing rules (kij = kji) can be used to match the asphaltene onset point for solvent diluted bitumen but severely under-predicts asphaltene yields at high dilution and over-predicts the pressure dependence of the onset point. Figure 6a confirms that the use of symmetric mixing rules (SvdW) under-predicts the asphaltene yields. Figure 6b shows that these mixing rules slightly under-predict the saturation pressures but are within the accuracy of the measurement.


Energy & Fuels

20

1,000

AvdW Model A AvdW Model B SvdW Model Yield Data

15

Saturation Pressure, kPa

Asphaltene Yield, wt%

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 37

10

5

AvdW Model A AvdW Model B SvdW Model Psat Data

800

600

400

200

(b)

(a) 0

0 20

40 60 80 Pentane Feed Content, wt%

100

0

10 20 30 40 50 60 Pentane Feed Content, wt%

70

Figure 6. Measured and modeled (symmetric model and asymmetric Models A and B) phase behavior data for n-pentane diluted bitumen: a) asphaltene yield at 21°C and 0.1 MPa; b) saturation pressures of n-pentane diluted bitumen at 90°C. Data from Johnston et al. (32). SvdW and AvdW indicate symmetric and asymmetric van der Waals mixing rules, respectively.

5.2 Asymmetric van der Waals Mixing Rules The use of asymmetric van der Waals mixing rules can improve the asphaltene yield match. Figure 6 shows the modeled yields and saturation pressures from the preliminary evaluation where the kij values for pseudo-component/pseudo-component binaries were set to zero and no temperature dependence was introduced into the other interaction parameters (kij = kijo). Figure 6a shows the modeled yields at 21°C for two different sets of tuning parameters that illustrate the benefits and limitations of the AvdW model (Models A and B). Both AvdW Model A and AvdW Model B have asymmetric interaction parameters for maltene-solvent and asphaltene-solvent binary pairs, Figure 7a. AvdW Model A has a step-change in ks-pc and kpc-s values at the transition from maltene to asphaltene pseudo-components, whereas the ks-pc and kpc-s values in AvdW Model B have a continuous trend from maltene to asphaltene pseudo components. AvdW Model A matches the shape of the asphaltene curve but does not predict the asphaltene components precipitating in the correct order; it predicts lighter pseudo-components precipitating before heavier pseudo-components. AvdW Model B predicts that the heaviest pseudo-components precipitate first and fits the ultimate yield, but over-predicts yields at intermediate dilutions. 20 ACS Paragon Plus Environment

Page 21 of 37

Figure 6b shows the predicted saturation pressures for n-pentane diluted bitumen at 90°C. While the AvdW model predicts saturation pressures within the error bars of the reported data, the shape of the curve is not representative of the observed pressure trends; the AvdW model overpredicts the saturation pressures at low dilutions. Only the saturation pressure curve at 90°C is shown here, but the model was found to over-predict saturation pressures at low dilutions at temperatures up to 180°C. This shortcoming was observed in every case where asymmetry was introduced into the maltene-solvent interaction parameters.

0.14

0.20

0.12 0.1

0.06

0.18

ks-pc

0.16

kpc-s

0.14

Model A: kpc-s Model A: ks-pc Model B: kpc-s Model B: ks-pc

0.12

kij

0.08

kij

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

0.10 0.08

0.04

0.06 0.04

0.02

(a)

0

(b)

0.02 0.00

Pseudo Component

Pseudo Component

Figure 7. Binary interaction parameters at 21°C for solvent/pseudo-component binaries for npentane diluted bitumen: a) AvdW Models A and B; b) AvdW Model C.

The AvdW model with a third set of tuning parameters, AvdW Model C, is shown in Figure 8. In this model the interaction parameters between maltene-solvent binary pairs were kept symmetrical, and were generated using the modified Gao et al.

(37)

correlation with n = 0.62.

Asymmetry was introduced into the interaction parameters between asphaltene-solvent binaries only, Figure 7b. Asphaltene/solvent interaction parameters were also calculated from the modified Gao et al. correlation, using a fitted value of n = 2 to calculate kij and n = 0.45 to calculate kji. This model is an improvement over the SvdW model for asphaltene yield, and matches saturation pressures generally to within the error of the measurements at temperatures 21 ACS Paragon Plus Environment

Energy & Fuels

up to 180°C (only the saturation pressures at 90°C are shown here). However, this model still under-predicts asphaltene yield at high dilutions because the AvdW model does not capture the asymmetry of system. Properties of the asphaltene components such as molecular mass and density change with the addition of n-pentane due to self-association. Simply setting kij ≠ kji in the AvdW does not fully account for this asymmetry.

A significant disadvantage of this approach is the arbitrary division between the asphaltenesolvent and maltene-solvent binary interaction parameters. This means that this model is not predictive for asphaltene yields for other solvents or other temperatures; that is, the interaction parameters must be tuned for each condition. Note, as with any cubic equation of state based on an oil characterization, the exact values of the binary interaction parameters are sensitive to the number of pseudo-components and must be slightly retuned for each characterization (see Appendix A). 20

1,000 AvdW Model C

AvdW Model C

Yield Data

Psat Data

15


Asphaltene Yield, wt%

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 22 of 37

10

5

800

600

400

200

(a)

(b)

0

0

20

40 60 80 Pentane Feed Content, wt%

100

0

10 20 30 40 50 60 Pentane Feed Content, wt%

70

Figure 8. Measured and modeled (AvdW Model C) phase behavior data for n-pentane diluted bitumen: a) asphaltene yield at 21°C and 0.1 MPa; b) saturation pressures of n-pentane diluted bitumen at 90°C. Data are from Johnston et al. (32).

5.3 Huron-Vidal Mixing Rules


Page 23 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

Since the asymmetric van der Waals mixing rules were inadequate to match asphaltene yields, compositionally dependent mixing rules were considered. One way to introduce a compositional dependence is via an excess energy mixing rule such as the Huron-Vidal mixing rules. In this case, the full bitumen characterization summarized in Table 1 and Figure 3 was used. Several sets of tuning parameters for the Huron-Vidal mixing rules were evaluated for matching the onset and yield curve for asphaltene precipitation from n-pentane diluted bitumen. Figures 9 and 10 show the results for two different sets of tuning parameters that illustrate the benefits and limitations of this approach (HV Models A and B). HV Model A matches both the shape of the yield curve and the ultimate yield, Figure 9a. However, this model was tuned by introducing a significant step-change in bij and bji for precipitating components, Figure 9b. A step change is arbitrary and not easily correlated for other solvents, a major shortcoming. HV Model B preserves smooth trends in the interaction parameters, Figure 10b, but under-predicts yields, Figure 10a. Both models predict saturation pressures up to 180°C to within the measurement error, Figure 11. Only the saturation pressures at 90°C are presented here.

The volume-based mixing rules were also examined (results not shown here) and yielded similar results. The model with volume-based HV mixing rules can better match asphaltene yields than the symmetric model, but it still under predicts asphaltene yield at high dilutions of n-pentane and requires tuning parameters that are not easily generalized. The HV models can better match asphaltene yields than the symmetric model because the effect of asphaltene association can somewhat be represented as a local variation in concentration. However, the local composition model does not fully capture the effect of asphaltene association.


Energy & Fuels

20

10000

(a)

0.2

(b)

8000 0.15

15

7000

10

5000

0.1

Alpha

6000

Bij


9000

4000 3000 5 Yield Data

Bij

2000

Bji

1000

Alpha

0

HV Model A

0.05

0

0 40


100

Pseudo Component

Figure 9. a) Measured and modeled asphaltene yield data from n-pentane diluted bitumen at 21°C and 0.1 MPa and b) tuning parameters used in HV Model A.

20 10000

(a)

9000

0.2

(b)

8000

15

0.15 7000

10

5000

0.1

Alpha

6000

Bij


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 37

4000 3000 2000

5 Yield Data

Bij Bji

0.05

Alpha

1000

HV Model B

0

0

0 40


100 Pseudo Component

Figure 10. a) Measured and modeled asphaltene yield data from n-pentane diluted bitumen at 21°C and 0.1 MPa and b) tuning parameters used in HV Model B.


Page 25 of 37

1000 Psat Data HV Model A HV Model B

800


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

600

400

200

0 0

10


70

Figure 11. Saturations pressures for n-pentane diluted bitumen at 90°C. Symbols are experimental data and lines are modeled using HV Model A and HV Model B.

5.4 Compositionally Dependent Mixing Rules Another way to introduce a compositional dependence is to relate the asymmetric van der Waals binary interaction parameters to composition. The Sandoval et al. mixing rules, the simplest of the compositionally dependent mixing rules, assume a linear relationship between the asymmetric binary interaction parameter and the mole fractions of the asymmetric pair of components. The Sandoval et al. mixing rules (SWV) were tested with the simplified bitumen characterization with a single pseudo-component for each of the maltenes and the asphaltenes. The binary interaction parameters were set to zero for maltene/solvent and maltene/asphaltene binaries and no temperature dependence was introduced. The solvent/asphaltene binary interaction parameters are calculated for each phase individually and were tuned to match the yield data. The SWV mixing rules are an improvement over the symmetric mixing rules, but cannot match ultimate yield, Figure 12. Therefore, no further tests were performed. It appears that a linear relationship between the binary interaction parameter and the mixture composition is not sufficient to match asphaltene yields at high dilutions.


Energy & Fuels

20 Yield Data SWV Model


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 26 of 37

15

10

5

0 40

50


100

Figure 12. Asphaltene yield data for n-pentane diluted bitumen at 21°C and 0.1 MPa. Symbols are data from Johnston et al. (32), solid line is the SWV model with kas = -0.05 and ksa = 0.0015.

Since the Sandoval et al.

(23)

mixing rules could not adequately match asphaltene yield data for

solvent diluted bitumen, the form of the compositional dependence necessary to capture asphaltene yield data was investigated. For the sake of simplicity, the asphaltene fraction of the oil was represented by only one asphaltene pseudo-component when determining the compositional dependence of the asphaltene/solvent binary interaction parameters, ka-s. The maltene fraction was represented by the ten pseudo-components summarized in Table 1. The binary interaction parameters for maltene/maltene and maltene/asphaltene pairs were set to zero. In order to accurately model saturation pressures at temperatures above 180°C, solvent/maltene binary interaction parameters were estimated using Equation 8 with k1ij = 2550 and k2ij = -1.50 as recommended by Agrawal et al. (20).

The value of ka-s was adjusted at each temperature, pressure, and composition to fit asphaltene yields at temperatures from 20 to 250°C and pressures up to 13.8 MPa. The fitted ka-s are shown in Figure 13 as a function of the feed composition. Since the ka-s appear to be insensitive to pressure, they were correlated to temperature and feed composition as follows:

7 = 7 = 0.00001& + 0.0345< =.>?.?@>AB C

(23) 26

ACS Paragon Plus Environment

Page 27 of 37

where ws is the mass fraction of solvent in the feed. The correlation is also shown in Figure 13.

0.45 250C 10.3 MPa

0.40

180C 13.8 MPa 180C 4.8 MPa

0.35

140C 4.8MPa 90C 4.8 MPa

0.30

ka-s

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

20C 0.1 MPa Correlation

0.25 0.20 0.15 0.10 0.05 0.00 40

60 80 Feed Pentane Content, wt%

100

Figure 13. Values of ka-s used to fit asphaltene yield. Symbols are fitted ka-s and the lines are Equation 23.

The model with van der Waals mixing rules and the compositionally dependent kij shown in Equation 23 (the CDvdW model) was used to generate PX phase diagrams and yield curves, Figures 4 and 5. The CDvdW model matches the saturation pressures and yield curves over the whole composition range. It also matches the pressure dependence of the LL boundary at high temperatures, Figure 4, which is an improvement over APR model with symmetric van der Waals interaction parameters. However, like all of the mixing rules evaluated in this study, the CDvdW model significantly under-predicted the solvent content in the asphaltene-rich phase at temperatures above 100°C, Table 4. The low predicted solvent contents are consistent with the values of less than 6 wt% reported at ambient conditions (42).

A similar fitting exercise was performed to match asphaltene yield data from n-heptane diluted WC-B-B2 bitumen data from Shafiee

(41)

, Figure 14. Similar trends are observed in the fitted n-

pentane/asphaltene and n-heptane/asphaltene binary interaction parameters at ambient conditions, Figure 15. The similarity of the trends is promising for the development of a 27 ACS Paragon Plus Environment

Energy & Fuels

correlation to relate ka-s to solvent properties; however, more data are required before proceeding further.


20

15

10

5 nC5 Yield Data nC7 Yield Data Model 0 40

60 80 Feed Diluent Content, wt%

100

Figure 14. Asphaltene yields from n-pentane diluted bitumen (32) and n-heptane diluted bitumen (41) . Symbols are experimental data and lines are the fitted CDvdW model.

0.12 n-C5 n-C7 Correlation

0.10

0.08

ka-s

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 37

0.06

0.04

0.02

0.00 40

60 80 Feed Diluent Mass Fraction

100

Figure 15. Values of ka-s used to fit asphaltene yield for n-pentane diluted bitumen and n-heptane diluted bitumen at ambient conditions. Symbols are fitted ka-s and the lines are lines of best fit.


Page 29 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

Note the compositional model presented above is not a full compositional implementation because ka-s was fit for the feed composition and not recalculated for each phase at its individual phase composition. This extra recalculation step was not taken because the beta version software used for this project did not have that capability. Some preliminary tests were performed with a full compositional implementation by translating the fitted kij’s into equivalent Huron Vidal inputs. The fully compositional Huron-Vidal model was run with these inputs and the results were similar to the CDvdW model. In addition, the methodology used here is similar to the methodology applied by Adachi and Sugie

(21)

to the a parameter of the equation of state. They

demonstrated that similar results were obtained when fitting to the feed composition or recalculating for each phase composition.

6. Conclusions Several forms of mixing rules were evaluated with the APR CEoS to model asphaltene yields and phase boundaries of n-alkane diluted bitumen. While the asymmetric mixing rules are an improvement over the APR with symmetric mixing rules, each of these models still has several limitations as noted in Table 6. None the models could accurately predict phase compositions at temperatures above 90°C.

Table 6. Summary of limitations for mixing rules examined in this study. Mixing Rule

Tunable to match

Tunable to match

Tuning parameters

asphaltene yield?

asphaltene yield and

easily correlated?

saturation pressures? AvdW

Yes

No

-

SWV

No

-

-

HV

Yes

Yes

No

CDvdW

Yes

Yes

No

The use of asymmetric van der Waals mixing rules (kij ≠ kji) adds more flexibility to the model and matches asphaltene yields. However, the tuning parameters that give the best match for 29 ACS Paragon Plus Environment

Energy & Fuels

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 30 of 37

asphaltene yields do not match VLE data with the same accuracy as symmetric mixing rules. The parameters that retain a good match to the saturation pressures require an arbitrary division between maltene/solvent and asphaltene/solvent binary interaction parameters and therefore are not predictive or easily correlated for other solvents. These models slightly over-predict asphaltene yields at intermediate dilutions and under-predict asphaltene yields at high dilutions.

The APR with Huron-Vidal excess energy based mixing rules can be tuned to match asphaltene yield. However, the tuning parameters that give the best match for asphaltene yield data require an arbitrary division between maltene/solvent and asphaltene/solvent binary interaction parameters and therefore are not predictive or easily correlated for other solvents. A model that can be predictive for other solvents cannot match the shape of the yield curve.

Mixing rules with compositionally dependent binary interaction parameters can successfully model asphaltene yield data but only with an appropriate mathematical form for the compositional dependence fits. For example, the Sandoval et al. mixing rules cannot be tuned to match asphaltene yields at high dilution, because the mixing rule includes a linear dependence of the binary interaction parameter on mole fraction. An exponential relationship to mass fraction provided better results and a correlation was found for the asphaltene/n-pentane binary interaction parameter as a function of the solvent mass fraction in the feed and absolute temperature. While the model is not very sensitive to the characterization methodology, the binary interaction parameters must be slightly retuned if the characterization is changed; for example, if the number of pseudo-components is altered.

Overall, CEoS are not well suited to modeling asphaltene precipitation. CEoS with asymmetric mixing rules cannot account for the changing properties of the asphaltene components as they associate in the presence of an alkane like n-pentane. The asphaltene nanoaggregates were represented with a fixed molecular weight distribution and there is no mechanism in a CEoS model to alter this distribution with changing conditions. Data fitting is possible over a limited range of conditions using asymmetric or compositionally dependent binary interaction parameters but a generalized correlation for a wide range of solvents and conditions remains to be achieved. 30 ACS Paragon Plus Environment

Page 31 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

Acknowledgements The authors are grateful for financial support from the sponsors of the NSERC Industrial Research Chair in Heavy Oil Properties and Processing, including NSERC, CNOOC Nexen, Petrobras, Schlumberger, Shell Canada Energy Ltd., Suncor Energy, and Virtual Materials Group.

References 1.

Dini, Y.; Becerra, M.; Shaw, J.M. Phase Behavior and Thermophysical Properties of Peace River Bitumen + Propane Mixtures from 303K to 393K. J. Chem. Eng. Data, 2016, 61, 2659-2668. 2. Mehrotra, A.K.; Svrcek, W.Y. Correlation and Prediction of Gas Solubility in Cold Lake Bitumen. Can. J. Chem. Eng., 1988, 66, 666-670. 3. Badamchi-Zadeh, A.; Yarranton, H.W.; Svrcek, W.Y.; Maini, B.B., Phase Behavior and Physical Property Measurements for VAPEX Solvents: Part 1. Propane and Athabasca Bitumen. JCPT, 2009, 48, 54-61. 4. Badamchi-Zadeh, A.; Yarranton, H.W.; Maini, B.B.; Satyro, M.A. Phase Behavior and Physical Property Measurements for VAPEX Solvents: Part 2. Propane, Carbon Dioxide and Athabasca Bitumen. JCPT, 2009, 48, 57-65. 5. Zou, X.; Zhang, X.; Shaw, J.M. The Phase Behavior of Athabasca Vacuum Bottoms + nAlkane Mixtures. SPE 97661, SPE International Thermal Operations and Heavy Oil Symposium, Calgary, Nov 1-3, 2005. 6. Tharanivasan, A.K.; Yarranton, H.W.; Taylor, S.D. Application of a Regular SolutionBased Model to Asphaltene Precipitation from Live Oils. Energy & Fuels, 2011, 25, 528538. 7. Alboudwarej, H.; Akbarzadeh, K.; Beck, J.; Svrcek, W.Y.; Yarranton, H.Y. Regular Solution Model for Asphaltene Precipitation from Bitumens and Solvents. AIChE J., 2003, 49, 2948-2956. 8. Wang, J.X.; Buckley, J.S. A Two-Component Solubility Model of the Onset of Asphaltene Flocculation in Crude Oils. Energy & Fuels, 2001, 15, 1004-1012. 9. Li, Z.; Firoozabadi, A. Cubic-Plus-Association Equation of State for Asphaltene Precipitation in Live Oil. Energy & Fuels, 2010, 24, 2956-2963. 10. Shirani, B.; Manouchehr, N.; Mousavi-Dehghani, S.A. Prediction of Asphaltene Phase Behavior in Live Oil with a CPA Equation of State. Fuel, 2012, 97, 89-96. 11. Panuganti, S.R.; Vargas, F.M.; Gonzalez, D.L.; Kurup, A.S.; Chapman, W.G. PC-SAFT Characterization of Crude Oils and Modeling of Asphaltene Phase Behavior. Fuel, 2012, 93, 658-669. 12. Ting, P.D.; Hirasaki, G.J.; Chapman, W.G. Modeling of Asphaltene Phase Behavior with the SAFT Equation of State. Petr. Sci. Technol., 2003, 21, 647-661. 31 ACS Paragon Plus Environment

Energy & Fuels

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 37

13. Tavakkoli, M.; Panuganti, S.R.; Taghikhani, V.; Pishvaie, M.R.; Chapman, W.G. Understanding the Polydiverse Behavior of Asphaltenes During Precipitation. Fuel, 2014, 117, 206-217. 14. AlHammadi, A.A.; Vargas, F.M.; Chapman, W.G. Comparison of Cubic-Plus-Association and Perturbed-Chain Statistical Associating Fluid Theory Methods for Modeling Asphaltene Phase Behavior and Pressure-Volume-Temperature Properties. Energy & Fuels 2015, 29, 2864-2875. 15. Jamaluddin, A.K.M.; Kalogerakis, N.E.; Chakma, A. Predictions of CO2 Solubility and CO2 Saturated Liquid Density of Heavy Oils and Bitumens using a Cubic Equation of State. Fluid Phase Equilibria, 1991, 64, 33-48. 16. Saber, N.; Shaw, J.M. Toward Multiphase Equilibrium Prediction for Ill-Defined Asymmetric Hydrocarbon Mixtures. Fluid Phase Equilibria, 2009, 285, 73-82. 17. Saber, N.; Shaw, J.M. On the Phase Behavior of Athabasca Vacuum Residue + n-Decane. Fluid Phase Equilibria, 2011, 302, 254-259. 18. Saber, N.; Zhang, X.; Zou, X.; Shaw, J.M. Simulation of the Phase Behavior of Athabasca Vacuum Residue + n-Alkane Mixtures. Fluid Phase Equilibria, 2012, 312, 25-31. 19. Castellanos-Díaz, O.; Modaresghazani, J.; Satyro, M.A.; Yarranton, H.W. Modeling the Phase Behavior of Heavy Oil and Solvent Mixtures. Fluid Phase Equilibria, 2011, 304, 7485. 20. Agrawal, P.; Schoeggl, F.F.; Satyro, M.A.; Taylor, S.D.; Yarranton, H.W. Measurement and Modeling of the Phase Behavior of Solvent Diluted Bitumens. Fluid Phase Equilibria, 2012, 334, 51-64. 21. Adachi, Y.; Sugie, H. A New Mixing Rule – Modified Conventional Mixing Rule. Fluid Phase Equilibria 1986, 28, 103-118. 22. Panagiotopoulos, A.Z.; Reid, M.C. A New Mixing Rule for Cubic Equations of State for Highly Polar, Asymmetric Mixtures. ACS Symposium Series, 1986, 300, 571-582. 23. Sandoval, R.; Wilczek-Vera, G.; Vera, J.H.Prediction of Ternary Vapor-Liquid Equilibria with the PRSV Equation of State. Fluid Phase Equilibria, 1989, 52, 119-126. 24. Stryjek, R.; Vera, J.H. PRSV – An Improved Peng-Robinson Equation of State with New Mixing Rules for Strongly Nonideal Mixtures. Can. J. Chem. Eng., 1986, 64, 334-340. 25. Huron, M.; Vidal, J., New Mixing Rules in Simple Equations of State for Representing Vapour-Liquid Equilibria of Strongly Non-Ideal Mixtures. Fluid Phase Equilibria, 1979, 3, 255-271. 26. Wong, D.S.H.; Sandler, S.I. A Theoretically Correct Mixing Rule for Cubic Equations of State. AIChE J., 1992, 38, 671-680. 27. Pedersen, K.S.; Christensen, P.L. Phase Behavior of Petroleum Reservoir Fluids; Taylor & Francis Group, 2006. 28. Kristensen, J.N.; Christensen, P.L. A Combined Soave-Redlich-Kwong and NRTL Equation for Calculating the Distribution of Methanol Between Water and Hydrocarbon Phases. Fluid Phase Equilibria, 1993, 82, 199-206. 29. Pedersen, K.S.; Michelsen, M.L.; Fredheim, A.O. Phase Equilibrium Calculations for Unprocessed Well Streams Containing Hydrate Inhibitors. Fluid Phase Equilibria, 1996, 126, 13-28. 30. Sørensen, H.; Pedersen, K.S.; Christensen, P.L. Modeling of Gas Solubility in Brine. Organic Geochemistry, 2002, 33, 635-642.


Page 33 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

31. Gregorowicz, J.; de Loos, T.W. Prediction of Liquid-Liquid-Vapor Equilibria in Asymmetric Hydrocarbon Mixtures. Ind. Eng. Chem. Res., 2001, 40, 444-451. 32. Johnston, K.A.; Schoeggl, F.F.; Satyro, M.A.; Taylor, S.D.; and Yarranton, H.W. Phase Behavior of Bitumen and n-Pentane. Fluid Phase Equilibria, 2017, 442, 1-19. 33. Peng, D.; Robinson, D.B. A New Two-Constant Equation of State. Ind. Eng. Chem. Fundam., 1976, 15, 59-64. 34. Péneloux, A.; Rauzy, E. A Consistent Correction for Redlich-Kwong-Soave Volumes. Fluid Phase Equilibria, 1982, 8, 7-23. 35. Jhaveri, B.S.; Youngren, G.K. Three Parameter Modification of the Peng Robinson Equation of State to Improve Volume Predictions. SPE Res. Eng., 1988, 3, 1033-1040. 36. VMG, VMGSim Version 9.5, VMGSim User’s Manual; Virtual Materials Group Inc.; Calgary, Canada, 2015. 37. Gao, G.; Dadiron, J.L.; Saint-Guirons, H.; Xans, P.; Montel, F.A. A Simple Correlation to Evaluate Binary Interaction Parameters of the Peng-Robinson Equation of State: Binary Light Hydrocarbon Systems. Fluid Phase Equilibria, 1992, 74, 85–93. 38. Søreide, I. Improved Phase Behavior Predictions of Petroleum Reservoir Fluids from a Cubic Equation of State; PhD Thesis, Norwegian Institute of Technology and Applied Geophysics; Trondheim, Norway, 1989. 39. Twu, C.H. An Internally Consistent Correlation for Predicting the Critical Properties and Molecular Weights of Petroleum and Coal-Tar Liquids. Fluid Phase Equilibria, 1984, 16, 137-150. 40. Barrera, D.M.; Ortiz, D.P.; Yarranton, H.W.; Molecular Weight and Density Distributions of Asphaltenes from Crude Oils. Energy & Fuels, 2013, 27, 2474-2487 41. Shafiee Neistanak, M. Kinetics of Asphaltene Precipitation and Flocculation from Diluted Bitumen; M.Sc. Thesis, University of Calgary; Calgary, Canada, 2014. 42. Yarranton, H.W.; Schoeggl, F.S.; George, S.; Taylor, S.D. Asphaltene-Rich Phase Compositions and Sediment Volumes from Drying Experiments. Energy & Fuels, 2011, 25, 3624-3633.


Energy & Fuels

Appendix A: Model Sensitivities Number of Components The predicted asphaltene yield is not very sensitive to the number of components in the characterization. With the AvdW model for example, representing the maltene and asphaltene fraction each with one pseudo-component results in a very similar asphaltene yield curve to the curve modeled with the full 16 pseudo-component characterization, Figure A1. The binary interaction parameters must be retuned for each characterization but have similar values for both cases, Table A1. While the one pseudo-component model is not recommended, (because saturation pressure behavior cannot be modeled with only one maltene pseudo-component), it is a straightforward starting point for tuning the model because the predicted yields and magnitudes of the kij are similar to the full characterization. Similar results were found with the other mixing rules. In all cases, when the characterization is changed, the model can be re-tuned to give similar results with small adjustments to the kij values.

25


1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 34 of 37

Yield Data Full Characterization Simplified Characterization

20

15

10

5

0 40


100

Figure A1. Measured and modeled asphaltene yield data for n-pentane diluted bitumen at 21°C and 0.1 MPa. Symbols are data from Johnston et al. (32), and lines are the AvdW model with a) bitumen characterized as only two components (simplified characterization) and b) bitumen characterized as described in Table 1 (full characterization).


Page 35 of 37

Table A1. Binary interaction parameters used in AvdW model with simplified and full characterizations. Component Pair Maltene-Solvent Asphaltene-solvent Solvent-Maltene Solvent-Asphaltene

Binary Interaction Parameter Simplified Characterization Full Characterization 0.02 0.030 - 0.046 0.04 0.047 – 0.055 0.08 0.095 – 0.114 0.11 0.116 – 0.125

Temperature Dependence Agrawal et al. (20) recommended a temperature dependence for kij to match saturation pressures at high temperatures. The temperature dependence has relatively little effect on the predicted saturation pressures up to 180°C but becomes significant at higher temperatures, Figure A2. For this reason and because the focus of this study was on modeling asphaltene yields, the temperature dependence for binary interaction parameters was neglected in all cases except the recommended CDvdW model.

3500

12000

Psat No Temperature Dependence Temperature Dependence

Psat No Temperature Dependence Temperature Dependence

10000


3000

Satiuration Pressure, MPa

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

2500 2000 1500 1000

8000

6000

4000

2000

500

(a)

(b)

0

0 0

10 20 30 40 50 60 Feed Pentane Content, wt%

70

0

10 20 30 40 50 60 Feed Pentane Content, wt%

70

Figure A2. Measured and modeled saturation pressure data for n-pentane diluted bitumen. Symbols are data from Johnston et al. (32) and lines are the SvdW model with and without temperature dependence.


Energy & Fuels

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 36 of 37

Appendix B: Additional Correlations

Twu method Reference Critical Temperature

&DE = & 0.533272 + 0.191017-10H & + 0.779681-10K & − 0.284376-10@ &H + 0.959468-10L ⁄&@H @

Reference Critical Pressure

D = 3.83354 + 1.19629N.> + 34.8888N + 36.1952N + 104.193NO N = 1−

& &D

Reference Critical Volume

PD = =1 − 0.419869 − 0.505839N − 1.56436NH − 9481.7N@O CL Reference Specific Gravity

QR = 0.843593 − 0.128624N − 3.36159NH − 13749.5N@ Deviation Parameters

SQR = + ,−0.182421 + 3.01721&.> .SQR Z


Page 37 of 37

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Energy & Fuels

WU = SQRU Y,2.53262 − 46.1955&.> − 0.00127885& .

+ ,−11.4277 + 252.140&.> + 0.00230535& .SQRU Z

Lee Kesler method Critical Pressure

0.0566 2.2898 0.11857 − 1-10H & 0.24244 + + QR QR QR 3.6480 0.47277 + 1-10K & 1.46850 + + QR QR 1.69770 − 1-10@ &H 0.42019 + QR

$%D = 8.3634 −

Critical Temperature

1-10> 0.4669 − 3.2623QR &D = 341.7 + 8.11QR + & 0.4244 + 0.1174QR + & Acentric Factor

[=

$% \

101.325 6.09648 _ D ] − 5.92714 + &^ + 1.28862$%&^ − 0.169347&^ &^ < 0.8 15.6875 _ 15.2518 − & − 13.4721$%&^ + 0.43577&^ ^

[ = −7.904 + 0.1352ab − 0.007465ab + 8.359&^ + @/H

1.8& ab = QR

1.408 − 0.01063ab &^ ≥ 0.8 &^

Søreide Correlation & = 1928.3 − 1.695 ∗ 10> fg .H> QR H.__ < =.O?∗hbO.K_L>∗ij.HO_∗hb∗ijC


Can a Cubic Equation of State Model BitumenâSolvent Phase Behavior?

Recommend Documents

Can a Cubic Equation of State Model BitumenâSolvent Phase Behavior?