Aspen Hydraulics Dynamics Reference

Aspen Dynamic Pipeline Solver Reference Guide

Who Should Read this Guide

This guide is intended as a reference aid for using Aspen Hydraulics functionality within the HYSYS 3.4 Oil & Gas Option and in particular the Aspen Dynamic Pipeline Solver used within Aspen Hydraulics.

Who Should Read this Guide

2

Contents INTRODUCING THE ASPEN DYNAMIC PIPELINE SOLVER................................. 4 ASPEN DYNAMIC PIPELINE SOLVER REFERENCE............................................. 5 Solution Procedure ................................................................................................... 5 Numerical Stability and the Courant Limit ............................................................... 6 Semi-implicit Methods.......................................................................................... 6 The SETS Method ................................................................................................ 6 Linearisation of the Finite Difference Equations...................................................... 10 Physical Properties ................................................................................................. 15 Required Properties ........................................................................................... 15 Closure Laws and Models ................................................................................... 17 Interfacial Friction ............................................................................................. 19 Wall Friction ..................................................................................................... 23 Heat Transfer Coefficients .................................................................................. 29 Cylindrical Wall Heat Conduction ......................................................................... 33 Interfacial Mass Transfer .................................................................................... 36

GENERAL INFORMATION............................................................................... 41 Copyright.............................................................................................................. 41

TECHNICAL SUPPORT.................................................................................... 43 Online Technical Support Center .............................................................................. 43 Phone and E-mail................................................................................................... 44

Contents

3

Introducing the Aspen Dynamic Pipeline Solver

The Aspen Dynamic Pipeline Solver is a code for modelling transient multiphase hydrocarbon flows in wells, pipelines and components. The Aspen Dynamic Pipeline Solver solves mass, momentum and energy equations for each phase using a one-dimensional finite difference scheme. Appropriate flow pattern maps and constitutive relationships are provided for wall and interfacial friction and heat transfer, and a model for multi-component phasechange is included.

Introducing the Aspen Dynamic Pipeline Solver

4

Aspen Dynamic Pipeline Solver Reference

The Aspen Dynamic Pipeline Solver is a transient multiphase flow model for oil and natural gas pipelines. The code is a six-equation, semi-implicit, finite difference computational fluid dynamics model. In order to close the equation set and obtain solutions, the code requires physical property data for the fluids, models for wall friction and interfacial friction based on mechanistic flow regime prediction, thermodynamic modelling of interfacial mass transfer, relationships for heat transfer between the fluids and the pipe walls and a model for two-phase critical flow.

Solution Procedure The finite difference method is very popular for numerically solving systems of partial differential equations, such as occur in single and multiphase flows. Its main feature is the replacement of the derivatives by finite differences of function values, ending up with equations, which have function evaluations at discrete values of the function variables (time and up to three space dimensions, in the case of fluid flow). This manual will not give an introduction to these methods, and familiarity with basic concepts such as implicit and explicit differencing is assumed. You should refer to standard texts such as Roache (1976). The most appropriate method for discretising a given set of fluid flow equations by finite difference methods will depend on various factors, including: •

Numerical stability - some ways of discretising the equations are numerically unstable: errors, however small, grow and the solution diverges. Numerical instability may be conditional, in which case the scheme may be used, but subject to certain restrictions (see Courant limit below) or unconditional, in which case the scheme cannot be used.

•

Computer time - the ideal is to find a scheme that will give accurate, stable results in the minimum of computer time.

•

Accurate representation of the continuous partial differential equations - for instance, some methods of Discretization may result in a finite difference equation, which tends, in the limit of, to a differential equation which is different to the original equation.


5

The methods used in many practical codes, therefore, are usually compromises between these various factors and, as such, can be expected to have ranges of invalidity, though these may be not too constraining.

Numerical Stability and the Courant Limit The question of numerical stability is a very important one and is closely connected with the type of differencing employed. The fully explicit methods, whilst easy to solve, are subject to strict limits on the size of the mesh used. If this condition is violated, the solution will become unstable. There is, in general, a condition, known as the Courant limit, such that (in a onedimensional system):

V∆t ≤1 ∆x (1)

The Courant limit was first identified by Courant, Friedrichs and Lewy (1928), who discussed the propagation of information in fluid flow simulations. Essentially, they realized that if the velocity of the fluid was fast enough, then information about momentum and other properties could be transported out of the current cell before the current timestep had finished, and this led to numerical instabilities. The Courant limit means practically that, given a system which has been discretized in space, there is thus a limit on the maximum size of the timestep. This means that any program using such a scheme to difference the flow equations may take a very long time to run.

Semi-implicit Methods In order to overcome the limitations of explicit methods, whilst retaining their computational simplicity, a class of finite difference methods has arisen known as the semi-implicit methods. These treat some of the terms in the differential equations explicitly and others, specifically those involved in the Information Propagation discussed above, implicitly. By doing this the idea is to eliminate, or at least relax, the Courant limit so that a larger timestep may be used, but without having the computational overhead of a fully implicit system each time. The SETS (Stability-Enhancing Two-Step) method, used in the Aspen Dynamic Pipeline Solver, is a semi-implicit method which treats the terms, and implicitly.

The SETS Method SETS is a two-step method, consisting of a basic step and a stabilizing step. The basic step is a semi-implicit equation set, and it provides information about pressure wave propagation. It treats the convective terms implicitly and this helps to relax the Courant limit on mesh size. However, studies by Mahaffy (1979, 1982) showed that in some circumstances numerical instabilities can arise and so the method is stability enhancing rather than totally stable. The second step is thus added as a stabilizing step, and it


6

provides information about the propagation of density, energy and momentum across cell boundaries.

A Simple Example - Single Phase Flow In order to illustrate the use of the SETS method, its application to onedimensional single-phase flow in a pipe will be considered, following Mahaffy. The equations he used for mass energy and momentum conservation respectively in single-phase flow were:

∂ρ + ∇⋅ρV = 0 ∂t ∂ρe + ∇ ⋅ ρ eV = − p∇ ⋅ V + h(Tw − T ) ∂t 1 ∂V + V ⋅ ∇V = − ∇ p − K V V ρ ∂t (2)

In these equations, K is a wall friction factor, h is a heat transfer coefficient for the heat-transfer area, and Tw is the temperature of the pipe wall. The pipeline is divided into discrete cells for the finite difference solution, and a staggered mesh is used for the discretized flow variables - that is, the velocities are defined on the cell edges and the bulk properties, such as density, energy, etc., are defined in the cell centers. In order to define values for cell-centered properties at the cell edges, the SETS method uses donor cell weighting. For any group of state variables Y, defined at the cell centers:

〈YV 〉 j + 1 = Y jV j + 1 , if V j + 1 ≥ 0 2

2

2

= Y j +1V j + 1 , if V j + 1 < 0 2

2

(3)

