Numerical Modeling of the Dynamic Elastic Modulus of Concrete

This article introduces simulations of theoretical material with controlled properties for the evaluation of the effect of key parameters, as volumetric fractions, elastic properties of each phase and transition zone on the effective dynamic elastic modulus. The accuracy level of classical homogenization models was checked regarding the prediction of dynamic elastic modulus. Numerical simulations were performed with the finite element method for evaluations of the natural frequencies and their correlation with Ed through frequency equations. An acoustic test validated the numerical results and obtained the elastic modulus of concretes and mortars at 0.3, 0.5 and 0.7 water–cement ratios. Hirsch calibrated according to the numerical simulation (x = 0.27) exhibited a realistic behavior for concretes of w/c = 0.3 and 0.5, with a 5% error. However, when the water-to-cement ratio (w/c) was set to 0.7, Young’s modulus displayed a resemblance to the Reuss model, akin to the simulated theoretical triphasic materials, considering matrix, coarse aggregate and a transition zone. Hashin-Shtrikman bounds is not perfectly applied to theoretical biphasic materials under dynamic situations.


Introduction
Dynamic modulus of elasticity (E d ) is related to a very small instantaneous strain, geometrically similar to the initial tangent modulus [1]. Mehta and Monteiro [1] reported E d higher than the initial tangent modulus in approximately 20%, 30% and 40%, for concretes of high, medium and low resistance. Nevertheless, the relation between the static and dynamic moduli is not precise, since no physics law correlates their values to the necessary precision [1][2][3]. The use of E d , which can be accurately determined by acoustic tests [4,5], is more appropriate for analyses of structures subjected to impact loading [1], e.g., vibration serviceability limit state. E d is strongly influenced by mixture parameters, such as the water-cement ratio of mortar (porosity), maturity, the elastic properties of coarse aggregate, the volumetric fraction of the interfacial transition zone (ITZ) and the aggregate-cement proportion [1,3,[6][7][8].
In 2013, Nemat-Nasser and Srivastava [18] reported another limitation of H-S bounds: the approach cannot apply to dynamic situations perfectly. This is a severe limitation, especially for the estimation of dynamic elastic modulus through modal analysis.
Other classical homogenization models have predicted the elastic properties of biphasic materials. Counto [1,11] assumed a prismatic coarse aggregate in the center of the mortar prism, whereas Hansen [1,11] adopted a spherical coarse aggregate in the same region. Hirsch [1,11] proposed a semi empirical model that relates the modulus of elasticity of the concrete to the Young's modulus' of the two phases (aggregate and matrix), their volume fractions and an empirical constant, x. The "x value" is a combination factor obtained by the Voigt and Reuss models and calibrated either experimentally, or via a numerical simulation. Topçu [11] obtained x = 0.3 experimentally. The present study aims to obtain x through numerical simulations.
Sophisticated models have been proposed towards predictions of effects, such as the size and shape of aggregate particles and the influence of porosity and transition zone [19][20][21][22][23][24][25][26]. Hashin and Monteiro [27] described concrete as a three-phase material formed by a matrix with embedded spherical particles and surrounded by a concentric spherical shell that represents an ITZ. According to the authors, the Young's modulus of ITZ is 50% of the elastic modulus of the cement mortar (i.e., E ITZ = E m /2). Duplan et al. [28] designed a very realistic model constituted by aggregates surrounded by a layer of ITZ and a layer of cement paste, while air bubbles are considered mono-sized inclusions with no elastic behavior, which result in a 5% error. Zheng, Zhou and Jin [22,23] proposed a model in which concrete is represented as a three-phase composite material composed of mortar, spherical or elliptical aggregate and an inhomogeneous ITZ, based on a semi-empirical gradient model, dependent on the water-cement ratio, the degree of hydration and the porosity at the ITZ. Pichler et al. [26] predict the dynamic modulus of cement pastes with w/c = 0.35 and 0.60, considering a material constituted by clinker, pores (with spherical shapes) and hydration products (with spherical and/or non-spherical shapes).
Although such models are complex for practical applications, those recent studies have opened up new possibilities for the micromechanical modeling of concrete through numerical simulations, especially simulations that allow the study of dynamic elastic modulus (E d ). The current research focuses on simulations of the concrete as a 2D theoretical material (biphasic or triphasic) with controlled elastic properties of each phase (mortar and aggregate). The purpose is to evaluate the behavior of dynamic effective elastic modulus (E d ) regarding porosity, individual elastic properties of coarse aggregate and/or mortar (i.e., E a , ν a , E m and ν m , respectively), volumetric fraction and the presence of ITZ with known properties (i.e., ITZ thickness and elastic modulus E ITZ ). Further, this approach with controlled properties was compared with classical models, such as Voigt, Reuss, Hirsch, H-S, Hansen and Counto, aiming to show the most accurate classical model for predicting the E d . An experimental validation showed the coherence and limitations of our numerical simulations. Therefore, this paper advances the methodology to predict the dynamic elastic modulus of concrete and contributes to the comprehension of its response under dynamical excitations.

Numerical Simulations
Numerical simulations were performed with the general-purpose finite element software Abaqus release 6.14. A theoretical material was modeled by a 150 mm × 150 mm × 50 mm square with fined mesh 1 mm maximum dimension and a quadrilateral and four-node bilinear element with reduced integration (CPS4R, Figure 1a). Plane stress was considered in this 2D shell structure. This model was chosen due to the low computational cost. The boundary conditions were not considered in the model, representing a free model subject to vibrations. model subject to vibrations.
In this study, modal analysis in Abaqus with the Lanczos solution method was utilized to investigate the dynamic behavior of a free shell structure. The Lanczos solution is an iterative method that enables the calculation of the eigenvectors and eigenvalues of a large sparse symmetric matrix, making it ideal for modal analysis. The shell structure was modeled using a linear elastic material with a Poisson ratio and elastic modulus, and a density of 2400 kg/m 3 was considered. The shell was allowed to vibrate freely to demonstrate all modes of vibration. Abaqus was employed to generate the mass and rigidity matrices and to determine the eigenvalues and eigenvectors for the shell structure.

Frequency Equations of an Isotropic and Homogeneous Material
The relationship between the first natural frequency (bending or shear mode of vibration) and effective elastic properties are established. The theoretical modal analysis determines the shear and flexural first natural frequency for a free vibration of a homogeneous isotropic square structure (see Figure 1). "Shear" frequency can be related to the Effective Dynamic Shear Modulus (Gd), whereas "flexural" frequency is associated with Effective Dynamic Elastic Modulus (Ed) [4,5,29] through "Frequency equations". Table 1 shows an analysis performed with a homogeneous isotropic material of ρ = 2400 kg/m 3 density and variations in the elastic and shear modulus. Simulations determined the frequency equation and confirmed the expected relationships of the "shear" frequency and the "shear" modulus (Equation (1)), as well as the "bending" frequency and the elastic modulus (Equation (2)), see Figure 2. Constants C1 and C2 are dependent on the geometry of the specimens and were determined by a linear regression analysis, C1 = 131.16 and C2 = 87.56 (with R 2 → 1), respectively. In this study, modal analysis in Abaqus with the Lanczos solution method was utilized to investigate the dynamic behavior of a free shell structure. The Lanczos solution is an iterative method that enables the calculation of the eigenvectors and eigenvalues of a large sparse symmetric matrix, making it ideal for modal analysis. The shell structure was modeled using a linear elastic material with a Poisson ratio and elastic modulus, and a density of 2400 kg/m 3 was considered. The shell was allowed to vibrate freely to demonstrate all modes of vibration. Abaqus was employed to generate the mass and rigidity matrices and to determine the eigenvalues and eigenvectors for the shell structure.

Frequency Equations of an Isotropic and Homogeneous Material
The relationship between the first natural frequency (bending or shear mode of vibration) and effective elastic properties are established. The theoretical modal analysis determines the shear and flexural first natural frequency for a free vibration of a homogeneous isotropic square structure (see Figure 1). "Shear" frequency can be related to the Effective Dynamic Shear Modulus (G d ), whereas "flexural" frequency is associated with Effective Dynamic Elastic Modulus (E d ) [4,5,29] through "Frequency equations". Table 1 shows an analysis performed with a homogeneous isotropic material of ρ = 2400 kg/m 3 density and variations in the elastic and shear modulus. Simulations determined the frequency equation and confirmed the expected relationships of the "shear" frequency and the "shear" modulus (Equation (1)), as well as the "bending" frequency and the elastic modulus (Equation (2)), see Figure 2. Constants C 1 and C 2 are dependent on the geometry of the specimens and were determined by a linear regression analysis, C 1 = 131.16 and C 2 = 87.56 (with R 2 → 1), respectively.

Effect of Porosity
The influence of porosity (P) was checked on effective elastic properties G d and E d . We suppose a theoretical material constituted by spherical voids of 3 mm diameter, distributed in a matrix (e.g., Figure 3). Voids were incrementally added to the matrix, i.e., four voids (2 by 2-0.13% porosity) until 1296 voids (36 by 36-40.72% porosity). For all cases, the cement paste had constant properties of ρ p = 2100 kg/m 3 , E p = 20,000 MPa and ν p = 0.2, whereas the voids had ρ v = 1 kg/m 3 and E v = 0,1 MPa (i.e., ρ v /ρ p = 0.0004 and E v /E p = 0.000005). An approximately 0.13% to 40.72% porosity increase decreased E d (Figure 4), since porosity is inversely linked to the material rigidity [8]. Boccacini and Fan [30] reviewed several theoretical models to correlate E d and the porosity of ceramic materials. Hanselmann-Hashin [31,32] (Equation (3)) and Mackenzie [33] (Equation (4)) achieved the best fit with numerical results (see Figure 4 and Table 2). The fit with R 2 close to 0.97 showed coherence of simulations.
where E is the effective elastic modulus, E 0 is the intact elastic modulus under porosity equal to zero (value recovered from E 0 ∼ = 30,000 MPa), P is porosity and A and K are empirical adjustments.
Materials 2023, 16, 3955 5 of 20 where E is the effective elastic modulus, E0 is the intact elastic modulus under porosity equal to zero (value recovered from E0 ≅ 30,000 MPa), P is porosity and A and K are empirical adjustments.

Biphasic Material
A perfectly biphasic material with coarse aggregate (of elastic modulus Ea) pe encrusted in a mortar phase (with elastic modulus Em) without transition zon simulated and the influence of individual elastic properties of each phase ( volumetric fraction of the aggregates and Poisson ratio of each phase on effective E checked. The accuracy of homogenization composite models was measured und above-mentioned conditions.

Biphasic Material
A perfectly biphasic material with coarse aggregate (of elastic modulus E a ) perfectly encrusted in a mortar phase (with elastic modulus E m ) without transition zone was simulated and the influence of individual elastic properties of each phase (E a /E m ), volumetric fraction of the aggregates and Poisson ratio of each phase on effective E d was checked. The accuracy of homogenization composite models was measured under the above-mentioned conditions.

Biphasic Material
A perfectly biphasic material with coarse aggregate (of elastic mod encrusted in a mortar phase (with elastic modulus Em) without tran simulated and the influence of individual elastic properties of eac volumetric fraction of the aggregates and Poisson ratio of each phase o checked. The accuracy of homogenization composite models was me above-mentioned conditions.   Figure 6 shows an increase in Ea increases Ed in the same volumetr constant Em, explained by the "Mixture Rule" [9,13]. The gray area in Fi a ± 5 % error, as of the effective model. This error is used in the verificatio of classical models and is similar to the admissible error adopted by D When the elastic modulus of the mortar (Em = 20,000 MPa) and aggre MPa) were relatively close (i.e., Ea/Em = 1.5), the prediction of Ed for al  Figure 6 shows an increase in E a increases E d in the same volumetric proportion and constant E m , explained by the "Mixture Rule" [9,13]. The gray area in Figure 6 represents a ± 5% error, as of the effective model. This error is used in the verification of the accuracy of classical models and is similar to the admissible error adopted by Duplan et al. [28]. When the elastic modulus of the mortar (E m = 20,000 MPa) and aggregates (E a = 30,000 MPa) were relatively close (i.e., E a /E m = 1.5), the prediction of E d for all classical models was accurate; however, when E a /E m > 2, almost all homogenization models overestimated effective E d (i.e., Voigt, H-S, Counto, mean of H-S and Voigt-Reuss). An exception to this was the Reuss model (series), which underestimated E d for all cases of E a /E m . The Hirsch model calibrated to x = 0.27 accurately represented the phenomenon. The x value is very similar to that obtained experimentally by Topçu for static elastic modulus [11] (i.e., x = 0.3). H-S-Low (and Hansen) also showed a good approximation level, with an almost 5% error. The data indicates that the H-S bounds include the effective biphasic material, while the estimation of the H-S lower bound demonstrates good accuracy in representing the increase of Ea/Em. According to Nilsen and Monteiro [16], the H-S limits diverge because concrete is not considered a biphasic material and a third phase (i.e., transition zone) should be included. This is partially true, since a portion of error due to biphasic H-S does not consider the ITZ phase; nevertheless, Figure 6 shows only cases of a perfectly biphasic material with no ITZ. This divergence of H-S and the biphasic material occurs because the Hashin-Shtrikman bounds were formulated from a static consideration [18], which is not valid here. Therefore, H-S bounds are not perfectly applied to dynamic situations and should be used with caution for modal analysis. Voigt and Reuss bounds always contain a biphasic material, and a simple average of those limits was inferior to 10%, when E a /E m < 3. The same tendency was observed for the H-S average. On the other hand, when E a extrapolates E a /E m > 3, Voigt-Reuss and H-S significantly diverge, which causes a very high error. because concrete is not considered a biphasic material and a third phase (i.e., transition zone) should be included. This is partially true, since a portion of error due to biphasic H-S does not consider the ITZ phase; nevertheless, Figure 6 shows only cases of a perfectly biphasic material with no ITZ. This divergence of H-S and the biphasic material occurs because the Hashin-Shtrikman bounds were formulated from a static consideration [18], which is not valid here. Therefore, H-S bounds are not perfectly applied to dynamic situations and should be used with caution for modal analysis. Voigt and Reuss bounds always contain a biphasic material, and a simple average of those limits was inferior to 10 %, when Ea/Em < 3. The same tendency was observed for the H-S average. On the other hand, when Ea extrapolates Ea/Em > 3, Voigt-Reuss and H-S significantly diverge, which causes a very high error.  In the same volumetric proportions and E a , a decrease in E m (and increase in ration E a /E m ) decreases E d , according to Figure 7. Classical biphasic models presented the same anteriorly postulated tendencies.

Effect of Volumetric Fraction
Variations in the volumetric portion can also influence Ed, with fixed Ea and Em [1,3,34]. Theoretical materials whose volumetric fraction increases progressively (i.e., 196 to 324 particles, according to Figure 8

Effect of Volumetric Fraction
Variations in the volumetric portion can also influence E d , with fixed E a and E m [1,3,34]. Theoretical materials whose volumetric fraction increases progressively (i.e., 196 to 324 particles, according to Figure 8 Figure 9 shows the correlation between volumetric fraction and E d . In this case, an H-S and/or Voigt-Reuss average did not describe the phenomenon successfully, since E a /E m > 3. Hirsch x = 0.27 and H-S-Low accurately captured the phenomenon with a 5% error. All tendencies already mentioned are observed.  Figure 9 shows the correlation between volumetric fraction and Ed. In this case, an H-S and/or Voigt-Reuss average did not describe the phenomenon successfully, since Ea/Em > 3. Hirsch x = 0.27 and H-S-Low accurately captured the phenomenon with a 5% error. All tendencies already mentioned are observed.

Triphasic Composite Material
The models were constituted by coarse aggregate, mortar and a transition zone with known properties, i.e., 0.05 mm [

Triphasic Composite Material
The models were constituted by coarse aggregate, mortar and a transition zone with known properties, i.e., 0.05 mm [35,36]-0.25 mm thickness and E m /2 (E ITZ = 10 GPa) elastic modulus [1,27,35,36]. Two diameters of aggregates, namely φ a = 4.8 mm (256 particles) and φ a = 7.5 mm (625 particles), were used at 50.27% constant volumetric fraction (c a ). For simplicity and reductions in the computational costs, the theoretical material contains no voids. phenomenon when ITZ increased to 0.25 mm and the behavior of the triphasic model was similar to the inferior biphasic model, i.e., Reuss (see Figure 12). Hashin-Shtrikman confirmed again its limitations, i.e., the influence of ITZ and the application of H-S under dynamic situations. Although H-S bounds showed limitations, H-S-Low revealed good accuracy when ITZ was 0.05 mm ( Figure 10) and Ea/Em < 3. A simple average of Voigt-Reuss or H-S showed errors higher than 5%.  According to Figures 10 and 11, the Hirsch model with x = 0.27 best predicted E d for ITZ of 0.05 mm and 0.10 mm, with a 5% error. However, it failed to capture the phenomenon when ITZ increased to 0.25 mm and the behavior of the triphasic model was similar to the inferior biphasic model, i.e., Reuss (see Figure 12).  Hashin-Shtrikman confirmed again its limitations, i.e., the influence of ITZ and the application of H-S under dynamic situations. Although H-S bounds showed limitations, H-S-Low revealed good accuracy when ITZ was 0.05 mm ( Figure 10) and E a /E m < 3. A simple average of Voigt-Reuss or H-S showed errors higher than 5%.
Classic biphasic models can consider ITZ, if applied twice, for the obtaining of a triphasic material. The methodology consists in the obtaining of a preliminary composite constituted by mortar and ITZ and associated with a coarse aggregate for the generation of a triphasic composite. Figure 13 shows the numerical data for φ a = 4.5 mm, E a /E m = 4.5 or 3, E m = 20,000 MPa and variable ITZ. Although ITZ is incorporated, the same tendencies above mentioned are observed. Classic biphasic models can consider ITZ, if applied twice, for the obtaining of triphasic material. The methodology consists in the obtaining of a preliminary composit constituted by mortar and ITZ and associated with a coarse aggregate for the generation of a triphasic composite. Figure 13 shows the numerical data for ϕa = 4.5 mm, Ea/Em = 4. or 3, Em = 20,000 MPa and variable ITZ. Although ITZ is incorporated, the same tendencie above mentioned are observed.

Materials and Mixes
Concrete and mortar samples were produced for the validation of the num study. Portland cement type III of 3.10 g/cm 3 density, specified according t classification of ASTM C150 [37] and ABNT NBR 5733:1991 [38], was used in the sam Natural sand of 2.19 fineness modulus [39] and 2.65 g/cm 3 density and crushed diab 3.30 g/cm 3 density and 12.5 mm nominal maximum sizes were used as aggregates. grading curves are shown in Figure 14.

Materials and Mixes
Concrete and mortar samples were produced for the validation of the numerical study. Portland cement type III of 3.10 g/cm 3 density, specified according to the classification of ASTM C150 [37] and ABNT NBR 5733:1991 [38], was used in the samples. Natural sand of 2.19 fineness modulus [39] and 2.65 g/cm 3 density and crushed diabase of 3.30 g/cm 3 density and 12.5 mm nominal maximum sizes were used as aggregates. Their grading curves are shown in Figure 14.
The proportions (in weight) of mortar and concrete mixtures are shown in Table 3. Water cement ratios, amounts of sand, coarse aggregates and superplasticizers were varied among the mixtures. The production of mortars and/or concrete samples consisted in (i) mix the cement and aggregates with no water addition, then (ii) water was added and the mixing procedure was performed until the material had shown a homogeneous appearance. After hardening, the molds were removed and immediately stored in a moist chamber for curing until the age of testing. The proportions (in weight) of mortar and concrete mixtures are shown in Table 3. Water cement ratios, amounts of sand, coarse aggregates and superplasticizers were varied among the mixtures. The production of mortars and/or concrete samples consisted in (i) mix the cement and aggregates with no water addition, then (ii) water was added and the mixing procedure was performed until the material had shown a homogeneous appearance. After hardening, the molds were removed and immediately stored in a moist chamber for curing until the age of testing.

Specimens
The following samples were produced:

Acoustic Tests and Evaluation of Experimental Ed
Acoustic tests (Figure 15a) were performed in mortars and concrete specimens, according to ASTM C215 [4]. The mass of each sample was measured prior to each acoustic test. The sample was then positioned on steel wires attached to a metal frame (Figure 15a) and an impact was manually applied with a hammer on the surface of the specimen, while

Specimens
The following samples were produced:

Acoustic Tests and Evaluation of Experimental E d
Acoustic tests (Figure 15a) were performed in mortars and concrete specimens, according to ASTM C215 [4]. The mass of each sample was measured prior to each acoustic test. The sample was then positioned on steel wires attached to a metal frame (Figure 15a) and an impact was manually applied with a hammer on the surface of the specimen, while a microphone captured the sound radiated by the specimen's surface in another position. Depending on the vibration mode evaluated (flexural, torsional or longitudinal), different impact and measurement positions were selected, as suggested by ASTM E1876-01 [5] and shown schematically in Figure 16. The first two modes of each type (flexural, torsional and longitudinal) (see Figure 17) were easily identified in the prismatic specimens. However, this study focused on 1st flexural and 1st longitudinal modes (see Figures 15c and 17), i.e., the necessary vibration modes for the obtaining of the dynamic elastic modulus (E d ). After excitation, an onboard sound card of a regular notebook captured the acoustic signal at a 96 kHz acquisition rate. Figure 15b shows a typical time vs. amplitude signal in a prismatic concrete sample tested at 28 days of age, obtained from an impact applied for activating the flexural and longitudinal modes (eccentric excitation). The first 1024-point block was selected, multiplied by a Hanning window and then zero-padded for the obtaining of an 8192-point vector. A fast Fourier transform was applied for the detection of peaks of natural frequencies.  Once experimental f1,flex and f1,long are obtained by acoustic test and exciting flexural and/or longitudinal vibrations modes, the elastic properties are evaluated through frequency equations, according to ASTM E1876 [5]. The sample's dynamic modulus of elasticity (Ed) was obtained from the first flexural mode of vibration (Ed,f) using Equations (5) and (6) for prismatic or cylindrical specimens, respectively: Once experimental f 1,flex and f 1,long are obtained by acoustic test and exciting flexural and/or longitudinal vibrations modes, the elastic properties are evaluated through frequency equations, according to ASTM E1876 [5]. The sample's dynamic modulus of elasticity (E d ) was obtained from the first flexural mode of vibration (E d,f ) using Equations (5) and (6) for prismatic or cylindrical specimens, respectively: where b and t are dimensions of the prismatic cross section, D is the diameter of the cylindrical cross section, L is the specimen length, m is the specimen mass (in kg) and A and B are a correction factor dependent on the Poisson coefficient (ν, generally 0.20, according to Mehta and Monteiro [1]) and dimensions of the specimens (b, t and D), computed by Equations (7) and (8): Similarly, the dynamic modulus of elasticity was obtained from the first longitudinal mode of vibration (E d,l ), in accordance with Equation (9) for prismatic specimens: In this case, K is a correction factor given by Equation (10) and D e is the effective diameter of the bar, equal to D for a cylinder and given by Equation (11) for a prismatic specimen.

Experimental Results and Discussion
The problem of evolution of elastic modulus over time was chosen to demonstrate the observed tendencies presented in the above numerical results. All curves were calibrated by CEB FIB [40] (Equations (12) and (13)) for ages of 1 to 60.
where E(t) is the function that represents the evolution of the elastic modulus over time, E 28 is the modulus of elasticity at 28 days and β E is the contribution of time, according to Equation (13): The input of the biphasic model was the evolution over time of the dynamic elastic modulus of mortar E m (see Figure 18a) and the dynamic elastic modulus of samples of diabase rock E a (see Figure 18b).

Experimental Results and Discussion
The problem of evolution of elastic modulus over time was chosen to demonstrate the observed tendencies presented in the above numerical results. All curves were calibrated by CEB FIB [40] (Equations (12) and (13)) for ages of 1 to 60. 28 ) where E(t) is the function that represents the evolution of the elastic modulus over time, E28 is the modulus of elasticity at 28 days and βE is the contribution of time, according to Equation (13): The input of the biphasic model was the evolution over time of the dynamic elastic modulus of mortar Em (see Figure 18a) and the dynamic elastic modulus of samples of diabase rock Ea (see Figure 18b).    The experimental investigation showed Hirsch (x = 0.27) achieved good accuracy when w/c = 0.3 and 0.5, with errors within the 5% limit and identical to the maximum error observed in numeric simulations. However, its absolute error was slightly higher than 5% for concretes with w/c = 0.7. This phenomenon is explained by a complex interaction of mortar, ITZ, water and coarse aggregate surface at early ages. According to Simeonov and Ahmad [17], a fresh concrete mixture undergoes a process of water migration to the surface of the coarse aggregates and this "water attraction" increases the ITZ thickness. Such a decrease is more pronounced at high water-cement ratios (e.g., w/c = 0.7); the E d (t) curve moves to the Reuss model, which also justified the same tendency of concretes with w/c = 0.5 (within 5% error). Therefore, numerical simulations revealed that, when ITZ increases, the best prediction is achieved by the Reuss model, which is similar to the experimental behavior. On the other hand, this phenomenon is reduced for concretes with little ITZ (e.g., w/c = 0.3) and dynamic elastic modulus is very close to Hirsch (x = 0.27).

Conclusions and Final Remarks
The research involved the development of numerical simulations for a porous material that is biphasic and triphasic, with a focus on investigating the dynamic elastic modulus.
The key findings of the study are: (i) Homogenization composite models such as Reuss and Hirsch (i.e., x = 0 and x = 0. 27) provide good accuracy to predict the dynamic elastic modulus, presenting a maximum error of 5%. Hirsch (x = 0.27) is the most reliable predictor of E d for ITZ values of 0.05 mm and 0.10 mm, while the Reuss model is the best fit for ITZ = 0.25 mm. In experimental studies, a similar trend was observed, where Hirsch (x = 0.27) successfully predicted the E d (t) phenomenon for concretes with low and moderate water-cement ratios (w/c = 0.3 and 0.5), while the Reuss was the only model to predict the mixture with w/c = 0.70 with an error consistently lower than 5%. (ii) In both numerical simulations and experimental validation, Voigt, H-S, and Hansen models consistently produce overestimated values for E d . (iii) Hashin-Shtrikman limits do not contain biphasic theoretical models and cannot be applied to dynamic situations perfectly [18]. (iv) The effective Poisson ratio determined by means of natural frequency measurements demonstrates a notable margin of error, owing to the anisotropy produced by the configuration and orientation of coarse aggregates within the mortar. Furthermore, the arrangement of these coarse aggregates may also affect the measurement, as varying arrangements with equivalent elastic properties and volumetric fractions can result in divergent Poisson ratio values.
The present study has improved our comprehension of how concrete responds to linear dynamic excitations. Through numerical simulations of the theoretical material, the classical models have been refined and evaluated.

Conflicts of Interest:
The authors declare no conflict of interest.