Ultrasonic Coda Wave Experiment and Simulation of Concrete Damage Process under Uniaxial Compression

Abstract

Using the coda wave interferometry (CWI) method to obtain the ultrasonic coda wave characteristics of loaded concrete is an important method to evaluate the mechanical response of concrete. In this paper, the ultrasonic coda wave characteristics of C40–C70 concrete specimens (four strengths of concrete) under uniaxial compression were tested by laboratory experiments. Furthermore, to clarify the relationship between the internal damage process of concrete and the change rate of coda wave velocity, an ultrasonic coda wave discrete element simulation model combined with digital image processing technology was established. The results show that the coda wave is very sensitive to small changes in the compressive damage to concrete, and the change in coda wave velocity can correspond to the development process of concrete damage. This research is conducive to a better understanding of the complex material behavior of compressive concrete and proves the feasibility of ultrasonic field simulation and processing by using numerical simulation images of concrete damage.

1. Introduction

Concrete ultrasonic nondestructive testing technology refers to the use of signal processing methods to analyze the changes in acoustic parameters such as wave speed, frequency, and amplitude of ultrasonic waves after passing through the concrete medium, and then to evaluate the mechanical properties and damage status of concrete [1]. However, concrete is a nonuniform and anisotropic composite material with a wide range of heterogeneous interfaces inside; complex reflection and refraction phenomena will occur when ultrasonic waves propagate through it [2], making it difficult to measure wave speed changes accurately. However, coda wave interferometry (CWI), developed by seismologists, provides a way to solve this problem [3]. A coda wave refers to the long tail part that extends behind the direct wave in the waveform [4]. Compared with direct waves, it is more sensitive to changes in the internal scattering of the medium, so it can be used to identify small changes in concrete [5,6]. Coda wave interferometry compares the waveforms at different times to determine the travel time offset change law between the waveforms and then obtains the relative change rate in the wave velocity [7]. Related research indicates that its wave velocity measurement accuracy can reach up to 0.001% [8]. Therefore, coda wave interferometry has been introduced into civil engineering as a high-sensitivity monitoring tool, and it has shown great potential and broad application prospects in the health monitoring of concrete structures [9,10]. Presently, many scholars have performed a series of concrete ultrasonic nondestructive testing studies using coda wave interference technology. The research directions include stress monitoring [11,12], damage estimation [13,14], crack positioning imaging [15,16,17], temperature influence, and correction [18,19,20].

This research focuses on the sound velocity response of concrete mechanical behavior under load [21]. Laboratory sample tests are a common method to obtain the ultrasonic coda wave characteristics of concrete under stress. In the process of uniaxial step-by-step compression, the coda wave velocity first slowly increases with the increase in the stress level and then decreases sharply due to the intensified damage [22]. For the multi-cycle loading test, after each unloading, the initial wave velocity of re-compression continues to decrease, reflecting the increase in the degree of internal damage to the concrete [23]. Considering the influence of the number of cycles, the slope of the wave velocity load (stress) function can be used to characterize fatigue damage, which increases with the number of cycles [24]. In the short-time fixed-loading experiment, the wave velocity change increases with increasing fixed-loading time [25]. On the other hand, under the action of tensile stress, the wave velocity decreases with the increase in stress; during the loading and unloading processes, the attenuation of the velocity continues to accumulate [26]. To realize the application of CWI technology in real-world engineering monitoring, the coda wave test was performed on real-size reinforced concrete components, and more attention was given to obtaining local acoustic velocity changes through a reasonable sensor layout. A four-point bending test of reinforced concrete T-beams demonstrated that the coda wave velocity can detect changes in compressive and tensile stresses and accurately locate the position of the neutral axis [27]. A 16-ton truck was used to carry out on-site load measurement of an old box girder bridge that had never been put into use. The results show that the coda wave velocity changes obtained by different sensor combinations can detect the position of the truck to varying degrees [28].

The above experimental results qualitatively explain the characteristics of concrete ultrasonic coda waves, but the experiment volume was limited by the research cost. Numerical simulations can perform a large number of simulation calculations on concrete load conditions and damage conditions by adjusting parameter settings, and can effectively eliminate the adverse effects of environmental factors. It is an effective tool for analyzing and interpreting the CWI features of concrete. Using the two-dimensional finite-difference time-domain model within the framework of the locadiff imaging technique, the extended changes in multiple scattering media characterized by coda wave decorrelation can be calculated [29]. A 2D spectral element method was applied to simulate complex coda signals to study the sensitivity of nonlinear coda interference parameters of samples with local damage, and the research scale was gradually refined: the effective elastic characteristics of the sample, the local crack area, and the effective length of the microcrack [30,31,32]. In order to explore the susceptibility of the coda wave to various damages in the concrete failure process at the micro-scale, a research system of concrete mesostructure generator discrete element calculation of concrete fracture finite-difference simulation of the ultrasonic field was constructed [33].

