Thermo-elastic-plastic porous material undergoing thermal loading

International Journal of Engineering Science 37 (1999) 1985±2005

www.elsevier.com/locate/ijengsci

Thermo-elastic-plastic porous material undergoing thermal loading Andrzej Søu_zalec

*

Technical University of Czestochowa, 42-201 Czestochowa, Poland Received 8 August 1997; received in revised form 2 December 1998

Abstract A model for predicting elastic±plastic stresses within a surface-heated porous structure has been developed. The relevant phenomena for the moisture, pressure, temperature and displacement ®elds in thermoelastic-plastic porous material are analysed. Considering mass and energy transfer processes, a set of governing di€erential equations is presented. The solution of the problem has been obtained with a ®nite di€erence scheme. The results demonstrate the in¯uence of the evaporation mechanism on pressure and thermal stresses within the porous material. Ó 1999 Elsevier Science Ltd. All rights reserved.

1. Introduction Transport phenomena and elastic±plastic stresses in porous media have recently received growing attention in light of the common usage of such media in various applications in the ®elds of energy technology. The physical phenomena of moisture transfer in porous media are usually explained by the following theories: di€usion theory [1]; capillary ¯ow theory [2]; and evaporation condensation theory [3]. Di€usion theory, which has been used extensively in the past does not give results in agreement with experimental data. Capillary ¯ow theory, in which the di€usivity depends on the pore water content, has been employed [4,5] to predict the drying rate in concrete slabs. However, the coecient of di€usivity is a complex function of pore moisture, temperature and other variables of the porous system. It cannot be de®ned as a simple function of pore moisture as done in Refs. [4,5]. Later, drying experiments at elevated temperatures revealed high pore pressure as a consequence of intense water vaporization [6]. The measured pore pressure, determined under equilibrium heating conditions, correspond well to the sum of the calculated

*

Corresponding author. Tel.: +48-34-3-250-966; fax: +48-34-3-250-920. E-mail address: [email protected] (A. Søu_zalec)

0020-7225/99/$ ± see front matter Ó 1999 Elsevier Science Ltd. All rights reserved. PII: S 0 0 2 0 - 7 2 2 5 ( 9 9 ) 0 0 0 1 1 - 7

1986

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

saturated water vapor pressure and ambient air pressure. As a consequence of the above models, in a previous work the evaporation mechanism was assumed with concentration and pressure gradient terms [7]. An analysis of poroelastic stresses was developed and improved by Biot [8±10]. In 1976 Rice and Cleary [11] gave a rational reformulation of the Biot theory, and their rationalised version has proved more convenient for solving relevent problems and in interpreting the solutions obtained (e.g. [12±15]). Modern theories of mixtures, taking into account ®nite deformations and thermal e€ects, for most cases have been developed by Crochet and Naghdi [16] and Bowen [17]. In 1971 Schi€man [18] developed an extended Biot theory including thermal e€ects. The present paper describes equations of heat and mass transfer in a thermo-elasto-plastic material. This theory shows that the displacement ®eld is, in general, completely coupled with the pore pressure and temperature ®elds.

2. Assumptions and equation 2.1. Theoretical assumptions and simpli®cations 1. 2. 3. 4. 5.

The following assumptions are made for the theoretical model: Heat transfer between the ¯uid and the solid is neglected. The response time for local heat transfer between the ¯uid and the solid is several orders of magnitude smaller than the times of interest. Liquid vapor equilibrium exists in the presence of free water, which makes the partial pressure of the vapor equal to the saturation pressure. Movement of the liquid is neglected. Darcy's law with a variable coecient holds for the gases. Air and water vapor are treated as ideal gases.

2.2. Conservation equations for mass The conservation equations for mass can be written as q_ i ˆ ÿr…qi wi † ‡ Wi ;

…1†

where qi is the density of species i, wi is the velocity of species i and Wi is the production rate of species i. Since no movement of the liquid (subscript c) is assumed wc ˆ 0. Also, W ˆ ÿWV ˆ ÿWm (v ˆ vapor, m ˆ air±vapor mixture), since the rate of liquid evaporation is the same as the rate of vapor production, and since the air does not change phase, the rate of mixture production equals the rate of vapor production. The above equation thus simpli®es to q_ c ˆ ÿWm :

