A model of incomplete phase transitions of gas hydrates in a porous medium

Journal of Applied Mathematics and Mechanics 75 (2011) 27–40

Contents lists available at ScienceDirect

Journal of Applied Mathematics and Mechanics journal homepage: www.elsevier.com/locate/jappmathmech

A model of incomplete phase transitions of gas hydrates in a porous medium夽 V.I. Kondaurov, A.V. Konyukhov Dolgoprudnyi, Russia

a r t i c l e

i n f o

Article history: Received 24 December 2009

a b s t r a c t The kinetics of incomplete phase transitions of gas hydrates in a porous medium is considered. In such transformations, the solid hydrate and its decomposition products coexist in extended regions. The conservation laws that take into account mass, momentum, and energy transfer between the components of the medium are formulated. The concept of phase transition dissipation is introduced. A general form of the constitutive relations, which is necessary and sufficient for the entropy inequality to be satisfied in any processes involving a change of state of the medium, is found. A potential of the skeleton that takes into account the surface energy, the latent energy of the phase transition and the temperature dependence is proposed. A thermodynamically consistent kinetic equation is formulated. The conditions under which the phase transition of a hydrate begins and is completed are found. The problem of isothermal dissociation of a hydrate when the pressure on part of the boundary of the body is reduced is examined. The influence of several parameters of the model on the seepage of the decomposition products is investigated. © 2011 Elsevier Ltd. All rights reserved.

1. Introduction Gas hydrates, which are solid compounds of gas and water molecules, are one of the promising sources of hydrocarbons. According to the estimates in Ref. 1, the world reserves of methane trapped in gas hydrates considerably exceed the known reserves of conventional natural gas. The bulk of all gas hydrates is trapped in porous subterranean sedimentary strata at a temperature of 7–12 ◦ C and a pressure of 107 Pa (100 atm). The behaviour of a porous medium saturated with gas hydrates has several special features. At a low temperature and a high pressure, the pore space is filled with hydrates in the solid state. Raising the temperature or reducing the compression level results in decomposition of the hydrate into the gas and water. When there is a pressure gradient in the pore channels, seepage of the gas and water occurs. At the present time, the most commonly used approach to describing gas hydrate dissociation is a generalization of the Stefan problem.2–5 A porous sedimentary stratum is divided into two parts by a dissociation front, on which complete decomposition of the hydrate occurs. The regions separated by the front have different permeabilities. A heat sink associated with the endothermic decomposition reaction acts on the front. The temperature and pore pressure are continuous on the front and satisfy the phase equilibrium condition, while the temperature gradient and the heat flux have a discontinuity on the front. A similar approach was also considered previously.6–9 The scheme works well when the permeabilities of the sedimentary stratum are low. At a high permeability, k ≥ 10−13 m2 , the removal of energy by the gas flow compensates the conductive supply of heat. As a result, the dissociation temperature and the velocity of the front decrease, and the temperature before the front exceeds the decomposition temperature, i.e., the hydrate becomes superheated. As was assumed in Ref. 10, this means that an extended phase transition region forms. In our opinion, a systematic description of the phase transition of a hydrate in a porous medium should differ from the classical models that rely on the Gibbs conditions11 for equilibrium between the phases of a homogeneous material, i.e., continuity of the pressure, temperature and chemical potentials of the phases on a macroscopic interfacial surface. The main reason is the difference between the capillary pressures in pore channels of different diameter. At equilibrium the narrow channels are filled with the hydrophilic liquid (water), and the wide channels are filled with the gas. The pressure of the water is lower than the pressure of the gas due to the surface tension. Therefore, when the Gibbs conditions are considered on the microscopic scale, they dictate that equilibrium is reached first in the narrow channels and, only after a further increase in temperature, in the wide channels. On the macroscopic level, where an element of a porous medium containing a large number of pore channels is considered, this leads to an incomplete, partial transition, whose extent is a function

夽 Prikl. Mat. Mekh. Vol. 75, No. 1, pp. 39–60, 2011. E-mail address: konyukhov [email protected] (A.V. Konyukhov). 0021-8928/$ – see front matter © 2011 Elsevier Ltd. All rights reserved. doi:10.1016/j.jappmathmech.2011.04.005

28

V.I. Kondaurov, A.V. Konyukhov / Journal of Applied Mathematics and Mechanics 75 (2011) 27–40

