Simulation of seismic response of passively stabilised sand Open Access

A_d

a positive variable entering the dilatancy function of the NTUA-Sand model

A_o

positive scalar multiplier of the dilatancy function of the NTUA-Sand model

a_max

peak value of acceleration: g

CS

percentage per weight of silica particles in the colloidal silica aqueous dispersion

CSR

cyclic stress ratio $τcyc/σvo′$

D

dilatancy function of the NTUA-Sand model

e

void ratio

e_cs

critical value of the void ratio

G

shear modulus: kPa

G_max

maximum (elastic) shear modulus: kPa

h_bd

positive variable entering the plastic modulus of the NTUA-Sand model

h_o

positive scalar multiplier of the plastic modulus of the NTUA-Sand model

K

bulk modulus of the pore fluid: kPa

K_p

plastic modulus of the NTUA-Sand model: kPa

K_w

bulk modulus of water: kPa

M

critical-state deviatoric stress ratio

M^b

peak deviatoric stress ratio

M^d

dilatancy deviatoric stress ratio

N_L

number of cycles for liquefaction

n

reduction denominator of the pore fluid modulus

P_f

pore fluid pressure: kPa

p′

mean effective stress: kPa

q

(triaxial) deviatoric stress: kPa

r_u

excess pore pressure ratio (⁠$= Δu/σvo′$⁠)

t

time: s

γ

single-amplitude cyclic shear strain

Δu

excess pore pressure: kPa

$ϵ˙qp$

plastic deviatoric strain rate

$ϵ˙volp$

plastic volumetric strain rate

ϵ_vol

volumetric strain of the soil element

ζ

fluid volume increase per soil element volume

η

deviatoric stress ratio

ξ

hysteretic damping ratio

$σvo′$

initial vertical effective stress: kPa

τ_cyc

single-amplitude cyclic shear stress: kPa

ψ

state parameter

For developed sites, the large majority of conventional liquefaction mitigation techniques cannot be applied, due to vibration (e.g. vibro replacement) or accessibility (e.g. deep soil mixing) restrictions. The few existing alternatives also have their own limitations (e.g. prefabricated drains cannot be easily installed under the whole structure). Passive site stabilisation (Gallagher, 2000) is a novel (and very promising) liquefaction mitigation technique for developed sites. It comprises the low-pressure injection of a stabiliser (colloidal silica grout) in the soil pores, which transforms into a firm gel after a predetermined (chemically controlled) time.

The concept is portrayed in Figure 1, where injection and pumping wells can be used around the perimeter of a structure at risk, making sure that the grout has replaced the water in the pores of the liquefiable soil under the structure. The novelty of the method is its requirement of low-pressure injection, since colloidal silica (an aqueous dispersion of silica particles) initially has a viscosity slightly higher than that of water (<1·5 cP), which also enables its injection for long distances (tens of metres; see Figure 1). After a predetermined time (named ‘gel time’, which may range from a few hours to several days, depending on the area to be treated), the viscosity increases rapidly, prohibiting further injection, and the colloidal silica gels in the pores of the non-cohesive soil. This change alters the mechanical response of the (stabilised) sand, making it less vulnerable to plastic strain accumulation related to liquefaction or cyclic mobility during earthquakes (e.g. see data from Diaz-Rodriguez et al. (2008), Gallagher and Mitchell (2002)).

This ground improvement technique is currently at the experimental stage. One important issue that needs to be addressed before it becomes practically usable is how to alter the colloidal silica properties to control the gel time, thus allowing effective injection at the desired distance within a soil layer with known properties. Hence, for the quality assurance of this ground improvement method, important issues are the control of gel time (e.g. Agapoulaki and Papadimitriou, 2015; Gallagher, 2000) and the injection potential for different soil types (e.g. Agapoulaki et al., 2015; Gallagher and Lin, 2009). However, both of these issues are beyond the scope of this paper.

Another significant issue is to establish the microscopic mechanism of improvement, as the first step for the formulation of reliable constitutive models for this new geomaterial named ‘stabilised sand’. This has not yet been established, and, thus, any attempt to formulate a numerical simulation approach will be necessarily phenomenological, on the basis of what is measured in element and physical modelling tests on such materials. Such an approach is performed herein, with the aid of an existing state-of-the-art constitutive model for the liquefaction response of sands (National Technical University of Athens-Sand model (NTUA-Sand); Andrianopoulos et al., 2010a, 2010b), but with an eye on the literature proposals for the microscopic mechanism of improvement. These simulation efforts, if successful, will enable the use of numerical techniques for analysing the response of stabilised geostructures, thus paving the way for the use of this liquefaction mitigation technique in practice.

