yade-dev team mailing list archive
-
yade-dev team
-
Mailing list archive
-
Message #07584
[Branch ~yade-dev/yade/trunk] Rev 2850: - some changes in normalInelasticity-test
------------------------------------------------------------
revno: 2850
committer: jduriez <jduriez@c1solimara-l>
branch nick: yade
timestamp: Tue 2011-05-03 17:04:37 +0200
message:
- some changes in normalInelasticity-test
- add of a check version of this script
- add of a sentence about these check tests in prog.rst
- changes of deprecated function in tutorial scripts
- typo in some other rst
- add of a note about plot.plot in plot.py
- beginning changes in SimpleShear preprocessor
added:
scripts/test/checks/checkTestNormalInelasticity.py
modified:
doc/sphinx/prog.rst
doc/sphinx/tutorial-hands-on.rst
doc/sphinx/tutorial/02-gravity-deposition.py
doc/sphinx/tutorial/03-oedometric-test.py
doc/sphinx/user.rst
pkg/dem/SimpleShear.cpp
pkg/dem/SimpleShear.hpp
py/plot.py
scripts/normalInelasticity-test.py
scripts/simpleShear.py
--
lp:yade
https://code.launchpad.net/~yade-dev/yade/trunk
Your team Yade developers is subscribed to branch lp:yade.
To unsubscribe from this branch go to https://code.launchpad.net/~yade-dev/yade/trunk/+edit-subscription
=== modified file 'doc/sphinx/prog.rst'
--- doc/sphinx/prog.rst 2011-04-28 14:20:15 +0000
+++ doc/sphinx/prog.rst 2011-05-03 15:04:37 +0000
@@ -243,7 +243,7 @@
Since the check tests history will be mostly based on standard output generated by "yade ---check", a meaningfull checkTest should include some "print" command telling if something went wrong. If the script itself fails for some reason and can't generate an output, the log will contain "scriptName failure". If the script defines differences on obtained and awaited data, it should print some useful information about the problem and increase the value of global variable resultStatus. After this occurs, the automatic test will stop the execution with error message.
An example check test can be found in checkTestTriax.py. It shows results comparison, output, and how to define the path to data files using "checksPath".
-Users are encouraged to add their own scripts into the scripts/test/checks/ folder. Discussion of some specific checktests design in users question is welcome.
+Users are encouraged to add their own scripts into the scripts/test/checks/ folder. Discussion of some specific checktests design in users question is welcome. Note that re-compiling is required before that added scripts can be launched by "yade ---check" (or direct changes have to be performed in "lib" subfolders).
A check test should never need more than a few seconds to run. If your typical script needs more, try and reduce the number of element or the number of steps.
Conventions
=== modified file 'doc/sphinx/tutorial-hands-on.rst'
--- doc/sphinx/tutorial-hands-on.rst 2011-01-20 16:35:46 +0000
+++ doc/sphinx/tutorial-hands-on.rst 2011-05-03 15:04:37 +0000
@@ -480,7 +480,7 @@
#. define particles as in the previous exercise (cloud of spheres inside a box open from the top)
#. define a simple simulation loop, as the one given above
- #. set $\Dt$ equal to the critial P-Wave $\Dt$
+ #. set $\Dt$ equal to the critical P-Wave $\Dt$
#. save the initial simulation state to memory
#. Run the previously-defined simulation multiple times, while changing the value of timestep (use the :guilabel:`⟳` button to reload the initial configuration).
=== modified file 'doc/sphinx/tutorial/02-gravity-deposition.py'
--- doc/sphinx/tutorial/02-gravity-deposition.py 2011-02-15 11:36:29 +0000
+++ doc/sphinx/tutorial/02-gravity-deposition.py 2011-05-03 15:04:37 +0000
@@ -6,7 +6,7 @@
from yade import pack, plot
# create rectangular box from facets
-O.bodies.append(utils.facetBox((.5,.5,.5),(.5,.5,.5),wallMask=31))
+O.bodies.append(utils.geom.facetBox((.5,.5,.5),(.5,.5,.5),wallMask=31))
# create empty sphere packing
# sphere packing is not equivalent to particles in simulation, it contains only the pure geometry
=== modified file 'doc/sphinx/tutorial/03-oedometric-test.py'
--- doc/sphinx/tutorial/03-oedometric-test.py 2011-01-20 16:35:46 +0000
+++ doc/sphinx/tutorial/03-oedometric-test.py 2011-05-03 15:04:37 +0000
@@ -17,7 +17,7 @@
# create box with free top, and ceate loose packing inside the box
from yade import pack, plot
-O.bodies.append(utils.facetBox((.5,.5,.5),(.5,.5,.5),wallMask=31))
+O.bodies.append(utils.geom.facetBox((.5,.5,.5),(.5,.5,.5),wallMask=31))
sp=pack.SpherePack()
sp.makeCloud((0,0,0),(1,1,1),rMean=rMean,rRelFuzz=rRelFuzz)
sp.toSimulation()
=== modified file 'doc/sphinx/user.rst'
--- doc/sphinx/user.rst 2011-04-29 13:28:28 +0000
+++ doc/sphinx/user.rst 2011-05-03 15:04:37 +0000
@@ -1412,7 +1412,7 @@
without further explanation. Frequent causes of such conditions are
* program error in Yade itself;
-* fatal condition in you particular simulation (such as impossible dispatch);
+* fatal condition in your particular simulation (such as impossible dispatch);
* problem with graphics card driver.
Try to reproduce the error (run the same script) with debug-enabled version of Yade. Debugger will be automatically launched at crash, showing backtrace of the code (in this case, we triggered crash by hand)::
=== modified file 'pkg/dem/SimpleShear.cpp'
--- pkg/dem/SimpleShear.cpp 2010-12-05 17:10:06 +0000
+++ pkg/dem/SimpleShear.cpp 2011-05-03 15:04:37 +0000
@@ -21,12 +21,13 @@
#include<yade/pkg/common/InsertionSortCollider.hpp>
#include<yade/core/Interaction.hpp>
#include<yade/pkg/common/Dispatching.hpp>
+#include<yade/pkg/common/InteractionLoop.hpp>
#include<yade/pkg/common/ForceResetter.hpp>
#include<yade/pkg/dem/NewtonIntegrator.hpp>
#include<yade/pkg/common/GravityEngines.hpp>
-#include<yade/pkg/dem/KinemCNDEngine.hpp>
+#include<yade/pkg/dem/KinemCTDEngine.hpp>
#include<yade/pkg/dem/Ig2_Sphere_Sphere_ScGeom.hpp>
#include<yade/pkg/dem/Ig2_Box_Sphere_ScGeom.hpp>
@@ -43,6 +44,7 @@
using namespace std;
YADE_PLUGIN((SimpleShear))
+CREATE_LOGGER(SimpleShear);
SimpleShear::~SimpleShear ()
{
@@ -76,6 +78,7 @@
shared_ptr<Body> w4; // The upper one
createBox(w4,Vector3r(length/2.0,height+thickness/2.0,0),Vector3r(length/2.0,thickness/2.0,width/2.0));
YADE_PTR_CAST<FrictMat> (w4->material)->frictionAngle = sphereFrictionDeg * Mathr::PI/180.0; // so that we have phi(spheres-superior wall)=phi(sphere-sphere)
+ scene->bodies->insert(w4);
// To close the front and the back of the box
shared_ptr<Body> w5; // behind
@@ -92,10 +95,11 @@
vector<BasicSphere> sphere_list;
// to use the TriaxialTest method :
-// GenerateCloud(sphere_list,Vector3r(0,0,-width/2.0),Vector3r(length,height,width/2.0),nBilles,0.3,porosite);
+ string out=GenerateCloud(sphere_list,Vector3r(0,0,-width/2.0),Vector3r(length,height,width/2.0),1000,0.3,0.7);// generates a sample of 1000 spheres, with a required porosity of 0.7
+ LOG_INFO(out);
// to use a text file :
- std::pair<string,bool> res=ImportCloud(sphere_list,filename);
+// std::pair<string,bool> res=ImportCloud(sphere_list,filename);
vector<BasicSphere>::iterator it = sphere_list.begin();
vector<BasicSphere>::iterator it_end = sphere_list.end();
@@ -107,8 +111,7 @@
scene->bodies->insert(body);
}
- message =res.first;
- return res.second;
+ return true;
}
void SimpleShear::createSphere(shared_ptr<Body>& body, Vector3r position, Real radius)
@@ -158,6 +161,7 @@
shared_ptr<Aabb> aabb(new Aabb);
// shared_ptr<BoxModel> gBox(new BoxModel);
shared_ptr<Box> iBox(new Box);
+ iBox->wire = true;
body->setDynamic(false);
@@ -210,19 +214,28 @@
globalStiffnessTimeStepper->timeStepUpdateInterval = timeStepUpdateInterval;
globalStiffnessTimeStepper->defaultDt=1e-5;
-
+ shared_ptr<KinemCTDEngine> kinemEngine (new KinemCTDEngine);
+ kinemEngine->compSpeed = 10.0;
+ kinemEngine->targetSigma=1000.0;
+
+
+ shared_ptr<InteractionLoop> ids(new InteractionLoop);
+ ids->geomDispatcher=interactionGeometryDispatcher;
+ ids->physDispatcher=interactionPhysicsDispatcher;
+ ids->lawDispatcher=shared_ptr<LawDispatcher>(new LawDispatcher);
+ shared_ptr<Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity> ldc(new Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity);
+ ids->lawDispatcher->add(ldc);
scene->engines.clear();
scene->engines.push_back(shared_ptr<Engine>(new ForceResetter));
scene->engines.push_back(globalStiffnessTimeStepper);
scene->engines.push_back(collider);
- scene->engines.push_back(interactionGeometryDispatcher);
- scene->engines.push_back(interactionPhysicsDispatcher);
-// scene->engines.push_back(shared_ptr<Engine>(new Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity));
+ scene->engines.push_back(ids);
if(gravApplied)
scene->engines.push_back(gravityCondition);
scene->engines.push_back(shared_ptr<Engine> (new NewtonIntegrator));
+ scene->engines.push_back(kinemEngine);
}
=== modified file 'pkg/dem/SimpleShear.hpp'
--- pkg/dem/SimpleShear.hpp 2011-01-21 08:14:28 +0000
+++ pkg/dem/SimpleShear.hpp 2011-05-03 15:04:37 +0000
@@ -33,7 +33,7 @@
bool generate(std::string& message);
- YADE_CLASS_BASE_DOC_ATTRS(SimpleShear,FileGenerator,"Preprocessor for creating a numerical model of a simple shear box.\n\n\t- Boxes (6) constitute the different sides of the box itself\n \t- Spheres are contained in the box. The sample could be generated via the same method used in TriaxialTest Preprocesor (=> see GenerateCloud) or by reading a text file containing positions and radii of a sample (=> see ImportCloud). This last one is the one by default used by this PreProcessor as it is written here => you need to have such a file.\n \n \t Thanks to the Engines (in pkg/common/Engine/PartialEngine) KinemCNDEngine, KinemCNSEngine and KinemCNLEngine, respectively constant normal displacement, constant normal rigidity and constant normal stress are possible to execute over such samples.\n \n \tNB about micro-parameters : their values correspond to those used in [Duriez2009a]_.",
+ YADE_CLASS_BASE_DOC_ATTRS(SimpleShear,FileGenerator,"Preprocessor for creating a numerical model of a simple shear box.\n\n\t- Boxes (6) constitute the different sides of the box itself\n \t- Spheres are contained in the box. The sample could be generated via the same method used in TriaxialTest Preprocesor (=> see GenerateCloud) or by reading a text file containing positions and radii of a sample (=> see ImportCloud). This last one is the one by default used by this PreProcessor as it is written here => you need to have such a file.\n \n \t Thanks to the Engines (in pkg/del) :yref:`KinemCNDEngine`, :yref:`KinemCNSEngine` and :yref:`KinemCNLEngine`, respectively constant normal displacement, constant normal rigidity and constant normal stress are possible to execute over such samples.\n \n \tNB about micro-parameters : their values correspond to those used in [Duriez2009a]_ and [Duriez2011]_.",
((string,filename,"../porosite0_44.txt",,"file with the list of spheres centers and radii"))
((Vector3r,gravity,Vector3r(0,-9.81,0),,"vector corresponding to used gravity [$m/s^2$]"))
((Real,thickness,0.001,,"thickness of the boxes constituting the shear box [$m$]"))
@@ -49,6 +49,7 @@
((bool,gravApplied,false,,"depending on this, :yref:`GravityEngine` is added or not to the scene to take into account the weight of particles"))
((int,timeStepUpdateInterval,50,,"value of :yref:`TimeStepper::timeStepUpdateInterval` for the :yref:`TimeStepper` used here"))
);
+ DECLARE_LOGGER;
};
REGISTER_SERIALIZABLE(SimpleShear);
=== modified file 'py/plot.py'
--- py/plot.py 2011-02-15 11:36:29 +0000
+++ py/plot.py 2011-05-03 15:04:37 +0000
@@ -566,7 +566,7 @@
def plot(noShow=False,subPlots=True):
- """Do the actual plot, which is either shown on screen (and nothing is returned: if *noShow* is ``False``) or, if *noShow* is ``True``, returned as matplotlib's Figure object or list of them.
+ """Do the actual plot, which is either shown on screen (and nothing is returned: if *noShow* is ``False`` - note that your yade compilation should present qt4 feature so that figures can be displayed) or, if *noShow* is ``True``, returned as matplotlib's Figure object or list of them.
You can use
=== modified file 'scripts/normalInelasticity-test.py'
--- scripts/normalInelasticity-test.py 2011-01-24 14:45:35 +0000
+++ scripts/normalInelasticity-test.py 2011-05-03 15:04:37 +0000
@@ -1,19 +1,20 @@
# -*- coding: utf-8 -*-
-# Script to test the constitutive law contained in NormalInelasticityLaw : consider two spheres moving one to each other (this script illustrates different ways of moving spheres)
+# Script to test the constitutive law contained in NormalInelasticityLaw : consider two spheres moving one to each other (this script illustrates few different ways of moving spheres)
#- first penetration of the contact evolves => Monitor of the normal force
#- then, test in tangential direction
#- finally for what concerns moment transfer, with relative orientation changes (use frictionAngle=0.0 in this case to study more easily what's going on)
#Different graphs illustrate the effects of the different loadings. The run is paused at each plot window (so that there is time to observe them). Push on "Return", while being in the Yade terminal, to go ahead.
-#No crash warranty with r2676
+#No crash warranty with r2676 and r2849
from yade import plot
#Def of the material which will be used
O.materials.append(NormalInelasticMat(density=2600,young=4.0e9,poisson=.04,frictionAngle=.6,coeff_dech=3.0,label='Materiau1'))
+#O.materials.append(FrictMat(young=10e9,poisson=.25,frictionAngle=0.5,density=1e3))
#Def of the bodies of the simulations : 2 spheres, with names which will be useful after
O.bodies.append(utils.sphere([0,0,0], 1, fixed=True, wire=False, color=None, highlight=False)) #'Materiau1', as the latest material defined, will be used
@@ -36,14 +37,16 @@
]
-#Def of the python command letMove() :
-# Will move "by hand" the upperSphere towards or away from the lower one.
+#Def of the python commands which will impose required displacements to the moving sphere
def letMove():#Load for the first 10 iterations, unload for the 7 following iterations, then reload
+ if mode=='normal':
vImposed=[0,-1,0]
if O.iter < 25 and O.iter>14:
vImposed=[0,1,0]
- upperSphere.state.vel=vImposed
- upperSphere.state.pos=upperSphere.state.pos+upperSphere.state.vel*O.dt
+ if mode=='tangential':
+ vImposed=[1,0,0]
+ upperSphere.state.vel=vImposed
+ upperSphere.state.pos=upperSphere.state.pos+upperSphere.state.vel*O.dt
#Def of the python command defData() : which will be used once the interaction will really exist
@@ -59,10 +62,15 @@
# ------ Test of the law in the normal direction, using python commands to let move ------ #
print 'Beginning of normal loading'
+mode='normal'
O.dt=1e-5
+#print 'Ok 1'
-yade.qt.View()
+#yade.qt.View()
+print 'Going to first cycles, for free'
+#stop
O.run(2,True) #cycles "for free", so that the interaction between spheres will be defined (with his physics and so on)
+print 'First cycles done'
O.engines=O.engines+[PyRunner(iterPeriod=1,command='defData()')]
O.run(40,True)
@@ -71,33 +79,30 @@
# define of the plots to be made : un(step), and Fn(un)
plot.plots={'step':('unTrue','torque',),'unPerso':('normFn',),'unTrue':('normFnBis',)}
-plot.plot()
-raw_input()
+plot.plot(noShow=False, subPlots=False)
print 'Type Return to go ahead'
print ''
+raw_input()
#NB : these different unTrue and unPerso illustrate the definition of penetrationDepth really used in the code (computed in Ig2_Sphere_Sphere_ScGeom) which is slightly different from R1 + R2 - Distance (see for example this "shift2"). According to the really used penetrationDepth, Fn evolves as it should
#O.saveTmp('EndComp')
#O.save('FinN_I_Test.xml')
-# ------ Test of the law in the tangential direction, using StepDisplacer ------ #
+# ------ Test of the law in the tangential direction, still with python function ------ #
print 'Beginning of tangential loading'
-
-dpos=Vector3.Zero
-Vector3.__init__(dpos,1*O.dt,0,0)
-O.engines=O.engines[:3]+[StepDisplacer(ids=[1],mov=dpos,setVelocities=True)]+O.engines[4:]
+mode='tangential'
O.run(1000)
print 'End of tangential loading'
plot.plots={'step':('gamma','torque',),'gamma':('fx',)}
plot.plot()
+print 'Type Return to go ahead'
raw_input()
-print 'Type Return to go ahead'
plot.plots={'normFn':('fx',)}
plot.plot()
+print 'Type Return to go ahead'
raw_input()
-print 'Type Return to go ahead'
#pylab.show() #to pause on the plot window. Effective only first time
print ''
@@ -109,13 +114,13 @@
# - during this phase O.forces.t(0).norm() / O.forces.f(0)[0] seems to increase between 0.502 and 0.507 (according to r=[T[i]/F[i] for i in range(50,T.__len__()-20) ])
# Please note, to explain this, that Fx = O.forces.f(0)[0] is more and more different from Ft, from which we can expect Ft = Torque /Radius
-## ------ Test of the law for the moment, using blockedDOF_s ------ #
+## ------ Test of the law for the moment, using blockedDOF_s and NewtonIntegrator ------ #
#O.loadTmp('EndComp')
print 'Beginning of rotationnal loading'
##To use blockedDOF_s, the body has to be dynamic....
upperSphere.dynamic=True
-upperSphere.state.blockedDOFs='x','rx','y','ry','z','rz'
+upperSphere.state.blockedDOFs='xyzXYZ'
upperSphere.state.angVel=Vector3(0,0,1)
upperSphere.state.vel=Vector3(0,0,0)
i=O.interactions[1,0]
@@ -136,3 +141,6 @@
print 'Type Return to go ahead'
#-- Comments : TO DO
+
+
+
=== modified file 'scripts/simpleShear.py'
--- scripts/simpleShear.py 2010-11-03 15:50:02 +0000
+++ scripts/simpleShear.py 2011-05-03 15:04:37 +0000
@@ -67,7 +67,7 @@
[Ip2_2xNormalInelasticMat_NormalInelasticityPhys()],
[Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity()]
),
- NewtonIntegrator(damping=.2)
+ NewtonIntegrator(damping=.1)
,PyRunner(iterPeriod=50,command='defData()')
]
=== added file 'scripts/test/checks/checkTestNormalInelasticity.py'
--- scripts/test/checks/checkTestNormalInelasticity.py 1970-01-01 00:00:00 +0000
+++ scripts/test/checks/checkTestNormalInelasticity.py 2011-05-03 15:04:37 +0000
@@ -0,0 +1,93 @@
+# -*- coding: utf-8 -*-
+
+# Check test version of normalInelasticityTest
+
+
+#Def of the material which will be used
+O.materials.append(NormalInelasticMat(density=2600,young=4.0e9,poisson=.04,frictionAngle=.6,coeff_dech=3.0,label='Materiau1'))
+
+#Def of the bodies of the simulations : 2 spheres, with names which will be useful after
+O.bodies.append(utils.sphere([0,0,0], 1, fixed=True, wire=False, color=None, highlight=False)) #'Materiau1', as the latest material defined, will be used
+O.bodies.append(utils.sphere([0,2,0], 1, fixed=True, wire=False, color=None, highlight=False))
+
+lowerSphere=O.bodies[0]
+upperSphere=O.bodies[1]
+
+
+#Def of the engines taking part to the simulation loop
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Sphere_Aabb()],verletDist=0),# use of verletDist>0 without NewtonIntegrator let crash ??
+ InteractionLoop(
+ [Ig2_Sphere_Sphere_ScGeom6D()],
+ [Ip2_2xNormalInelasticMat_NormalInelasticityPhys(betaR=0.24)],
+ [Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity()]
+ ),
+ PyRunner(iterPeriod=1,command='letMove()')
+ ]
+
+
+#Def of the python commands which will impose required displacements to the moving sphere
+def letMove():#Load for the first 10 iterations, unload for the 7 following iterations, then reload
+ if mode=='normal':
+ vImposed=[0,-1,0]
+ if O.iter < 25 and O.iter>14:
+ vImposed=[0,1,0]
+ if mode=='tangential':
+ vImposed=[1,0,0]
+ upperSphere.state.vel=vImposed
+ upperSphere.state.pos=upperSphere.state.pos+upperSphere.state.vel*O.dt
+
+
+
+
+# ------ Test of the law in the normal direction, using python commands to let move ------ #
+mode='normal'
+O.dt=1e-5
+O.run(40,True)
+
+
+
+# ------ Test of the law in the tangential direction, still with python function ------ #
+
+mode='tangential'
+O.run(1000)
+
+
+## ------ Test of the law for the moment, using blockedDOF_s and NewtonIntegrator ------ #
+
+##To use blockedDOF_s, the body has to be dynamic....
+upperSphere.dynamic=True
+upperSphere.state.blockedDOFs='xyzXYZ'
+upperSphere.state.angVel=Vector3(0,0,1)
+upperSphere.state.vel=Vector3(0,0,0)
+
+O.engines=O.engines[:3]+[NewtonIntegrator()]
+
+
+
+O.run(8000,True)
+
+
+
+
+# Reference value of force acting on upperSphere, (r2849)
+fRef=Vector3(-519943.97434309247,760000.00000497897,0)
+# Reference value of torque acting on upperSphere, (r2849)
+tRef=Vector3(0,0,-1279894.5796705084)
+
+#Actual values:
+f=O.forces.f(1)
+t=O.forces.t(1)
+
+dF=fRef-f
+dT=tRef-t
+
+tolerance=0.01
+if (dF.norm()>tolerance or dT.norm()>tolerance):
+ print "Regression concerning normalInelasticity-test. We compare loads acting on one sphere, moving with respect to another. At the end of the script (after some relative movements) we have a difference, compared to reference values, of"
+ print "- force (the norm) equal to ", dF.norm()
+ print "- torque (the norm) equal to ", dT.norm()
+ print "Whereas tolerance is ", tolerance
+ resultStatus +=1
+