…2†

The continuity equation for the air (subscript a) takes the form q_ a ˆ ÿr…qa wa †

…3†

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

1987

and for the vapor, q_ v ˆ ÿr…qv wv † ‡ Wv :

…4†

The conservation of gas phase mass gives q_ m ˆ ÿr…qm wm † ‡ Wm :

…5†

Fick's law allows the ¯uxes to be presented in the forms of Eqs. (6) and (7): ja ˆ qa …wa ÿ wm † ˆ ÿqm Drqba ;

…6†

where qbi ˆ qi =qm is the mass fraction of species i with respect to the density of the air±vapor mixture, and D is the di€usion coecient for Fick's law for the air±vapor mixture; and jv ˆ qv …wv ÿ wm † ˆ ÿqm Drqbv :

…7†

Finally, we get the following species equations: qmba ‡ qm wm rqba ˆ r…qm Drqba † ÿ qba Wm

…8†

qmbv ‡ qm wm rqbv ˆ r…qm Drqbv † ‡ …1 ÿ qbv †Wm :

…9†

and

2.3. Thermal equations The ¯uxes of heat q and ¯owing gases r can be expressed as q ˆ ÿkrh

…10†

r ˆ qa wa ha ‡ qv wv hv ;

…11†

and

where k is the thermal conductivity, h is the temperature and hi is the enthalpy of component i per unit mass of component i. Eq. (11) can be transformed [19] to the form r ˆ qm wm hm ÿ qm Dha rqba ÿ qm Dhv rqbv

…12†

r ˆ qm wm hm ‡ qm D…hv ÿ ha †rqba :

…13†

or

Assuming that qe ˆ qs hs ‡ qc hc ‡ qm hm ÿ qm Rm h;

…14†

where e is the thermal energy and R is the gas constant, the thermal equations can be expressed as

1988

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

qcp h_ ˆ rkrh ÿ qm ‰wm cpm ‡ D…cpv ÿ cpa †rqba Šrh ÿ ‰…hv ÿ ha †Wm ÿ qm rm_ hŠ;

…15†

where cp is the speci®c heat at constant pressure. 2.4. Darcy's law The velocity of the air±vapor mixture is given by wm ˆ ÿkD rp ;

…16†

where kD is Darcy's coecient and p is the pressure. 2.5. Thermodynamic relations Assuming that the vapor and air are ideal gases we have the following relations. 2.5.1. Ideal gas equation for the vapor pv V v ˆ qv Rv h;

…17†

where V i ˆ qi =qai represents the volume occupied by component i per unit total volume. 2.5.2. Ideal gas equation for the air pa V a ˆ qa Ra h:

…18†

2.6. Clausius±Clapeyron equation Since the liquid and vapor are assumed to be in equilibrium, pv ˆ psat …h† in the presence of liquid water. An analytic expression for psat is   A ÿ…B=Rv † exp ÿ : psat …h† ˆ Ch Rv T

…19†

…20†

2.7. State equation Using the notations qbi and V i , the state equation can be presented as pv …/ ÿ V c † ˆ qm qbv Rv h ˆ …1 ÿ qba †qm Rv h and

…21†

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

pa …/ ÿ V c † ˆ qm qba Ra h;

1989

…22†

where V v ˆ V a ˆ …/ ÿ V c † and / is the porosity. Combining the above equations we get p…/ ÿ V c † ˆ qm Rm h;

…23†

where Rm ˆ qbv Rv ‡ qba Ra :

…24†

2.8. Porosity In a linearized approach, the current porosity / depends on the current values of ¯uid mass density, strain and increase in ¯uid mass content according to relation   qc ÿ qc;o m ‡ tr e …25† ÿ /o / ˆ /o ‡ qc;o qc;o where m is the ¯uid mass content, qc;o is the ¯uid mass density in the reference state, /o is the porosity in the reference state. In in®nitesimal transformation, the trace tr e ˆ eii of the linearized stress tensor e represents the volume change per unit volume in the deformation. It is called volume dilatation. 2.9. Thermo-elastic-plastic constitutive equations 2.9.1. Isotropic hardening The basic hypothesis of isotropic hardening is the assumption that the shape of the yield surface is unchanging, and its growth can be described by one scalar parameter which is a function of plastic deformation. We can describe the plastic potential function F of the porous material by the porous material by the equation [24] F ˆ F …rij ; K; h; p†