To date, relevant scholars have achieved good results in research on establishing the relationship between ultrasonic coda wave characteristics and concrete stress or damage. Larose performed linear regression on the data of coda wave relative velocity change and stress change and defined the estimated value of the nonlinear parameter with an estimation accuracy of 5% for this parameter [8]. Schurr used the acoustoelasticity constant and nonlinear elastic parameters to estimate the damage degree and proved the two parameters of the damaged samples are typically multiple times higher than those of the damaged samples [24]. Hafiz used magnitude-squared coherence to evaluate stress change versus ultrasonic waveform change and obtained a test constant as a measure of the sensitivity by linear regression [34]. Clauß derived the correlation between ultrasonic coda wave relative velocity changes and reinforcement strains, covering 90% of the complex load-bearing behavior in reinforced concrete members, and demonstrated a feasible method for evaluating structural state objectives using ultrasonic coda wave characteristic parameters [35].

Existing studies have qualitatively analyzed the characteristics of ultrasonic coda waves in the development process of loaded concrete damage. The research objects are mainly concrete strength (C30–C50) samples commonly used in practical engineering [11,22,23,25]. In this paper, through laboratory experiments and numerical simulations, the correlation between the damage degree of concrete and the ultrasonic coda wave velocity under uniaxial load was studied. A high-strength concrete (C60–C70) ultrasonic coda wave test group was set up to broaden the research scope. Based on the variation curve of concrete coda wave velocity with the percentage of maximum load, a mathematical formula for evaluating concrete damage was fitted. The uniaxial compression simulation model established by this discrete element method can clarify the damage and failure mechanism of concrete. The ultrasonic field simulation model based on digital image processing technology, as an attempt to quantify the relationship between the damage degree and coda wave velocity change rate, provides a feasible platform for further research.

2. Experimental Study

2.1. Test Specimen

This test requires concrete test blocks of 100 mm and 150 mm cubes. The cement material is ordinary Portland cement. The fine aggregate is the finished river sand (2.8 mm), which has not been screened again, and the gradation curve cannot be provided. The coarse aggregate is continuous graded basalt crushed stone (5–10 mm), and the proportion of each particle size will be illustrated in Section 4.1.1. Based on relevant engineering practical experience and specifications, the material mix ratio of the concrete selected in this paper is shown in

Table 1. Material mix ratios.

Concrete Strength Grade	Material Consumption (Weight Ratio)
	Cement	Water	Sand	Pebble	Mineral Powder	Fly Ash	Water Reducing Agent
C40	1.000 (P.O42.5)	0.450	1.360	3.030	/	/	/
C50	1.000 (P.O42.5)	0.350	0.910	2.380	/	/	/
C60	1.000 (P.O52.5)	0.359	1.663	2.298	0.138	0.069	0.032
C70	1.000 (P.O52.5)	0.380	1.620	2.650	0.235	0.117	0.054

.

2.2. Test Device and Coda Wave Test Parameters

The equipment required for the test includes a pressure loading system and an ultrasonic coda wave test system, as shown in Figure 1. The pressure loading system is a YAN-2000 electrohydraulic servo pressure testing machine (Jinan Dongtest Testing Machine Technology Co., Ltd., Jinian, China). The ultrasonic coda wave test system consists of a signal generator, sensor, oscilloscope, signal amplifier, and corresponding signal lines (33522B arbitrary function generator, DSOX2014G dual-channel color digital oscilloscope, Keysight Technologies, USA. Sensor, 25K-P28F, 50K-P28F, 100K-P40F, 0.25P-40F, Guangdong Goworld Co., Ltd., Shantou, China).

During the implementation of the test, in order to reduce the friction between the concrete and the steel plate of the press, lubricating oil was applied to the end face of the test sample. The contact surface of the sample was polished and then coated with couplant. A rubber sleeve was used to fix the ultrasonic probe. Finally, the probe and the sample were squeezed to remove air, so that the two items were in full contact to ensure the testing effect.