Essentially, this ensures that the value used for any bulk fluid property is that from the cell which is ‘providing’ the material: if the velocity on the right hand cell edge is positive, then the flow is out of this cell in the direction of increasing j, and the current cell's fluid properties are used. If the velocity is negative, then fluid is flowing in from the next `higher' cell (j+1), and the density, energy etc. of that cell's fluid are used instead. This helps stability and is also more physically realistic. The one-dimensional equations, when discretized by the SETS method, are as follows:


7

Stabilizer Momentum Equation

~ n +11 − n ) (V V j+ 1 j+ 2

+

2

∆t 1 ρ

~ ~ ~ + V jn+ 1 ∇ j + 1 V n +1 + ξ(V jn++11 − V jn+ 1 )∇ j + 1 V n 2

2

2

(

2

2

)

~ ( p nj+ 1 − p nj ) + K nj+ 1 2 V jn++11 − V jn+ 1 V jn+ 1 = 0

n +1

∆x j + 1 j+ 1

2

2

2

2

2

2

2

(4)

Basic Mass Equation

(~ ρ nj +1 − ρ nj ) ∆t

+ ∇ j ⋅ (ρ n +1 V n +1 ) = 0 (5)

Basic Momentum Equation

~ (V jn++11 − V jn+ 1 ) 2

2

∆t

+

2

1 ρ

n +1 j + 12

~ ~ + V jn+ 1 ∇ j + 1 V n +1 + ξ(V jn++11 − V jn+ 1 )∇ j + 1 V n

∆x j + 1

2

2

2

2

(~ p nj++11 − ~ p nj +1 ) + K nj+ 1 (2 V jn++11 − V jn+ 1 ) V jn+ 1 = 0 2

2

2

2

2

(6)

Basic Energy Equation

e jn +1 − ρ nj e nj ) (~ ρ nj +1 ~

+ ∇ j ⋅ (ρ n +1 e n +1 V n +1 ) + ~ p nj +1∇ j ⋅ (V n +1 ) ∆t ~ − h jn (Twn, j − T jn +1 ) = 0 (7)

Stabilizer Mass Equation

(ρ nj +1 + ρ nj ) ∆t

+ ∇ j ⋅ (ρ n +1V n +1 ) = 0 (8)


8

Stabilizer Energy Equation

(ρ nj +1 e nj +1 − ρ nj e nj )

+ ∇ j ⋅ (ρ n +1 e n +1 V n +1 ) + ~ p nj +1∇ j ⋅ (V n +1 ) ∆t ~ − h jn (Twn, j − T jn +1 ) = 0 (9)

The operators in these equations are defined, using the donor-cell notation, as follows:

∇ j ⋅ (YV ) =

( A j + 1 YV

j + 12

2

− A j − 1 YV 2

j − 12

)

vol j (10)

and

V j + 1 ∇ j + 1V = 2

V j + 1 (V j + 1 − V j − 1 ) 2

2

2

∆x j + 1

2

, if V j + 1 ≥ 0 2

2

=

V j + 1 (V j + 3 − V j + 1 ) 2

2

2

∆x j + 1

, if V j + 1 < 0 2

2

(11)

Also

~ ξ = 0 if ∇ j + 1 V n < 0 2 ~ ξ = 1 if ∇ j + 1 V n > 0 2

(12)

These equations may be solved in any order, although Mahaffy found superior stability properties if they are solved in the above order, i.e. the stabilizer momentum equations are solved first, followed by a set of basic equations and finally the stabilizer mass and energy equations. The stabilizer momentum equation is tridiagonal in V n +1 and is solved by standard methods. The basic equations are a nonlinear, coupled set of equations which are solved by an iterative method, as described below. Finally, once V n +1 is known, the stabilizer mass and energy systems can be solved as they are reduced to tridiagonal systems in the variables ρ (mass equation) and e ρ (energy equation).


9

Linearisation of the Finite Difference Equations The basic set of finite difference equations is non-linear and hence to solve them at each timestep it is necessary to use an iterative method, based on Newton’s method. Starting with some estimated values for the independent variables at the new timestep, the derivatives of the equations with respect to those variables are used to give the next best estimates - based on linear extrapolation from the last value - continuing until the latest estimates are equal (within prescribed tolerance) to the previous ones. The following is an extension of the summary given in Appendix C of Liles et al (1984). At a timestep n+1, given an initial guess of the independent variables (in this case p and T), the values on the next iteration (variables without primes) are assumed to be related to those at the last (with primes) by the relations:

p nj +1 = p 'nj +1 + δ p (13)

T jn +1 = T 'nj +1 + δ T (14)

Since the finite-difference equations are functions of T and p, a Taylor expansion about the last iteration’s value, retaining only the terms linear in dpj and dTj gives:

f ( p, T ) = f ( p '+ δ p, T '+ δ T ) = f ( p' ,T ' ) +

∂f ∂p

δp+ p= p'

∂f δ T + higher order terms ∂p T =T ' (15)

In practice, an expansion is performed only on the mass and energy equations, and the momentum equation is treated differently. First, it is rearranged to yield as a function of the pressures in the current and adjacent cells. This gives:

V jn+ 1 2 V jn++11 = 2

  ~ p nj++11 − ~ p nj  ~ n +1 ~n n n n n  − ∆t V j + 1 ∇ j + 1 (V − ξ V ) − K j + 1 V j + 1 V j + 1 + n 2 2 2 2  2 ρ j + 1 ∆x j + 1  2  2  ~n n n 1 + ∆t (2 K j + 1 V j + 1 + ξ ∇ j + 1 V j + 1 ) 2

2

2

2

(16)


10

From this equation it can be seen that:

∂V jn++11 2

∂p nj++11

=−

∂V jn++11 2

∂p nj +1 (17)

Also from Equation (16), by replacing the index j with j - 1:

∂V jn−+11 2

∂p nj−+11

=−

∂V jn−+11 2

∂p nj +1 (18)

Performing the differentiations indicated in Equation (15), using the mass equation as an example, leads to:

~ ~ ∆t f ( ~ p , T ) = ∆t f ( ~ p',T ' )  ∂ ρ' n +1 ∂V jn++11 ∂V jn−+11  n n 2 2 + δ pj + ∆t A j + 1 ρ j + 1 − ∆t A j − 1 ρ j − 1 n +1 n +1  2 2 2 2 ∂ p ∂p j ∂p j   j   ∂V jn++11 ∂V jn−+11  n n +1 2 2 + δ p j +1 ∆t A j + 1 ρ j + 1 − ∆t A j − 1 ρ j − 1 n +1 n +1  2 2 2 2 ∂ ∂ p p j +1 j +1     ∂V jn++11 ∂V jn−+11  n n 2 2 + δ p j −1 ∆t A j + 1 ρ j + 1 − ∆t A j − 1 ρ j − 1 n +1 n +1  2 2 2 2 ∂ ∂ p p j −1 j −1     ∂ ρ' n +1  + δT    ∂T j  (19)

The ρ nj + 1 and ρ nj − 1 terms are the donor-cell weighted averages of the density 2

2

at the right and left cell boundaries, respectively. Derivatives of velocity with respect to pressure occur in these equations, and Equations (17) and Equation (18) above may now be used to eliminate the

∂V jn++11 2

∂p nj++11

and

∂V jn−+11 2

∂p nj−+11

terms.


11

From Equation (16), the following derivatives are zero:

∂V jn++11 ∂V jn−+11 2 2 , ∂p nj−+11 ∂p nj++11 (20)

Rearranging Equation (19) gives:

 A j + 1 ρ nj + 1 ∆t ∂V jn++11  2 2 2  (δ p j +1 − δ p j ) δ pj −  n +1 vol ∂ p   j j j j   n n +1  A j − 1 ρ j − 1 ∆t ∂V j − 1  2 2 2  (δ p j − δ p j −1 ) − n +1 vol ∂ p   j j   ∂ ρ' ∂T

n +1

δT +

∂ ρ' ∂p

n +1

= ρ nj − ρ'nj +1 −∆t∇ j ⋅ (ρ n V n +1 ) (21)

A similar procedure is applied to the energy equation, resulting in:

n +1  n +1 ∂e' n +1  ∂ ρ' + e'nj +1 + ∆t ⋅ h jn +1  δ T ρ' j ∂T j ∂T j   n +1 n +1   ∂e' ∂ ρ' + ρ'nj +1 + e'nj +1 + ∆t∇ ⋅ V 'nj +1  δ p j ∂p j ∂Tp j   n +1  ρ nj +1   A j − 12 ∆t ∂V j + 12  n +1  (δ p j +1 − δ p j ) − ρ' j + n   2e j + 1   vol j ∂p nj +1   2    + 1 n  ρ nj −1   A j − 12 ∆t ∂V j − 12   (δ p j − δ p j −1 ) − ρ'nj +1 + n   2e j − 1   vol j ∂p nj +1   2   

= ρ nj e nj − ρ'nj +1 e'nj +1

[

− ∆t ∇ j ⋅ ρ n e nV 'n +1 −h jn +1 (Twn, j − T 'nj +1 )

] (22)

