[Yade-users] [Question #702407]: CGAL ERROR: assertion violation! when running HM coupling
New question #702407 on Yade:
https://answers.launchpad.net/yade/+question/702407
Hi,
I am trying to do HM coupling with cylinder sample.Under the confining pressure
of 30MPa,I have get the result of which the flow pressure is 10,20,30MPa.
However,when I add the flow pressure to 40 and 50MPa,I get the following error
prompt and my simulation stop:
###
AREA <= 0!!
AREA <= 0!!
AREA <= 0!!
AREA <= 0!!
AREA <= 0!!
AREA <= 0!!
AREA <= 0!!
ThreadRunner:35 void yade::ThreadRunner::run(): Exception
occured:
CGAL ERROR: assertion violation!
Expr: dexp != 2047
File: /usr/include/CGAL/Mpzf.h
Line: 411
Explanation: Creating an Mpzf from infinity or NaN.
###
I do not know what happened.I will attach my code in the end.Thanks for anyone
help!
###code###
from yade import pack, ymport, plot, utils, export, timing
import numpy as np
import sys
damp=0.4
young=66.2e9
name='JCFPM_triax'
poisson=0.522
name='cylinder'
preStress=-30e6
strainRate = -0.1
OUT=str(preStress)+'_JCFPM_triax'
r1=0.002
r2=0.01
nw=24
nh=15
rParticle=0.000731723
bcCoeff = 5
width = 0.025
height = 0.05
allx,ally,allz,r=np.genfromtxt('packing-cylinder.spheres', unpack=True)
mnx=min(allx)*1.01 #same as aabbExtrema, get the xyz of the lowest corner,
multiply by factor to reduce it to let the wall surround all the spheres.
mny=min(ally)*1.01
mnz=min(allz)*0.99
mxx=max(allx)*1.01 #same as aabbExtrema, get the xyz of the top corner,
multiply by factor to increase it to let the wall surround all the spheres.
mxy=max(ally)*1.01
mxz=max(allz)*1.01
Sy=50e6
O.materials.append(JCFpmMat(type=1,density=2640,young=young,poisson=poisson,tensileStrength=8e6,cohesion=40e6,frictionAngle=radians(25),label='sphere'))
O.materials.append(JCFpmMat(type=1,label='wall1',young=young,poisson=poisson,frictionAngle=radians(25)))
O.materials.append(JCFpmMat(type=1,frictionAngle=0,density=0,label='wall2'))
mn,mx=Vector3(mnx,mny,mnz),Vector3(mxx,mxy,mxz)
walls=aabbWalls([mn,mx],thickness=0,material='wall2')
wallIds=O.bodies.append(walls)
bottom,top=Vector3(0,0,0),Vector3(0,0,0.05)
sp=pack.SpherePack()
pred=pack.inCylinder(bottom,top,radius)
sp=pack.randomDensePack(pred,radius=0.0005,rRelFuzz=0,returnSpherePack=True,memoizeDb='/tmp/triax.sqlite')
spheres=sp.toSimulation(color=(0.9,0.8,0.6),material='sphere')
bot = [O.bodies[s] for s in spheres if
O.bodies[s].state.pos[2]height-rParticle*bcCoeff]
vel = strainRate*(height-rParticle*2*bcCoeff)
for s in bot:
s.shape.color = (1,0,0)
s.state.blockedDOFs = 'xyzXYZ'
#s.state.vel = (0,0,-vel)
for s in top:
s.shape.color = Vector3(0,1,0)
s.state.blockedDOFs = 'xyzXYZ'
s.state.vel = (0,0,vel)
facets = []
rCyl2 = .5*width / cos(pi/float(nw))
for r in range(nw):
for h in range(nh):
v1 = Vector3( rCyl2*cos(2*pi*(r+0)/float(nw)),
rCyl2*sin(2*pi*(r+0)/float(nw)), height*(h+0)/float(nh) )
v2 = Vector3( rCyl2*cos(2*pi*(r+1)/float(nw)),
rCyl2*sin(2*pi*(r+1)/float(nw)), height*(h+0)/float(nh) )
v3 = Vector3( rCyl2*cos(2*pi*(r+1)/float(nw)),
rCyl2*sin(2*pi*(r+1)/float(nw)), height*(h+1)/float(nh) )
v4 = Vector3( rCyl2*cos(2*pi*(r+0)/float(nw)),
rCyl2*sin(2*pi*(r+0)/float(nw)), height*(h+1)/float(nh) )
f1 = facet((v1,v2,v3),color=(0,0,1),material='wall1')
f2 = facet((v1,v3,v4),color=(0,0,1),material='wall1')
facets.extend((f1,f2))
O.bodies.append(facets)
mass = O.bodies[7].state.mass
for f in facets:
f.state.mass = mass
f.state.blockedDOFs = 'XYZz'
def lateral():
elatTot=0.0
nTot=0
for b in O.bodies:
x=b.state.refPos[0]
y=b.state.refPos[1]
d=sqrt(pow(x,2)+pow(y,2))
if d > r1 and d < r2:
b.shape.color=(0.6,0.5,0.15)
xnew=b.state.pos[0]
ynew=b.state.pos[1]
dnew=sqrt(pow(xnew,2)+pow(ynew,2))
elat=(dnew-d)/d
elatTot+=elat
nTot+=1
elat_avg=elatTot/nTot
return elat_avg
def addForces():
for f in facets:
n = f.shape.normal
a = f.shape.area
O.forces.addF(f.id,preStress*a*n)
def stopIfDamaged(maxEps=0.001):
extremum = max(abs(s) for s in plot.data['s'])
s = abs(plot.data['s'][-1])
e = abs(plot.data['e'][-1])
if O.iter < 1000 or e < maxEps:
return
if abs(s)/abs(extremum) < 0.5 :
print('Simulation finished')
presentcohesive_count = 0
for i in O.interactions:
if hasattr(i.phys, 'isCohesive'):
if i.phys.isCohesive == True:
presentcohesive_count+=1
print('the number of cohesive bond now
is:',presentcohesive_count)
print('Max stress and strain:',extremum,e)
O.pause()
O.dt=.5*utils.PWaveTimeStep()
newton=NewtonIntegrator(d
Re: [Yade-users] [Question #702407]: CGAL ERROR: assertion violation! when running HM coupling
Question #702407 on Yade changed: https://answers.launchpad.net/yade/+question/702407 Status: Open => Answered Bruno Chareyre proposed the following answer: Hi, This error suggests huge overlaps between some spheres, such that in some pore throats the intersected solid area is more than the area of the triangle connecting sphere centers (hence a negative area left for the fluid). It may come from a small stress to stiffness ratio, or more simply a numerical instability. Try with a smaller time-step? Bruno -- You received this question notification because your team yade-users is an answer contact for Yade. ___ Mailing list: https://launchpad.net/~yade-users Post to : yade-users@lists.launchpad.net Unsubscribe : https://launchpad.net/~yade-users More help : https://help.launchpad.net/ListHelp
Re: [Yade-users] [Question #702344]: Permeability of non-spherical (e.g. polyhedral) particles packing
Question #702344 on Yade changed: https://answers.launchpad.net/yade/+question/702344 Status: Answered => Solved Soheil Safari confirmed that the question is solved: Dear Robert, Thank you very much for your kind reply and useful link . Best regards, Soheil -- You received this question notification because your team yade-users is an answer contact for Yade. ___ Mailing list: https://launchpad.net/~yade-users Post to : yade-users@lists.launchpad.net Unsubscribe : https://launchpad.net/~yade-users More help : https://help.launchpad.net/ListHelp
Re: [Yade-users] [Question #702405]: This is a test question for checking email notifications from Launchpad
Question #702405 on Yade changed: https://answers.launchpad.net/yade/+question/702405 Status: Open => Answered Janek Kozicki proposed the following answer: Looks like it works! I received the email notification :) -- You received this question notification because your team yade-users is an answer contact for Yade. ___ Mailing list: https://launchpad.net/~yade-users Post to : yade-users@lists.launchpad.net Unsubscribe : https://launchpad.net/~yade-users More help : https://help.launchpad.net/ListHelp
[Yade-users] [Question #702411]: Calculation of a time averaging of a physical quantity
New question #702411 on Yade: https://answers.launchpad.net/yade/+question/702411 Hello everyone, I would like to calculate a time averaging of a physical quantity in yade (the average of the velocity of a group of particles made over a period of rotation of the drum). Below is a part of my code that allows me to calculate the average velocity of several groups of particles, but this average is at a given time. What I would like is the average made over one revolution. Any suggestion will be appreciated. Part of my code: for i in range(0,19): teta0=i*dteta global listcase listcase=[] #global listcase for b in O.bodies: if isinstance(b.shape,Sphere): rayon=sqrt((b.state.pos[0])**2 +(b.state.pos[1])**2) dr=d teta=atan(b.state.pos[1]/b.state.pos[0]) if (b.state.pos[1]>0 and b.state.pos[0]<0) or (b.state.pos[1]<0 and b.state.pos[0]<0): teta=pi-teta if (b.state.pos[1]<0 and b.state.pos[0]>0): teta=2*pi+teta if (rayonrcylint) and (teta>teta0 and teta
Re: [Yade-users] [Question #702407]: CGAL ERROR: assertion violation! when running HM coupling
Question #702407 on Yade changed: https://answers.launchpad.net/yade/+question/702407 Ziyu Wang posted a new comment: Hi,Bruno. I change the time-step from .5*PWaveTimeStep() to .1*PWaveTimeStep(),and it seems fine right now(I haven't got the final result yet).I will give feedback when I get the final result. A small question about this:Whether changing the time step will affect the final result (I am also changing the time step of the previous simulation to verify the result, and will give feedback in time). Thanks again! -- You received this question notification because your team yade-users is an answer contact for Yade. ___ Mailing list: https://launchpad.net/~yade-users Post to : yade-users@lists.launchpad.net Unsubscribe : https://launchpad.net/~yade-users More help : https://help.launchpad.net/ListHelp
[Yade-users] [Question #702412]: O.interactions[0, 1] is unable to detect the interactions between sphere and wall, throwing an error stating "No such interaction"
New question #702412 on Yade:
https://answers.launchpad.net/yade/+question/702412
Hi everyone,
I was doing a simulation of a sphere hitting on a wall at some distance in an
acceleration field. My code is given bellow:
## code starts here
## SCRIPT TO TEST IMPACT BETWEEN A SPHERE AND A WALL (MINDLIN - nonlinear
elastic model)
## list of engines
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_FrictMat_FrictMat_MindlinPhys()],
[Law2_ScGeom_MindlinPhys_Mindlin()]
),
NewtonIntegrator(damping=0.0,gravity=(10,0,0)),
###
### NOTE this extra engine:
###
PyRunner(iterPeriod=1,command='myAddPlotData()')
]
## define and append material
mat=FrictMat(young=7.0e10,poisson=0.3,density=2699,frictionAngle=radians(0),label='Al.alloy')
O.materials.append(mat)
## create a sphere and a wall and append them
s0=sphere([0,0,0],0.1,color=[0,1,0],fixed=False,wire=False,material='Al.alloy')
w=wall([0.2,0,0],axis=0,sense=0,color=[1,0,0])
O.bodies.append(s0)
O.bodies.append(w)
## time step
O.dt=.2*PWaveTimeStep()
O.saveTmp('Mindlin')
from yade import qt
qt.View()
qt.Controller()
# now the part pertaining to plots #
from yade import plot
def myAddPlotData():
if O.interactions[0,1].isReal:
i=O.interactions[0,1]
## O.interactions[0,1] is not able to detect the interactions
between sphere and wall
## error is throwing here stating "No such interaction"
## but while running the test, it is clearly visible that the
sphere is colliding with wall and bouncing back
## rest of the code is working good. PLEASE HELP!!!
plot.addData(fn=i.phys.normalForce[0],step=O.iter,un=s0.shape.radius-w.state.pos[0]+s0.state.pos[0],kn=i.phys.kn)
O.run(100,True)
## make one plot: step as function of fn
plot.plots={'un':('fn')}
## this function is called by plotDataCollector
## it should add data with the labels that we will plot
## if a datum is not specified (but exists), it will be NaN and will not be
plotted
plot.plot()
## code ends here
In the above code, an error is throwing in the 46th line. The error message is
given below:
IndexErrorTraceback (most recent call last)
/usr/bin/yade in
/usr/bin/yade in myAddPlotData()
44
45 def myAddPlotData():
---> 46 if O.interactions[0,1].isReal:
47 i=O.interactions[0,1]
48 ## O.interactions[0,1] is not able to detect the
interactions between sphere and wall
IndexError: No such interaction
The exact same code runs good when I am trying with two spheres instead of
sphere and a wall.
Please help me to figure this out. I have to track contact forces while they
are impacting.
Please also tell me if there is any other way to track contact force rather
using O.interactions[0,1].
Thanks in advance.
--
You received this question notification because your team yade-users is
an answer contact for Yade.
___
Mailing list: https://launchpad.net/~yade-users
Post to : yade-users@lists.launchpad.net
Unsubscribe : https://launchpad.net/~yade-users
More help : https://help.launchpad.net/ListHelp
Re: [Yade-users] [Question #702411]: Calculation of a time averaging of a physical quantity
Question #702411 on Yade changed: https://answers.launchpad.net/yade/+question/702411 Status: Open => Answered Bruno Chareyre proposed the following answer: Hi, Something like this? timeSeries = [] while totalRotationhttps://launchpad.net/~yade-users Post to : yade-users@lists.launchpad.net Unsubscribe : https://launchpad.net/~yade-users More help : https://help.launchpad.net/ListHelp
Re: [Yade-users] [Question #702412]: O.interactions[0, 1] is unable to detect the interactions between sphere and wall, throwing an error stating "No such interaction"
Question #702412 on Yade changed: https://answers.launchpad.net/yade/+question/702412 Status: Open => Answered Bruno Chareyre proposed the following answer: Hi, If there are only two bodies O.forces.f(0) and O.interactions[0,1] should be the same. Except that the first will not give an error when there is no interaction. Alternatively, this boolean could help: O.interactions.has(0,1) Cheers Bruno -- You received this question notification because your team yade-users is an answer contact for Yade. ___ Mailing list: https://launchpad.net/~yade-users Post to : yade-users@lists.launchpad.net Unsubscribe : https://launchpad.net/~yade-users More help : https://help.launchpad.net/ListHelp
Re: [Yade-users] [Question #702407]: CGAL ERROR: assertion violation! when running HM coupling
Question #702407 on Yade changed: https://answers.launchpad.net/yade/+question/702407 Status: Answered => Open Ziyu Wang is still having a problem: Hi, For the case of flow pressure is 50MPa,there was no problem most of the time, but in the end, another error was reported: # ThreadRunner:35 void yade::ThreadRunner::run(): Exception occured: /builds/yade-dev/trunk/deb/yadedaily/pkg/pfv/Network.ipp : switch default case error. # For my small question,I studied the case where the flow pressure was 10MPa, and found that the curve became smoother, which had no great impact on the overall result. -- You received this question notification because your team yade-users is an answer contact for Yade. ___ Mailing list: https://launchpad.net/~yade-users Post to : yade-users@lists.launchpad.net Unsubscribe : https://launchpad.net/~yade-users More help : https://help.launchpad.net/ListHelp
Re: [Yade-users] [Question #702316]: how to get the change of permeability during triaxial compression in flowengine
Question #702316 on Yade changed: https://answers.launchpad.net/yade/+question/702316 Ziyu Wang confirmed that the question is solved: Thanks Robert Caulk, that solved my question. -- You received this question notification because your team yade-users is an answer contact for Yade. ___ Mailing list: https://launchpad.net/~yade-users Post to : yade-users@lists.launchpad.net Unsubscribe : https://launchpad.net/~yade-users More help : https://help.launchpad.net/ListHelp
Re: [Yade-users] [Question #702316]: how to get the change of permeability during triaxial compression in flowengine
Question #702316 on Yade changed: https://answers.launchpad.net/yade/+question/702316 Status: Open => Solved Ziyu Wang confirmed that the question is solved: Hi, I have known how to measure the permeabilty about my question,but still not understand the bndcondvalue you provided..I will open an new thread asking about the bndcondValue. Thanks again Robert! -- You received this question notification because your team yade-users is an answer contact for Yade. ___ Mailing list: https://launchpad.net/~yade-users Post to : yade-users@lists.launchpad.net Unsubscribe : https://launchpad.net/~yade-users More help : https://help.launchpad.net/ListHelp
[Yade-users] [Question #702418]: About flow.bndCondIsPressure and flow.bndCondValue in flowengine
New question #702418 on Yade: https://answers.launchpad.net/yade/+question/702418 Hi, I am doing fluid-solid coupling simulation of a cylinder.In the lab,I first apply a fixed confining pressure, and then carry out axial loading while maintaining the osmotic pressure in the longitudinal direction (z direction) until failure.So I want to do the same thing in simulation. I read the explanation in [1] about flow.bndCondIsPressure and flow.bndCondValue,but I have some doubts about the setting of these two values. 1)flow.bndCondIsPressure:defines the type of boundary condition for each side. True if pressure is imposed, False for no-flux. As I described above,there is osmotic pressure only in the Z direction. More specifically, it is only applied in the positive direction of the Z axis. Measure the permeability from the other end. Can it be set to =[0,0,0,0,0,1] or [1,1,1,1,1,1] 2)flow.bndCondValue:Imposed value of a boundary condition. Only applies if the boundary condition is imposed pressure, else the imposed flux is always zero presently.Ditto, I set it as =[0,0,0,0,0,P], is it reasonable? Cheers everyone. [1]https://yade-dem.org/doc/yade.wrapper.html?highlight=flowengine#yade.wrapper.FlowEngine -- You received this question notification because your team yade-users is an answer contact for Yade. ___ Mailing list: https://launchpad.net/~yade-users Post to : yade-users@lists.launchpad.net Unsubscribe : https://launchpad.net/~yade-users More help : https://help.launchpad.net/ListHelp
