Interaggregate Forces and Energy Potential Effect on Clay Deformation

This study aims to propose a new model for clay behavior by incorporating physical–chemical effects acting between clay clusters using the Chang and Hicher micromechanical approach. Local mechanisms are thus introduced through repulsive and attractive forces similar to double-layer van der Waals forces, which are obtained from the derivation of energy potentials. A specific study of local parameters and their evolution is conducted using experimental data from the mineralogy variation provoked by the variation of a remolded saturated reconstituted clayey mixture consisting of kaolinite and montmorillonite ranging from 0 to 100% montmorillonite. Using scale transition, the results of micro–macrocalculations under an isotropic loading path show a very good agreement with experimental results.


Introduction
The aggregate structure of saturated clayey sediment makes this medium a mechanically complex system. The aggregates, i.e., the assembly of packing of clay particles, are the result of aggregation phenomenon (Johnson et al. 1996;Veerapaneni and Wiesner 1996) due to physical-chemical reactions between particles and surrounding water (Derjaguin and Landau 1941;Verwey and Overbeek 1948;van Damme 2000). The properties of the water (more or less loaded with salts, for instance) will thus influence in a significant way the microstructure of the aggregates, and the manner in which particles are associated (Collins and McGrown 1974;Mitchell 1976). During the sedimentation and consolidation process (Blewett et al. 2001;Been and Sills 1981) the aggregates will arrange themselves to form a given geometric configuration. The consolidation state notion, defined in the framework of the continuum medium, is directly related to the arrangement in an elementary volume (Biarez and Hicher 1994). Other phenomena that are not discussed in this article, such as cementation, time effect (Mesri and Hayat 1993), etc., can equally develop over time making the structure and the mechanical behavior of the sediment even more complex (Skempton and Northey 1952;Mesri et al. 1975;Mitchell 1976;Leroueil et al. 1979;Burland 1990;Yin et al. 2011, etc.).
The mechanical behavior of this kind of medium has been often approached using continuum media principles in the formal framework of elastoplasticity (Biarez and Hicher 1994). This leads to the development of constitutive models arising from numerous works, among them Roscoe et al. (1958), Schofield and Wroth (1968), Rowe (1962), Dafalias (1987), and others. These works reproduce rather well the experimental reality at the macroscopic scale. However, the microstructural properties (such as arrangement and interparticle forces) are not directly taken into account in these models, and the link to the local mechanisms remains unknown.
Partly because of the lack of experimental data, especially in clayey materials, incorporating in continuous models the local mechanisms that occur at the microscopic scale remains a complex approach. In this way, work can be found in the literature, in particular for granular material, where fabric is regarded as a tensor representing the anisotropy of the granular packing, e.g., Oda (1993); Emeriault and Cambou (1996); Nemat-Nasser and Zhang (2002); Wan and Al-Mamun (2005); Chang and Hicher (2005).
Based on Chang and Hicher (2005), Chang et al. (2009) proposed that clay can be regarded as an assembly of clusters (or aggregates). The deformation of the assembly can then be obtained by integrating the movement of the intercluster contacts in all orientations. The effect of the intercluster interaction was explicitly approached by a contact law. Nonetheless, these models do not consider the mechanisms related to physical-chemical interactions that become significant when the clay's behavior is viewed in terms of their microstructure.
This study proposes a new model for clay behavior by incorporating into the Chang and Hicher (2005) micromechanical approach both an attractive and repulsive force, similar to van der Waals (1873) and double-layer forces (Gouy 1910;Chapman 1913;van Olphen 1977), acting between clay clusters.
The calculations were carried out on an isotropic loading path to obtain the results of two mineralogically different clays, kaolinite and montmorillonite, and a mix between them at 30 and 60% montmorillonite. A specific study of local parameters and their evolution is conducted using experimental data from the mineralogy variation provoked by the variation of the mix-clay ranging from 0 to100% montmorillonite.

