Capillary bridges between spherical particles under suction control: Rupture distances and capillary forces

Journal Pre-proof Capillary bridges between spherical particles under suction control: rupture distances and capillary forces Chao-Fa Zhao, Niels P. Kruyt, Olivier Millet PII:

S0032-5910(19)30802-2

DOI:

https://doi.org/10.1016/j.powtec.2019.09.093

Reference:

PTEC 14807

To appear in:

Powder Technology

Received Date: 5 August 2019 Revised Date:

23 September 2019

Accepted Date: 29 September 2019

Please cite this article as: C.-F. Zhao, N.P. Kruyt, O. Millet, Capillary bridges between spherical particles under suction control: rupture distances and capillary forces, Powder Technology (2019), doi: https:// doi.org/10.1016/j.powtec.2019.09.093. This is a PDF file of an article that has undergone enhancements after acceptance, such as the addition of a cover page and metadata, and formatting for readability, but it is not yet the definitive version of record. This version will undergo additional copyediting, typesetting and review before it is published in its final form, but we are providing this version to give early visibility of the article. Please note that, during the production process, errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain. © 2019 Published by Elsevier B.V.

Graphical Abstract Capillary forces between two spherical particles: influence of suction

Large particle Contact angle

Small particle

RL Capillary bridge

Separation distance

Rs

2

Capillary bridges between spherical particles under suction control: rupture distances and capillary forces

3

Chao-Fa Zhao1,2 , Niels P. Kruyt1∗ , Olivier Millet2

1

4

1

Department of Mechanical Engineering, University of Twente, P.O. Box 217, Enschede 7500 AE, the Netherlands 2

5

6

LaSIE-UMR CNRS, Universit´e de La Rochelle, 23 Avenue Albert Einstein, La Rochelle 17000, France

Abstract

7

This study focuses on properties of axisymmetric capillary bridges between spherical particles in

8

the pendular regime under suction control. Using the toroidal approximation in combination with the

9

governing Young-Laplace equation, analytical expressions for the rupture distance (dependent on the

10

suction) and for the capillary force (dependent on the suction and the interparticle separation distance)

11

have been obtained that do not involve any calibrated coefficient. The developed analytical expres-

12

sions are effective for values of the dimensionless suction larger than 104 . To predict capillary forces

13

and rupture distances for a wider range of values of the dimensionless suction, closed-form expressions

14

for rupture distances and capillary forces have been obtained by curve-fitting to a large dataset of

15

numerical solutions of the Young-Laplace equation. It is shown that these expressions can also be em-

16

ployed for capillary bridges between spheres with unequal radii when the Derjaguin radius is employed.

17 18

Keywords:

19

capillary bridge; suction control; rupture distance; capillary force; Young-Laplace equation; toroidal

20

approximation

∗

Corresponding author; e-mail addresses: [email protected] (C.F. Zhao); [email protected] (N.P. Kruyt); [email protected] (O. Millet); tel: +31-53-489 2528 Preprint submitted to Powder Technology

October 11, 2019

21

1. Introduction

22

Different types of cohesive forces, such as van der Waals forces [1], electrostatic forces [2] and

23

capillary forces [3] occur in granular materials. For wet granular materials, capillary forces are larger

24

than the other cohesive forces.

25

Hence, cohesion in wet granular materials originates from capillary forces between the constituting

26

particles [4–8]. With decreasing liquid content, various regimes are encountered [7, 9–12]: (1) the

27

capillary regime where the particles are fully immersed in liquid, (2) the funicular regime with coalesced

28

capillary bridges, with liquid clusters that are simultaneously in contact with at least three particles

29

and (3) the pendular regime where a capillary bridge is associated with exactly two particles.

30

The current study focuses on capillary bridges in the pendular regime and situations where the

31

surface of the capillary bridge is axisymmetric. The effect of viscosity is not considered here, since the

32

capillary bridges are assumed to be in quasi-static equilibrium. The influence of gravity on capillary

33

bridges that is investigated in [13] is neglected in the present study.

34

Important properties of capillary bridges between two spheres in the pendular regime are the

35

rupture distance (i.e. the maximum separation distance between the particles for which a capillary

36

bridge exists) and the capillary force. The dependence of the rupture distance on the volume of the

37

capillary bridge and of the capillary force on the volume of the capillary bridge and on the separation

38

distance has been extensively investigated [3, 14–18]. This type of problem is called volume control

39

here.

40

Besides properties of capillary bridges under volume control (where the volume is kept constant),

41

such properties under suction control (where the pressure difference between the surrounding gas and

42

the liquid inside the capillary bridge is kept constant) are also important. For instance, in Discrete

43

Element Method [19] simulations of the behaviour of wet granular materials, expressions for given

44

capillary bridge volume are appropriate for undrained tests, whereas expressions for given suction

45

pressure are appropriate for drained tests [20–23]. In practice, the behaviour of wet granular materials

46

under suction control is of great interest when the humidity of the surrounding gas is controlled.

47

The relationship between suction and humidity under conditions of thermodynamic equilibrium is 2

48

discussed in detail in Section 2.5. The current study focuses on properties of capillary bridges under

49

suction control.

50

Properties of capillary bridges are determined by the capillary bridge surface geometry, which

51

is described by the well-known Young-Laplace equation [3, 14–17, 24–26], according to which the

52

mean curvature of the surface is constant. It has been shown experimentally that the Young-Laplace

53

equation accurately describes the properties of capillary bridges [27–35]. According to the Young-

54

Laplace equation [14–18, 24], the mean curvature of the capillary bridge surface H is proportional to

55

the value of suction ∆p

H=

∆p γ

=

pext −pint γ

H∗ = H · R ,

(1)

56

where pext and pint are the pressures outside and inside the capillary bridge; γ is the surface tension

57

between water and gas (γ = 72.8 × 10−3 N/m for water-air interfaces at a temperature of 20◦ C); H is

58

made dimensionless using the Derjaguin radius R [18, 36], giving the dimensionless suction H ∗ . The

59

(effective) Derjaguin radius R corresponding to a small spherical particle with radius Rs and a large

60

spherical particle with radius RL is defined by

R=

2Rs RL . Rs + RL

(2)

61

Properties of capillary bridges under suction control have been studied using the following com-

62

plimentary approaches: numerically [20, 22, 30, 37–39], analytically [37, 38] and experimentally [30].

63

Numerical solutions of the nonlinear Young-Laplace equation are effectively exact, but are cumbersome

64

and time-consuming when used a large number of times. Analytical methods have a clear physical

65

meaning, but they have a limited range of validity as they involve some approximations. Experiments

66

directly provide information on the capillary bridge properties, but they are time-consuming and dif-

67

ficult to perform due to the effects of gravity and particle surface roughness. To the authors’ best

68

knowledge, explicit expressions for the capillary force under suction control have not been reported in

69

the literature.

3

70

According to analyses of numerical solutions to the Young-Laplace equation, the relation between

71

rupture distance and capillary bridge volume under suction control is very different from that under

72

volume control [22, 37–40]. Numerical solutions also show that the capillary force decreases with

73

increasing value of suction (for the same separation distance).

74

Analytical theories are generally based on the toroidal approximation [3, 16, 37, 38, 40, 41]. In

75

this approximation the meridional profile of capillary bridges is represented by part of a circle. Under

76

volume control, expressions for the capillary force have been reported [37, 38, 40, 41]. Recently,

77

analytical theories based on the ellipse approximation [14, 15, 18] have been developed for volume

78

control. In this approximation, the meridional profile of capillary bridges is represented by part of

79

an ellipse, and the Young-Laplace equation has been taken into account, contrary to the toroidal

80

approximation. For given capillary bridge volume and separation distance, expressions for the rupture

81

distance and the capillary force have been derived, in which no calibrated coefficients are required.

82

Gras et al. [30] have experimentally measured capillary forces and rupture distances of capillary

83

bridges between two equal-sized spherical particles under suction control. Their experimental findings

84

are consistent with those obtained from numerical solutions to the Young-Laplace equation [20, 22,

85

30, 37–39]. They have also suggested a closed-form expression for the rupture distance, based on

86

their experimental results and numerical solutions to the Young-Laplace equation. However, the

87

dimensionless suction in their closed-form expression for the rupture distance is limited to small values

88

(0 ≤ H ∗ ≤ 15).

89

For capillary bridges under volume control or suction control, it has been numerically shown that

90

two solutions to the Young-Laplace equation may exist [14–16, 18, 22, 37]. To select the physically-

91

relevant branch (also referred to as the stable solution), various criteria have been developed in [22,

92

25, 42]. However, the effectiveness of these criteria has not been verified for both volume and suction

93

control conditions.

94

In experimental studies with wet granular materials, it is not possible to continually increase the

95

liquid pressure to large values while keeping the gas pressure at atmospheric level, since the liquid may

96

vaporize under high (negative) pressure levels [43, 44]. Therefore, the axis translation technique that

4

97