The pre-experiment was carried out under no-pressure conditions to determine the relevant parameters of the ultrasonic coda wave test, namely the ultrasonic excitation frequency (50 kHz, 100 kHz, 250 kHz), the ultrasonic excitation voltage (1 V, 5 V, and 10 V), and the ultrasonic test path length (100 mm and 150 mm, which is the size of the concrete sample). Using the control single-variable method, the original waveforms of the corresponding parameter conditions are shown in Figure 2. The degree of ultrasound scattering in concrete is related to its frequency; usually, lower frequencies correspond to weak scattering. Therefore, the original waveform corresponding to the 50 kHz excitation frequency has a larger amplitude and a longer waveform duration, which is convenient for data processing and analysis. Similarly, the excitation voltage and sample size were determined to be 10 V and 100 mm, respectively. The optimal parameters determined above were used in subsequent ultrasonic wake wave tests and numerical simulation studies: ultrasonic excitation frequency (50 kHz), the ultrasonic excitation voltage (10 V), and the ultrasonic test path length (100 mm, i.e., the concrete sample is a 100 × 100 × 100 mm cube).

2.3. Loading Plan and Coda Wave Test Process Control

Uniaxial loading was performed at a rate of 2 kN/s and in steps of 20 kN. The load was fixed for 40 s for every 20 kN increase, and the process continued until the design load of the concrete specimen was reached. To prevent the rupture of the concrete from damaging the sensor probe, 75% of the average axial compressive strength of the concrete of each strength was taken as the maximum load of the corresponding test, as seen in

Table 2. Maximum load of each strength of concrete.

Concrete Strength	Axial Compressive Strength/MPa	Maximum Load/kN
C40	44.1	460
C50	53.87	520
C60	63.8	620
C70	73.6	740

.

The ultrasonic coda wave test used pulse waves as the source signal. Before the start of the test, the parameters of instruments were adjusted. Before loading, the ultrasonic coda wave test system was started, and the waveform $u_0(t)$ at the moment of 0 kN was recorded. When the load reached 20 kN, the pressure remained constant for 40 s, and signal $u_1(t)$ was recorded at 30 s. Then, the ultrasonic signal $u_i(t)$$(i\ge 1)$ was recorded for every 20 kN increase in the load. The loading step and data acquisition control are shown in Figure 3 (taking the C40 concrete specimen as an example), and the loading and ultrasonic test process is shown in Figure 4.

Concrete loading and ultrasonic coda wave tests were performed at a stable room temperature (26 °C), and the experimental time was relatively short; thus, it was considered that environmental factors had little influence on the experimental results [8,18,19]. It is worth noting that the creep effect of concrete is unavoidable during the experiment. However, compared with concrete stress damage, the degree of change in ultrasonic coda wave velocity caused by creep is not significant [26]. Therefore, this study did not analyze this aspect.

2.4. Coda Wave Test Data Processing

The existing CWI data processing methods include the doublet technique (DT) and the stretching technique (ST). The calculation results of the DT method vary with the selected time window. In contrast, the ST method can process the entire waveform to obtain more consistent results that are also closer to the actual situation [9]. Therefore, this paper adopts the ST method to analyze the signal of the ultrasonic coda wave. Although the stretching technique eliminates the time dependence, the corresponding velocity change is the average of the entire propagation period, so the wave velocity sensitivity is slightly reduced. To solve this problem, a stepwise CWI method was introduced; that is, only the previous coda wave measurement signal was used as a reference to calculate the stepwise change in acoustic wave velocity [11].

$$u_i^{\epsilon }(t)=u_i(t(1+\epsilon ))$$

(1)

$$C{C}_i(\epsilon )=\frac{{\displaystyle \int u_i^{\epsilon }(t)u_{i-1}(t)dt}}{\sqrt{{\displaystyle \int {(u_i^{\epsilon }(t))}^{2}dt}{\displaystyle \int {(u_{i-1}(t))}^{2}dt}}}$$

(2)

$${\epsilon }_i=\arg \max (C{C}_i(\epsilon ))=\frac{d{v}_i}{v_i}$$

(3)

$$\frac{dv}{v}\cong (1+{\epsilon }_i)(1+{\epsilon }_{i-1})\cdots (1+{\epsilon }_1)-1$$

(4)

The collected waveforms of the two different loads before and after are $u_{i-1}(t)$$(i\ge 1)$ and $u_i(t)$. First, the stretching processing is performed in the time domain of $u_i(t)$ to obtain $u_i^{\epsilon }(t)$, and the stretching coefficient is $\epsilon$. Then, when the cross-correlation coefficient $C{C}_i(\epsilon )$ calculated by the stretched waveform $u_i^{\epsilon }(t)$ and the reference waveform $u_{i-1}(t)$ is largest, the corresponding stretching coefficient ${\epsilon }_i$ is the change rate of the i-th coda wave velocity $d{v}_i/v_i$. Finally, the stepwise wave velocity change value is calculated by the corresponding multiplication to obtain the accumulated relative change of coda wave velocity $dv/v$.