…26†

where rij is the stress tensor, and K is the work ± hardening parameter. According to the Quinney±Taylor (QT) hypothesis [23] K_ ˆ rij e_ pij :

…27†

If we di€erentiate F by using the chain rule for partial di€erentiation we obtain oF oF _ oF _ oF _ F_ ˆ r_ ij ‡ p: K‡ h‡ orij oK oh op

…28†

Using the QT hypothesis, it is seen that the second term in the right-hand side of Eq. (28) can be expressed in terms of plastic strains rates e_ pij as

1990

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

oF _ oF oK p e_ : Kˆ oK oK oepij ij

…29†

Combining Eqs. (28) and (29) gives oF oF oK p oF _ oF _ F_ ˆ r_ ij ‡ e_ ‡ h‡ p: orij oK oepij ij oh op

…30†

Equilibrium conditions require that variation of the plastic energy be stationary oF oF oK p oF _ oF F_ ˆ r_ ij ‡ e_ ‡ h‡ p_ ˆ 0: orij oK oepij ij oh op

…31†

In the approach presented we assume the small strain thermo-elastic-plasticity theory in which the total small strain rate is a sum of elastic, thermal, plastic and strain rate due to pressure changes. By this assumption strain rate tensor in a thermo-elasto-plastic process in a solid body is of the form e_ ij ˆ e_ eij ‡ e_ Tij ‡ e_ bij ‡ e_ pij ;

…32†

where e_ ij ; e_ eij ; e_ Tij ; e_ bij ; e_ pij are the rates of the total, elastic, thermal, pressure and plastic strain tensor, respectively. After rearranging Eq. (32) we can obtain the components of the elastic strain rate tensor e_ eij ˆ e_ ij ÿ e_ Tij ÿ e_ bij ÿ e_ pij :

…33†

Making use of Hooke's law the rates of change of the total stress are given as e …_eekl ÿ e_ Tkl ÿ e_ bkl ÿ e_ pkl †; r_ ij ˆ Cijkl e where Cijkl are components of the elasticity tensor   1 2 e Cijkl ˆ dik djl ‡ dil djk ÿ dij dkl ; 4G 1‡m

…34†

…35†

and G is the shear modulus, m is Poisson's ratio, dij is the Kronecker delta. Combining Eqs. (31), (33) and (34) and considering the ¯ow rule e_ pij ˆ k

oF ; orij

…36†

we get oF e oF oK oF oF _ oF _ Cijkl …_ekl ÿ e_ Tkl ÿ e_ bkl ÿ e_ pkl † ‡ k‡ h‡ p: p orij oK oeij orij oh op

…37†

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

After some calculations we obtain the relation for the proportionality factor k   T b oF oF _ oF e _ _ _ h ‡ p C …_ e ÿ e ÿ e † ‡ kl kl orij ijkl kl oh op   kˆ : oF oF oF oK oF e C ÿ p orpq pqrs orrs oK oe opq

1991

…38†

pq

Introduce the notations Sˆ

oF e oF oF oK oF C ÿ ; orpq pqrs orrs oK oeppq orpq

…39†

_ e_ Tij ˆ aij h;

…40†

aij ˆ ah dij

…41†

_ e_ bij ˆ Bij p;

…42†

Bij ˆ

3…mu ÿ m† dij ; 2GBu …1 ‡ m†…1 ‡ mu †

…43†

where mu is the undrained Poisson ratio, ah is the thermal expansion coecient and Bu is the induced pore pressure parameter. In general, 0 < Bu 6 1;

1 0 < m < mu 6 : 2

From Eqs. (38)±(42) we get   1 oF e oF _ oF _ _ ‡ C …_ekl ÿ akl h ÿ Bkl p† kˆ h‡ p_ : S orij ijkl oh op

…44†

…45†

Next combining Eqs. (34) and (45) we have e Cijkl oF e e _ _ekl ÿ Cijkl _rij ˆ Cijkl _ …akl h ‡ Bkl p† ‡ S orkl   oF e oF _ oF _ _ ‡  C …_ekl ÿ akl h ÿ Bkl p† h‡ p_ : orij ijkl oh op