of the current thermodynamic state. The surface separating the phases of the material is replaced by an extended region, in which a gradual transition from the solid to a gas-liquid state occurs. Another special feature of the behaviour of a porous medium saturated with hydrates is the kinetics of the phase transition of the hydrates in it. This is indicated, for example, by the temperature difference between the forward and reverse transitions.12 One of the first kinetic models of gas hydrate dissociation was proposed by Kim et al.13 It is based on an equation that relates the hydrate decomposition rate to the difference between the gas volatilities in the current state and the equilibrium state. The concept of the gas volatility over a rigid pore surface is fairly obscure, and its experimental determination presents certain difficulties. In our opinion, the kinetics of the transition is related not to nucleation (the formation and growth of nuclei of the new phase), but to another factor, namely, the finite time for the establishment of capillary equilibrium.14,15 The redistribution of fluids in pore channels of different diameter in an elementary volume under the action of capillary forces is a considerably slower process than nucleation. The traditional kinetics of hydrate dissociation is only significant for fast mass-transfer processes. The fact is that the characteristic times of the decomposition of gas hydrates and of their seepage differ by an order of magnitude. Hence it follows that taking into account the kinetics of the transition results in significant differences only in a small vicinity of the melting curve, allowing the structure of the decomposition front to be revealed. States associated with the conditions of the Gibbs phase equilibrium adjoin this narrow region. For hydrates in a porous medium, where a partial, rather than a complete, transition occurs, the situation is different. Kinetics with a short transition time, during which a state close to the equilibrium state is produced, creates an extended region, in which the solid hydrate coexists with the gas and water, instead of an interfacial surface with a strong discontinuity. The dimensions of this region depend on the structure of the porous medium, i.e., the size distribution of the pores, as well as of the skeletal material, which determines the surface tension of the water and the gas. Below, the principle of thermodynamic consistency,16 which requires the entropy inequality to be satisfied for any history of variation of the parameters of state, is used to construct a model. This principle, which has been extended to the case of a porous medium,17 is used to find the general form of the constitutive equations of a porous medium, whose non-deformable skeleton is saturated with water and a gas and is capable of exchanging mass, momentum and energy with them. One of the important consequences of the general form of the constitutive equations is the representation of the total dissipation of the medium in the form of a sum of three components: the thermal dissipation due to the thermal conductivity of the medium, the seepage dissipation associated with the viscous friction of the water and the gas on the skeleton and the phase transition dissipation. The last component is equal to the product of the transition rate (the kinetic function) and the difference between the chemical potentials of the hydrate in the solid and liquid states. It is shown that, unlike the chemical potentials of the water and the gas, the chemical potential of the non-deformable skeleton is defined as the difference between its free energy and a quantity that is proportional to the derivative of the free energy with respect to the extent of the transition. Since the difference between the chemical potentials is a sign-variable function, in order for the product of the transition rate and the difference between the chemical potentials to be non-negative, it is sufficient that the transition rate be proportional to the difference indicated with a positive proportionality factor. The approach used directly relates the thermodynamic potential and the kinetics of the phase transition. The model developed shows that in the case of a linear dependence of the latent energy of the phase transition on the extent of the transition, the region of the transition is pressed against the melting curve in the pressure–temperature plane with the resultant appearance of classical phase equilibrium conditions. When the dependence is non-linear, the model of an incomplete, partial phase transition is realized. Phase transitions in a porous medium also lead to some special features of the seepage process. The conductivity of the medium in the initial state is low because the pores in the skeleton are partially or completely occupied by solid hydrates. When the hydrate dissociates, connected porous channels form, in which motion of the decomposition products occurs, i.e., multiphase seepage takes place under the conditions of variable permeability. When the temperature is reduced or the pressure is increased, the reverse process occurs, i.e., the mixture of the gas and water is converted into the hydrate, the porosity decreases, and seepage ceases. 2. Kinematics and conservation laws The model is constructed using an approach based on the hypothesis of interpenetrating and interacting continua.18,19 A continuum with the properties of the solid is defined to correspond to the rock and hydrates forming the skeleton, and continua with the properties of fluids are defined to correspond to the decomposition products. The quantities corresponding to the skeleton are marked with a subscript s, and the quantities referring to the water and the gas are denoted by the subscripts w and g, respectively. For brevity, all three continua are also marked with the subscript A, which runs through the values s, w and g, and the fluids (the water and gas) are marked with the Greek subscript ␣ = w, g. It is assumed that the skeleton is non-deformable and fixed in the reference system under consideration. This assumption, of course, is invalid if its cementing properties are lost as a result of the hydrate dissociation and the skeleton acquires viscous or viscoelastic properties. In this case, the strains of the skeleton should be taken into account by an element of the model. It is presumed that the components of the medium are in local thermal equilibrium, i.e., the skeleton, water and gas at the point x at the time t are characterized by one absolute temperature ␪(x, t) > 0. We assign to the material element of each of the continua comprising the porous medium the volume fractions ␾A , which are such that

where ␾ is the porosity of the medium. In the initial state the porosity is equal to ␾0 . The degree of saturation of the medium with fluid ␣ is defined as S␣ = ␾␣ /␾; therefore, the following relations hold for the volume fractions:

In the initial state the mass of the hydrate per unit volume equals ␳h ␾0 h (x), where ␳h is the true density of the undeformed hydrate and ␾0 h (x) is the volume fraction of the hydrate in the initial state, which is different in the general case for different points of the skeleton.

V.I. Kondaurov, A.V. Konyukhov / Journal of Applied Mathematics and Mechanics 75 (2011) 27–40

29

In the current state the mass of the hydrate per unit volume equals ␳h ␾h , where ␾h (x, t) is the current volume fraction of the hydrate. Suppose (2.1) is the extent of the phase transition, which is equal to the difference between the volume fractions of the hydrate in the initial and current states. The parameter m has an upper limit m0 ␾0 h . The value m = m0 corresponds to complete decomposition of the hydrate at a given point of the skeleton. Negative values of this parameter indicate formation of the hydrate. We will confine ourselves below to the case when mass transfer occurs between the skeleton and the fluids. There is no mass transfer between the water and the gas, i.e., the variation of the moisture content of the gas is neglected. Then the mass conservation laws of each continuum have the form

(2.2) Here and below rs is the averaged mass density of the skeleton, including both the rock and the hydrate, r␣ = ␾␣ ␳␣ is the averaged mass density of the fluid ␣ in the current state, ␳␣ is the true mass density, and ␰␣ = const denotes the stoichiometric coefficients (the mass fractions of water and gas in the hydrate), which are such that ␰1 + ␰2 = 1. Regarding the force fA acting on the continuum A enclosed in the region V, it is assumed that it is equal to the sum of three terms: the mass gravitational force fb A , the interaction force fint A of continuum A with other components of the medium, and the contact force fc A acting on continuum A on the boundary of the region under consideration, i.e.,

(2.3) where g is the gravitational force density, tA is the stress vector, and bA int is the interaction force density. In the case of slow, quasistatic flows, the interaction force is determined by the friction of the other components of the medium on the continuum under consideration. When dynamic phenomena are studied, the interaction force must also take into account the momentum exchange between the continua. It is assumed that the sum of the interaction forces is equal to zero: (2.4) As in the case in which there is no mass exchange between the continua,18,19 this hypothesis can be substantiated by averaging the momentum conservation laws for the skeleton and the fluids on the microscopic level. All the continua are assumed to be non-polar. It is also assumed that Cauchy’s postulate tA = tA (x, n) holds, i.e., the stress vector depends on the place x where the surface element is located and on the outward normal n, which specifies the orientation of the element. In this case the following representation holds: (2.5) where TA = TT A is the symmetric Cauchy partial stress tensor of continuum A. Henceforth the flow of the fluids is assumed to be slow, i.e., Euler’s number is