These two equations can be arranged into the form:

δ p j  Bj  = b j + c j (δ p j +1 − δ p j ) − d j (δ p j − δ p j −1 )  δT j  (23)


12

Bj is a 2 x 2 matrix containing the coefficients of dpj and dTj from Equations (21) and Equations (22). If this equation is now multiplied by B-1, two equations result giving dpj and dTj in terms of the variations in pressure:

δ p j = b'1 +c'1 (δ p j +1 + δ p j ) − d '1 (δ p j − δ p j −1 ) (24)

and

δ T j = b'2 +c'2 (δ p j +1 − δ p j ) − d '2 (δ p j − δ p j −1 ) (25)

where:

 b'  b' =  1  = B −1b b'2  (26)

 c'  c' =  1  = B −1c c ' 2  (27)

and

d '  d' =  1  = B −1d d '2  (28)

The first of these equations has pressure only as an unknown, and is a tridiagonal system which can easily be solved for the dpjs. Once these are known, the dTjs can be calculated from the second equation. This gives the current iteration’s new p and T, and once these are known, the densities and energies can be found (from the thermodynamics) and the velocities calculated. Finally, the solution at this iteration is compared to the solution at the last to see whether it has converged within the prescribed tolerance. If not, the new p and T values from this iteration are taken to be the old ones for next iteration. This cycle continues until either a solution is reached or a maximum number of iterations have taken place, in which case no solution can be found.


13

Extension of SETS to Two-phase Flow It is relatively straightforward to extend the above analysis to two-phase flow, for which the governing equations, as implemented in the Aspen Dynamic Pipeline Solver, are: Gas Mass Conservation

∂ αgρg

+ ∇ ⋅ (α ρ g Vg ) = Γ

∂t

(29)

Liquid Mass Conservation

∂ (1 − α) ρ l + ∇ ⋅ ((1 − α) ρ l Vl ) = −Γ ∂t (30)

Gas Momentum Conservation

∂Vg ∂t −

+ V g ⋅ ∇V g = −

c 1 ∇p − i (Vg − Vl ) Vg − Vl ρg αρ g

cwg Γ+ Vg Vg + g sin β (Vg − Vl ) − αρ g αρ g (31)

Liquid Momentum Conservation

ci ∂Vl 1 + Vl ⋅ ∇Vl = − ∇p + (Vg − Vl ) Vg − Vl ∂t ρl (1 − α) ρ l −

cw l Γ− (Vg − Vl ) − Vl Vl + g sin β (1 − α) ρ l (1 − α) ρ l (32)

Total Energy Conservation

[

∂ (1 − α) ρ l el + αρ g e g

] + ∇ ⋅ [(1 − α) ρ

∂t = p∇ ⋅ (1 − α)Vl + α V g + q wl + q wg

[

]

l

el Vl + αρ g e g V g

] (33)


14

Gas Energy Conservation

∂ (α ρ g eg )

+ ∇ ⋅ (α ρ g e gVg ) = − p

∂t + qwg + qig + Γhsg

∂α − p∇ ⋅ (α Vg ) ∂t (34)

A total energy equation is used instead of a liquid energy equation because this gives the opportunity to force thermal equilibrium by using one equation rather than two; this implementation is easier if a total energy equation is included from the start. In the case of more complex piping, terms must also be added to deal with tee junctions. Mass, momentum and energy source terms need to be added to these equations at the tee. Another modification found to be beneficial to stability is the addition of a pair of explicit momentum equations (one per phase) to predict the initial new timestep velocity on the first iteration when solving the basic equation set. The SETS equations for this model, are as implemented in the Aspen Dynamic Pipeline Solver, but assuming no tee components, are given in Sets equations for two-phase flow. For each timestep: First, the stabilizer Momentum equations are solved, yielding values for V jn +1 . The main iteration begins, solving the basic equations. An outer loop controls the iterations as follows: 1

On the first iteration only, the explicit predictor momentum equations are solved to provide a good initial estimate for the new timestep velocities.

2

Then the whole set of basic equations are solved. As for the single-phase case, an iterative solution must be used, linearising the equations at each stage.

Finally, the stabilizer mass and energy equations are solved, giving the final solution for this timestep.

Physical Properties Required Properties Certain physical properties are required. These include: •

Vapor mass fraction

•

Liquid density

•

Vapor density

•

Liquid viscosity

•

Vapor viscosity

•

Liquid internal energy

•

Vapor internal energy

•

Surface tension


15

The following derivatives are also required:

 ∂ ρl   ∂p

  ∂ ρl   ,   T  ∂T  p

 ∂ρg   ∂p

  ∂ ρg  ,  T  ∂T

  p

 ∂el   ∂el    ,    ∂p T  ∂T  p  ∂e g   ∂p

  ∂eg  ,  T  ∂T

  p (35)

These derivatives are required at constant composition of the appropriate phase for correct application in the Aspen Dynamic Pipeline Solver. The derivatives are, therefore, different to those at constant overall composition.