…46†

After rearranging the terms

1 e oF oF e e e _ ÿ Cijkl …akl h_ ‡ Bkl p† C e_ kl e_ kl ÿ Cijkl r_ ij ˆ Cijkl S orkl orij ijkl   1 e oF oF e 1 e oF oF _ oF _ _ ÿ Cijkl ‡ Cijkl C …akl h ‡ Bkl p† h‡ p_ : S orkl orij ijkl S orkl oh op

…47†

1992

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

Assuming the Huber±Mises (HM) condition and making use of the fact oF ˆ sij orij

…48†

we can de®ne the so-called plasticity tensor 1 e oF oF e 1 p e ˆ Cijkl C ˆ C e skl sij Cijkl Cijkl S orkl orij ijkl S ijkl

…49†

where sij are components of the deviatoric stress tensor. 1 sij ˆ rij ÿ rkk dij : 3

…50†

By substituting Eq. (49) into Eq. (47) we get p e e e _ ÿ Cijkl _ r_ ij ˆ Cijkl …akl h_ ‡ Bkl p† …akl h_ ‡ Bkl p† e_ kl ÿ Cijkl e_ kl ‡ Cijkl   1 e oF oF _ oF ÿ Cijkl h‡ p_ : S orkl oh op

…51†

De®ning the so-called elasto-plasticity tensor ep p e ˆ Cijkl ÿ Cijkl Cijkl

…52†

the thermo-elasto-plastic constitutive equation can be expressed as r_ ij ˆ

ep Cijkl e_ kl

ÿ

ep Cijkl …akl h_

e Cijkl skl _ ÿ ‡ Bkl p† S



 oF _ oF h‡ p_ : oh op

…53†

One may ®nd that during plastic deformation of a solid e skl ˆ 2Gsij Cijkl

…54†

due to the fact that the ®rst deviatoric strain invariant vanishes. It leads to the following relation p on Cijkl : e e Cijkl skl sij Cijkl ˆ 2G sij skl 2G ˆ 4G2 sij skl :

…55†

e srs ˆ 2Gspq the ®rst term in Eq. (39) can be Now we evaluate the S given by Eq. (39). Since Cpqrs expressed as

oF e oF 2 r2 4 2 e Cpqrs ˆ spq Cpqrs srs ˆ spq 2Gspq ˆ 2G ˆ G r; 3 orpq orrs 3

…56†

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

1993

where 1 2 1  ˆ spq spq : r 3 2

…57†

By Eq. (27) and Eq. (36) we have oF K_ ˆ krij : orij

…58†

From Eq. (27) we get rij ˆ

K_ : e_ pij

…59†

During plastic deformation of work-hardening materials, the yield strength increases with load. The de®nitions of the HM plastic potential function can be used to derive the following useful relation:   oF o 1 2 2 o r  ˆÿ r  ˆ J2s ÿ r ; oK oK 3 3 oK

…60†

where J2s is the second invariant of the devatoric stress. If we refer to a typical ¯ow curve, we have dep ; dK ˆ r

…61†

dep 1 ˆ :  dK r

…62†

or

Since o r o r dep 1 H0 ˆ H0 ˆ ; ˆ p   r oK oe dK r

…63†

where H 0 is the plastic modulus of the material in a multiaxial stress state (see Fig. 1) H0 ˆ

o r : oep

Eq. (60) will have the form  0 oF 2 H 2H 0  ; ˆ ˆÿ r  3 r oK 3 which leads to the following relation:

…64†

…65†

1994

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

Fig. 1. Slope H 0 as a plastic modulus of the material during plastic deformation.

oF oK oF 2H 0 2 2 2 4 2 0  ˆÿ r  H: ˆÿ rpq ssq ˆ ÿ H 0 r p 3 oK oepq orpq 3 3 9

…66†

Then S can be expressed in terms of the material properties and the state of stress as follows:   4 4 2 0 4 2 2  H ˆÿ r  …3G ‡ H 0 †: …67† r ÿ ÿ r S ˆ G 3 9 9 2.9.2. Kinematic hardening In processes of kinematic hardening, in order to describe the motion of the initial yield surface, one introduces the translation tensor aij , whose components determine a new position of the yield surface centre. Ziegler [20] assumes a motion of the surface in the direction of the di€erence of r and a: _ ij ÿ aij †; a_ ij ˆ l…r

