Analysis of Eshelby-Cheng’s model in anisotropic porous cracked medium: An ultrasonic physical modeling approach

Journal Pre-proofs Analysis of Eshelby-Cheng’s Model in Anisotropic Porous Cracked Medium: an Ultrasonic Physical Modeling Approach Murillo J.S. Nascimento, J.J.S. de Figueiredo, C. Barros da Silva, Bruce F. Chiba PII: DOI: Reference:

S0041-624X(19)30313-0 https://doi.org/10.1016/j.ultras.2019.106037 ULTRAS 106037

To appear in:

Ultrasonics

Received Date: Revised Date: Accepted Date:

19 May 2019 12 September 2019 26 September 2019

Please cite this article as: M.J.S. Nascimento, J.J.S. de Figueiredo, C.B. da Silva, B.F. Chiba, Analysis of Eshelby-Cheng’s Model in Anisotropic Porous Cracked Medium: an Ultrasonic Physical Modeling Approach, Ultrasonics (2019), doi: https://doi.org/10.1016/j.ultras.2019.106037

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.

Analysis of Eshelby-Cheng’s Model in Anisotropic Porous Cracked Medium: an Ultrasonic Physical Modeling

1

Approach 2

Murillo J. S. Nascimento† , J. J. S. de Figueiredo∗ , C. Barros da Silva∗ and

3

Bruce F. Chiba∗

4

∗

UFPA, Faculty of Geophysics, Petrophysics and Rock Physics Laboratory - Prof. Dr.

5

Om Prakash

6

Verma, Bel´em, PA, Brazil.

7

†

National Institute of Petroleum Geophysics (INCT-GP), Brazil.

8

(October 11, 2019)

9

Running head: Analysis of the Eshelby-Cheng effective model

ABSTRACT 10

Many effective medium theories are designed to describe the macroscopic properties of

11

a medium (the rock, or reservoir in this case) in terms of the properties of its constituents

12

(the background matrix of the rock and the inclusions, for our scenario). A very well known

13

effective medium theory is the Eshelby-Cheng model, which was studied by us in previous

14

work, being tested for the case where the background medium was weakly-anisotropic and

15

porous. The analysis was done testing elastic velocities and Thomsen parameters - as a

16

function of crack density for fixed values of aspect ratio - predicted by the model with data

17

acquired from synthetic rock samples. In this work, we aim to complete the analysis of the

18

Eshelby-Cheng model capabilities when applied to rocks with porous and vertical trans-

19

versely isotropic (VTI) backgrounds, testing the model for the elastic velocities as functions

20

of aspect ratio - for fixed values of crack density - against experimental data. The data

21

used to test the model were obtained from 17 synthetic rock samples, one uncracked and

22

16 cracked, the latter divided into four groups of four samples each, each group with cracks

23

having the same aspect ratio, but with the samples having different crack densities. In these

24

samples, ultrasonic pulse transmission measurements were performed to obtain the experi1

25

mental velocities used to test the model. As was not possible to acquire data for velocity as

26

a function of aspect ratio for fixed values of crack density, we performed interpolations of

27

the experimental data to estimate these velocities. Eshelby-Cheng model effective velocities

28

and Thomsen parameters were calculated using three formulations proposed for the crack

29

porosity: one proposed by Thomsen, the second one proposed in our previous work (which

30

depends only on the crack density) and the third one proposed in this work (which depends

31

on the crack porosity and the aspect ratio, just like Thomsen’s proposal). The comparisons

32

between elastic velocities and Thomsen parameters - as function of crack aspect ratio, for

33

fixed values of crack density - predicted by the model and estimated from the data via

34

interpolation showed that the third formulation produced better fittings (lower root-mean-

35

square errors) between model and experimental data for all ranges of aspect ratio and crack

36

density.

2

INTRODUCTION 37

The analysis of the physical characteristics of cracked rocks has attracted consider-

38

able interest from the academy and the industry over the years. That happens because

39

the knowledge concerning fractures/cracks physical properties (e.g., aperture, orientation,

40

and aspect ratio) and concerning properties related to them (e.g., crack porosity and den-

41

sity) is important for hydrocarbon recovery (Crampin and Peacock, 2005; Nishizawa and

42

Yoshino, 2001) and is used for reservoir characterization (Far et al., 2013). The contribu-

43

tion to the rock’s physical properties associated to its cracks can be modeled using effective

44

medium theories, which are mathematical descriptions that relate the macroscopic (effec-

45

tive) physical properties of a medium to the physical properties and volumetric fractions

46

of its constituents. Some, very well known in the literature, effective medium theories -

47

constructed for describing a cracked rock as inclusions into a homogeneous background -

48

are the Hudson-Crampin Hudson (1981) and Eshelby-Cheng models Cheng (1993). Some,

49

more complex models, were constructed to describe more complicated kinds of media, as for

50

example: transversally cracked media with low crack density (Eshelby, 1957; Hudson, 1981)

51

and high crack density (Cheng, 1993; M. Kachanov et al., 2003; Grechka and Kachanov,

52

2006).

53

The construction of synthetic rock samples (physical modeling) followed by the per-

54

forming of ultrasonic measurements is a tool used to experimentally test the capabilities

55

and limitations of the effective medium theories already described. De Figueiredo et al.

56

(2018) investigated the application of Eshelby-Cheng’s model to cracked samples whose

57

backgrounds are VTI and porous. The model proposed by De Figueiredo et al. (2018) is an

58

extension of Eshelby-Cheng’s model (Cheng, 1993), but as it is not constructed from first

59

principles, we don’t have (from a mathematical analysis) an estimate of how much error it

60

may produce; therefore we investigate its limitations via experimental study. This investi-

61

gation was done through comparisons between experimental velocities - as a function of the

62

samples cracks densities, for fixed values of aspect ratio - and velocities predicted by the

63

effective medium model (testing two alternatives for the calculation of the crack porosity

64

used in the model). Due to technical difficulties in the experiment, it was not possible to 3

65

construct samples with different aspect ratios but same crack densities, therefore it was not

66

possible to acquire the data necessary to obtain the experimental velocities as a function of

67

aspect ratio. In this work we aim to fill this gap, by estimating velocities as function aspect

68

ratio - for fixed crack densities - via interpolation of the original data, and then to perform

69

the comparisons between data and model predictions (again testing the alternatives for the

70

model porosity) in order to complete the previously-mentioned investigation of the model.

71

Furthermore, as the velocities predicted by the model fed with a crack porosity - depen-

72

dent only on crack density (and not on crack aspect ratio) -showed some mathematical

73

inconsistencies in certain regions. For parameter-space constituted by crack density and

74

crack aspect ratio (the two main variables of the problem), we propose in this work a third

75

formulation (beyond that two shown in De Figueiredo et al. (2018) ) for the crack porosity,

76

adapted from Thomsen’s formulation.

METHODOLOGY 77

Sample Preparation

78

Seventeen porous synthetic rock samples were constructed under controlled conditions.

79

These were constructed, based on the method developed by Santos et al. (2017), from a

80

mixture between water and a matrix of 65 % sand and 35% cement. All seventeen samples

81

have a VTI background matrix but, among them, sixteen are cracked, and there is one

82

uncracked purely VTI. Figure 1 (c) shows the complete set of samples. The inclusions -

83

used to form the cracks - were arranged in planes parallel to the layering of the samples.

84

These sixteen anisotropic cracked samples were divided into four groups, each with four

85

samples and one crack aspect ratio for the samples within its group.

86

[Figure 1 about here.]

87

[Table 1 about here.]

88

To create the VTI background, we deposited our mixture layer by layer into a rectangular

89

mold. To create (approximately) penny-shaped voids inside the samples, we placed low 4

90

aspect-ratio-cylindrical styrofoam cuts along with the layers, between each layer deposition,

91

as shown in Figure 1 (a). After dried, the mixture formed a rectangular synthetic rock,

92

from it we extracted cylindrical plugs (as shown in Figure 1 (b)) using a plugging machine

93

from CoreLab instruments. To form the penny-shaped voids, a chemical leaching using

94

paint thinner was performed to the samples (immersing them in paint thinner for around

95

24 hours), dissolving the styrofoam and leaving the desired voids. Test samples were cut to

