A comparative study of Biot’s theory and the linear Theory of Porous Media for wave propagation problems

Wave propagation in porous media is an important topic for example in geomechanics or oilindustry. Especially due to the interplay of the solid skeleton with the fluid the so-called second compressional wave appears. The existence of this wave is reported in the literature not only for Biot’s theory (BT) but also for theoretical approaches based on the Theory of Porous Media (TPM – mixture theory extended by the concept of volume fractions). Assuming a geometrically linear description (small displacements and small deformation gradients) and linear constitutive equations (Hooke’s law) the governing equations are derived for both theories, BT and the TPM, respectively. In both cases, the solid displacements and the pore pressure are the primary unknowns. Note that this is only possible in the Laplace domain leading to the same structure of the coupled differential equations for both approaches. But the differential equations arising in BT and TPM possess different coefficients with different physical interpretations. Correlating these coefficients to each other leads to the well-known problem of Biot’s ‘‘apparent mass density’’. Furthermore, some inconsistencies are observed if Biot’s stress coefficient is correlated to the structure arising in TPM. In addition to the comparison of the governing equations and the identification of the model parameters, the displacement and pressure solutions of both theories are presented for a one-dimensional column. The results show good agreement between both approaches in case of incompressible constituents whereas in case of compressible constituents large differences appear.


Introduction
For a wide range of fluid infiltrated materials, such as water saturated soils, oil impregnated rocks, or air filled foams, the elastic as well as a viscoelastic description of the material behavior is a crude approximation for the investigation of wave propagation in such media. Due to their porosity and due to the interaction of the skeleton and the pore content, a different theoretical approach is necessary to describe the observed effects like the second compressional wave.
A historical review on the subject of multiphase continuum mechanics identifies two theories which have been developed and are used nowadays, namely the BT and the TPM. For more details, the reader is directed to the work of de Boer and Ehlers [1], [2] or to the recently published monograph [3]. The early works on porous media are attributed to Fillunger in 1913 [4]. In this paper and in subsequent ones, Fillunger was concerned with the question of buoyancy of barrages. At the same time, a more intuitively theory has been developed by von Terzaghi [5]. These two basic works form the basis for the two different theories used up to day.
In a two-phase material not only each constituent, the solid and the fluid, may be compressible on a microscopic level but also the skeleton itself possesses a structural compressibility. If the compression modulus of one constituent is much larger on the microscale than the compression modulus of the bulk material this constituent is assumed to be materially incompressible. A common example for materially incompressible solid constituent is soil. In this case, the individual grains are much stiffer than the skeleton itself. In the following, the governing equations are given for materially compressible and incompressible constituents, respectively, cf. [15], [16], [19], [20], [21]. Beside these two extreme cases, there exist so-called hybrid models where only one of the constituents is modelled as incompressible and the other one as compressible [19]- [21]. As these intermediate cases can be simply deduced from the equations given next they are not considered here in detail.
Throughout this paper, the summation convention is applied over repeated indices, and Latin indices receive the values 1, 2 and 1, 2, 3 in two-dimensions (2-d) and three-dimensions (3-d), respectively. Commas () ,i denote spatial derivatives and primes ðÞ 0 denote the material time derivative with respect to the moving skeleton.

Theory of porous media
As the TPM is mostly presented in a general nonlinear fashion [15]- [18], here the focus is given on the linearization process. All given nonlinear equations are formulated with respect to the reference configuration of the solid skeleton, therefore, special indication of the reference coordinate system is skipped. Furthermore, time derivatives are given as material derivatives with respect to the moving skeleton. But according to the subsequent linearization no distinction will be made between the material time derivative and the partial time derivative.