…68†

where l_ is the multiplier. Melan [21] has proposed de®nition of the tensor aij as a_ ij ˆ C e_ pij :

…69†

This de®nition was also used by Prager [22]. Considering the modi®cation of the yield surface with respect to temperature and pressure we assume C ˆ C…h; p†:

…70†

The plastic potential function F can be expressed as F ˆ F …rij ; aij ; h; p†:

…71†

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

1995

The total di€erential of F is oF oF oF _ oF F_ ˆ r_ ij ‡ a_ ij ‡ h‡ p_ ˆ 0: orij oaij oh op

…72†

The condition, Eq. (72), is sometimes called the compatibility equation for the plastic yielding condition. On the basic of Eq. (68) we can obtain the following relations: oF oF ˆÿ ; oaij orij orij ˆ ÿ1: oaij

…73† …74†

If one substitutes the above expressions into Eq. (72) the following relations is obtained: oF oF _ oF …r_ ij ÿ a_ ij † ‡ h‡ p_ ˆ 0 orij oh op

…75†

oF oF oF _ oF _ a_ ij ˆ r_ ij ‡ h‡ p: orij orij oh op

…76†

or

Further substitution of Eq. (68) into the above expression yields   oF oF _ oF _ _ h ‡ p ‡ r ij orij oh op  : l_ ˆ  oF …rkl ÿ akl † orkl

…77†

In kinematic hardening theory we de®ne the so-called translated stress tensor rij ˆ rij ÿ aij ;

…78†

and translated stress deviators 1 sij ˆ rij ÿ rkk dij : 3

…79†

The yield function becomes 2 F ˆ sij sij ; 3

…80†

which leads to the expression osij oF ˆ 3sij ˆ 6sij : orij orij

…81†

1996

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

By substituting Eq. (81) into Eqs. (77) and (78) we get u_ ˆ

kˆ

h_ ‡ oF p_ 6sij r_ ij ‡ oF oh op 6rpq spq

;

…82†

h_ ‡ oF p_ 6sij r_ ij ‡ oF 1 oh op ;   36spq spq C…h; p†

…83†

e_ pij ˆ 6sij k:

…84†

Assuming the decomposition (32) e_ ij ˆ e_ eij ‡ e_ Tij ‡ e_ bij ‡ e_ pij ;

…85†

and using Hooke's law one obtains e r_ ij ˆ Cijkl …_ekl ÿ e_ Tkl ÿ e_ bkl ÿ e_ pkl †

…86†

e e r_ ij ˆ Cijkl …akl h_ ‡ Bkl p_ ‡ 6skl k†: e_ kl ÿ Cijkl

…87†

or

Substituting Eq. (69) and Eq. (81) into Eq. (36) gives a_ ij ˆ 6C…h; p†sij k:

…88†

The above expression, after substituting into Eq. (72), gives oF oF _ oF F_ ˆ …r_ ij ÿ 6C…h; p†sij k† ‡ h‡ p_ ˆ 0: orij oh op

…89†

By Eqs. (89) and (83) we get oF oF _ dF e e _ ‡ …Cijkl …akl h_ ‡ Bkl p†† e_ kl ÿ Cijkl h‡ p_ orij oh op : kˆ oF e s † 6…C…h; p†spq ‡ Cpqrs rs orpq

…90†

Substituting Eq. (90) into Eq. (87) we get e e _ e_ kl ÿ Cijkl …akl h_ ‡ Bkl p† r_ ij ˆ Cijkl

ÿ

oF e  orij Cijkl skl

e e _ ‡ oF h_ ‡ oF p_ …Cijkl …akl h_ ‡ Bkl p†† e_ kl ÿ Cijkl oh op  s † oF 6…C…h; p†spq ‡ Cpqrs rs orpq

:

…91†

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

1997

Let e srs † S ˆ …C…h; p†spq ‡ Cpqrs

oF e ˆ 6spq spq C…h; p† ‡ 6spq Cpqrs srs : orpq

…92†

