**3.3 Synthetic textile and tyres**

Browne et al. [25] assorted that process of laundries in household and industries leads to mass production of primary microplastics via scraping and dispersion of fibers, which was then emitted out in sewage water and culminates in the ocean [31]. Tyres contain profuse mixture of several synthetic polymers in addition with natural rubber. Tyres get deteriorated when used and tyre dust that contains synthetic rubber circulated by wind or swiped away by rain. Large segments of such particles were reported to congregate in the sea [19, 23].


[163] Evangelista J.P.C.; Chellappa T.; Coriolano A.C.F.; Fernandes Jr. V.J.; Souza L.D.; Araujo A.S.; *Synthesis of Alumina Impregnated with Potassium Iodide Catalyst for Biodie‐*

[164] Arzamendi G.; Campo I.; Arguiñarena E.; Sánchez M.; Montes M.; Gandía L.M. *Syn‐ thesis of Biodiesel with Heterogeneous NaOH/Alumina Catalysts: Comparison with Homoge‐*

[165] Baroutian S.; Aroua M.K.; Raman A.A.; Sulaiman N.M.N. *Methanol Recovery During Transesterification of Palm Oil in a TiO2/Al2O3 Membrane Reactor: Experimental Study and*

[166] Wu H.; Zhang J.; Wei Q.; Zheng J.; Zhang J. *Transesterification of Soybean Oil to Biodie‐ sel Using Zeolite Supported CaO as Strong Base Catalysts.* Fuel Process. Technol. 2013,

[167] Babajide O.; Musyoka N.; Petrik L.; Ameer F. *Novel Zeolite Na-X Synthesized from Fly Ash as a Heterogeneous Catalyst in Biodiesel Production.* Catal. Today 2012, 190, 54-60.

[168] Alves C.T.; Oliveira A.; Carneiro S.A.V; Silva A.G.; Andrade H.M.C.; Vieira de Melo S.A.B.; Torres E.A. *Transesterification of Waste Frying Oil Using a Zinc Aluminate Cata‐*

[169] Borges M.E.; Díaz L. *Recent Developments on Heterogeneous Catalysts for Biodiesel Pro‐ duction by Oil Esterification and Transesterification Reactions: A Review.* Renew. Sust. En‐

[170] Di Serio M.; Tesser R.; Pengmei L.; Santacesaria E. *Heterogeneous Catalysts for Biodiesel*

[171] Wen Z., Yu X., Tu S.T., Yan J., Dahlquist E. *Biodiesel Production from Waste Cooking Oil Catalyzed by TiO2–MgO Mixed Oxides.* Bioresource Technol. 2010, 101, 9570–9576. [172] Guo F.; Peng Z.G.; Dai J.Y.; Xiu Z.L. *Calcined Sodium Silicate as Solid Base Catalyst for*

[173] Singh A.K.; Fernando S.D.; Hernandez R. *Base Catalyzed Fast Transesterification of Soy‐*

[174] Kalva A.; Sivasankar T.; Moholkar V.S. *Physical Mechanism of Ultrasound-assisted Syn‐*

[175] Nishimura CSMV; Maeda Y.R. *Conversion of Vegetable Oil to Biodiesel Using Ultrasonic*

[176] Salamatinia B., Mootabadi H., Bhatia S., Abdullah A.Z. *Optimization of Ultrasonic-As‐ sisted Heterogeneous Biodiesel Production from Palm Oil: A Response Surface Methodology*

*Biodiesel Production*. Fuel Process. Technol. 2010, 991, 322–328.

*bean Oil Using Ultrasonication.* Energ. Fuel 2007, 21, 1161–1164.

*thesis of Biodiesel*. Ind. Eng. Chem. Res. 2008, 48, 534–544.

*Approach.* Fuel Process. Technol. 2010, 91, 441–448.

*sel Production from Rice Bran Oil.* Fuel Process. Technol. 2012, 104, 90–95.

*Neutral Network Modeling*. Sep. Purif. Technol. 2010, 76, 58–63.

*neous NaOH.* Chem. Eng. J. 2007, 134, 123–130.

*lyst*. Fuel Process. Technol. 2013, 106, 102–107.

*Production.* Energ. Fuels 2008, 22, 207–17.

*Irradiation.* Chem. Lett. 2003, 32, 716–717.

erg. Rev. 2012, 6, 2839–2849.

109, 13-18.

440 CO2 Sequestration and Valorization


**4. Conclusions**

This paper describes a WRMS-based ultrasonic non-destructive testing approach based on the scanning laser Doppler vibrometry technique. The damage identification and imaging in composite steel–concrete beams were performed. Based on the results of the experimental investigations evaluated by numerical FEM simulations and supported by theoretical analysis, the following final remarks can be formulated.

The leading conclusion of the presented research is that the guided wave propagation method supported with WRMS calculation can be successfully applied for damaged imaging in steel–concrete composite beams. However, considering the complexity of guided wave propagation phenomenon, it is important to perform some initial analysis, including dispersion curves determination and FEM numerical calculations. Further works will be directed to develop the proposed method to be applicable for non-destructive diagnostics of real-scale engineering structures.

**Author Contributions:** Conceptualization and Methodology, E.W. and M.R.; Experimental Investigations, M.R. and E.W.; FEM Calculations, E.W. and M.K.; Formal Analysis, E.W.; Visualization, E.W.; Writing—Original Draft Preparation, E.W. and M.K.; Writing—Review and Editing, MR.; Supervision, M.R.; All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Acknowledgments:** Abaqus calculations were carried out at the Academic Computer Centre in Gda ´nsk.

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


**2. Case Study**

A problematic reach of the upper Hungarian Danube reach (river km (rkm) 1796–rkm 1794, Figure 1) has undergone major morphological changes during the last decades. Many studies presented in the literature [13–19] that because of installations of many river regulation measures (e.g., groin fields, ripraps and hydropower plant at rkm 1819) in the last decades, intensive gravel formations [19], important bed level incision [20] and bed armoring processes [16,19] could be detectable, mainly in the main channel [21]. In contrast, the bed content is much finer in the groin fields, causing siltation and erosion of the finer sediments during flood waves [19,22].

**Figure 1.** The sketch of the investigated Danube study reach (**left**) and grain size distributions taken from the investigated reach (**right**). The characteristic water discharges are mean flow *Qm* = 2000 m3/s, range of bankfull discharge *Qbf* = 4300–4500 m3/s, 2, 10 and 100-year flood event *Q2* = 5950 m3/s, *Q*10 = 7950 m3/s and *Q*100 = 10,400 m3/s [23].

The river reach can be characterized by the following parameters: the channel bed width at mean water-stage ranges between 150 m and 350 m [19] with the average water surface gradient of 0.0002–0.00025. As Figure 1 shows, the river section is regulated by conventional structures, such as groins and the banks are protected against erosion by ripraps. Also, sidearms, islands, gravel bars, confluence zone can be observed, which refer to the complex topography of the river section. Bed material samples were taken from the main stream, groin fields and gravel bars. Some of the grain size distributions can be seen in Figure 1 right. The figures refer to very diverse spatial bed contents (0.32 mm < *d*50 < 70.5 mm, where *d*50 is the median grain size). Such a wide dispersion of the bed content is a unique feature of the Danube River (~rkm 1600–rkm 1800); at the Lower Austrian Danube, (~rkm 1885, 90 km upstream), the Danube flows through a gravel bed, where *d*50 is 21.1 mm without any finer fractions [24]. In turn, the middle Hungarian Danube (~200 km downstream) has a typical sandy bed with *d*50 < 0.05 mm [25]. The complexity of the topography and bed content suggest spatially and temporally varied sediment transport nature [22]. That is in some places the gravel, elsewhere the sand transport dominates [26]. This kind of individual complexity was presented e.g., by the field measurements of Török and Baranya [19], or in [27,28].

Based on field measurements, the researchers could make approximate estimates regarding the ongoing local- and reach-scale morphodynamic processes. For instance, the essential bed changes in the last decades caused important water management related problems and also difficulties in navigation, which could be detected with the field measurements [19,20,27–29]. However, the studies emphasize that additional studies are required for even more reliable morphodynamic investigation. Computational modeling offers an alternative for analyzing and predicting the recent and expected processes [30,31]. That is, the numerical modeling based investigation could be an important complementary manner with meaningful added value.

For this reason, the reliable calculation of the morphological changes is a major interest to researchers and the application of a 3D computational fluid dynamics (CFD) sediment transport model became justified. However, the choice of the applied sediment transport model was not obvious. Many formulas can be found in the literature (e.g., [1–4,32]). Most of them are developed focusing on a given morphodynamic process (e.g., Wilcock and Crowe model: bed armoring [8]; van Rijn model: sand erosion and deposition [33] etc.). However, in the case of the examined river section, spatially and temporally varied sediment transport nature occurs. That is, none of the existing sediment transport formula are expected to operate reliably for both the sand and coarse bed material, within a given river reach.

#### **3. Materials and Methods**

#### *3.1. Introducing the Combining Sediment Transport Calculation Method*

A novel combined approach [12] of the van Rijn and the Wilcock and Crowe bed load sediment transport formulas were applied because of the complex and spatially varied bed material and dominant sediment transport nature. From now, the Wilcock and Crowe formula will be indicated with *W&C*, while the van Rijn will be with *vR*. A brief summary of these models can be read in the Appendix A. The combined manner was already presented and validated against laboratory measurements [12]. In that study, laboratory experiments were used, which were carried out in an 11 m long, 1 m wide straight flume with non-uniform bed material. A groin was installed 5 m far from the outlet, which caused scouring, bed armoring and sand aggradation processes in the channel [34].

Török and Baranya [22,26] pointed out a novel decision criterion, which is a suitable method for indicating whether the sand or rather the sand transport dominates locally. The herein presented combining method was based on this statement. Namely, if the shear Reynolds number (*Re\**) is below 300, the sand transport is prevalent. Otherwise, if the *Re\** occurs above 400 the gravel transport dominates. Based on these, the combined calculation method said that the local bed load sediment transport rate was calculated as the following (*qbi*,*W*&*C* is the sediment transport rate calculated by *W&C* and *qbi*,*vR* is the rate by *vR*):