Compressible constituents
In order to describe the two different phases of the material the concept of volume fractions is introduced [l5], [l6], [18]. Therefore, the given volume element V is divided in two fractions V S and V F occupied by the solid skeleton (index S) and the interstitial fluid (index F), respectively. If the whole space is filled with matter, the saturation condition requires V¼V S + V F . The volume fraction of each constituent is defined by The partial densities q k of both constituents relate the mass element of the constituents to the volume element V of the mixture while the effective densities q kR relate the same element of mass to the volume element occupied by the constituent. Therefore, the partial densities are obtained by the product of the volume fraction and the respective effective density, Changes of the partial density are therefore possible due to changes of the effective density and of the volume fraction, i.e., the material itself as well as its porous structure allow for a compressibility.
Within the general framework of compressible constituents, the effective densities are state variables and the volume fractions are internal variables [15], [l6], [l8] which may be transformed to state variables under certain conditions. For an elastic and materially compressible solid phase, a nonlinear representation of the solid volume fraction is given by Diebels [21] as The current value of the volume fraction n S depends on the solid deformation gradient F S and the initial solid volume fraction n S 0 . Expression (3) is derived from an evolution equation for the volume fraction by appropriate assumptions. Due to the assumption of a geometrically linear description the determinant of the deformation gradient is approximated by as shown in [26]. The divergence of the solid displacement u i;i gives the linear expression for the volumetric strain. Inserting the linearized format of the deformation gradient (4) into the expression for the volume fraction (3) and a subsequent Taylor series expansion yields The balance equations of momentum of a two-phase continuum give the basis for the theoretical description within the TPM. They can either be given for both constituents separately or one of the individual balances may be replaced by the balance of momentum of the mixture as discussed in detail in [27]. In the present contribution, the mixture balance of momentum is used in combination with the fluid momentum balance. The balance equations of momentum for the two-phase mixture read (i) for the mixture (ii) for the fluid In Eqs. (6) and (7), w i denotes the seepage velocity defined as the relative velocity of the fluid with respect to the deforming solid skeleton. The stress tensor is given by T k ij with k ¼ S for the solid skeleton and k¼F for the fluid, respectively. The body force density in the fluid and in the solid is n F q FR b F i , and n S q SR b S i , respectively. The force density p Ã i results from a momentum production representing the interaction between both constituents. Therefore, it is obviously not present in the equation for the mixture (6).
Keeping in mind a linear version of the theory, the convective terms on the left hand sides of Eqs. (6) and (7) are of second order and will consequently be neglected. Furthermore, inserting the series expansion of the solid volume fraction (5) into the balance of momentum (6), in a consequent linearization remains only the constant part n S 0 in the final linear equation. Subsequently, due to the saturation condition 1 ¼ n S þ n F , both volume fractions are assumed to be constant within the balances of momenta This corresponds to the case called ''frozen volume fractions'' by Bowen [16].
Gathering all the linearizations formulated above the linear balances of momentum are (i) for the mixture (ii) for the fluid where no distinction between the partial time derivative and the material time derivative has to be made. In Eq. (9), the bulk body force qb i with q ¼ n S 0 q SR þ n F 0 q FR is introduced as an abbreviation for the sum of the solid and fluid body force. Additionally, the balance of moment of momentum is fulfilled if the stress tensors are symmetric.
Furthermore, constitutive assumptions must be specified which link the stress tensors and the momentum production term to kinematic quantities. Neglecting the fluid extra-stress [28], the stress tensor of the fluid is governed by the pore pressure p, where d ij denotes the Kronecker delta. Furthermore, the viscosity of the fluid is taken into account by the momentum production or by the interaction force between the solid and the fluid which is given by the linear relation [18] with the permeability j T ; ð Þ T1 ¼ TPM. This permeability depends on the intrinsic permeability k S and on the fluid viscosity l F according to the relation j T ¼ k S =l F (see, e.g., [28]). In the balance of momentum for the fluid (7) or (10) the stresses and the interaction forces combine to According to the choice (11) and (12), the viscous properties of the fluid are modelled by the momentum exchange term (12) while the fluid extra stress is neglected [28]. For the solid skeleton Hooke's law is taken into account assuming a linear elastic behavior. Hence, with the extra stress the stress tensor of the solid skeleton is given by [21] if a linear strain-displacement relation e ij ¼ 1=2ðu i;j þ u j;i Þ holds. The shear modulus G and the compression modulus K are introduced in the constitutive equations. These material constants refer to the bulk material and, therefore, the compression modulus includes also the compressibility of the skeleton structure. Furthermore, the state variable z S was introduced in [21] to separate effects related to material and structural compressibilities, respectively. In formulating (15), it is additionally assumed that the free Helmholtz energy is independent of the volume fractions. Subsequently, due to this assumption, the configuration pressure defined in [21] vanishes. In the balance of momentum for the mixture the divergence of the total stress tensor is needed, i.e., the combination of the solid and fluid stress tensor. Under the assumptions made above the divergence of the total stress is obtained as assuming constant volume fractions according to the linearization (8).
The balance of mass of the solid is formulated for the partial density q S ¼ n S q SR and is split into two parts by the introduction of the arbitrary function 0 z S 1 [21] For arbitrary values of z S , (17) is fulfilled if each part of the sum is equal to zero, Based on a micro mechanical investigation, Diebels [21] has proposed the dependence z S ¼ 1 À K S =K SR relating z S to the compression modulus of the structure K S and the compression modulus of the solid grains K SR . He showed that this choice is thermodynamically admissible and that in the case of an incompressible solid skeleton the limit z S ¼ 1 transforms (17) into the well known volume balance n S ¼ n S 0 det F À1 S . Finally, an equation of state for the fluid must be prescribed because in (11) for the fluid extra stress tensor no constitutive assumption was given. Within the framework of a linear theory, the simplest case of the ideal gas equation [29] is applied, with the reference density q FR 0 and the reference pore pressure (static pressure) p 0 . The alternative second expression uses the absolute temperature # and the specific gas constant R. More complex laws to describe the volumetric behavior of the fluid could be included here, however, the linearization neglects additional effects.
With these preliminaries the continuity equation for the fluid can be formulated. In general, this equation reads using the seepage velocity w i ¼ v i À u 0 i instead of the fluid velocity v i . Introducing the material time derivative and combining (20) with the definition of the partial density (2) and with the saturation condition in the form ðn F Þ 0 ¼ Àðn S Þ 0 ¼ z S n S u 0 i;i yields the following representation [21]: Equation (22) is the nonlinear form of the continuity equation of the fluid with respect to the moving solid reference system. This equation is linearized by a formal Taylor series, where in the first term on the left-hand side the gas equation (19) is substituted. In the second and third term, the density is multiplicated with the divergence of the solid velocity or seepage velocity, respectively. Consequently, only the constant factor q FR 0 of the series expansion of the density is used. Additionally, according to Eq. (8), constant volume fractions are introduced leading to the following linearized form of the continuity equation of the fluid: Gathering all above given linearizations the following set of coupled differential equations is obtained from the balance equations: The primary variables in (24) are the solid displacement u i , the seepage velocity w i , and the pore pressure p. Note that in Eqs. (24), due to the linearization, constant densities q SR 0 ; q FR 0 and are used with the exception of the body force terms (Boussinesq approximation), where a linear approximation of the density is inserted.
From a physical point of view it is sufficient to describe the problem with only two primary variables, namely the solid displacement u i and pore pressure p instead of three variables as discussed in [24], [25].
In the quasi-static case, i.e., u 00 i % 0; w 00 i % 0, the balance of momentum of the fluid (24.2) can be rearranged to express the seepage velocity in terms of the pore pressure gradient. In this case, Darcy's law is obtained. Inserting this expression into Eqs. (24.1) and (24.3) eliminates the seepage velocity as primary variable from the set of the governing equations. Since in the dynamic case w i is given as time derivative in (24.2), this procedure is only possible in the Laplace domain. Before the Laplace transformation can be performed the following assumptions are made: -All initial conditions vanish, i.e., -The pore pressure p is assumed to be the excess pressure relative to the static pressure p 0 . Therefore, the initial conditions for the pore pressure also vanish, Taking these assumptions into account, the transformed Eqs. (24) are whereð Þ indicates the Laplace transform and s is the complex Laplace variable.
Rearranging the Laplace transformed balance of momentum for the fluid (27.2) the seepage velocity is obtained, In (28) the complex valued abbreviation is introduced. Eliminating the seepage velocityŵ w i from the remaining balances (27.1) and (27.3) by use of (28), finally, the balance of momentum for the mixture and the mass balance of the fluid are achieved. These operations establish a system of coupled partial differential equations for the unknowns solid displacementû u i and pore pressurep p, An analytical representation of Eqs. (32) and (33) in time domain is only possible for a constant value b T . This is only achieved in the limit j T ! 1, i.e., l F ! 0. Consequently, the interaction force p Ã i between the solid and the fluid is proportional to the pore pressure p Ã i % pn F ;i and the influence of the seepage velocity on the momentum exchange vanishes. Evidently, this is only valid under equilibrium conditions where no fluid motion takes place.