where ␳0 , 0 and ␴0 are the characteristic density, velocity and stress. In this case the forces of inertia and the forces appearing when mass is transferred can be neglected. The equations of motion of the water and the gas in this case become the equilibrium equations (2.6) Later an important role will be played by the decomposition of the thermodynamic parameters into equilibrium and dissipative components. We will call the state of a material element an equilibrium state if the temperature gradient and the mass velocities are equal to zero in a certain vicinity of the current time. The interaction forces b␣ int , the stresses T␣ and the volume fractions ␾␣ can be represented in the form of the sums (2.7) where b␣ eq , T␣ eq and ␾␣ eq are the equilibrium components, and b␣ dis , T␣ dis and ␾␣ dis are the dissipative components, which are equal to zero in the equilibrium state. The equilibrium interaction force of the fluid ␣ with the other components of the medium equals (2.8) Here and below ␴␣ is the true stress in the fluid ␣, which is such that T␣ = ␾␣ ␴␣ . To prove this, we write equilibrium equation (2.6) in the form

Taking into account that in a state of rest ·␴␣ eq + ␳␣ eq g = 0, we hence obtain formula (2.8).

30

V.I. Kondaurov, A.V. Konyukhov / Journal of Applied Mathematics and Mechanics 75 (2011) 27–40

The first principle of equilibrium thermodynamics in the quasistatic approximation under consideration is the equation of total internal energy balance for the porous medium enclosed in the region V

Here uA is the internal energy of unit mass of continuum A, q is the total heat flux and V is the fixed region in space occupied simultaneously by all the continua (the skeleton, the water and the gas). The three-dimensional heat sources are neglected. For smooth motions the energy balance equation is written in the form

Using continuity equation (2.2), equilibrium equation (2.6) and relation (2.4), we can reduce the local energy conservation law to the form

(2.9) Here and below d␣ f/dt = ∂f/∂t + v␣ · f is the derivative with respect to time along the trajectory of a particle of continuum ␣ and

˙ Qm = ␳h m

␰␣ (us − u␣ ) is the distributed energy source (sink) associated with the difference between the internal energies of the gas

␣

hydrate and its decomposition products. The second principle of non-equilibrium thermodynamics for a porous medium, saturated with water and gas and characterized by a common temperature for all the continua, is taken in the form of the integral inequality

where ␩A is the specific entropy of continuum A. Hence, in the case of smooth motions, the following local inequality follows for the rate of the entropy production

When continuity equation (2.2) is taken into account, this inequality can be written in the form

The entropy inequality plays an important role in the mechanics of continuous media, being a condition19 which must be satisfied by the constitutive relations for any smooth histories of the variation of state. To derive the constraints on the constitutive relations, it is next convenient to use the so-called reduced form of the entropy inequality, which does not contain the rate of heating

(2.10) Here

(2.11)

(2.12) where ␺A is the free-energy density of continuum A, ␺m is the distributed energy source (sink) associated with the difference between the internal energies and entropies of the gas hydrate and its decomposition products, ␦f is the seepage dissipation, and ␦T is the thermal dissipation. We stress that in the general case the seepage dissipation is equal to the power of the dissipative stresses and the forces of interaction of the water and gas with the skeleton. Only when there are no dissipative stresses is ␦f equal to the power of the interaction forces.

V.I. Kondaurov, A.V. Konyukhov / Journal of Applied Mathematics and Mechanics 75 (2011) 27–40

31

To derive of inequality (2.10) we use expansion (2.7) and expression (2.8) for the equilibrium interaction force. When they are taken into account, we arrive at the relation

Therefore,

(2.13) The seepage dissipation ␦f is specified by the first expression in (2.12). We also note that for a non-deformable skeleton saturated with water and gas, the divergence of the partial seepage rate ␾␣ v␣ can be expressed from continuity equation (2.2) in terms of the derivatives of the true density, volume fraction and extent of the transition with respect to time

(2.14) Substituting expression (2.13) into equation (2.9) of the local energy balance and taking into account equality (2.14), we arrive at the expression for the divergence of the heat flux

where s␣ eq = ␴␣ eq + p␣ eq I is the deviator of the true equilibrium stresses of fluid ␣. When formulas (2.11) and (2.12) are taken into account, substituting this expression into the entropy inequality gives relation (2.10). 3. The constitutive equations The following assumptions underlie the constitutive relations of a porous medium saturated with gas hydrates that undergo phase transition. The non-deformable skeleton is capable of storing (or releasing) energy when there is a change in temperature, a phase transition and a change in the degree of saturation. The change in the energy of the non-deformable skeleton when there is a change in temperature is related to the specific heat of the solid. Apart from the thermal form of energy, the skeleton is also characterized by surface energy, which is associated with the different wettabilities of the skeleton with water and the gas. On the macroscopic level the surface effects are described by the degree of saturation, i.e., the fraction of the pore space occupied by water, and by the capillary pressure, i.e., the difference between the pore pressures of the gas and the water. The energy change when there is a phase transition is specified by the latent energy for restructuring the hydrates. The dissociation products (water and gas) have a thermoelastic response, which means a dependence of the thermodynamic variables (the free energy, entropy, pressure, etc.) on the current values of the temperature and the mass densities. Unlike a thermoelastic porous medium, in which there are no structural changes, the state of the material under consideration is also characterized by the extent of the transition. The rate of transition is specified by the kinetic equation. It is assumed that the skeleton and the decomposition products are capable of conducting heat. The interaction of the skeleton and the fluids depends on their rates of relative motion. The response of the medium also includes the dependence of the degree of saturation with water and the gas on the state of the material. These qualitative statements are formalized as follows. The state of an elementary volume of the porous medium at the point x at the time t is given by the set of parameters (3.1) where ␳␣ is the true density, v␣ is the velocity of the particles of continuum ␣, ␪ is the temperature, m is the extent of the transition, and ␥ = ␪ is the temperature gradient. The response (reaction) of the material at the point (x, t) is characterized by the relations