The next section presents a brief overview of the mechanical response of stabilised sands based on element test results. Then, comparison is made of data from element test results on both untreated and stabilised sand samples with data from simulations performed with different phenomenological approaches. Finally, the comparisons are extended to the system level, by the use of dynamic centrifuge test data.

To date, the mechanical response of stabilised sands has only been studied macroscopically, by (a) limited on-site testing (Gallagher et al., 2007a), (b) few physical modelling attempts (dynamic centrifuge tests: Conlee et al. (2012), Gallagher et al. (2007b), Pamuk et al. (2007)) and (c) a number of laboratory efforts outlined later in this paper. However, each of these laboratory efforts has shed light on partial aspects of the response, since they have been performed on different sands by employing different types of tests and sometimes with different types of colloidal silica (although the majority have been performed with Ludox SM). Given these variations, it is very difficult to establish a clear framework for the mechanical response of stabilised sand.

However, there are some indisputable traits of stabilised sand response. For example, the stabilisation offers unconfined compression strength in stabilised sands (e.g. Persoff et al., 1999), the value of which increases with CS(%) – that is, the percentage per weight of colloidal silica particles in the grout (e.g. Gallagher, 2000; Papadimitriou and Agapoulaki, 2013). Focusing on the dynamic response at small cyclic shear strains γ (e.g. smaller than 10⁻⁵), published data show that stabilisation with colloidal silica leads to a small increase in the elastic shear modulus G_max, equal to 10–25% on average, for relatively low confining pressures of up to 100 kPa (e.g. Conlee, 2010; Papadimitriou and Agapoulaki, 2013; Spencer et al., 2008). In contrast, at medium cyclic shear strains γ (e.g. up to 5 × 10⁻⁴), published data (Spencer et al., 2008) show no substantial effect on the normalised shear modulus G/G_max degradation and hysteretic damping ξ increase curves with cyclic shear strain γ.

At large cyclic shear strains, element tests reveal a much more stable behaviour of stabilised sands, in comparison with that of their untreated counterparts. For example Figure 2 presents the evolution of the excess pore pressure ratio r_u build-up and the corresponding (double-amplitude (DA)) cyclic shear strain accumulation with cycles for untreated (CS = 0%) and stabilised sand (CS = 7·25%, 14·5%), on the basis of cyclic simple shear tests with the same initial and loading conditions (data from Diaz-Rodriguez et al. (2008)). Observe that after stabilisation with colloidal silica, the same level of r_u is reached after a much increased number of cycles, and this beneficial effect increases with increasing CS(%) per weight in the grout. From a practical point of view, it is important to note in Figure 2 that the soil skeleton–pore fluid response of the stabilised samples show increased cyclic shear resistance, since they reach the same level of DA cyclic shear strains after a much increased number of cycles, in comparison with that of the untreated sand.

This increase in cyclic shear resistance has been observed repeatedly in different types of cyclic tests on stabilised sands (e.g. triaxial in the paper by Gallagher and Mitchell (2002), torsional shear in the article by Kodaka et al. (2005), simple shear in the paper by Diaz-Rodriguez et al. (2008)) and can be quantified better in terms of liquefaction resistance curves, where liquefaction is depicted when the DA cyclic shear strain exceeds a set level (e.g. of 1, 2 or 5%). For example, Figure 3 shows the significant increase in liquefaction resistance quantified in terms of higher cyclic stress ratios (⁠$CSR = τcyc/σvo′$⁠) for a given number of cycles to liquefaction N_L due to stabilisation, based on cyclic simple shear tests performed on stabilised and untreated sand samples initiating from the same initial conditions of stress and relative density D_r = 40% (data from Diaz-Rodriguez et al. (2008); also including the data for CS = 0 and 14·5% from Figure 2). In particular, these data show the clear beneficial effect of stabilisation for sand samples, irrespective of the stress level. The same paper also shows similar data for D_r = 60%, leading to the overall conclusion that quantitatively the beneficial effect of stabilisation is more pronounced when it is most needed – that is, for the lower-density sands shown in Figure 3.