With the aid of Eq. (92) the relation Eq. (91) becomes oF 1 e e e e _ ÿ Cijkl _ …akl h_ ‡ Bkl p† skl e_ kl ÿ Cijkl …C e_ kl …akl h_ ‡ Bkl p†† r_ ij ˆ Cijkl orij S ijkl      1 e  oF _ oF 1 e e  oF e ÿ Cijkl skl e_ kl C h‡ p_ ˆ Cijkl ÿ Cijkl skl S oh op S orij ijkl      1 1 e  oF _ oF e e  oF e _ _ ÿ Cijkl skl …akl h ‡ Bkl p† ÿ Cijkl ÿ C Cijkl skl h‡ p_ : S orij ijkl S oh op

…93†

Assuming that   1 1 p e  oF e e e Cijkl ˆ …6Cijkl skl sij Cijkl † ˆ Cijkl Cijkl skl s orij S

…94†

and ep p e Cijkl ˆ Cijkl ÿ Cijkl ;

…95†

the constitutive equation for kinematic hardening material subject to thermo-elasto-plastic deformation takes the form r_ ij ˆ

ep Cijkl e_ kl

ÿ

ep Cijkl …akl h_

  1 e  oF _ oF _ ÿ Cijkl skl ‡ Bkl p† h‡ p_ : S oh op

…96†

2.10. Momentum balance If no body force exists, the momentum balance in the quasi-static nonisothermal context is rij;j ˆ 0:

…97†

2.11. Strain±displacement relation The deformation of the material is described by the strain tensor eij , which is de®ned in terms of the displacement ui of the solid constituent as follows: 1 eij ˆ …ui;j ‡ uj;i †: 2

…98†

1998

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

3. Example 3.1. Simpli®ed equations for a 1D axisymmetrical problem One of the simple examples of the thermo-poro-elasticity theory presented in this paper is 1D axisymmetrical problem. In such a case the set of governing equations and suitable boundary conditions in a cylindrical coordinate system take the form as follows. 3.1.1. Thermal equations

qcp

    oh o2 h 1 o D…cpr ÿ cpa † oqba oh ˆk 2‡ r ÿ qm cpm wm ‡ or or ot or r or cpm   o ÿ …hv ÿ ha †Wm ÿ …qm Rm h† ; ot

…99†

h…r; 0† ˆ h0 …r†;

…100†

oh …0; t† ˆ 0; or

…101†

and ÿk

oh …R; t† ˆ h‰h…R; t† ÿ hf …t†Š ‡ f r‰h4 …R; t† ÿ h4f …t†Š; or

…102†

where R is the radius of the element, hf is the environmental temperature, h is the convective heat transfer coecient, f is the e€ective shape factor for radiation and r is the Stefan±Boltzman constant. 3.1.2. Species equations   oqba o2 qba oqba qba Wm 1 o ˆD ÿ ‡ ; …rqm D† ÿ wm 2 ot or or qm qm r or

…103†

qba …r; 0† ˆ qba;0 …r†;

…104†

oqba …0; t† ˆ 0: or

…105†

and

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

1999

3.1.3. Continuity equations oqm 1 o ‡ …rqm wm † ˆ Wm ; ot r or

…106†

wm …0; t† ˆ 0;

…107†

oqc ˆ ÿWm ; ot

…108†

qc …r; 0† ˆ qc;0 …t†:

…109†

and

3.1.4. Field equations (

drrr

)

" ˆ 2G

drhh

1ÿm 1ÿ2m

m 1ÿ2m

m 1ÿ2m

1ÿm 1ÿ2m

err ehh

derr

) ÿ

3…ml ÿ m† dp Bl …1 ‡ ml †…1 ÿ 2m† #( ) d 2rr srr shh

dehh " 2 2G…1 ‡ m† 2G srr ÿ adh ÿ 1 ÿ 2m So srr shh s2hh d 2hh " 2 #  srr shh 2G srr 3…ml ÿ m† adh ‡ dp ‡ So srr s2hh s2hh 2GBl …1 ‡ m†…1 ‡ ml † " 1ÿm #( )  m srr 2G 1ÿ2m 1ÿ2m oF oF dh ‡ dp ÿ m 1ÿm So oh op shh 1ÿ2m



#(

 ˆ

…110†

1ÿ2m

 ou  r

or ur r