(3.2) The rate of the phase transition is specified by the kinetic equation (3.3) In the material definition (3.1)–(3.3) we stress several special features. First, the system of constitutive relations (3.2) does not contain the stress tensor of the skeleton, since the skeleton is non-deformable. Second, the constitutive relations contain only one equation for the volume fractions. As such an equation, it is convenient to use the equation for the degree of saturation S = ␾w /␾, because the volume fractions are not independent quantities, since the sum of the volume fractions ␾␣ of the water and the gas is equal to the porosity ␾ = ␾0 + m, which depends linearly on the extent m of the phase transition. In addition, the response of the material under consideration includes only the dissipative part of the force of interaction of the fluids with the skeleton b␣ dis , rather than the total interaction force. This

32

V.I. Kondaurov, A.V. Konyukhov / Journal of Applied Mathematics and Mechanics 75 (2011) 27–40

is associated with the presence of equilibrium component (2.8), which has a universal character that does not depend on the properties of the saturating water and the gas. We also note that the fact that the constitutive relations contain kinetic equation (3.3) for the extent of the transition m means that the material under consideration has a memory of its previous states. This is associated with the fact that the solution of kinetic equation (3.3) is a functional defined on the prehistory of the parameters of state (␳␣ , ␪, ␪, v␣ ) at the point x. The assumption that the state of the medium is specified by the true densities of the water and the gas actually means that the constitutive equations of the water and the gas are invariant to any transformations of their initial state that do not alter the volume. The fact that the parameters of state do not include distortion rates (or gradients of the mass velocities that are related to them by compatibility equations17 ) means that the macroscopic viscous stresses of the water and gas are neglected. However, the viscous interaction of the fluids with the skeleton is believed to be a very important factor. In other words, each fluid is an “ideal gas”, whose particles do not interact with one another, but exchange momentum and energy with the skeleton. As evaluations19 and direct numerical modelling of flows of a viscous fluid in a pore space have shown, this assumption is valid for most models of a porous medium saturated with viscous fluids. The temperature gradient and the relative velocities, being “natural” for the heat flux and the interaction forces, are also included among the arguments of all the constitutive functions, including the thermodynamic potentials, by virtue of the “principle of equipresence”.16 The condition of thermodynamic consistency leads to partial decoupling, i.e., the potentials, the stress, the entropy, and the degree of saturation do not depend on the temperature gradient and the relative velocities. A dependence of the relative velocities on the temperature gradient remains only for the heat flux, the interaction forces, and the kinetic function. The principle of thermodynamic consistency16 means that inequality (2.10) for the rate of the entropy production should hold for any prehistory of states. For this requirement to be satisfied, it is necessary and sufficient that the constitutive equations of the porous medium under consideration should have the form

(3.4)

(3.5)

(3.6) (3.7)

(3.8) eq where the volume fractions ␾␣ , the degree of saturation S and the pore pressures p␣ are equal to their equilibrium values ␾␣ , S␣ eq

and p␣ eq . The quantity F(S, ␪, m) is the free energy of the skeleton expressed in terms of parameters of state other than those in the function ␺s (␳␤ , ␪, m), and ␦m is the phase transition dissipation, which is associated with the difference between the chemical potential of the skeleton ␹s and the chemical potential of the fluids ␹␣ . Formulae (3.4) mean that the free energies of each fluid, i.e., the functions of the true density of a given fluid and the temperature of the porous medium, do not depend on the density of the other fluid, the temperature gradient, the mass velocities and, most importantly, the phase transition parameter. The stress tensor of the fluid is spherical, i.e., the deviator of the stresses is identically equal to zero. The pore pressure and the entropy of the fluid are specified by partial derivatives of its free energy. Relations (3.5) show that the free energy of the non-deformable skeleton depends on the degree of saturation, the temperature and the extent of the phase transition. Like the free energies of the water and the gas, the function F does not depend on the temperature gradient and the velocities of the liquid and gaseous phases. The specific entropy of the skeleton and the capillary pressure, which is equal to the difference between the pore pressures of the gas and the water, are specified by the derivatives of the free energy F and do not depend on the temperature gradient and the velocities of the water and the gas. The equality (S, ␳␤ , ␪, m) = 0, which is formulated using the free energies of the components of the medium, specifies the degree of saturation S as a function of the true densities of the water and the gas, the temperature and the extent of the transition. Formulae (3.7) show that the heat flux, the dissipative interaction forces and the kinetic function depend, in general, on the temperature gradient and the velocities of the water and the gas, i.e., there is a relation between the heat transfer, the mass and the phase transitions. It follows from inequality (3.8) that the total dissipation is equal to the sum of the thermal dissipation associated with the propagation of heat, the seepage dissipation associated with the friction of the fluids on the skeleton and the phase transition dissipation, which depends on the kinetic function and the difference between the chemical potentials of the fluids and the skeleton. In the general case, the sum of all the forms of dissipation is non-negative. Proof of relations (3.4)–(3.8). We now turn to inequality (2.10) for the rate of entropy production. Suppose 0 = {␳␣ , ␪, m, ␥, v␣ } is an arbitrary state of the medium at the point (x0 , t0 ). We specify the state of the medium in a small vicinity of this point in the form of the linear function

where a and B are the rates of change and the gradients of the parameters of state of the porous medium at the point (x0 , t0 ). These ˙ which is specified by the kinetic function M(0 ), the temperature gradient ␥ = ␪, which quantities, except the phase transition rate m,