Furthermore, observe in Figure 3 that the increase in liquefaction resistance in terms of CSR for any given number of cycles N_L is typically of the order of 0·08–0·13, and this is more or less independent of the chosen N_L value for each case – that is, the data show an approximately parallel translation of the so-called liquefaction resistance curve CSR–N_L to higher CSR values due to stabilisation. A survey of related results from the literature shows that this is a common trend for all data, depicting that the stabilisation leads to an increase in CSR ranging from as low as 0·01 (for dense sands under low confinement; Diaz-Rodriguez et al., 2008) up to a maximum of 0·2 (for loose sands; Kodaka et al., 2005). It should be mentioned here that the amount of CSR increase due to stabilisation described earlier is also a function of the type of colloidal silica (not all tests have been performed with Ludox SM) and the percentage of CS(%) in the grout (Gallagher and Mitchell, 2002). In practice, for the mentioned type of stabiliser, the tests show that CS = 5% leads to imperceptible benefits, while a value of CS = 10% provides significant increase in CSR and poses as an upper limit for a safe and economically viable ground improvement (e.g. Gallagher and Mitchell, 2002). Finally, it is noted here that for DA cyclic shear strain values smaller than 5%, the benefits of treatment in terms of CSR continue to appear, but are relatively less pronounced.

These beneficial effects at the element level are also demonstrated at the system level, at least in the few pertinent physical modelling attempts in the literature, namely horizontal stabilised sand layers do not exhibit the liquefaction-related deamplification of seismic ground motion (e.g. Gallagher et al., 2007b), while their settlements are also significantly reduced (e.g. Gallagher et al., 2007a, 2007b). Similarly, stabilisation essentially averts horizontal displacements (Conlee et al., 2012; Pamuk et al., 2007) and drastically reduces bending moments of piles in laterally spreading sand layers (Pamuk et al., 2007).

At this time, there is no dedicated constitutive model for the simulation of this new geomaterial (stabilised sand) other than the preliminary proposal of Kodaka et al. (2005), who slightly differentiated a pre-existing model for sands (Oka et al., 1999) to account for less intense excess pore pressure development after stabilisation – for example, after load reversal from above the phase transformation line (Ishihara et al., 1975). It is believed that proposing a dedicated constitutive model for stabilised sands requires first understanding the mechanism of improvement from a micromechanical point of view and then establishing a clear framework of mechanical response from a macromechanical point of view. An outline of the latter has been briefly introduced in the previous section. However, the former is yet to be established, since there are two main approaches in the related literature for the increased deformation resistance after stabilisation. Specifically, this increased deformation resistance is considered to originate (a) from bonding between the gel and the grains, as well as gel encapsulation of individual grains (Gallagher et al., 2007a), and (b) from the reduced excess pore pressure development during shearing, due to the seemingly compressible gel in the soil pores, in comparison with the incompressible water of saturated deposits (Towhata, 2008).

In acknowledgement of the premature nature of constructing a dedicated constitutive model, this paper explores both these approaches from a phenomenological point of view, by way of the use of an existing state-of-the-art constitutive model for the liquefaction response of sands (NTUA-Sand; Andrianopoulos et al., 2010a, 2010b), without changes in its constitutive equations. For better comprehension of the simulation process, the employed model is briefly described in the next section.

The simulations are performed with the finite-difference method, by employing Flac (Itasca Inc., 2005), a two-dimensional explicit finite-difference program, where NTUA-Sand has been implemented by employing the user-defined-model capability, thus creating a subroutine that is readily available for potential users (Itasca Consultants GmbH, 2016).

NTUA-Sand is a bounding surface critical-state plasticity model, which does not consider a purely elastic region, but retains the elasto-plastic strain rate deconvolution of classical elasto-plasticity. The elastic strain rate is estimated on the basis of a Ramberg–Osgood-type non-linear hysteretic formulation, which governs the response for small and medium cyclic shear strains. Large cyclic shear strain response is mainly a function of the plastic strain rate of the model, which is explained with the aid of Figure 4. Observe that the model assumes three conical surfaces in the triaxial stress space (p′–q), where p′ is the mean effective stress and q is the stress deviator (difference between the major and minor principal stresses), while their ratio q/p′ depicts the deviatoric stress ratio η. These are the critical-state, the bounding and the dilatancy surface, which are defined fully by their apertures in terms of the homonymous deviatoric stress ratios M, M ^b and M ^d, respectively. The model calculates the deviatoric and volumetric plastic strain rates, $ϵ˙qp$ and $ϵ˙volp$⁠, respectively, in terms of the difference of the ever-current η from these ratios (which are Lode angle dependent). In particular, the strain rates are assumed functions of a plastic (hardening) modulus K _p and a dilatancy function D (implying non-associative flow rule), in the following form

Note that the model predicts plastic strains that are scaled by the value of K _p and are non-zero only for shear paths where the η changes in value (Equation 1a). Note also that the magnitude of the volumetric strain rate is governed by the value of D, while its sign is solely a function of the sign of D (Equation 1b). These key model ingredients are calculated by