96

verify the absence of styrofoam at the cracks after the leaching. As styrofoam is a polymer

97

made by 98 % air and 2 % polystyrene (Poletto et al., 2014), after leaching a practically

98

nil amount of polystyrene remains within the fractures. The complete description of the

99

samples and a Table 2 the values of all of its physical properties (e.g., porosity, mass density,

100

crack density) can be found in De Figueiredo et al. (2018). The values of crack density in

101

the samples constructed range from 0 to 0.12 while the values of aspect ratio range from

102

0.08 to 0.52. Notice that the aspect ratios are very large (typical values for natural rocks

103

are lower 0.1), this is mainly due to difficulties in cutting very thin styrofoam slices (which

104

are used as inclusions), besides, the goal is to test the limitations of the model, therefore to

105

use high values of aspect ratio is desired.

[Table 2 about here.]

106

107

Ultrasonics measurements

108

We performed the transmission measurements using a 500 kHz S-wave transducer that

109

allows recording both P- and S-waveforms in a single trace. The sampling rate per channel

110

for all measurements was 0.01 µs. We also used a pulse-receiver 5072PR, a pre-amplifier

111

5660B from Olympus and a USB oscilloscope of 200 Ms/s from Hantek (reference DSO

112

3064). Figure 2 (a) shows the picture of the experimental setup used in this work. The

113

complete description of ultrasonic experimental setup including the frequency source spectra

114

can be found at Santos et al. (2016) and Henriques et al. (2018).

115

[Figure 2 about here.]

5

116

From the first arrival pickings of the acquired P- and S-waveforms, we calculated the P-

117

and S-wave velocities. We obtained five experimental velocities using three P-wave arrival

118

times and two S-arrival times. From these velocities, we calculated the Thomsen parameters

119

(Thomsen, 1986), which quantify the anisotropy of a VTI medium, associated with each

120

sample. We define our Z-axis as being perpendicular to the crack planes, while the Y-axis is

121

defined as the central axis of the cylinder and the X-axis is perpendicular to both of them.

122

Following these definitions, the XY plane coincides with the crack planes - as can be seen in

123

Figure 1 (a) - and the system is orthogonal and right-handed. The three P-wave velocities

124

were measured (1) in the Z direction, (2) making 45 degrees with the Z-axis in the ZX plane

125

and (3) in the X-direction. Both S-wave velocities were measured in the Y direction. Here,

126

the velocity called SH was measured for the S-wave polarized parallel to the crack planes,

127

while the velocity called SV was measured for the S-wave polarized perpendicular to the

128

crack planes. The measurement configurations for both P- and S-wave pulses can be seen

129

in Figure 2 (b). These five velocities - for each sample - are our experimental data, that

130

will be used for comparisons against the effective medium theory proposed.

131

Data interpolation

132

As there was not possible to acquire experimental data for velocity as a function of aspect

133

ratio (for fixed values of crack density), we performed a linear interpolation of the velocities

134

as a function of crack density and ”resampled” the linear function for fixed values of crack

135

density. As this process is repeated - using the same crack density values for different curves -

136

for the four crack aspect ratios studied, it produces values that approximate how the velocity

137

would change in a sample for varying aspect ratios, maintaining the crack density constant.

138

The linear interpolation was chosen, instead of a higher polynomial fitting, because the

139

velocity values had a small relative rate of change inside the interval considered, as can be

140

seen in Figure 3.

141

[Figure 3 about here.]

6

142

The same methodology was applied to the Thomsen parameter’s data, in order to pro-

143

duce an estimation of the Thomsen parameters as a function of aspect ratio, estimation

144

that would not be possible to obtain experimentally, because the Thomsen parameters are

145

calculated from the five elastic velocities measured, which in turn could be not measured

146

for varying aspect ratios with fixed crack densities, as described earlier. Figure 4 shows the

147

Thomsen parameters calculated from the measured velocities as a function of crack density

148

for fixed aspect ratios, as well as the linear regression of the data corresponding to each as-

149

pect ratio. Sampling this curves in the same crack density values - noting that each column

150

corresponds to an aspect ratio value - produces the desired estimate: Thomsen parameters

151

as a function of cracks aspect ratio with fixed crack densities.

[Figure 4 about here.]

152

THEORETICAL BACKGROUND 153

Eshelby-Cheng model

154

Based on the Pad´e expansion and Eshelby (1957)’s solution - for the elastic field of

155

inclusions into an isotropic background - Cheng (1993) developed a mathematical formula-

156

tion for describing transversely isotropic cracked media. With this expansion Cheng (1993)

157

solved the problem, present in the description of Hudson (1981), of divergence at high crack

158

ef f densities. In Eshelby-Cheng’s model, the effective VTI elastic stiffness coefficients Cij are

159

given by (0)

(1)

ef f Cij = Cij − φCij ,

(1)

(0)

160

where Cij are the isotropic elastic coefficients associated with the host material, φ is the

161

medium’s porosity and Cij is the first-order correction - which carries the information

162

about the influence of the cracks in the stiffness of the rock.

163

The first-order correction term Cij is a function of both crack density  and crack aspect

164

ratio α, which are the two parameters describing the set of cracks present in rocks, the

165

first one associated to the relative volume occupied by cracks and the second one related

(1)

(1)

7

166

to the shape of the cracks. The samples crack density is estimated by the modified Hudson

167

(1981)’s equation, given by  = N πhr2 /Vm ,

(2)

168

where N is the total number of penny-shaped inclusions, h is the aperture (thickness) of

169

each inclusion, r is the radius of each inclusion, and Vm is the model volume (volume of the

170

sample in our case). The aspect ratio of a shape is defined as the ratio between two of its

171

dimensions. For our cracks, considering its cylindrical form, the aspect ratio is calculated

172

from α = h/D.

(3)

173

where h and D = 2r are the cracks aperture and diameter, respectively. For values of

174

aspect ratio much smaller than one, the cylindrical cracks look like discs or low aspect ratio

175

penny-shaped cracks, and this is our case. For accentuated values of aspect ratio, higher

176

than one, the cylindrical cracks are considered as needle-shaped. As spheres have aspect

177

ratio equal to one, cracks that exhibit aspect ratio close to one are well approximated by

178

spheres.

179

Adaptation for anisotropic porous samples

180

In the Eshelby-Cheng model, the input parameters used to feed the model are (1) the

181

background information used to calculate Cij , i.e., the background isotropic moduli, which

182

can be obtained from measurements of background density, P-wave velocity and S-wave

183

velocity, and (2) the cracks’ information - crack aspect ratio and crack density - which are

184

used to calculate φ and, combined to the background information, are used to calculate the

185

first-order correction terms.

(0)

186

In this model, porosity is exclusively associated with the cracks. To adapt this model

187

to porous rocks De Figueiredo et al. (2018) introduced the background porosity contribu-

188

tion to the total porosity φ. De Figueiredo et al. (2018) also tested the equation 4 with

189

background VTI stiffness coefficients instead of background isotropic stiffness coefficients

190

to try to estimate the effective elastic coefficients of cracked samples with VTI background 8

191

anisotropy. So for De Figueiredo et al. (2018)’s adaptation of Eshelby-Cheng’s model we

192

can write (1)

ef f V TI Cij = Cij − φT Cij ,

(4)

193

where the total porosity is obtained by adding up the crack porosity to the background

194

porosity φo (which can be estimated from a reference (uncracked) sample using a porosimeter

195

in practical cases). In other words, the total porosity is mathematically given by

φT = φo + φc .

(5)

196

where φc is the crack porosity. Another modification proposed was in the estimation of the

197

crack density. Cheng (1993) estimates the crack density via 6, while De Figueiredo et al.

198

(2018) tests two possible formulations:

199

4 φc = πα, 3

(6)

φc = (1 − φo ).

(7)

and

200

In this work, we propose an additional formulation, which combines both of the crack

201

porosities previously shown: 4 φc = π(1 − φo )α, 3

(8)

202

which is a linear function of α and , just like the first formulation, but tries to incorporate

203

the effect of the background porosity - to be more precise, to exclude the influence of the

204

background porosity from the crack porosity. The applications of the adapted models are

205

summarized in the following workflow, which also describes the measurements needed to