Microscopic Structure and Intercluster Interaction
Viewed from a microstructural scale, it appears justified to consider clay as an assembly of aggregates. clay structures can be visually observed from scanning electron microscope (SEM) images as a packing of assembled particles ], the mercury porosimetry test results, especially in the change of the incremental curves during loading, can provide more insights. The studies undertaken on sensitive clays and collected in Delage (2010), show that microstructural variation, during a loading process (oedometric), could be precisely described by the variation of pore size distribution, in which the pores gradually collapse by a decreasing order from the largest pores to smallest (Fig. 5). A similar phenomenon is also shown through the variation of pore size distribution obtained after isotropic loading for intact and remolded reconstituted marine sediments (Fig. 6). The experimental study by Hattab et al. (2013) showed that the variation of porosity under the common levels of stresses was directly related to the pore interaggregates' closure while at the same time the intraaggregates' mean pore diameter remained constant. That can be clearly observed on the parts of the Fig. 6 curves where pores are smaller than 0.05 μm for the remolded sediment, and than 0.03 μm for the intact sediment. This phenomenon permits to conclude that the mechanisms of volumetric strain (for low level of stresses) mainly take place between aggregates. Under the common stress level, the pores within the individual aggregate remain not affected by the loading.
It is customary to model the behavior of soil by treating the medium as an assemblage of grains interacting with each other. The elementary grain, termed as the cluster corresponding to the aggregate in the experimental descriptions, can be deformable. Thus, the macroscopic volume deformation is due to the integration of relative movements of grains and their deformability. Although this view seems to be academic and well accepted for the modeling of a granular medium such as sand, it has not been commonly adapted for the modeling of a clay medium because in a clay medium, the water, which partially or fully saturates the pore space, interacts with the clay particles and plays an important role in clay behavior (van Damme 2000).
The approach proposed in this article considers the clay medium as an assemblage of clusters ( Fig. 7) whose porous space is fully    (Hattab et al. 2010(Hattab et al. , © 2008 Canadian Science Publishing or its licensors. Reproduced with permission) saturated with water. As discussed above, there are two families of pores: (1) pores within the cluster (intracluster pores), and (2) pores between the clusters (intercluster pores). Deformation of a cluster is contributed by the change of intracluster pores, whereas the deformation of an assembly is primarily caused by the change of intercluster pores. As discussed earlier, under the common stress level, the pores within individual aggregate remain not affected by the loading. Therefore, for simplicity, the authors assume that the clusters are not deformable. Now the authors consider the interacting forces between two clusters, which is generated by the physicochemical potential at the surface of clusters. There are two types of interacting forces: repulsive and attractive, which are analogous to the double-layer and van der Waals forces, respectively.
The micromechanical model is based on the conceptual schematic plot shown in Fig. 8, where the two clusters are taken as a spherical shape; the water between clusters has a complicated interaction with the charged surface of the two neighboring clusters.

Equation and Parameters
In the saturated medium, the model considers two neighboring clusters (Cl 1 and Cl 2 ; see Fig. 8) of R radius, with the length connecting the two centers denoted as l c . Under an isotropic loading, the deformation of the assembly is assumed to caused by the change of l c due to the change of intercluster pores. Under the common levels of stresses, the change of intracluster pores is negligible and thus is neglected in the present paper.
The interaction of two neighboring clusters is governed by two types of potential: (1) repulsive type, denoted as w R , analogous to the double-layer potential, and (2) attractive type, denoted as w A , analogous to van der Waals potential. The theories of double layers and van der Waals interaction can easily be found in numerous works, for instance the calculation of the double-layer potential between particles by Verwey and Overbeek (1948), and the detailed calculation of the van der Waals potential between two spherical particles by Hamaker (1937). Based on these works the authors propose simplified expressions for the repulsive potential w R and the attractive potential w A as follows: B andÃ = parameters related to the surface potential of clusters, which depend on the mineralogical nature of the clay and the properties of water between the clusters. It is necessary to introduce a minimum distance (of 2d min ) possible between the two clusters, so that they became able to maintain their boundary surface. Then the following condition must be satisfied: l c > 2R þ 2d min . It is noted that similar types of expressions have been used by Misra and Yang (2010) to describe the behavior in the case of stronger cohesive system in a pseudogranular approach.
From the derivation of the two potentials, the intercluster force is obtained in Eq. (3) which represents the sum of the repulsive and attractive forces Eq.
(3) shows that the intercluster force is a function of distance between the two clusters. At a specific distance ðl c Þ 0 -2R, termed as at-rest distance, the corresponding force is zero. When the distance of the two clusters is greater than this distance, the intercluster force is attractive. Whereas, when the distance of the two clusters is smaller than this distance, the intercluster force becomes repulsive. However, the two clusters in contact do not have true physical contact. The limiting condition is that they are separated by a distance of 2d min , which is measured in micrometers and cannot perceived by naked eyes. Eq.
(3) can be written as a function of the variable l c normalized by d min , the minimum distance. Hence the following expression of f is obtained, consideringÃ n ¼Ã=d min ; R n ¼ R=d min ; and ) Eq. (5) shows that the relationship between the intercluster force and the intercluster distance is influenced by the two ratiosÃ N =B and R=d min . A simulation for the force-displacement relationship is represented in Fig. 9 for a givenÃ N =B value (= 5 in this example). For convenience, the authors define a net distance ratio d=d min , in which d ¼ l c − 2R (see illustration in Fig. 8). When the relative distance between two clusters (d=d min ) is large (>20), the interacting force become null. As the distance decreases, the model first gives an attractive behavior before becoming repulsive. Fig. 9 highlights the influence on the attractive force between two clusters, which increases when the R=d min ratio increases, showing a more and more pronounced peak that occurs at a smaller and smaller distance d=d min . For lower values of R=d min , the attractive part of the curve disappears and the interaction between clusters becomes solely repulsive.
This interacting force depends not only on the variation R=d min ratio, but also on the ratioÃ N =B, which is related to the charge potential at the surface of the clusters. Simulation results by using Eq. (5) are presented Fig. 10 to show the effect ofÃ N =B with R=d min ¼ 10. Fig. 10 highlights the increasingly marked prevalence of the attractive part of the force, especially when thẽ A N =B ratio ≥5 (continuous line in Fig. 10). TheÃ N =B ratio affects the peak only moderately, but controls the slope of the curve so that the latter can switch from the repulsion domain towards the attraction domain (Fig. 10). It can be noticed that the range of the curves for this simulation is relatively small. As an example, for the minimum possible value d=d min ¼ 2 (d ¼ l c − 2R; see illustration in Fig. 8), the value of f B ranges in the interval ½0.9; −0.2. This range can be affected by the value of R=d min , as shown in Fig. 11. The example of simulation exhibited in Fig. 10 with R=d min ¼ 10 is plotted in Fig. 11 (shaded area). The simulations show that, when R=d min ¼ 30, the range of the curves is much larger, and the  Additionally, it is necessary to identify the range of values for R=d min andÃ N =B that is physically possible (i.e., the mechanisms) and that can lead to realistic experimental mechanical behavior of clays. A discussion of this issue follows.

