yade-users team mailing list archive
-
yade-users team
-
Mailing list archive
-
Message #24334
[Question #693898]: Bonding particles with JCFpm yields unexpected forces
New question #693898 on Yade:
https://answers.launchpad.net/yade/+question/693898
Dear Yade community,
being new to yade I lack some basic understanding about the constitutitve law ‘Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM’ and related JCFpmMat, Ip2_JCFpmMat_JCFpmMat_JCFpmPhys classes:
Problem Formulation:
Question 635871 [1] states that the JCFpm law is an implementation of Potyondy & Cundalls ‘A bonded particle model for rock’ [3]. I would therefore expect a Hertz Mindlin contact force [2. eq: 33-36], updated by a force proportional to a set bond normal (and shear) stiffness and the cross section A [3. eq: 14], as proposed by Potyondy & Cundall.
The following minimal test example with two spheres, each customly set ‘onJoint=True’ and the interaction set ‘isOnJoint=True’, shows that the normal force is proportional to: 2*young*(r1*r2/(r1+r2))*relativeOverlap. The stiffness kn is not calculated with the JCFpmMat parameter jointNormStiffness as expected from the source code JointedCohesiveFrictionalPM.ccp - line: 560-580.
Minimal example:
mat=JCFpmMat(young=3.e10,poisson=0.25,density=5000.,frictionAngle=radians(10.),
tensileStrength=1.e8,cohesion=1.e8,jointNormalStiffness=2.e8,jointShearStiffness=2.e8,
jointCohesion=1.e8,jointFrictionAngle=radians(10.), label='Friction')
#Creating Spheres
s1 = utils.sphere([0.,0.,0.],radius=1,color=(1,0,1),material=mat,fixed=True)
s2 = utils.sphere([0.,0.,2.5],radius=1,color=(1,0,1),material=mat)
s2.state.vel = (0,0,-1.)
O.bodies.append([s1,s2])
O.engines=[
ForceResetter(),
InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=1.5,label='bo1s'),]),
InteractionLoop(
[Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=1.5,label='ig2ss'),],
[Ip2_JCFpmMat_JCFpmMat_JCFpmPhys(label='interactionPhys')],
[Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM(label='interactionLaw')]
),
NewtonIntegrator(damping=0.),
]
O.dt = .05*utils.PWaveTimeStep()
O.step()
bo1s.aabbEnlargeFactor = 1.
ig2ss.interactionDetectionFactor = 1.
#Set sphere state and interactions manually 'onJoint'
for i in O.interactions:
O.bodies[i.id1].state.joint += 1
O.bodies[i.id2].state.joint += 1
O.bodies[i.id1].state.onJoint, O.bodies[i.id2].state.onJoint = True, True
i.phys.isOnJoint = True
#Check if force is related to jointNormalStiffness:
O.run(10)
i = O.interactions[0,1]
print('isOnJoint = ',i.phys.isOnJoint)
print('normalStifness=',i.phys.kn, ' normalForce=',i.phys.normalForce, 'ref force=', (i.geom.penetrationDepth-i.phys.initD)*mat.young)
Output:
>>> isOnJoint = ', True, normalStifness=, 30000000000.0, normalForce=, Vector3(0,0,6063592.219396341), ref force=, 6063592.219396341
Questions:
1) Is it necessary to custom set every interaction ‘isOnJoint’ when forces according to jointNormalStiffness / jointShearStiffness are desired?
2) How does one achieve forces that are guided by the cross-section A and a defined jointNormalStiffness?
3) Similar to the Mindlin Physics Model, is it possible to integrate viscous damping with a coefficient of restitution en
(Ip2_FrictMat_FrictMat_MindlinPhys(en=en))?
[1]: [https://answers.launchpad.net/yade/+question/635871]
[2]: [https://onlinelibrary.wiley.com/doi/epdf/10.1002/ceat.201200116]
[3]: [https://www.sciencedirect.com/science/article/pii/S1365160904002874]
--
You received this question notification because your team yade-users is
an answer contact for Yade.