where the terms in parentheses quantify the distances of η from the M ^b and M ^d stress ratios and h _bf and A _d are positive model variables (whose expressions are of no importance here), while h _o and A _o are the two positive model constants that scale K _p and D, respectively. Observe that Equation 2a, in conjunction with Equation 1a, yields a hardening response (K _p > 0) for η < Μ ^b and a softening response (K _p < 0) for η > Μ ^b. Also, note that D takes positive and negative values in Equation 2b for η < Μ ^d and η > Μ ^d, respectively, thus depicting volume reduction and increase by way of Equation 1b.

The critical-state M, the bounding M ^b and the dilatancy M ^d deviatoric stress ratios are interrelated on the basis of the state parameter ψ (Been and Jefferies, 1985), a concept first proposed by the two-surface model of Manzari and Dafalias (1997). Since then, this concept has been adopted by many constitutive models, which may be collectively referred to as SANISAND models (simple anisotropic sand model, proposed by Taiebat and Dafalias (2008)). While the critical-state M is constant, the M ^b and M ^d stress ratios change with the state parameter ψ (which measures the distance from the unique critical-state line in the void ratio e–p′ space) and coincide with M for ψ = 0 – that is, when the void ratio e takes its critical value e _cs. Hence, at the critical state, when e = e _cs and η = Μ, the model predicts K _p = D = 0, and, hence, plastic deviatoric strains develop at a constant volume, due to Equations 1a and 1b, as required by the classical critical-state theory.

NTUA-Sand has a total of 13 model constants, which are outlined in Table 1, along with their values for the Nevada sand, after calibration on the data of Arulmoli et al. (1992). Figure 5 shows an example of the model accuracy for a cyclic undrained (direct) simple shear test on the Nevada sand with a relative density D _r = 60%. Observe how the model predicts well both the (shear stress) τ–(effective vertical stress) $σv′$ path and the (shear stress) τ–(shear strain) γ response. Of more interest is not the accuracy for a single element test, but how well the model predicts the liquefaction resistance of the Nevada sand in (direct) simple shear conditions (which resemble the seismic shaking response). Hence, Figure 6 compares the liquefaction resistance of the Nevada sand based on various (direct) simple shear element tests (hollow symbols) under different relative densities, D _r = 40 and 60% (Figures 6(a) and 6(b) respectively) with the pertinent simulations with NTUA-Sand (solid symbols) and the well-established liquefaction resistance curves of DeAlba et al. (1976). An overall satisfactory agreement is observed, which underlines the appropriateness of the NTUA-Sand model for liquefaction analyses. This is the reason that this model has been systematically used for analysis of various geotechnical systems in a liquefaction regime, such as footings (e.g. Karamitros et al., 2012), piles (e.g. Chaloulos et al., 2013) and lateral spreading (e.g. Valsamis et al., 2010). Further details on the model specifics, including the exact definition of all model constants in Table 1 and not only the h _o and A _o that enter Equation 2, may be found in the papers by Andrianopoulos et al. (2010a, 2010b).

NTUA-Sand is used for the simulation of both the stabilised sand and the untreated sand, without change in its constitutive equations. Therefore, to simulate the effects of gel bonding of the grains and/or the gel encapsulation of grains, one has to resort to targeted recalibration of its model constants. As outlined earlier, the elastic strain rate governs the response under small and medium cyclic shear strains, whereas the response under large cyclic shear strains is primarily a function of the plastic strain rate (and hence of K _p and D). This deconvolution of the simulated response under different strain levels makes the selection of NTUA-Sand very appealing, since the stabilisation affects mainly the response under large cyclic shear strains (see the section headed ‘Mechanical response of stabilised sands’), and, therefore, the target of recalibration is set on K _p and D. In any case, the reference sand here is Nevada, and, hence, the (reference) values for the plastic and dilatancy constants are h _o = 15 000 and A _o = 0·8, respectively, for the untreated sand (see Table 1).

In particular, based on the section headed ‘Mechanical response of stabilised sands’, stabilisation produces a stiffer and more dilative material, whose liquefaction resistance is significantly increased. Hence, the required recalibration of K _p and/or D for simulating stabilisation effects entails the increase in h _o and/or the decrease in A _o relative to their values for the untreated sand. The accuracy of this simulation approach is evaluated in the next section.