V.I. Kondaurov, A.V. Konyukhov / Journal of Applied Mathematics and Mechanics 75 (2011) 27–40

33

is one of the parameters of state, and the divergence of the partial seepage rate ·(␾␣ w␣ ), which can be expressed by virtue of (2.14), in terms of the parameters of state and their rates of variation, can take any values that do not depend on the state 0 . Since

the derivative of the equilibrium volume fraction has the form

(3.9) Substituting the expressions for the derivative of the equilibrium volume fraction (3.9) and the derivatives of the free energies of the water and the gas into inequality (2.10) and taking into account that d␣ m/dt = M + v␣ · m, we obtain

(3.10) eq

The coefficients in front of the gradients ⊗v␣ , m, the deviator of the tensor ⊗(␾␣ v␣ ) and the derivatives with respect to time at the point (x, t) → (x0 , t0 ) are functions of the state 0 , by virtue of definition (3.1)–(3.3) of the material under consideration. If the response of the material together with the derivatives of the functions ␺A and ␾␣ eq with respect to state is continuous at the point 0 , then, when (x, t) → (x0 , t0 ), from inequality (3.10) we obtain

(3.11)

(3.12)

(3.13)

(3.14) Relations (3.11) show that the deviators of the equilibrium stresses and the dissipative stresses of the water and the gas are equal to zero; therefore, the stresses T␣ = –␾␣ eq p␣ eq I. The free energies ␺␣ do not depend on the extent of the phase transition, the temperature gradient and the velocities of motion. The free energy of fluid ␣ does not depend on the density ␳␤ . Therefore, the thermodynamic potentials ␺␣ = ␺␣ (␳␣ , ␪). It follows from relations (3.12) that the entropy of fluid ␣ is specified by its free energy and that the volume fraction ␾␣ obeys the equality

Hence it follows that ␾␣ = ␾␣ (␳␣ , ␪, m). This means that the volume fractions, the degree of saturation and the pore pressures are identical to their equilibrium values and that the pore pressures p␣ = ␳2 ␣ ∂␺␣ /∂␳␣ . Thus, as in the classical case of a porous medium,17 we arrive at relations (3.4), by virtue of which the free energy of continuum ␣ is a function of the true density of this continuum and the temperature, and the pore pressure and specific entropy are specified by the corresponding partial derivatives. The first two relations in (3.13) show that the free energy of the skeleton does not depend on the temperature gradient and the velocities of motion of the water and the gas.

34

V.I. Kondaurov, A.V. Konyukhov / Journal of Applied Mathematics and Mechanics 75 (2011) 27–40

We introduce the two independent functions f␣ = f␣ (␳␤ , ␪, m), which are such that det ∂f␣ /∂␳␤ = / 0. We change from the arguments ␳␣ of the thermodynamic potential of the skeleton ␺s to the variables f␣ . Using the notation

and taking into account the rule for the differentiation of a composite function, we write the third relation in (3.13) in the form

We have taken into account here that the volume fractions and the pore pressures are equal to their equilibrium values. Suppose fw = S(␳␤ , ␪, m). For the volume fraction of the water ␾w = S␾(m) and the volume fraction of the gas ␾g = (1 – S)␾(m) we find

We reduce the third relation in (3.13) to the form

Since the determinant of the matrix ∂f␣ /∂␳␤ is non-zero, it follows that (3.15) When the results obtained are taken into account, the last equality in (3.13) is written in the form

Thus, the potential of the skeleton F does not depend on the temperature gradient and the mass velocities, and its derivatives with respect to the temperature and the degree of saturation are equal to the specific entropy of the skeleton and the capillary pressure. Using the definition of the capillary pressure pcap = pg – pw , the expression for the pore pressures p␣ = ␳2 ␣ ∂␺␣ /∂␳␣ and formula (3.15), we arrive at Eq. (3.6). The rule for the differentiation of a composite function together with the first formula in (3.15) gives

which leads to expression (3.8) for the phase transition dissipation. Inequality (3.8) follows from relation (3.14). Thus, the form (3.4)–(3.8) of the constitutive relations is necessary for the principle of thermodynamic consistency to hold. It is not difficult to verify that relations (3.4)–(3.8) are sufficient for inequality (2.12) to hold for any history of the variation of state (t). 4. The quasilinear approximation We will next examine a special case which is based on the following assumptions. The variation of the absolute temperature relative to its initial value and the extent of the phase transition are assumed to be small together with their derivatives with respect to space and time. The kinetic equation of the phase transition is taken in the form

(4.1) where K = const > 0 is the characteristic modulus, i.e., a quantity which has the dimensions of pressure, ␶ = const > 0 is the characteristic time of the transition, m0 is the volume fraction of the hydrate in the initial state, H(x) is the Heaviside function, which is equal to zero for x ≤ 0 and to unity for x > 0. The parameter ␬ vanishes when the hydrate decomposes completely (the extent of the transition reaches the value m0 ), as well as in the case of a shortage of water or gas for forming the hydrate when the degree of saturation S becomes equal to unity or zero. For a material with kinetics (4.1), the phase transition dissipation is non-negative for any states, since

The seepage dissipation and the thermal dissipation are assumed to be non-negative for any velocity of the water and the gas and any temperature field:

(4.2)

V.I. Kondaurov, A.V. Konyukhov / Journal of Applied Mathematics and Mechanics 75 (2011) 27–40

35

Conditions (4.1) and (4.2) are sufficient for the total dissipation to be non-negative, i.e., for condition (3.8) to hold, but, of course, they are not necessary conditions. The water is assumed to be a weakly compressible fluid. In this case the free energy in the vicinity of the initial state (␳0 w , ␪0 ) with the initial pressure p0 can be represented in the form of the quadratic expansion

where ␺0 w and ␩0 w are the energy and entropy in the initial state, ϑ = ␪ – ␪0 is the temperature change, ␳w = ␳w – ␳0 w is the mass density change, Kw is the bulk modulus, ␣w is the coefficient of thermal expansion, and cw is the specific heat of the water at constant volume. Hence we have the linear expressions for the pore pressure and the entropy