…111†

where   2 2 H0  1‡ So ˆ r : 3G 3

…112†

The plastic potential function is assumed in the form r 1 2 F ˆ …rrr ÿ rhh † ‡ …Bp ‡ ch ÿ k† ˆ 0 2 where b; c and k are material parameters.

…113†

2000

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

The thermodynamics relations are given by Eq. (17) and Eq. (18), the Clausius±Clapeyron equation by Eq. (20) and Darcy's law by Eq. (16). 3.2. The solution The solution is obtained by an implicit ®nite di€erence technique with a constant grid and variable time step sizes. Since the equations describing heat and mass transfer in porous materials are partial di€erential equations of parabolic type, the solutions of such equations are well known and can be found in many textbooks on numerical analysis. These data do not provide any special information regarding this paper and have been omitted. For most heating situations with rapid heating a time step size of 1 s is found to be satisfactory for the calculations. A mesh size of slightly less than 2  10ÿ3 m was found to be adequate for the example. 3.3. Temperature, pressure and stresses in a 1D axisymmetrical element The speci®c case of a 1D axisymmetrical structural element with a uniform initial temperature is considered to illustrate the results of the analysis (Fig. 2). The radius of the considered cylindrical element is 0.12 m. The initial moisture content is from 0 to 0.108 m and the remaining 0.108 to 0.12 m is supposed to be dry. The assumption of the existence of a dry region close to the surface is often validated. The thermal and mechanical parameters used are presented in Table 1. Fig. 3 presents the temperature in the porous element as a function of time for the assumed heating curve for a ˆ 10 J=s m2 K, kD ˆ 10ÿ11 m3 s=kg and qc;o ˆ 70 kg=m3 . In Fig. 4 one can observe maximum pressures for various values of kD and in Fig. 5 changes of pressure with time.

Fig. 2. Cylinder under consideration.

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

2001

Table 1 Thermal and mechanical parameters used in this study cp ˆ 1040 J=kg K

kD ˆ 5  10ÿ12 to 1 m3 s=kg

D ˆ 2:142  10ÿ5 m2 =s f ˆ 0:9 a ˆ 0 to A aD ˆ 1 A ˆ 3:18  106 J=kg B ˆ 2470J=kg K C ˆ 6:05  1026 N=m2 ah ˆ 10  10ÿ6 1=K b ˆ 0:1 c ˆ 0:1 K ˆ 5  106 N=m2

qba1 ˆ 1:0 a ˆ 5:22  10ÿ7 m2 =s / ˆ 0:2 q ˆ 2400 kg=m3 qc;o ˆ 0:200 kg=m3 Bu ˆ 0:7 m ˆ 0:16 mu ˆ 0:17 E ˆ 2  1010 N=m2

Fig. 6 shows the distribution of stresses for maximum pressures (kD ˆ 10ÿ12 m3 s=kg and t ˆ 3  103 s). Figs. 4 and 5 show the pressure pro®les and pressure histories respectively for the several values of Darcy's coecient. Fig. 4 indicates that for lower values of Darcy's coecients pressures in elements increases. So for very porous material such as sand the pressures are very low and we do not observe signi®cant deformations due to pressure e€ect during heating. Increasing of Darcy's coecient shifts the point of maximum pressure toward the inside. This is consistent with the preceding argument since to have reasonable ¯ow resistance the point should

Fig. 3. Temperatures in the structural element for a ˆ 10 =sm2 K, kD ˆ 10ÿ11 m3 s=kg and qc;0 ˆ 70 kg=m3 : (a) r ˆ 120 mm; (b) r ˆ 100 mm; (c) r ˆ 80 mm; (d) r ˆ 60 mm.

2002

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

Fig. 4. Maximum pressures as a function of the distance from the heating surface for various kD values: (a) kD ˆ 10ÿ12 m3 s=kg; (b) kD ˆ 10ÿ10 m3 s=kg.

Fig. 5. Changes of pressure as a function of time for various kD values: (a) kD ˆ 10ÿ12 m3 s=kg; (b) kD ˆ 10ÿ10 m3 s=kg.

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

Fig. 6. Radial (a) and circumferential (b) stresses in an element for kD ˆ 10ÿ12 m3 s=kg and t ˆ 3  103 s.

