yade-users team mailing list archive

[James X <question305894@xxxxxxxxxxxxxxxxxxxxx>]

Question #305894 on Yade changed:
https://answers.launchpad.net/yade/+question/305894

James X posted a new comment:
Hi Jan,

Thank you very much. I have included a simplified version of what I am
trying to do.  It is currently a pseudo-2D simulation, though I would
ideally like to be able to handle both 3D and 2D. What I would hope for
is to be able to export, at any user-defined timestep, the fracture
event coordinates, the coordinates of the element pair that broke, and
the material of the two elements.

However, having just the coordinates of the fracture event themselves
itself would be great.

Just as a FYI, the first half is just a 2D gravity-settling script to
generate a 2D mesh. While I would like to do 3D simulations, the
limitations of my experimental data does not allow me to go 3D at the
moment.


import numpy as np
import time

from yade import pack
from yade import geom
from yade import polyhedra_utils
from yade import utils
from yade import linterpolation
from yade import export
from yade import ymport

#other random parameters
start = time.time()

#this function just makes a 2D sphere mesh and exports it as <FILENAME>.sphere
def Make_Packed_file(filename, xSize, ySize, radius, rRelFuzz_value, nbIter=3000):
    #filename = filename #output a msh file
    sp=pack.SpherePack()
    sp.makeCloud(minCorner=(0,0,0),maxCorner=(xSize,ySize*2.5,0),rMean=radius, rRelFuzz = 0)
    sp.toSimulation()
    for b in O.bodies:
        b.state.blockedDOFs = 'zXY'
    box=O.bodies.append(utils.geom.facetBox((xSize/2,1.25*ySize,0),(xSize/2,1.25*ySize,radius),wallMask=63)) # the same box of makeCloud

    O.engines=[
        ForceResetter(),
        InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()]),
        InteractionLoop(
        [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom(),Ig2_Wall_Sphere_ScGeom()],
        [Ip2_FrictMat_FrictMat_FrictPhys()],
        [Law2_ScGeom_FrictPhys_CundallStrack()]
        ),
        NewtonIntegrator(gravity=(0,-2.0,0),damping=.1),
        PyRunner(iterPeriod=500,initRun=False,command='print("elapsed time is: " + str(time.time() - start) ) '),
        DomainLimiter(lo=(-0.01,-0.01,-0.01),hi=(xSize,ySize,1.1*radius),iterPeriod=nbIter-1)
    ]
    O.dt=utils.PWaveTimeStep()
    O.stopAtIter=nbIter
    O.run()
    O.wait()
    export.text(filename+'.spheres')


###random parameters
radius = 3
smoothContact=True
jointFrict=radians(20)
jointDil=radians(0)
new_scale = 1.0 
xSize = 250;
ySize = 150;
rRelFuzz_value = 0;

#actually make this packed mesh
Make_Packed_file("simplified_packed_2D_mesh", xSize, ySize, radius, rRelFuzz_value)

O.reset() #reset everything, just in case


#defin the material
def mat_1(): return JCFpmMat(type=1,young=30.0e9,poisson=0.3,frictionAngle=radians(30),density=3000,tensileStrength=1.23e8,cohesion=1.23e8,jointNormalStiffness=1e7,jointShearStiffness=1e7,jointCohesion=1e6,jointFrictionAngle=radians(20),jointDilationAngle=0.0) 


O.materials.append((mat_1()))

#import the mesh created in the previous simulation
import_spheres = ymport.text("simplified_packed_2D_mesh.spheres",scale=new_scale,shift=Vector3(0,0,0),material=mat_1)
O.bodies.append(import_spheres)


#add in wall information, from uniax.py
young = 90e9
mn,mx=Vector3(0,0,0),Vector3(xSize*new_scale,1.01*ySize*new_scale,2.00001*radius) # corners of the initial packing
O.materials.append(FrictMat(young=young,poisson=0.5,frictionAngle=0,density=0,label='walls'))
walls=aabbWalls([mn,mx],thickness=0,material='walls')
wallIds=O.bodies.append(walls)


#from uniax.py
bb=uniaxialTestFeatures()
negIds,posIds,axis,crossSectionArea=bb['negIds'],bb['posIds'],bb['axis'],bb['area']

strainRateTension=.05
setSpeeds=True


interactionRadius=1.60 # to set initial contacts to larger neighbours and a bit further

O.engines=[

	ForceResetter(),
	InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=interactionRadius,label='is2aabb'),]),
	InteractionLoop(
		[Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=interactionRadius,label='ss2d3dg')],
		[Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(cohesiveTresholdIteration=1,label='interactionPhys')],
		[Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(smoothJoint=smoothContact,label='interactionLaw')]
	),
	GlobalStiffnessTimeStepper(timestepSafetyCoefficient=0.8),
	NewtonIntegrator(damping=0.2,gravity=(0.,0,0.)),
        UniaxialStrainer(strainRate=strainRateTension,axis=axis,asymmetry=0,posIds=posIds,negIds=negIds,crossSectionArea=crossSectionArea,blockDisplacements=False,blockRotations=False,setSpeeds=setSpeeds,label='strainer'),
]

#### time step definition (low here to create cohesive links without big changes in the assembly)
O.dt=0.05*utils.PWaveTimeStep()

#### set cohesive links with interaction radius>=1
O.step();

#### initializes now the interaction detection factor to strictly 1
ss2d3dg.interactionDetectionFactor=-1.0
is2aabb.aabbEnlargeFactor=-1.0

-- 
You received this question notification because your team yade-users is
an answer contact for Yade.