206

calculate the input parameters that, in turn, are used in equation 4 in order to get the

207

effective elastic stiffness coefficients.

9

208

Anisotropic parameters and velocities

209

It is assumed that the velocities obtained via the pulse transmission method performed

210

in this work (described in ultrasonic measurements in the Methodology section) are phase

211

velocities because the transducers are wide and their diameters are comparable with sample

212

sizes (Dellinger and Vernick, 1994). As we are dealing with phase velocities, they can be

213

easily inverted to obtain the coefficients of the stiffness tensor and vice-versa (given the

214

medium’s density), both for isotropic and VTI media.For the elastic isotropic medium the P-

215

and S-wave velocities depend on two independent elastic moduli (the Lam´e parameters) and

216

the medium’s density, while for the elastic VTI medium the P-, SV- and SH-wave velocities

217

depend on five independent elastic parameters, the medium’s density and the propagation

218

angle with respect to the symmetry axis. The Thomsen parameters are three dimensionless

219

parameters introduced in the Thomsen (1986) to quantify the anisotropy (initially) on VTI

220

media. They are defined in terms of the five elastic stiffness tensor components of the

221

VTI medium. As the stiffness coefficients can be written in terms of elastic velocities, it

222

is also possible to write the Thomsen parameters as a function of five independent elastic

223

velocities. Based on the dependence between the anisotropic parameters with velocities, is

224

possible to see, that  quantifies the P-wave anisotropy, γ quantifies the S-wave anisotropy

225

and δ is related to both P- and S-waves.

[Figure 5 about here.]

226

227

Measurements performed in this work give us elastic wave velocities and densities for

228

the samples backgrounds, which are converted to background stiffness coefficients, which

229

in turn are used as input (together with crack information) in the adapted Eshelby-Cheng

230

model equation (equation 4), the equation outputs effective stiffness coefficients, which in

231

turn are converted to effective elastic velocities. As argued, the velocities can be used to

232

calculate the Thomsen parameters. That is the workflow used to produce the ”theoretical

233

elastic velocities and Thomsen parameters” proposed by De Figueiredo et al. (2018), and it

234

is well summarized in Figure 5.

10

RESULTS 235

Elastic velocities

236

The Eshelby-Cheng model, following the workflow adapted for rocks with porous VTI

237

background proposed by De Figueiredo et al. (2018), was used for the calculations of the

238

curves in Figure 6. Differently. from the previous work, the theoretical curves for the veloc-

239

ities (see Table 3 to see the experimental values) were calculated as functions of the crack

240

aspect ratio. The purple circles represent the experimental data (or our better estimate of

241

it, given our experimental limitations), that is, the results of the interpolation described in

242

the data interpolation subsection of the methodology section.

243

[Figure 6 about here.]

244

[Figure 7 about here.]

245

[Table 3 about here.]

246

The workflow shown in Figure 5 was applied for calculating the theoretical elastic ve-

247

locities, using -as estimation of the crack porosity- equation 6 to generate the continuous

248

blue curves, equation 7 to generate the dashed red curves and equation 8 to produce the

249

yellow curves made of dashes and points. It is possible to observe in Figure 6 that, in

250

general, the yellow curves follow closer to the circles (i.e., that they have a better fitting

251

with the data) than the other two curves and that the blue curves are just a little further

252

away from the points. The red curves have very different behavior from the other two and

253

are much further away from the data circles in some graphs - e.g., V p0 , V p45 and Vsv plots.

254

The fact that the velocity curves calculated using crack porosity equations 6 and 8 show

255

very similar behavior, while the velocities calculated using crack porosity equation 7 show

256

a different pattern. This can be associated to the fact that in the first pair of equations

257

the crack porosity has an explicit dependence on both crack density, and crack aspect ratio

258

(note that the expressions are basically the same function, differing only by a proportion-

11

259

ality factor of 1 − φo ), while the third crack porosity equation mentioned has an explicit

260

dependence only on crack density.

261

Thomsen parameters

262

From the theoretical velocity values (i.e., calculated with the effective medium theory)

263

the Thomsen Parameters , δ and γ were calculated, as functions of aspect ratio. As de-

264

scribed in the data interpolation subsection in the methodology section, the experimental

265

velocity values were also used to calculate the Thomsen parameters, and then these values

266

were interpolated and resampled to produce an approximation of the experimental data

267

that would be measured for the Thomsen parameters as a function of aspect ratio. The

268

comparisons between the theoretical (calculated from the theoretical velocities, using the

269

three formulations) and the experimental (or an approximation of what would be the ex-

270

perimental values) Thomsen parameters are shown in Figure 7.

271

[Figure 8 about here.]

272

In Figure 7 we observe basically the same features observed in Figure 6, that is: the

273

curves associated to the crack porosity equation 8 present the best fittings with the data,

274

while the curves associated to the crack porosity 6 shows the same pattern following the

275

trend of the circles, while being just a little further away from the data points. The curves

276

associated to the crack porosity equation 7 show a completely different behavior from the

277

other two, being much further away from the points -and not following the trend of the

278

data- at least for low aspect ratio values. Another problem concerning the red curves

279

(related to the crack porosity 7) is that its values seem not to respect the restriction values

280

of the Thomsen parameters, the curves assume negative values, then increase very fast

281

surpassing one and then come back to reasonable values. Visually, the misfit between the

282

latter mentioned curves and the circles looks greater for low aspect ratio values.

12

DISCUSSION 283

As this work tries to extend an analysis initiated in De Figueiredo et al. (2018) -which

284

studies velocities and Thomsen parameters as function of crack density- by the velocities and

285

Thomsen parameters calculated with the Eshelby-Cheng model as function of aspect ratio,

286

then it is useful to interpret the presented results and compare them to the ones shown in the

287

previous work from a perspective which treats elastic velocities and anisotropic parameters

288

as multivariable functions of aspect ratio and crack density. Therefore, we calculated and

289

presented in Figure 8 the velocity maps (we do not show the Thomsen parameters maps

290

because the analysis would be redundant) as a function of the two variables discussed here.

291

[Figure 9 about here.]

292

It is possible to observe in Figure 8 that the velocity maps produced by feeding the

293

Eshelby-Cheng model with crack porosity 6 and 8 show the same contours. That is ex-

294

pected, as the curves (which are essentially profiles of these maps) show very similar be-

295

havior and because the porosity equations 6 and 8 are the same functions, differing by a

296

proportionality factor. The velocity maps associated to the crack porosity equation 7 has

297

different contours from the other two maps in each elastic velocity, as expected, because its

298

porosity dependency is different from the other two, it depends only on aspect ratio. The

299

velocity maps associated to equation 7 show some problems: the V p0 and V sv maps have

300

some regions -low aspect ratio and high crack density regions- in which the model breaks

301

down, that is, produces unphysical results, in this case, it provides complex-valued velocity

302

values. That is, probably, because the crack porosity lacks the dependence on the aspect

303

ratio variable (and that analysis was one of our motivations to propose crack porosity equa-

304

tion 8, an equation which tries to take into account the background porosity and has the

305

same linear dependencies as equation 6).

306

The curves of velocity as a function of crack density (presented in De Figueiredo et al.

307

(2018)) can be seen as vertical profiles of the maps displayed in Figure 8, while the curves of

308

velocity as a function of aspect ratio can be seen as horizontal profiles of these maps. There-

309

fore both works compare profiles of these maps to experimental data, but De Figueiredo 13

310

et al. (2018) compares crack porosities 6 and 7, while this work includes a third formulation

311

(extending the previous work of De Figueiredo et al. (2018)) and compares it to the last

312

two. De Figueiredo et al. (2018) results shows that the fitting between the model and the

313

data -for velocities as a function of crack density- is, in general, better using crack porosity

314

equation 6 for lower aspect ratio values, while it is better using crack porosity equation 7 for

315

higher aspect ratios. A similar effect can be seen in our results mainly in Figure 6 (b), but

316

(as will be demonstrated) the overall errors in these windows related to the crack porosity 7

317

are still more significant than the ones related to the crack porosity equation 7. In order to

318

ensure a more rigid analysis in terms of comparisons between the results, we calculated the

319