In the previous section, NTUA-Sand is recalibrated by considering that the stabilisation affects the sand skeleton response. The alternative considers that the sand skeleton remains unaffected (i.e. NTUA-Sand retains the model constant values for the untreated sand; see Table 1), but it is the pore fluid of the saturated medium that is altered. In particular, Towhata (2008) showed results of unconfined compression tests on gelled samples of pure colloidal silica, which exhibited minimal strength (on the order of a few kilopascals) and volume reduction, while it is well known that water is practically incompressible. In other words, while any tendency for volume reduction (e.g. due to seismic shaking) translates to excess pore pressure development in untreated sand, in stabilised sand any similar loading is not expected to generate (significant) excess pore pressures. These characteristics of the gelled pore fluid may enter the fully coupled formulation by way of the equations of Biot’s consolidation theory (for incompressible grains)

where K is the bulk modulus of the pore fluid, ζ is the fluid volume increase per soil element volume (due to seepage or pumping), ϵ _vol is the volumetric strain of the soil element and P _f is the pressure of the pore fluid. For undrained conditions, which essentially prevail during seismic shaking, ζ = 0, and, hence, the pore fluid pressure P _f is analogous to the volumetric strain rate of the soil skeleton and the pore fluid bulk modulus K. Now, if the gelled colloidal silica is considered more compressible than water (whose bulk modulus is denoted by K _w) by n times, then

The use of Equation 4 with n > 1 in conjunction with Equation 3 leads to smaller P _f values for the same volumetric strain rate of the soil skeleton. This translates to smaller excess pore pressures under the same loading – that is, relatively higher strength and stiffness, which are both effective stress level dependent in NTUA-Sand, in comparison with what is expected for the untreated sand. Again, the accuracy of this simulation approach is evaluated in the next section.

Numerical simulations are performed to compare the accuracy of the two aforementioned simulation approaches against laboratory tests on stabilised sands. The emphasis of the simulation effort is on quantifying the effect of stabilisation on the various aspects of response, and, hence, simulations are also performed for the untreated sand under the same loading conditions.

Figure 7 compares the effect of stabilisation on the small strain (elastic) modulus G _max values from data and simulations. The data originate from the resonant column tests of Spencer et al. (2008), which show a small increase of 10% for p′ = 50 kPa. The effect of stabilisation is quantified for a wider range of p′ = 30–90 kPa approximately, by employing three simulation efforts: (a) an increase in the plastic modulus constant h _o by ten times, (b) a decrease in the dilatancy constant A _o by ten times and (c) a reduction in the fluid bulk K modulus by n = 50 times in comparison with that of water. The simulations show that the last two approaches produce imperceptible effects on the value of G _max due to stabilisation, whereas the (large) tenfold increase in h _o leads to an increase of 15–22% – that is, higher than the 10% of the shown data and lower than 25%, which is the average increase reported by Conlee (2010) on the basis of in-flight measurements of shear wave velocity in centrifuge tests. Hence, overall, it may be concluded that all three simulation approaches for stabilisation lead to relatively small effects on the G _max value, in general accordance with the test data for the small cyclic shear strain regime.

Figure 8 focuses on the normalised shear modulus G/G _max degradation and hysteretic damping ratio ξ increase curves with cyclic shear strain level γ, on the basis of tests and simulations. The data originate from the resonant column tests of Spencer et al. (2008) and show in Figures 8(a) and 8(c) that a stabilisation with CS = 9% (a relatively high value) does not alter the measurements, in comparison with the untreated sand (remember that G _max is higher in the stabilised sand by only 10%; see Figure 7). The pertinent numerical results are presented in Figures 8(b) and 8(d) and include the simulation for untreated sand and the three simulation efforts for stabilised sand, 10h _o, A _o/10 and K _w/50, whose G _max values are also presented in Figure 7. Observe that, similarly to the test data, the simulations for stabilisation do not show any significant difference from that of untreated sand. Specifically, the simulations for A _o/10 and K _w/50 show imperceptible differences from the untreated response, while only the (large) tenfold increase in h _o leads to slightly reduced non-linearity, particularly in terms of the damping ratio. Thus, overall, it may be concluded that all three simulation approaches for stabilisation lead to relatively very small (if existing) effects on the G/G _max–γ and ξ–γ curves, in general accordance with the test data for the medium cyclic shear strain regime.