Incompressible constituents
Naturally, the balances of momentum (6) and (7) are not changed due to the assumption of incompressible constituents. So, the linearization process is performed as shown in the previous section. Also, caused by linearization, the volume fractions n S and n F are assumed to be constant within the balance equations. So, the linearized balances of momentum for the mixture (9) and for the fluid (10) are valid also in case of incompressible constituents.
On the other hand, the continuity equation of the solid (17) reduces to a balance of volume. As stated above, the incompressible case is included in the general framework by the choice z S ¼ 1. The physical interpretation is obviously a constant density q SR resulting in the wellknown balance of volume [19]- [21] n S À Á 0 þn S u i;i ¼ 0: Assuming both constituents as materially incompressible and inserting the assumptions q FR ¼ const. and z S ¼ 1 into the nonlinear form of the continuity equation of the fluid (22) yields The constitutive equations for the incompressible solid and incompressible fluid can also easily be achieved. The stress tensor of the fluid (11) and the interaction force (12) are not changed yielding the well known principle of effective stress, but note that the pore pressure becomes a Lagrangian multiplier in this case which ensures the assumption of constant density. There is no longer an equation of state linking the density to the pressure. Finally, the divergence of the total stress is obtained by these assumptions in combination with the saturation condition n S þ n F ¼ 1, As in the compressible case, the incompressible model results in three equations for the three variables solid displacement u i , pore pressure p, and the seepage velocity w i , 2) leads to the same expression as given in (28). Eliminating the seepage velocity from the balance of momentum (37.1) and from the balance of volume (37.3) results in the set of coupled differential equations for the unknowns solid displacementû u i and pore pressurep p As in the compressible case, an analytical representation in time domain is only possible for j T ! 1.