$$q\_{hi} = \begin{cases} \begin{array}{c} q\_{hi, \mathcal{W} \& \mathcal{C}} \text{ if } \mathcal{R} \mathbf{c}^\* > 400 \\ q\_{hi, \mathcal{rR}} \text{ if } \mathcal{R} \mathbf{c}^\* \le 300 \\ \end{array} \\ \begin{array}{c} \text{else} \end{array} \end{cases} \tag{1}$$

where

$$f = \frac{1}{100}(R\varepsilon^\* - 300). \tag{2}$$

These models do not contain any variable to calibrate, the constants (e.g., the constant von Karman, or the lift coefficient) and base equations (e.g., the original Shields curve [35] was used in the van Rijn model) were defined as they are published and recommended in the original papers [8,33].

In addition to the bed load transport estimation, the suspended sediment transport was calculated in each computational grid according to the suspended *vR* formula [36]. Thus, the bed load was calculated according to Equations (1) and (2), while the suspended load was estimated by the van Rijn equation.

### *3.2. Applied 3D Flow Model*

The numerical model used in this study [37,38] solved the 3D Reynolds-averaged Navier–Stokes (RANS) equations with the *k-*ε turbulence closure (see e.g., [39]) by using a finite-volume method and the (semi-implicit method for pressure linked equations) SIMPLE algorithm [40,41] on a 3D non-orthogonal grid. At the boundaries, where the fluid flow cannot be considered as a free turbulence zone, the wall law was applied for the velocity profile calculation [42]. The momentum equations were in the complete form, describing the hydrodynamic effects in all directions. The roughening impact of the vegetation in the flood plain area was described as an energy loss term in the Navier–Stokes equations [43], and could be specified for each cell. Using this option, the effect of the vegetation was taken into account as a drag-effect.

In order to eliminate the boundary effect, the computational grid was longer in both upstream (rkm 1801) and downstream (rkm 1793.5) direction than the investigated ~ rkm 1795 and rkm 1799 river reach. The applied grid can be seen in Figure 2.

**Figure 2.** The computational grid of the investigated Danube study reach.

The study side was discretized with 355 cells in the streamwise direction and 150 cells in the lateral direction, respectively, resulting in the streamwise direction in an average resolution of 18 m and transversely 5 m in the main channel, while 13 m in the floodplain area. Vertically 11 layers were defined, at every ten percent of the depth, resulting in a maximum of 585,750 computational cells. The resolution of the grid was set based on previous 3D model studies [37,44–48] and Section 3.2.2 demonstrates its applicability.

The bed material was discretized by five fractions, which were: *d*1 = 0.3125 mm, *d*2 = 1.25 mm, *d*3 = 5.7 mm, *d*4 = 16.2 mm and *d*5 = 56.57 mm. The ripraps and groins were characterized by *d* = 300 mm. According to field measurement based considerations, the active layer thickness was set to 0.5 m.

#### 3.2.1. Parameterization

For instance, Parker mentioned [11] that the bed material of most river reaches is less complex and the grain sizes happen in a narrow range. In turn, in rare cases—e.g., the herein studied river section—the occurring grain sizes cover a significantly wider range (silt–gravel), resulting in a very complex spatial distribution. The bed material cannot be supposed as spatially uniform because of this (as many studies do at most river sections [44,46]), which makes the allocation of the bed material a less obvious method. According to Baranya [49], a relation can be stated between the calculated local bed shear stress value by 3D flow model at the mean-water stage and the local *d50*. Thus, applying the fitted function, a transitional and continuous *d50* map can be estimated based on the calculated bed shear stress distributions. This method was used based on 33 bed material samples [19,26]. The standard deviation of the calculated *d50* to the measured *d50* was 3.2 mm.

As the boundary conditions for the RANS equations, at the inflow boundaries the water discharge, at the outflow boundary, the water level was set. The discharge time series (Figure 3) and the water levels of the Danube were defined based on the measured time series. Additionally, the flow discharge time series of the Mosoni-Danube was set based on a 1D numerical Danube model [50].

**Figure 3.** Discharge time series at rkm 1801 for the period October 2012–October 2014. The different colored periods marked with letters a–h display the flood waves (including the 100-year historical flood wave from 2013, e) which exceed the bed-forming flow discharge (*Q* > 2100 m3/s).

An essential part of the model setup is the correct set of the inflow sediment rate. For this purpose, the flow discharge dependence of the suspended load [51] and the bed load [19] functions were used.

The riverbed topography of the main river channel was available from 2012 and 2014. The initial bed geometry was set according to the map from 2012. The calculated bed change map could be prepared for this two-year-long period, which included the historical flood wave from 2013 (Figure 3). This bed change map was used as a benchmark for the validation purpose (Figure 4).

**Figure 4.** Measured bed changes for the period October 2012–October 2014.

Four regions in the river reach were highlighted by green rectangles and an ellipse, marked by A, B, C and D. In these places, the following bedforms and morphodynamic processes were detected by field measurements [16,19]. At region A and D, the blue spots refer to a pronounced scouring downstream of the groins. In region B a groin field could be found. Here, local bed changes took place, both scouring processes (blue spots) and sediment depositions (brown spots). The ellipse (C) and the brown spot in region D show the places where gravel bars were located. As these phenomena basically determine the reach-scale morphodynamic processes, a key question is whether the novel sediment transport calculation manner introduced a more reliable estimation of them.

The numerical simulation of the 2 years, 722 days long period demands very large computational capacity. According to the preliminary estimation calculations, the simulation of one model variant for such a long time period would take around half a year. Therefore, to reduce the duration of the simulation, only the periods that exceeded the bed-forming flow discharge (*Q* > 2100 m3/s) [19,24] were simulated. It meant that the inflow boundary condition was set to a continuous flow discharge series, including the highlighted periods in Figure 3a–h, respectively and omitting the intermediate periods. In this range, 66% of the annual bed load amount passed [23]. In turn, it is emphasized that the ignored 34% of the annual bed load yield was still significant. That is, the numerical model neglects the simulation of bed changes that take place during the lower water regime. However, many studies (including papers regarding the investigated river reach) display that the major bed changes such as scouring, bar formations and flushing of the groin fields are expected rather during floods [16,19,24,52]. That is, the mean and lower water stages are less important in the view of bed changes. Based on the above, the bed changes caused by the eight flood waves (Figure 3a–h, a total of 211 days) approach well the real two years changes. Hence, the calculated bed changes were comparable to the measured.

#### 3.2.2. 3D Flow Model Validation

The herein applied 3D CFD model was already adapted and validated for the investigated Hungarian reach of River Danube, which was published in previous research works, e.g., [45,46,53,54]. Those studies have already demonstrated the reliable application of the 3D flow model. Regardless of these, the flow model validation was elaborated for the peak of the historical flood wave in 2013. The cross-sectional acoustic Doppler current profiler (ADCP) flow measurements by the North-Transdanubian Water Directorate regarding the peak stage of the flood wave were used as benchmark flow values. Figure 5 shows the measured (left) and the calculated (right) cross-sectional velocity distributions, from exactly the same points. However, as the ADCP instrument is not able to measure at the direct water - and bed surfaces, the measured cross-sections are smaller slices, respectively. Furthermore, regarding the comparison it is important to mention that the moving boat ADCP measurements record the momentary flow conditions, which include the effect of turbulence, resulting in velocity pulsation in 0.1 m/s order of magnitude [25,53,55]. In turn, the RANS model calculates the time-averaged velocity values.

Comparing the measured and calculated cross-sections, the following remarks can be stated. The velocity values were in the same ranges (0–3 m/s). The highest velocities (yellow and red spots) were calculated at the same place as the cross-sections than in the real case, which underlined the reliable estimation of the main stream. The locations of the lower velocities (blue and light green spots) were calculated as trustworthy also. That is the calculated flow pattern could be realized reliable at the groin fields (e.g., at the right-bank sides of cross-section III, IV and V, Figure 5) and gravel bars also (e.g., at the right-bank sides of cross-section VII). Finally, the root-mean-square deviation (RMSD) between the measured and calculated velocities could be calculated for all the seven cross-sections. These values can be seen in Table 1.

**Figure 5.** Measured and calculated cross-sectional horizontal velocity distributions.

**Table 1.** The root-mean-square deviation (RMSD) between the measured and calculated velocities, the maximum cross-sectional horizontal velocity and the percentage deviation of them for the seven cross-sections.

Considering the maximum cross-sectional horizontal velocities (Table 1, third row), the RMSD values displayed that the deviation of the flow calculation was below ~12% of the maximum values. The only two sections, where the inaccuracy happens ~15% were cross-sections V and VI. In these cross-sections, the higher deviation suggested that the numerical model estimated a narrower main stream as the measurements display. Around these cross-sections, the varying geometry increasec the transverse flow structure. However, a known limitation of the time-averaged RANS modeling was just this phenomenon, which is the underestimation of the secondary flow strength [56–58]. Next to this section, the flow model could be considered as validated for higher flood waves too, taking into account the available measurement data. For a more comprehensive validation, fixed boat ADCP measurements would have been expedient, which could be used to obtain time-averaged velocity values and also to estimate local bed shear stresses. The model validation for the 100-year flood wave should be evaluated in this light.

#### **4. Results**

#### *4.1. Comparison of the Calculation Methods*

In order to confirm the better operation of the novel *Re\** dependent combined method, the bed change calculation by the van Rijn, Wilcock and Crowe and the combined method were compared.

The simulations were performed only for the d and e flood waves because of the significant computational time (see Figure 3d,e). The initial model setup (the flow field, water levels and the bed material) was given for each run from the results of the model runs for the first three flood waves by the *Re\** dependent combined method. The flood wave d was a relatively low (peak is around 3460 m3/s), but durable (~2 months long) flood wave, while the e was the historical one with a peak higher than 10,000 m3/s. Thus, the comparative analysis presented the operational characteristics of the sediment transport models for both the durable lower and also for the extreme water regimes. As a benchmark, the measured bed change map indicated the extent of the possible bed changes. However, the measured and calculated maps could not be compared directly, because the measured belonged to the whole two-year-long period (Figure 3).

The calculated bed change maps for the flood wave d (Figure 3) are presented in Figures 6–8.