and the expression for the chemical potential of the water

(4.3) where ␹0 w = ␺0 w + p0 /␳0 w is the chemical potential in the initial state and cw (p) is the specific heat at constant pressure. The gas formed as a result of dissociation of the hydrate is assumed to be a perfect gas with a constant specific heat. The equations of state of such a gas have the form

where ␺0 g and ␩0 g are the free energy and entropy in the initial state, Mg is the molar mass, R is the universal gas constant, and cg = –␪∂2 ␺g /∂␪2 is the specific heat at constant volume. Hence we have the following expression for the chemical potential of the gas

(4.4) It is assumed that in the initial state the extent of the phase transition m = 0. Thermal and mechanical actions cause small changes in the temperature ␪ and the parameter m, but the final change in the degree of saturation satisfies the condition 0 < S < 1. It is assumed that the potential of the skeleton has the form

(4.5) Here Fc (S) is the surface energy associated with the difference between the surface tension of the water and the gas. If an element of the porous body is completely saturated with the former, more strongly wetting fluid in the initial state (S = 1), the value of the potential in this state is Fc = 0. Values in the range 0 < S < 1 correspond to Fc (S) > 0. The function Fm (S, m) is the latent energy of the phase transition per unit mass of the hydrate. This function is specified by the non-linear relation (4.6) It is assumed that

where K is the characteristic modulus. In the general case, ␥ and ␤ are functions of the degree of saturation S. The quantity F␪ (S, ϑ, m) is the thermal component of the potential of the skeleton, which depends on the degree of saturation, the temperature of the medium and the extent of the phase transition. The following approximation is taken for it (4.7) where ␩0 s is the entropy of the skeleton when there is no temperature change or phase transition, cs > 0 is the specific heat of the skeleton, and ␣ > 0 is the coefficient of mutual influence of the temperature and the phase transition, which are such that

These relations together with the assumptions

lead to the fact that all the terms in the expression for the potential of the skeleton F(S, ␪, m) are of the same order of smallness.

36

V.I. Kondaurov, A.V. Konyukhov / Journal of Applied Mathematics and Mechanics 75 (2011) 27–40

The potential (4.5)–(4.7) and formulae (3.5) generate linear expressions of the entropy and the capillary pressure

If the functions Fc (S), ␥(S) and ␩0 s (S) depend weakly on the degree of saturation, i.e., their derivatives with respect to the degree of saturation S are of the order of ␦, the following well-known relation19 holds in the linear approximation

By virtue of formulae (4.6) and (4.7), the chemical potential of the skeleton ␹s = F−(rs /␳h )∂F/∂m is written in the linear approximation in the form (4.8) (S), ␩0

1+␻ , ␣ϑm and c ϑ 2 /␪ , which have a higher order of smallness than the first, were discarded. Assuming s ϑ, ␥m, ␤m s 0 the chemical potentials ␹0 ␣ are of the order of ␦2 and taking into account formulae (4.3), (4.4) and (4.8), we arrive

where the terms Fc that the initial values of at the equation for the kinetics of the phase transition

(4.9) where

It follows from Eq. (4.9) that the phase transition begins (m = 0) under the condition that the state (S, ϑ, p␣ ) of an element of the medium reaches the surface

Within this surface (f0 (S, ϑ, p␣ ) ≤ 0) there are no phase transitions. It is seen that gas hydrate dissociation can be initiated by raising the temperature (ϑ > 0) or by lowering the pore pressure of the water and the gas (p␣ < p0 ). This intuitively obvious assertion holds under the condition ␥(S) > 0. The function ␥(S) > 0 plays the role of a threshold value, which must be exceeded as a result of an increase in temperature and/or a decrease in pressure in order for dissociation of the hydrate to begin. When the parameter m reaches the value m0 , dissociation is complete. This occurs when the state (S, ϑ, p␣ ) of an element of the medium reaches the surface

Outside this surface in the space of states, there are no phase transitions. It is not difficult to see that in the approximation under consideration the initial surface differs from the limiting shift by a quantity equal to ␤m0 ␻ . If the parameter ␶ in kinetic equation (4.9) is small compared with the characteristic time t0 of the change in state, it follows from Eq. (4.9) that the transition rate is constrained as ␶ → 0 if the quantity in the square brackets is equal to zero. Hence we obtain the following expression for the parameter m in the range 0 < m < m0

(4.10) This means that in such a process, which we will call an equilibrium phase transition, the extent of the phase transition is a function of the current thermodynamic state. The temperature change ϑ and the pore pressures p␣ can serve as independent variables, since the degree of saturation S can be expressed in terms of the capillary pressure pcap (S) = pg – pw . To formulate the constitutive equations for the heat flux and the dissipative components of the interaction forces, we will consider the equilibrium state (S, ϑ, m), in which the temperature gradient and the velocities of the fluids are equal to zero. The heat flux and the dissipative forces in such a state are also equal to zero: (4.11) Also, the symmetric second-rank equilibrium thermal conductivity and resistance tensors, which are defined by the expressions

have the properties for being “positive-definite”

(4.12)

V.I. Kondaurov, A.V. Konyukhov / Journal of Applied Mathematics and Mechanics 75 (2011) 27–40

37

but there is no mutual influence of the temperature gradient and the velocities, i.e., (4.13) In fact, it follows from the definitions of the seepage dissipation and the thermal dissipation that these quantities are equal to zero in the equilibrium state and are positive in all other states, i.e., the dissipation ␦f has a minimum when v␣ = 0, and the dissipation ␦T has a minimum when ␥ = 0. The necessary condition for an extremum

when v␣ = 0 and ␥ = 0 leads to equalities (4.11). Therefore, in the linear approximation

Consequently, the expressions for the dissipations have the form