3. Experimental Results and Discussion

3.1. Waveform Cross-Correlation Analysis

The changes in internal damage of concrete were detected through the coda wave interferometry method, mainly to determine the small difference between the coda wave waveforms before and after the concrete is loaded. A comparison of the received signal waveforms before and at the beginning of loading is shown in Figure 5. In time window 1, the two received signals are completely coincident, which ensures the repeatability of the wave source. In time windows 2 and 3, there is a phase difference between the two signals that gradually increases. Comparing windows 1, 2, and 3 proves that the coda wave is very sensitive to small changes in the random medium of concrete.

The value of the cross-correlation coefficient $C{C}_i(\epsilon )$ ranges from 0 to 1. A value of 0 means that the two waveforms are completely uncorrelated; 1 means that the two waveforms are completely the same [23]. Figure 6 shows the correlation between the coda wave cross-correlation coefficient and the maximum load percentage of the concretes. Within 50% of the maximum load percentage (A–B–C), the correlation coefficient of each concrete is approximately equal to 1, so the corresponding wave velocity change is linear and uniform [36]. When the maximum load percentage is between 50–80% (C–D), the cross-correlation coefficient decreases continuously, indicating that coda wave propagation is strongly affected by the internal damage to the concrete, which causes the waveform to be distorted. As the load exceeds 80% of the maximum load percentage (D–E), the cross-correlation coefficient exhibits a larger range of fluctuations based on further reduction, indicating that the typical causing of scatter inside the concrete, that is, the cracks, has further widened or linked up.

3.2. Coda Wave Velocity Change Rate Analysis

The curves of the maximum load percentage–wave velocity change rate of the four strengths of concrete samples are shown in Figure 7. According to the trend of the changes, the curve can be divided into four parts: strengthening (A–B), stabilizing (B–C), weakening (C–D), instability (D–E).

The strengthening section (A–B) of the curve ranges from 0% to 30% of the maximum load percentage. The coda wave velocity of the concrete sample shows a slowly increasing trend with increasing the load percentage, and the increased speed is negatively correlated with the load increase. This is because, in the early stage of loading, the internal cracks slowly close, making the internal microstructure denser. The stabilizing section (B–C) of the curve ranges from 30% to 50% of the maximum load percentage. The rate of change of the coda wave velocity of the C40 and C50 concrete specimens remained stable. The coda wave velocity changes of C60 and C70 concrete specimens first increase and then tend to be generally stable. It is speculated that within this range, the original microcracks in the concrete are further closed, and at the same time, the material particles are squeezed and elastically deformed. However, with the increase of pressure, the degree of development of internal microcracks may increase slightly. The range of the weakening section (C–D) is approximately 50% to 80% of the maximum load percentage. The wave velocity of the coda wave decreases with the increase of external load. At this time, the closed cracks inside the concrete begin to expand and secondary cracks are generated, causing the wave velocity to drop. The instability section (D–E) ranges from 80% to 100% of the maximum load percentage. In this interval, the coda wave velocities of the four strengths of concrete showed a rapid downward trend with the increase of the load, indicating that the scale of the internal cracks in the concrete continued to increase and gradually connected. Eventually, a macroscopic fracture surface will be formed.

This experiment comprised a coda wave test on a limited number of concrete samples. Due to the complexity and randomness of concrete materials, the relationship between the rate of change of coda wave velocity of uniaxially compressed concrete and its strength grade has not yet been clarified. To complete the research for this aspect, many tests of concretes of various strengths are still needed.

3.3. Coda Wave Velocity Change Rate Fitting

Piecewise linear fitting was performed on the coda wave velocity change curve of the abovementioned concrete uniaxial compression process. The fitting function expression is $y=ax+b$, where $x=F/F_{\max }$. $F$ is the concrete load value, $F_{\max }$ is the concrete ultimate load value, and $F/F_{\max }$ is the ratio of the concrete load value to its ultimate load value; that is, the maximum load percentage. The notations a and b are constants. The curve fitting results of four kinds of concrete coda wave velocity change rates are shown in

Table 3. Fitting curve parameters of concrete with different strengths.

