yade-users team mailing list archive

[Bruno Chareyre <question247806@xxxxxxxxxxxxxxxxxxxxx>]

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

Bruno Chareyre proposed the following answer:
@Ning,
I see two possible simple solutions to your problem:
One needs cpu time (see example below):
dtheta=0.01 #something small (aka rotation per timestep)
thetaDot=dtheta/O.dt
velGrad=Matrix3(0,thetaDot,0, -thetaDot,0,0, 0,0,0)
N=theta/dtheta
O.run(N) #now everything is rotated by theta

The other is instantaneaous: instead of rotating the packing, rotate the imposed velocity gradient in order to define it in the cells reference frame. For theta=45°:
q=Quaternion().setFromTwoVectors((1,0,0),(1,1,0))
velGrad=Matrix3(1,0,0,0,-1,0,0,0,0) #pure shear, written in the fixed cartesian frame
q.toRotationMatrix()*velGrad*q.inverse.toRotationMatrix()
Matrix3(-2.220446049250313e-16,1,0, 1,2.220446049250313e-16,0, 0,0,0)  #pure shear, in the cells reference frame

Actually, you may not have to rotate anything if the strain from FEM is
defined in an axisymmetric frame, since r,theta are always aligned on
the periodic cell axis...

@Vaclav
My impression is that W/2 comes out in simple shear (and so it is visible in many books - and correct), as the skew-symmetric part of L when L is:
0 W 
0 0
However, when the gradient is antisymmetric in the first place: 1/2(L+LT) = L, with L=
 0  W
-W 0
The 1/2 disappears, and the spin is just W.

Whatever it is, it will work when you look at a simulation. Even with
yade, without defining any spin, they rotate correctly in order to
maintain equilibrium (exception for a few spheres without contacts,
which clearly do not rotate correctly).

Example of rotating (modified from test/periodic-simple-shear.py):
============

from yade import pack,qt
O.periodic=True
O.cell.hSize=Matrix3(0.1,0,0, 0,0.1,0, 0,0,0.1)
radius=5e-3
sp=pack.SpherePack().makeCloud((0,0,0),(.1,.1,.1),5e-3,.2,500,periodic=True) 

num=sp.makeCloud((0,0,0),(.1,.1,.1),radius,.2,500,periodic=True) # min,max,radius,rRelFuzz,spheresInCell,periodic
O.bodies.append([sphere(s[0],s[1]) for s in sp])
O.engines=[
	ForceResetter(),
	InsertionSortCollider([Bo1_Sphere_Aabb()],verletDist=.05*radius),
	InteractionLoop(
		[Ig2_Sphere_Sphere_ScGeom()],
		[Ip2_FrictMat_FrictMat_FrictPhys()],
		[Law2_ScGeom_FrictPhys_CundallStrack()]
	),	PeriTriaxController(dynCell=True,mass=0.2,maxUnbalanced=0.01,relStressTol=0.02,goal=[-1e4,-1e4,0],stressMask=3,globUpdate=5,maxStrainRate=[1.,1.,1.],doneHook='triaxDone()',label='triax'),
	NewtonIntegrator(damping=.2),
]

phase=0
def triaxDone():
	global phase
	if phase==0:
		print 'Here we are: stress',triax.stress,'strain',triax.strain,'stiffness',triax.stiff
		print 'Now rotating.'
		O.cell.velGrad=Matrix3(0,1,0, -1,0,0, 0,0,0)
		triax.dead=1
		phase+=1

O.dt=PWaveTimeStep()
Gl1_Sphere.stripes=1

-- 
You received this question notification because you are a member of
yade-users, which is an answer contact for Yade.