Biot's theory
In this section, Biot's model of a poroelastic continuum is presented. Using different notation for all variables which are not definitely the same as in the TPM approach allows for an a posteriori comparison between both approaches including an identification of the individual terms and of their physical meaning.

Compressible constituents
Following Biot's approach to model the behavior of porous media, an elastic skeleton with a statistical distribution of interconnected pores is considered [7]. The porosity is denoted by where V F is the volume of the interconnected pore space contained in a sample of bulk volume V. Contrary to these connected pores the sealed pores are considered as parts of the solid. Therefore, / ¼ n F is only valid if all pores are interconnected. As mentioned above, full saturation is assumed leading to V ¼ V F þ V S with V S the volume of the solid including sealed pores.
If the constitutive equations are formulated for the elastic solid and the viscous interstitial fluid, the following partial stress formulation is obtained [7]: using the solid strain e S ij and the volumetric strain of the fluid e F kk . The elastic skeleton is assumed to be isotropic and homogeneous. Its elastic behavior is governed by the two material constants compression modulus K and shear modulus G. The coupling between the solid and the fluid is characterized by two additional parameters Q and R B . In these equations, the sign conventions for stress and strain follow that of elasticity, namely, tensile stress and strain are denoted positive. Therefore, in Eq. (41.2) the pore pressure p is the negative hydrostatic stress in the fluid r F .
An alternative representation of the constitutive equations (41) is used in Biot's earlier work [6]. There, the total stress is obtained. In Eq. (42), the solid strain is replaced by the common linear strain-displacement relation e S ij ¼ 1=2ðu i;j þ u j;i Þ. In addition to the total stress r ij , the variation of fluid volume f per unit reference volume is introduced as a second constitutive equation, The variation of fluid content f is governed by the mass balance, i.e., by the continuity equation with the specific flux q i . This flux is identified with the filter velocity q i ¼ /w i . A time integrated form of (44) identifies f as a kind of volumetric strain describing the motion of the fluid relative to the solid as discussed in [30]. Additional to (44), the balance of momentum for the bulk material must be fulfilled. The dynamic equilibrium is given by with the bulk body force per unit volume f i ¼ ð1 À /Þf S i þ /f F i , and the bulk density q ¼ q S ð1 À /Þ þ /q F . For these densities the subscript is used instead of the superscript in order to distinct them from the quantities arising in TPM.
Furthermore, the fluid transport in the interstitial space in terms of the specific flux q i is modelled with a generalized Darcy's law which is given constitutively. Here, j B ¼ k S =l F denotes the permeability defined by the intrinsic permeability k S and the viscosity of the fluid l F . The superscript B is chosen to distinct this permeability from the permeability in the TPM. Furthermore, f F i is the fluid body force per unit volume. In Eq. (46), an additional density, the so-called apparent mass density q a was introduced by Biot [9]. The apparent mass describes the dynamic interaction between fluid phase and solid skeleton. Typically it is written as q a ¼ C/q F where C is a factor depending on the geometry of the pores and the frequency of excitation. At low frequencies, Bonnet and Auriault [31] measured C¼0.66 for a sphere assembly of glass beads. In higher frequency ranges, a certain functional dependence of C on the frequency has been proposed based on conceptual porosity structures [10], [31]. In the following, C¼0.66 is assumed.
The equation of motion for the poroelastic model is obtained from the above balance laws and constitutive equations. As shown in [24], it is sufficient to use the solid displacements u i and the pore pressure p as basic variables to describe a poroelastic continuum. Therefore, the above equations are reduced to these two primary variables as was already done for the TPM model. First, Darcy's law (46) is transformed into Laplace domain and rearranged to obtain In Eq. (47), the abbreviation b B is defined for further usage, and the superscript B is chosen to distinct from the similar abbreviation in TPM. Moreover, as in TPM vanishing initial conditions for u i ; p, and w i are assumed. Now, the final set of differential equations for the displacementû u i and the pore pressurep p is obtained by inserting the constitutive eqs. (42) and (43) in the Laplace transformed dynamic equilibrium (45) and into the continuity equation (44). Taking into accountŵ w i according to (47) leads to This set of coupled differential equations describes the behavior of a poroelastic continuum on the basis of BT. As in TPM, an analytical representation in time domain is only possible for j B ! 1. This case represents negligible friction between solid and interstitial fluid.