root-mean-squared (RMS) relative error between theoretical values (elastic velocities and

320

Thomsen parameters calculated using the Eshelby-Cheng model) and experimental data,

321

for each graph shown in Figures 6 and 7. As the RMS error is calculated adding differences

322

(the squares of the relative differences and then taking the square root to be more precise)

323

between the curve and the data points along the entire range of aspect ratio, then it is a

324

measure related to the entire curve, so it is useful to solve the issue raised in the previous

325

paragraph (comparing curves that fit better the data in different regions). Each plot is

326

associated with a crack density; therefore the RMS error shown in Figure 9 is a function of

327

crack density.

328

[Figure 10 about here.]

329

In Figure 9, it is possible to observe that the error increases with crack density for all

330

formulations. The distance between data points and calculated curves increases with crack

331

density in the comparison plots shown both in this work and in De Figueiredo et al. (2018),

332

so this is further demonstration that all formulations proposed have a limitation concerning

333

high crack density values, differently from their behavior relative to crack aspect ratio,

334

which is different for each formulation: equation 7 produces velocities that behave better

335

for high crack aspect ratios and gets worse until it breaks down for low aspect ratios (as

336

discussed earlier), while values related to equation 6, as in Figure 6 for example, seem to

337

get further away from data points as aspect ratio increases.

14

338

In all graphs of Figure 9 it is clear that the errors associated to the velocities calculated

339

from the Eshelby-Cheng model with crack porosity equation 7 are the greatest, much higher

340

than the other ones, as anticipated by the analysis of Figure 6. As shown from the previous

341

observation, it is noting that the smaller error in all graphs is the one related to the use of

342

crack porosity equation 8, also as anticipated by the analysis of Figure 6. The RMS-errors

343

were also calculated for the Thomsen parameters and are shown in Figure 10. As the errors

344

associated with the crack porosity equation 7 surpassed 100%, they were normalized based

345

on the greatest error in each graph. Therefore it is not possible to draw comparisons of the

346

errors between different parameters, but that is not our goal here. [Figure 11 about here.]

347

348

Figure 9 shows an, even more, accentuated difference between errors, that may be due

349

to the fact that the additional mathematical operations between the calculations of the

350

elastic velocities and the Thomsen parameters may have amplified the errors and, as the

351

error associated with 7 was initially bigger, the differences between the errors also were

352

amplified. The same patterns from Figure 9 can be seen in Figure 10: The errors increase

353

with crack density for all formulations of crack porosity used and the errors associated to

354

the crack porosity equation 7 are the highest, while the errors associated to the use of the

355

crack porosity equation 8 are the lower.

CONCLUSIONS 356

357

The conclusions drawn from the results and the analysis made can be summarized in three points:

358

1. The workflow for applying the Eshelby-Cheng model adapted to porous anisotropic

359

cracked samples works best (produces the overall best fitting with the data) when fed

360

with the crack porosity calculated from equation 8.

361

2. When fed with crack porosities calculated from equation 7, the model produces the

362

worst fittings, mainly for low aspect ratio regions, and is unstable -producing un15

363

physical results- in regions of low aspect ratio and high crack density for the elastic

364

velocities V p0 and V sv. This behavior may have a relationship to the fact that equa-

365

tion 7 does not have the original dependence between crack porosity and cracks aspect

366

ratio.

367

3. For all three formulations tested, the RMS error grows with crack density, but for the

368

model fed with the crack porosity equation 8 the error keeps bellow 10%, except for

369

V p0 when the crack density is the highest tested (0.08).

370

Based on the points presented above, we conclude that the proposal which works best for the

371

adapted Eshelby-Cheng model proposed in De Figueiredo et al. (2018) is the crack porosity

372

equation 8, which present acceptable, considerably low in fact, values of RMS error. We

373

propose, for further works, testing this adapted model against models constructed from

374

first principles, like Guo et al. (2019), to describe cracked rocks with anisotropic porous

375

background.

ACKNOWLEDGEMENTS 376

We would like to thank the Brazilian agencies PET/GEOF´ISICA-UFPA-Minist´erio da

377

Educa¸c˜ao, INCT-GP and CNPq (Grant No.: CNPQ 459063/2014-6 and Grant No.: CNPQ

378

140174/2016-8) from Brazil and the Geophysics Graduate Program at Federal University

379

of Par´a for the financial support in this research.

APPENDIX A ESHELBY-CHENG EFFECTIVE MODEL 380

Here, in this work, for limitations related to saturation equipment in the laboratory,

381

we compared the Eshelby-Cheng’s effective with experimental results only in dry condi-

382

tion. However, for an extension, in this section, we show the mathematical formulation of

383

the Eshelby-Cheng’s effective model for both dry and saturated condition. According to

384

ef f Eshelby-Cheng’s effective model, the Cij for a rock containing weak-fluid filled ellipsoidal

16

385

cracks, can be evaluated by: dry,sat(ef f )

Cij

386

dry,sat(0)

= Cij

(dry,sat)−(1)

− φCij

,

(A-1)

where

387

c111 = λ(dry,sat) (S31 − S33 + 1) + c133 = 388

c113 =

2µ(dry,sat) (S33 S11 − S31 S13 − (S33 + S11 − 2C − 1) + C(S31 + S13 − S11 − S33 )) , D(S12 − S11 + 1)

(λ(dry,sat) + 2µ(dry,sat) )(−S12 − S11 + 1) + 2λ(dry,sat) S13 + 4µ(dry,sat) C , D

(λ(dry,sat) + 2µ(dry,sat) )(S13 + S31 ) − 4µ(dry,sat) C + λ(dry,sat) (S13 − S12 − S11 − S33 + 2) , 2D

S33

c144 =

µ(dry,sat) , 1 − 2S1313

c166 =

µ(dry,sat) , 1 − 2S1212

S11 = QIαα + RIα ,   4π 2 =Q − 2Iαc α + Ic R, 3 S12 = QIαb − RIα , S13 = QIαc α2 − RIα ,

389

S31 = QIαc − RIc , S1212 = QIαb + RIα , S1313 =

Q(1 + α2 )Iαc R(Iα + Ic ) + , 2 2

17

Iα =

2πα(cos−1 α − αSα ) , Sα3 Ic = 4π − 2Iα ,

Iαα = π − (3/4)Iαc , 390

Ic − Iα , 3Sα3 Iαα = , 3

Iαc = Iab

Sα =

p 1 − α2 ,

3K (dry,sat) − 2µ(dry,sat) , 6K + 2µ(dry,sat) 1 − 2σ R= , 8π(1 − σ) 3R Q= , 1 − 2σ kf C= , (dry,sat) 3(K − kf )

σ=

391

I = K (dry,sat) −

2µ(dry,sat) , 3

392

where α is the crack aspect ratio, K (dry,sat) is the bulk modulus of solid matrix, kf (0.0

393

GPa and 2.2 GPa for air and water, respectively) is the fluid bulk modulus. The λ(dry,sat)

394

and µ(dry,sat) are the Lam´e parameters for dry or water saturated reference sample (ma-

395

trix without inclusions), respectively. According to the flowchart showed in Figure 5, the

396

reference bulk and shear modulus are

K 397

dry,sat

dry,sat

=ρ

  4 dry,sat 2 dry,sat 2 (VP 0 ) + (VS0 ) , 3

dry,sat 2 µ(dry,sat) = ρdry,sat (VS0 ) .

398

18

(A-2)

(A-3)

APPENDIX B EFFECTIVE THOMSEN PARAMETER 399

Figure 5 shows a flowchart with the ”ingredients” used to feed the modified Eshelby-

400

Cheng’s model. As our reference sample is not isotropic, to feed the Eshelby-Cheng’s model,

401

we used ρref (see Table 1), VP 0 = (VPref (θ = 0o ) + VPref (θ = 90o ))/2, VS0 = (VSref (ϕ =

402

0o ) + VSref (ϕ = 90o ))/2 (see Table 3). The output of model (see branch (2) of flowchart)

403

are elastic coefficient constant as function of crack parameter (, α). These coefficients can

404

be converted in velocities using equations such as

ef f (dry,sat)

VP,S

v u ef f (dry,sat)−mod u Cij ∝t , (dry,sat) ρef f

(B-1)

405

