This paper is the third installment from a series of publications pertaining to modeling of porous medium. The first paper, Teagle et al. , derived a coupled set of fluid/structure equations for a porous medium applying asymptotic and homogenization techniques. In  it is shown that there are mainly 3 modes of energy transformation: 1) The first mode is through the connection between the micro and macro structural framework of the porous skeleton, 2) The second is via the viscous boundary layer, and 3) The third interaction is through thermal (entropy) boundary layer. The combination of the viscous boundary layer and how tortuous the porous material is, results in the encapsulation of the fluid medium. This, in turn, changes the apparent mass of the structural medium. Details pertaining to these encapsulating phenomena can be found in the work by Johnson, et al. . In  the concepts of tortuosity, viscous length, and viscous permeability are explained. The thermal interchange is described in Teagle, et al. . In , the mathematical description of the thermal energy dissipated by the thermal boundary layer is explained along with the relaxation process. This thermal exchange changes the acoustic bulk modulus of the porous medium and thus the speed of sound inside the porous layer. These thermal phenomena are then represented by the parameters of thermal length and thermal permeability.
This paper is a continuation of . The equations of motion derived in the aforementioned publication will be presented in a form that is used for calculation of TL for porous material used in the automobile industry. The combination of 1) the change of density due to the viscous effect and 2) the change of acoustic bulk modulus due the thermal exchange will be applied. Additionally, a matrix representation of the wave propagation (for both forward and reflected propagating wave) that incorporates the aforementioned viscous and thermal effects is developed. Each porous layer makes contact with another type of porous medium, elastic panel, or air. This study will establish the correct boundary conditions between mediums in order to run numerical simulations that lead to stable and unique solutions. In order to translate these boundary conditions to the numerical model, interface matrices are developed. In continuation, the global matrix that represents the multi-layered material is assembled and applying the definition of impedance the diffuse field TL is calculated. In the TMM formulation it is assumed that each layer is of infinite extent. A correction applying Green’s function technique is applied.
2.1. Basic Tensor Calculus and Notation
Tensor calculus and concepts from differential geometry are used extensively in this paper. This notation affords a level of abstraction that leads to an efficient explanation of the stresses, strains, and their relationships. Here, the usual notation and represent contravariant and covariant vectors (tensors of order 2). Physical tensors like stress, that are neither contravariant or covariant, are designated here as or (with over bar and void of any ij subscript/superscript or ij subscripts in parenthesis). The relationship between physical tensors, like stress, and their cotravariant and covariant counterparts are the following
where are scaling factors that satisfy the following
It should be noted that is the components of the metric tensor. The Einstein summation convention is used in which covariant index followed by the identical contravariant index is implicitly summed over, thus contracting the order of the tensor. In this paper, if 2 tensors of different order are shown next to each other, this contraction rule is followed. To obtain expression for strains, the gradient of the first term of the asymptotic expansion of the displacement vector, , is required. The expression for the components of the tensor of order 2 is the following
and the expression for strain is given as
The elastic constant, , is a fourth rank tensor. Due to tensor contraction rules the stress tensor, , is of second rank. Additionally, the divergence of the stress tensor, , results in the force field. From Equation (1), the contravariant component of the physical stress tensor, , is
Thus, the physical component of the force field in the jth direction is given by
The use of orthogonal coordinates results in
is a tensor operator that transforms the stress data into a hydrostatic form.
2.2. Fluid-Structure Interaction: Dynamic Equations
This paper applies the definitions and results obtained in . In that publication, the fluid/structure interaction equations are derived via asymptotic and homogenization techniques. A slight modification is applied to the averaged relative displacement of the fluid with respect to the skeleton,. In this paper this equation is represented as
In , Equations (54) and (55) are derived and represent the structural stress, , and internal fluid pressure,. These are rewritten below as Equations (9) and (10)
where and are elastic constant tensor (fourth rank) per their explanations given in .
is the interstitial fluid bulk modulus, details on its derivation can be found in Teagle et al. . Additionally, in  it is shown that
They act as scalar multiples of the identity operator.
Set and therefore
Analyzing the derivations in , represents the proportion of fluid pressure that produces the same strains as the total stress. Additionally, it is important to note that Equations (12a, b) (see Equations (21, 22)) are conditions that make it possible for there to be a strain-energy term for the porous material. This justifies Biot’s assumption of the existence of a potential when he derived his landmark equation, pertaining to porous rocks, for oil exploration.
Equations (44) and (45) from  represent the inertial forces for the structural and fluid components and are rewritten below
Subtracting Equation (13) from (12c) and substituting the definition of, Equation (8), the following expression for the structural portion is obtained
and applying the definition for, an efficient representation of Equation (13) is
Substituting for, in Equation (9), the expression in (10), and in turn taking the divergence the dynamical equation for the structural portion of the foam can be expressed as
In order to derive (17), the definition for, Equation (12b), was used. When, it is easy to verify that and are new Lame’ constants of the elastic portion when the porous material is drained.
Additionally, taking the gradient of Equation (10) followed by a multiplication by, the following expression dealing with pressure is obtained
Combining Equations (14) and (17), the dynamic equations pertaining to the structural/skeleton portion of the porous medium is obtained
Similarly, combining Equations (16) and (18) results in the macro level equations describing the fluid motion
Equation (19) can be expressed in a compact form
and by applying the coupling definition of Q, Equation (20) becomes
2.3. Wave Participation Factor: Eigenvalue Problem
In this section, the study adapts the concepts from Brouard, et al. . The motion of the structural/skeleton and the fluid portion are described by introducing the potential scalar functions and and the potential vector functions, and, such that
Here, the basic definition from mechanics is used where the gradient represent portion of the displacement vector that is purely dilatational or in compression. The curl of and represent any shear motion associated with the displacement. Substitute these definitions into Equations (21) and (22)
Recognize that the operator commutes with i.e. and for any general vector.
Using these relationships, Equations (25) and (26) are rewritten in the following form
Gathering terms corresponding to the gradient operator results in
These equations will give the formulation to solve for the compression waves.
Equations (29) and (30) are represented in matrix form
Applying (34), Equation (33) is turned into an eigenvalue problem
The eigenvalues produce the complex wave number for the compression waves and they will have the following expression
Two sets of eigenvectors are generated the relationship is known as the participation factor. Substituting the eigenvectors into Equation (31) ones gets
These results indicate that there exist two compression waves traveling the porous medium, one is fast and the other one is slower. is important since the number will indicate which wave, whether acoustic or solid, has the most contribution at that particular frequency. For the shear wave we accumulate all terms in Equation (27) and (28) that contain the curl operator
Equation (40) results in, insert this in (39)
Setting, the complex wave number for the shear wave can be given as
and the participation factor is
2.4. Wave Propagation: Porous Medium
This portion will analyze wave propagation in a semi infinite porous medium. Consider Figure 1.
The disturbance wave is traveling down at an angle of incidence of θ in the general xz plane and will impinge the porous medium at the z = 0 level. A portion of this wave will be reflected and some of it will be transmitted to travel through the porous medium until it impinges the next layer, at z = L, of different impedance. A portion of the wave will again be reflected back into the porous medium and the rest will be transmitted to the next layer.
Since there exist 3 wave types, 2 compression (dilatation) waves and 1 shear wave, and since each wave type has 2 waves (forward moving and the reflected wave) there will be 6 variables. These variables are:
structural velocity in the x and z direction
acoustic velocity in the x direction and z direction
= structural stress component perpendicular to z face in the z direction
structural shear component
acoustic stress component
An array of these components will be represented by . The continuity of these variables at the interface between 2 layers will become the boundary conditions that the traveling wave has to satisfy. Using the definition of the participation factor above, the potential functions can be expressed as
Figure 1. Porous medium element: Thickness L, angle of incidence θ.
where, indicates the matching of acoustic and porous waves along the x direction (k is the acoustic wave number and k sin(θ) is the projected acoustic wave number in the x direction). projects the porous wave in the z direction. Note: 1) if is real, there is a propagating wave, 2) if is complex there is a decaying wave in the z direction. The following notations will be used, and . An expression for can now be derived
The six variables can be written as a function of the potential functions and using the above notation:
The following three equations are essential to define the stress components
Utilizing the stress definitions established in section 2.2, the stress components can be written as
Applying the corresponding expressions in (44) into Equations (45)-(54) a set of equations involving, and complex trigonometry equations will be obtained. For example if (44) is applied to (45) the following expression in sues
Similar operations are done to derive expressions for and. The final results are represented in matrix form
The idea behind a transfer matrix is to relate the 6 variable, V(z), at z = L to the conditions at z = 0, i.e. relate the matrix with respect V(L) . The conditions at z = 0 can be written as
Here, is known as the transfer matrix (58)
2.5. Interface Conditions: Uniqueness
A general form for the required boundary conditions between the layers of multi-layered materials is developed. These boundary conditions are developed to ensure uniqueness and numerical stability. Special focus is on deriving interface conditions for the porous material. This section puts concepts introduced by Deresiewicz, et al. , under the context and language set forth in this paper.
The kinetic energy, per unit volume for the two phase system is expressed as
, , and are, , void the portion with tortuosity (without viscous effect). The dynamic equations in Equation (17) assume that the porous material is statistically isotropic. Recall that is the stress in the skeleton and is the stress in the acoustic medium, , is the macro-level identity (hydrostatic) tensor. The strain energy, per unit volume is given by
, = Strain tensor for liquid (solid) (60)
Taking the time derivative of the kinetic and strain energies results in
The following tensorial relationship holds
Combining Equations (61, 62) and the divergence theorem the expression for power is obtained
The portion inside the volume integral represents the traction and inertial work applied to the skeleton and acoustic medium respectively. Equations (21) and (22) are written in shorthand form
The term is the viscous force term due to the interstitial fluid, Β is a viscous transfer function replicating the viscous effects. When (64) and (65) are plugged into Equation (63), the volume portion of the equation becomes
is the dissipation function per unit volume. Taking (66) into consideration the power expression is expressed as
The right hand side of (67a) expresses the rate at which work is done on the material by surface forces, is the normal vector to the surface. Now consider 2 different porous medium with volumes and, boundaries and and where represents the common boundary between porous medium 1 and 2. To establish uniqueness, for each of the two mediums the field quantities will be replaced by a difference representing the possibility of two different solutions but with the same boundary conditions, making the right hand side of Equation (67a) equal to zero. is positive definite, therefore the only way that the difference version will satisfy the zero condition is that everywhere, indicating that the solution is unique. The surface integral in Equation (67a) is partitioned in the following way,
is the element number (67b)
In the common portion,. Combining the two surface integrals results in
For the non-intersecting boundaries (S1 and S2) the boundary conditions will have to be given. At the Sc boundary, continuity is required across the interface to be able to maintain uniqueness. If the skeleton phase is to remain in contact with each other and to maintain the principle of conservation of mass of the acoustic medium it is continuity of the normal relative velocity of the acoustic medium with respect to the skeleton, i.e.
Applying (69) into (68), the Sc integrand becomes
To assure continuity of condition (70), continuity in the following quantities has to be maintained
The non-alignment of the pores can produce a pressure drop across the interface so the continuity condition for pressure p is modeled as, is a coefficient of resistance. The materials used in this study are highly porous so was set equal to zero. This also seemed to be the coefficient that gave the best results after a quick parameter study. If the intersection or interface is between a porous medium and a plate then, or. The continuity conditions for this case (porous-plate) will be the following
The over-bar in conditions (72) and (73) represent quantities in the plate. In case the interface is between porous medium and air, for the layer pertaining of air will be set equal to one, the stresses and will be set equal to zero and. The interface conditions for porous and air becomes
The continuity of was applied to obtain Equation (75).
3. Results and Conclusion
A Simple application of TMM is applied to the simple layer configuration shown in Figure 2. The figure shows a multi layer system where there is a plate that is glued to foam 1 which in turn is connected to foam 2. The parameters pertaining to plate 1, foam 1, and foam 2 are listed in Table 1. This simulation also applied the following parameters that are not listed in the table: Viscous length = 40 μm, Thermal length = 80 μm. The goal of this simulation is to calculate Random Incidence Transmission Loss (TL). This is achieved by applying Finite Size correction Transfer Matrix Method (FTMM)  . The Finite Size Correction is achieved by incorporating a Green’s Function integration to the radiation efficiency. The dynamic equations pertaining to the skeleton and fluid along with the boundary conditions depicted in Equations (71)-(75) are applied in order to simulate how the structural and fluid stresses will change as the acoustic wave travels through the multi-layered material.
Equations (9) and (10) represent the structural equations in its most general form. A detailed Finite Element (FE) can be developed in order to analyze the interplay between macro and micro levels and in turn obtain expressions for P, Q, R, and N; parameters required in Equations (21) and (22). In this study there is not enough information to construct an FE model, the only structural information available are the experimentally obtained bulk modulus, , and the
Figure 2. Multi layer configuration.
Table 1. Parameters pertaining to plate and foam.
shear modulus, N, both pertaining to the porous material. Materials used in the acoustic field usually consist of a highly porous material whose structural portions are considerably stiffer. Due to one of Biot’s  experiment, a jacketed porous medium experiences a hydrostatic pressure pj, but the air inside the porous medium experiences no change in pressure. The definition for porosity
before and after deformation respectively are and.
Assuming highly stiff solid frame will lead to small changes in solid volume, therefore and, is jacketed frame dilation. This information leads to the following equation for
Similarly, a relationship for the fluid dilatation is obtained,. Combining this last expression with that of Equation (76) the following important relationship which relates the fluid and frame dilatation is obtained
Analyzing Equations (17) and (18), the following expressions for structural stress is obtained
For the porous material with relatively very stiff frame, Equation (11) shows that M can be approximated by. In turn, taking the trace of Equation
(78), dividing that result by 3, and setting in Equation (79), the jacketed experiment is simulated
= porous bulk modulus (80)
Equation (77) is applied to obtain (80) and (81). Firstly, from Equation (81), it is deduced that
Plugging Equation (82a, b) into (80) and solving for P results in
The results incorporating Equations (82abc) (indicated as “case 9”) into the simulation are shown in Figure 3. The results are compared to measured
Figure 3. Top: TL comparison (dB); Bottom: TL difference (dB) measured data is the reference.
(“meas”) results and also to the empirical formulas of Delaney and Bazley (“DB”)  . The graph shows that the calculated results come within 0.7 dB.
tr = trace of a tensor, = averaged at the microscopic level
v = Fluid Velocity, p = Fluid Pressure
u = Structural displacement, U = Fluid displacement
= Porosity, n = unit normal pointing into the solid
= Viscosity, = Second Viscosity
= Domain Occupied by Fluid, = Domain Occupied by Structure
= Coefficient of thermal Conductivity, Pr = Prandtl Number =
= Stress Tensor in the Fluid, = Stress Tensor in the structure
= density of the Fluid, = density of the Structure
= Lame Parameters, n = unit normal pointing into the solid
T = Temperature deviation, = Coefficient of thermal Conductivity
= Dynamic Viscous Permeability, = Dynamic Thermal Permeability
Cp = specific heat at constant pressure, Cv = Specific heat at constant volume
= Elastic fourth ranked Contra variant tensor, operates on
= Acoustic bulk Modulus
= asterisk superscript means scaled variable (dimensionless)
 Teagle-Hernandez, A., Ohtmer, O. and Nguyen, D. (2018) Mathematical Modeling of Porous Medium for Sound Absorption Simulations: Application of Multi-Scales and Homogenization. Journal of Applied Mathematics and Physcis, 6, 2705-2717. https://doi.org/10.4236/jamp.2018.612224
 Teagle-Hernandez, A., Ohtmer, O. and Nguyen, D. (2019) Model for Frequency Dependence of Thermal Permeability in Order to Quantify the Effects of Thermal Exchange on Wave Propagation in Multi Layered Porous Medium. Journal of Applied Mathematics and Physics, 7, 1440-1462. https://doi.org/10.4236/jamp.2019.77097
 Nguyen, D. and Teagle-Hernandez, A. (2019) Application of Transfer Matrix Method, Green’s Functions, and Variational Techniques to Predict Diffuse Absorption Coefficients. CalTech-SOCAMS, Pasadena, CA, 27 April 2019.