Figure 9 focuses on the large cyclic shear strain regime and specifically on the stabilisation-induced increase in cyclic stress ratios (CSR) required for liquefaction for any given number of cycles (e.g. see Figure 3). It emphasises the values of load cycles usually related to earthquake loading (N _L = 5–15) and quantifies the difference in CSR based on the tests of Diaz-Rodriguez et al. (2008), which are presented in Figure 3. Figure 9 also includes numerical predictions for the gain in CSR due to stabilisation following the three simulation efforts for stabilised sand, 4h _o, A _o/2·67 and K _w/50. Observe that the data present a gain in CSR of the order of 0·1 (from 0·08 to 0·125) with a slight reduction in the number of N _L and the three simulation efforts compare satisfactorily with the data, both quantitatively and qualitatively. Of interest here is that this good agreement in terms of the gain in CSR is obtained with less intense increase in h _o (four times, against ten times in Figures 7 and 8), as well as less intense decrease in A _o (2·67 times, against 10 times in Figures 7 and 8). Hence, if one repeated the pertinent simulations in Figures 7 and 8 with these reduced changes in the values of h _o and A _o, the thereby simulated small effects of stabilisation would be reduced even further. This does not hold for the pore fluid modulus approach, whose value of K _w/50 is used in all simulations, showing good simulative accuracy for all cyclic shear strain regimes.

The previous section presents two approaches for simulating passive stabilisation, one by focusing on its (possible) effects on the sand skeleton (by way of recalibration of the NTUA-Sand model in terms of plastic modulus and dilatancy) and another by accounting for changes in the pore fluid modulus (by way of reduction of the bulk modulus of the pore fluid in the coupled analysis). It also proceeds in their calibration on the basis of element test data. The current section explores the potential of both approaches in a boundary value problem involving stabilised sand. For this purpose, a dynamic centrifuge test is selected (Gallagher et al., 2007b) that refers to the one-dimensional (1D) dynamic response of a uniform sand layer under sinusoidal motion. The test essentially replicates the 1D test performed during the Verification of Liquefaction Analyses by Centrifuge Studies (VELACS) project (model 1; Taboada, 1995) by using loose saturated Nevada sand (D _r = 40% here, against 45% in VELACS), with the basic difference being that the Nevada sand here is treated with Ludox SM colloidal silica grout with CS = 6% by weight before shaking. In prototype scale, both experiments refer to a layer 10 m deep, which is excited by 20 cycles of a 2 Hz sinusoidal horizontal input, with a uniform peak base acceleration of 0·2g (against 0·235g in VELACS).

Figure 10 illustrates the test configuration of the dynamic centrifuge test on stabilised sand – that is, its dimensions (23 × 10 m in prototype scale), its ten accelerometers and its five linear variable differential transformers (LVDTs) measuring displacements. No pore pressure transducers were installed, since it was still uncertain whether the pertinent measurements were reliable or not, given the extremely high viscosity of the pore fluid once gelled (Gallagher et al., 2007b). The figure also includes the finite-difference grid used for the analysis of the prototype sand layer, which consisted of 1·0 × 1·0 m quadrilateral zones in Flac (Itasca Inc., 2005). The bottom nodes of this grid were set to follow the prescribed input motion in both directions (zero acceleration in the vertical direction). The lateral boundaries were free to move in the horizontal and vertical directions, but nodes at the same elevation on both boundaries were tied to one another to enforce the same horizontal and vertical displacements of the two boundaries, given the use of a laminar box in the centrifuge test. After stabilisation, the permeability of the sand decreases dramatically, and a value of 10⁻⁹ m/s is considered appropriate for CS = 6% on the basis of Persoff et al. (1999).

Before proceeding to the comparison of measurements and simulations for the centrifuge test on stabilised Nevada sand, it is instructive to focus on the pertinent comparison for the corresponding centrifuge test on the untreated Nevada sand performed for the VELACS project (model 1). Hence, Figure 11 compares the time histories of horizontal acceleration at the ground surface (Figure 11(a)) and excess pore pressure ratio $ru= Δu/σvo′$ at the mid-depth of the 10 m thick layer (Figure 11(b)) from the measurements and the simulation by using NTUA-Sand (whose simulation details may be found in the paper by Andrianopoulos et al. (2010a)). What is important to observe here is that the loose saturated Nevada sand (under this sinusoidal motion) liquefies after t = 4 s at the ground surface and slightly later at greater depths (e.g. t = 6 s), and this leads to intense deamplification (almost nullification) of the horizontal acceleration. Furthermore, note that NTUA-Sand simulates very accurately this liquefaction response, both in terms of acceleration and in terms of pore pressures.