where ρef f the the theoretical effective density shown in the flowchart of Figure 5.e Specif-

406

ically, the effective velocities (dry and/or sat) as well as the effective Thomsen parameters

407

(also dry and/or sat) can be estimated (from effective stiffenes coefficient) by the following

408

equations: ef f C11

= ρef f (VP2 (θ = 90o ))ef f ,

(B-2)

ef f C33

= ρef f (VP2 (θ = 0o ))ef f ,

(B-3)

ef f C44

2 = ρef f (VS2 (ϕ = 90o ))ef f ,

(B-4)

ef f C66

2 = ρef f (VS1 (ϕ = 0o ))ef f , (B-5) 1   2  2 2 ef f ef f ef f (V 2 (θ = 45o ))ef f − C ef f − C ef f − 2C ef f 4ρ − C − C 11 33 44 11 33 P   ef f =   − C(B-6) 44 , , 4

ef f C13

409

with ef f ef f ef f C12 = C11 − 2C66 ,

19

(B-7)

410

and phic rampinγ ef f

=

ef f ef f C66 − C44

,

(B-8)

εef f

=

ef f ef f C11 − C33

,

(B-9)

δ ef f

=

ef f ef f 2 ef f ef f 2 (C13 + C44 ) − (C33 − C44 )

ef f 2C44 ef f 2C33

ef f ef f ef f 2C33 (C33 − C44 )

20

.

(B-10)

REFERENCES 411

412

413

414

Cheng, C. H., 1993, Crack models for a transversely isotropic medium: Journal of Geophysical Research: Solid Earth, 98, 675–684. Crampin, S., and S. Peacock, 2005, A review of shear-wave splitting in the compliant crackcritical anisotropic Earth: Wave Motion, 41, 59–77.

415

De Figueiredo, J. J. S., M. J. do Nascimento, E. Hartmann, B. F. Chiba, C. B. da Silva,

416

M. C. De Sousa, C. Silva, and L. K. Santos, 2018, On the application of the eshelby-cheng

417

effective model in a porous cracked medium with background anisotropy: An experimental

418

approach: Geophysics, 83, C209–C220.

419

420

421

422

Dellinger, J., and L. Vernick, 1994, Do traveltimes in pulse-transmission experiments yield anisotropic group or phase velocities?: Geophysics, 59, 1774–1779. Eshelby, J. D., 1957, The determination of the elastic field of an ellipsoidal inclusion, and related problems: Proc. Roy. Soc. London, 241, 376–396.

423

Far, M. E., C. M. Sayers, L. Thomsen, D.-h. Han, and J. P. Castagna, 2013, Seismic char-

424

acterization of naturally fractured reservoirs using amplitude versus offset and azimuth

425

analysis: Geophysical Prospecting, 61, 427–447.

426

427

Grechka, V., and M. Kachanov, 2006, Effective elasticity of rocks with closely spaced and intersecting cracks: GEOPHYSICS, 71, D85–D91.

428

Guo, J., T. Han, L.-Y. Fu, D. Xu, and X. Fang, 2019, Effective elastic properties of rocks

429

with transversely isotropic background permeated by aligned penny-shaped cracks: Jour-

430

nal of Geophysical Research: Solid Earth, 124, 400–424.

431

Henriques, J. P., J. J. S. de Figueiredo, L. Kirchhof, D. L. Macedo, I. Coutinho, C. B.

432

da Silva, A. L. Melo, and M. Sousa, 2018, Experimental verification of effective anisotropic

433

crack theories in variable crack aspect ratio medium: Geophysical Prospecting, 66, 141–

434

156.

435

436

437

438

Hudson, J. A., 1981, Wave speeds and attenuation of elastic waves in material containing cracks: Geophysical Journal of the Royal Astronomical Society, 64, 133–150. M. Kachanov, B. Shafiro, and Igor Tsukrov, 2003, Handbook of Elasticity Solutions: Kluwer Academic Publishers.

21

439

440

441

442

443

444

Nishizawa, O., and T. Yoshino, 2001, Seismic velocity anisotropy in mica-rich rocks: an inclusion model: Geophysical Journal International, 145, 19–32. Poletto, M., H. L. Junior, and A. J. Zaterra, 2014, Polystyrene : synthesis, characteristics, and applications: Nova Science Publishers. Santos, L. K., J. J. S. de Figueiredo, and C. B. da Silva, 2016, A study of ultrasonic physical modeling of isotropic media based on dynamic similitude: Ultrasonics, 70, 227–237.

445

Santos, L. K., J. J. S. de Figueiredo, C. B. Silva, D. L. Macedo, and A. M. Leandro, 2017, A

446

new way to construct synthetic porous fractured medium: Journal of Petroleum Science

447

and Eng., 156, 763–768.

448

Thomsen, L., 1986, Weak elastic anisotropy: GEOPHYSICS, 51, 1954–1966.

22

LIST OF FIGURES 449

1

450 451 452 453 454 455 456 457

458

2

459 460 461 462 463 464 465 466 467 468 469 470

471

3

472 473 474 475 476

477

4

478 479 480 481 482 483

484 485 486 487 488 489 490 491 492

5

(a) Sample construction of a group of samples with inclusion’s aspect ratio equal to 0.08. (b) Cylindrical samples (plugs) are obtained by a plugging machine. c) Photograph of all samples used in this work. The samples are divided by the aspect ratio of their inclusions, forming four groups (α1 = 0.08, α2 = 0.2, α3 = 0.32, α4 = 0.52) of four samples each, plus a reference -with no inclusions- sample identified as REF in the image. The experimental crack density (), aspect ratio (α) as well as the dry sample density-ρb and porosity (φ) are shown in the Table 1. These pictures were modified from in De Figueiredo et al. (2018). . . . . . . . . . . . . . . . . . . . . . . . . . . .

25

(a) Picture of the ultrasonic system setup used in this work. (b) Experimental arrangement of the transducers, in the samples, for P- and S- wave acquisitions. S1 and S2 indicate the faster and slower S-wave polarization, respectively. PZ indicates P-wave propagation in the Z direction. c) Experimental sketch of wave propagation direction and polarization to obtain the P- and S-waveform records. The angle θ is the angle between the direction of propagation of the P-wave and the Z axis, while ϕ is the angle between the direction of propagation of the S-wave and the Z axis. For the P-wave, the propagation is radial changing from θ = 0o (Z direction) to θ = 90o (X direction) and for S-wave, the propagation is in the Y direction and the slow and fast polarization correspond to ϕ = 0o (X direction) and ϕ = 90o (Z direction), respectively. These pictures were modified from in De Figueiredo et al. (2018). . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

26

Interpolation performed to the velocities as function of crack density in order to produce an estimation of velocities as function of crack aspect ratio. Each column is associated with a fixed value of crack aspect ratio, while each graph shows velocities as function of crack density. The blue circles represent the experimental velocity data. The red lines represent the linear interpolation of the experimental data. . . . . . . . . . . . . . . . . . . . . . . . . . . . .

27

Linear regression performed to the Thomsen parameters as a function of crack density in order to produce an estimation of Thomsen parameters as a function of crack aspect ratio. Each column is associated with a fixed value of crack aspect ratio, while each line of graphs shows a Thomsen parameter as a function of crack density. The blue circles represent the Thomsen parameters calculated directly from experimental data. The red lines represent the linear interpolation of the data. . . . . . . . . . . . . . . . . . . . . . . . . . . . .

28

Workflow used to calculate the effective theoretical velocities and Thomsen parameters. Block (1) shows the input parameters, obtained from elastic parameters (VP , VS and ρ) measured in the background reference sample as well as the physical parameters of inclusions (α, , φref and φc ). Block (2) shows the output obtained from equation 1, i.e., the effective elastic stiffness parameters. These coefficients can be used to evaluate the Thomsen parameters, as represented by block (3), or combined with the theoretical effective density, they can be inverted on theoretical elastic velocities, as represented by block (4). . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

29

23

493

6

494 495 496 497 498 499 500 501 502 503 504

505

7

506 507 508 509 510 511 512 513 514 515 516 517 518 519 520

521

8

522 523 524 525 526 527 528

529

538