Since these expressions are positive for ␥ = / 0 and v␣ = / 0, and since ␥ and v␣ are independent, properties (4.12) and (4.13) hence follow. Thus, in the linear approximation we have

(4.14) int

Since the interaction force b␣ = –p␣ ␾␣ + b␣ of the generalized multiphase seepage law

dis ,

using the second formula in (4.14), we can write equation of motion (2.6) in the form

(4.15) where W␥ = ␾␥ v␥ is the partial seepage vector. If instead of the resistance tensors Y␤␥ (S), we use the true permeability coefficients K␣␤ , which are such that

where ␮␣ denotes the dynamic viscosities of the liquid and gas phases, we arrive at the multiphase seepage law

(4.16) which is analogous to Darcy’s law, but is distinguished from it by the presence of cross terms. 5. Isothermal seepage of the dissociation products of a gas hydrate Consider the case in which the variation of the temperature can be neglected. It is assumed that the viscosities of the water and the gas are constant and that ␮g  ␮w . The skeleton is assumed to be isotropic. The gravitational force and the motion of the water in the pore space are neglected. The absolute permeability is assumed to be a function of the porosity and is specified by the formula (5.1) where k0 is the permeability in the initial state. The cross terms and the capillary pressure are neglected. Consequently, the seepage law in an isotropic sedimentary stratum is written in the form

(5.2) where fg (S) is the relative permeability of the gas and p pg is the pore pressure of the gas. The thermal component of the potential F␪ (S, ϑ, m) in expansion (4.5) is equal to zero. The quantities ␤, ␥ = const. The water is assumed to be incompressible. In the linear isothermal approximation the chemical potential of the water is (5.3)

38

V.I. Kondaurov, A.V. Konyukhov / Journal of Applied Mathematics and Mechanics 75 (2011) 27–40

The pressure of a perfect gas in an isothermal process is proportional to the density, that is, (5.4) and the chemical potential of the gas (4.4) equals

(5.5) An equilibrium phase transition is considered. By virtue of expression (4.10) for the extent of the transition m in the range 0 < m < m0 , the extent is given by the relation

(5.6) Hence it follows that the isothermal dissociation of the hydrate begins (m = 0) when the pressure of the gas is p* and ends when the pressure is p*, so that

(5.7) These two relations provide a way to determine the parameters ␥ and ␤ from data on the pressure at the beginning and end of dissociation. To formulate the differential equations that describe the variation of the pressure of the gas p and the degree of saturation S, we turn to continuity equations (2.2). Taking into account seepage law (5.2) and the incompressibility of water, in the approximation of stationary water we obtain

(5.8) where ␾(m) = m + ␾0 by virtue of formula (2.1) and m = meq (p) is the function of the pressure defined by expression (5.6). Next consider the one-dimensional problem of seepage of the gas formed as a result of dissociation of a gas hydrate in a porous medium occupying the half-space x ≥ 0. At the initial instant of time, the pore space is occupied by water with initial pressure p0 . On the boundary x = 0, a pressure pb < p0 is maintained. At infinity (x → ∞) the pore pressure is equal to the initial pressure p0 , and the degree of saturation is equal to S0 = 1. Therefore, the boundary and initial conditions for system (5.8) have the form

(5.9) We convert Eqs (5.8) and boundary conditions (5.9) to dimensionless form using the variables

As a result, we obtain the system of equations (here and below, the bars over the dimensionless variables are omitted)

(5.10) where

(5.11) and all the quantities on the left-hand sides of equalities (5.11) are of the order of unity. System (5.10) is closed by conditions (5.9). Taking the initial conditions into account, we have from the first equation in system (5.10)

Substituting this expression into the second equation in system (5.10), we obtain

(5.12)

V.I. Kondaurov, A.V. Konyukhov / Journal of Applied Mathematics and Mechanics 75 (2011) 27–40

39

Fig. 1.

Hence it is seen that in the region where the extent of transition is equal to zero (m = 0), the pressure distribution has a stationary character, since in this case H(p) = 0. In the region of complete transition (m = m0 ) we have the function H(p) = (1 – aw )m0 , i.e., the pressure varies with time. In the dissociation region the function takes the form

Suppose the pressure p tends to the pressure p*, at which dissociation begins. Then, when ␻ < 1, the functions (p) and H(p) tend to zero. Taking into account that H’(p) < 0 in the vicinity of p = p*, we find that the boundary-value problem for Eq. (5.12) is proper at least in the vicinity of the point p = p* if the parameter ␻, which characterizes the non-linearity of the energy of the phase transition, is less than unity. When ␻ = 1, we have

The necessary condition for the initial-value problem for Eq. (5.15) to be well-posed in this case has the form (5.13) In fact, let the pressure be equal to the pressure at the beginning of dissociation p = p*. Then (p*) = 0. Condition (5.13) follows from the inequality H(p*) > 0. To be specific, we will take a hydrate in which there are six molecules of water for every molecule of methane. Suppose the initial pressure is p0 = 100 atm ≈ 107 Pa, the initial temperature is ␪0 = 280 K, the molar mass of the gas is Mg = 16 × 10−3 kg/mole, the stoichiometric coefficients are ␰w = 0.87 and ␰g = 0.13, and the hydrate density is ␳h = 0.9 × 103 kg/m3 . Then the initial gas density is ␳0 g = 102 kg/m3 , and the coefficients are aw = 0.78 and ag = 1.16, i.e., condition (5.13) holds in this √ case. We will seek a self-similar solution that depends on the variable ␰ = x/ t. In this case Eq. (5.12) is written in the form of the ordinary second-order differential equation (5.14) with boundary conditions (5.15) Figure 1 shows self-similar solutions of boundary-value problem (5.14), (5.15). The upper curves represent the dependences p(␰), and the lower curves represent the dependences m(␰). The calculation was performed for the set of parameters