Discussion
The behavior of remolded saturated reconstituted clay in isotropic loading is now quite well known especially on the macroscopic level; see for instance experimental works on kaolinite clay by Biarez and Hicher (1994) or by Hattab and Hicher (2004). As shown in Fig. 20, for different remolded clay behavior on an isotropic path, during isotropic loading the force between clay clusters is repulsive, and increases with a decrease of l value, which corresponds to the decrease of the void ratio (or water content). Thus the force should be repulsive during the compression stage, although during the slurry state before the clay is consolidated, a slight attraction force can be justified between the clay clusters.
Therefore it can be assumed that as long as the clusters are defined by the same boundary surface (so damage or blend of clusters), the attractive total potential between clusters is physically unacceptable during compression stage. For instance, the case in Fig. 11 for R=d min ¼ 30,Ã N =B ≥ 5 is considered to be physically impossible.
Making this assumption, and considering the simulation results for different ratios R=d min and different values ofÃ N =B, an available domain of parameters for clays can be proposed on the normalized plane (Ã N =B; R=d min ). This domain is represented in Fig. 12 by the red circles.
It can be noticed in Fig. 12 that, with high values of R=d min and low values ofÃ N =B, the interaction between clusters tends toward a contact (densely packed cluster) interaction problem. With high values ofÃ N =B and low values of R=d min , the interaction between clusters is dominated by attractive potential interaction. With both R=d min andÃ N =B values approaching zero, the interaction between clusters is dominated by repulsive potential interaction.

