Dynamic replacement is one of the methods of soil improvement [1], consisting of forming columns of aggregate in low-bearing soil by blowing a pounder of a large mass (usually approximately 10-15 t) from a certain height (up to 15-25 m). During the hammering of the column material, large displacements occur in the ground. Under field conditions, it is possible to form columns with lengths of up to approximately 5 m and diameters of approximately 1.5-2 times greater than the width (diameter) of the pounder used (approximately 1.8 m for a pounder with a diameter of 1 m) [2]. During dynamic replacement, the contractor has no direct control over the length of the column and its diameter. An indirect influence on the geometric parameters of the columns can be obtained by selecting the appropriate impact energy and the amount of material filled in the crater. Therefore, the field trial tests performed at the beginning of the work site are crucial, allowing the clarification of technological parameters. An alternative to trial field tests may be numerical analyses enabling the modelling of column driving with variable technological parameters (pounder blow height, pounder mass and shape, amount of material inserted in the crater and its grain size) and variable soil and water conditions. However, this issue is related to the need to model the aforementioned large displacements and the propagation of vibrations in the ground.
In numerical modelling of this type of problem, mesh models are used [3,4,5], a combination of mesh and non-mesh models [6,7] and meshless models [8].
Work on modelling the hammering process in dynamic replacement technology dates back to the 1990s. The first attempt to model pounder drops with crater charge and the comparison of the numerical results with the results of field tests was carried out by Gunaratne et al. [3]. The calculations were carried out using the finite element method (FEM) in a twodimensional model, with automatic mesh change due to large soil displacements and with the use of constitutive models. CAP (modelling of the work platform and column material) and Modified Cam-Clay (modelling of soft soil behaviour). The calculations were carried out in stages, during which the program assumed changes that occurred in the ground during previous pounder drops. The proposed model estimated the depth of the first crater and the pressure of the pore water with an accuracy of about 90%, overestimating by 50% the values of horizontal ground displacements in the vicinity of the first drop. In subsequent pounder drops, the estimated horizontal displacements ranged from 75 to more than 200% with good convergence of soil uplift (up to 10% difference). Despite the differences obtained, the presented model should be evaluated very well, except that it did not include the prediction of the shape of the columns and the propagation of their material during hammering. Some of the authors continued to create a numerical FEM model [4] that estimates the length of the formed column, the effect of the shape of the pounder on its length, and the soil water pressure in the pores being strengthened based on the results in the research plot on a natural scale. The calculations were carried out in ANSYS Mechanical in three-dimensional analysis using the MO Granular constitutive model for backfill and the Multilinear Isotropic Strain Hardening model for cohesive soils. The model parameters were determined by field tests (CPT) and laboratory tests (triaxial compression tests). Satisfactory convergences of the results were obtained.
Danilewicz and Sikora [6] and Danilewicz [7] modelled the formation of the crater [6] and the dynamic replacement column process [7] using the ANSYS LS-Dyna numerical program combining the following methods: meshless Smoothed Particle Hydrodynamics (SPH) and finite element (FEM). The first one (SPH) allows one to model the behaviour of the area at the pounder drop site, where large displacements occur. The second (FEM) is used to model the area where there are smaller displacements, significantly reducing the calculation time for the entire model. The basis for the calculations is research on a research plot carried out on a natural scale on the occasion of a road investment. The columns were formed in silt and peat, layered with sand. The parameters of the constitutive model (MAT_SOIL_AND_FOAM) adopted in the calculations were determined in the triaxial compression tests and static soundings CPTU and dilatometer soundings DMT. Noteworthy is the fact that a threedimensional model was used (a repeatable quadrant of the full issue), which was also associated with the adoption of certain simplifying assumptions, e.g., omitting the extraction of the pounder and replacing it with a backfill. Unfortunately, the shape obtained in the calculations was not compared with the field results (no excavation was carried out), but a very good convergence of horizontal displacements of the soil next to the column was obtained.
Sołowski et al. [8] proposed a numerical model using the method of material points in their own program, calibrating it on the basis of laboratory tests in which the process of forming columns was observed. Soft soil in the form of sawdust was simulated with the Modified Cam-Clay model and the aggregate of the columns with the elastic-perfectly plastic model with the Mohr-Coulomb failure criterion, respectively. The results of calculations in the form of aggregate penetration of columns, horizontal and vertical displacements of the improved medium, were compared with laboratory results. A good convergence of the results was obtained, although difficulties in modelling related to the use of more complex numerical models, the need to improve modelling with the material point method and the difficulties of numerical modelling in the presence of a larger number of pounder discharges and aggregate backfills were emphasised. However, this type of approach requires, in the first place, the preparation of dedicated software.
In this paper, the approach of combining meshed (FEM) and meshless (SPH) methods using the ANSYS LS-Dyna program with model calibration based on authors’ laboratory test results was used [9]. First, an attempt was made to model the propagation process of the hammered aggregate into the soil.
The laboratory test stand allowed one to form a column in dynamic replacement technology by freely dropping the pounder into the soft soil layer lined with the bearing layer.
Research was carried out in a box with dimensions of 1x1x1 m, in which the walls were made of acrylic glass with a thickness of 20 mm, circumferentially improved and in height with steel sections (Fig. 1a). To dampen vibrations generated during the formation of the column and transmitted through the ground, 50 mm thick mineral wool (Fig. 1a), protected against moisture (Fig. 1b), was placed on the walls of the chamber. First, compacted medium sand (relative density ID = 88%) with a thickness of 400 mm was formed in layers in the box, and then a separation geotextile was laid on it and a soft soil was formed from silty clay with a thickness of 400 mm (Fig. 1b). The clay used was taken from the excavation of one of the road investments. After being brought to the laboratory, the clay was dried to a water content of about 1% and very carefully crushed. In this form, it was placed in a box in layers with a thickness of about 20 mm, evenly poured with enough water to achieve the assumed liquidity index (IL = 0.7), consistency index (IC = 0.3).
Test bench: a) damping of the side walls with mineral wool, b) view of the clay layer with foil protection
Columns were formed from basalt aggregate with a fraction of 5-25 mm.
A total of two laboratory tests were carried out, during which the following were performed:
- –
measurements of pounder penetrations and masses of aggregate used (test Nos. 1 and 2),
- –
measurements of vertical and horizontal vibration accelerations of the pounder and the soil (tests Nos. 1 and 2),
- –
measurements of horizontal displacements of the ground taking place next to the column (test No. 1),
- –
measurements of vertical displacements of the soil (uplift, test No. 1),
- –
CPT (cone penetration test) probing before and after amplification as a function of time and distance from the column (test 1 and partially test 2), – tests of changes in the water content of the improved soil as a function of time and distance from the column (test No. 1 and partially test No. 2)
- –
consolidation and settlement of the column-soft subsoil system (test No. 2),
- –
open-pit mines and column inventory (study Nos. 1 and 2).
Additionally, the physical and mechanical parameters of the soil used were determined. They were triaxial compression tests (in conditions without consolidation and without drainage (UU)) and oedometer tests of soft plastic clays, as well as direct shear box tests of basalt aggregate used for columns and the bottom bearing layer.
In both studies, dynamic replacement columns were formed in the same way, by dropping a barrel-shaped pounder with a base diameter of 90 mm, a maximum diameter of 100 mm, a height of 200 mm and a weight of 10.55 kg, which forced crushed basalt aggregate with a fraction of 5-25 mm into soft plastic clay. The blows were divided into three series. The first one, during which the crater was formed and the first aggregate backfills were applied, consisted of discharges from the height of 0.3, 0.5, 0.7 and 3x0.85 m (discharges Nos. 1-6). The second (blows Nos. 7-32) was the main series, in which the pounder fell from a height of 0.85 m. The last series is to finish off the column head with blows from the height of 0.85, 0.5, 0.3 and 0.15 m. The crater backfills followed discharges No. 3, 4, 5, 6, 33, 34, 35, 36 and after every two blows in the main series. In total, the pounder was dropped 36 times. After individual pounder discharges, the depth of the resulting crater was measured (Fig. 2), and the masses and heights of aggregate backfills were recorded.
Crater depths after individual pounder blows
In total, 22.80 kg of aggregate was used to form the stone column in test No. 1, while in test No. 2, 16.65 kg was used. The columns were 350–390 mm long and had diameters of 150–165 mm (heads) and 220–280 mm (the largest diameter) – Fig. 3. Despite repeatable ground conditions and the same pounder discharge heights, the formed columns differ significantly, which is quite characteristic of this method of improvement.
Inventoried shape of the columns: a) No. 1, b) No. 2
In this paper, a two-dimensional (2D) model in an axisymmetric state of deformation is proposed using the LS-Dyna program for calculations (Fig. 4) [10]. This type of approach significantly speeds up the calculation time and allows for the initial calibration of the model, which can then be prepared three-dimensionally after a positive evaluation. Due to the fact that the tests are carried out in a cuboid box, the radius of the 2D model (approx. 540 mm) was selected so that the horizontal cross-sectional area of the cylindrical cell was equal to the horizontal cross-sectional area of the laboratory box. The top layer, in which the column was formed, was divided into two parts: in the place of the examined largest displacements, the substrate is modelled by SPH elements (with a total width of approximately 220 mm) and in the remaining part by FEM elements (Fig. 4). At the junction of these two zones, contact elements (CONSTRAINED_TIE-BREAK) were used [10]. The load-bearing layer and the pounder are modelled with FEM elements (Fig. 4), while the aggregates are modelled with SPH elements.
Numerical model
The SPH (Smoothed Particle Hydrodynamics) method belongs to the family of gridless methods, with the Lagrange approach, in which the area is divided into particles, representing parts of the body, remaining in the vicinity during the calculations [7]. For each particle, the so-called nucleus function is defined, with a given distribution decreasing as it moves away from the particle's axis, allowing for the interpolation of physical quantities over adjacent particles. The introduction of this function allows us to present a partial differential equation in the form of a system of ordinary differential equations in a discrete form, with a time variable, due to which such a system can be solved using one of many numerical integration schemes [7]. The main advantages of the method include the ability to describe large deformations and fragmentations [7].
The thickness of the layers in the model is consistent with that used in laboratory research.
Standard boundary conditions have been adopted in the form of blocking vertical and horizontal movements on the bottom edge and horizontal movement on the side edges of the model. In addition, due to the damping of vibrations on the walls of the site (use of mineral wool), an edge condition was applied on its numerical side, which did not allow for the reflection of waves reaching it (BOUNDARY_NON_REFLECTING_2D).
For the modelling of clay, sand and grit layers, an elastic-perfectly plastic model with the Mohr-Coulomb failure criterion (MAT_173) was used, and its parameters such as the angle of internal friction (ϕ), cohesion (c), Poisson’s ratio (ν), the shear modulus (G) were determined in laboratory tests (ϕ, c, G) and based on the literature (ν). The model assumes conditions without drainage, which is characteristic of the issues of very fast loading of cohesive soil [11]. The pounder was modelled using the RIGID model (MAT_20), assuming the parameters of linear elasticity (modulus of elasticity (E) and Poisson’s ratio (ν)) as for steel.
The calculations were carried out in stages, using the so-called small restart, which allows for the continuity of the model's calculations. The analysis modelled pounder discharges, its extraction, insertion of the charge, and then its free discharge into the crater, assuming the height of the charge as in laboratory tests. It is planned that 36-pounder discharges and 20 portions of aggregates will be divided into 110 stages.
Contact elements (2D_NODE_TO_SOLID) were used between the pounder and the SPH elements, for which the static coefficient of friction (FS) and the dynamic coefficient of friction (FD) were defined with values adopted in this type of issue [12, 13].
All analyses were performed using the explicit method, i.e., methods based on the explicit scheme are used to integrate the equations of motion, searching for a solution for the next step using the current step, and the state of primary stress was achieved by applying the dynamic relaxation phase.
The model uses 4140 FEM elements and 26837 SPH elements. The dimensions of the FEM elements were 10x10 mm, and the spacing of the SPH elements was 2.5 mm.
The values of the parameters of the models used are presented in Table 1.
Adopted constitutive models and their parameters
| Lp. | Kind | Constitutive model adopted | Model parameters |
|---|---|---|---|
| 1. | Pounder | 020-Rigid | E=210 GPa, ν=0.3, γ=74 kN/m3, |
| 2. | Base layer | 173-Mohr-Coulomb | G=41.7 MPa, ϕ=33, c=0 kPa, ν=0.2, γ=18 kN/m3 |
| 3. | Improved layer | 173-Mohr-Coulomb (Tresca criterion) | G=1.39 MPa, Su=10.2 kPa, ν=0.31, γ=18.7 kN/m3 |
| 4. | DR column aggregate | 173-Mohr-Coulomb | G=6.86 MPa, ϕ=47, c=0 kPa, ν=0.2, γ=18 kN/m3 |
| 5. | Pounder contact elements | Nd | FS=0.5 |
The values of the primary stresses in the soil before the start of dynamic replacement are shown in Figure 5. Assuming the weights of the soil layers (Table 1) and their thickness (2x0.4 m), the values presented are correct.
Primary stress values before modelling dynamic replacement
Fig. 6 presents selected stages related to the first crater backfill in the form of a view: model after the 3rd pounder drop (complete pounder plunge, Fig. 6a), free discharge of the aggregate (Fig. 6b) and after falling into the crater (Fig. 6c) and after the pounder discharge into the aggregate (Fig. 6d).
View of the model in the characteristic stages associated with the first fill: a) after the 3rd pounder drop, b) free-falling aggregate into the crater, c) backfill after pouring into the crater, d) after pounder drop into the aggregate
Fig. 7 shows the view of the model after driving a different number of aggregates (4-11). The proposed model simulated 22 out of 36 pounder discharges (11 out of 20 charges), and in the next stage, there was a problem with the convergence of calculations. However, qualitative results were obtained following laboratory tests conducted by one of the authors [14]. It is visible that the first two backfills are plunged practically vertically into the ground and the next ones are compacted and spread to the sides. The reason for this is the increasing rigidity of the crater floor. With each subsequent drop, the thrust of the column base decreases (Fig. 8), which causes the column to widen edgeways. At the interface between the SPH elements and the FEM elements, there is a discontinuity of the nuclear function (Figs. 7c-d), which was not found in laboratory tests. This may be due to a zone that is too small, modelled with SPH elements or a numerical formulation of SPH particles.
View of the model after the pounder blow in: a) 4 charges, b) 6 charges, c) 8 charges, d) 11 charges
Diagram of vertical displacements of the SPH element occurring directly under the base of the column (the first discharge into the aggregate after 2 s, preceded by 3 discharges forming the crater)
After the 11th charge, the diameter of the column head was approximately 115 mm, the maximum diameter was approximately 178 mm and the diameter of the base was 98 mm, respectively (Fig. 9). They are smaller than those obtained in laboratory tests, but numerical modelling of the column-drive process has not been completed. After about 60% of the calculations, the column obtained in the calculations was characterised by a much smaller diameter in its lower parts, with similar dimensions of the upper part (Fig. 9). The maximum diameter was higher than in laboratory tests.
Shapes of the obtained columns: a) numerical column, b) comparison of shapes: numerical and laboratory
As a result of the calculations, the uplifts of the soil at the column were obtained equal to approximately 45 mm (Fig. 10) and are much larger than those obtained in the laboratory (10-14 mm). The numerical model did not fully take into account the preservation of soil with a rather unusual structure. During soil sampling for laboratory tests, it was noticed that the soft-plastic cohesive soil is characterised by high porosity, visible to the naked eye. The reason for this is the method of formation in the laboratory (pouring layers of shredded soil and pouring water). Therefore, during the hammering of the column, the soil was compacted, which made its uplifts smaller. This was not taken into account by the adopted constitutive model of the soil.
Diagram of vertical displacements of the SPH element occurring in the upper part of the model directly next to the dropped pounder
Fig. 11 shows a comparison of the depths formed after the discharge of craters obtained in calculations and laboratory tests. Up to the 5th discharge, a satisfactory convergence of results was achieved (in the laboratory interval or up to 10% higher). In the next stages in the proposed model, the crater depths are up to 70% lower than those observed in the laboratory. This may be determined by the choice of a too simple constitutive model of the ground. While in the first stage the plunging of the pounder and aggregate is mainly determined by the failure of the soil, well simulated by the elastic-perfectly plastic model with the Mohr-Coulomb failure criterion, in the further part of the modelling of soil compaction is also important.
Comparison of the depth of the crater obtained in laboratory tests and SPH+FEM calculations
Work on the proposed model allowed us to recognise and determine the method of modelling the column insertion by the dynamic replacement method in the LS-Dyna numerical program from the technical side. This allowed for full modelling of pounder blows, its free discharge and extraction, and modelling of the aggregate. The authors are not aware of works and calculations that present the hammering process in such detail.
Modelling of dynamic replacement using an elastic-perfectly plastic model with the Mohr-Coulomb failure criterion for the improved layer allowed for a good reflection of the initial stage of column driving, associated mainly with soil failure and greater discrepancies with laboratory results for some of the remaining stages related to soil compaction. It was also not possible to model the full driving process due to the lack of convergence of the calculations.
Undoubtedly, the presented model requires changes, which in the first place may concern the modelling of the upper improved layer only with SPH elements (due to the small width of the site, this is possible and should allow one to avoid the formation of the observed discontinuity and increase the stability of calculations) and to check the modelling of density changes in the scope of the constitutive model used and settings related to SPH elements. In the event of a further negative assessment, the next step should be the use of constitutive models modelling failure and volume changes in the improved soil and the transition to 3D modelling using symmetry (analysis of a quarter of the defined problem).