If the pressure pb = 0.2 on the boundary at ␰ = 0 is less than the pressure p* = 0.5 for the completion of dissociation, the solution represented by the solid curves holds. Complete dissociation of the gas hydrate occurs in a region of finite width. The pore pressure curve has two points of inflection. In the region adjacent to the boundary at ␰ = 0 the hydrate does not exist, and the extent of the transition m(␰) is equal to the limiting value m0 . The distributions of the pore pressure p(␰) and of the extent of the phase transition m(␰) for a small decrease in pressure (p* < pb < p*) are indicated by the dashed curves. In this case the hydrate decomposes incompletely, and the extent of the transition m(␰) in the region adjacent to the boundary at ␰ = 0 is less than the limiting value m0 . It is seen that unlike the model of complete decomposition of the hydrate on the dissociation front, the model √ under consideration leads to a broad dissociation region for any decrease in pressure. The width of this region increases with time as t. The limiting case in which the dissociation region contracts to a point is only obtained when p* → p* and corresponds to a solution with a kink on the pressure

40

V.I. Kondaurov, A.V. Konyukhov / Journal of Applied Mathematics and Mechanics 75 (2011) 27–40

Fig. 2.

and extent of the transition curves at the point of the phase transition. The derivatives of the pressure, the degree of saturation and the parameter m have a discontinuity. The slope of the curve p(␰) at ␰ = 0, which does not depend on the exponent n in formula (5.1), indicates that the yield of the gas decreases with time as t–1/2 . This phenomenon, which is associated with the increasing distance of the decomposition region from the extraction site, can become stronger if decomposition of the gas hydrate is accompanied by degradation of the strength properties of the skeleton, the appearance of viscous strains, filling of the pore channels and a reduction in the permeability. The study of these phenomena necessities taking into account the strains of the skeleton and the influence of dissociation of the hydrate on the rheological properties of the skeleton. The influence of the exponent n in the definition of the permeability and of the exponent l in the definition of the relative phase permeability of the gas on the solution was investigated. It was found that all other conditions being equal, an increase in n or l results in a decrease in the size of the decomposition region. The curves in Fig. 2 show the dependence of the pore pressure on the self-similar coordinate ␰ for various values of the initial porosity ␾0 . A significant influence of the value of ␾0 on the solution obtained can be seen. As ␾0 → 0, the coordinate of the dissociation front increases. The pore pressure distribution in the vicinity of the front becomes very abrupt, i.e., on the boundary of the dissociation region the pore pressure undergoes a jump with a height equal to p0 – p*. In the case in which the pores are filled with the gas in the initial state, an increase in the initial porosity causes the opposite effect. Acknowledgements This research was financed by the Programme of the Department of Power Engineering, Machine Building, Mechanics and Control Processes of the Russian Academy of Sciences for Controlling the Tribological and Strength Properties of Materials and Parts through Physicomechanical and Chemical Actions, the Russian Foundation for Basic Research (09-05-00542) and the Programme of the Russian Ministry of Education and Science for the Development of the Scientific Potential of Higher Schools (2009-2010). References 1. 2. 3. 4. 5. 6. 7. 8. 9. 10. 11. 12. 13. 14. 15. 16. 17. 18. 19.

Sloan ED. Clathrate Hydrates of Natural Gases. New York: Dekker; 1998. Verigin NN, Khabibullin IL, Khalikov GA. Linear problem of the decomposition of gas hydrates in a porous medium. Izv Akad Nauk SSSR MZhG 1980;(1):174–7. Cherskii NV, Bondarev EA. On the theory of exploiting gas hydrate deposits. Dokl Akad Nauk SSSR 1984;275(1):158–62. Bondarev EA, Maksimov AM, Tsypkin GG. On the mathematical modelling of gas hydrate dissociation. Dokl Akad Nauk SSSR 1989;308(3):575–8. Tsypkin GG. Analytical solution of the nonlinear problem of gas hydrate dissociation in a formation. Izv Ross Akad Nauk MZhG 2007;(5):133–42. Kamath VA, Godbole SP. Evaluation of hot-brine stimulation technique for gas production from natural gas hydrates. J Petroleum Technol 1987;39(11):1379–88. Selim MS, Sloan ED. Heat and mass transfer during the dissociation of hydrates in porous media. AIChE J 1989;35(6):1049–52. Nong H, Pooladi-Barvish M, Bishnoi PR. Analytical modeling of gas production from hydrates in porous media. J Can Petrol Technol 2003;42(11):45–51. Chuang J, Goodarz A, Duane HS. Constant rate natural gas production from a well in a hydrate reservoir. Energ Convers Manag 2003;44(15):2403–23. Vasilev VI, Popov VV, Tsypkin GG. Numerical investigation of the decomposition of gas hydrates co-existing with gas in natural sedimentary strata. Izv Ross Akad Nauk MZhG 2006;(4):127–34. Landau LD, Lifshitz EM. Statistical Physics. Oxford: Pergamon; 1969. Makogon YF, Holditch SA. Experiments illustrate hydrate morphology and kinetics. Gas Oil J 2001;99(7):45–50. Kim HC, Bishnoi PR, Heidemann RA, Rizvi SS. Kinetics of methane hydrate decomposition. Chem Eng Sci 1987;42(7):1645–53. Barenblatt GI, Yentov VM, Ryzhik VM. Motion of Liquids and Gases in Natural Sedimentary Strata. Moscow: Nedra; 1984. Kondaurov VI. A non-equilibrium model of a porous medium saturated with immiscible fluids. Prikl Mat Mekh 2009;73(1):121–42. Truesdell C. A First Course in Rational Continuum Mechanics. Baltimore: The Johns Hopkins University; 1972. Kondaurov VI. Mechanics and Thermodynamics of a Saturated Porous Medium. Moscow: MFTI; 2007. Nigmatulin RI. Principles of the Mechanics of Heterogeneous Media. Moscow: Nauka; 1978. Coussy O. Poromechanics. New York: Wiley; 2004.

Translated by P.S.