Principle of Integration and Basic Equations
The stress-strain relation at the continuous medium scale can be determined by the integration of local relations obtained at the elementary system, which is composed of cluster/water/cluster, as illustrated in Fig. 8. The micro-and macrovariables can thus be related.
The approach used in the transition from the micro-to the macrovariables is based on the micromechanical model of Chang (1988) Yin et al. 2011), in which the clay is assimilated into an assembly of clusters. The orientation of two interacting clusters is along the direction of a length vector joining the centers of two clusters. The macroscopic behavior law can then be approximated by averaging the behavior of two interacting clusters in all orientations. Taking two interacting clusters oriented in the direction α, the local forces f α i and the local displacement δ α i can be expressed as f α i ¼ ff α n ; f α s ; f α t g and δ α i ¼ fδ α n ; δ α s ; δ α t g, where n, s, and t = the unit vectors in the three directions of the local coordinate.
In the current model, the two-cluster systems will be integrated to the elementary systems (Fig. 13), in which the total local forces are deduced from the sum of the potentials activated within elementary systems. The sum of the local forces in the assembly may be calculated as follows: Using the equilibrium principle, the stress tensor of an assembly can be obtained from the local forces in all directions. Further details, particularly for the definition of the stress in a grain assembly, can be found in works such as those of Christofferson et al. (1981), Rothenburg and Selvadurai (1981), and Cambou and Jean (2001). It is generally accepted that where l a = branch vector connecting the center of two neighboring clusters, V = elementary volume, and N c = total number of branches of the assembly.

Relation Between Macrostress and Local Forces under Isotropic Loading Condition
Considering the isotropic loading case, Eqs. (6) and (7) permit to deduce the expression of the stress as follows: where f = magnitude of intercluster force and δ ij = Kronecker symbol. Hence, the mean stress

Relation between Macrostrain and Local Displacements under Isotropic Loading Condition
Assuming that the work done at the representative elementary volume (REV) level is equal to the work done by all the elementary systems, we obtain the following equation: Strain under isotropic loading can thus be related to the displacement as follows: In Eq. (12), δ = displacement between clusters belonging to the same elementary system and δ ij = Kronecker symbol. δ is easily deduced as a function of the volumetric deformation

Micro-Macro Transition and the Case of Isotropic Loading
The relationship can also be established to link the local displacement to the void ratio e through Eq. (14), or to the water content w (expressed in percentage) through Eq. (15) e ¼ e 0 þ 3δ l 0 ð1 þ e 0 Þ ð 14Þ Fig. 13. Local coordinate system and potential between clusters into the elementary system Following the kinematic approach represented in Fig. 14, the mean stress can be directly linked to the volumetric strain in an isotropic loading case. Using Eq. (14) or Eq. (15), the mean stress can be directly linked to the void ratio or water content.

Discussion: Calculation and Meaning of Each Parameter
The parameters of the model are: 1. For the microstructural description • R 0 = mean radius of clusters that depends on the material. R 0 can be estimated from SEM photoanalysis. • N=V = contact number per unit volume with N=V ¼ 12=½ðπ=3Þð2RÞ 3 ð1 þ eÞ (e = void ratio and R = mean size of clusters), which makes the number of contacts per unit volume vary during the deformation. • l 0 = distance joining the center of two neighboring clusters at the liquid limit state of the material. This local parameter can be calculated as a function of the equivalent void ratio e L with l 0 ¼ 2R 0 ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi ð1 þ e L Þ=ð1 þ 0.35Þ 3 p .

For the intercluster properties
•B andÃ = parameters related to the surface potential of clusters and of the properties of the intercluster water. • d min = minimum distance possible between the surfaces of two neighboring clusters, also depending on the surface potential of clusters. This parameter imposes the following condition: l 0 =d min > 2R 0 =d min þ 2. • The identification and calculation methods forB,Ã, and d min parameters are presented in the section "Discussion on the calculated model parameters." 3. For the properties of the assembly • γ s =γ w , weight density of the solid matter.

Verification with the Experimental Results
The micro-macro calculations, using the framework of this work, were first performed for isotropic loading tests on two kinds of clay with different mineralogy: kaolin (kaolinite P300) and montmorillonite (Greek clay). The simulation was then performed on mixtures from these two clays with different proportions (see Table 1 for the geotechnical properties of the mix-clays used). One can also refer to Hammad et al. (2013), where experimental results of kaolin/montmorillonite mixtures, consisting of a detailed analysis of geotechnical properties and the behavior of a onedimensional compression path, can be found. In the Appendix, a description is given of the change of properties for clays mixed with different fractions of kaolinite and montmorillonite.

