# Assessing aperture size while injectin in a joint plane

Asked by Yousef Golabchi on 2020-04-10

Hello Everyone,

I doing a simulation using Yade. In my case the sample has a joint plane (with a dip angle). I'm doing injection in the middle of the joint plane. I would like to assess how the aperture will change during injection.

Please let me know if there is any code to understand and monitor the aperture change during injection.

Best regards,
Yousef

## Question information

Language:
English Edit question
Status:
Solved
For:
Assignee:
No assignee Edit question
Solved by:
Luc Scholtès
Solved:
2020-04-16
Last query:
2020-04-16
2020-04-16

## This question was reopened

 Luc Scholtès (luc) said on 2020-04-10: #1

Hi Yousef,

2 ways:

1) record the positions of particles located on each side of the joint

I would suggest option 2) if you use the neverErase=True which keep track of the interactions during the duration of the simulation. To do so, you could (this is just one possibility from someone who is not very good at programming):

1) create a list containing all interactions "onJoint". You can do that right after the firts timestep byt looping over all interactions, checking for those which are "onJoint":

intOnFracPlane=[]
for i in O.interactions:
if i.phys.isOnJoint:
intOnFracPlane.append(i)

2) record their aperture "periodically " during the simulation with a dedicated fonction that you would call every n iterations like:

def recordInteractionStatus():
global intOnFracPlane
inFile=open('interactions_timePositionsAperturePressure_'+str(O.iter)+'_'+str(round(O.time,3)),'a')
for i in intOnFracPlane:
inFile.write( str(i.geom.contactPoint[0]) + '\t' + str(i.geom.contactPoint[1]) + '\t' + str(i.geom.contactPoint[2]) + '\t' + str(i.phys.crackJointAperture) + '\t' + str(flow.getPorePressure((i.geom.contactPoint[0],i.geom.contactPoint[1],i.geom.contactPoint[2]))) + '\n')
inFile.close()

Doing so, you will get several text files with the information you want and you'll then be able to map the evolution of the aperture over all the joint surface.

I guess you could also adapt these lines to get some sort of average value (probably not meaningful though) or, another option, the aperture at a given location (e.g., injection point).

Luc

 Yousef Golabchi (yousef-golabchi) said on 2020-04-12: #2

Hi Luc,

Thanks a lot. The code works fine and I am really appreciated you help me with that. It just needed a tiny change in the lines.

#####

def recordInteractionStatus():
global intOnFracPlane
intOnFracPlane=[]
for i in O.interactions:
if i.phys.isOnJoint:
intOnFracPlane.append(i)
inFile=open('interactions_timePositionsAperturePressure_'+str(O.iter)+'_'+str(round(O.time,3)),'a')
for i in intOnFracPlane:
inFile.write( str(i.geom.contactPoint[0]) + '\t' + str(i.geom.contactPoint[1]) + '\t' + str(i.geom.contactPoint[2]) + '\t' + str(i.phys.crackJointAperture) + '\t' + str(flow.getPorePressure((i.geom.contactPoint[0],i.geom.contactPoint[1],i.geom.contactPoint[2]))) + '\n')
inFile.close()

#####

Best regards,
Yousef

 Yousef Golabchi (yousef-golabchi) said on 2020-04-12: #3

Thanks Luc Scholtès, that solved my question.

 Yousef Golabchi (yousef-golabchi) said on 2020-04-16: #4

Hi Luc,

Could you please inform me how to obtain aperture at a given point (for instance injection point).

Should I have the cell IDs of that point or is it possible to do it by defining the coordinates of that point?

Best regards,
Yousef

 Luc Scholtès (luc) said on 2020-04-16: #5

Hi Youssef,

There is not direct solution. You need to identify the interaction which is the closest to the injection point (you can adapt the loop in the above script to identify the contactPoint coordinates x,y,z that are approximately equal to the injectionPoint coordinates x,y,z). Keep in mind that you won't get the aperture at the exact same position where the injection takes place since the injection is done in a pore and not at the interaction position.

Luc

 Yousef Golabchi (yousef-golabchi) said on 2020-04-16: #6

Many thanks for your swift response Luc.

Cheers,
Yousef

 Yousef Golabchi (yousef-golabchi) said on 2020-04-16: #7

Thanks Luc Scholtès, that solved my question.

 Luc Scholtès (luc) said on 2020-04-17: #8

Something like this should work to identify the interactions on the whole joint and the one at the injection point Xp:

intOnFracPlane=[]
intOnInjPoint=[]
for i in O.interactions:
if i.phys.isOnJoint:
intOnFracPlane.append(i)
if (i.geom.contactPoint[0]>(Xp[0]-Rmin)) and (i.geom.contactPoint[0]<(Xp[0]+Rmin)) and (i.geom.contactPoint[1]>(Xp[1]-Rmin)) and (i.geom.contactPoint[1]<(Xp[1]+Rmin)):
intOnInjPoint=i

Of course, you need to define Xp (the injection point) and to compute the minimum radius Rmin before that.

Luc

 Yousef Golabchi (yousef-golabchi) said on 2020-04-17: #9

Hi and good morning Luc,

Thank you so much for providing the code for identifying the interaction.

Please set me straight if I am wrong, what I understood yesterday from your command was that I can open one of the text files (that is produced by the command) and check first three coloums. If the coordinates are close to the injection point (as you pointed out it can never be exactly same) then the fourth coloums is the crack joint aperture of the injection point.

So how the code should be written is to check every row to see if X,Y, and Z are in the vicinity of the injection point coordinates (imposing a range) then if so, add the crack joint aperture (fourth coloum) and finally make an average.

Do I understand the general concept, right?

Thanks for your time and consideration.

Best regards,
Yousef

 Luc Scholtès (luc) said on 2020-04-17: #10

No, no need to work on the crack file. Just implement the lines I gave you after the first iteration of your simulation (where the bonds are created). Then, you can record the aperture at the injection point (and the corresponding pressure) inside a dedicated function called every n-iteration.

Like this one that I use to record "meta data" in a textfile that I can then process with a dedicated plotting script (in my case, a Python script):

intOnInjPoint=[]
def recordData():
global e10,e20,e30,intOnInjPoint
e1=triax.strain[0]-e10,e2=triax.strain[1]-e20,e3=triax.strain[2]-e30,
s1=triax.stress(triax.wall_right_id)[0],s2=triax.stress(triax.wall_top_id)[1],s3=triax.stress(triax.wall_front_id)[2],
p=flow.getPorePressure(Xp),
a=intOnInjPoint.phys.crackJointAperture,
tc=interactionLaw.nbTensCracks,
sc=interactionLaw.nbShearCracks,
unbF=utils.unbalancedForce()
)
plot.saveDataTxt(OUTPUT)

Just add the recordData() function in the engines list such as:

PyRunner(iterPeriod=10,initRun=False,command='recordData()',label='data'),

Luc