Velocity maps produced using Eshelby-Cheng effective medium theory. The colors represent the velocity values in m/s. the vertical axis shows values of crack density, while the horizontal axis shows values of aspect ratio. Each line of maps shows one of the five elastic velocities studied in this work. Each column of maps is associated with one of the three crack porosity equations tested in this work. There are some regions in the maps (first line, second column and fourth line, second column) in which the model, using 5, produced complex valued velocities, that is, unphysical results. . . . . . . . . . . . . .

33

Root-mean-squared error between Thomsen parameters calculated from the Eshelby-Cheng model, using the three crack porosity formulations, and experimental elastic velocities. The errors are shown for each fixed crack density associated with the graph from which the anisotropic parameter errors were calculated. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

533

537

32

10

532

536

The vertical axis shows values of the Thomsen parameters -which are dimensionless and always vary between zero and one- while the horizontal axis shows values of crack aspect ratio. The first line of graphs shows the values for the  parameter (also known as P-wave anisotropy), the second line of graphs shows the values for the γ parameter (also known as P-wave anisotropy), and the third line shows the values for the δ parameter, which is related both to P- and S-waves, while each column is associated to a fixed value of crack density. The purple circles represent the values of the parameters obtained via interpolation of the experimental data. All three curves represent parameters obtained from velocities calculated using the Eshelby-Cheng model workflow shown in Figure 5. The continuous blue curves represent the velocity values calculated using the crack porosity estimated by equation 7, while the red dashed curves represent the velocity values calculated using the crack porosity estimated by equation 5 and the yellow curves made of dashes and points represent the velocities calculated using the model fed by crack porosity equation 8. . . . . . . . . . . . . . . . . . . . . . . . . . . . .

Root-mean-squared error between elastic velocities calculated from the EshelbyCheng model, using the three crack porosity formulations, and experimental elastic velocities. The errors are shown for each fixed crack density associated to the graph from which the velocity errors were calculated. (a) P-wave velocity errors. (b) S-wave velocity errors. . . . . . . . . . . . . . . . . . . . 34

531

535

31

9

530

534

The vertical axis shows values of velocity, while the horizontal axis shows values of crack aspect ratio. Each line of graphs shows an elastic velocity with a specific direction of propagation, while each column of graphs is related to a fixed value of crack density. The purple circles represent the velocity values obtained via interpolation of the experimental data. All three curves represent velocity values calculated using the Eshelby-Cheng model workflow shown in Figure 5. The blue continuous curves represent the velocity values calculated using the crack porosity estimated by equation 7, while the red dashed curves represent the velocity values calculated using the crack porosity estimated by equation 5. The yellow curves made of dashes and points represent the velocities calculated using the model fed by crack porosity equation 8. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

24

35

c)

Figure 1: (a) Sample construction of a group of samples with inclusion’s aspect ratio equal to 0.08. (b) Cylindrical samples (plugs) are obtained by a plugging machine. c) Photograph of all samples used in this work. The samples are divided by the aspect ratio of their inclusions, forming four groups (α1 = 0.08, α2 = 0.2, α3 = 0.32, α4 = 0.52) of four samples each, plus a reference -with no inclusions- sample identified as REF in the image. The experimental crack density (), aspect ratio (α) as well as the dry sample density-ρb and porosity (φ) are shown in the Table 1. These pictures were modified from in De Figueiredo et al. (2018).

25

(a)

(b)

(c )

Figure 2: (a) Picture of the ultrasonic system setup used in this work. (b) Experimental arrangement of the transducers, in the samples, for P- and S- wave acquisitions. S1 and S2 indicate the faster and slower S-wave polarization, respectively. PZ indicates P-wave propagation in the Z direction. c) Experimental sketch of wave propagation direction and polarization to obtain the P- and S-waveform records. The angle θ is the angle between the direction of propagation of the P-wave and the Z axis, while ϕ is the angle between the direction of propagation of the S-wave and the Z axis. For the P-wave, the propagation is radial changing from θ = 0o (Z direction) to θ = 90o (X direction) and for S-wave, the propagation is in the Y direction and the slow and fast polarization correspond to ϕ = 0o (X direction) and ϕ = 90o (Z direction), respectively. These pictures were modified from in De Figueiredo et al. (2018). 26

(a) α = 0.08

Vp0 (m/s)

4500

4000

4000

3500

3500

3500

3500

3000 0.05

0.1

3000 0

0.05

0.1

3000 0

0.05

0.1

4500

4500

4500

4500

4000

4000

4000

4000

3500

3500

3500

3500

3000

3000 0

0.05

0.1

3000 0

0.05

0.1

0.05

0.1

4500

4500

4500

4000

4000

4000

4000

3500

3500

3500

3500

3000 0

0.05

0.1

3000 0

crack density

0.05

0.1

0

0.05

0.1

0

0.05

0.1

0.05

0.1

3000 0

4500

3000

α = 0.52

4500

4000

0

Vp45 (m/s)

α = 0.32

4500

4000

3000

3000 0

crack density

0.05

0.1

0

crack density

crack density

(b)

Vsv (m/s)

α = 0.08

α = 0.2

α = 0.32

α = 0.52

2400

2400

2400

2400

2200

2200

2200

2200

2000

2000

2000

2000

0

Vsh (m/s)

Vp90 (m/s)

α = 0.2

4500

0.05

0.1

0

0.05

0.1

0

0.05

0.1

0

2400

2400

2400

2400

2200

2200

2200

2200

2000

2000

2000

2000

0

0.05

0.1

crack density

0

0.05

0.1

0

crack density

0.05

0.1

0

crack density

0.05

0.1

0.05

0.1

crack density

Figure 3: Interpolation performed to the velocities as function of crack density in order to produce an estimation of velocities as function of crack aspect ratio. Each column is associated with a fixed value of crack aspect ratio, while each graph shows velocities as function of crack density. The blue circles represent the experimental velocity data. The red lines represent the linear interpolation of the experimental data.

27

α = 0.08

ǫ

1

0.5

γ

0.05

0.1

0.5

0 0

0.05

0.1

0 0

0.05

0.1

1

1

1

1

0.5

0.5

0.5

0.5

0

0 0

0.05

0.1

0 0

0.05

0.1

0.05

0.1

1

1

1

0.5

0.5

0.5

0.5

0 0

0.05

crack density

0.1

0 0

0.05

0.1

0

0.05

0.1

0

0.05

0.1

0.05

0.1

0 0

1

0

α = 0.52

1

0.5

0 0

α = 0.32

1

0.5

0

δ

α = 0.2

1

0 0

crack density

0.05

crack density

0.1

0

crack density

Figure 4: Linear regression performed to the Thomsen parameters as a function of crack density in order to produce an estimation of Thomsen parameters as a function of crack aspect ratio. Each column is associated with a fixed value of crack aspect ratio, while each line of graphs shows a Thomsen parameter as a function of crack density. The blue circles represent the Thomsen parameters calculated directly from experimental data. The red lines represent the linear interpolation of the data.

28