Discussion on the Calculated Model Parameters
Method of Determining Intercluster ParametersB,Ã, and d min Parameter determination is based on the experimental data for the isotropic loading of the two basis clays, kaolin and montmorillonite, and for the mixtures between these two kinds of clay. The clay named M 35% contains 35% montmorillonite and M 65% contains 65% montmorillonite. The results are presented in terms of the water content w with respect to isotropic effective stress p 0 0 [ Fig. 15(a)]. Following Biarez's correlations (Biarez and Favre 1975;Favre Hattab 2008), one can assume that the liquid limit w L and plastic limit w P correspond to p 0 0 of 6.5 kPa and p 0 0 of 1,000 kPa, respectively. Those two limits are thus introduced in the (w-p 0 0 ) plane for all materials. The discrepancy observed between the plastic limit measurement and the water content obtained from the isotropic tests at p 0 0 ¼ 1,000 kPa [ Fig. 15(a)] is due to experimental difficulties in the w P measurement and in the real possibilities of obtaining a good value with correct precision. Fig. 15(b) is obtained by using global w-p 0 0 data. The obtained experimental curves express the relations between local distance l c (branch length joining the centers of two neighboring clusters) and f, the local force that develops between two neighboring clusters.

Calculation of l c and f from Experimental Global Data
Considering spheres of equal radius R, one can express the porosity n as Let R 0 = radius of clusters, e = void ratio, N sph = number of spheres, and V = total volume where N sph is counted. Hence Assuming the value of the radius ¼ l c =2 in a closely packed assembly, and assuming the closely packed assembly has a void ratio 0.35 (corresponding to hexagonal packing), then Thus, the branch length can be related to the void ratio as Fig. 14. Micro-to macrotransition kinematic method by Chang and Liao (1994)  The total number of contact N c per unit volume can also be calculated using the void ratio parameter (assuming coordination number = 12 for the assembly), hence and then using the equations in Section 3, the force-displacement relationship becomes Local Properties and Analysis of Experimental Results (l c -f ) The microrelation between f and l c is expressed by Eq.
(3) [also referenced by Eq. (22)]. This relationship contains three parameters B,Ã, and d min , which represent the physical-chemical properties of the clays Three particular points on the experimental (l c -f) curves of Fig. 15(b) will permit to deduceB,Ã, and d min by resolving the system Eq. (23). Thus, for a given clayey material, the three parameters can be calibrated from the isotropic compression curve obtained in experiments.
For a given clay, the following approach can be proposed: two specific points corresponding to liquid limit and the plastic limit, i.e., w L -p 0 0 of 6.5 kPa and w P -p 0 0 of 1,000 kPa; and one a priori selected point on the curve, for instance corresponding to a p 0 0 of 300 kPa. The measured water content w associated with p 0 0 ¼ 300 kPa is termed w 300 . Hence the equation system to resolve is in the following form: where the three parametersB,Ã, and d min can be determined (see Table 2). The calculations were performed using the model parameters given in Table 2. Those were calculated for each clay mixture, following their definition and relations described in the previous section. In Fig. 12, the normalized parameters are calculated for the mix-clays, and give Table 3, permitting location of the valid domain on the (Ã N =B; R=d min ) plane for this kind of material (Fig. 16). In this plane, the kaolin parameters suggest a cluster interaction much closer to the contact (closely packed cluster) type of interaction than the other mixtures are. On the other hand, those  calculated parameters (whereB,Ã, and d min vary significantly from K to M 35% ) reflect the experimental observations and the strong influence of the montmorillonite fraction quite well, especially from 35% in the clay mixture, in both mineralogy and mechanical behavior . For instance, from the mineralogical point of view, one can observe the change in clay activity when the montmorillonite fraction increases, which is directly related to the complex physical-chemical interaction acting at the local scale. Fig. 21, presented in the Appendix, shows that normal activity appears around 35% montmorillonite in the mixture according to the Skempton (1953) classification. Before that the clay remains inactive.