2003

2004

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

be farther inside the medium for a higher value of Darcy's coecient. The pressure term in¯uences on the stresses pro®les in the element.

4. Concluding remarks An analysis is developed for heat and mass transfer in a wet porous medium subject to unsteady, nonlinear boundary conditions. The simpli®ed equations have been solved simultaneously by an implicit ®nite di€erence technique. The assumption of no liquid movement was made in the development of the present theory. This assumption is probably valid for low and medium pressures. However, for high pressures it is questionable.

References [1] V.K. Sherwood, Application of the theoretical di€usion equations to the drying of solids, Trans Am. Inst. Mech. Eng. 27 (1931) 190±202. [2] R. Buckingham, Studies in the movement of soil moisture, US Department Agr. Bur. Soils Bull. 38 (1907) 29±61. [3] C.G. Gurr, T.J. Marshall, J.T. Hutton, Movement of water in soil due to a temperature gradient, Soil Science 74 (1952) 335±345. [4] S.E. Pihlajavaara, Introductory bibliography for research on drying of concrete. The state Institute for Technical Research, Helsinki, 1964. [5] Z.P. Bazant, L.J. Najjar, Nonlinear water di€usion in nonsaturated concrete, Materials Constructions 5 (25) (1972) 3±20. [6] S. Whitaker, A theory of drying in porous media, Adr. Heat Transfer 12 (1977) 34. [7] Z.P. Bazant, Constitutive equation for concrete and shrinkage based on thermodynamics of multiphase systems, Materials Constructions 3 (1970) 13. [8] M.A. Biot, General theory of three-dimensional consolidation, J. Appl. Phys. 12 (1941) 155. [9] M.A. Biot, Theory of elasticity and consolidation for a porous anisotropic solid, J. Appl. Phys. 26 (1955) 182. [10] M.A. Biot, General solutions of the equations of elasticity and consolidation for a porous material, J. Appl. Phys. 78 (1956) 91. [11] J.R. Rice, M.P. Cleary, Some basic stress di€usion solutions for ¯uid-saturated elastic porous media with compressible constituents, Rev. Geophys. Space Phys. 14 (1976) 227. [12] M.P. Cleary, Fundamental solutions for ¯uid-saturated porous media and application to localized rupture phenomena, Ph.D. thesis, Univ. Micro®lms Int., Ann Arbor, MI, 1976. [13] M.P. Cleary, Fundamental solutions for a ¯uid-saturated porous solid, Int. J. Sol. Struc. 13 (1977) 785. [14] M.P. Cleary, Moving singularities in elasto-di€usive solids with applications to fracture propagation, Int. J. Sol. Struc. 14 (1978) 81. [15] J.W. Rudnicki, Plane strain dislocations in linear elastic di€usive solids, J. Appl. Mech. 54 (1987) 545. [16] M.J. Crochet, P.M. Naghdi, On constitutive equations for ¯ow of ¯uid through an elastic solid, Int. J. Eng. Sci. 4 (1966) 383. [17] R.M. Bowen, Incompressible porous media models by use of the theory of mixtures, Int. J. Eng. Sci. 18 (1980) 1129. [18] R.L. Schi€man, A thermoelastic theory of consolidation, Envir. Geophys. Heat Transfer 4 (1971) 78. [19] A. Søu_zalec, J. Paluszy nski, Thermal and moisture phenomena in heat-resisting concrete, Archiv of Thermodynamics 11 (3/4) (1990) 135±157. [20] H. Ziegler, A modi®cation of Prager's hardening rule, Q. Appl. Math 17 (1) (1959) 55±56. [21] E. Melan, Zur Plastizit at der r aumliichen Kontinuums, Ing. Archiv 9 (1938) 116±126.

A. Søu_zalec / International Journal of Engineering Science 37 (1999) 1985±2005

2005

[22] Prager, W. The general theory of limit design, Proceedings of the Eighth International Congress on Appl. Mech, Istanbul, 2 (1955) 65±72. [23] Taylor, G., Quinney, H., The plastic distortion of metals, Phil. Trans. R. Soc. Ser. A 230 (1931) 323±362. [24] A. Søu_zalec, Introduction to Nonlinear Thermomechanics, Springer, Berlin, 1992.