(1) dry,sat dry,sat o = 90 ), Vref −P (θ = 0 ), Vref −P (θ dry,sat dry,sat o Vref −S (φ = 90 ), Vref −S (φ = 0o )

Input Parameters

dry,sat Vref −P (θ

o

ρref {z

|

}

(0)−V T I Cij

dry,sat dry,sat o o Vref −P (θ = 90 ) + Vref −P (θ = 0 ) 2 dry,sat dry,sat o o Vref −S (φ = 90 ) + Vref −S (φ = 0 ) = 2 , α, φref , φc {z } (1)−mod (2)−mod , Cij Cij

VPdry,sat = 0 dry,sat VS0

|

= 45o )

4 φc = πα 3 or φc = (1 − φref ) or 4 φc = π(1 − φref )α 3

ef f (dry,sat)−mod

Cij Hudson-Crampin’s (Modified) Model ef f (dry,sat)

C11

ef f (dry,sat)

, C33

ef f (dry,sat)

C44

ρdry ef f

(dry,sat)

ef f

ef f (dry,sat)

, C13

(dry,sat)

, γef f

(dry,sat)

, δef f

ef f (dry,sat)

, C66

ef f (dry,sat)

VP,S

(3) Theoretical Thomsen Parameters

(4) v u ef f (dry,sat)−mod u Cij ∝t (dry,sat) ρef f

= (1 −

f φef total )ρmin

or

ef f water ef f ρsat ef f = (1 − φtotal )ρmin + φtotal ρf l

ρref = (1 − φref )ρmin ρref ρmin = (1 − φref ) ρmin = 2186.5 kg/m3 and ρwater = 1500 kg/m3 fl

Theoretical Velocities

Output Parameters

(2)

ef f (dry,sat)

VP

(θ = 0o )

ef f (dry,sat) VP (θ = 45o ) ef f (dry,sat) VP (θ = 90o ) ef f (dry,sat) VS1 (φ = 0o ) ef f (dry,sat) VS2 (φ = 90o )

4 f φef total = φref + πα 3 or f φef total = φref + (1 − φref ) or 4 f φef = φ + π(1 − φref )α ref total 3

Figure 5: Workflow used to calculate the effective theoretical velocities and Thomsen parameters. Block (1) shows the input parameters, obtained from elastic parameters (VP , VS and ρ) measured in the background reference sample as well as the physical parameters of inclusions (α, , φref and φc ). Block (2) shows the output obtained from equation 1, i.e., the effective elastic stiffness parameters. These coefficients can be used to evaluate the Thomsen parameters, as represented by block (3), or combined with the theoretical effective density, they can be inverted on theoretical elastic velocities, as represented by block (4).

29

(a) ǫ = 0.02

Vp0 (m/s)

4500

ǫ =0.04

4500 4000

4000

4000

3500

3500

3500

3500

3000

3000

3000

3000

2500

2500

2500

0.2

0.4

0.6

0

Vp45 (m/s)

aspect ratio

0.2

0.4

0.6

2500 0

aspect ratio

0.2

0.4

0.6

0

aspect ratio

4500

4500

4500

4000

4000

4000

4000

3500

3500

3500

3500

3000

3000

3000

3000

2500 0

0.2

0.4

0.6

2500 0

aspect ratio

0.2

0.4

0.6

aspect ratio

0.2

0.4

0.6

0

aspect ratio 4500

4500

4000

4000

4000

4000

3500

3500

3500

3500

3000

3000

3000

3000

2500 0.2

0.4

0.6

aspect ratio

2500 0

0.2

0.4

0.6

0.2

0.4

0.6

0

aspect ratio

φ = φ 0 + 4 π α ǫ/3

30

0.2

0.4

0.6

2500 0

aspect ratio

0.6

aspect ratio

4500

0

0.4

2500 0

4500

2500

0.2

aspect ratio

4500

2500

ǫ =0.08

4500

4000

0

Vp90 (m/s)

ǫ =0.06

4500

φ = ǫ (1- φ 0 )

φ = φ 0 + 4 π αǫ (1- φ 0 )/3

0.2

0.4

0.6

aspect ratio experimental

(b)

ǫ =0.02

Vsv (m/s)

2500

ǫ =0.04

2500 2400

2400

2400

2300

2300

2300

2300

2200

2200

2200

2200

2100

2100

2100

2100

2000

2000

2000

2000

1900 0

0.5

1900 0

aspect ratio

0.5

1900 0

aspect ratio

0.5

0

aspect ratio

2500

2500

2500

2400

2400

2400

2400

2300

2300

2300

2300

2200

2200

2200

2200

2100

2100

2100

2100

2000

2000

2000

2000

1900 0

0.5

aspect ratio

1900 0

0.5

1900 0

aspect ratio

0.5

aspect ratio

2500

1900

ǫ =0.08

2500

2400

1900

Vsh (m/s)

ǫ =0.06

2500

0.5

aspect ratio

0

0.5

aspect ratio

Figure 6: The vertical axis shows values of velocity, while the horizontal axis shows values of crack aspect ratio. Each line of graphs shows an elastic velocity with a specific direction of propagation, while each column of graphs is related to a fixed value of crack density. The purple circles represent the velocity values obtained via interpolation of the experimental data. All three curves represent velocity values calculated using the Eshelby-Cheng model workflow shown in Figure 5. The blue continuous curves represent the velocity values calculated using the crack porosity estimated by equation 7, while the red dashed curves represent the velocity values calculated using the crack porosity estimated by equation 5. The yellow curves made of dashes and points represent the velocities calculated using the model fed by crack porosity equation 8.

31

ǫ = 0.02

ǫ

1

ǫ = 0.04

1

0.5

0.5

0 0.5

γ

aspect ratio

0.5

0 0

aspect ratio

0.5

0

aspect ratio

1

1

1

0.5

0.5

0.5

0.5

0 0

0.5

0 0

aspect ratio

0.5

0 0

aspect ratio

0.5

0

aspect ratio

1

1

1

0.5

0.5

0.5

0.5

0 0

0.5

aspect ratio

0 0

0.5

0 0

aspect ratio

0.5

0

aspect ratio

φ c=4 παǫ/3

0.5

aspect ratio

1

0

0.5

aspect ratio

1

0

δ

0.5

0 0

ǫ = 0.8

1

0.5

0 0

ǫ = 0.06

1

φ c= ǫ(1- φ 0 )

φ c=4 παǫ(1- φ 0 )/3

0.5

aspect ratio experimental

Figure 7: The vertical axis shows values of the Thomsen parameters -which are dimensionless and always vary between zero and one- while the horizontal axis shows values of crack aspect ratio. The first line of graphs shows the values for the  parameter (also known as P-wave anisotropy), the second line of graphs shows the values for the γ parameter (also known as P-wave anisotropy), and the third line shows the values for the δ parameter, which is related both to P- and S-waves, while each column is associated to a fixed value of crack density. The purple circles represent the values of the parameters obtained via interpolation of the experimental data. All three curves represent parameters obtained from velocities calculated using the Eshelby-Cheng model workflow shown in Figure 5. The continuous blue curves represent the velocity values calculated using the crack porosity estimated by equation 7, while the red dashed curves represent the velocity values calculated using the crack porosity estimated by equation 5 and the yellow curves made of dashes and points represent the velocities calculated using the model fed by crack porosity equation 8. 32

(a) φc = ǫ (1- φ0 )

φc =4 παǫ /3

Vp 0 (m/s)

φc =4 παǫ (1- φ0 )/3

Vp 0 (m/s)

Vp 0 (m/s)

4000

0.02

4000

0.04

3000

0.04

3000

0.04

3000

0.06

2000

0.06

2000

0.06

2000

1000

0.08

1000

0.08

0.08 0

0.2

0.4

ǫ

0.02

ǫ

4000

ǫ

0.02

0.6

0

0.2

complex-valued velocities

α Vp 45 (m/s)

0.4

0.6

1000

0

0.2

0.4

α

0.6

α

Vp 45 (m/s)

Vp 45 (m/s)

5000

0.02

5000

0.04

4500

0.04

4500

0.04

4500

0.06

4000

0.06

4000

0.06

4000

3500

0.08

3500

0.08

0.08 0

0.2

0.4

ǫ

0.02

ǫ

5000

ǫ

0.02

0.6

0

0.2

α

0.4

0.6

3500

0

0.2

0.4

α

Vp 90 (m/s)

0.6

α

Vp 90 (m/s)

Vp 90 (m/s)

4320

0.02

4320

0.04

4300

0.04

4300

0.04

4300

0.06

4280

0.06

4280

0.06

4280

4260

0.08

4260

0.08

0.08 0

0.2

0.4

ǫ

0.02

ǫ

4320

ǫ

0.02

0.6

0

0.2

α

0.4

0.6

4260

0

0.2

0.4

α

0.6

α

(b) φc = 4 παǫ /3

φc = ǫ (1- φ0 )

Vsv (m/s) 0.02

2000

0.03

1800

φc =4 παǫ (1- φ0 )/3

Vsv (m/s)

0.04

0.02

2000

0.03

1800

0.03

1800

0.04

0.04 1600

0.05

1600

ǫ

0.05

ǫ 1400

1400

0.06

0.07 0.08 0

0.2

0.4

1200

0.07

1000

0.08

0.6

α

0.06

0

0.2

1000

0.08

2270

ǫ

0.04 2260

0.05

2240

0.2

0.4

α

0.6

2260

0.05 0.06

2250

0.07

0

2280

2270

2250

0.08

0.6

0.03

0.06

0.07

0.4

0.02

0.04 2260

0.2

Vsh (m/s)

2280

2270

0.04

0.06

1000

0

0.03

0.05

1200

α

0.02

0.03

0.07

0.6

Vsh (m/s)

2280

0.02

0.4

1200

α

complex-valued velocities

Vsh (m/s)

0.05

1400

0.06

ǫ

ǫ

1600

ǫ

Vsv (m/s)

0.02

2000

2250

0.07

0.08

2240

0

0.2

0.4

α

0.6

0.08

2240

0

0.2

0.4

0.6

α

Figure 8: Velocity maps produced using Eshelby-Cheng effective medium theory. The colors represent the velocity values in m/s. the vertical axis shows values of crack density, while the horizontal axis shows values of aspect ratio. Each line of maps shows one of the five elastic velocities studied in this work. Each column of maps is associated with one of the three crack porosity equations tested in this work. There are some regions in the maps (first line, second column and fourth line, second column) in which the model, using 5, produced 33 complex valued velocities, that is, unphysical results.

(a) Vp0-RMS-error

0.8

Vp45-RMS-error

0.8

Vp90-RMS-error

0.1

φ = 4 παǫ/3 c

φ c= ǫ(1- φ 0 )

0.7

0.09

0.7

φ c=4 παǫ(1- φ 0 )/3

0.08

0.6

0.6

0.5

0.5

0.06

error

error

0.4

0.4

0.3

0.3

0.2

0.2

0.1

0.1

0.05 0.04 0.03 0.02

0

0.01

0 0.02

0.04

0.06

0.08

0.1

0 0

0.02

crack density

0.04

0.06

0.08

0.1

0

0.02

crack density

0.04

0.06

0.08

crack density

(b) Vsv-RMS-error

0.5

Vsh-RMS-error

0.1

0.45

0.09

0.4

0.08

0.35

0.07

0.3

0.06

error

0

error

error

0.07

0.25

0.05

0.2

0.04

0.15

0.03

0.1

0.02

0.05

0.01

0

0 0

0.02

0.04

0.06

0.08

0.1

0

crack density

0.02

0.04

0.06

0.08

0.1

crack density

Figure 9: Root-mean-squared error between elastic velocities calculated from the EshelbyCheng model, using the three crack porosity formulations, and experimental elastic velocities. The errors are shown for each fixed crack density associated to the graph from which the velocity errors were calculated. (a) P-wave velocity errors. (b) S-wave velocity errors.

34

0.1

epsilon-RMS-error

1

gamma-RMS-error

0.9

delta-RMS-error

1

φ = 4 παǫ/3 c

0.9

0.8

0.8

0.7

0.6

0.6

0.5

0.5

error

error

0.6

error

c

0.4

0.5 0.4

0.4 0.3

0.3

0.3

0.2

0.2

0.2

0.1

0.1 0

0.1

0 0

0.02

0.04

0.06

crack density

0.08

0.1

0 0

0.02

0.04

0.06

0.08

0.1

crack density

0

0.02

0.04

0.06

0.08

0.1

crack density

Figure 10: Root-mean-squared error between Thomsen parameters calculated from the Eshelby-Cheng model, using the three crack porosity formulations, and experimental elastic velocities. The errors are shown for each fixed crack density associated with the graph from which the anisotropic parameter errors were calculated.

35

0

φ c=4 παǫ(1- φ 0 )/3

0.8

0.7

0.7

φ = ǫ(1- φ )

0.9

LIST OF TABLES 539

1

540 541

542

2

543 544

545

3

Sample parameters description (part I): D = sample diameter; H = sample height; h = crack aperture; N = number of inclusions. The separation between crack planes was 5 mm. . . . . . . . . . . . . . . . . . . . . . . . .

37

Sample parameters description (part II): α = aspect ratio;  = fracture density; ρdry = dry density sample (g/cm3 ); ρsat = saturated density sample, φ = total sample porosity. . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

38

P- and S-wave velocities for dry condition. . . . . . . . . . . . . . . . . . . .

39

36

Table 1: Sample parameters description (part I): D = sample diameter; H = sample height; h = crack aperture; N = number of inclusions. The separation between crack planes was 5 mm. Sample Reference α1 1 α1 2 α1 3 α1 4 α2 1 α2 2 α2 3 α2 4 α3 1 α3 2 α3 3 α3 4 α4 1 α4 2 α4 3 α4 4

D (mm) 38.00 38.00 38.00 38.00 38.00 38.00 38.00 38.00 38.00 38.00 38.00 38.00 38.00 38.00 38.00 38.00 38.00

H (mm) 59.00 53.10 53.00 53.90 51.10 55.50 54.10 55.70 52.50 52.50 51.31 51.83 50.81 61.15 60.60 61.00 61.60

37

h (mm) 0 0.50 0.50 0.50 0.50 1.25 1.25 1.25 1.25 2.00 2.00 2.00 2.00 3.75 3.75 3.75 3.75

N 0 36 54 72 90 18 30 66 96 18 36 54 72 36 48 60 72

Table 2: Sample parameters description (part II): α = aspect ratio;  = fracture density; ρdry = dry density sample (g/cm3 ); ρsat = saturated density sample, φ = total sample porosity. Sample No cracks α1 1 α1 2 α1 3 α1 4 α2 1 α2 2 α2 3 α2 4 α3 1 α3 2 α3 3 α3 4 α4 1 α4 2 α4 3 α4 4

α 0.00 0.08 0.08 0.08 0.08 0.20 0.20 0.20 0.20 0.32 0.32 0.32 0.32 0.52 0.52 0.52 0.52

 0.0000 0.0092 0.0138 0.0181 0.0238 0.0110 0.0188 0.0401 0.0618 0.0185 0.0380 0.0564 0.0767 0.0518 0.0696 0.0865 0.1028

38

ρdry 1.9460 1.9260 1.9159 1.9065 1.8939 1.9220 1.9050 1.8584 1.8108 1.9054 1.8630 1.8228 1.7784 1.8328 1.7937 1.7569 1.7213

ρsat 2.0560 2.0451 2.0396 2.0346 2.0277 2.0430 2.0338 2.0085 1.9826 2.0340 2.0110 1.9891 1.9650 1.9946 1.9734 1.9534 1.9341

φ 0.1100 0.1191 0.1237 0.1281 0.1338 0.1210 0.1288 0.1501 0.1718 0.1286 0.1480 0.1663 0.1866 0.1618 0.1797 0.1965 0.2128

Table 3: P- and S-wave velocities for dry condition. Sample Reference α1 1 α1 2 α1 3 α1 4 α2 1 α2 2 α2 3 α2 4 α3 1 α3 2 α3 3 α3 4 α4 1 α4 2 α4 3 α4 4

VP (0o ) 4194.3 4102.8 4032.7 3997.1 3938.9 4079.1 3983.6 3762.9 3574.8 3985.7 3771.2 3606.5 3431.2 3634.2 3453.3 3305.2 3163.7

VP (45o ) 4231.6 4171.2 4155.8 4134.2 4102.8 4165.6 4116.8 3986.5 3847.7 4123.3 3992.8 3883.9 3769.4 3908.3 3819.6 3689.8 3652.2

VP (90o ) 4260.1 4257.4 4254.0 4251.6 4250.8 4256.0 4253.4 4248.7 4236.4 4259.2 4249.8 4243.2 4237.5 4242.4 4238.5 4231.5 4225.4

39

VSH (90o ) 2267.5 2265.4 2264.1 2262.3 2258.6 2261.7 2258.4 2256.9 2255.2 2256.4 2253.8 2252.9 2252.6 2254.8 2252.0 2251.3 2250.9

VSV (90o ) 2159.6 2146.1 2137.8 2134.6 2122.8 2143.2 2131.3 2094.6 2065.9 2128.6 2098.7 2067.7 2031.3 2073.7 2045.0 2019.3 1992.8

1. Ultrasonic analysis on cracked medium; 2. Effective medium is important to theoretical characterization of cracked medium; 3. Unconventional shale reservoir are an isotropic; 4. Characterization of fractured medium is important fro hydraulic fracturing.