Question #665857 on Yade changed:
https://answers.launchpad.net/yade/+question/665857
Status: Answered => Open
Seti is still having a problem:
Hi Robert,
thanks for your reply, no error just below message:
init done, will now run.
Damaged, stopping.
Compressive strength fc=691231
gnuplot 20180330T170719p8266.gnuplot
Bye.
here is the code:
#!/usr/bin/python
# -*- coding: utf-8 -*-
from __future__ import division
from yade import plot,pack,timing
import time, sys, os, copy
#import matplotlib
#matplotlib.rc('text',usetex=True)
#matplotlib.rc('text.latex',preamble=r'\usepackage{concrete}\usepackage{euler}')
"""
A fairly complex script performing uniaxial tension-compression test on
hyperboloid-shaped specimen.
Most parameters of the model (and of the setup) can be read from table using
yade-multi.
After the simulation setup, tension loading is run and stresses are
periodically saved for plotting
as well as checked for getting below the maximum value so far. This indicates
failure (see stopIfDamaged
function). After failure in tension, the original setup is loaded anew and the
sense of loading reversed.
After failure in compression, strain-stress curves are saved via
plot.saveGnuplot and we exit,
giving some useful information like peak stresses in tension/compression.
Running this script for the first time can take long time, as the specimen is
prepared using triaxial
compression. Next time, however, an attempt is made to load
previously-generated packing
(from /tmp/triaxPackCache.sqlite) and this expensive procedure is avoided.
The specimen length can be specified, its diameter is half of the length and
skirt of the hyperboloid is
4/5 of the width.
The particle size is constant and can be specified using the sphereRadius
parameter.
The 3d display has displacement scaling applied, so that the fracture looks
more spectacular. The scale
is 1000 for tension and 100 for compression.
"""
# default parameters or from table
readParamsFromTable(noTableOk=True, # unknownOk=True,
young=24e9,
poisson=.2,
sigmaT=3.5e6,
frictionAngle=atan(0.8),
epsCrackOnset=1e-4,
relDuctility=30,
intRadius=1.5,
dtSafety=.8,
damping=0.4,
strainRateTension=.5,
strainRateCompression=.5,
setSpeeds=True,
# 1=tension, 2=compression (ANDed; 3=both)
doModes=2,
specimenLength=.15,
sphereRadius=3.5e-3,
# isotropic confinement (should be negative)
isoPrestress=0,
)
from yade.params.table import *
from yade.params.table import *
if 'sigmaT=3.5e6, compression' in O.tags.keys():
O.tags['id']=O.tags['id']+O.tags['sigmaT=3.5e6, compression']
# make geom; the dimensions are hard-coded here; could be in param table if
desired
# z-oriented hyperboloid, length 20cm, diameter 10cm, skirt 8cm
# using spheres 7mm of diameter
#
mat1=CpmMat(young=young,frictionAngle=frictionAngle,poisson=poisson,density=4800,sigmaT=sigmaT,relDuctility=relDuctility,epsCrackOnset=epsCrackOnset,isoPrestress=isoPrestress)
concreteId1=O.materials.append(mat1)
#sp=pack.randomDensePack(pack.inHyperboloid((0,0,-.5*specimenLength),(0,0,.5*specimenLength),.25*specimenLength,.17*specimenLength),spheresInCell=2000,radius=sphereRadius,memoizeDb='/tmp/triaxPackCache.sqlite',returnSpherePack=True)
#
sp=pack.SpherePack()
pred=pack.inCylinder((0,0,0.0002),(0,0,0.3),0.05)
O.bodies.append(pack.randomDensePack(pred,radius=0.007))
bb=uniaxialTestFeatures()
negIds,posIds,axis,crossSectionArea=bb['negIds'],bb['posIds'],bb['axis'],bb['area']
O.dt=dtSafety*PWaveTimeStep()
print 'Timestep',O.dt
mm,mx=[pt[axis] for pt in aabbExtrema()]
coord_25,coord_50,coord_75=mm+.25*(mx-mm),mm+.5*(mx-mm),mm+.75*(mx-mm)
area_25,area_50,area_75=approxSectionArea(coord_25,axis),approxSectionArea(coord_50,axis),approxSectionArea(coord_75,axis)
O.engines=[
ForceResetter(),
InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=intRadius,label='is2aabb'),],verletDist=.05*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'),
]