Calculation Results
The calculation results are then compared with all experimental curves. As shown in Fig. 17 on the local scale, the obtained simulation curves fit well the experimental curves previously deduced from the global experimental results. For all materials, whatever the percentage of montmorillonite, isotropic compression leads to a decrease of l c , the distance between the centers of two neighboring clusters, and an increase of intercluster force. One can also note that for a given macroisotropic load, the corresponding local intercluster forces became larger when the percentage of montmorillonite increased. Additionally, the value of l c over the cluster radius increased with the percentage of montmorillonite. These properties can be seen in Fig. 18, which describes the intercluster force (and l c =R 0 ) change, corresponding to the isotropic stress of 1,000 kPa, with respect to the montmorillonite fraction. The observed variation in the local properties, highlighted by changes in f and l c , is due to the more significant interactions under the elementary system cluster/water/cluster provoked by increases of the more reactive montmorillonite fraction in the material. This montmorillonite effect is also given by the simulations on large l c values, i.e., values ≥l 0 (the clay being the slurry state), where negative forces developed. As one can observe in Fig. 17, the model predicts a very weak attraction phenomenon for the kaolin that increases when the montmorillonite fraction increases in the mixture. Good agreement between experimental results and the simulations is achieved from the calculation of the micro-macro transition. This agreement can be observed in Fig. 19, which exhibits the results for the mix-clays starting from the pure kaolin clay to the 100% montmorillonite clay.

Conclusion
Based on the Chang and Hicher micromechanical approach, a new model for clay behavior is developed in this study, considering energy potentials in the local mechanisms between clay aggregates. This approach allows taking into account physical-chemical effects between clusters through repulsive and attractive forces, whose relations are similar to double layers and van der Waals forces. The parameter choice of the local law permits the estimation of these forces between clusters, depending on the energy potential of the mineralogical nature of the clay (more or less active in a global experimental approach), and on the properties of the water occupying the pore space. The local law consists of four parameters (B,Ã, d min , and R), which permit estimation of the influence of repulsive and attractive energy potential. Simulation results of the local law highlighted an available domain of possible parameters for the clayey material. Densely packed clusters appear as one of the boundaries of this domain.
Calculations using micro-macroscale transition were thus performed for isotropic loading on the mixture clay made of kaolinite and montmorillonite, so that mineralogy variation of the tested clay is assured. A specific method was also developed for determining the values ofB,Ã, and d min local parameters based on a calibration from the experimental isotropic compression curve.
A good agreement between experimental results and simulations is then highlighted, showing the behavior variation related to the variation of the montmorillonite fraction under an isotropic loading path.

Appendix. Tested Clayey Material Properties
The two types of tested natural clays were a yellow kaolinite marketed under the name of Kaolin P300, and a montmorillonite (smectite) designated as Greek clay. Both are provided as dried powder. Mix-materials are made by mixing the two kinds of clays; M x = percentage x of montmorillonite. For example, the Greek clay M 100 has 100% montmorillonite, and the mixture with 35% montmorillonite is called M 35 . Kaolin P300 (M 0 ) is noted K.

Kaolin P300
Kaolin P300 is marketed by Dousselin (Fontaines sur Saône, Rhône, France). Analysis by X ray diffraction shows that the clay consists mainly of kaolin traces, illite, and some quartz. Observations by SEM  showed that the kaolinite particles present themselves as rigid hexagonal plates whose size varies between 0.5 and 3 μm. The particle itself consists of a set of stacked sheets.

Greek Clay
The Greek clay is a calcium montmorillonite, taken from Milos Island in the Aegean Sea. Observation by SEM ) shows a high porosity fabric presenting a complex honeycomb network with overlapping particles.

Atterberg's Consistency Limits and Activity Variation according to the Percentage of Montmorillonite
Consistency limits (w L and w P ) of the mix-material have been determined by Hammad et al. (2013) from two different methods of    Fig. 20 shows that w L and w P limits vary linearly with montmorillonite percentage. This could be expected, considering the activity values (A a ) of the mix-material (Fig. 21). According to the Skempton (1953) classification, the mix-clay evolves from inactive for K (A a ¼ 0.25) to active for M 100% (A a ¼ 1.31), with a normal activity appearing between 35 and 50% montmorillonite. This demonstrates the appearance of a more complex local behavior when the montmorillonite fraction increases in the mix-clay.