Figure 12 compares the time histories of horizontal acceleration near the ground surface (depth of 2 m) from the centrifuge test on stabilised sand (Gallagher et al., 2007b) with pertinent simulations performed with NTUA-Sand after increasing the plastic modulus constant h _o of the untreated sand, by four times (Figure 12(a)) and ten times (Figure 12(b)) respectively to account for the effect of stabilisation. First, observe that the test data (grey lines) do not depict deamplification of the ground motion; rather, they show significant amplification after the three to four cycles, with peak accelerations reaching 0·55g (in comparison with 0·2g at the base). The measured response includes significant (dilation) spikes, which are essentially the cause of the foregoing amplification. By comparing this acceleration record with the corresponding one in Figure 11(a), one may conclude indirectly that the stabilisation has effectively mitigated the extensive liquefaction that would have occurred without treatment and that this stabilisation may lead to significant amplifications of the seismic ground motion.

The reason for this significant difference in the seismic amplification pattern due to stabilisation may be explained if one focuses on what happens in the two cases. In the untreated sand, liquefaction leads to very small values of the shear wave velocity (10–20% of its initial small strain value) and very large values of damping (around 20%), according to Miwa and Ikeda (2006). Hence, heavily damped out-of-phase vibration is observed generally, which is consistent with the acceleration record in Figure 11(a). On the other hand, stabilisation leads to a relatively stiffer and more dilative response, seemingly due to no (significant) build-up of excess pore pressures (as outlined in the section headed ‘Mechanical response of stabilised sands’). Hence, the stabilised sand layer undergoing a moderate intensity excitation (of 0·2g) exhibits a lightly damped vibration, whose amplification (or not) depends essentially on the ratio of the fundamental soil period over the predominant excitation period. For this case of a 10 m thick layer undergoing an excitation with predominant period of 0·5 s, this ratio is expected to be slightly smaller than 1·0, leading to amplification of the motion, which is consistent with the acceleration record in Figure 12.

Focusing on the comparison of data with simulations (black lines), it is clear that the simulations do not agree with the data. Particularly, the simulation with 4h _o that agrees well with the cyclic element tests (see Figure 9) depicts intense deamplification of the motion after two cycles (Figure 12(a)), and even the simulation with 10h _o (Figure 12(b)) produces minimal benefit in predictive accuracy (the intense deamplification initiates after three cycles) of the centrifuge test. From a performance-based design point of view, in addition to the horizontal accelerations, the seismic settlements are also of interest. Hence, Figure 13 compares the time histories of surface settlements from the centrifuge test on stabilised sand (grey lines; Gallagher et al., 2007b) with pertinent simulations (black lines) performed with NTUA-Sand, after increasing the plastic modulus constant h _o of the untreated sand by four times and ten times to account for the effect of stabilisation (i.e. the same simulations included in Figure 12). Observe that the simulations cannot capture the small settlements (approximately 2 cm) of the stabilised sand layer, regardless of the amount of increase in h _o.

Based on the preceding discussion, the approach of increasing the plastic modulus constant h _o to simulate the effects of gel bonding and encapsulation of the grains does not seem successful on a boundary value problem level. However, these effects have also been alternatively simulated by a decrease in the dilatancy constant A _o. Thus, Figure 14 compares the time histories of horizontal acceleration near the ground surface (depth of 2 m) from the centrifuge test on stabilised sand (Gallagher et al., 2007b) with pertinent simulations performed with NTUA-Sand after decreasing the dilatancy constant A _o of the untreated sand by 2·67 times (Figure 14(a)) and 10 times (Figure 14(b)), to account for the effect of stabilisation. Similar to what is observed in Figure 12, the comparison of data (grey lines) with simulations (black lines) shows a non-satisfactory accuracy. Specifically, the simulation with A _o/2·67 that agrees well with the cyclic element tests (see Figure 9) depicts intense deamplification of the motion after three cycles (Figure 14(a)), while the simulation withA _o/10 (Figure 14(b)) produces some benefit in predictive accuracy (the intense deamplification initiates after eight cycles) of the centrifuge test, but remains far from capturing the true response.

Overall, the simulations of horizontal acceleration with decreased values of the dilatancy constant A _o (Figure 14) are slightly more accurate than those with increased values of the plastic modulus constant h _o (Figure 12). However, a complete comparison requires the evaluation of simulated seismic settlements. Consequently, Figure 15 compares the time histories of surface settlements from the centrifuge test on stabilised sand (grey lines; Gallagher et al., 2007b) with pertinent simulations (black lines) performed with NTUA-Sand, after decreasing the dilatancy constant A _o of the untreated sand by 2·67 times and 10 times to account for the effect of stabilisation (i.e. the same simulations included in Figure 14). Again, the simulations cannot capture the small settlements (approximately 2 cm) of the stabilised sand layer, regardless of the amount of decrease in A _o.

