yade-users team mailing list archive
-
yade-users team
-
Mailing list archive
-
Message #27078
[Question #700561]: Convert script to Openmpi
New question #700561 on Yade:
https://answers.launchpad.net/yade/+question/700561
Dear all,
I am new on running Yade on multiple nodes using openmpi.
Still checking how to convert the code below to run in mpi.
Runnning the code below I get the messsage:
Friction: 1 porosity: 1.0 None
[95mMaster: single-core, fall back to O.run() [0m
Probably I am not using O.run and mpirun correctly.
If anyone with more experimence see what I am duing wrong, please let me know.
Cheers,
#!/usr/bin/python
# -*- encoding=utf-8 -*-
#*************************************************************************
# Copyright (C) 2010 by Bruno Chareyre *
# bruno.chareyre_at_grenoble-inp.fr *
# *
# This program is free software; it is licensed under the terms of the *
# GNU General Public License v2 or later. See file LICENSE for details. *
#*************************************************************************/
import os
from yade import mpy as mp
from yade import pack
from yade import bodiesHandling
from yade import export
from yade import utils
from yade import ymport
import math
############################################
### DEFINING VARIABLES AND MATERIALS ###
############################################
# The following 5 lines will be used later for batch execution
nRead=readParamsFromTable(
num_spheres=3000,# number of spheres
compFricDegree = 1, # contact friction during the confining phase
key='_triax_base_', # put you simulation's name here
unknownOk=True
)
from yade.params import table
num_spheres=table.num_spheres# number of spheres
key=table.key
targetPorosity = 0.35 #the porosity we want for the packing
compFricDegree = table.compFricDegree # initial contact friction during the confining phase (will be decreased during the REFD compaction process)
finalFricDegree = 30 # contact friction during the deviatoric loading
rate=0 # loading rate (strain rate)
damp=0.8 # damping coefficient
stabilityThreshold=0.01 # we test unbalancedForce against this value in different loops (see below)
#2e4+70e4medio 1e4+70e4bom 1e4+60e4bom 3e4+90e4+w3,1,-1-the best
young=20e5 # contact stiffness200e4
young2=20e5
youngcoat=20e5
bondstr=0.3e3#2e7
bondstr2=0.3e3
bondstrcoat=1e3
## create materials for spheres and plates
mat=O.materials.append(JCFpmMat(type=1,young=young,poisson=0.3,frictionAngle=radians(compFricDegree),density=2000,tensileStrength=bondstr,cohesion=bondstr,jointNormalStiffness=0,jointShearStiffness=0,jointCohesion=bondstr,jointFrictionAngle=radians(0),jointDilationAngle=0.0,label='spheres'))
O.materials.append(JCFpmMat(type=1,young=20e7,poisson=0.3,frictionAngle=radians(0),density=2600,tensileStrength=0,cohesion=0,jointNormalStiffness=0,jointShearStiffness=0,jointCohesion=0,jointFrictionAngle=radians(0),jointDilationAngle=0.0,label='walls'))
O.materials.append(JCFpmMat(type=1,young=youngcoat,poisson=0.3,frictionAngle=radians(1),density=1500,tensileStrength=bondstrcoat,cohesion=bondstrcoat,jointNormalStiffness=0,jointShearStiffness=0,jointCohesion=bondstrcoat,jointFrictionAngle=radians(0),jointDilationAngle=0.0,label='spherescoat'))
## create walls around the packing
mn,mx=Vector3(0,0,0),Vector3(0.0101,0.0101,0.0101)
mnbox,mxbox=Vector3(0,0,0),Vector3(0.0101,0.0115,0.0101)
walls=aabbWalls([mnbox,mxbox],thickness=0,material='walls')
wallIds=O.bodies.append(walls)
O.bodies.append(ymport.textExt("matrix_vtest.txt", format='x_y_z_r', shift=Vector3(0,0,0), scale=1.0,material='spherescoat',color=(0,1,1)))
################Particle substitution by large aggregate######################################################################
bodid=[]
a=[]
for b in O.bodies:
if b and isinstance(b.shape,Sphere):
# print (b.shape.radius)
if b.shape.radius==0.0005:
bodid.append(b.id)
a.append(b.state.pos)
i=0
for p in bodid:
t=a[i]
f1=O.bodies.append(ymport.textExt("agg5e4_5e5.txt", format='x_y_z_r', shift=t-Vector3(0,0,0), scale=1.0,material='spheres',color=(0,1,1)))
O.bodies.erase(bodid[i])
i=i+1
bodidd=[]
aa=[]
for bb in O.bodies:# in sp:
if bb and isinstance(bb.shape,Sphere):
# print (bb.shape.radius)
if bb.shape.radius==0.0003:
bodidd.append(bb.id)
aa.append(bb.state.pos)
ii=0
for pp in bodidd:
tt=aa[ii]
f2=O.bodies.append(ymport.textExt("agg3e4_3e5.txt", format='x_y_z_r', shift=tt-Vector3(0,0,0), scale=1.0,material='spheres',color=(0,1,1)))
O.bodies.erase(bodidd[ii])
ii=ii+1
bodiddd=[]
aaa=[]
for bbb in O.bodies:# in sp:
if bbb and isinstance(bbb.shape,Sphere):
# print (bbb.shape.radius)
if bbb.shape.radius==0.0002:
bodiddd.append(bbb.id)
aaa.append(bbb.state.pos)
iii=0
for ppp in bodiddd:
ttt=aaa[iii]
f3=O.bodies.append(ymport.textExt("agg2e4_2e5.txt", format='x_y_z_r', shift=ttt-Vector3(0,0,0), scale=1.0,material='spheres',color=(0,1,1)))
O.bodies.erase(bodiddd[iii])
iii=iii+1
##############################################################################################################################
#====================================================================================================================================================================
############################
### DEFINING ENGINES ###
############################
triax=TriaxialStressController(
## TriaxialStressController will be used to control stress and strain. It controls particles size and plates positions.
## this control of boundary conditions was used for instance in http://dx.doi.org/10.1016/j.ijengsci.2008.07.002
maxMultiplier=1.+2e4/young, # spheres growing factor (fast growth)
finalMaxMultiplier=1.+2e3/young, # spheres growing factor (slow growth)
thickness = 0,
## switch stress/strain control using a bitmask. What is a bitmask, huh?!
## Say x=1 if stess is controlled on x, else x=0. Same for for y and z, which are 1 or 0.
## Then an integer uniquely defining the combination of all these tests is: mask = x*1 + y*2 + z*4
## to put it differently, the mask is the integer whose binary representation is xyz, i.e.
## "100" (1) means "x", "110" (3) means "x and y", "111" (7) means "x and y and z", etc.
stressMask = 0,
internalCompaction=False, # If true the confining pressure is generated by growing particles
wall_front_activated=True,
wall_back_activated=True,
wall_top_activated=True,
wall_bottom_activated=True,
wall_left_activated=True,
wall_right_activated=True,
goal1=-5,
goal2=-10,
goal3=-5,
)
newton=NewtonIntegrator(damping=damp)
O.engines=[
ForceResetter(),
InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()]),
InteractionLoop(
[Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
[Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(cohesiveTresholdIteration=-1,label='Physicspheres')],#,xSectionWeibullShapeParameter=1.5, xSectionWeibullScaleParameter=1
[Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(smoothJoint=False)]
),
GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.3),
triax,
# VTKRecorder(iterPeriod=2000,recorders=['all','cracks'],fileName="/home/user/Ãrea de Trabalho/Paper_Biopore/SandBox_Laser/Video/Videotest"),
# newton
# NewtonIntegrator(damping=0.4,gravity=[0,0,0]),
newton
]
#Display spheres with 2 colors for seeing rotations better
#Gl1_Sphere.stripes=0
#if nRead==0: yade.qt.Controller(), yade.qt.View()
###################################################
### REACHING A SPECIFIED POROSITY PRECISELY ###
###################################################
## We will reach a prescribed value of porosity with the REFD algorithm
## (see http://dx.doi.org/10.2516/ogst/2012032 and
## http://www.geosyntheticssociety.org/Resources/Archive/GI/src/V9I2/GI-V9-N2-Paper1.pdf)
import sys #this is only for the flush() below
#poro=utils.voxelPorosityTriaxial(triax, 200, offset=0.0018)
while triax.porosity>targetPorosity:
# we decrease friction value and apply it to all the bodies and contacts
# compFricDegree = 0.95*compFricDegree
setContactFriction(radians(compFricDegree))
print ("\r Friction: ",compFricDegree," porosity:",triax.porosity,
sys.stdout.flush())
# while we run steps, triax will tend to grow particles as the packing
# keeps shrinking as a consequence of decreasing friction. Consequently
# porosity will decrease
mp.mpirun(500,1)
#O.save('compactedState'+key+'.yade.gz')
mp.mprint("### Compacted state saved ###")
#print(triax.stress(3)[1])
##############################
### DEVIATORIC LOADING ###
##############################
#We move to deviatoric loading, let us turn internal compaction off to keep particles sizes constant
#triax.internalCompaction=False
# Change contact friction (remember that decreasing it would generate instantaneous instabilities)
setContactFriction(radians(finalFricDegree))
#O.materials.append(JCFpmMat(type=1,young=young,poisson=0.3,frictionAngle=radians(finalFricDegree),density=2000,tensileStrength=bondstr2,cohesion=bondstr2,jointNormalStiffness=0,jointShearStiffness=0,jointCohesion=bondstr2,jointFrictionAngle=radians(0),jointDilationAngle=0.0,label='spheres'))
#=======================================
#set stress control on x and z, we will impose strain rate on y
triax.stressMask = 2
triax.wall_bottom_activated=0
#now goal2 is the target strain rate
triax.goal1=rate#triax.stress(1)[0]
triax.goal3=rate#triax.stress(5)[2]
triax.goal2=triax.stress(3)[1]
##Save temporary state in live memory. This state will be reloaded from the interface with the "reload" button.
#O.saveTmp()
#####################################################
### Example of how to record and plot data ###
#####################################################
#from yade import plot
from yade import plot
mp.mpirun(10,True)
#strain is logarithmic strain or true strain which is ls=(ln1+e) where e=dl/L (strain)
ei0=-triax.strain[0];ei1=-triax.strain[1];ei2=-triax.strain[2]
si0=-triax.stress(0)[0];si1=-triax.stress(2)[1];si2=-triax.stress(4)[2]
## a function saving variables
def history():
plot.addData(e11=-triax.strain[0]-ei0, e22=-triax.strain[1]-ei1, e33=-triax.strain[2]-ei2,
ev=-triax.strain[0]-triax.strain[1]-triax.strain[2],
s11=-triax.stress(triax.wall_right_id)[0]-si0,
s22=-triax.stress(triax.wall_top_id)[1]-si1,
s33=-triax.stress(triax.wall_front_id)[2]-si2,
e=math.exp(-triax.strain[1]-ei1)-1,
pc=-unsat.bndCondValue[2],
sw=unsat.getSaturation(False),
z1=O.bodies[3].state.pos[1],
i=O.iter)
if 1:
## include a periodic engine calling that function in the simulation loop
O.engines=O.engines[0:5]+[PyRunner(iterPeriod=20,command='history()',label='recorder')]+O.engines[5:7]
#plot.plots={'pc':('sw',None,'e22')}
#plot.plot()
#print("test",math.exp(2))
#######################################################
## Drainage Test under oedometer conditions ###
#######################################################
##Instantiate a two-phase engine
unsat=TwoPhaseFlowEngine()
#meanDiameter=(O.bodies[-1].shape.radius + O.bodies[6].shape.radius) / 2.
##set boundary conditions, the drainage is controlled by decreasing W-phase pressure and keeping NW-phase pressure constant
unsat.bndCondIsPressure=[0,0,1,1,0,0]
unsat.bndCondIsWaterReservoir=[0,0,1,0,0,0]
unsat.bndCondValue=[0,0,-1e8,0,0,0]
unsat.isPhaseTrapped=True #the W-phase can be disconnected from its reservoir
unsat.initialization()
unsat.surfaceTension = 0.0728
#start invasion, the data of normalized pc-sw-strain will be written into pcSwStrain.txt
f5=open('SwPc3.txt',"w")
ts=O.dt
pgstep= 45000000*ts #30Pa/s
print (pgstep)
pgmax= 9000#9316 #Pa
for pg in arange(1.0e-8,pgmax,pgstep):
unsat.bndCondValue=[0,0,(-1.0)*pg,0,0,0]
unsat.invasion()
unsat.computeCapillaryForce()
unsat.meshUpdateInterval=500
unsat.defTolerance=-1
unsat.updateTriangulation=True
print(unsat.getSaturation(False),pg,-triax.strain[1])
for b in O.bodies:
O.forces.setPermF(b.id, unsat.fluidForce(b.id))
while 1:
mp.mpirun(100,True)
unb=unbalancedForce()
if unb<0.1:
break
f5.write(str(pg)+" "+str(unsat.getSaturation(False))+" "+str(triax.strain[1])+"\n")
f5.close()
plot.plots={'pc':('sw',None,'e22')}
plot.plot(noShow=True).savefig('Fig.png')
mp.ERASE_REMOTE_MASTER = True #keep remote bodies in master?
mp.DOMAIN_DECOMPOSITION= True #automatic splitting/domain decomposition
mp.mpirun(20000) #passive mode run
mp.MERGE_W_INTERACTIONS = False
#mp.mpirun(NSTEPS,numThreads,withMerge=True) # interactive run, numThreads is the number of workers to be initialized, see below for withMerge explanation.
mp.mergeScene() #merge scene after run.
--
You received this question notification because your team yade-users is
an answer contact for Yade.