Interpolation of Property Data Physical properties at a given pressure and temperature are interpolated. In previous distributed versions of the software, two interpolation methods were provided, linear and bicubic. In the Aspen Dynamic Pipeline Solver only linear interpolation of property data is available. Previous experience has shown that the bicubic interpolation did not greatly enhance the performance of the software when compared to the linear interpolation.

Linear Interpolation The pressure, p, lies between p1 and p2; the temperature T, between T1 and T2. A physical property φ for each of the four points surrounding (p, T) is:

φ11 = φ ( p1 , T1 ) φ12 = φ ( p1 , T2 ) φ 21 = φ ( p2 , T1 ) φ 22 = φ ( p2 , T2 ) (36)

By interpolating in the pressure and temperature planes, we obtain:

 p − p1   T − T1  (φ 21 − φ11 ) +  (φ12 − φ11 )  p2 − p1   T2 − T1 

φ ( p, T ) = φ11 + 

 T − T1  p − p1   (φ 22 + φ11 − φ12 − φ 21 ) +   T2 − T1  p2 − p1  (37)


16

Closure Laws and Models The Solution Procedure section show in detail the six equation set used in the Aspen Dynamic Pipeline Solver. In order to close this set of equations and provide all of the necessary information for them to be solved, a set of additional models are required. These closure models are required for wall and interfacial friction, wall and interfacial heat transfer and interfacial mass transfer. In addition to these models the software also includes models for choking flow and for pigging of a pipeline. These are also discussed in this section.

Flow Regimes The Aspen Dynamic Pipeline Solverhas different flow regime maps for vertical and horizontal flow. These are needed because the friction behavior for two-phase flow is dependent on the orientation of the flow. The flow regime is predicted for each cell edge based on the velocities of each phase, void fraction and other parameters. The vertical flow pattern map is used if the angle of inclination is above 10° and the horizontal flow pattern map at inclinations below 10°.

Vertical Flow Regimes The flow regime boundaries in vertical flow are mainly based on void fraction. There is, however, an additional transition from slug flow to bubble flow for a mass flux above 2000 kg/m2s. This flow regime map was developed by the authors of the TRAC code. They based this on physical intuition and it has been found over many years of use in TRAC, PLAC and ProFES Transient, to work well for a wide variety of transients.


17

Horizontal Flow Regimes The horizontal flow map was similarly developed by the TRAC developers. This includes a stratified flow model, the basis for which is a model by Taitel and Dukler (1976) based on a modified Kelvin-Helmholz model for a circular pipe.

The model defines a critical gas velocity, Ucrit, above which stratification is impossible.

U crit

 h  (ρ1 − ρ g ) g cos β Ag = 1 − l  ρ g dAl / dhl  D 

1

2    (38)

where: 1 dAl = ( D 2 − (2hl − D) 2 ) 2 dhl

(39)

β is the angle formed by the pipe axis with the horizontal and hl is the liquid height. Stratified flow is only possible if hl > D/1000. For gas velocities between Ucrit and 2Ucrit the flow regime is assumed to be in transition between stratified and other horizontal flow patterns. The transitions between slug flow and annular flow are based simply on void fraction. A transition region is assumed between these two regimes.


18

Interfacial Friction Models for interfacial friction are used for each of the vertical and horizontal flow regimes. There is a relaxation between the previous timestep friction factor and the new friction factor for stability. The transitions between flow regimes are also relaxed.

Interfacial Friction for Vertical Flow Bubble Flow The interfacial friction is calculated by predicting the bubble size and shear coefficient. The bubble diameter is calculated as follows:

Db =

Webσ ρ l Vr2 (40)

where Db is the bubble diameter and Web is the Weber number, assumed to be 7.5 from Crowley et al (1977). The bubble diameter must lie between the cell hydraulic diameter and 10-4 meters. The interfacial shear coefficient is provided by a standard set of formulae for a sphere, from Govier and Aziz (1972):

ci =

3cb α ρ l 4 Db (41)

where the shear coefficient is based on the bubble Reynolds number as follows: Reb

Cb

< 0.1

180

0.1 < = Reb <= 989

(18/Reb) (1 + 0.15 Re0.687)

> 989

0.33

Plug Flow The plug flow regime is treated in the same manner as bubble flow. The plug diameter is calculated by interpolation between the bubble size calculated above and the pipe diameter, based on the void fraction. A similar interpolation is done between 2000 and 2700 kg/m2s.


19

Annular Flow The entrainment fraction is calculated based on Kataoka and Ishii (1982) as follows:

[

E = 1 − exp − 0.23(Vg − VE )

] (42)

where:

1

 (ρ l − ρ g )σWed  4 VE = 2.33  ρ l2   (43)

The remainder of the liquid is in a film or sheet. The interfacial shear is a volume average of the film and droplet relations in the annular-mist regime. The wetted surface area of the cell is determined from the portion of the geometric flow area that is blocked. The total interfacial surface area is determined by the sum of the areas contained in the wetted film and droplets. A critical Weber number, equal to 4 for the drops, is used with a calculation procedure similar to that for bubbly flow. This value of the Weber number is appropriate for accelerating drops. The interfacial drag coefficient for the annular-droplet regime combines the droplet drag (see above) and the Wallis (1969), correlation for annular flow:

Ci =

0.01ρ g Dh

(1.0 + 75.0(1 − α)(1 − E )) (44)

where E is given by Equation (42). To avoid a singularity in the liquid acceleration, a void fraction, αcut, is calculated that corresponds to the minimum allowed film thickness. Above this void fraction value, the above equation is multiplied by (1 - α)/(1 - αcut). To obtain the interfacial drag coefficient, droplet drag is weighted by the liquid fraction that is entrained, and Wallis’ annular flow is weighted by the fraction remaining as a film.


20

Churn Flow For the regime between the bubbly/slug flow and annular flow, a cubic spline interpolation in the void fraction is made between the conditions that would exist if the void fraction were 0.75, in the annular or annular-mist regime, and the conditions that would exist if the void fraction were 0.5, in the bubbly/slug regime. If 0.5 < a < 0.75, then a weight factor, W, is calculated from:

W = (4 α − 2) 2 × (7 − 8 α) (45)

This interpolation assures that the correlation for the interfacial friction is a continuous function of the void fraction, the relative velocity, the mass flux and the various fluid thermodynamic and transport properties.

Interfacial Friction for Horizontal Flow Interfacial friction is calculated for each flow regime and relaxed over the timestep. Smoothing occurs between slug flow and annular flow, based on void fraction, for 0.4 < a < 0.6.

W = (5 α − 2) 2 × (7 − 10 α) (46)

and between stratified flow and other horizontal flow regimes, based on the ratio of gas velocity to the critical velocity from Equation (38), for Ucrit < Vg < 2Ucrit :

 Vg  W =  3 − 2 2 −  U crit  

V     2 − g  U crit  

  

2

(47)

Stratified Flow The interfacial friction in stratified flow is assumed to be a multiple of the gas wall friction factor. This multiple is obtained using the Sinai (1983) model or alternatively may be set to a constant value by the user in the input to the Solver. The Sinai model is summarized by the following equations. The interfacial friction factor is calculated from the ratio of the friction velocity to the gas velocity, thus:

U f i = 2 τ Ug 

   

2

(48)


21

where:

 Dg   + 4.73 = 5.75 log10  Uτ  2 εi 

Ug

(49)

and

ε i = 180

ρg S i U τ2 S g + S i g (ρ l − ρ g ) (50)

The interfacial friction factor is limited to a maximum of 6 times the gas phase wall friction factor. This is to give realistic values for interfacial friction factor in high pressure, large diameter pipelines, based on results of studies by Kawaji et al (1987), Oliemans (1987), Crowley and Rothe (1988) and Spedding and Hand (1990). The drag coefficient is then given by:

ci = 0.5 f i ρ g

Si FA (51)

where Si is the width of the interface.

Slug Flow The horizontal slug flow model is identical to the vertical plug flow model, except that the slug size is assumed to equal the pipe diameter.

Annular Flow The droplet entrainment fraction, E, is calculated in the same manner as for vertical annular flow. The drag coefficient is based on Whalley (1987). The droplet core void fraction is calculated from:

 V  α c = (1 − α)1 − l E   Vg    (52)

The film radius is calculated from:

r film = 0.5 D α c (53)


22

and the film thickness from:

d f = 0.5D − r film (54)

The ratio of interfacial friction factor to gas-wall friction factor is then calculated from:

fi f wg

= 1+ 360

df D (55)

The interfacial drag coefficient is then given by:

ci = 0.5(1 − α) f i ρ g

P FA (56)

where the length of the interface is:

P = 2πr film (57)

Wall Friction The total pressure gradient calculated in the momentum equations is expressed as the sum of the fractional dissipation, acceleration head and potential head terms. The Aspen Dynamic Pipeline Solver calculates coefficients for the frictional dissipation terms and for losses associated with abrupt area changes. Under single-phase flow conditions, pressure drops associated with frictional losses are correlated as functions of fluid velocity, fluid density, fluid viscosity, channel hydraulic diameter and surface roughness of the channel wall. When a two-phase mixture is flowing in a channel, a correction to the single-phase frictional loss is necessary to account for added dissipation between phases and interactions with the channel walls. This correction factor is the two-phase flow multiplier. The wall shear coefficients cwg and cwl are defined as:

cwg = αρ g

c fg Dh (58)

and


23

c wl = (1 − α) ρ l

c fl Dh (59)

where cfg and cfl are the gas and liquid friction factors. The options available to calculate the wall friction are: 1

Constant friction factor (user input).

2

Homogeneous model for smooth pipe walls.

3

Homogeneous model for rough pipe walls.

4

Annular flow model.

5

Smooth + form loss.

6

Rough + form loss.

7

Rough + annular flow model + form loss.

8

Form loss only.

The first option allows you to specify a two-phase friction factor in order to model pressure drop data. Since this is a two-phase friction factor, a value of 0.01 generally gives similar results to using the homogeneous friction factor (for smooth pipes). Roughness for the third option is specified as a relative roughness (absolute roughness height/pipe internal diameter). For the fourth option then, if annular flow is detected, the homogeneous friction factor (for rough pipe) will be replaced by one calculated from an annular flow model. Options 5, 6 and 7 are the same as 1, 2 and 3, with the addition of an automatic calculation of an appropriate form loss coefficient, if there are abrupt area changes. The final option just calculates form loss coefficients.