**Figure 6.** Calculated bed changes by the van Rijn (*vR*) formula for a 2.5 month-long period (Figure 3).

**Figure 7.** Calculated bed changes by the Wilcock and Crowe (*W&C*) formula for a 2.5 month-long period (Figure 3).

**Figure 8.** Calculated bed changes by the combined method for a 2.5 month-long period (Figure 3).

As the results show, the *vR* model estimated unrealistic changes both spatially and in magnitude. The unrealistically huge bed changes suggest that the *vR* model did not seem to be an appropriate model for the given Danube reach, particularly not for bed change calculation in the main channel.

The *W&C* sediment transport model estimated a more stable bed surface than the *vR* (Figure 7). In this particular case, the bed surface seemed so resistant that the mean flow field was too weak to cause any significant bed changes. The motion of the very fine, basically suspended inlet load was calculated by the *W&*C model as bed load. Therefore, that part of the inlet sediments settled progressively along the channel. However, because of the quite low suspended load [51], the bed level rise caused by sedimentation was negligible (<0.005 m).

Figure 8 shows the bed changes calculated by the combined method. The red lines illustrate the borderline which separated the areas where the *vR* or the *W&C* model was activated in the initial moment of the model run. Accordingly, it can be seen that the *vR* formula was invoked at the near-bank areas, at the groin fields and also at a smaller part of the Vének lower gravel bar. At these regions, more significant (~0.05 m) sedimentation was estimated. That is, at these less hydraulically rough parts of the river bed, the deposition of both the finer bed load and suspended load were expected, which areas can be detected by the *Re\** [22]. In turn, in the navigation channel, no considerable bed change happened. According to the suspended form of the *vR* formula, the finer suspended load passes over the calculation domain, while the bed surface remains still, calculated by the *W&C* formula. This assumption is consistent with the conclusions of the field measurements [19]; the main channel seems to be armored enough to be resistant at the mean water regime.

The other flood wave, for which the comparative analysis of the three calculation methods (*vR*: Figure 9, *W&C*: Figure 10 and combined: Figure 11) was established was the historical flood wave from 2013 (Figure 3). Regarding this hydrological case, the *vR* model estimated also an unrealistic bed change map (Figure 9.). This was mainly true for the main channel. There, the motion of the coarser grains was probably overestimated, resulting in huge erosions and depositions. In turn, at the near-bank regions, at gravel bars and at the groin fields, the changes seemed to be partly in the expectable order of magnitude. However, it is clearly visible that the *vR* formula is not an appropriate choice for the morphological change calculation of such a complex river reach.

**Figure 9.** Calculated bed changes by the *vR* formula for the historical flood wave (Figure 3).

**Figure 10.** Calculated bed changes by the *W&*C formula for the historical flood wave (Figure 3).

**Figure 11.** Calculated bed changes by the combined method for the historical flood wave (Figure 3).

The *W&C* model calculated more realistic bed changes (Figure 10), especially in the main channel. However, the measurements at the lower gravel bar and also at the whole main channel showed significant (Δ*z* > ±0.2 m) changes. Based on these, the *W&C* formula likely overestimated the stability

of the channel. A remarkable bed level increase can be pointed out, which moderated towards the downstream direction. These bed level changes could be explained by the settling of the inlet finer load, which could not be taken into account as suspended load. The publication of Török et al. [59] pointed out that in the case of mixed bed content, the Shields diagram predicts lower critical bed shear stress for the bed load of the finer, sand particles, than the reference shear stress of the *W&C* model. Accordingly, the *W&C* model estimates respectively higher stability for the sand particles, than the Shields diagram and thus the *vR* model. This leads to such kind of deposition pattern in the main stream which is unrealistic compared to the measured bed changes map [19].

The combined method predicts (Figure 11) more significant bed changes, compared to the one resulted by the *W&C* model. The red borderline suggested that as at the flood wave d, the combined method calculated the sediment transport by the *vR* formula in the near-bank regions. During the flood wave e, the remarkable erosion at the groin field B meant that the groin field got flushed and the earlier deposited finer sands got eroded. Also, notable changes took place at the vicinity of the gravel bar at region D. Here, the widening of the downstream sides of the gravel bar can be seen, in accordance with the measured bed change map in Figure 4. Around the groin pair at the left bank, on the opposite side of the gravel bar, the blue spots refer to erosion. This process was probably the result of a similar, flushing process like the one which happened at the upstream groin field. In turn, in the main and navigation channel, no considerable bed change happened. According to the suspended form of the *vR* formula, the finer suspended load passed over the calculation domain, while the bed surface remained still, calculated by the *W&C* formula. The conclusions of the field measurements [19] referred also to the resistant main channel.

The results indicate that the interactions between the sediment transport at different channel sections (groin fields, gravel bars and main channel) could not be estimated by the *vR* or *W&C* formulas. However, the expedient combination of them gives an opportunity to deal with the interaction-mechanism between the local- and reach-scale processes.

Analyzing the calculated bed changes in the marked boxes, the following assumptions can be stated. In this part, the results of the *vR* model were skipped because of the unrealistic bed change calculations. At region A, the real bed level deepening (Figure 4) could not be reproduced by any method. In region B, only the combining method was able to predict significant depositions and erosions. Accordingly, the depositions probably took place during the lower water regime, while the bed level incision occurred during the flood waves. Thus, the measured bed level changes during the two years were formed most likely indeed during the whole two-year-long period. In region C, all two model results suggested that the gravel bar was in a stable state. Besides region A, the bed level increase in the main channel between region C and D (Figure 4) could not be pointed out by any sediment transport formula. Finally, in region D the model results showed that the widening of the gravel bar and also the erosion at the vicinity of the groin pair occurred rather during the higher flood wave than during the slighter flood waves, or mean water regime.

#### *4.2. Measured Data-Based Verification of the Combined Method*

A quantitative assessment of the tested sediment transport formulas was performed based on the results of the comparative analysis. First, in accordance with the conclusions of the field measurement based investigation, it was assumed that the measured erosion and deposition took place mainly during higher water regimes in region D, when the flow discharge was higher than the bed forming discharge [19,23]. Thus, the measured and calculated data could be considered as indirectly comparable in this region. Therefore, from the measured bed level change maps, the total volume of the erosion and deposition could be calculated. Furthermore, counting the number of days of the higher water levels, the average daily rate of the volume of both the erosion and deposition could be estimated. Likewise, based on the calculated bed change maps and knowing the duration of the historical flood wave, the daily average volume changes of the deposition and erosion could be estimated.

Table 2 shows data about these volumes. The table presents the ratio of the calculated volume to the measured, regarding the deposition and erosion separately, for each sediment transport formula. A value of 1 would indicate a perfect match to the measured volume change. Table 2 shows that the *vR* model was the one which overestimates most both the deposition and erosion volumes. The *W&C* model calculated the deposition rate more accurately, but still indicated more than the measured. This was partly explained by the lack of the suspended sediment calculation. In turn, the *W&C* model estimated negligible low erosion. In total, the more reliable results were provided by the combined method. With this, the erosions at the near-bank parts were calculated more accurately by the *vR* model, resulting in sediment feed for region D. Thus, because of the capturing of the coming sediments, the widening of the deposition could be better represented. As the *Re\** dependent criterion activated the *vR* formula at the groin pair, the bed level incision at its vicinity was also estimated better.

**Table 2.** The average daily volume changes for calculated (Δ*Vc*) and measured (Δ*Vm*) volume ratio values Δ*Vc*/Δ*Vm* for region *D*. A value of 1 would indicate the perfect match.

Even though the measured and calculated bed changes could not be compared directly, the nature, the magnitudes and the locations of the remarkable bed changes suggested the greater aptitude of the combined method.

As it was discussed in Section 3.2.1, the calculated bed changes were elaborated only for the higher flow discharges (>2100 m3/s, Figure 3), so the results could not be compared with the measured changes directly. Therefore, to achieve a notionally common scale, the bed change values were normalized. That is both the measured and calculated bed changes got divided by the highest bed level decrease or increase the value of the main channel:

$$
\Delta z\_{\text{norm}} = \left\{ \begin{array}{c} if \, \Delta z > 0 \to \frac{\Delta z}{\Delta z\_{\text{max}}} \\ \epsilon lse \to \frac{\Delta z}{|\Delta z\_{\text{min}}|} \end{array} \right\}. \tag{3}
$$

Thus, the occurring values develop between −1 and 1, where 1 indicates the maximum deposition height along the river reach, while −1 belongs to the biggest erosion (Δ*zmax*, *meas* = 1.5 *m*, Δ*zmin*, *meas* = 1.8 *m*; Δ*zmax*, *calc* = 0.5 *m*, Δ*zmin*, *meas* = 0.25 *m*). The measured bed changes of both the erosions and depositions were consequently higher than the calculated. That is, the numerical model underestimated the magnitudes of the bed changes. This could be partly explained by the ignoring of the third part of the annual bed load in the numerical model estimation.

By the comparison of the measured (Figure 12) and calculated normalized bed changes (Figure 13) the following remarks can be stated. The modeled bed changes did not represent the remarkable erosion at region A. It is noted that this difference also meant a sediment supply loss for the downstream in the model calculation. Since the model did not manifest any bed level decrease, the bed material was probably finer here than it was set in the model. At region B, the measured scours appeared in the calculated results (blue spots). However, not as concentric scours, but rather as lengthwise formations. The numerical model represented depositions close to the measured magnitude at region B also. However, their location iwa not accurate; the brown spots occur between the groins, instead of in the front of the groins, like in Figure 12.

**Figure 12.** Normalized bed changes of the measured values (Figure 4) regarding the whole period October 2012–October 2014.

**Figure 13.** Normalized bed changes of the calculated values. The calculation was elaborated for the eight flood waves in the period October 2012–October 2014 (Figure 3).

At the gravel bar in region C, negligible bed changes were measured. Likewise, the numerical model predicted a stable bed surface. Finally, the combined method pointed out the growing of the lower part of the Vének lower gravel bar (region D, Figure 13). The measured bed changes indicated two, separable depositions. These double depositions were also represented by the model results. Like the measured changes, the model also calculated scouring in front of the left bank groin pair. However, the extension of it was considerably lower than the scale of the measured deepening. In turn, the lengthwise deposition between the groins (between regions C and D) was also indicated in both measured and calculated maps. However, the calculated deposition occurred in a significantly lower range and formed at the right bank side and not in the main stream. Also, an important match was that neither the measured nor the calculated values suggested any essential large-scale bed changes in the main channel.

