yade-users team mailing list archive

[Xin Li <question701938@xxxxxxxxxxxxxxxxxxxxx>]

New question #701938 on Yade:
https://answers.launchpad.net/yade/+question/701938

Hi everyone,

I am facing some issues when using cylinderConnection with large size (radius > width of periodic cell)  under periodic conditions.
1.Axial force is very high or 'Nan' when the length of segment >= 1/2 cell size (sometimes, other small length has the same issue). Is there a way to delete the contact geom and phys between nodes and gridconnections as I don't need the inner interactions of cylinder in my case. 
2.Losing interactions with particles or some particles cross the cylinder and then interactions are detected after a period of time, this will cause explosion. Some specific number and length of segments of the cylinder can work normally (for example, in my MWE, length of each segment of the cylinder  is cell width/2.1 and repeat 20 times). I have not tried the same size as cell, because the error 'float division by zero' from gridPfacet.py line140-i.phys.kn=E*math.pi*(radius**2)/L (two nodes located in the boundaries of different cells and their distance is zero?).
My MWE as below:
# encoding: utf-8

from builtins import range
import math
import glob, os
import numpy as np
from yade import pack, plot
from yade import geom
from yade import export
from yade import polyhedra_utils
from yade.gridpfacet import *

# Soil sphere parameters
E_m=8.1654644182
E=pow(10,E_m) # micro Young's modulus
v=1.1131826506e-01 # micro Poisson's ratio
mu=2.6666064103e+01 # interparticle friction angle
beta=7.2094747073e-01 # rolling stiffness
eta=6.9562665265e-02 # plastic limit of rolling
rho = 2650 # soil density
         
spMat = O.materials.append(CohFrictMat(young=E,poisson=v,frictionAngle=radians(mu),density=rho,isCohesive=False,alphaKr=beta,alphaKtw=beta,momentRotationLaw=True,etaRoll=eta,etaTwist=eta))
wallMat = O.materials.append(CohFrictMat(young=E,poisson=v,frictionAngle=radians(0),density=rho,isCohesive=False,alphaKr=0,alphaKtw=0,momentRotationLaw=False,etaRoll=0,etaTwist=0)) 
extcylMat = O.materials.append(CohFrictMat(young=E,poisson=v,frictionAngle=radians(0),density=1190,isCohesive=False,alphaKr=0,alphaKtw=0,momentRotationLaw=True,etaRoll=0,etaTwist=0)) 
intcylMat = O.materials.append(CohFrictMat(young=0,poisson=0,frictionAngle=radians(0),density=1190,isCohesive=False,alphaKr=0,alphaKtw=0,momentRotationLaw=True,etaRoll=0,etaTwist=0)) 

O.periodic = True
x=0.024 # 0.02
y=0.1
z=0.1
O.cell.setBox(x,y*3,z*3)

D=0.075 # 37.5
r=0.00225 # 0.0015
d=r*2
sp=pack.SpherePack()
mn,mx=Vector3(0,y,z+D),Vector3(x,y*2,z*2)                                                                                     
sp.makeCloud(mn,mx,rMean=r,rRelFuzz=0,seed=1)
sphereIds1=sp.toSimulation(material=spMat)
print(len(sp))

for i in range(0,3):
    O.bodies.append(utils.box((x/6+x*i/3,y,z*1.5),(x/6,0,z/2),wire=True,fixed=True,color=[1,0,0],material=wallMat))
    O.bodies.append(utils.box((x/6+x*i/3,y*2,z*1.5),(x/6,0,z/2),wire=True,fixed=True,color=[1,0,0],material=wallMat))
    O.bodies.append(utils.box((x/6+x*i/3,y*1.5,z),(x/6,y/2,0),wire=True,fixed=True,color=[1,0,0],material=wallMat))

#define engines
O.engines=[
   ForceResetter(),
   InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb(),Bo1_GridConnection_Aabb()],allowBiggerThanPeriod=True),
   InteractionLoop(
        [Ig2_Sphere_Sphere_ScGeom6D(),Ig2_Box_Sphere_ScGeom6D(),
         Ig2_Sphere_GridConnection_ScGridCoGeom(),
        ],
		[Ip2_CohFrictMat_CohFrictMat_CohFrictPhys(setCohesionNow=True, setCohesionOnNewContacts=False)],
		[Law2_ScGridCoGeom_CohFrictPhys_CundallStrack(),
         Law2_ScGeom6D_CohFrictPhys_CohesionMoment(
		   always_use_moment_law=True,
		   useIncrementalForm=True
		)],
   ),
   GlobalStiffnessTimeStepper(timestepSafetyCoefficient=0.8),
   NewtonIntegrator(gravity=(0,0,-9.81),damping=0.2,label='newton'),
   PyRunner(initRun=False,iterPeriod=20000,command='check()'),
]

nodesIds=[]
cylIds=[]
vertices=[]
L=x/2.1 
for i in range(0,21):
    vertices.append([i*L,y*1.5,z+0.5*D])
cylinderConnection(vertices,radius=0.5*D,fixed = True,dynamic=False,nodesIds=nodesIds, cylIds=cylIds, intMaterial=intcylMat, extMaterial=extcylMat)
print(nodesIds,cylIds)

def check():
    for i in O.interactions:
        if i.id1 in nodesIds and i.id2 in nodesIds:
            print(i.id1,i.id2,i.phys.normalForce.norm(),i.phys.shearForce.norm())
        if i.id1 in cylIds and i.id2 in cylIds:
            print(i.id1,i.id2,i.phys.normalForce.norm(),i.phys.shearForce.norm())
    print('unb=',unbalancedForce())
    if unbalancedForce()<0.01:
        print('stable')
        O.pause()

Take my MWE for example,settings of cylinderConnection remain the same and I tried several tests as below:
1.cell width=0.024m, cylinder radius=0.075m, minimum particle radius=0.00225mm-work normally
2.cell width=0.024m, cylinder radius=0.075m, minimum particle radius=0.00150mm-explosion
3.cell width=0.024m, cylinder radius=0.0375m, minimum particle radius=0.00225mm-explosion
4.cell width=0.024m, cylinder radius=0.0375m, minimum particle radius=0.00150mm-explosion
5.cell width=0.02m, cylinder radius=0.075m, minimum particle radius=0.00225mm-work normally
6.cell width=0.02m, cylinder radius=0.075m, minimum particle radius=0.00150mm-explosion
7.cell width=0.02m, cylinder radius=0.0375m, minimum particle radius=0.00225mm-work normally
8.cell width=0.02m, cylinder radius=0.0375m, minimum particle radius=0.00150mm-work normally
I am confused about how to set the cylinder especially when I need larger cylinder and smaller particle.
My yade version is Yade 20220218-6333~d331682~focal1 and system is Ubuntu 20.04.3 LTS
Dose anyone has the same problem when using large size cylinderConnection under periodic condition? 
Thanks for any help.
Relative questions: https://answers.launchpad.net/yade/+question/181411; https://answers.launchpad.net/yade/+question/689446

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