Concrete Strength	x = F/Fmax	a	b	R2
C40	0–30%	0.0526	0.00153	0.94931
	30–50%	−0.00707	0.01709	1
	50–80%	−0.21354	0.12101	0.99099
	80–100%	−1.06024	0.78913	0.99304
C50	0–30%	0.05787	0.00301	0.9146
	30–50%	−0.00801	0.01954	1
	50–80%	−0.20256	0.11869	0.97374
	80–100%	−1.03074	0.78799	0.99365
C60	0–30%	0.08213	0.00686	0.84953
	30–50%	0.02818	0.01879	1
	50–80%	−0.15884	0.11397	0.98072
	80–100%	−1.23438	0.98672	0.98572
C70	0–30%	0.08939	0.00305	0.96862
	30–50%	−0.00745	0.03426	−0.10323
	50–80%	−0.16207	0.11139	0.95158
	80–100%	−1.93605	1.53482	0.99379

. The fitting effect of C70 concrete is poor when the load level is in the range of 30–50% of the maximum load (R² = −0.10323). The data show a trend of increasing first and then decreasing. The poor fitting effect is mainly due to the large difference between the maximum value and the minimum value of the data in this group (the difference is equal to 0.00469). However, as shown in Figure 7, the general trend in this section is the same as that of C60 concrete. It is speculated that the dispersion of the existing data is caused by the complexity of the test material and the deviations in the operating process.

Combined with Table 3, a mathematical formula for concrete damage under uniaxial compression based on the rate of change in coda wave velocity is constructed:

$$\frac{dv}{v}=\left\{\begin{array}{ll}{a}_1\frac{F}{F_{\max }}+b_1& 0<\frac{F}{F_{\max }}\le 30\%\hfill \\ a_2\frac{F}{F_{\max }}+b_2& 30\%<\frac{F}{F_{\max }}\le 50\%\hfill \\ a_3\frac{F}{F_{\max }}+b_3& 50\%<\frac{F}{F_{\max }}\le 80\%\hfill \\ a_4\frac{F}{F_{\max }}+b_4& 80\%<\frac{F}{F_{\max }}\le 100\%\hfill \end{array}\right.$$

(5)

In the above formula, $a_i$ is the slope of the fitted straight line for each interval. Considering the discreteness of individual data, the overall law of the experiment is $0<a_1<0.1$, $a_2\approx 0$, $-1\le a_3<0$, and $a_4<-1$. $b_i$ is a constant and depends on the strength of the concrete.

4. Numerical Simulation

4.1. Numerical Simulation of Damaged Concrete under Uniaxial Compression

4.1.1. Concrete Model Establishment

In the particle discrete element program, clump theory is used to overlap most particles to form an irregularly shaped aggregate. It is a rigid body that will not break under the action of external force and always maintains a complete shape [37,38]. Therefore, it is feasible to use irregular clump bodies to simulate the real coarse aggregates in concrete when establishing the PFC2D model (PFC software: Particle Flow Code, Copyright Itasca Consulting Group, Inc. PFC Version 5.0, last updated 24 June 2016.)

First, a picture of the coarse basalt aggregate was taken (Figure 8a), closed multiline segments were used to outline the boundary of the aggregate, and 21 types of coarse bone particle patterns were finally selected. The abovementioned particle geometry model was imported into the software to generate a clump template, which represents the real coarse aggregate. The contact model between particles was defined as a parallel bonding model in which there are forces and moments between the particles. This feature is consistent with the material stress on the inside of the structure when the actual concrete is under pressure. The boundary wall was generated according to the size of the concrete sample, and then the “distribute” command was used in the wall area to randomly generate the set size particle unit to represent the mortar. Subsequently, clump bodies are generated randomly in the corresponding particle size gradation range, as shown in Figure 8b. In this paper, the “solve aratio” command in PFC2D is used to make the initial sample reach the equilibrium state. The aratio (or average ratio) is the ratio of the average value of the unbalanced force magnitude (i.e., the magnitude of the sum of the contact forces, body forces, and applied forces) on overall bodies to the average value of the sum of the magnitudes of the contact forces, body forces and applied forces on overall bodies. In this simulation, the value of aratio is set to 1e-4. In the enclosed space surrounded by the surrounding walls, the particles are calculated to an equilibrium state, and the unbalanced force gradually decreases, making aratio less than 1e-4, which has no impact on the subsequent calculation. The specimen is not consolidated; in other words, the initial pressure is 0. The voids in the concrete c in the equilibrium state are all distributed, and there is a certain contact between the particles, as shown in Figure 8c.

4.1.2. Calibration of the Meso-Parameters of the Concrete Model

In this paper, five meso-parameters, namely the elastic modulus $E_c$, stiffness ratio $k_n/k_s$, friction coefficient $\mu$, parallel bond normal strength ${\overline{\sigma }}_c$, and parallel bond tangential strength ${\overline{\tau }}_c$, are calibrated. Using the control single variable method, the changes in the peak strength and damage failure mode of concrete caused by the change of meso-parameters of each particle are discussed. After the values of each meso-parameter are tested, the parameter set that best matches the stress–strain curve obtained from the experiment is selected. This chapter presents the established numerical model of C40 concrete specimen, and the discussion values and final simulation values of each particle meso-parameter are shown in

Table 4. Mesoscopic parameter values of C40 concrete.

Meso-Parameters	${\mathit{E}}_{\mathit{c}}$/GPa	${\mathit{k}}_{\mathit{n}}/{\mathit{k}}_{\mathit{s}}$	$\mathbf{\mu }$	${\overline{\mathbf{\sigma }}}_c$/MPa	${\overline{\mathbf{\tau }}}_{\mathit{c}}$/MPa
Discussed Value	27.9	0.38	0	8.5	33
	32.9	1.38	0.3	13.5	38
	37.9	2.38	0.6	18.5	43
	42.9	3.38	0.9	23.5	48
	47.9	4.38	1.2	28.5	53
	52.9	5.38	1.5	32.5	58
The best value for simulation	42.9	3.38	0.9	23.5	48

.

The uniaxial compression simulation results of the C40 concrete sample were compared with the test results. The main mechanical parameters are shown in

Table 5. Comparison of the simulation results and test results of C40 concrete uniaxial compression.

Comparison	Elastic Modulus/GPa	Peak Stress/MPa	Peak Strain
Simulation	33.43	43.85	0.002603
Experiment	33.32	43.35	0.002593
Deviation percentage	0.33%	1.12%	0.38%

. Figure 9 shows that before the concrete reaches the peak stress, the experimental and simulated stress–strain curves are almost the same, but there is a clear difference in the stress softening stage. Compared with the experimental value, the simulated curve yields and decreases earlier. The main reason is that the two-dimensional model used in the simulation cannot truly reflect the three-dimensional compression situation. The experimental and simulated values are very close, and the deviation is small. Figure 10 shows that the concrete failure modes obtained by the uniaxial compression experiment and PFC2D simulation are the same; both are “C”-shaped crack fracture zones.

4.1.3. Numerical Simulation Results of Concrete Uniaxial Compression

Using the abovementioned concrete uniaxial compression simulation model, the development process of concrete internal cracks is obtained, as shown in Figure 11. At the initial stage of loading, the particle stress level in the model is low, and a small number of particles gradually generate microcracks. Early damage occurs between the coarse aggregate and small particle mortar on the upper bearing surface and the middle of the model. Subsequently, stress redistribution occurs in the concrete sample, and the stress concentration phenomenon is easily formed in the place where the microcrack occurs, which intensifies its continuous expansion. As the load further increases, the original cracks at the interface between the aggregate and the mortar, as well as inside the mortar, will extend, while many new cracks will appear at the interface of other aggregates. When peak stress is reached, displacement of the block under the action of external force is increased so that the local cracks become fully expanded and linked up, forming evident macroscopic cracks.

4.2. Ultrasonic Field Simulation Based on Damaged Concrete Images

4.2.1. Ultrasonic Field Model of Damaged Concrete

• Ultrasonic field setting

Based on the existing elastic wave equation, the first-order velocity–stress equation of the two-dimensional acoustic wave equation is derived, namely:

$$\left\{\begin{array}{l}\frac{\partial P}{\partial t}=\lambda \left(\frac{\partial v_x}{\partial x}+\frac{\partial v_z}{\partial z}\right)\\ \rho \frac{\partial v_x}{\partial t}=\frac{\partial P}{\partial x}\\ \rho \frac{\partial v_z}{\partial t}=\frac{\partial P}{\partial z}\end{array}\right.$$

(6)

where $P$ is the sound pressure, $\rho$ is the density of the particle, t is the time, $\lambda$ is the Lame constant, and $v_x$ and $v_z$ are the velocity components of the particle vibration in the x and z directions, respectively.

The finite-difference time-domain method is used to solve the above equation to simulate the propagation of ultrasonic waves in concrete. The receiving boundary of the concrete model is set as the total reflection boundary, and the other boundaries are set as the perfect matching layer (PML) absorbing boundary [39,40]. The source excitation waveform is the Ricker wavelet, and the seismic source is loaded on the center of the upper surface of the concrete model as a point source.

• Digital image processing

Digital image processing [41,42] was performed on concrete images with different damage degrees obtained in the numerical simulation of concrete damage under uniaxial compression in Section 4.1. In a grayscale image, the brightness of different colors represents a grayscale value, which can be divided into values from 0 to 255. The higher the grayscale value, the lighter the color. Figure 12a is the damage image of 100% concrete ultimate strength calculated by uniaxial compression simulation and Figure 12b is the corresponding grayscale damage image.

• Threshold segmentation and assignment

The image is segmented according to the difference in gray value, where gray 0 represents aggregate, gray 255 represents fracture, and other gray values represent mortar. Image units of the same brightness are assigned the same number, and the wave velocity and density are assigned to the corresponding numbered image areas. The actual physical parameters of the different grid materials are assigned in the model; that is, the wave speed of mortar is set to 4000 m/s and the density is 2073 kg/m³; the wave speed of aggregate is set to 5500 m/s and the density is 2760 kg/m³; the wave speed of air is 340 m/s and the density is 1200 kg/m³. After multi-threshold segmentation is performed on the grayscale damage grayscale images of 0% and 100% concrete ultimate strength, the velocity assignment statistics can be obtained, as shown in Figure 13. The concrete numerical model after wave velocity assignment is shown in Figure 14.

4.2.2. Analysis of the Ultrasonic Field of Damaged Concrete

Twenty-four concrete images with different damage degrees were selected to carry out the numerical simulation of the ultrasonic field of the damaged concrete. Since the lateral boundaries of the 2D concrete model are simplified as PML absorption boundaries, the simulated coda signal is very weak. Figure 15 shows the comparison diagram of the ultrasonic signal of the concrete non-damage and the initial damage models; the corresponding load value is about 10% of the maximum load. Comparing the waveforms in time windows 1, 2, and 3, the two waveform signals in window 3 have an obvious phase difference, indicating that the coda wave is very sensitive to the slight changes in the concrete random medium. The waveform signal comparison chart obtained by simulation and experiment showed similar characteristics, which proved the credibility of ultrasonic field simulation.

The coda wave signal of 24 concrete models with different damage degrees was processed by the stretching technique, and the corresponding coda velocity change rate was calculated, as shown in Figure 16. The change rate of the coda wave velocity measured in the experiment is roughly the same as the change rate of the coda wave velocity obtained by the ultrasonic field simulation. However, the numerical simulation result of the ultrasonic field is slightly larger than the actual test result. It is speculated that the reason may be that the concrete simulation model is simplified to two dimensions, and the boundary condition settings are also different from the real boundary.

The load percentage is in the range of 30–100%, and the wave velocity change trend of simulation and experiment is roughly the same. However, there is a large difference between the two in the load percentage of 0–30%. The main reason is that the uniaxial compression damage simulation of concrete based on the PFC2D discrete element method cannot simulate the process of concrete compaction during the initial compression stage in the test.

Table 6. Comparison of the simulation results and test results of C40 concrete uniaxial compression.

Slope	0–30%	30–50%	50–80%	80–100%
Simulation	−0.03845	−0.06839	−0.27649	−1.48015
Experiment	0.052	−0.00707	−0.21354	−1.06024
Fitting formula	$0<a_1<0.1$	$a_2\approx 0$	$-1\le a_3<0$	$a_4<-1$

compares the fitting slope of the wave velocity curve calculated by the simulation with the test results. In the 0–30% interval, as mentioned earlier, the difference between the two is mainly caused by the deficiencies of the simulation calculation software. In the load range of 30–50%, the slope of the two still has a large deviation. However, the curve corresponding to the range of 30–40% tends to be horizontal, and the calculated slope is −0.00821, which conforms to the law of the fitting formula. In the 50–80% and 80–100% sections, the simulation, experiment, and fitting formulas all show a high degree of agreement. In general, the numerical simulation results are consistent with the experimental results, which proves the feasibility of monitoring concrete damage using the fitting formula of the coda wave velocity change rate and the maximum load percentage of concrete under uniaxial compression.

5. Discussion

In this paper, four types of concrete test blocks, namely C40, C50, C60, and C70, were tested under uniaxial compression to study the coda wave velocity changes during the damage process in concrete of different strengths. The stretching technique was used to analyze the ultrasonic coda signal of concrete with different strengths of concrete under uniaxial stepwise compression. Both the correlation coefficient and the rate of change of coda wave velocity can reflect the internal damage degree of concrete under load. According to the slope of the rate of change of coda wave velocity, the relationship curve between the rate of change of coda wave velocity and the percentage of maximum load during uniaxial compression can be divided into four sections: strengthening, stabilizing, weakening, and instability. The turning points of adjacent sections are 30%, 50%, and 80% of the maximum load percentage. A concrete damage fitting formula based on the relationship between the rate of change of the coda wave velocity and the percentage of the maximum load under uniaxial compression is established. The slope of the coda wave velocity change rate is greater than or close to 0, the internal microcracks in the concrete are closed, and the particles are squeezed and deformed. The slope of the coda wave velocity rate change between −1 and 0 indicates that the initial concrete cracks are propagating and new cracks are initiating and developing. When the slope of the coda wave velocity change rate is less than −1, the internal cracks of the concrete will be connected and link up, then the block will slide under the action of external force, causing macroscopic rupture. It is worth noting that within the scope of this study, no clear relationship between the ultrasonic coda wave signal characteristics of uniaxially compressed concrete and its intensity level was found. Due to the discrete nature of concrete materials, a lot of experimental research is still needed on this problem.

The simulation method is used to capture the whole process of uniaxial compression concrete damage and destruction. Damage cracks should first appear at the weak interface between aggregate and mortar. As the load increases, the original cracks and the new cracks that are generated continue to expand and connect. In the end, the concrete was severely damaged and its bearing capacity was lost, and the damage was in the form of double “C”-shaped macroscopic cracks. This method provides images of different degrees of damage for concrete ultrasonic field simulation. The digital image processing technology was used to extract aggregate, cement mortar, and crack information from concrete images of different degrees of damage. The concrete damage model established by processing grayscale images by threshold segmentation technology simulated the propagation of ultrasonic coda waves in concrete with different damage levels and quantitatively calculated the coda wave velocity change rate. The results show that the concrete wave velocity change rate obtained by the numerical simulation is consistent with the overall change trend obtained by the experiment, and the expected effect is achieved. Because the numerical modeling requires a lot of trial calculations and adjustment, within this research period, only the uniaxial compression simulation and the ultrasonic field simulation of the C40 concrete was completed. In the later stage, it is still necessary to optimize the numerical simulation, so as to realize the ultrasonic field numerical simulation of a large number of concrete damage processes using the above method and provide a database for related research.

Based on the perspective of correlating the characteristics of ultrasonic coda waves with the damage and failure process of concrete structures, only the compressive failure of pure concrete was studied in this work. However, in practical engineering, structural components are usually made of reinforced concrete, and the failure modes of components under different mechanical effects are complex. This study simplifies the research problem but can provide a basis for further research on the ultrasonic coda wave characteristics of reinforced concrete members under complex mechanics.

6. Conclusions

In this paper, a combination of tests and simulations is used to explore the ultrasonic coda wave characteristics of uniaxially compressed concrete. The main conclusions are as follows:

The ultrasonic coda signal is processed by the stretching technique. The cross-correlation coefficient and the relative change in coda wave velocity are very sensitive to small changes in the concrete. The regularity of the relationship between the change rate of the coda wave velocity and the percentage of the ultimate load is more significant than the cross-correlation coefficient.

According to the slope of the curve, the relationship curve between the concrete coda wave velocity change rate and the maximum load percentage is divided into four sections, namely strengthening (the slope is greater than 0), stabilizing (the slope is approximately equal to 0), weakening (the slope is between −1 and 0), and instability (the slope is less than −1). The internal damage changes of the concrete corresponding to the curves of the four sections are pore compaction, grain deformation, crack propagation, and block slip in sequence.

Using the PFC uniaxial compression simulation, the image of the internal damage development process of the concrete can be obtained, and the ultrasonic field simulation of the concrete damage and rupture process can be completed by further combining the image processing technology. The calculation results obtained by the above methods are in line with expectations. However, in the process of PFC modeling, it is necessary to perform trial calculation and calibration of multi-parameters to match the experimental results. Due to the uniqueness of concrete compression test results, the obtained damage rupture model (C40) is not suitable for simulating the compressive failure mode of concrete with various strengths. The modeling method of concrete compressive failure still needs to be further optimized.

The numerical simulation method of concrete ultrasonic coda wave characteristics established in this paper is conducive to the simulation calculations for many concretes and provides a feasible method for further obtaining a more accurate quantitative relationship between concrete damage degree and coda wave characteristics.