This section evaluates the simulation of stabilisation, by enforcing a reduction in the fluid bulk modulus K, as per Equation 4. This is performed in terms of horizontal accelerations in Figure 16, where the recording near the ground surface (depth of 2 m) from the centrifuge test on stabilised sand (Gallagher et al., 2007b) is compared with pertinent simulations performed with NTUA-Sand after adopting different values of the denominator n in Equation 4, namely n = 50 (Figure 16(a)), n = 500 (Figure 16(b)), n = 1000 (Figure 16(c)) and n → ∞ (Figure 16(d)) to simulate the fully drained conditions – that is, conditions that do not appear in nature, but are considered here only for comparison purposes. Observe in Figure 16(a) that the calibration (of n = 50) on the basis of element tests (Figure 9) fails to predict the system response (it depicts intense deamplification of the motion after four cycles), similar to what is shown in the previous section for the other element test-calibrated approaches. However, further increase in n seems to ameliorate considerably the predicted time histories of horizontal acceleration. For example the use of n = 500 in Figure 16(b) does not allow for intense deamplification of the motion (related to liquefaction) throughout the 20 cycles, while further increase to n = 1000 in Figure 16(c) leads to slight amplification. Finally, the (reference) analysis for fully drained (i.e. no excess pore pressures or almost zero pore fluid bulk modulus) in Figure 16(d) reproduces essentially the recorded amplification of the motion.

Figure 17 compares the time histories of surface settlements from the centrifuge test on stabilised sand (grey lines; Gallagher et al., 2007b) with pertinent simulations (black lines) performed with NTUA-Sand, after decreasing the pore fluid modulus by n = 50, 500, 1000 and n → ∞ (simulating fully drained conditions). As shown in this figure, the value of n = 50 leads to underprediction of the surface settlements, but the simulated settlements increase significantly with an increase in denominator n. Practically speaking, a value of n = 500 provides very good agreement of measurements to simulations, while larger values lead to overpredictions, which reach the unrealistic value of 13 cm for fully drained conditions.

This paper explores the possibility of using existing constitutive models for simulating the seismic response of sands stabilised passively with colloidal silica. Two alternative phenomenological approaches were followed for simulating the stabilisation and both were calibrated on the basis of element test data with considerable accuracy. However, when both thus calibrated methodologies were used in a boundary value problem involving seismic response of a stabilised sand layer, the beneficial effect of stabilisation was generally underpredicted.

This can be considered an artefact of not using a dedicated constitutive model for stabilised sands and can be attributed to the non-uniformity of the seismic response of a soil medium of considerable volume (whose boundaries are far from the area of interest), as compared with the response of a single soil element (under fully prescribed boundary conditions), where the slightest change in model constants (e.g. h_o, A_o) or simulation ingredients (e.g. K) has a clear effect on its response. This is particularly true for the latter approach of reducing the bulk modulus of the pore fluid K, where different areas of the medium may be deforming quite differently given that allowance for pore fluid compressibility. On the other hand, this latter approach is the one that gives the most promising simulation potential, since it is the only one that may accurately predict both horizontal accelerations and seismic settlements, with a value for the denominator n between 500 and 1000 being the optimal choice.

Based on this computational study for the seismic response of passively stabilised sands, the following may be concluded.

• Stabilisation with colloidal silica mainly affects the response of sands at large cyclic shear strains (significant increase in liquefaction resistance), while its effect on the stiffness G and damping ratio ξ values at small and medium strains is relatively small (e.g. 10–25% increase in the value of G_max, but no effect on the G/G_max–γ degradation).

• Existing state-of-the art constitutive models (such as the NTUA-Sand used here) can simulate macroscopically the response of stabilised sands at the element level by mere recalibration of model constants (e.g. increase in plastic modulus, decrease in dilatancy), or by introducing a decrease in the pore fluid modulus in the coupled analysis.

• Successful simulations of element tests on stabilised sands do not necessarily guarantee successful simulations of the system response in boundary value problems involving stabilised sand. However, adopting a decrease in the pore fluid modulus K (= K_w/n) in comparison with that of water (K_w) proves the most promising simulation approach for boundary value problems, with the denominator n = 500–1000 proving the optimal choice.

• The compressibility of colloidal silica after gelation needs to be measured experimentally. Existing evidence shows that it is more compressible compared with water, but the exact value of its bulk modulus, if measured directly, will ascertain the validity of the phenomenological approach proposed herein.

This research was co-financed by the EU (European Social Fund) and Greek national funds through the Operational Programme ‘Education and Lifelong Learning’ of the National Strategic Reference Framework–Research Funding Program: Thales – Investing in Knowledge Society through the European Social Fund.