is suitable for lower values of suction and the vapour equilibrium method that is applicable for higher

98

values of suction have been developed. In the vapour equilibrium method the relative humidity of the

99

gas is controlled. As described in [44, 45], the vapour equilibrium technique cannot accurately control

100

the dimensionless suction when the relative humidity h is close to 100%, due to the strongly nonlinear

101

h−H ∗ curve that requires a high-precision humidity sensor. The relation between the relative humidity

102

h and the dimensionless suction H ∗ is generally given by the Kelvin equation [37, 38, 46]. Therefore,

103

it is valuable to (also) formulate the developed expressions for the rupture distance and the capillary

104

force in terms of humidity and separation distance.

105

The range of dimensionless suction H ∗ of interest in the current study is H ∗ ≥ 1. For H ∗ < 1, the

106

capillary force is relatively small and hence it is not important. For granular soils, the mean particle

107

size is about 7.4 × 10−5 m [44]; accordingly, the dimensionless suction H ∗ is approximately larger than

108

102 (H ∗ ' 102 ). The values of H ∗ considered here have a much larger range than that in [30], in

109

which values 0 ≤ H ∗ ≤ 15 have been considered.

110

In the present study, properties of capillary bridges under suction control are determined using

111

the complimentary approaches of an analytical theory as well as analyses of numerical solutions to the

112

governing Young-Laplace equation. The analytical theory is valid for large values of the dimensionless

113

suction (H ∗ ≥ 104 ) that correspond to small capillary bridge volumes. A large dataset of numerical

114

solutions to the Young-Laplace equation is used to formulate closed-form expressions (by curve-fitting)

115

that are valid for a wider range of dimensionless suction H ∗ ≥ 1, and thus both small and large capillary

116

bridge volumes have been taken into account. Specific objectives of this study are to:

117

• demonstrate the effectiveness of the surface free-energy based criterion for selecting the physically-

118

relevant solutions to the governing Young-Laplace equation (obtained by a high-resolution inte-

119

gration method) under both suction and volume control conditions;

120

121

• develop an analytical theory for the rupture distance and the capillary force under suction control, by employing the toroidal approximation;

122

• formulate closed-form expressions for the rupture distance and the capillary force under suction

123

control by curve-fitting to a large dataset of numerical solutions of the governing Young-Laplace

5

124

equation;

125

• consider the influence of the particle size ratio on the rupture distance and the capillary force

126

under suction control by using the Derjaguin radius (which is valid for small capillary bridge

127

volumes).

128

The outline of this study is as follows. Section 2 gives the description of the geometry of the cap-

129

illary bridge, the numerical solution method for the governing Young-Laplace equation, the criterion

130

for selecting the physically-relevant solutions and the vapour equilibrium method for suction control.

131

In Section 3 analytical expressions are derived for the rupture distance and the capillary force for

132

given suction and separation distance, based on the toroidal approximation. In Section 4 closed-form

133

expressions are given for the capillary force and the rupture distance by curve-fitting to the large

134

dataset of numerical solutions to the Young-Laplace equation. Section 5 develops the obtained ex-

135

pressions for unequal-sized spherical particles by using the Derjaguin radius. Section 6 summarizes

136

the contributions of this study and discusses further developments.

137

2. Capillary bridge geometry and numerical solution of the Young-Laplace equation

138

In this Section, the geometry of the capillary bridge is described. The Young-Laplace equation is

139

given that describes the geometry of the capillary bridge surface. The numerical solution method for

140

the Young-Laplace equation is presented. The surface free-energy criterion for selecting physically-

141

relevant branches of the numerical solutions is then discussed, for both suction and volume control

142

conditions. The method for controlling suction by the vapour equilibrium method using the Kelvin

143

equation is also described.

144

2.1. Geometry of the capillary bridge

145

An axially symmetric capillary bridge between two spherical particles with radii Rs and RL (RL ≥

146

Rs ) is considered. The subscripts s and L refer to properties of the small and large particle, respectively

147

(conform the notation in [18]).

148

As shown in Fig. 1, the axial coordinate is denoted by X, and the meridional profile of the capillary

149

bridge is described by Y (X). The neck where the radius of the capillary bridge is smallest is taken at 6

150

X = 0. The neck radius is denoted by Y0 . Hence

(3)

Y (0) = Y0 Y 0 (0) = 0 ,

151

where the slope of the meridional profile Y (X) is denoted by Y 0 (X) = dY /dX.

152

The meridional coordinates (Xcs , Ycs ) and (XcL , YcL ) of the three-phase contact circles at the small

153

and large spheres are related to separation distances (Ss and SL ), and half-filling angles (δs and δL )

154

(see Fig. 1) by [18]

Xcs = Ss + Rs (1 − cos δs ) Ycs = Rs sin δs XcL = −SL − RL (1 − cos δL ) YcL = RL sin δL .

155

156

(4)

For small capillary bridge volumes, the contact radii Ycs and YcL are very small. Therefore, Eq.(4) can be approximated by

Xcs ∼ = Ss +

2 1 Ycs 2 Rs

Ycs ∼ = Rs δs XcL ∼ = −SL −

2

1 YcL 2 RL

YcL ∼ = RL δL .

(5)

157

In addition, the slopes of the meridional profiles at the contact circles involve the contact angle θ (see

158

also Fig. 1). These slopes are determined by [18] √ Y 0 (Xcs ) = cot(δs + θ) =

159

2 −T Y Rs2 −Ycs cs

Ycs +T

√

2 Rs2 −Ycs

√ Y 0 (XcL ) = − cot(δL + θ) = −

0 ≡ Ycs

2 −Y 2 −T Y RL cL cL

YcL +T

√

2 −Y 2 RL cL

0 , ≡ YcL

(6)

where T is a short-hand notation for the term

T = tan θ .

160

Here the contact angle θ is in radians.

161

2.2. Young-Laplace equation and capillary force

(7)

162

It is well-known that the geometry of the capillary bridge surface is described by the Young-

163

Laplace equation [14–18, 24], according to which the mean curvature of the surface is constant. For

164

the axially-symmetric capillary bridges considered in this study, the Young-Laplace equation can be 7

165

expressed as

K[Y ](X) ≡

Y 00 1 1 − =H, (1 + Y 02 )3/2 Y (1 + Y 02 )1/2

(8)

166

where H is given by Eq.(1) and K[Y ] is the mean curvature operator; Y 0 = dY /dX and Y 00 =

167

d2 Y /dX 2 . For this nonlinear second-order ordinary differential equation, boundary conditions are

168

given by Eqs.(3), (4) and (6).

169

The Young-Laplace equation (8) has a first integral

Λ[Y ](X) ≡ √

1 Y + HY 2 = λ , 02 2 1+Y

(9)

170

where λ is a constant. This expression corresponds to axial force equilibrium along the capillary bridge

171

[14, 16, 25]. By substitution of Eqs.(3) and (6) into Eq.(9), it follows that the mean curvature H and

172

the constant λ are given by

Y0 − √

YcL 1+Y 02 cL 2 −Y 2 YcL 0

173

Y0 − √ Ycs

=

1 2H

=

02 1+Ycs

2 −Y 2 Ycs 0

λ = Y0 + 12 HY02 .

(10)

The capillary force Fcap between two spherical particles is then determined by

Fcap = λ (2πγ) ,

174

where Fcap and λ are both dimensional (contrary to in [14, 15]).

175

2.3. Numerical solutions of the Young-Laplace equation

(11)

176

When values of the contact radius of the small particle Ycs and the neck radius Y0 are specified,

177

the meridional profile X(Y ), the dual of Y (X), is obtained by using a high-resolution numerical

178

integration method, as described in [18]. The capillary bridge volume (for each particle, i.e. the part

179

up to the neck located at X = 0) is evaluated by numerical integration from the meridional profile

180

Y (X). The total capillary bridge volume V and total separation distance 2S are determined by the

8

181

capillary bridge volumes and the separation distances of each particle

V = Vs + VL 2S = Ss + SL ,

(12)

182

where Vs and VL are the capillary bridge volumes for the small and large particle respectively; Ss and

183

SL are the separation distances for the small and large particle respectively, as shown in Fig. 1.

184

Zhao et al. [18] compiled a database consisting of a large number of numerical solutions to the

185

Young-Laplace equation, for various values of half-filling angles δs ≤ 60◦ , ratios of neck radius to

186

contact radius 0 ≤ Y0 /Ycs ≤ 1 and size ratios of RL /Rs =1, 2, 8 and 128 and contact angles of θ =

187

0◦ , 20◦ and 40◦ . For each half-filling angle δs and ratio Y0 /Ycs , the resulting separation distances (Ss ,

188

SL and S), capillary bridge volumes (Vs , VL and V ), capillary force (λ) and mean curvature (H) have

189

been stored. For any combination of mean curvature and separation distance (H, S) of interest, the

190

(scaled) capillary force λ is obtained by interpolating the dataset, using the method described in [14].

191

The same database and interpolation method have been employed here.

192

2.4. Stable branches of solutions to the Young-Laplace equation

193

For specified values of the neck radius Y0 and contact radius Yc , the solution for the meridional

194

profile Y (X), as determined by direction integration (see for instance Eq.(13) in [18]), is unique.

195

However, this may not be the case with specified volume V and separation S (so under volume

196

control) [14–16, 18] or with specified suction H and separation S (so under suction control) [22].

197

In general, only a single branch of such solutions is stable and is observed experimentally. Criteria

198

based on the minimum of some energy have been proposed to identify the stable branch of solutions

199

[14–16, 22, 25, 47, 48].

200

The energy considered in [16, 47, 48] is the total surface energy of the system Esurf , given by the

201

sum of surface areas weighted by the corresponding interfacial energies per area (i.e. surface tensions)

Esurf = γAlg + γls Als + γsg Asg .

202

(13)

Here γ, γls and γsg are the surface tensions for liquid-gas, liquid-solid and solid-gas interfaces, and the 9

203

204

corresponding interfacial surface areas are Alg , Als and Asg respectively. An energy that also involves other sources of free energy has been considered in [22]

Esurf+∆pV = Esurf + ∆pV .

205

(14)

where ∆p is the pressure difference between gas and liquid; V is the capillary bridge volume.

206

Besides the energies in Eqs.(13) and (14), the energy considered by [25] is equal to that in Eq.(13)

207

but shifted by the constant 4πγR and the energy Esurf+pext V = Esurf + pa V in [42] does not involve

208

the pressure difference ∆p that determines the mean curvature in the Young-Laplace equation (8).

2

209

Under volume control (with fixed capillary bridge volume V ), Lian et al. [16] have demonstrated

210

that the criterion of minimal Esurf is effective in identifying the stable branch of solutions of the Young-

211

Laplace equation. They also found that the minimum free surface energy criterion is equivalent to the

212

minimum surface energy of the liquid-gas interface, i.e. γAlg . However, for capillary bridges under

213

suction control, the validity of the criterion of a minimum of the energy defined in Eq.(13) has not

214

been investigated in literature.

215

Under suction control (with fixed suction H), Duriez and Wan [22] have shown that the criterion

216

of minimal Esurf+∆pV is effective in identifying the stable branch of solutions of the Young-Laplace

217

equation. However, the appropriateness of this criterion for capillary bridges under volume control

218

has not been investigated.

219

The suitability of criteria based on Esurf and Esurf+∆pV for capillary bridges under volume as well

220

as suction control, respectively, has been investigated in Appendix A. There it is shown that the

221

criterion of a minimum of the surface free energy defined by Eq.(13) is appropriate for identifying

222

the stable branch of solutions to the Young-Laplace equation, under both volume and suction control

223

conditions.

224

2.5. Relationship between dimensionless suction and relative humidity

225

The suction in wet granular materials can be controlled through several techniques. Among these

226

methods, the vapour equilibrium method has been extensively employed for unsaturated soils [44, 45].

10

227

With this method, the suction is controlled by varying the vapour pressure, which is equivalent to

228

changing of the relative humidity. The theoretical basis of this method is provided by the Kelvin

229

equation [37, 38, 46, 49] that relates vapour pressure pv and saturation vapour pressure psat to the

230

dimensionless suction H ∗ under conditions of thermodynamic equilibrium

ln(

pv −∆p −γH ∗ )= = , psat ρl Rv Tt ρl Rv Tt R

(15)

231

Here the ratio pv /psat is defined as the relative humidity h; ρl is the density of liquid in kg/m3 ; Rv

232

is the (specific) gas constant of liquid vapour (Rv = 461.5 J/kgK for water vapour) and Tt is the

233

absolute temperature in K. Accordingly, the dimensionless suction H ∗ can be related to the relative

234

humidity h by H∗ =

−ρl Rv Tt R ln h . γ

(16)

235

The relationship between dimensionless suction H ∗ and relative humidity h is illustrated in Fig. 10

236

of [37].

237

The large dataset of numerical solutions to the Young-Laplace equation described in Section 2.3

238

can also be used to study properties of capillary bridges under humidity control, by employing Eq.(16).

239

In addition, the effect of temperature Tt on the dimensionless suction H ∗ at a prescribed humidity

240

level can also be investigated. By using Eq.(16), the dimensionless suction H ∗ present in closed-form

241

expressions in the following Sections can be converted to terms involving the relative humidity h.

242

3. Analytical expressions for rupture distances and capillary forces

243

In this Section, capillary bridges between equal-sized spherical particles under suction control are

244

considered. The meridional profile of the capillary bridge is represented by part of a circle, i.e.

245

the toroidal approximation is employed [3, 37, 38, 40, 41]. Based on the toroidal approximation,

246

analytical expressions for the rupture distance and the capillary force are derived in Sections 3.2 and

247

3.3, respectively. The accuracy of these expressions is demonstrated by comparisons with numerical

248

solutions of the Young-Laplace equation.

11

249

3.1. Toroidal approximation

250

In the toroidal approximation, part of a circle is used to represent the meridional profile of the

251

capillary bridge, see for instance [3, 16, 37, 38, 40, 41]. In earlier developments [3, 16, 37, 38, 40, 41],

252

capillary forces were expressed as functions of separation distance and half-filling angle. It is therefore

253

valuable to further develop this approximation in order to directly calculate rupture distances and

254

capillary forces for given values of suction and separation distance.

255

Numerical solutions to the Young-Laplace equation can be represented in dimensionless form using

256

the scaling method suggested in [14, 15, 18]. For equal-sized spherical particles, the capillary bridge

257

is symmetric with respect to the neck located at X = 0, thus Xcs = −XcL and Ycs = YcL . Therefore,

258

the scaling method is expressed as

P =

|Xci | Yci

Q=

Y0 Yci

.

(17)

259

As discussed in [14, 15, 18], the dimensionless coordinate P represents the scaled axial position of the

260

contact circle and the dimensionless coordinate Q represents the scaled neck radius.

261

262

The ’closure’ relationship P = F (Q) between the scaled coordinates P and Q in Eq.(17) according to the toroidal approximation is given by [15, 18]

√ P = ( T 2 + 1 + T )(1 − Q) .

(18)

263

Under volume control (so for fixed value of the capillary bridge volume), it has been found that the

264

toroidal approximation Eq.(18) is accurate for small separation distances, but not for large separation

265

distances [14, 15, 18, 28]. For suction control, Fig. 2 shows the comparison between the relationship

266

between the scaled coordinates P and Q from numerical solutions of the Young-Laplace equation and

267

according to the toroidal approximation, for contact angles θ of 0◦ and 40◦ . This Figure shows that,

268

under suction control, the toroidal approximation well represents numerical solutions to the Young-

269

Laplace equation, both for small and large separation distances, for the stable branch of solutions.

270

Under volume control, the critical point (Q∗ , P ∗ ) (separating the stable and the unstable branch)

12

271

of the relation P = F (Q) corresponds to the maximum value of P [14, 15, 18]. However, the P − Q

272

relation for suction control does not have such a property. Therefore, it is not necessary to employ

273

the ellipse approximation [14, 15, 18] that has a peak value of P for capillary bridges under suction

274

control.

275

For capillary bridges between equal-sized spherical particles, Rs = RL = R, and hence

S∗ =

S R

=

Ss R

=

SL R

Yc∗ =

Ycs R

=

YcL R

Y0∗ =

Y0 R

H ∗ = H · R λ∗ =

λ R

V∗ =

V R3

.

(19)

276

These dimensionless variables have been adopted in the following derivations for the rupture distance

277

and the capillary force.

278

3.2. Rupture distance

279

Here it is shown that the critical point (Q∗ , P ∗ ) (indicated in Fig. 2) that separates the stable

280

and unstable branch of solutions can be well predicted by the toroidal approximation, and hence the

281

rupture distance under suction control can be estimated by the toroidal approximation.

282

283

For capillary bridges between equal-sized particles with small volumes, it follows from Eqs.(5) and (17) that P can be expressed in terms of separation distance S ∗ and contact radius Yc∗

P =

S∗ 1 + Yc∗ . ∗ Yc 2

(20)

284

Substituting Eqs.(17) and (20) into Eq.(18), it follows that the separation distance S ∗ can be expressed

285

in terms of the contact radius Yc∗ and the neck radius Y0∗

p p 1 S ∗ = ( T 2 + 1 + T )Yc∗ − Yc∗2 − ( T 2 + 1 + T )Y0∗ . 2

(21)

286

Using Eq.(10), the dimensionless neck radius Y0∗ can be expressed in terms of the dimensionless contact

287

radius Yc∗ and the dimensionless suction H ∗

Y0∗ = − H1∗ +

q

2Yc∗ (T +Yc∗ ) √ H ∗ T 2 +1

13

+

1 H ∗2

+ Yc∗2 .

(22)

288

Substituting Eq.(22) into Eq.(21), the separation distance S ∗ is obtained as a function of contact

289

radius Yc∗ and dimensionless suction H ∗ :

 √ √ S ∗ = ( T 2 + 1 + T )Yc∗ − 12 Yc∗2 − ( T 2 + 1 + T ) −

290

1 H∗

+

q

2Yc∗ (T +Yc∗ ) √ H ∗ T 2 +1

+

1 H ∗2

 + Yc∗2 .

(23)

Eq.(23) can be rewritten as a quartic equation in Yc∗

K4 Yc∗4 + K3 Yc∗3 + K2 Yc∗2 + K1 Yc∗ + K0 = 0 ,

(24)

291

where the coefficients K4 , K3 , K2 , K1 and K0 (functions of the dimensionless suction H ∗ , the dimen-

292

sionless separation distance S ∗ and the contact angle θ) are given in Appendix B. For given suction

293

and separation distance, Eq.(24) has two real and positive solutions for Yc∗ that are physically rele-

294

∗ , these two solutions are the same. This occurs when vant. At the critical separation distance 2Scrit

295

the discriminant of the quartic equation (24) equals zero [50]. Based on this property, the rupture

296

∗ distance 2Scrit can be expressed as

∗ 2Scrit

  h 27 27 2 1 81 2  1 i1/3 = ∗√ · 1− + T+ T . H 32 16 64 H∗ T2 + 1

(25)

297

The detailed derivation of Eq.(25) is given in Appendix B. With relationship Eq.(25), the critical

298

∗ separation distance 2Scrit of capillary bridges can be determined for specified suction values H ∗ (and

299

contact angle θ).

300

∗ extracted from the large database Fig. 3 shows comparisons of the critical separation distance 2Scrit

301

of numerical solutions of the Young-Laplace equation and that predicted by the rupture criterion of

302

Eq.(25). Note that the product of rupture distance and dimensionless suction is compared, as that

303

more clearly shows its accuracy (as the rupture distance itself is proportional to 1/H ∗ , which varies by

304

multiple orders of magnitude). The analytical rupture criterion agrees very well with the numerical

305

solutions of the Young-Laplace equation, for contact angles 0◦ ≤ θ ≤ 40◦ . When H ∗ ∼ = 104 , the

306

deviations are smaller than 2% for θ = 0◦ , 1% for θ = 20◦ and 0.5% for θ = 40◦ . For larger values of

14

307

H ∗ , the deviations are much smaller, not more than 0.5% when H ∗ ≥ 106 . Note that this analytical

308

result does not involve any calibrated fitting coefficient.

309

3.3. Capillary force

310

The dimensionless capillary force λ∗ can be determined, for given values of mean curvature H ∗

311

and neck radius Y0∗ , through the second expression of Eq.(10). For capillary bridges under suction

312

control, the neck radius Y0∗ is a function of the dimensionless suction H ∗ and the dimensionless

313

separation distance S ∗ , i.e. Y0∗ = Y0∗ (H ∗ , S ∗ ). Therefore, for given values of dimensionless suction

314

H ∗ and separation distance S ∗ , the capillary force between two equal-sized spherical particles can be

315

determined by 1 λ∗ = Y0∗ (H ∗ , S ∗ ) + H ∗ Y0∗2 (H ∗ , S ∗ ) , 2

316

(26)

where Y0∗ (H ∗ , S ∗ ) is analytically derived by using the toroidal approximation in the following manner.

317

For 0◦ ≤ θ ≤ 40◦ , it is observed from the numerical solutions of the Young-Laplace equation that

318

0 < K4 /H ∗ ≤ 1 and −4 ≤ K3 /H ∗ < −2. Since the contact radius is small, Yc∗  1, the quartic

319

term in Yc∗ in Eq.(24) can be neglected in comparison to the cubic term. Accordingly, Eq.(24) can be

320

reduced to

K3 Yc∗3 + K2 Yc∗2 + K1 Yc∗ + K0 = 0 .

(27)

321

The discriminant of Eq.(27) was found to be negative by numerical verification. Thus, Eq.(27) has

322

three unequal real solutions (see also [50]). These solutions can be determined by the trigonometric

323

formula [50], of which the physically-relevant one is given by √

Yc∗ = −

2

−3K1 K3 +K22 3K3

cos φ −

K2 3K3

where cos 3φ =

27K0 K32 −9K1 K2 K3 +2K23

√

2

−(3K1 K3 −K22 )3

,

(28)

324

and K3 , K2 , K1 and K0 are given in Eq.(34). The accuracy of this solution for Yc∗ has been verified by

325

comparisons of numerical solutions of Eq.(24) and numerical solutions of the Young-Laplace equation.

326

Substitution of this solution for Yc∗ into Eq.(22) gives an expression for the dimensionless neck radius

15

327

Y0∗ as a function of H ∗ and S ∗ , i.e. Y0∗ (H ∗ , S ∗ ).

328

Given the dimensionless suction H ∗ and the dimensionless separation distance S ∗ , the procedure for

329

determining the capillary force between equal-sized spherical particles using the toroidal approximation

330

∗ consists of: (1) checking that the separation distance S ∗ is smaller than the critical separation Scrit

331

determined by Eq.(25); (2) calculating the dimensionless contact radius Yc∗ through Eq.(28) from the

332

given dimensionless suction H ∗ and dimensionless separation distance S ∗ ; (3) substituting Yc∗ and H ∗

333

into Eq.(22) to obtain the dimensionless neck radius Y0∗ ; (4) with the obtained value of Y0∗ , determining

334

λ∗ from Eq.(26) and the capillary force Fcap from Eq.(11).

335

Fig. 4 shows the comparison of the capillary forces from the numerical solutions of the Young-

336

Laplace equation and those predicted by the toroidal approximation. The analytical theory agrees very

337

well with the numerical solutions of the Young-Laplace equation when the mean curvature H ∗ ≥ 104

338

and the contact angle 0◦ ≤ θ ≤ 40◦ .

339

4. Curve-fitting expressions for rupture distances and capillary forces

340

The toroidal approximation described in Section 3 is limited to the dimensionless suction H ∗ ≥ 104 .

341

In this Section, computationally-fast closed-form expressions for the rupture distance and the capil-

342

lary force of capillary bridges between equal-sized spherical particles under suction control have been

343

obtained by curve-fitting to the large dataset of numerical solutions of the Young-Laplace equation.

344

These expressions are applicable over a wider range of dimensionless suction than the analytical ex-

345

pressions that have been formulated in Section 3.

346

4.1. Rupture distance

347

∗ for given values of Gras et al. [30] obtained a closed-form expression for the rupture distance 2Scrit

348

dimensionless suction H ∗ , by curve-fitting to their limited dataset of numerical solutions of the Young-

349

Laplace equation. As described in [30], such an expression is applicable only for the dimensionless

350

suction in the range of 0 ≤ H ∗ ≤ 15. This range of dimensionless suction corresponds to large capillary

351

bridge volumes V ∗1/3 ≥ 0.1 and is very limited, compared to the whole range of dimensionless suction

352

indicated in Fig. 5. Therefore, it is necessary to formulate a closed-form expression that is capable of 16

353

predicting rupture distances for a wider range of dimensionless suction than expression (25) and that

354

in [30].

355

By curve-fitting to a large dataset of numerical solutions of the Young-Laplace equation that is

356

obtained by a high-resolution integration method in Section 2, a more accurate closed-form expression

357

∗ for determining the rupture distance 2Scrit has been obtained, given by

∗ 2Scrit =

358

  2 1 ∗−1/3 ∗−2/3 √ · 1 + c H + c H , 1 2 H∗ T 2 + 1

(29)

where the coefficients c1 and c2 are dependent upon the contact angle θ, given in Appendix B.

359

As shown in Fig. 5, the obtained closed-form expression Eq.(29) agrees very well with the numerical

360

solutions of the Young-Laplace equation for the very large range of dimensionless suction considered

361

with H ∗ ≥ 1. For H ∗ ≥ 102 , the deviations between the fit and the numerical solutions are smaller

362

than 0.5%, much better than that of Gras et al. [30] which gives deviations up to 20%. For smaller

363

values of the dimensionless suction, the deviations are less than 2.5%.

364

4.2. Capillary force

365

The analytical expression for the capillary force in Eq.(26) is not accurate for 1 ≤ H ∗ ≤ 104 (data

366

not shown). Therefore, it is desirable to formulate a closed-form expression for the capillary force

367

for a wider range of dimensionless suction, H ∗ ≥ 1. Note that the dimensionless capillary force λ∗ is

368

determined by the dimensionless neck radius Y0∗ .

369

Based on observations on numerical solutions of the Young-Laplace equation, the neck radius

370

Y0∗ for capillary bridges between equal-sized spherical particles is dependent on the separation dis-

371

tance S ∗ , the dimensionless suction H ∗ and the contact angle θ, i.e. Y0∗ = Y0∗ (S ∗ , H ∗ ; θ). Given

372

the number of dependent variables, the steps for curve-fitting to the neck radius from numerical

373

solutions of the Young-Laplace equation are: (1) curve-fitting the neck radius at zero separation dis-

374

∗ ∗ ≥ 1 and for contact angles θ of 0◦ , 20◦ and 40◦ , i.e. tance Y0{S ∗ =0} for dimensionless suction H

375

∗ ∗ ∗ ∗ Y0{S ∗ =0} = Y0 (S = 0, H ; θ); (2) curve-fitting the neck radius normalized by its value at zero sepa-

376

∗ ∗ ◦ ◦ ration distance Y0∗ /Y0{S ∗ =0} for dimensionless suction H ≥ 1 and for contact angles θ of 0 , 20 and

17

Y0∗ ∗ , H ∗ ; θ); (S ∗ /Scrit ∗ Y0{S ∗ =0}

377

∗ 40◦ , i.e. Y0∗ /Y0{S ∗ =0} =

378

∗ ∗ Y0∗ /Y0{S ∗ =0} and Y0{S ∗ =0} .

379

380

and (3) obtaining the neck radius Y0∗ by multiplying

By curve-fitting to the large dataset extracted from numerical solutions of the Young-Laplace ∗ equation, an expression for the neck radius at zero separation distance, Y0{S ∗ =0} , has been obtained

r ∗ Y0{S ∗ =0}

=

  i 2 h p2 √ · p exp − + p 1 3 , H∗ H∗

(30)

381

Here the coefficients p1 , p2 and p3 are functions of T (T = tan θ), given in Appendix C. The accuracy

382

of this expression is demonstrated in Fig. 6.

383

By curve-fitting to the large dataset, the neck radius normalized by the neck radius at zero sepa-

384

∗ ∗ ration distance Y0∗ /Y0{S ∗ =0} for the selected dimensionless suction H ≥ 1 and for contact angles θ of

385

∗ 0◦ , 20◦ and 40◦ , the ratio of Y0∗ /Y0{S ∗ =0} is obtained as

Y0∗ ∗ Y0{S ∗ =0}

386

s  S∗ S∗ 1− ∗ −1 +B ∗ , =1+A Scrit Scrit

(31)

where A and B are functions of H ∗ and θ, given by

A = a1 ln2 (H ∗ ) + a2 ln(H ∗ ) + a3 . 2

B = b1 ln

(H ∗ )

+ b2

ln(H ∗ )

(32)

+ b3

387

Here the coefficients ai and bi only depend on the contact angle θ, as given in Appendix C. The

388

accuracy of this expression is demonstrated in Fig. 7.

389

Fig. 8 shows the comparison for the dimensionless capillary forces λ∗ between those from numerical

390

solutions to the Young-Laplace equation and the predictions by Eqs.(26), (30) and (31), for the

391

dimensionless suction H ∗ of 1, 25, 103 and 105 and the contact angle θ of 0◦ and 40◦ . The deviations

392

between the predictions and the numerical solutions are smaller than 1% for both θ = 0◦ and θ = 40◦ .

393

With these curve-fitting expressions, the steps for determining the capillary force of capillary

394

bridges between equal-sized spherical particles for given dimensionless suction H ∗ ≥ 1 and separation

395

distance S ∗ are: 18

396

397

398

399

∗ determined • checking that the separation distance S ∗ is smaller than the critical separation Scrit

by Eq.(29); ∗ • calculating the neck radius at zero separation distance Y0{S ∗ =0} by Eq.(30) from the given

dimensionless suction H ∗ ;

400

∗ ∗ • substituting Y0{S ∗ =0} into Eq.(31) to obtain Y0 ;

401

• with the obtained values of Y0∗ and specified H ∗ , determining λ∗ from Eq.(26) and the capillary

402

403

force Fcap from Eq.(11).

5. Capillary bridges between unequal-sized spherical particles under suction control

404

In this Section, the expressions for capillary force and rupture distance for the case of equal-sized

405

spherical particles under suction control are developed to the case of unequal-sized spherical particles.

406

5.1. Considering size effect by Derjaguin radius

407

The effect of size ratio on capillary force and rupture distance has been studied in [17, 18, 28, 51]

408

(under volume control). It has been analytically and numerically proven in [18] that, for small capillary

409

bridge volumes, the expressions for capillary force and rupture distance for equal-sized particles can be

410

employed for unequal-sized particles by using the Derjaguin radius (also called the harmonic radius)

411

that is defined in Eq.(2). The accuracy of the Derjaguin radius R has been demonstrated for small

412

capillary bridge volumes (V ∗1/3 ≤ 0.1). Such a volume range corresponds to the dimensionless suction

413

H ∗ ' 102 .

414

415

Using the Derjaguin radius R, the separation distance 2S, the capillary bridge volume V , the capillary force λ and the mean curvature H can be made dimensionless

2S ∗ =

2S R

V∗ =

V 3 R

λ∗ =

λ R

H∗ = H · R .

(33)

416

Substituting the dimensionless S ∗ , λ∗ and H ∗ into the expressions for rupture distances (Eqs.(25)

417

and (29)) and capillary forces (Eqs.(26), (30) and (31)) for equal-sized spherical particles, rupture

418

distances and capillary forces for unequal-sized spherical particles can be determined. The accuracy

419

of these predictions is investigated in the following subsections. 19

420

5.2. Rupture distance

421

∗ Fig. 9a shows the comparison of rupture distances 2Scrit extracted from numerical solutions to

422

the Young-Laplace equation and the predictions according to the toroidal approximation, for particle

423

size ratios of RL /Rs =1, 2, 8 and 128, contact angle of 40◦ and dimensionless suction H ∗ ≥ 104 . As

424

expected, the toroidal approximation agrees very well with numerical solutions to the Young-Laplace

425

equation, with deviations smaller than 1%.

426

∗ Fig. 9b compares the rupture distance 2Scrit extracted from numerical solutions of the Young-

427

Laplace equation with the predictions of the closed-form expression in Eq.(29), for particle size ratios

428

of RL /Rs = 1, 2, 8 and 128 and contact angle of 40◦ . The range of dimensionless suction is much

429

larger than that for the toroidal approximation in Fig. 9a. The closed-form expression agrees well

430

with the numerical solutions to the Young-Laplace equation, with deviations smaller than 1% for

431

H ∗ ≥ 102 . For 1 ≤ H ∗ ≤ 102 (corresponding to capillary bridge volume V ∗1/3 ≥ 0.1), the use of the

432

Derjaguin radius gives larger deviations [18]. Hence the deviations between the numerical solutions of

433

the Young-Laplace equation and the predictions of Eq.(29) are larger.

434

5.3. Capillary force

435

Fig. 10a shows the comparison between capillary forces λ∗ extracted from numerical solutions

436

to the Young-Laplace equation and the predictions of the toroidal approximation, for particle size

437

ratios of RL /Rs = 1, 2, 8 and 128, contact angles of θ = 0◦ and 40◦ when the dimensionless suction

438

H∗ ∼ = 105 . As expected, the toroidal approximation agrees very well with numerical solutions to the

439

Young-Laplace equation, with deviations smaller than 1%.

440

Fig. 10b compares the capillary force λ∗ extracted from numerical solutions to the Young-Laplace

441

equation with the predictions of the closed-form, curve fitting expression in Eqs.(26), (30) and (31),

442

for particle size ratios of RL /Rs = 1, 2, 8 and 128 and contact angles of θ = 0◦ and 40◦ when the

443

dimensionless suction H ∗ ∼ = 102 . It is shown that the predictions are in good agreement with the

444

numerical solutions of the Young-Laplace equation.

445

From these analyses, it is concluded that the closed-form expressions are effective for capillary

446

bridges between unequal-sized particles by using the Derjaguin radius, when the dimensionless suction 20

447

H ∗ ≥ 102 (corresponding to V ∗1/3 ≤ 0.1, and thus the Derjaguin radius is valid [18]). The procedure

448

for determining rupture distance and capillary force of capillary bridges between unequal-sized particles

449

is then the same as that for equal-sized particles described in Sections 3 and 4.

450

6. Conclusions

451

This study focuses on capillary forces and rupture distances for capillary bridges between spherical

452

particles under suction control, using three complementary different approaches: (1) analyses of nu-

453

merical solutions of the governing Young-Laplace equation using a high-resolution integration method,

454

(2) an analytical theory based on the toroidal approximation, that does not require any calibrated

455

parameters and (3) closed-form expressions with coefficients determined by curve fitting to the large

456

dataset of numerical solutions of the Young-Laplace equation. Specific contributions of this study are:

457

• the effectiveness of the surface free-energy criterion has been demonstrated, under suction control

458

and volume control, for identifying the physically-relevant solutions to the governing Young-

459

Laplace equation.

460

• based on the toroidal approximation, an analytical theory has been developed for the rupture

461

distance and the capillary force (given values of the suction), that does not involve any calibrated

462

parameters. The analytical expressions are valid when the dimensionless suction H ∗ ≥ 104 and

463

the contact angle θ ≤ 40◦ .

464

• closed-form expressions have been formulated for the rupture distance and the capillary force

465

(given values of the suction) by curve-fitting to the large dataset of numerical solutions of the

466

Young-Laplace equation. These expressions are in good agreement with the numerical dataset

467

when the dimensionless suction H ∗ ≥ 1 and the contact angle θ ≤ 40◦ .

468

• the influence of particle size ratio on the rupture distance and the capillary force has been

469

considered by using the Derjaguin radius that is valid when the dimensionless suction H ∗ ≥ 102

470

(corresponding to small capillary bridge volumes V ∗1/3 ≤ 0.1).

471

The numerical solutions to the Young-Laplace equation and the closed-form expressions for cap-

472

illary forces and rupture distances have not been compared with the experimental results in [30] for

21

473

capillary bridges between spherical particles under suction control, since the values of the dimension-

474

less suction H ∗ in [30] are very small (in the range of 0 ≤ H ∗ ≤ 15) and the contact angle θ is not

475

constant in these experiments when the separation between the spherical particles is varied.

476

For future studies, it is recommended to implement the developed closed-form expressions into

477

Discrete Element Method simulation software or into micromechanics-based constitutive relations (see

478

for instance [52–54]) for granular materials in order to investigate the effect of capillary cohesion under

479

suction control. In addition, considering the different behaviour of capillary bridges under volume and

480

suction control, the mechanical response of wet granular assemblies is highly loading-path dependent

481

[55]. Therefore, it is important to considering such properties in effective stress tensors [56–58] and

482

constitutive relations [52, 53, 59] for wet granular materials.

483

Acknowledgements

484

The authors acknowledge support from the international research network GDRI-GeoMech (Multi-

485

Physics and Multi-Scale Couplings in Geo-Environmental Mechanics). The research described was

486

funded in part by the European Union’s Horizon 2020 research and innovation program under the

487

Marie Sklodowska-Curie grant agreement No. 832405 (ICARUS).

488

Appendix A: stable branches of numerical solutions to the Young-Laplace equation

489

As discussed in Section 2.4, branches of (multiple) solutions of the Young-Laplace equation may

490

exist under volume as well as suction control. To identify the physically-relevant stable branch of

491

solutions, criteria have been proposed based on the minimum of the energy in Eq.(13) for volume

492

control and on the minimum of the energy in Eq.(14) for suction control. In this Appendix the

493

appropriateness is investigated of (”opposite”) criteria based on the energy in Eq.(13) for suction

494

control and on the energy in Eq.(14) for volume control. To this end, the large dataset of numerical

495

solutions to the Young-Laplace equation (see Section 2.3) is employed.

496

Under volume control (so with fixed capillary bridge volume V ), the dependence of the capillary

497

force on the separation S is shown in Fig. 11a. Two branches of solutions exist: the capillary forces

22

498

larger than the capillary force at critical separation are shown by (many) black dots, whereas the

499

capillary forces smaller than the capillary force at critical separation are indicated by red triangles.

500

The energies Esurf and Esurf+∆pV in Eqs.(13) and (14) are shown in Figs. 11b and 11c, respectively.

501

According to [16, 47, 48], the stable branch under volume control corresponds to a minimum in the

502

energy Esurf in Eq.(13). For this stable branch, the energy Esurf+∆pV in Eq.(14) is larger than for the

503

unstable branch. Therefore, the criterion of a minimum of the energy Esurf+∆pV in Eq.(14) (proposed

504

in [22] for suction control) can not be used to identify the stable branch of solutions under volume

505

control.

506

Under suction control (so with fixed suction H), the dependence of the capillary force on the

507

separation S is shown in Fig. 12a. The two branches of solutions exist: the capillary forces larger

508

than the capillary force at critical separation are shown by (many) black dots, whereas the capillary

509

forces smaller than the capillary force at critical separation are indicated by red triangles. The energies

510

Esurf and Esurf+∆pV in Eqs. (13) and (14) are shown in Figs. 12b and 12c, respectively. According to

511

[22], the stable branch under suction control corresponds to a minimum in the energy Esurf+∆pV in

512

Eq.(14). For this stable branch, the energy Esurf in Eq.(14) is smaller than for the unstable branch.

513

Therefore, the criterion of a minimum of the energy Esurf in Eq.(13) (proposed in [22] for volume

514

control) can also be used to identify the stable branch of solutions under volume control. Note that

515

under suction control the stable branch has larger volume and surface area than the unstable branch

516

of solutions [22].

517

518

Similar analyses have been performed (data not shown) for other values of the contact angle θ, particle size ratio RL /Rs , capillary bridge volume V and suction H, yielding the same conclusions.

519

From these analyses, it is concluded that the criterion of a minimum of the surface free energy Esurf

520

defined in Eq.(13) is appropriate to identify the stable branch of the solutions to the Young-Laplace

521

equation, under both volume and suction control conditions.

23

522

Appendix B: rupture distance for given suction value H ∗

523

In Eq.(24) for the dimensionless contact radius Yc∗ , the coefficients K4 , K3 , K2 , K1 and K0 are

524

dependent upon dimensionless suction H ∗ , dimensionless separation distance S ∗ and contact angle θ

525

(via T defined in Eq.(7)). These coefficients are given by

K4

 √  √ = H ∗ 2 T 2 + 1T 2 − 2T 3 + T 2 + 1 − 2T

K3

  √ = −4H ∗ T 2 − T 2 + 1T + 1

K2

=



8H ∗ S ∗ T 2 + 4H ∗ S ∗ + 4T

K1

=



8H ∗ S ∗ T + 8

K0

= −4S ∗

h

√

T2 + 1 +

√



T 2 + 1 − 8H ∗ S ∗ T 3 − 8H ∗ S ∗ T − 4T 2 − 12 .

(34)

 − 8T 2 − 8 S ∗ H ∗ − 8T

√

i
− 2H ∗ T 2 − H ∗ S ∗ − 2T T 2 + 1 + 2 T 2 + 1 H ∗S∗T + 1

526

The quartic equation (24) has four real solutions, of which two are (potentially) physically relevant.

527

At the critical separation distance Scrit , there is only a single solution, and hence the discriminant ∆

528

of the quartic equation (24) (see for example [50]) becomes zero. For values of the separation distance

529

larger than the critical separation, no real solution with physical meaning exists.

530

Upon expansion of the equation ∆ = 0 in the small parameter 1/H ∗ , it follows after some algebra

531

(performed with the symbolic mathematics program Maple) that the discriminant is zero (up to second

532

order in the small parameter 1/H ∗ ) if

∗ Scrit

533

=

1 √ α T 2 +1 H ∗

α=1−

h

β(T ) H∗

i1

3

,

(35)

where β is a (complicated) function of T , whose third-order Taylor expansion in β is given by

β=

27 27 81 + T + T2 . 32 16 64

(36)

534

In the expression for the rupture distance obtained from curve-fitting to the large dataset of

535

numerical solutions of the Young-Laplace equation, Eq.(29), the coefficients ci only depend on the

536

contact angle θ. The resulting expressions for these coefficients are given by

24

c1 = −0.236T 2 − 0.120T − 1.358 ,

(37)

c2 = 0.239T 2 + 0.111T + 0.458 537

in which the calibrated coefficients have been incorporated; T = tan θ, where θ is the contact angle in

538

radians.

539

Appendix C: constants in the closed-form expression for neck radius

540

∗ In the expression for the neck radius at zero separation distance Y0{S ∗ =0} obtained from curve-

541

fitting to the large dataset of numerical solutions of the Young-Laplace equation, Eq.(30), the coeffi-

542

cients pi only depend on the contact angle θ. The resulting expressions for the functions pi (expressed

543

concisely as a matrix-vector multiplication) are given by  







2 p1   −0.151 −0.015 0.850  T            p  =  0.302  T  , 0.038 1.217  2              p3 0.012 −0.017 0.152 1

544

(38)

where the calibrated parameters are incorporated, and T = tan θ as defined in Eq.(7).

545

In the same manner, in the expression for the ratio of neck radius to its value at zero separation

546

∗ distance Y0∗ /Y0{S ∗ =0} obtained by curve-fitting to the large dataset of numerical solutions of the

547

Young-Laplace equation, Eq.(31), the coefficients ai and bi present in Eq.(32) only depend on the

548

contact angle θ (also through T = tan θ). The resulting expressions for the functions ai (expressed

549

concisely as a matrix-vector multiplication) are given by

 





a1   −0.0015 0.0001 0.0019       a  =  0.0347 −0.0142 0.0021  2        −0.1165 0.0311 0.6087 a3 550

and for the functions bi (in terms of T ) are

25

 2

 T    T    1

   ,   

(39)

 



 T2

0.0026   b1   −0.0016 0.0006         b  =  0.0330 −0.0201 −0.0365   T  2           b3 −0.0949 0.0674 0.0088 1 551

References

552

References

553

554

    .   

(40)

[1] H. C. Hamaker, The London–van der Waals attraction between spherical particles, Physica 4 (10) (1937) 1058–1072.

555

[2] I. Aranson, D. Blair, V. Kalatsky, G. Crabtree, W.-K. Kwok, V. Vinokur, U. Welp, Electrostat-

556

ically driven granular media: phase transitions and coarsening, Physical Review Letters 84 (15)

557

(2000) 3306–3309.

558

559

[3] R. Fisher, On the capillary forces in an ideal soil; correction of formulae given by W.B. Haines, The Journal of Agricultural Science 16 (3) (1926) 492–505.

560

[4] F. Souli´e, F. Cherblanc, M. El Youssoufi, C. Saix, Influence of liquid bridges on the mechanical

561

behaviour of polydisperse granular materials, International Journal for Numerical and Analytical

562

Methods in Geomechanics 30 (3) (2006) 213–228.

563

564

565

566

567

568

[5] V. Richefeu, M. El Youssoufi, E. Az´ema, F. Radja¨ı, Force transmission in dry and wet granular media, Powder Technology 190 (1–2) (2009) 258–263. [6] L. Scholt`es, B. Chareyre, F. Nicot, F. Darve, Micromechanics of granular materials with capillary effects, International Journal of Engineering Science 47 (1) (2009) 64–75. [7] J.-Y. Delenne, V. Richefeu, F. Radja¨ı, Liquid clustering and capillary pressure in granular media, Journal of Fluid Mechanics 762 (R5) (2015) 1–10.

569

[8] M. Pouragha, R. Wan, J. Duriez, N. H. Sultan, Statistical analysis of stress transmission in wet

570

granular materials, International Journal for Numerical and Analytical Methods in Geomechanics

571

42 (16) (2018) 1935–1956. 26

572

573

574

575

576

577

578

579

580

581

582

583

[9] M. E. D. Urso, C. J. Lawrence, M. J. Adams, Pendular, funicular, and capillary bridges: Results for two dimensions, Journal of Colloid and Interface Science 220 (1) (1999) 42–56. [10] J.-P. Wang, E. Gallo, B. Fran¸cois, F. Gabrieli, P. Lambert, Capillary force and rupture of funicular liquid bridges between three spherical bodies, Powder Technology 305 (2017) 89–98. [11] K. Murase, T. Mochida, Y. Sagawa, H. Sugama, Estimation on the strength of a liquid bridge adhered to three spheres, Advanced Powder Technology 19 (2008) 349–367. [12] D. Lievano, S. Velankar, J. J. McCarthy, The rupture force of liquid bridges in two and three particle systems, Powder Technology 313 (2017) 18–26. [13] M. J. Adams, S. A. Johnson, J. P. Seville, C. D. Willett, Mapping the influence of gravity on pendular liquid bridges between rigid spheres, Langmuir 18 (16) (2002) 6180–6184. [14] N. P. Kruyt, O. Millet, An analytical theory for the capillary bridge force between spheres, Journal of Fluid Mechanics 812 (2017) 129–151.

584

[15] C.-F. Zhao, N. P. Kruyt, O. Millet, Capillary bridge force between non-perfectly wettable spherical

585

particles: An analytical theory for the pendular regime, Powder Technology 339 (2018) 827–837.

586

[16] G. Lian, C. Thornton, M. J. Adams, A theoretical study of the liquid bridge forces between two

587

588

589

590

591

592

593

594

rigid spherical bodies, Journal of Colloid and Interface Science 161 (1) (1993) 138–147. [17] G. Lian, J. Seville, The capillary bridge between two spheres: new closed-form equations in a two century old problem, Advances in Colloid and Interface Science 227 (2016) 53–62. [18] C.-F. Zhao, N. P. Kruyt, O. Millet, Capillary bridges between unequal-sized spherical particles: Rupture distances and capillary forces, Powder Technology 346 (2019) 462–476. [19] P. A. Cundall, O. D. L. Strack, A discrete numerical model for granular assemblies, G´eotechnique 29 (1) (1979) 47–65. [20] L. Scholt`es, P.-Y. Hicher, F. Nicot, B. Chareyre, F. Darve, On the capillary stress tensor in wet

27

595

granular materials, International Journal for Numerical and Analytical Methods in Geomechanics

596

33 (10) (2009) 1289–1313.

597

598

599

600

601

602

603

604

605

606

607

608

609

610

611

612

613

614

615

616

617

[21] J.-P. Wang, X. Li, H.-S. Yu, A micro–macro investigation of the capillary strengthening effect in wet granular materials, Acta Geotechnica 13 (3) (2018) 513–533. [22] J. Duriez, R. Wan, Contact angle mechanical influence in wet granular soils, Acta Geotechnica 12 (1) (2017) 67–83. [23] J. Duriez, R. Wan, Subtleties in discrete-element modelling of wet granular soils, G´eotechnique 67 (4) (2017) 365–370. [24] C. Delaunay, Sur la surface de r´evolution dont la courbure moyenne est constante, Journal de Math´ematiques Pures et Appliqu´ees (1841) 309–314. [25] M. A. Erle, D. Dyson, N. R. Morrow, Liquid bridges between cylinders, in a torus, and between spheres, AIChE Journal 17 (1) (1971) 115–121. [26] F. Orr, L. Scriven, A. P. Rivas, Pendular rings between solids: meniscus properties and capillary force, Journal of Fluid Mechanics 67 (4) (1975) 723–742. [27] G. Mason, W. Clark, Liquid bridges between spheres, Chemical Engineering Science 20 (10) (1965) 859–866. [28] C. D. Willett, M. J. Adams, S. A. Johnson, J. P. Seville, Capillary bridges between two spherical bodies, Langmuir 16 (24) (2000) 9396–9405. [29] Y. I. Rabinovich, M. S. Esayanur, B. M. Moudgil, Capillary forces between two spheres with a fixed volume liquid bridge: theory and experiment, Langmuir 21 (24) (2005) 10992–10997. [30] J.-P. Gras, J.-Y. Delenne, M. S. El Youssoufi, Study of capillary interaction between two grains: a new experimental device with suction control, Granular Matter 15 (1) (2013) 49–56. [31] H. N. G. Nguyen, O. Millet, G. Gagneux, On the capillary bridge between spherical particles

28

618

of unequal size: analytical and experimental approaches, Continuum Mechanics and Thermody-

619

namics (2018) 1–13.

620

[32] H. N. G. Nguyen, O. Millet, G. Gagneux, Exact calculation of axisymmetric capillary bridge

621

properties between two unequal-sized spherical particles, Mathematics and Mechanics of Solids

622

(2018). doi:10.1177/1081286518787842.

623

624

[33] G. Gagneux, O. Millet, An analytical framework for evaluating the cohesion effects of coalescence between capillary bridges, Granular Matter 18 (2) (2016) 16.

625

[34] B. Mielniczuk, O. Millet, G. Gagneux, M. S. El Youssoufi, Characterisation of pendular capillary

626

bridges derived from experimental data using inverse problem method, Granular Matter 20 (1)

627

(2018) 14.

628

[35] G. Gagneux, O. Millet, B. Mielniczuk, M. S. El Youssoufi, Theoretical and experimental study of

629

pendular regime in unsaturated granular media, European Journal of Environmental and Civil

630

Engineering 21 (7-8) (2017) 840–853.

631

632

633

634

[36] B. V. Derjaguin, Untersuchungen u ¨ber die Reibung und Adh¨asion, IV (in German), Colloid & Polymer Science 69 (2) (1934) 155–164. [37] F. Molenkamp, A. Nazemi, Interactions between two rough spheres, water bridge and water vapour, G´eotechnique 53 (2) (2003) 255–264.

635

[38] A. Nazemi, A. Majnooni-Heris, A mathematical model for the interactions between non-identical

636

rough spheres, liquid bridge and liquid vapor, Journal of Colloid and Interface Science 369 (1)

637

(2012) 402–410.

638

[39] J.-P. Gras, J.-Y. Delenne, F. Souli´e, M. S. El Youssoufi, DEM and experimental analysis of the

639

water retention curve in polydisperse granular media, Powder Technology 208 (2) (2011) 296–300.

640

[40] J. Lechman, N. Lu, Capillary force and water retention between two uneven-sized particles,

641

Journal of Engineering Mechanics 134 (5) (2008) 374–384.

29

642

643

[41] O. Harireche, A. Faramarzi, A. M. Alani, A toroidal approximation of capillary forces in polydisperse granular assemblies, Granular Matter 15 (5) (2013) 573–581.

644

[42] F. R. De Bisschop, W. J. Rigole, A physical model for liquid capillary bridges between adsorptive

645

solid spheres: The nodoid of Plateau, Journal of Colloid and Interface Science 88 (1) (1982)

646

117–128.

647

[43] N. Lu, W. J. Likos, Unsaturated Soil Mechanics, Wiley, 2004.

648

[44] D. G. Fredlund, H. Rahardjo, M. D. Fredlund, Unsaturated Soil Mechanics in Engineering Prac-

649

650

651

652

653

tice, John Wiley & Sons, 2012. [45] J. A. Blatz, Y.-J. Cui, L. Oldecop, Vapour equilibrium and osmotic technique for suction control, in: Laboratory and Field Testing of Unsaturated Soils, Springer, 2008, pp. 49–61. [46] M. D¨ormann, H.-J. Schmid, Distance-dependency of capillary bridges in thermodynamic equilibrium, Powder Technology 312 (2017) 175–183.

654

[47] B. J. Lowry, P. H. Steen, Capillary surfaces: stability from families of equilibria with application

655

to the liquid bridge, Proceedings of the Royal Society of London. Series A: Mathematical and

656

Physical Sciences 449 (1937) (1995) 411–439.

657

658

659

660

661

662

663

664

[48] J. Bostwick, P. Steen, Stability of constrained capillary surfaces, Annual Review of Fluid Mechanics 47 (2015) 539–568. [49] K. Galvin, A conceptually simple derivation of the Kelvin equation, Chemical Engineering Science 60 (2005) 4659–4660. [50] G. A. Korn, T. M. Korn, Mathematical Handbook for Scientists and Engineers, McGraw-Hill, 2000. [51] V. Mehrotra, K. Sastry, Pendular bond strength between unequal-sized spherical particles, Powder Technology 25 (2) (1980) 203–214.

30

665

666

[52] P.-Y. Hicher, C. S. Chang, A microstructural elastoplastic model for unsaturated granular materials, International Journal of Solids and Structures 44 (7-8) (2007) 2304–2323.

667

[53] C.-F. Zhao, Y. Salami, P.-Y. Hicher, Z.-Y. Yin, Multiscale modeling of unsaturated granular ma-

668

terials based on thermodynamic principles, Continuum Mechanics and Thermodynamics (2018)

669

1–19.

670

[54] C.-F. Zhao, Z.-Y. Yin, A. Misra, P.-Y. Hicher, Thermomechanical formulation for micromechan-

671

ical elasto-plasticity in granular materials, International Journal of Solids and Structures 138

672

(2018) 64–75.

673

[55] H. Rahardjo, O. B. Heng, L. E. Choon, Shear strength of a compacted residual soil from consol-

674

idated drained and constant water content triaxial tests, Canadian Geotechnical Journal 41 (3)

675

(2004) 421–436.

676

677

678

679

680

681

682

683

[56] A. W. Bishop, G. Blight, Some aspects of effective stress in saturated and partly saturated soils, G´eotechnique 13 (3) (1963) 177–197. [57] E. E. Alonso, J.-M. Pereira, J. Vaunat, S. Olivella, A microstructurally based effective stress for unsaturated soils, G´eotechnique 60 (12) (2010) 913–925. [58] J. Duriez, R. Wan, Stress in wet granular media with interfaces via homogenization and discrete element approaches, Journal of Engineering Mechanics 142 (12) (2016) 04016099. [59] A. Gens, M. S´ anchez, D. Sheng, On constitutive modelling of unsaturated soils, Acta Geotechnica 1 (3) (2006) 137–147.

31

684

Figures

Y YcL



RL

Y(X )

Ycs Y0



L O

X cL

s X cs

Rs

X

Ss

SL

2S = S L + S s Figure 1: Meridional geometry of a capillary bridge between two unequal-sized spheres (RL ≥ Rs ) with contact angle θ. The meridional profile is described by Y (X); the neck radius is denoted as Y0 ; the meridional coordinates of the contact circles of the small and the large particles are (Xcs , Ycs ) and (XcL , YcL ), respectively; the half-filling angles are δs and δL , respectively; the separation distances, with respect to the neck located at X = 0, for the small and the large particle are Ss and SL , respectively. The total separation distance between the two particles is 2S. (Figure from [18]).

Figure 2: Relation between scaled coordinates P and Q (defined in Eq.(17)) from numerical solutions to the YoungLaplace equation and predictions according to the toroidal approximation for dimensionless suction H ∗ ≥ 104 . The critical point (Q∗ , P ∗ ) is indicated that separates the stable and unstable branches. (a) contact angle θ = 0◦ ; (b) contact angle θ = 40◦ .

32

∗ H ∗ , for various values of the dimensionless Figure 3: Product of rupture distance and dimensionless suction, 2Scrit suction H ∗ : numerical solutions of the Young-Laplace equation and predictions by the analytical theory; contact angles θ = 0◦ , 20◦ and 40◦ .

Figure 4: Dependence of the capillary force on the separation distance, obtained from the numerical solutions of the Young-Laplace equation and predictions by the toroidal approximation, for contact angle θ of 0◦ and 40◦ . (a) dimensionless suction H ∗ = 104 (b) dimensionless suction H ∗ = 105 .

33

∗ H ∗ , for various values of dimensionless suction Figure 5: Rupture distances multiplied by dimensionless suction, 2Scrit ∗ H of numerical solutions to the Young-Laplace equation and the predictions by Eq.(29).

∗ Figure 6: Product of (H ∗ /2)1/2 and the scaled neck radius at zero separation distance, Y0{S ∗ =0} , extracted from numerical solutions of the Young-Laplace equation and predicted by the closed-form expression in Eq.(30), for contact angles of θ = 0◦ , θ = 20◦ and θ = 40◦ .

34

∗ of the neck radius Y0∗ normalized by that at zero Figure 7: Dependence on normalized separation distance S ∗ /Scrit ∗ separation distance, Y0{S ∗ =0} , extracted from numerical solutions of the Young-Laplace equation and that predicted by the closed-form expression in Eq.(31), for dimensionless suctions H ∗ = 1, 25, 103 and 105 : (a) contact angle θ = 0◦ , (b) contact angle θ = 40◦ .

Figure 8: Capillary forces extracted from numerical solutions of the Young-Laplace equation and those predicted by the closed-form expressions Eqs.(26), (30) and (31), for dimensionless suctions H ∗ = 1, 25, 103 and 105 : (a) contact angle θ = 0◦ , (b) contact angle θ = 40◦ .

35

∗ H ∗ , as extracted from numerical solutions to Figure 9: Product of rupture distance and dimensionless suction, 2Scrit the Young-Laplace equation and that predicted by the closed-form expressions for particle size ratio RL /Rs = 1, 2, 8 and 128 and contact angle θ = 40◦ : (a) toroidal approximation for dimensionless suction H ∗ ≥ 104 ; (b) predictions by Eq.(29) for dimensionless suction H ∗ ≥ 1.

Figure 10: Capillary forces from numerical solutions to the Young-Laplace equation and according to the predictions by (a) the toroidal approximation of Eq.(26) when the dimensionless suction of H ∗ ∼ = 105 and (b) curve-fitting expression ∗ ∼ 2 of Eqs.(26), (30) and (31) when the dimensionless suction of H = 10 , for particle size ratios of RL /Rs = 1, 2, 8 and 128, contact angles of θ = 0◦ and θ = 40◦ .

36

Figure 11: Numerical solutions to the Young-Laplace equation under volume control, with capillary bridge volume V∗ ∼ = 10−6 , particle size ratio RL /Rs = 1 and contact angle θ = 0◦ : (a) dimensionless capillary force λ∗ as a function of dimensionless separation distance S ∗ ; note that two branches of solutions are present: the branch where the capillary forces is larger than the capillary force at critical separation is indicated by black dots, while the branch where the capillary forces smaller than the capillary force at critical separation is indicated by red triangles; (b) dimensionless free energy Esurf+∆pV ∗ ∗ E ∗ = Esurf 2 as a function of dimensionless separation distance S ; (c) dimensionless free energy Esurf+∆pV = 2 γR γR as a function of dimensionless separation distance S ∗ .

∼ Figure 12: Numerical solutions to the Young-Laplace equation under suction control, with dimensionless suction H ∗ = 103 , particle size ratio RL /Rs = 1 and contact angle θ = 0◦ : (a) dimensionless capillary force λ∗ as a function of dimensionless separation distance S ∗ ; note that two branches of solutions are present: the branch where the capillary forces is larger than the capillary force at critical separation is indicated by black dots, while the branch where the capillary forces smaller than the capillary force at critical separation is indicated by red triangles; (b) dimensionless free energy Esurf+∆pV ∗ ∗ E ∗ = Esurf 2 as a function of dimensionless separation distance S ; (c) dimensionless free energy Esurf+∆pV = 2 γR γR as a function of dimensionless separation distance S ∗ .

37

Highlights • effectiveness of surface free-energy based criteria has been investigated • the toroidal approximation method has been developed for suction control • expressions for rupture distance and capillary force have been formulated • particle size effect has been considered by using the Derjaguin radius