Incompressible constituents
To find the respective governing equations for incompressible constituents the material parameters a and R B have to be rewritten in a different way. Considerations of constitutive relations at micro mechanical level as given in [30] lead to a more rational model for this purpose, a ¼ 1 À K K S and ð50:1Þ where K S denotes the compression modulus of the solid grains and K F the compression modulus of the fluid. Based on these expressions materially incompressible behavior of the constituents may be described by the above given constitutive assumptions. Note that material incompressibility means that the compression modulus of each individual constituent is much larger than the one of the bulk material. The respective conditions are [30] The corresponding limit process shows [30] a % 1 and R B ! 1: ð52Þ According to (43) f % u i;i , i.e., the fluid is influenced only by the solid volumetric strain. With these considerations at hand the set of governing differential equations reduces to

Analytical solution
For the subsequent comparison of both theories not only the governing equations but also numerical results will be compared. For this comparison, a semi-analytical solution for a onedimensional fluid-saturated porous continuum is considered. Several of such solutions are found in the literature, e.g., for an infinitely long column consisting of incompressible constituents modelled with an incompressible version of TPM in [32] or for a finite column consisting of compressible constituents modelled with BT [33]. The latter is also applied here, in addition a solution for incompressible constituents is given in the following. Formally, Biot's theory is used as the starting point, however, because the governing equations differ only in the coefficients but not in the type of the differential operator the same solution is valid for the TPM.
A one dimensional column of length ' as sketched in Fig. 1 is considered. The side walls and the bottom are assumed to be rigid, frictionless, and impermeable. Hence, the displacements normal to the surface are blocked and, on the other hand, the column is free to slide parallel to the wall. At the top, the stress r y ðy ¼ 'Þ ¼ ÀP 0 f ðtÞ is prescribed as a function of time while the pressure vanishes, i.e., pðy ¼ 'Þ ¼ 0 N/m 2 . Therefore, the surface is ideally drained. At the bottom the column is fixed and impermeable, i.e., u y ( y¼0)¼0 m and q y ( y¼0)¼0 m/s. Due to these restrictions only the displacement component u y in vertical direction and the pore pressure p remain as degrees of freedom.
Hence, the governing set of differential eqs. (53) and (54) is reduced to two scalar-valued coupled ordinary differential equations in the Laplace domain, Eû u y;yy À ð1 À b B Þp p ;y À s 2 ðq À bq F Þû u y ¼ 0; ð55Þ b B sq Fp p ;yy À ð1 À b B Þsû u y;y ¼ 0; with Young's modulus E ¼ K þ 4 3 G. The boundary conditions are û u y ðy ¼ 0Þ ¼ 0 m;q q y ðy ¼ 0Þ ¼ 0 m s and The applied load on the top of the column is an impulse function f ðtÞ ¼ dðtÞ in time, with dðtÞ denoting the Dirac distribution. In addition, zero initial conditions are assumed. Due to the neglected body forces this is a system of homogeneous ordinary differential equations with inhomogeneous boundary conditions which can be solved by the exponential ansatẑ u u y ðyÞ ¼ Ue ksy ;p pðyÞ ¼ Pe ksy : Inserting the ansatz functions (58) in eqs. (55) and (56) yields an eigenvalue problem for k, with the characteristic equation The characteristic Eq. (60) has the following four complex roots: leading to the complete solution of the homogeneous problem in the form u u y ðyÞ ¼ U 1 e k 1 sy þ U 2 e Àk 1 sy þ U 3 þ yU 4 ;p pðyÞ ¼ P 1 e k 1 sy þ P 2 e Àk 1 sy þ P 3 þ yP 4 : ð62Þ The eight unknown constants U i and P i ; i ¼ 1; . . . ; 4, cannot be determined by the four boundary conditions (58) alone. Also none of the complex roots can be excluded due to physical reasons. But the eigenvector of the system (59) gives the relation and into the one-dimensional form of Darcy's law (46), q q y ðs; The remaining four constants P i can be fit to the four boundary conditions (57). This leads to four equations for four unknowns which can be solved, preferably, with the aid of computer algebra.
Consequently, the solution for the displacement and the pore pressure is achieved by inserting these coefficients in the ansatz functions (62). This results in e Àk 1 sð'þyÞ À e Àk 1 sð'ÀyÞ 1 þ e À2k 1 s' ; ð66Þ The corresponding stress and flux are calculated from the constitutive Eq. (64) and from Darcy's law (65), respectively. As we are dealing with a linear problem the superposition principle is valid. Therefore, solutions for different load cases, e.g., a pressure load or a prescribed displacement, can be achieved by the same procedure. Note, due to the dependence of b B on the Laplace parameter s, the roots k i , are dependent on s. Therefore, an analytical inverse Laplace transform of the above given solutions is in general not possible. However, if the viscosity l F of the fluid tends to zero the damping due to the relative motion of the fluid and the solid may be neglected and, therefore, In this case, an analytical inverse Laplace transform can be found (for the case of compressible constituents see [33]). The same limit for j T can be taken in case of TPM.
For an arbitrary value of j B a numerical inverse Laplace transformation is necessary. A number of methods is available in the literature, and the advantages and disadvantages have been studied, e.g., in [34] or [35]. In the present case, the response in the time domain can be calculated with the convolution integral where L À1 is the inverse Laplace transform operator. Because one function in the convolution integral (69) is only available in the Laplace domain and the other function in the time domain, it is preferable to take the 'Convolution Quadrature Method' proposed by Lubich [36]. This method approximates the convolution integral (69) numerically by the quadrature formula x nÀkû u y ; Dt ð Þf ðkDtÞ; n ¼ 0; 1; . . . ; N: Details concerning the integration weights x nÀk ðû u y ; DtÞ and the used parameters can be found in Appendix A. In the following, the time dependent responses are evaluated with this method, choosing a backward differentiation formula of order 2 (BDF 2) as the underlying multi-step method.

