yade-users team mailing list archive
-
yade-users team
-
Mailing list archive
-
Message #08490
Re: [Question #238289]: Use Law2_ScGeom_CpmPhys_Cpm() with two spheres
Question #238289 on Yade changed:
https://answers.launchpad.net/yade/+question/238289
Hicham BENNIOU posted a new comment:
Thanks for the quick answer Jan !
->Could you try to "debug" the four values values from UX (just 'print axis'
is enough)? This is the first suspicious place..
negIds=[0], posIds=[1], axis=2, area=0.0
for further info, here is the complete script (might be familiar to you
as it is a modified version of uniax.py :-) )
from yade import plot,pack,timing,utils
import time, sys, os, copy
readParamsFromTable(
noTableOk=True,
young=24e9,
poisson=.2,
sigmaT=3.5e6,
frictionAngle=0,
epsCrackOnset=1e-4,
relDuctility=30,
intRadius=1,
dtSafety=1,
damping=0.4,
strainRateTension=.05,
strainRateCompression=.05,
setSpeeds=True,
doModes=3, # 1=tension, 2=compression , 3=both
sphereRadius=3.5e-3,
isoPrestress=0, # isotropic confinement
)
from yade.params.table import *
if 'description' in O.tags.keys():
O.tags['id']=O.tags['id']+O.tags['description']
concreteId=O.materials.append(CpmMat(young=young,frictionAngle=frictionAngle,poisson=poisson,density=4800,sigmaT=sigmaT,relDuctility=relDuctility,epsCrackOnset=epsCrackOnset,isoPrestress=isoPrestress))
O.bodies.append([
utils.sphere(center=(0,0,-1.5*sphereRadius),radius=sphereRadius,material=concreteId,fixed=True),
utils.sphere((0,0,1.5*sphereRadius),radius=sphereRadius,material=concreteId)
])
yade.qt.View()
UX=uniaxialTestFeatures()
negIds,posIds,axis,crossSectionArea=UX['negIds'],UX['posIds'],UX['axis'],UX['area']
O.dt=dtSafety*PWaveTimeStep()
print 'Timestep',O.dt
O.engines=[
ForceResetter(),
InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=intRadius,label='is2aabb')],verletDist=.5*sphereRadius),
InteractionLoop(
[Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=intRadius,label='ss2sc')],
[Ip2_CpmMat_CpmMat_CpmPhys()],
[Law2_ScGeom_CpmPhys_Cpm()],
),
NewtonIntegrator(damping=damping,label='damper'),
CpmStateUpdater(realPeriod=.5),
UniaxialStrainer(strainRate=strainRateTension,axis=axis,asymmetry=0,posIds=posIds,negIds=negIds,crossSectionArea=crossSectionArea,blockDisplacements=False,blockRotations=False,setSpeeds=setSpeeds,label='strainer'),
PyRunner(virtPeriod=1e-6/strainRateTension,realPeriod=1,command= 'addPlotData()',label='plotDataCollector',initRun=True),
PyRunner(realPeriod=4,command='stopIfDamaged()',label='damageChecker'),
]
plot.plots={'eps':('sigma',)}
O.saveTmp('initial');
O.timingEnabled=False
global mode
mode='tension' if doModes & 1 else 'compression'
def initTest():
global mode
print "init"
if O.iter>0:
O.wait();
O.loadTmp('initial')
print "Reversing plot data"; plot.reverseData()
else: plot.plot(subPlots=False)
strainer.strainRate=abs(strainRateTension) if mode=='tension' else -abs(strainRateCompression)
try:
from yade import qt
renderer=qt.Renderer()
renderer.dispScale=(1000,1000,1000) if mode=='tension' else (100,100,100)
except ImportError: pass
print "init done, will now run."
O.step();
ss2sc.interactionDetectionFactor=1.
is2aabb.aabbEnlargeFactor=1.
O.run()
def stopIfDamaged():
global mode
if O.iter<2 or not plot.data.has_key('sigma'): return
sigma,eps=plot.data['sigma'],plot.data['eps']
extremum=max(sigma) if (strainer.strainRate>0) else min(sigma)
minMaxRatio=0.5 if mode=='tension' else 1
if extremum==0: return
import sys; sys.stdout.flush()
if abs(sigma[-1]/extremum)<minMaxRatio or abs(strainer.strain)>(5e-3 if isoPrestress==0 else 5e-2):
if mode=='tension' and doModes & 2:
mode='compression'
print "Damaged, switching to compression... "; O.pause()
import thread; thread.start_new_thread(initTest,())
return
else:
print "Damaged, stopping."
ft,fc=max(sigma),min(sigma)
print 'Strengths fc=%g, ft=%g, |fc/ft|=%g'%(fc,ft,abs(fc/ft))
title=O.tags['description'] if 'description' in O.tags.keys() else O.tags['params']
print 'gnuplot',plot.saveGnuplot(O.tags['id'],title=title)
print 'Bye.'
O.pause()
def addPlotData():
yade.plot.addData({'t':O.time,'i':O.iter,'eps':strainer.strain,'sigma':strainer.avgStress+isoPrestress,
})
plot.plot(subPlots=False)
initTest()
waitIfBatch()
Cheers !
--
You received this question notification because you are a member of
yade-users, which is an answer contact for Yade.