A significant difference between the measured and calculated bed changes happened in region A. Here, the effect of the potential error in the initial bed material was further examined. An investigation

was performed, which was based on the assumption that the initial bed material was set inaccurately around region A. The bed material samples around *d*50 ≈ 0.01 m. Considering the grain-size distributions [19], a still realistic, but considerably finer bed material was presupposed. Therefore, the model was set up by 30% lower *d*50, which was *d*50 ≈ 0.007 m. With this only one difference, the model was run for the historical flood wave. Figure 14 presents the bed changes at region A, and at the downstream of it.

**Figure 14.** Calculated bed changes by the combining method for the historical flood wave. The model in the middle figure was set up with the initial, while in the right figure with finer bed material.

The right side of Figure 14 shows the bed changes in the case of finer bed material. It can be seen that the 30% decrease of the *d50* resulted in major erosion at region A. Considering the measured changes in the left figure it is obvious that the decreasing of the *d50* led to a better match to the real bed changes. The lower row of Figure 14 represents the bed changes at the downstream. Here, important deposition formations could be measured (left Figure 14) in front of the right bank groin pair (region B) and also in the main stream, between the two gravel bars (between region B and D). These changes could not be represented by the original model setup (middle Figure 14). In turn, in the case of finer bed material (right Figure 14), the model predicted important depositions in these regions. The extension of the deposition in front of the groin pair (region B) was very similar to the measured. Also, the lengthwise extension of the deposition downstream (between regions B and D), in the main stream were also reproduced. However, the location of it was not correct. It seems that the model underestimated the crosswise sediment transport, which is a known limitation of the Reynolds averaged description of the flow field [56–58,60]. Concluding, the herein presented investigation suggested that the bed material at region A was finer than the predicted *d50* allocation for the original model setup. With this assumption, the deposition nature at the downstream has become also detectable.

#### **5. Discussion**

There are existing proven models in the literature that work reliably for given morphodynamic cases (e.g., the van Rijn model for hydraulic smooth regimes, which mainly occurs in clear sand bed; or the Wilcock and Crowe model for hydraulic rougher regimes, which develops in coarser bed surfaces). However, there are river reaches, within which hydraulic smoother and rougher zones can form, resulting in a completely different type of sediment transport, such as bed armoring and silt aggradation parallelly. For such cases, Török et al. developed a novel sediment transport calculation method and validated it against laboratory measurements.

In this study, its validation against field measurements was introduced. The results showed that the combined application of the Wilcock and Crowe and the van Rijn models could significantly increase the precision of numerical bed change calculations. The results showed that the combined method calculated both erosion and deposition the most reliably. It estimated the volume of the scouring by order of magnitude more precisely than the other two models. Additionally, the results also evinced that the location, extension and shape of the bed formations (e.g., scours, deposition, bar evolvement) could be calculated more reliably by the novel combined approach. These results were consistent with the laboratory measurements based on validation of the calculation method.

The sediment transport interaction between the sand and gravel dominated parts of a river reach played an essential role in the bed changes. For instance, the flushing of the silt dominated groin fields can feed the downstream gravel bars, leading to the growth of the bars. The application of the combined approach really makes sense, within river reaches, where such well-separated morphodynamic situations prevail (e.g., sand aggradation in groin fields; bed armoring in the main stream). In such cases, the local-scale morphodynamic processes are calculated by the proven models, while the interaction between them leads to a more accurate reach-scale calculation. The comparison of the numerical simulations to the measured data presented that such reach-scale interactions can be taken into account by the novel combining method.

#### **6. Conclusions**

The validation against both laboratory and field measurements suggests that not necessarily the development of completely new models can lead to the evolution of the sediment transport calculation. This study highlights that the optimally combined use of sediment transport models can be a promising alternative in the sediment transport calculation. In this context, it is important to state that regardless of whether a sediment transport model is applied in combination or simply, the determination of its applicability limits by any morphodynamic variables in a well-defined way is vital. In this study, a shear Reynolds number based description was applied which allowed for the authors an easily described combination criteria. However, the essence of the combining theory does not insist on any variable for the applicability limit description. However, the herein presented example highlights the need for such a kind of description method and also the potential of the parallel use of any sediment transport models.

**Author Contributions:** Conceptualization, methodology, software, validation, writing—original draft preparation, visualization, G.T.T., resources, data curation, G.T.T. and S.B.; formal analysis, writing—review and editing, S.B.; supervision, project administration, funding acquisition, S.B. and J.J.

**Funding:** This research received no external funding.

**Acknowledgments:** This research was supported by MTA TKI (Office for Research Groupsof the Hungarian Academy of Sciences). Support of grant BME FIKP-VÍZ by EMMI (Ministry of Human Capacities) is also kindly acknowledged. The authors acknowledge the funding of the OTKA (Hungarian Scientific Research Fund) FK 128429 grant. The research was partly supported by the ÚNKP-19-4 New National Excellence Program of the Ministry for Innovation and Technology. The last author acknowledges the funding from the János Bolyai fellowship of the Hungarian Academy of Sciences.

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

#### **Appendix A**

#### *Appendix A.1 Wilcock and Crowe Bed Load Transport Formula*

The well-known surface-based bed load formula of *W&C* [8] was developed for mixed gravel/sand sediments. The gravel ranged in size from 2.0 to 64 mm, the sand from 0.5 to 2.0 mm. They conducted laboratory experiments in a tilting laboratory flume, in which both the water and the sediments were recirculated. Five different sediment mixtures were investigated, varying the sand content of the sediment mixtures from 6.2% up to 34.3% [60]. The mixture grain-size distribution used in the experiments of Török et al. [34] and of Wilcock et al. [60] shows a good match. In case of the experiments of Wilcock et al. [60], *dm* grain size of the initial bed materials was in the range of 4.1–10.5 mm, whereas *dm* was 8.7 mm the measurements of Török et al. [34]. The *W&C* formula is therefore expected to be capable of describing the gravel/sand bed load transport processes in our experiments. In the case of such non-uniform bed material, the interaction between (hiding and exposure) the particles of different sizes plays an essential role in the stability of the sediments and can result in bed armoring. The *W&C* formula is reported to be capable of predicting the transient conditions of bed armoring or scouring processes. The preliminary numerical simulations of Török et al. [46] also supported this statement.

The transport formula uses the collapse similarity hypothesis, which estimates the fractional sediment transport rate (*wi\**) as the function of the ratio of the shear stress (τ) and the reference shear stress (τ*ri*). The reference shear stress τ*ri* is the value of τ at which *wi\** is equal to a small but already perceptible value *wi\** = 0.002. That is the reference shear stress has a similar physical meaning to the critical shear stress (τc), which indicates the initiation of motion of particles. The dimensionless transport rate is obtained according to the following equation:

$$\mathcal{W}\_i^\* = \begin{cases} 0.002 \theta^{7.5} & \text{if } \theta < 1.35\\ 14 \Big( 1 - \frac{0.894}{\theta^{0.5}} \Big)^{4.5} & \text{if } \theta \ge 1.35 \end{cases} \tag{A1}$$

where

$$
\theta = \pi/\tau\_{\dot{n}}.\tag{A2}
$$

The volumetric transport rate per unit width *qbi* can be expressed from the following formula:

$$\mathcal{W}\_i^\* = \frac{(s-1)\mathcal{g}q\_{hi}}{\mathcal{F}\_i \mu\_\*^3},\tag{A3}$$

where *Fs* is the proportion of sand in surface size distribution.

For the estimating of τ*ri*, the *W&C* model suggests that the stability of the sediment particles depends on the sand content. Thus τ*ri* can be calculated by the following steps:

$$
\tau\_{rm}^\* = 0.021 + 0.015e^{-20F\_S},\tag{A4}
$$

where τ\* *rm* is the dimensionless reference shear stress for the mean size of the bed surface. The reference shear stress for the mean size can be obtained from the particle Froude number equation:

$$
\tau\_{rm}^\* = \frac{\tau\_{rm}}{\mathcal{g}(\rho\_S - \rho\_W)d}. \tag{A5}
$$

Then, the last step is the estimation of the reference shear stress for the particle sizes separately:

$$\frac{\pi\_{ri}}{\pi\_{rm}} = \left(\frac{D\_i}{D\_{sm}}\right)^b,\tag{A6}$$

where

$$b = \frac{0.67}{1 + e^{(1.5 - \frac{D\_i}{D\_{\rm sun}})}}.\tag{A7}$$

#### *Appendix A.2 Van Rijn Bed Load Transport Formula*

The development of the *vR* formula is based on laboratory and field measurements. Van Rijn presented a mathematical method for the computation of bed load transport rate. The formula is basically the numerical solution of the equations of motions for a solitary particle. The mathematical model was calibrated based only on the experiments of Fernandez Luque [61,62], where the motions of uniform particles (*d* = 0.18 mm) under steady flow (*u\** = 0.04 m/s) were monitored.

The lifting term of the equation of motion for a particle is derived as [63]:

$${}\_{L}F\_{L}(shar) = \alpha\_{L} \rho \nu^{0.5} D^{2} v\_{r} \left(\frac{\partial u}{\partial z}\right)^{0.5},\tag{A8}$$

where *FL* is the lift force, α*L* is the lift coefficient, ν is the kinematic viscosity coefficient, ρ is the density of the fluid, *vr* is the relative particle velocity and ∂u/∂z is the velocity gradient. The mathematical method uses the turbulent wall law Equation (A9) for calculation of the vertical flow velocity distribution:

$$
\Delta U(z) = \frac{u\_\*}{\kappa} \ln \left( \frac{z}{z\_0} \right) \tag{A9}
$$

where *U(z)* is the horizontal velocity as the function of the *z*, *z* is the level above the bed, *u\** is the bed shear velocity, κ is constant von Karman (κ = 0.41), and *z0* is the zero-velocity level above the bed.

The volumetric bed load transport rate is defined as:

$$\frac{q\_b}{[(s-1)g]^{0.5}D\_{50}^{1.5}} = 0.053 \frac{T^{2.1}}{D\_\*^{0.3}}\,\text{\AA}\tag{A10}$$

where *qb* is the bed load transport per unit width, *D50* is the particle size, *s* is the specific density, *g* is the acceleration of gravity, *T* is the transport stage parameter, can be represented as:

$$T = \frac{\left(u\_\star'\right)^2 - \left(u\_{\star,cr}\right)^2}{\left(u\_{\star,cr}\right)^2},\tag{A11}$$

where *u\*,cr* is the critical bed shear velocity according to Shields [35], *u* ∗ = *g*0.5/*C u* is the bed shear velocity related to grains, where *C´* is the Chézy-coefficient related to grains and *u* is the mean flow velocity.

In the Equation (A10), *D\** is the particle parameter defined as:

$$D\_\* = D \text{50} \left[ \frac{(s-1)g}{\nu^2} \right]^{1/3} \tag{A12}$$

Finally, the *vR* formula was verified based on 56 field measurements and 524 flume data. The bed load formula is recommended to use for particles in the range of 200–2000 μm (i.e., below 2 mm). As such the formula is not expected to reliably calculate the motion of gravels and the interaction between particles of different sizes and to describe such a process, like bed armoring. However, the model can presumably describe the transport of finer grains.

#### **References**

© 2019 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).


*Article* **Retrospective Analysis of FED Method Treatment Results in 11–17-Year-Old Children with Idiopathic Scoliosis**

**Sandra Trzci ´nska 1 and Kamil Koszela 2,\***

**Abstract:** (1) Background: Idiopathic scoliosis is a major treatment problem due to its unknown origin and its three-dimensional nature. Attempts to cure it and search for new methods of physiotherapeutic treatment that would lead to its correction are one of the key issues of modern medicine. One of them is the fixation, elongation, de-rotation method (FED), used in the conservative treatment of idiopathic scoliosis. The aim of the study was evaluation of the short-term effectiveness of the FED method in the treatment of patients with idiopathic scoliosis. (2) Methods: Each patient underwent therapy based on the guidelines of the FED method. Patients were tested with the Bunnell scoliometer and the Zebris computer system. The treatment period was three weeks, after which the examinations were repeated. (3) Results: The results appeared to be statistically significant for all tested variables. (4) Conclusions: The examinations showed that the FED method had a statistically significant effect on the improvement of all parameters of posture examination, regardless of the size of the scoliotic deformation angle and bone maturity.

**Keywords:** scoliosis; musculoskeletal disorders; spine; conservative treatment; rehabilitation

### **1. Introduction**

Idiopathic scoliosis is a major treatment problem due to its unknown origin and its three-dimensional nature. Attempts to cure it and search for new methods of physiotherapeutic treatment that would lead to its correction are one of the key issues of modern medicine.

One of the less known methods used in Poland for the treatment of patients with idiopathic scoliosis is the FED method (fixation, elongation, de-rotation) established in Spain. The FED treatment method uses a device where corrective forces act on the curvature. The strength of the device is focused on stabilizing, stretching, and de-rotating the spine, under the control of an innovative computer program [1]. As a method that uses a special apparatus for the correction of scoliotic deformity, it has been of great interest to medics in recent years. However, reports that assess the effects of its application are still insufficient [2].

There are several ways to test the effectiveness of the FED method in treating patients with idiopathic scoliosis. The most popular, but harmful due to radiation exposure, is a radiological examination, which is still the standard for diagnosing scoliosis [3–5]. Scoliosis patients can get 10 to 25 spinal X-rays over several years equating to a maximum 10 to 25 mGy of cumulative exposure. Patients who were diagnosed at a younger age and received early and ongoing treatment may be subjected to up to 40 to 50 X-rays, 50 mGy in total [6]. Other non-invasive devices have been used more and more frequently to assess the effects of therapeutic activities undertaken, including the assessment and comparison of various treatments and methods. Such a test can be both a clinical examination, including a scoliometer test as a simple and reliable device for measuring the transverse plane of the

**Citation:** Trzci ´nska, S.; Koszela, K. Retrospective Analysis of FED Method Treatment Results in 11–17-Year-Old Children with Idiopathic Scoliosis. *Children* **2022**, *9*, 1513. https://doi.org/10.3390/ children9101513

Academic Editors: Axel A. Horsch, Maher A. Ghandour and Matthias Christoph M. Klotz

Received: 28 August 2022 Accepted: 30 September 2022 Published: 3 October 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

spine-trunk rotation, and a more modern three-dimensional posture test, such as the Zebris system. The Zebris system is a modern and specialized apparatus that uses ultrasound, enabling non-invasive examination of the body posture by creating a three-dimensional image of the patient's figure. The system performs a computer analysis to create a report [7]. While the scoliometer is traditionally used to measure trunk rotation at the apex of the curvature, it can also be used to assess global trunk rotation, which evaluates the overall impact of therapy on the spine. A similar test assessing the overall effect of the therapy, but in the frontal plane, is the examination of scoliotic deformity with the Zebris system.

The aim of the study was to evaluate the results of the FED method in patients with idiopathic scoliosis in the short term.

#### **2. Materials and Methods**

#### *2.1. Study Population*

The study included 81 subjects, 72 girls and 9 boys, aged from 11 to 17 (mean 14.28 ± 1.63). Each of the subjects had a current radiograph, which assessed the following: Risser test, the size of the Cobb angle for individual curvatures and the type of scoliosis was marked according to the King-Moe classification. All patients had idiopathic double-curve scoliosis of type I or II, characterized by the presence of a sigmoid curve in the thoracic (type I) or in the lumbar (type II) segments greater than the other one.

#### Inclusion Criteria

#### *2.2. Study Protocol*

Each patient underwent therapy based on the guidelines of the FED method, which consisted of three basic elements: physical therapy as well as analytical and instrumental kinesiotherapy. The main component of the treatment was a special device that corrected the spine in 3 planes. With the help of a special vest, the patient was suspended in the device. Elongation was performed by a computer-controlled hoist, which, at the same time, regulated the pressure of the mobile arm correcting the apex of the scoliotic curve. Other arms stabilized the scoliotic curve at its ends. The pressure force was determined depending on the patient's ability, up to a maximum of 100 kg. The time of the procedure was 30 min, the time of corrective pressure performed by the pneumatic movable arm was 20 s, and the break was 10 s. The arm corrected the curve both in the frontal and rotational planes, owing to the possibility of its angular positioning. In order to prepare the patient for the therapy in the device, the tissues were made more flexible and blood supply improved in the places to be subjected to the therapy, so, in this study, electrostimulation of the muscles on the convex side of the curve and thermal treatment, in which warm compresses were placed in the deformation concavity, were used. Both procedures lasted 15 min. Then the patient performed exercises for about 20–30 min, individually selected in accordance with the guidelines of the FED method. The selection of exercises was based on the King-Moe classification, which divides scoliosis into 5 types in terms of the

location of scoliotic curvatures. In addition, patients wore the Boston brace every day for approximately 21–22 h a day, except for FED therapy (up to 3 h) and personal care. The Boston brace was made on the basis of a plaster cast and self-report was used.

Each participant was tested with the Bunnell scoliometer and the Zebris computer system on the day before therapy began. The treatment period was 3 weeks, after which the examinations were repeated on the day the therapy completed. In the study, the scoliometer assessed both the trunk rotation angle at the apex of both scoliotic curves for the thoracic (ATR Th) and lumbar (ATR L) spine, and the total spine rotation using the SDR summing parameter, which consisted in summing the values of the rotation on both curves as positive values regardless of their direction. The computer examination assessed overall scoliotic deformation (SD) in the frontal plane. This parameter was the sum of the angles of tangents from the seventh cervical to the fifth lumbar vertebrae (C7-L5).

This study was conducted according to the guidelines of the Declaration of Helsinki, and approved by the Bioethics Committee for Scientific Research at the College of Rehabilitation in Warsaw, number 100/2022.

#### *2.3. Data Analysis*

In order to answer the research questions, statistical analyses were performed with the IBM SPSS Statistics 27 (Armonk, NY, USA). This was used to analyze basic descriptive statistics, with Shapiro–Wilk test, and Student's t-test for dependent samples, and Wilcoxon test and Spearman's rho correlation analysis. The level of significance was considered to be α = 0.05.

#### **3. Results**

#### *3.1. Data Analysis*

Most patients (66.7%) had type II scoliosis according to the King-Moe classification. The mean Cobb angle at the thoracic level [◦ ] was 35.91 ± 10.43, and at the lumbar level it was 33.54 ±10.94. Detailed results are presented in Tables 1 and 2.

**Table 1.** Gender and type of scoliosis.

*n*—number.

**Table 2.** Age and X-ray analysis of indicators.

*n*—number, *x*—mean, *SD*—standard deviation, *Min*—the lowest value, *Max*—the highest value.

#### *3.2. Basic Descriptive Statistics with the Shapiro-Wilk Test*

In the first step of the analysis, the distributions of quantitative variables were checked. For this purpose, basic descriptive statistics were calculated together with the Shapiro-Wilk test examining the normality of the distribution. The results of the analysis are presented in Table 3.

**Table 3.** Basic descriptive statistics of the studied variables together with the Shapiro-Wilk test for normality.

*n*—number, *x*—mean, *SD*—standard deviation, *Me*—median, *Sk*—skewness, Kurt—kurtosis, *Min*—the lowest value, *Max*—the highest value, *W*—Shapiro–Wilk test, *p*—level of significance, \*—statistical significance.

The results of the Shapiro-Wilk test appeared to be statistically significant for the lumbar rotation angle-ATR L [◦ ] and for the scoliotic deformation angle-SD [◦ ] both before and after the therapy. This meant that distributions of these variables differed from the normal distribution. However, in the case of the lumbar rotation angle-ATR L [◦ ], the value of the skewness did not exceed the absolute value of one, which indicated that the asymmetry was insignificant. Therefore, the analyses for the trunk and lumbar rotation angles-ATR Th [◦ ] and AR L [◦ ], and for the sum of two rotations-SDR [◦ ], were performed based on parametric tests. Yet, for the scoliotic deformation angle-SD [◦ ], the skewness exceeded the value of one. A detailed analysis showed that this value of skewness resulted from the presence of two outliers (+ 3SD). Thus, for these variables, the analyses were based on non-parametric tests.

#### *3.3. Comparison of the Value of the Trunk Rotation Angle-ATR Th [*◦ *], the Lumbar Rotation Angle-ATR L [*◦ *], the Sum of Two Rotations SDR [*◦ *] and the Angle of Scoliotic Deformation-SD [* ◦ *], before and after the Therapy*

In the next step, it was checked whether the applied therapy influenced values of the trunk rotation angle-ATR Th [◦ ], the lumbar rotation angle-ATR L [◦ ], the sum of two rotations-SDR [◦ ] and the angle of scoliotic deformation-SD [◦ ]. For this purpose, the Student's *t*-test was performed for dependent samples and, in the case of the scoliotic deformation angle, the non-parametric Wilcoxon test was performed. The analyses were performed for all patients, taking into account the division into gender and type of scoliosis, based on the King-Moe classification (Table 4).

The results appeared to be statistically significant for all tested variables. The values of individual parameters, the trunk and lumbar rotation angles, the sum of two rotations and the scoliotic deformation angle, were significantly lower after the therapy than before. Differences in measurements were observed in both girls and boys, and in subjects with scoliosis types I and II. Each of the differences between the measurements before and after the therapy was significant.

**Table 4.** Results of Student's *t*-test for dependent samples and Wilcoxon's test comparing individual parameters measured with a scoliometer and Zebris, before and after therapy.

*n*—number, *x*—mean, *SD*—standard deviation, *t*—*t*-test, *Z*—Wilcoxon test, *p*—level of significance, 95% *CI*—95% Confidence Interval, *LL*—Lower Limit, *UL*—Upper Limit, *d Cohena*—effect size for Student's *t* test for dependent samples, *r*—effect size for Wilcoxon test, \*—statistical significance.

#### *3.4. Correlations between Parameters Measured with X-ray, Scoliometer and Scolioscan*

In the next step, it was verified whether the measurements made with the use of a scoliometer and Zebris computer system correlated with each other. For this purpose, Spearman's rho correlation analyses were performed. The parameters were compared before and after therapy and for all observations in general. The results are presented in Table 5.

The analysis showed six statistically significant correlations. For the measurements performed before the therapy, a significant correlation was observed only in the case of the sum of two rotations-SDR [◦ ]. In the case of measurements performed after the therapy, statistically significant correlations concerned the trunk rotation angle-ATR Th [◦ ] and again the sum of two rotations-SDR [◦ ]. On the other hand, the analysis performed for all measurements in total showed that the angle of scoliotic deformation-SD [◦ ] was significantly related to each of the results obtained with the use of the scoliometer. Parameters that were

measured with the scoliometer were positively related to the scoliotic deformation angle, which was measured with the Zebris computer system. Three of the statistically significant correlations were weak and three were moderately strong.

**Table 5.** Correlation analysis for measurements made with a scoliometer and Zebris computer system.

*n*—number, \*—statistical significance.

#### *3.5. Comparison of the Value of the Trunk Rotation Angle-ATR Th [*◦ *], the Lumbar Rotation Angle-ATR L [*◦ *], the Sum of Two Rotations SDR [*◦ *] and the Angle of Scoliotic Deformation-SD [* ◦ *], before and after the Therapy Depending on the Dimensions of the Larger Curve*

In the next step, it was verified whether, as a result of the therapy, in the groups distinguished by the dimensions of the larger curvature, the individual measured parameters changed, i.e., the trunk rotation angle-ATR Th [◦ ], the lumbar rotation angle-ATR L [◦ ], the sum of two rotations-SDR [◦ ] and the angle of scoliotic deformation-SD [◦ ]. The analyses were carried out in groups distinguished according to the value of the larger Cobb angle (at the level of the thoracic or lumbar spine). It was assumed that I◦ results were up to 24◦ , II◦ , 25–45◦ , and III◦ , above 45◦ . Due to the very small number of individuals with I◦ (*n* = 4), the analyses were performed only for those with II◦ and III◦ . The Student's *t*-test was used for dependent samples, and, in the case of the scoliotic deformation angle, the non-parametric Wilcoxon test was used (Table 6).

The results were statistically significant for all variables in every group. The values of trunk and lumbar rotation angles [◦ ] (ATR Th and ATR L), sum of two rotations [◦ ] and scoliotic deformation angle-SD [◦ ] were lower after the therapy compared to those before the therapy. In both groups and for all variables, the observed differences were strong.

**Table 6.** Results of Student's *t*-test for dependent samples and Wilcoxon's test-comparison of individual parameters measured with the scoliometer and Zebris, before and after therapy, with division into groups distinguished on the basis of the value of greater Cobb angle.

**Table 6.** *Cont.*

*n*—number, *x*—mean, SD—standard deviation, *t*—*t*-test, *Z*—Wilcoxon test, *p*—level of significance, 95% *CI*—95% Confidence Interval, *LL*—Lower Limit, *UL*—Upper Limit, *d Cohena*—effect size for Student's *t* test for dependent samples, *r*—effect size for Wilcoxon test, \*—statistical significance.

*3.6. Comparison of the Value of the Trunk Rotation Angle-ATR Th [*◦ *], the Lumbar Rotation Angle-ATR L [*◦ *], the Sum of Two Rotations SDR [*◦ *] and the Angle of Scoliotic Deformation-SD [* ◦ *], before and after the Therapy Depending on the Risser Sign Grade*

The results were statistically significant for all variables in every group. The values of trunk and lumbar rotation angles [◦ ] (ATR Th and ATR L), sum of two rotations [◦ ] and scoliotic deformation angle-SD [◦ ] were lower after the therapy compared to those before the therapy. In both groups and for all variables, the values were lower after the therapy compared to those before the therapy independently of the Risser Sign Grade, the observed differences were strong (Table 7).

**Table 7.** Results of Student's *t*-test for dependent samples and Wilcoxon's test-comparison of individual parameters measured with the scoliometer and Zebris, before and after therapy, with division into groups distinguished on the basis of Risser Sign Grade.

**Table 7.** *Cont.*

*n*—number, *x*—mean, *SD*—standard deviation, *t*—*t*-test, *Z*—Wilcoxon test, *p*—level of significance, 95% *CI*—95% Confidence Interval, *LL*—Lower Limit, *UL*—Upper Limit, *d Cohena*—effect size for Student's *t* test for dependent samples, *r*—effect size for Wilcoxon test, \*—statistical significance.

#### **4. Discussion**

The FED method is a relatively little-known method of treatment. It originated in Spain, where most of the scientific reports on its effectiveness come from [8–11]. Conducted studies have concerned mainly single cases or studies in smaller groups of patients. At the time when the method became fairly common in Poland, projects concerning larger research groups began to appear [12–14]. Today, many reports assessing its effectiveness in the treatment of patients with idiopathic scoliosis come from Poland. The number of centers using this type of treatment is constantly growing and the method has a large number of supporters. The reports that appear, especially in recent years, prove its effectiveness in the treatment of patients with idiopathic scoliosis [15].

Studies on the effectiveness of the FED method in the treatment of patients with idiopathic scoliosis show a significant improvement in all measured parameters, both in the examination with the scoliometer and with the Zebris system. These changes occurred both in the transverse and frontal planes. It turned out to be important to correct not only individual curves, but also the entire area of the spine in both tested planes. These studies indicated a positive effect on the entire spine, in contrast to the existing therapeutic view on the correction of one curve at the expense of deteriorating the other [16,17].

The use of traditional scoliometer examination and other posture parameters enable the assessment of the effect of the undertaken treatment procedures, mainly focusing not only on the diagnosis, but also on the effectiveness of individual methods and ways of the therapy. However, there are more and more reports indicating that basic examination (e.g., with a scoliometer) is insufficient. The impact of the therapy not only on the rotation of individual curves, but also on the global rotation, should also be evaluated. For this reason, summing parameters were used in the assessment of body posture. These are even less well known, but more and more often they are willingly used to evaluate therapy and its effect on the transverse plane [18–20]. A similar study, assessing the overall effect of therapy on the entire spine, in terms of the frontal plane, involved the scoliotic deformation parameter tested with the Zebris system. The conducted research demonstrated a correlation between these two parameters, despite the fact that they relate to two different transverse and frontal planes.

The project also investigated the impact of FED therapy on scoliotic deformations in terms of its size. The modified SOSORT angular division was used, which divides scoliosis on the basis of size assessed by the Cobb method on the radiograph. The studies showed a significant improvement in all parameters, regardless of the size of the scoliosis.

The examinations showed an improvement in all the tested parameters at various stages of bone maturity (Risser Sign Grade 1–4). Improvement was observed in all patients after the use of the FED method, regardless of the Risser Sign Grade.

The harmfulness of radiological examinations has more and more frequently led to the use of modern computer diagnostics to assess the posture [21–24]. These computer diagnostics allow for a harmless assessment of treatment stages, in terms of posture, especially since they correlate with radiographic images [25]. However, a distinction should be made here regarding when to use these two different examinations: radiological examinations should be used to diagnose and assess the progression of scoliotic deformity, and computer-based monitoring and evaluation should be used to diagnose and assess the treatment course. Computer-based methods are probably an alternative to radiology, but not in all respects. These two examinations differ from each other but should complement each other. Carrying out diagnostic activities on various grounds results in better control of scoliotic deformation progressing to impaired quality of life or to surgery, which often causes complications, and which we would like to avoid in the treatment of patients with idiopathic scoliosis [26].

Due to the nature of the therapy, the patients' follow-up was short-term. However, it would be worth performing repeated examinations at a longer time interval. The FED method is relatively new in Poland, so such studies should be carried out, especially related to the assessment of the impact of the method on the sagittal plane of the posture. It requires further analyses supplemented with long-term follow-up in a larger group of patients. A 3-month and 6-month follow-up examination is planned.

#### **5. Conclusions**

All assessed parameters of posture examination, both with the scoliometer and the Zebris system, showed a statistically significant improvement in patients treated with the FED method. Statistical improvement occurred both in boys and girls, as well as in all types of scoliosis. The analysis showed the occurrence of statistically significant correlations between the parameters of the posture examination with a scoliometer and the computer examination. The examinations showed that the FED method had a statistically significant effect on the improvement of all parameters of posture examination, regardless of the size of the scoliotic deformation angle and bone maturity. However, due to the short observation time, this method requires further research with a long follow-up period.

**Author Contributions:** Conceptualization, S.T.; methodology, S.T.; software, S.T.; validation, S.T., and K.K.; formal analysis, S.T. and K.K.; investigation, S.T.; resources, S.T. and K.K.; data curation, S.T.; writing—original draft preparation, S.T. and K.K.; writing—review and editing, S.T. and K.K.; visualization, S.T. and K.K.; supervision, S.T. and K.K.; project administration, S.T.; funding acquisition, S.T. and K.K. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Institutional Review Board Statement:** This study was conducted according to the guidelines of the Declaration of Helsinki, and approved by the Bioethics Committee for Scientific Research at College of Rehabilitation in Warsaw, number 100/2022.

**Informed Consent Statement:** Informed consent was obtained from all subjects involved in the study.

**Data Availability Statement:** The datasets analyzed during the current study are available from the corresponding author upon reasonable request.

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

#### **References**


**2. Traditional Monitoring Techniques Applied to Industrial Processes**

The use of sensors in the monitoring of manufacturing processes is consolidated in the manufacturing industry. Some monitoring techniques, such as the passive and active-passive, are well known. In this section, a brief review of researches that made use of such techniques is presented.


**1. Introduction**

Conventional flexible pavement, namely asphalt concrete (AC) pavement, suffers from many distresses including rutting, moisture damage, and fatigue failure after longterm operation, because of weather wearing and traffic volume growing worldwide [1,2]. With soaring demand of transportation, the costs of highway construction and pavement maintenance have increased dramatically over time. As a solution for these problems, Grouted Semi-flexible Pavement (GSP), as known as Resin-Modified Pavement, is a type of high-performance pavement that was derived in France in the 1950s. GSP consists of an open-graded matrix asphalt (OMA) mattress with a 20–35% void rate, grouted by high-fluidity cement mortar (HCM) [3–8]. GSP has great resistance to rutting, top-down cracking, oil corrosion, and fatigue damage, combining advantages of both asphalt and concrete [9–11]. It is cost-effective and has been used in many fields, which may be under unfavorable or complex traffic environments that are slow-speed, heavy-duty, or hightemperature, such as airport runways, factory field pavement, bus rapid transit lanes, and high-performance reclaimed asphalt pavement [12–14].

Construction of GSP is normally a two-phase operation. First, the OMA mattress is prepared and paved by equipment lighter than or equal to that used for AC. After the asphalt is cooled, HCM can be spread on the surface. Due to the good connectivity of voids in OMA, HCM penetrates the whole layer to obtain a very low residential void rate for GSP by rubber scrapers and light vibratory rollers. The GSP eventually forms the required strength for traffic through curing for a few days. The beneficial properties of this material rely on careful construction process control, and recent related challenges are problems for mechanical mechanism investigation and distress prevention, such as cracking.

In Europe, GSP was originally invented as a patent of heat-resistant pavement, namely "*Salviacim*", by French company Jean Lefebvre at Cognac airport in 1954 [3,15]. In 1987, two

**Citation:** Guo, X.; Hao, P. Influential Factors and Evaluation Methods of the Performance of Grouted Semi-Flexible Pavement (GSP)—A Review. *Appl. Sci.* **2021**, *11*, 6700. https:// doi.org/10.3390/app11156700

Academic Editor: Amir Tabakovic

Received: 22 June 2021 Accepted: 19 July 2021 Published: 21 July 2021

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2021 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

School of Highway, Chang'an University, Xi'an 710064, China<EMAIL_ADDRESS><EMAIL_ADDRESS>
companies, Densit and Phoenix, jointly developed and promoted the second generation of GSP, called "*Densiphalt*", of which the residual void rate was reduced after being grouted to improve in-field performance [15–17]. In 1989, Roffe, the researcher belonging to the Jean Lefebvre Company in France, published a monograph on the "*Salviacim*", which was still considered a special construction method by the public [3,18].

The U.S. Army issued a series of evaluation reports of GSP and practiced this technology on airport aprons and tank runways from 1976 [19–21]. The Departments of Transportation of states and the U.S. were also concerned for the performance of GSP and conducted many projects to prove its performance [22–24]. Consequently, the specification of GSP by the U.S. Army was established and applied in many fields, and were updated several times [25–28].

Wang et al., conducted a trial road of GSP on Huishen pavement in the Guangdong province in China in 1986 [29]. In 1995, Xu et al., introduced laboratory tests for GSP and presented the bending tensile strength as an indicator to control its quality [30]. Subsequently, Zhang and Pan adapted the Main Aggregate Filling (MAF) method and the orthogonal tests to design GSP in 2000 [31,32]. Further, Hao et al., proposed that the optimal asphalt content of OMA should be determined by the Cantabro test and the Schellenberg Binder Drainage test, while the low-temperature cracking of GSP was found to be the main problem affecting its operational life [4–6].

Although GSP was just used originally as a construction technique for special pavements such as airport aprons, it is now considered as a rather new technology because it is found to be strongly different from rigid pavement and flexible pavement for their working behaviors. With the theory of perpetual pavement arousing and traffic load increasing in the recent years, GSP has been welcomed into a rapidly developing period and has become a potential high-performance and cost-effective pavement which has a higher complex shear modulus than AC. However, standards of this technology are seldom published in the world. There are still many unsolved problems in GSP design, such as low-temperature cracking, unreasonable testing parameters, and undeveloped laboratory methods, etc., which need intensive research [6,33]. The poor understanding of these problems limited further development and application of GSP.

The objectives of this study are to summarize the findings in GSP over the past two decades and analyze the relationships of various factors and the methods of laboratory tests on the performance to identify knowledge gaps, promote this technology and deliver recommendations for future research. 143 papers have been selected and reviewed in total on three principles: first, the study should illustrate a qualitative and quantitative evaluation for GSP; second, the research should come from a new perspective for factors and properties of GSP; third, new findings and different conclusions are included for more discussion. Finally, relevant recommendations and future studies are proposed.

#### **2. Influential Factors on the Performance of GSP**

### *2.1. Microstructure*

Raw materials and structures of GSP are the essential elements affecting laboratory performance and also yielding technical problems. In a microscopic view, the relationships among asphalt, aggregate, and cement particles build the macro-strength of GSP. Indeed, the micro-analysis method is most commonly used in laboratory tests to demonstrate the intuitive influential factors of micromorphology on the behavior of GSP.

In recent years, resolution of micromorphology techniques, especially the Scanning Electron Microscope (SEM), has reached to the nanometer level, which can capture the microstructure views between different particles to analyze their relationships [34,35]. The GSP samples are observed in various instruments, resolutions, and observing ranges. The methods and results are shown in Table 1 and Figure 1.

**Table 1.** Methods of microscopic observation for GSP.

**Figure 1.** The GSP microscopic images observed by SEM: (**a**) Microstructure at 90-day age; (**b**) Interface and fibrous structures; (**c**) Wrapping phenomenon of asphalt film; (**d**) Generation of micro-cracking and bridging effect of the fiberlikenetwork.

GSF has higher density and solidity than traditional AC, which contributes to great complex shear modulus and compressive strength [36]. The strength of GSP is reinforced with increase of age (at ages of 7 d, 28 d, and 90 d), and the inner hydrated cement forms a fibrous-like network that penetrates through asphalt films, as shown in Figure 1 [39]. This exclusive phenomenon, refered to as three-dimensional lattices or a mosaic spatial network, bridges asphalt and cement, which increases the thickness of interface transition zones and enhances the bonding force between asphalt and cement [36]. Meanwhile, the un-hydrated cement plays a role as mineral powder to strengthen the adhesion between asphalt and aggregate.

The propagation of micro-cracks in GSP can be also investigated by SEM (Figure 1). The observation results indicate cracking emerges from the bottoms of samples and expands along with the interface between asphalt and cement [38,40]. This reveals that the asphalt-cement interface is the weak interface existing through the whole depth of GSP structure, approved by the Heavy Vehicle Simulator test and the Dynamic Cone Penetrometer test [41,42]. Additionally, the reaction of internal among raw materials is mainly physical, because no chemical or asphalt aging processes are identified at the interfaces by the Fourier Infrared Spectrum Analyzer (FISA) (200SXVFT-IR) [37].

In general, these characteristics of micro-structure determine the macro properties of GSP, which has high rutting resistance and low cracking resistance. Moreover, the micro-analysis method can effectively and intuitively gain the relationships of each of these components and raw materials. On the other hand, this method is also limited, owing to the subjectivity of corresponding researchers and the complexity of GSP's structure. For example, some different fibrous structures can be observed from GSP samples, but which has a better coherent strength between materials is unknown because it is hard to get quantitative analysis from SEM images.

### *2.2. Raw Materials and Admixtures*

GSP is a complex composite, combining multiple materials, such as asphalt binder, aggregate, cement, and admixtures. Generally, the effects of them can be separated into two independent steps: intermediate mixes (OMA and HCM) and final grouted or combined mix (GSP). In other words, the properties and the formation of raw materials and admixtures can directly determine the performance of intermediate mixes, which further affects the quality of the following GSP.

OAC is an asphalt mattress with high void rate similar to the open grade friction course (OGFC). According to design principles of OGFC, aggregate gradation can contribute to the volume of void and binder is able to provide inner cohesive strength in OAC. Likewise, strength and fluidity are two main properties of HCM, and are controlled by two factors: formation and types of raw materials.

### 2.2.1. Aggregate Gradation

Open-grade aggregate is designed primarily as a skeleton of OAC for GSP. It mainly maintains compressive pressures from traffic loads and forms inter-air voids of OAC. Therefore, a good gradation can not only provide a good compressive strength but also construct an even and interconnective space of void.

An optimal gradation mostly can be realized by two design methods: an experiential method and a volumetric method called the main aggregate filling (MAF) method. The latter assumes coarse aggregate is supporting the main structure of OAC and can achieve the required void volume by filling fine aggregate and binder into the compacted coarse aggregate. Due to back-calculation of the void volume, this method is concise to control the final void rate of OMA [43,44]. However, interconnectivity of void is overlooked in this method, which is another important factor for permeability of OAC equally to void volume [45]. Therefore, factors affecting particle sizes and consecutiveness of aggregate gradation are considered to establish the relationships between the given characters above, illustrated in Figure 2, Tables 2–4.

**Table 3.** Effects of gradations on physical properties of OAC [47].

**Figure 2.** Some used gradations with different particle sizes and proportions.

The main coarse aggregate used in GSP normally falls in a range of 4.75 mm–13.2 mm particle size. Several combinations of different tiers in this range can form OMA with the requested void rate. Under the similar void rate, gradation with a large particle size will benefit the whole structure. Moreover, the water permeability of samples which have a large proportion in larger particles (14~10 mm) is significantly higher than that of average samples (Table 2) [46]. A lower residual void rate of GSP will be obtained after grouting [50,51]. Overall, a coarser aggregate gradation can be more conducive to forming a larger void space, which has better interconnection, to make the grouted cement mortar more easily able to penetrate the OMA mat [46,52]. Subsequently, flexural tensile strength and compressive strength of the following GSP are also be improved, which are positively linked with low-temperature performance and high-temperature performance, respectively [49]. On the contrary, some other properties may be slightly harmed with the

increase of particle size. Cantabro loss rises to a certain value, which indicates that the cohesion between aggregate in GSP is weakened.

From a different perspective, Ding et al., proposed that uniform gradation could conduct smoother void spaces and less stress concentration than continuous gradation. He compared two types of gradations: uniform gradation (one tier: 10~5 mm) and continuous gradation (three tiers: 13~10 mm, 10~5 mm, and 5~3 mm) [53]. It was found that GSP with uniform gradation was superior to that with continuous gradation in the properties of Marshall stability, splitting strength, compressive resilience, and low-temperature bending capacity [53].

Saboo et al., evaluated seven types of gradations using the hierarchical ranking strategy considering parameters as void content, permeability, abrasion resistance, and tensile strength of OAC; BSI-4% and Densiphalt12-4.5% were selected as the optimal gradation which had a main tier of 4.75 mm~12.5 mm [47]. In general, gradation with the main particle sizes in the scope of 4.75 mm~12.5 mm is appropriate to improve the performance of GSP to form a better structure (Figure 2).

However, a conflict exists in previous hypotheses that mainly concerns coarse particles. For example, a given continuous gradation may also contain larger-tier aggregate. The reason may be that the influence of gravel morphology and fine aggregate is overlooked which needs further research [52]. In our opinion, the two viewpoints are all correct, though only the applying scope of these principles is different. The definitions for describing gradation types are difficult to quantify, due to the existence of few samples with untested properties, such as residual void rate and water permeability. Large particle size and uniform gradation can both contribute to the formation of evener and larger void space.

### 2.2.2. Asphalt Binder

The types of asphalt binder in OMA include base asphalt, rubber asphalt, SBSmodified asphalt, and high-viscosity modified asphalt. Due to various characteristics at different temperatures, the optimal type and content of asphalt binder as two factors varies in OMA design, as shown in Table 5. Orthogonal experiment uncovers the effects of binder on the relative OMA or GSP sample, as illustrated in Tables 5 and 6.

**Table 5.** Optimal binder and high-temperature performance of OMA.

**Table 6.** Optimal binder and low temperature performance of OMA.

Asphalt binder with high viscosity and low penetration, such as 50-pen straight-run asphalt, SBS asphalt, rubber-asphalt, and high-viscosity asphalt, can increase strength and reduce Cantabro loss for OMA which will result in good high-temperature performance of GSP. However, binder types have little impact on void space structure in OMA, as well as strength and moisture resistance of GSP [46,54].

OMA samples using SBS asphalt as binder can obtain favorable water permeability. Setyawan et al., explained that modified asphalt protected OMA from drainage that caused void-blocking [46]. Therefore, the performance of following GSP samples was also enhanced, especially in low-temperature cracking resistance [46,54,55].

Rubber asphalt also promotes high-temperature performance of GSP. It is considered as a potential candidate due to its excellent cost-effectiveness and environmentfriendliness [48].,High-viscosity asphalt also shows advantages of flexural tensile strain and stiffness modulus at a low temperature, which is associated with better anti-cracking performance [54].

Optimal binder content of OMA can be determined by three main parameters in orthogonal testing: compressive strength, tensile strength, and Cantabro loss. The value of it mostly falls in a range of 3%~5%, due to different test methods. However, it can be argued that a high content rate of binder in OMA can increase toughness of OMA, leading to good anti-cracking performance of GSP [57], because the Marshall stability of OMA continues to rise with the increase of the binder content rate, even up to 9% [58].

Besides these two factors: binder type and content in OMA, bonding force for the asphalt-cement interface is overlooked in steps of influential factors analysis, which can contribute to low-temperature performance of GSP. To achieve a criterion for AC (the flexural tensile strain in small beam bending test should exceed 2000 με at −10 ◦C [59]), the factors for anti-cracking characteristics of GSP should be considered of both OMA and HCM to bridge two materials in whole structure in design.

### 2.2.3. High-Fluidity Cement Mortar

High-fluidity cement mortar (HCM) is a specialized grouting material which can penetrate OMA due to its high fluidity to construct a fibrous-like network structure inside the asphalt mattress. Strength and shrinkage rate are the two main factors of HCM controlling the performance of GSP. Therefore, to investigate their effects, various formulas of raw materials and different additives are compared, as shown in Table 7.

**Table 7.** Formulas of raw materials and performance of HCM.

The water-cement ratio is the critical factor affecting both fluidity and strength of cement mortar. The fluidity of HCM is positively associated with the water-cement ratio and negatively related to the sand-cement ratio [65]. However, an excess water-cement ratio or a short sand-cement ratio will lead to high dry shrinkage [60,61]. In addition, with the increase of water-cement ratio, the bleeding rate of HCM shows an upward trend, which may cause void-blocking and slurry leakage [66,67]. Therefore, Cheng et al., recommends that the water-cement ratio should be less than or equal to 0.55 for HCM, while Saboo et al., advises a ratio scope of 0.4~0.6 [45,66]. Additionally, 10~14 s is a recommended range for HCM fluidity [45,49,67].

The additive type is another factor in determining the performance of HCM. Superplasticizer and fly ash are commonly used to increase fluidity of HCM and minimize residual void after grouting. Some studies indicates that silica fume, mineral powder, and ultrafine sand also have good effects on fluidity [43,64,68,69]. Expansion agent UEA can reduce the shrinkage of HCM and furtherly enhance the anti-cracking ability of following GSP [60–63,70]. However, the strength of cement mortar is inevitably underestimated by fluidity improving in the methods above. Therefore, to achieve the balance design between fluidity and strength, a range of strength value is proposed as 10~30 Mpa, and furtherly narrowed to 20~25 Mpa for HCM [45,49]. Although other factors such as the shrinkage and flexural strength of HCM have important effects on the low-temperature performance of GSP, their values are not considered as parameters by recent studies comprehensively.

Cai et al., used ABAQUS software to simulate the shrinkage and expansion of cement mortar in GSP to estimate the inner stresses [71]. Index values of the deformation were calculated and determined to protect GSP from cracking that is caused by stress concentration of shrinkage, as shown in Table 8.

**Table 8.** Recommended deformation indexes of cement mortar in GSP [71].

### *2.3. Admixtures*

Admixtures for GSP refer to specific materials adding to OMA or HCM to improve workability or strength of following GSP. Resin is the admixture early used in GSP which can be traced back to 1976 in the U.S. Amy [72], and polymers have been used to enhance the anti-cracking performance of GSP from the 1980s on the Huishen highway, China [73]. Over the past 40 years, more and more types of admixtures have been developed into a big family for GSP.

According to different ingredients, GSP admixtures can be roughly classified into five categories: polymers, fibers, interface modifiers, emulsified asphalt, and other new functional admixtures. New functional admixtures have been recently explored for extension in road functions such as emission-reducing, water conservation, and weather cooling. In addition, two blending methods are usually used for admixtures of GSP: firstly, admixture is blended into cement mortar to reduce rigidity and dry shrinkage of HCM, which furtherly improve cracking resistance for following GSP; secondly, it is mixed within asphalt to enhance the strength of OMA, which achieves a high rutting resistance for the whole structure. Admixtures and their effects on GSP are illustrated in Table 9.

**Table 9.** Admixtures and their effects on GSP.

**Table 9.** *Cont.*

1 A = Asphalt; C = Cement; P = Polymer; E = Epoxy; F = Fibre; I = Interface modifier; M = Modifier.


**Celebrating** *Applied Sciences* **Reaches 20,000 Articles Milestone: Invited Papers in "Applied Dentistry and Oral Sciences" Section**

Editor

**Mary Anne Melo**

MDPI • Basel • Beijing • Wuhan • Barcelona • Belgrade • Manchester • Tokyo • Cluj • Tianjin

*Editor* Mary Anne Melo General Dentistry University of Maryland School of Dentistry Baltimore United States

*Editorial Office* MDPI St. Alban-Anlage 66 4052 Basel, Switzerland

This is a reprint of articles from the Special Issue published online in the open access journal *Applied Sciences* (ISSN 2076-3417) (available at: www.mdpi.com/journal/applsci/special issues/ Milestone Invited Papers).

For citation purposes, cite each article independently as indicated on the article page online and as indicated below:

LastName, A.A.; LastName, B.B.; LastName, C.C. Article Title. *Journal Name* **Year**, *Volume Number*, Page Range.

**ISBN 978-3-0365-7727-2 (Hbk) ISBN 978-3-0365-7726-5 (PDF)**

Cover image courtesy of Mary Anne Melo

© 2023 by the authors. Articles in this book are Open Access and distributed under the Creative Commons Attribution (CC BY) license, which allows users to download, copy and build upon published articles, as long as the author and publisher are properly credited, which ensures maximum dissemination and a wider impact of our publications.

The book as a whole is distributed by MDPI under the terms and conditions of the Creative Commons license CC BY-NC-ND.