Comparison of both theories
Both theories, namely BT and TPM, describe the same physical behavior of fluid saturated porous media. In the Sects. 2 and 3, the governing equations for each theory were given for compressible as well as incompressible constituents. Now, a comparison of these equations is performed in order to identify the physical interpretation of the parameters of both approaches and to show whether there are discrepancies between the theories even if the underlying structure of the governing equations is the same. Preliminary for this comparison is the (evident) assumption that the independent variables solid displacement u i and pore pressure p have the same physical meaning in both approaches. Furthermore, from the constitutive equations, Eqs. (14) and (41.1), of the solid extra stress it can be concluded that the shear modulus G and the compression modulus K are macroscopic moduli valid for the porous skeleton, i.e., solid material including its structure. In order to compare the other parameters arising in the governing equations the model equations are recalled, first, for incompressible constituents: -Displacement, TPM (38) and BT (53) (39) and BT (54) Comparing the densities it is found that i.e., in BT the densities correspond to the effective densities introduced in TPM. Furthermore, the porosity can be identified with the initial fluid volume fraction n F 0 ¼ / which is also indicated by their initial definitions (1) and (40). However, as mentioned in Sect. 3, this conclusion requires that all pores are interconnected.
With theses identifications in mind the body forces can be compared. While in TPM the body force density b k i is defined as force per mass, in BT f k i is defined as force per volume (k¼S, F). According to the definition of the densities the identification f k i ¼ q kR b k i is obtained. Clearly, this difference is due to the different underlying definitions of the body forces.
Finally, the governing equations of both theories become identical if b B ¼ b T holds. Comparing the definitions (47) and (29) the following identity has to hold: Before evaluating (72) the different expressions for the permeabilities have to be discussed. In BT, the permeability is defined by the quotient of the intrinsic permeability k S and the viscosity l F of the fluid, i.e., j B ¼ k S =l F . Note that the intrinsic permeability k S describes only the pore structure. As discussed in Sect. 2, the permeability j T in TPM has the same physical meaning as j B . Taking this equality into account from Eq. (72) it follows that However, this can only be achieved if the apparent mass density q a vanishes, which is in accordance with Ehlers and Kubik [22]. Therefore, for the incompressible case it can be concluded that the linearized governing equations of both theories are identical if condition (73) holds. This equivalence of both approaches is also verified by the one-dimensional example. In Fig. 2, the displacement at the top of the 1-d column (see Fig. 1) is plotted versus time for both theories. The loading is assumed to be r y ðy ¼ 'Þ ¼ À1 N/m 2 HðtÞ, i.e., it is kept constant over the complete observation time t ! 0. The material data are those of a rock as given in Table 1. The agreement of both results is perfect as expected due to the identification of the parameters.
As stated above, the perfect agreement of the results shown in Fig. 2 can only be achieved for vanishing apparent mass density. So, the question is which influence has the apparent mass density in BT and what differences appear if condition (73) is not fulfilled? To answer this, numerical tests have shown that the apparent mass density has no influence on the results for the given set of material data according to Table 1. However, if the permeability j ¼ j B ¼ j T is increased or if the viscosity of the fluid is decreased, differences appear depending on the apparent mass. As already reported by Schanz and Cheng [33], the second slow compressional wave becomes visible for increased permeabilities.  To study this effect, the pressure 5 m behind excitation of a long column (' ¼ 1000 m) is depicted versus time for different values of j in Fig. 3. The smallest value of the permeability j ¼ 1:9 Á 10 À10 m 4 /(Ns) represents the realistic case. In this case, both graphs show no difference. Furthermore, two effects are observed. Firstly, an initial jump indicates that the fast compressional wave travels with infinite wave speed. This is due to the incompressible model and is mathematically found in the analytical solution by the zeros of the eigenvalues (61). Secondly, by increasing j a second jump becomes visible corresponding to the highly damped second compressional wave [33]. This slow compressional wave is of negative amplitude because it represents the out-of-phase movement of fluid and solid. The arrival time and, therefore, the wave speed of this wave depend on the chosen value of the apparent mass density. Additionally, the pressure level in case of q a 0 is smaller than that of q a ¼ 0:66/q F . These results show that the interaction between solid and fluid is influenced by the apparent mass density according to the structure of the BT, where q a was introduced as 'mass coupling parameter' [9]. Therefore, it is concluded that TPM and BT model the interaction between solid and fluid in a different way, however, for the test data this has no significant influence.
In the next step, the comparison of the compressible models is performed. For this purpose, the governing equations are recalled: F i Àf f i ; -for the pore pressure, TPM (33) and Biot (49) Obviously, if the equivalences found in case of incompressible constituents are taken into account only two additional parameters have to be identified. Firstly, the comparison of the pressure term yields the not expected result On the left-hand side of Eq. (74), the material parameter R B (50.2) depending on the fluid as well as on the solid properties is compared with the gas constant R and the absolute temperature # on the right hand side characterizing only the fluid. The same appears in identifying the last parameter a. If both expressions are the same in both theories it must hold that Inserting the identification given by Diebels [21] into (75), leads to a contradiction. It should be kept in mind that K and K S have the same physical interpretation in Biot's theory as K S and K SR have in TPM, respectively. Therefore, the constitutive relations derived on micro mechanical models for both theories are not in agreement to each other. This contradiction becomes obvious in Fig. 4 where the displacement at the top of the column is given versus time. While in the incompressible case no differences between the displacement solutions of both approaches are visible, here, large differences appear. To compute these results the specific gas constant of water R¼461,61 Nm/(kgK) at the absolute temperature # ¼ 293 K is used.
Finally, to close the comparison the influence of the apparent mass density in case of compressible constituents is studied. As before in the incompressible case, the pressure in a long column is presented for both theories. The pressure result 5 m behind excitation is depicted versus time for different values of j ¼ j B ¼ j T in Fig. 5. Additionally to the observations in the incompressible case, in the compressible case the fast compressional wave is observed as a first jump. This appears in BT at t % 0:0017 s and in TPM at t % 0:0022 s. The arrival time of the slow compressional wave in BT is, as in the incompressible case, different depending on whether the apparent mass density is zero or not, i.e., at t % 0:004 s or at t % 0:005 s. The arrival time of the slow compressional wave in TPM is t % 0:013 s and does not coincide with one of Biot's models. Furthermore, the pressure level in TPM is much smaller than in BT. This is in accordance with the different displacement levels as shown in Fig. 4.

Conclusions
In the present article, Biot's theory for both compressible and incompressible constituents is recalled. Additionally, under the assumption of a linear theory, the dynamic equations for the mixture theory based Theory of Porous Media (TPM) are presented. Both theories model a two-phase continuum consisting of a porous solid skeleton saturated with an interstitial pore fluid. A comparison of the governing equations as well as wave propagation results for a onedimensional poroelastic column are presented.
Summarizing the results of the comparison, the structure of the governing differential equations in BT and in TPM is the same. So, the wave forms predicted by both theories are equal. As a side effect, it was shown that Darcy's law results naturally from the balance of momentum in the fluid. Due to this, also the generalised version of Darcy's law, which takes the inertia terms into account, is a consequence of the fluid balance of momentum.
In case of incompressible constituents, the model parameters are identified in a way that the governing equations are the same in both theories if the apparent mass density is set to zero. The equivalence between both approaches is also verified numerically. On the other hand, in case of compressible constituents neither the identification procedure nor the numerical results match. This is related to the definition of Biot's stress coefficients a and R B and the identification of the state variable z S in TPM, respectively. These parameters are motivated by micro mechanical considerations in both theoretical approaches. The micro mechanical motivation for the parameters is not mandatory and should be critically surveyed. Furthermore, comparing the constitutive equations for the solid stress tensors of both theories shows that the configuration pressure should be included into the TPM formulation in order to achieve comparable equations.

Convolution Quadrature Method
The ''Convolution Quadrature Method'' developed by Lubich numerically approximates a convolution integral x nÀk ðf f ; DtÞgðkDtÞ; n ¼ 0; 1; . . . ; N; ð77Þ by a quadrature rule whose weights are determined by the Laplace transformed functionf f and a linear multi-step method. This method was originally published in [36] and [37]. An application to the boundary element method may be found in [38]. Here, a brief overview of the method is given. In formula (77), the time t is divided in N equal steps Dt. The weights x n ðDtÞ are the coefficients of the power serieŝ with the complex variable z. The coefficients of a power series are usually calculated with Cauchy's integral formula. After a polar coordinate transformation, this integral is approximated by a trapezoidal rule with L equal steps 2p L . This leads to x n ðf f ; DtÞ ¼ 1 2pi where R is the radius of a circle in the domain of analyticity off f ðzÞ. The function cðzÞ is the quotient of the characteristic polynomials of the underlying multistep method, e.g., for a BDF 2, cðzÞ ¼ 3 2 À 2z þ 1 2 z 2 . The used linear multi-step method must be A (a)-stable and stable at infinity [37]. Experience shows that the BDF 2 is the best choice [39]. Therefore, it is used in all calculations in this paper.
If one assumes that the values off f ðzÞ in (79) are computed with an error bounded by e, then the choice L¼N and R N ¼ ffiffi e p yields an error in x n of size Oð ffiffi e p Þ [36]. Several tests conducted by the first author lead to the conclusion that the parameter e ¼ 10 À10 is the best choice for the kind of functions dealt with in this paper [40]. The assumption L¼N results in N 2 coefficients x n ðDtÞ to be calculated. Due to the exponential function at the end of formula (79) this can be done very fast using the technique of the Fast Fourier Transformation (FFT).