Question #226352 on Yade changed:
https://answers.launchpad.net/yade/+question/226352
Description changed to:
Dear Jan Stránský:
Thanks for helping me to finish this task and I can accomplish all the
process of 2D biaxial compression simulation, from preprocessing to
postprocessing, using this script that is shown as follows for others to
reference.
### fundamental details of application ###
# unicode: UTF-8
Filename='2D-simulation'
from yade import pack,os
################################# preprocessing for simulation
##########################################
### prescribing variables and functions for simulation controller ###
# material defination
spheremat =
O.materials.append(ViscElMat(kn=4e6,ks=4e6,cn=.0,cs=.0,density=1500,frictionAngle=25.565))
wallmat =
O.materials.append(ViscElMat(kn=1e8,ks=1e8,cn=.0,cs=.0,density=2600,frictionAngle=25.565))
# walls defination
mn,mx=Vector3(0,0,0),Vector3(7,7,1)
wallIds=O.bodies.append(utils.aabbWalls([mn,mx],thickness=.01,material=wallmat))
# ThreeDTriaxialEngine defination for initial-state determination(the first
calculation step)
triax01=ThreeDTriaxialEngine(
wall_bottom_id=wallIds[2],wall_top_id=wallIds[3],
wall_left_id=wallIds[0],wall_right_id=wallIds[1],
wall_back_id=wallIds[4],wall_front_id=wallIds[5],
wall_front_activated = False,wall_back_activated = False,
internalCompaction=False,
stressControl_1 = True, stressControl_2 = True,stressControl_3 = True,
computeStressStrainInterval =10,
sigma_iso = 1.25e5,
sigma1 = 1.25e5,
sigma2 = 1.25e5,
sigma3 = 1.25e5,
strainRate1 = 0.01,strainRate2 = 0.01,
)
# ThreeDTriaxialEngine defination for biaxial compression(the second
calculation step)
triax02=ThreeDTriaxialEngine(
wall_bottom_id=wallIds[2],wall_top_id=wallIds[3],
wall_left_id=wallIds[0],wall_right_id=wallIds[1],
wall_back_id=wallIds[4],wall_front_id=wallIds[5],
wall_front_activated = False,wall_back_activated = False,
internalCompaction=False,
stressControl_1 = True, stressControl_2 = False,stressControl_3 = True,
computeStressStrainInterval =10,
sigma_iso = 1.25e5,
sigma1 = 1.25e5,
sigma2 = 1.25e5,
# sigma3 = 1.25e5,
strainRate1 = 0.001,strainRate2 = 0.15,
)
# Simulation stop conditions defination
def checkUnbalanced():
unb=unbalancedForce()
meanS=(triax01.stress(triax01.wall_right_id)[0]+triax01.stress(triax01.wall_top_id)[1])/2
q=unb
r=abs(meanS-triax01.sigma_iso)/triax01.sigma_iso
if q<0.01 and r<1e-5:
O.pause()
O.save('first-step.xml'.format(Filename))
################################# control flow for simulation
##########################################
# particles generation
O.periodic=1
O.cell.setBox(7,7,7)
sp=pack.SpherePack()
sp.makeCloud((0,0,.5),(7,7,.5),rMean=-1,rRelFuzz=0,num=800,periodic=True)
sp.toSimulation(material=spheremat)
# determining colors for particles in different aeras of the cell
for b in O.bodies:
if isinstance(b.shape,Sphere):
pos = b.state.pos
if pos[0] <3.5 and pos[1] < 3.5: b.shape.color = (1,0,0) # area 1
elif pos[0] >= 3.5 and pos[1] <3.5: b.shape.color = (0,1,0) # area 2
elif pos[0] >= 3.5 and pos[1] >= 3.5: b.shape.color = (0,0,1) # area 3
else: b.shape.pos = (1,1,0) # area 4
O.periodic=0
# blockedDOFs
for b in O.bodies:
if isinstance(b.shape,Sphere):
b.state.blockedDOFs='zXY'
# Simulation assembly for the first step
O.engines=[
ForceResetter(),
InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb(),Bo1_Wall_Aabb()]),
InteractionLoop(
[Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom(),Ig2_Wall_Sphere_ScGeom()],
[Ip2_ViscElMat_ViscElMat_ViscElPhys()],
[Law2_ScGeom_ViscElPhys_Basic()]
),
GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.8),
triax01,
NewtonIntegrator(damping=.1),
PyRunner(command='checkUnbalanced()',iterPeriod=200)
]
# first step of simulation startting with a correct inheriting for the next step
O.dt = 2e-4
O.run(); O.wait()
# loading inheriting
O.load('first-step.xml')
# Simulation assembly for the second step
O.engines=[
ForceResetter(),
InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb(),Bo1_Wall_Aabb()]),
InteractionLoop(
[Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom(),Ig2_Wall_Sphere_ScGeom()],
[Ip2_ViscElMat_ViscElMat_ViscElPhys()],
[Law2_ScGeom_ViscElPhys_Basic()]
),
GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.8),
triax02,
NewtonIntegrator(damping=.1),
]
# second step of simulation startting
O.dt = 2e-4
O.run(20000,True);
# whole task over
O.save('final-step.xml'.format(Filename));
O.wait()
################################## postprocessing for simulation
######################################################
f = file("/home/fzg/fu/result.dat",'w')
f.write('# This is the result data of 2D simulation\n\n')
f.write('# There are 8 types of varibles in this data as follows:\n\n')
f.write("varibles='X-cordinate','Y-cordinate','Z-cordinate','Radius','X-displacement','Y-displacement','Z-displacement','Ids'\n\n")
f.write('%-16s %-16s %-16s %-16s %-16s %-16s %-16s %-16s\n'%
('X-cordinate','Y-cordinate','Z-cordinate','Radius','X-displacement','Y-displacement','Z-displacement','Ids'))
for b in O.bodies:
if isinstance(b.shape,Sphere):
pos = b.state.pos
rad = b.shape.radius
displ = b.state.displ()
f.write('%-16g %-16g %-16g %-16g %-16g %-16g %-16g
%-16d\n'%(pos[0],pos[1],pos[2],rad,displ[0],displ[1],displ[2],b.id))
f.write('################################ ends
##########################################')
f.close()
def rename():
global Filename
os.rename("/home/fzg/fu/result.dat","/home/fzg/fu/{0}.plt".format(Filename))
rename()
After this script showing, I have two problems to ask you as that:
(1). This script can almost do all the things the simulation need to do except
the module "determining colors for particles in different aeras of the cell".
After using these cammands, the colors of area 1 to area 3(only one color for
each area) are differently assigned but area 4 has no reaction to this
requirement that it comes in a variety of colors. I don not know why and desire
your suggestion.
(2).The script can run well but is too long and inconvenient, but now I have
almost no experience to optimize it using Python, please give me some advises
for this purpose.
SEEKING YOUR HELP!
--
You received this question notification because you are a member of
yade-users, which is an answer contact for Yade.
_______________________________________________
Mailing list: https://launchpad.net/~yade-users
Post to : [email protected]
Unsubscribe : https://launchpad.net/~yade-users
More help : https://help.launchpad.net/ListHelp