Homogeneous Model The homogeneous friction factor model alters the single phase value by using a two-phase viscosity defined in terms of the flow quality (x) (Collier (1972):

1 x (1 − x) = + µ µg µ1 (60)

The homogeneous friction factor (Rohsenow and Choi (1971), is then given by:

Re ≤ 500, f = 0.032 (61) −6

500 ≤ Re ≤ 5000, f = 0.032 − 5.25 × 10 ( Re − 500) (62)


24

Re ≤ 5000, f =

0.0625   15 ε / D   +  log10     Re 3.715   

2

(63)

where: Re = GDh /µ and G = ρ m Vm . The three equations above (Equations (61), Equations (62), Equations (63)) represent a constant friction factor for Re <= 500, a linear interpolation between 500 < Re < 5000, to the value given by Equation (63). The coefficient of friction for the liquid is then given by:

c fl = 2 fΦ lo2 (64)

where, for α < 0.9:

2  µ1

 (1 − x)   µ Φ lo2 = (1 − α)

−0.2

(65)

and, for α > 0.9:

Φ lo2

 (1 − x) 2 (10 α − 9) 2 (21 − 20 α) x 2 ρ l  µ l = +  αρ g  (1 − α)  µ

  

−0.2

(66)

If the void fraction is greater than 0.9, the coefficient of friction for the gas phase is:

c fg = (10 α − 9) 2 (21 − 20 α)c fl (67)

The purpose of this function is to ensure a smooth transition from zero gas wall friction at α = 0.9 to the single phase gas value at α = 1.0.


25

Stratified Flow Model If the flow regime is stratified, then the homogeneous friction factors are replaced by those for stratified flow, which are based on the wall area with which each phase is in contact. Given the void fraction, the liquid height is determined from the following relationship:

α=

θ 1  2hl  + 1 −  sinθ D  π π (68)

together with: −1

 2h  θ = cos 1 − l  D   (69)

Then the following geometric parameters are evaluated:

S g = (π − θ) D Sl = θ D S i = D sinθ π 2 D 4 Ag = A α Al = A(1 − α) A=

Dhg =

4 Ag S g + Si

Dhl =

4 Al Sl + Si (70)


26

The Reynolds numbers for the gas and liquid phases are then calculated as:

Reg =

ρ g V g Dg µg (71)

and

Rel =

ρ l Vl Dl µl (72)

Note: The hydraulic diameter for the liquid phase is calculated using the sum of the liquid and interfacial lengths (Sl+ Si). This has the effect of decreasing the hydraulic diameter compared to the usual definition involving only Si and thereby increasing the friction factor. The physical justification for such a change of definition is that the interface is usually rough and the interfacial shear stress is comparable to the liquid wall shear stress. Thus the interface acts more as a rough wall than a free surface. This modification has been tested against experimental data and shown to give better reproduction of the data.

The friction factors fwg and fwl are then calculated using the relationship cited under the homogeneous friction factor section:

f =

0.0625   15 ε / D   +  log10     Re 3.715   

2

(73)

Drag coefficients are then calculated as:

c fg =

Dh S g DS 1 1 f wg ; c fl = f wl h l Ag Al 2 2 (74)

Form Losses The finite-difference equations yield the correct pressure loss for an abrupt expansion. However, this is not true for an abrupt contraction or an orifice. For one-dimensional components, a form loss option (see above) should be specified for the input friction at the location of any abrupt area change. An


27

appropriate loss correction is calculated by including an extra term in the Bernoulli equation of the form:

∆p =

k ρV 2 2 (75)

where k is a form loss coefficient. For an abrupt expansion or zero length orifice:

 A  k = 1 − 1  A2  

2

(76)

and for an abrupt contraction:

A  A  k = 0.5 − 0.7 1  + 0.2 1   A2   A2 

2

(77)

where A1 and A2 are the smaller and the larger flow areas, respectively. The above equation is a curve that was fitted to the values reported in Massey (1968).

Annular Flow If the user selects the annular flow option for wall friction (see above) then, if the flow regime is annular, the homogeneous wall friction factors are replaced by wall friction factors calculated on the basis of an annular flow model. Since, in annular flow, the gas phase flows in the core surrounded by a liquid film on the walls, there is no gas in contact with the pipe wall and hence the gas-wall friction factor is set to zero. The liquid-wall friction factor is calculated as follows:

Re ≤ 200, f = 16 /( Re) (78)

200 ≤ Re ≤ 8000, f = 0.001069 + 143.38 (ln Re) −4.5 (79)

Re ≤ 8000, f =

0.0625   15 ε / D   +  log10     Re 3.715   

2

(80)


28

where the thickness of the liquid annulus is taken as the characteristic diameter in the Reynolds number.

Heat Transfer Coefficients The heat fluxes per unit volume required for the basic equations may be obtained from the heat transfer coefficients:

qi = hi Ai

(Tg − Tl vol (81)

qwg = hwg Awg

(Tw − Tg ) vol (82)

and

qwl = hwl Awl

(Tw − Tl ) vol (83)

Simple correlations are used to predict the heat transfer coefficients as discussed below.

Interfacial Heat Transfer Currently, if the default 2 energy equations are selected, the interfacial heat transfer coefficient, hi, is fixed such that:

hi Ai = 10 7 W / K (84)

with the result that the liquid and gas temperatures are almost always equal. The user can force the liquid and gas temperatures to vary independently by selecting the option to have separate energy equations with no interfacial heat transfer.

Fluids to Wall Heat Transfer The heat transfer coefficients between the liquid or gas phases and the pipe wall are calculated according to flow pattern. The determination of flow patterns and the interpolation between them is identical to that used for calculation of friction factors. Values of the product of heat transfer coefficient and contact area are determined by flow regime as follows.


29

Bubbly Flow There is assumed to be contact only between the liquid and the wall. For calculation of the liquid to wall heat transfer coefficient, the velocity is the mean fluid velocity, the viscosity is the homogeneous viscosity and the diameter is the hydraulic diameter of the channel. The gas to wall heat transfer coefficient is taken to be zero.

Stratified Flow The hydraulic diameters appropriate to the gas and liquid layers in stratified flow are calculated as discussed earlier for calculation of stratified flow wall friction factors. The liquid and gas to wall heat transfer coefficients are then calculated separately using the velocities and physical properties for the gas and liquid phases. These heat transfer coefficients are then weighted by the wetted perimeters of the two phases.

Slug Flow Here, the liquid velocity is taken as the slug velocity and the diameter is the hydraulic diameter of the pipe. The liquid and gas to wall heat transfer coefficients are then calculated separately and weighted by the volume fraction of each phase in the pipe.

Annular Flow In annular flow, it is assumed that there is heat transfer only between the liquid film and the wall. Hence the liquid film velocity and film thickness are used in the determination of liquid to wall heat transfer coefficient. The film is assumed to be uniform circumferentially.

Single Phase For single-phase gas or liquid flow, the heat transfer coefficient is calculated based on velocity and physical properties of the phase and the channel hydraulic diameter. For void fraction greater than 0.98, single phase gas is assumed for calculation of heat transfer coefficient.

Gas to Wall Heat Transfer The gas to wall heat transfer coefficient is taken as the maximum of the turbulent natural convection equation (McAdams (1954):

 ρ 2g g Tw − Tg hnc = 0.13k g   µ 2 Tg 

1

3 1  Pr 3 g   (85)


30

and the turbulent Dittus-Boelter equation:

hturb = 0.023

kg Dg

Reg0.8 Prgn (86)

where the gas Reynolds number is:

Reg =

ρ g V g Dg µg (87)

and the gas Prandtl number is:

 µ cp   Prg =   k g (88)

The Prandtl number exponent, n, is 0.3 for cooling (heat transfer from fluids to pipe wall) and 0.4 for heating (heat transfer from pipe wall to fluids).

Liquid to Wall Heat Transfer The liquid-wall heat transfer coefficient is taken as the maximum of the laminar heat transfer coefficient:

hlam = 4.0

kl Dl (89)

and the Dittus-Boelter equation for turbulent flow:

hl = 0.023

kl Rel0.8 Prl n Dl (90)

where the liquid Reynolds number is:

Rel =

ρ l Vl Dl µl (91)


31

and the liquid Prandtl number is:

 µ cp   Prl =   k l (92)

The Prandtl number exponent, n, is 0.3 for cooling (heat transfer from fluids to pipe wall) and 0.4 for heating (heat transfer from pipe wall to fluids).

Enhancement Due to Roughness The above calculations for liquid or gas to wall heat transfer coefficients are based on flow in a smooth pipe. For turbulent flow in a rough pipe, the heat transfer coefficient is enhanced in the same way as the friction factor. The friction factor is calculated using the relationship:

f =

0.0625   15 ε / D   +  log10     Re 3.715   

2

(93)

and hence the heat transfer coefficient enhancement factor, from the ratio of friction factor for a rough pipe to that for a smooth pipe, is given by:

   15  log10      Re   E=  15 ε / D     log10  Re + 3.715    

2

(94)


32

Cylindrical Wall Heat Conduction The Aspen Dynamic Pipeline Solver solves the radial conduction equation (below) through the pipe wall to update the wall temperatures radially:

ρ cp

∂T 1  ∂  ∂T  =  rk  ∂t r  ∂r  ∂r  (95)

The wall temperatures are updated after all the fluid mechanics calculations are completed for that timestep. The wall is split into elemental volumes, as shown below:

The conduction equation is re-cast in finite difference form, with boundary conditions applied at the inner and outer walls (i=1,N) such that:

−k

∂T ∂r

= hl (Tl − Twl ) + hg (Tg − Twl ) i =1

(96)

and

−k

∂T ∂r

= − hout (Tout − TwN ) i=N

(97)

The heat transfer coefficient from the last node to the surrounding heat sink, hout, and the ambient temperature, Tout, must be specified by the user. If warm up of the pipe wall is not important, for example if the material has a low thermal inertia (such as steel) or is uninsulated and exposed, the number of heat transfer nodes can be set to 1. In this case, the heat conduction is


33

based on the fluid temperature, the external temperature and the overall heat transfer coefficient, hout. For thermal transient calculations, such as pipeline warm up, the surrounding media (such as the soil) should be simulated as additional layers of insulation material. The number of nodes is one more than the number of material layers. If correct warm up of a single material layer is required, nodes = 2 (for example). Where thick insulation layers are present, these need to be subdivided in order to accurately represent the heat transfer through them. In previous versions, this had to be done by the user, it is now automatic. The method used is as follows: Each sub-layer has an inner radius, Ri, and outer radius, Si, where i is the number of the sub-layer. The ratio of the sub-layer thickness dRi is a maximum fraction x of the inner radius, Ri. Layer

R1

dR1

S1

1

R

xR

(R+xR=(1+x)R

2

(1+x)R

x(1+x)R

i-1

i-1

I

(1+x) R

(1+x)R+x(1+x)R=(1+x) 2 R (1+x)iR

x(1+x) R

n

n

The total thickness of n sub-layers is (1+x) R-R = ((1+x) -1)R. Thus, for a given layer:

[

]

dR = (1 + x) n − 1 R (98)

Rearranging then leads to the following expression for the number of sublayers:

 dR  + 1 log R   n= log(1 + x) (99)

where must be rounded up to the next integer. The actual value of x to use for the local sub-layer is then given by:

dR   +1  log −1 R  −1 x = log  N       (100)

where N is the integer number of sub-layers.


34

Burial Calculations Another feature of the Aspen Dynamic Pipeline Solver, is the automatic calculation of the equivalent soil layer thickness for a buried pipe (BD), i.e. the distance from the top of the pipe (plus any insulation layers) to the soil surface. Define:

z = BD +

D 2 (101)

where D is the outer diameter of the pipe (including insulation). Making dimensionless:

z' =

z BD 1 = + D D 2 (102)

The dimensionless equivalent outer diameter is then given by:

D'o = exp(cosh −1 (2 z '−1)) : z ' <

3 2 (103)

and

D ' o = 4 z ': z ' ≥

3 2 (104)

and the thickness of the soil layer by:

∆soil =

D'o D − D  D'o −1  =  D. 2  2  (105)


35

Interfacial Mass Transfer Interfacial mass transfer occurs in multi-component hydrocarbon pipelines due to evaporation or condensation. This arises due to changing pressure and temperature along the pipe (convective phase change), or with time (transient phase change). The rate of evaporation of the liquid phase is given by:

 ∂ ζ  ∂p  ∂ ζ  ∂T   ∂ ζ  ∂p ∂x  ∂ ζ  ∂T ∂x   Γ =   + +    +    (m g + ml )  ∂p  T ∂t  ∂T  p ∂t   ∂p  T ∂x ∂t  ∂T  p ∂x ∂t   (106)

where: Γ = Mass of liquid evaporated in cell per unit time ζ = Equilibrium liquid mass fraction

m g , ml = Mass of gas, liquid in cell The equilibrium liquid fraction is a function of the pressure and temperature:

ζ = ζ( p, T ) (107)

and so the differential change in equilibrium liquid fraction is given by:

 ∂ζ   ∂ζ  d ζ =   dp +   dt p ∂  ∂T  p  T (108)

Hence the derivatives of with respect to time and to distance are:

∂ ζ  ∂ ζ  ∂p  ∂ ζ  ∂T + =   ∂t  ∂p T ∂t  ∂T  p ∂t (109)

and

∂ ζ  ∂ ζ  ∂p  ∂ ζ  ∂T + =   ∂x  ∂p T ∂x  ∂T  p ∂x (110)


36

Substituting Equations (109) and Equations (110) into Equation (106) gives:

 ∂ ζ ∂ ζ ∂x  Γ = − +  ( m g + ml )  ∂t ∂x ∂t  ∂ ζ ∂x  ∂ζ = −  ( m g + ml ) + (m g + ml ) ∂x ∂t   ∂t (111)

The first term on the right hand side of represents the transient contribution to the phase change and can be expressed in finite difference form as:

 ζ nj − ζ nj −1  ∂ζ (m g + ml ) = (m g + ml )   ∂t  ∆t  (112)

The second term on the right hand side of Equation (111) represents the convective contribution to the phase change. To express in finite difference form, it is first assumed that the rate of convection into and out of the cell is much greater than the rate of change of total mass in the cell, i.e.,

[

∂x ∂ (m g + ml ) ≈ x(m g + ml ) ∂t ∂t

] (113)

Hence the convective phase change component may be expressed as:

 ζ nj − ζ nj −1  ∂ ζ ∂x (mg + ml ) = (αα g Vg + (1 − x) ρl Vl ) A∆x   ∂x ∂t  ∆t  (114)

In the Aspen Dynamic Pipeline Solver the conservation equations are multiplied by the timestep, thus giving, for the phase change from Equations (111), Equations (112) and Equations (114):

Γ∆t = (m g + ml )(ζ nj−1 − ζ nj ) + A(α ρ gV g + (1 − α) ρ l Vl )∆t (ζ nj−1 − ζ nj ) (115)

In the above equations, the subscript j refers to the cell position, and the superscript n to the time. Equation (115) is implemented in the Aspen Dynamic Pipeline Solver, subject to the following: •

The transient and convective components are separated and may be suppressed if required for stability.


37

•

The convective component is aligned with the direction of flow. If the flow is from cell j+1 to cell j then the difference

is used.

•

The transient term is approximated using a backward-time difference, as shown in Equation (112). This may lead to numerical instability, and the transient component is therefore multiplied by an under-relaxation factor.

•

The under-relaxation factor also serves as a limited control over the rate of transient phase change during rapid depressurization. In these cases, the equilibrium liquid mass fraction can change significantly between timesteps, but thermodynamic equilibrium is not attained over the timestep. This rather crude method allows the simulation of blow down, but it is recognised that a more formal calculation method to allow for nonequilibrium effects is required. This work is under development.

The amount of phase change in a cell at the current timestep is modified by a factor dependent on the amount of liquid in the cell, PGRAT. If the liquid fraction is less than 0.9, then PGRAT=1.0; if the liquid fraction is greater than 0.95 then PGRAT=0.0; if the liquid fraction lies between these limits then PGRAT is linearly interpolated. This has the effect of smoothly reducing the amount of phase change to zero as a cell fills with liquid.

Choked Flow Model Critical flow can occur either at area restrictions in the pipe, such as valves, or at boundaries. The conditions in the fluid change very rapidly close to the choke plane, and hence it is not possible to predict critical flow directly. Special models have been written which calculate the velocities at the choke plane based on the conditions immediately upstream.

Single Phase Choking No model has been implemented for pure liquid flow, since this condition is unlikely to occur. The critical velocity for pure gas flow is calculated using the standard formula below:

q = γRTg (116)

where R is the gas constant, is the ratio of specific heats and q is the critical velocity; and R are calculated from the physical property data provided. If the velocity calculated at the new timestep is above the critical velocity, the new velocity is reduced to the critical velocity over several timesteps.


38

Two-phase Choking The two-phase choked flow model is based on the work by Henry and Fauske (1971). In this paper the critical mass flux for a general two-phase choked system is derived as: −1

dv  G2   = (1 + x(k − 1)) x g + (v g (1 + 2 x(k − 1) + kvl 2 −  dp  k 

( x − 1) + k (1 − 2 x))

dx dv + k (1 + x(k − 2) − x 2 (k − 1)) l dp dp

v g  dk  + x(1 − x) kvl −  k  dp  (117)

where k is the slip ratio, Vg/Vl ; x is the quality (gas mass fraction); vg and vl are the specific volumes of the gas and liquid (1/density).

Homogeneous Frozen Model The homogeneous frozen model, which you can select from the input interface, is derived from Equation (117) by making the assumptions: • Quality, χ, is constant at the choked plane. • The flow is homogeneous at the choked plane, i.e.,

νg = νl and k = 1.

The quality is evaluated from the void fraction, using:

x=

αρ g αρ g + (1 − α) ρ l (118)

Then the critical mass flux is calculated by:

G=

αρ g + (1 − α) ρ l  α   ρg 

 ∂ρg   ∂p

 (1 − α)  ∂ ρ l  +  ρ l  ∂p s

    s  (119)

Once the total mass flux across the choke plane increases above the critical value, the phase velocities are set to the homogeneous critical velocity over several timesteps. There is likely to be some slip between the phases as the


39

mass flux approaches the critical value. Hence, the gas velocity will initially increase beyond the critical velocity before it is reduced, and conversely for the liquid velocity. Care must be taken in evaluating the derivatives at constant entropy required in Equation (119). The derivatives available from the PVT package MULTIFLASH are at constant pressure or constant temperature. The required derivatives are given by the relationship: 2

 ∂ρ   ∂ρ  T  ∂ρ    =   − 2    ∂p  s  ∂p T ρ c p  ∂T  p (120)

and where the heat capacity at constant pressure, cp, is given by:

p  ∂ρ   ∂U  cp =    − 2 T ∂  p ρ  ∂T  p  (121)

Homogeneous Equilibrium Model In the homogeneous equilibrium model, which the user can also select from the input interface, the quality at the choke plane is assumed to be the thermodynamic equilibrium quality, which is related to the equilibrium liquid . The liquid and gas velocities are still assumed to be fraction by equal at the choked plane. This leads to the following expression for the critical mass flux:

G=

1  (1 − ζ)  ∂ ρ g  2   ρ g  ∂p

  1 1  α  ∂ρg  +  −      s  ρ g ρ l  ρ m  ∂p

 ζ  ∂ρ  + 2  l  s ρ l  ∂p

     s  (122)

After the choked mass flux is calculated, implementation is similar to the homogeneous frozen model.


40

General Information

This section provides Copyright details and lists any other documentation related to this release.

Copyright Version Number: HYSYS 3.4 March 2004 Copyright  1981 - 2004 Aspen Technology, Inc. All rights reserved. HYSYS®, HYSYS Optimizer™, ACM Model Export™, HYSYS Amines™, HYSYS Crude Module™, HYSYS Data Rec™, HYSYS DMC+ Link™, HYSYS Dynamics™, HYSYS Electrolytes™, HYSYS Lumper™, HYSYS Neural Net™, HYSYS Olga Transient™, HYSYS OLGAS 3-Phase™, HYSYS OLGAS™, HYSYS PIPESIM Link™, HYSYS Pipesim Net™, HYSYS PIPESYS™, HYSYS RTO™, HYSYS Sizing™, HYSYS Synetix Reactor Models™, HYSYS Tacite™, HYSYS Upstream™, HYSYS for Ammonia Plants™, MUSE™, PIPE™, Polymers Plus®, Process Manuals™, Process Tools™, ProFES 2P Tran™, ProFES 2P Wax™, ProFES 3P Tran™, ProFES Tranflo™, STX™, TASC-Thermal™, TASCMechanical™, Aspen Plus®, ACOL™, ACX™, APLE™, Aspen Adsim™, Aspen Aerotran™, Aspen CatRef®, Aspen Chromatography®, Aspen Custom Modeler®, Aspen Decision Analyzer™, Aspen Dynamics®, Aspen Enterprise Engineering™, Aspen FCC®, Aspen Hetran™, Aspen Hydrocracker®, Aspen Hydrotreater™, Aspen Icarus Process Evaluator™, Aspen Icarus Project Manager™, Aspen Kbase™, Aspen Plus® HTRI®, Aspen OLI™, Aspen OnLine®, Aspen PEP Process Library™, Aspen Plus BatchFrac™, Aspen Plus Optimizer™, Aspen Plus RateFrac™, Aspen Plus SPYRO®, Aspen Plus TSWEET®, Aspen Split™, Aspen WebModels™, Aspen Pinch®, Aspen Properties™, Aspen SEM™, Aspen Teams™, Aspen Utilities™, Aspen Water™, Aspen Zyqad™, COMThermo®, COMThermo TRC Database™, DISTIL™, DISTIL Complex Columns Module™, FIHR™, FLARENET™, FRAN™, HX-Net®, HX-Net Assisted Design Module™, Hyprotech Server™, the aspen leaf logo and Enterprise Optimization are trademarks or registered trademarks of Aspen Technology, Inc., Cambridge, MA. All other brand and product names are trademarks or registered trademarks of their respective companies.

General Information

41

This document is intended as a guide to using AspenTech's software. This documentation contains AspenTech proprietary and confidential information and may not be disclosed, used, or copied without the prior consent of AspenTech or as set forth in the applicable license agreement. Users are solely responsible for the proper use of the software and the application of the results obtained. Although AspenTech has tested the software and reviewed the documentation, the sole warranty for the software may be found in the applicable license agreement between AspenTech and the user. ASPENTECH MAKES NO WARRANTY OR REPRESENTATION, EITHER EXPRESSED OR IMPLIED, WITH RESPECT TO THIS DOCUMENTATION, ITS QUALITY, PERFORMANCE, MERCHANTABILITY, OR FITNESS FOR A PARTICULAR PURPOSE.

Corporate Aspen Technology, Inc. Ten Canal Park Cambridge, MA 02141-2201 USA Phone: (1) (617) 949-1000 Toll Free: (1) (888) 996-7001 Fax: (1) (617) 949-1030 URL: http://www.aspentech.com/

General Information

42

Technical Support

Online Technical Support Center AspenTech customers with a valid license and software maintenance agreement can register to access the Online Technical Support Center at: http://support.aspentech.com You use the Online Technical Support Center to: • Access current product documentation. • Search for technical tips, solutions, and frequently asked questions (FAQs). • Search for and download application examples. • Search for and download service packs and product updates. • Submit and track technical issues. • Search for and review known limitations. • Send suggestions. Registered users can also subscribe to our Technical Support e-Bulletins. These e-Bulletins proactively alert you to important technical support information such as: • Technical advisories. • Product updates. • Service Pack announcements. • Product release announcements.

Technical Support

43

Phone and E-mail Customer support is also available by phone, fax, and e-mail for customers who have a current support contract for their product(s). Toll-free charges are listed where available; otherwise local and international rates apply. For the most up-to-date phone listings, please see the Online Technical Support Center at: http://support.aspentech.com Support Centers

Operating Hours

North America

8:00 – 20:00 Eastern time

South America

9:00 – 17:00 Local time

Europe

8:30 – 18:00 Central European time

Asia and Pacific Region

9:00 – 17:30 Local time

Technical Support

44

Aspen Hydraulics Dynamics Reference

Recommend Documents