The speed of run is too slow ! What are the ways to run the code faster?
Hello everyone,
I will be really appreciated it if anyone can help me with this question.
I am using Ubuntu 18.04, and Yade 2019-08-
I coded a Triaxial test based on [1] and developed it to store the micro + macro variables. The problem is it takes sooo long to be run when I define my main sample for the triaxial test. From last week, I run the code and it has been still running and I believe it would be running for more 2 weeks with respect to the finishing criterion. Do you know any way to make the running faster? because it is not useful in the way it is.
[1] https:/
My code:
%%%%%%%
#######
######### TRIAXIAL PROBLEM, Y IS THE VERTICAL AXIS, X IS THE RIGHT AXIS, Z IS THE FRONT AXIS #########
#######
# =======
# Please copy the file to the Desktop and then run the code.
# - to record the micro in the same iterations as the history() function, only modify the pyRunner iterPeriod : "PyRunner(
# - to change the micro iterPeriod, change the "micro_
# - to record the micro at each "ORN" iteration, basically remove the pyRunner (for the record micro function), and call "record_
# =======
print ('************** START **************')
import numpy as np
import time
import datetime, os
start_time=
from datetime import datetime
import math
from yade import qt, export, utils
from yade import pack
#######
######### DEFINING VARIABLES #########
print ('============ DEFINING VARIABLES ============')
nRead=readParam
num_spheres=20000,
compFricDegree = 29,
key='_triax_',
unknownOk=True
)
from yade.params import table
num_spheres=
key=table.key
targetPorosity = 0.4
compFricDegree = table.compFricD
finalFricDegree = 29
IP=100 # iteration period to record data and stuff
micro_record_
ORN=3000 # O.Run Number of iterations
micro_record_
micro_record_
damp=0.2
thick=0
stabilityThresh
PCPC=0.0001 # Precision of Confining Pressure Convergence
r_min=0.1*1e-3 # m
d_min=2*r_min # m
r_max=0.3*1e-3 # m
d_max=2*r_max # m
r_avr=(
d_avr=2*r_avr # m
r_fuz=(
Kn=10e8*(d_avr) ### FIXME
Kt=10e8*(d_avr) ### FIXME
young=Kn/r_avr # 2 (E r1 E r2 / E r1 + E r2) >>> E = Kn/r_avr
#young=5e6 # contact stiffness
poisson=Kn/Kt # Kt/Kn
#poisson=0.3 # Kt/Kn
Ls=0.02 # m length of specimen ### FIXME
L_REV=7*(d_avr) # m
if Ls < L_REV:
sys.exit("*** ERROR! The specimen's dimension is too samll! ***")
elif Ls==L_REV:
print ("*** This is the minimum specimen's dimension you can take! ***")
else:
print ("*** The specimen's dimension is good enough! ***")
mn,mx=Vector3(
Vt=-1*1e-3 # m/s # negative sign describes the compression direction
strainRate=Vt/Ls # %/sec
#strainRate=-0.01 # %/sec ### FIXME
target_strain=0.25 ### FIXME %
print ("The target strain has been set to:", target_strain)
sigmaIso=-5e5 # Pa ### FIXME
particleDensity
#particleDensit
#######
################# DEFINING MATERIALS #################
print ('============ DEFINING MATERIALS ============')
O.materials.
O.materials.
#######
################# DEFINING PACKING #################
print ('============ DEFINING PACKING ============')
walls=aabbWalls
for w in walls:w.
wallIds=
sp=pack.
clumps=False
sp.makeCloud(
#sp.makeCloud(
O.bodies.
from yade import export
os.mkdir(
os.chdir(
export.
#######
################# DEFINING TRIAXIAL TEST #################
print ('============ DEFINING TRIAXIAL TEST ============')
triax=TriaxialS
maxMultiplier=
finalMaxMultip
thickness = thick,
stressMask = 7,
internalCompac
)
#######
################# DEFINING FUNCTIONS #################
print ('============ DEFINING FUNCTIONS ============')
from yade import plot
def history():
plot.addData(
e11 = -triax.strain[0],
e22 = -triax.strain[1],
e33 = -triax.strain[2],
ev = -triax.
s11 = -triax.
s22 = -triax.
s33 = -triax.
i = O.iter,
t = O.time, # virtual (yade) time --- time of simulation
fab = utils.fabricTen
#FK: Micro data recording
#NOTE: Data structure in hdf5 is:
# hdf5-file/
# ├── bodies/
# | ├── positions(3)
# | ├── radius(1)
# └── interactions/
# ├── ids(2)
# ├── normalForce(3)
# └── shearForce(3)
#IMPORTANT NOTE: the array indexation in bodies/positions and bodies/radius is the same as in interactions/ids : it is the Yade indexation (O.bodies[
#It means that the values in interactions/ids hdf5 array can be used as index in bodies/positions and bodies/radius.
import h5py # works with ORN and macro period of recording
def record_
if not os.path.
filename=
with h5py.File(
#Record time:
#Create a "folder" for bodies:
#Add bodies positions
#Add bodies radius. Radius=0 means its a wall
#Create a "folder" for interactions:
#Prepare the output array, performance is way better when you set its size a-priori.
out_array = np.empty(
#Fill the output_array
for i,I in enumerate(
if micro_record_
with h5py.File(
for i,I in enumerate(
#######
################# DEFINING ENGINES #################
print ('============ DEFINING ENGINES ============')
newton=
O.engines=[
ForceResetter(),
#GravityEng
InsertionSortC
InteractionLoop(
[Ig2_
[Ip2_
[Law2_
),
## We will use the global stiffness of each body to determine an optimal timestep (see https:/
GlobalStiffnes
triax,
PyRunner(
#PyRunner(
TriaxialStateR
newton
]
Gl1_Sphere.
if nRead==0: yade.qt.
#######
################# APPLYING CONFINING PRESSURE #################
print ('============ APPLYING CONFINING PRESSURE ============')
triax.goal1 = sigmaIso
triax.goal2 = sigmaIso
triax.goal3 = sigmaIso
while 1:
O.run(ORN,True)
unb = unbalancedForce()
meanS=
ConfStressR
print ('~~~~~~~~~~~~~ Phase_01: Converging to Isotropic Compression, 50kPa ~~~~~~~~~~~~~')
print ('mean stress engine:
print ('porosity:
if unb<stabilityTh
break
export.
e22Check=
print ('Axial Strain',e22Check)
print ('Mean stress engine: ',triax.meanStress)
print ('Mean stress (Calculated): ',meanS)
print ('#####
#######
################# REACHING TARGET POROSITY ##################
print ('============ REACHING TARGET POROSITY ============')
import sys
while triax.porosity>
compFricDegree = 0.95*compFricDegree
setContactFric
print ('\r Friction: ',compFricDegree,' porosity:
sys.stdout.flush()
O.run(ORN,True)
print ('#####
#######
################ PRINT SOME CHECK VARIABLES #################
RRmax=max(
RRmin=min(
print('Maximum Radius:',RRmax)
print('Minimum Radius:',RRmin)
print ('Number of elements:', len(O.bodies))
print ('Box Volume engine:', triax.boxVolume)
Vt = (mx[0]-
print ('Box Volume calculated:', Vt)
if Vt == triax.boxVolume:
print ("*** Volume calculation is Correct. ***")
else:
sys.exit("*** ERROR! Volume calculation is WRONG. ***")
Vs=triax.
print('Total volume of particles (Vs):',Vs)
Vv=Vt-Vs
print('Total volume of voids Calculated (Vv):',Vv)
print ('porosity:
n=Vv/Vt
print ('porosity Calculated (n):',n)
if n == triax.porosity:
print ("*** Porosity calculation is Correct. ***")
e=n/(1-n)
print ('Void ratio Calculated (e):',e)
else:
sys.exit("*** ERROR! Porosity calculation is WRONG. ***")
print ('step that starts the deviatoric loading ', O.iter/ORN)
#######
################# DEVIATORIC LOADING #################
print ('============ APPLYING DEVIATORIC LOADING ============')
triax.internalC
setContactFrict
# Allow walls to move. >>> False means fixed , True means free to move
triax.wall_
triax.wall_
triax.wall_
triax.wall_
triax.wall_
triax.wall_
# Make Stress Control on X,Y,Z directions by setting True, False means Strain Control
triax.stressCon
triax.stressCon
triax.stressCon
# set the raite of strain to the desirable axis
#triax.
triax.strainRat
#triax.
# set the deviatoric conditions
triax.stressMask = 5
triax.goal1 = sigmaIso
triax.goal2 = strainRate
triax.goal3 = sigmaIso
newton.damping=0.1 ### FIXME
os.mkdir(
os.chdir(
#micro_
while 1:
O.run(ORN,True)
export.
record_
unb=
axialS=
eps2=
eps2cal=
print ('~~~~~~~~~~~~~ Phase_02: Converging to Deviatoric Compression, Strain Rate ~~~~~~~~~~~~~')
print ('sigma2: ',axialS, ' q: ', axialS-sigmaIso,' step= ', O.iter/ORN,' Time:',O.time, ' TimeStep',O.dt)
print ('Strain2 Calculated: ',eps2cal, ' Strain2 engine: ',eps2,'axial deformation (%)', (eps2-e22Check)
print('A:', abs(eps2-e22Check), ' B:',target_strain)
if abs(eps2-
break
print ('#####
#######
################# END #################
print ('************** END **************')
O.realtime
print ('Analysis has been taken for',O.realtime, 'seconds or', O.realtime/60, 'minutes')
print('Real time of run:',(time.time() - start_time), 'seconds or',(time.time() - start_time)/60, 'minutes')
#######
################# RECORD Macro DATA and Plot #################
print ('============ RECORD AND PLOT DATA ============')
O.run(ORN,True) # Show the deformed shape of specimen, otherwise the specimen shows no deformation schematically
# =======
# plot.plots=
# plot.labels=
# plot.plot()
# =======
os.chdir(
plot.saveDataTx
from matplotlib import pyplot as plt
e22 = plot.data["e22"]
ev = plot.data["ev"]
s22 = plot.data["s22"]
s11 = plot.data["s11"]
s33 = plot.data["s33"]
i = plot.data["i"]
t = plot.data["t"]
fig1, ax1=plt.
fig2, ax2=plt.
fig3, ax3=plt.
fig4, ax4=plt.
ax1.plot(
ax1.legend()
#ax1.set_
ax1.set_
ax1.set_
ax1.grid(True)
fig1.tight_layout()
fig1.savefig(
ax2.plot(
ax2.legend()
#ax2.set_
ax2.set_
ax2.set_
ax2.grid(True)
fig2.tight_layout()
fig2.savefig(
ax3.plot(
ax3.legend()
#ax3.set_
ax3.set_
ax3.set_
ax3.grid(True)
fig3.tight_layout()
fig3.savefig(
ax4.plot(
ax4.plot(
ax4.legend()
#ax4.set_
ax4.set_
ax4.set_
ax4.grid(True)
fig4.tight_layout()
fig4.savefig(
plt.show()
O.pause()
%%%%%%%
Thank you,
Ehsan
Question information
- Language:
- English Edit question
- Status:
- Solved
- For:
- Yade Edit question
- Assignee:
- No assignee Edit question
- Solved by:
- ehsan benabbas
- Solved:
- Last query:
- Last reply: