yade-dev team mailing list archive
-
yade-dev team
-
Mailing list archive
-
Message #11451
[Branch ~yade-pkg/yade/git-trunk] Rev 3447: Merge github.com:yade/trunk into chaoUnsat
Merge authors:
Alexander Eulitz [Eugen] (kubeu)
Alexander Eulitz [Eugen] (kubeu)
Anton Gladky (gladky-anton)
Anton Gladky (gladky-anton)
Bruno Chareyre (bruno-chareyre)...
------------------------------------------------------------
revno: 3447 [merge]
committer: cyuan <chaoyuan2012@xxxxxxxxx>
timestamp: Sat 2014-08-02 09:15:35 +0800
message:
Merge github.com:yade/trunk into chaoUnsat
removed:
core/Bound.cpp
core/Dispatcher.cpp
core/EnergyTracker.cpp
core/Engine.cpp
core/FrontEnd.cpp
core/PartialEngine.cpp
core/Shape.cpp
core/TimeStepper.cpp
lib/multimethods/Indexable.cpp
pkg/common/Aabb.cpp
pkg/common/Bo1_Box_Aabb.cpp
pkg/common/Bo1_Box_Aabb.hpp
pkg/common/Bo1_Facet_Aabb.cpp
pkg/common/Bo1_Facet_Aabb.hpp
pkg/common/Bo1_Sphere_Aabb.cpp
pkg/common/Bo1_Sphere_Aabb.hpp
pkg/common/BoundaryController.cpp
pkg/common/Box.cpp
pkg/common/Callbacks.cpp
pkg/common/CylScGeom6D.cpp
pkg/common/ElastMat.cpp
pkg/common/FieldApplier.cpp
pkg/common/ForceResetter.cpp
pkg/common/GLDrawFunctors.cpp
pkg/common/Gl1_Aabb.cpp
pkg/common/Gl1_Aabb.hpp
pkg/common/Gl1_Box.cpp
pkg/common/Gl1_Box.hpp
pkg/common/Gl1_Facet.cpp
pkg/common/Gl1_Facet.hpp
pkg/common/Gl1_Sphere.cpp
pkg/common/Gl1_Sphere.hpp
pkg/common/NormShearPhys.cpp
pkg/common/PeriodicEngines.cpp
pkg/common/PyRunner.cpp
pkg/common/Recorder.cpp
pkg/common/Sphere.cpp
pkg/common/StepDisplacer.cpp
pkg/common/TorqueEngine.cpp
pkg/dem/DemXDofGeom.cpp
pkg/lbm/LBMbody.cpp
pkg/lbm/LBMlink.cpp
added:
examples/bouncingbubble.py
pkg/common/Bo1_Aabb.cpp
pkg/common/Bo1_Aabb.hpp
pkg/common/Gl1_Primitives.cpp
pkg/common/Gl1_Primitives.hpp
pkg/common/common.cpp
modified:
CMakeLists.txt
core/BodyContainer.cpp
core/Dispatcher.hpp
core/EnergyTracker.hpp
core/Engine.hpp
core/FrontEnd.hpp
core/Functor.hpp
core/Material.cpp
core/Material.hpp
core/Shape.hpp
core/TimeStepper.hpp
core/main/main.py.in
doc/references.bib
doc/sphinx/conf.py
doc/sphinx/installation.rst
doc/yade-conferences.bib
examples/gts-horse/gts-horse.py
examples/stl-gts/gts-stl.py
examples/triax-tutorial/script-session1.py
lib/base/Math.cpp
lib/multimethods/Indexable.hpp
pkg/common/Aabb.hpp
pkg/common/BoundaryController.hpp
pkg/common/Box.hpp
pkg/common/Callbacks.hpp
pkg/common/CylScGeom6D.hpp
pkg/common/Cylinder.cpp
pkg/common/Cylinder.hpp
pkg/common/Dispatching.hpp
pkg/common/ElastMat.hpp
pkg/common/FieldApplier.hpp
pkg/common/ForceResetter.hpp
pkg/common/Grid.cpp
pkg/common/Grid.hpp
pkg/common/InteractionLoop.cpp
pkg/common/MatchMaker.cpp
pkg/common/MatchMaker.hpp
pkg/common/NormShearPhys.hpp
pkg/common/OpenGLRenderer.cpp
pkg/common/ParallelEngine.cpp
pkg/common/PeriodicEngines.hpp
pkg/common/PyRunner.hpp
pkg/common/Recorder.hpp
pkg/common/Sphere.hpp
pkg/common/StepDisplacer.hpp
pkg/common/TorqueEngine.hpp
pkg/dem/BubbleMat.cpp
pkg/dem/BubbleMat.hpp
pkg/dem/CapillaryTriaxialTest.cpp
pkg/dem/CohesiveFrictionalContactLaw.cpp
pkg/dem/CohesiveFrictionalContactLaw.hpp
pkg/dem/CohesiveTriaxialTest.cpp
pkg/dem/ConcretePM.cpp
pkg/dem/ConcretePM.hpp
pkg/dem/DemXDofGeom.hpp
pkg/dem/DomainLimiter.cpp
pkg/dem/ElasticContactLaw.cpp
pkg/dem/ElasticContactLaw.hpp
pkg/dem/FrictViscoPM.cpp
pkg/dem/FrictViscoPM.hpp
pkg/dem/HertzMindlin.cpp
pkg/dem/HertzMindlin.hpp
pkg/dem/InelastCohFrictPM.cpp
pkg/dem/InelastCohFrictPM.hpp
pkg/dem/JointedCohesiveFrictionalPM.cpp
pkg/dem/JointedCohesiveFrictionalPM.hpp
pkg/dem/L3Geom.cpp
pkg/dem/L3Geom.hpp
pkg/dem/LudingPM.cpp
pkg/dem/LudingPM.hpp
pkg/dem/NormalInelasticPM.cpp
pkg/dem/NormalInelasticPM.hpp
pkg/dem/PeriIsoCompressor.hpp
pkg/dem/Polyhedra.cpp
pkg/dem/Polyhedra.hpp
pkg/dem/Shop.hpp
pkg/dem/Shop_01.cpp
pkg/dem/Shop_02.cpp
pkg/dem/SimpleShear.cpp
pkg/dem/Tetra.cpp
pkg/dem/Tetra.hpp
pkg/dem/TriaxialStressController.cpp
pkg/dem/TriaxialStressController.hpp
pkg/dem/TriaxialTest.cpp
pkg/dem/VTKRecorder.cpp
pkg/dem/VTKRecorder.hpp
pkg/dem/ViscoelasticCapillarPM.cpp
pkg/dem/ViscoelasticCapillarPM.hpp
pkg/dem/ViscoelasticPM.cpp
pkg/dem/ViscoelasticPM.hpp
pkg/dem/WirePM.cpp
pkg/dem/WirePM.hpp
pkg/lbm/LBMbody.hpp
pkg/lbm/LBMlink.hpp
pkg/lbm/LBMnode.cpp
pkg/pfv/FlowEngine.hpp.in
pkg/pfv/FlowEngine.ipp.in
py/_utils.cpp
py/tests/__init__.py
py/utils.py
scripts/checks-and-tests/checks/DEM-PFV-check.py
scripts/checks-and-tests/checks/checkList.py
--
lp:yade
https://code.launchpad.net/~yade-pkg/yade/git-trunk
Your team Yade developers is subscribed to branch lp:yade.
To unsubscribe from this branch go to https://code.launchpad.net/~yade-pkg/yade/git-trunk/+edit-subscription
=== modified file 'CMakeLists.txt'
--- CMakeLists.txt 2014-07-17 04:21:07 +0000
+++ CMakeLists.txt 2014-08-02 01:15:35 +0000
@@ -433,7 +433,8 @@
FILE(GLOB_RECURSE SRC_PKG "pkg/*.cpp")
FILE(GLOB SRC_LIB "lib/*.cpp")
-SET(SRC_LIB "${SRC_LIB};lib/base/Math.cpp;lib/factory/ClassFactory.cpp;lib/factory/DynLibManager.cpp;lib/multimethods/Indexable.cpp;lib/serialization/Serializable.cpp;lib/pyutil/gil.cpp;core/main/pyboot.cpp;${GUI_SRC_LIB};${CGAL_SRC_LIB}")
+SET(SRC_LIB "${SRC_LIB};lib/base/Math.cpp;lib/factory/ClassFactory.cpp;lib/factory/DynLibManager.cpp")
+SET(SRC_LIB "${SRC_LIB};lib/serialization/Serializable.cpp;lib/pyutil/gil.cpp;core/main/pyboot.cpp;${GUI_SRC_LIB};${CGAL_SRC_LIB}")
#===========================================================
=== modified file 'core/BodyContainer.cpp'
--- core/BodyContainer.cpp 2014-07-17 06:34:00 +0000
+++ core/BodyContainer.cpp 2014-07-17 08:05:16 +0000
@@ -55,7 +55,8 @@
for(Body::MapId2IntrT::iterator it=b->intrs.begin(),end=b->intrs.end(); it!=end; ++it) { //Iterate over all body's interactions
scene->interactions->requestErase((*it).second);
}
- body[id]=nullptr;
+
+ body[id].reset();
return true;
}
=== removed file 'core/Bound.cpp'
--- core/Bound.cpp 2009-12-04 23:07:34 +0000
+++ core/Bound.cpp 1970-01-01 00:00:00 +0000
@@ -1,12 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Olivier Galizzi *
-* olivier.galizzi@xxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#include "Bound.hpp"
-
-
-
=== removed file 'core/Dispatcher.cpp'
--- core/Dispatcher.cpp 2014-07-03 17:20:40 +0000
+++ core/Dispatcher.cpp 1970-01-01 00:00:00 +0000
@@ -1,15 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Olivier Galizzi *
-* olivier.galizzi@xxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-
-#include<yade/core/Functor.hpp>
-Functor::~Functor(){}; // vtable
-
-#include "Dispatcher.hpp"
-
-Dispatcher::~Dispatcher(){}
=== modified file 'core/Dispatcher.hpp'
--- core/Dispatcher.hpp 2014-07-03 17:20:40 +0000
+++ core/Dispatcher.hpp 2014-07-19 19:52:41 +0000
@@ -24,7 +24,7 @@
virtual int getDimension() { throw; };
virtual string getBaseClassType(unsigned int ) { throw; };
//
- virtual ~Dispatcher();
+ virtual ~Dispatcher() {};
YADE_CLASS_BASE_DOC(Dispatcher,Engine,"Engine dispatching control to its associated functors, based on types of argument it receives. This abstract base class provides no functionality in itself.")
};
REGISTER_SERIALIZABLE(Dispatcher);
=== removed file 'core/EnergyTracker.cpp'
--- core/EnergyTracker.cpp 2011-02-15 11:36:29 +0000
+++ core/EnergyTracker.cpp 1970-01-01 00:00:00 +0000
@@ -1,12 +0,0 @@
-#include<yade/core/EnergyTracker.hpp>
-
-Real EnergyTracker::total() const { Real ret=0; size_t sz=energies.size(); for(size_t id=0; id<sz; id++) ret+=energies.get(id); return ret; }
-py::list EnergyTracker::keys_py() const { py::list ret; FOREACH(pairStringInt p, names) ret.append(p.first); return ret; }
-py::list EnergyTracker::items_py() const { py::list ret; FOREACH(pairStringInt p, names) ret.append(py::make_tuple(p.first,energies.get(p.second))); return ret; }
-py::dict EnergyTracker::perThreadData() const {
- py::dict ret;
- std::vector<std::vector<Real> > dta=energies.getPerThreadData();
- FOREACH(pairStringInt p,names) ret[p.first]=dta[p.second];
- return ret;
-}
-
=== modified file 'core/EnergyTracker.hpp'
--- core/EnergyTracker.hpp 2014-07-03 17:20:40 +0000
+++ core/EnergyTracker.hpp 2014-07-19 19:52:41 +0000
@@ -38,10 +38,15 @@
void clear(){ energies.clear(); names.clear(); resetStep.clear();}
void resetResettables(){ size_t sz=energies.size(); for(size_t id=0; id<sz; id++){ if(resetStep[id]) energies.reset(id); } }
- Real total() const;
- py::list keys_py() const;
- py::list items_py() const;
- py::dict perThreadData() const;
+ Real total() const { Real ret=0; size_t sz=energies.size(); for(size_t id=0; id<sz; id++) ret+=energies.get(id); return ret; };
+ py::list keys_py() const { py::list ret; FOREACH(pairStringInt p, names) ret.append(p.first); return ret; };
+ py::list items_py() const { py::list ret; FOREACH(pairStringInt p, names) ret.append(py::make_tuple(p.first,energies.get(p.second))); return ret; };
+ py::dict perThreadData() const {
+ py::dict ret;
+ std::vector<std::vector<Real> > dta=energies.getPerThreadData();
+ FOREACH(pairStringInt p,names) ret[p.first]=dta[p.second];
+ return ret;
+ };
typedef std::map<std::string,int> mapStringInt;
typedef std::pair<std::string,int> pairStringInt;
=== removed file 'core/Engine.cpp'
--- core/Engine.cpp 2010-09-02 09:11:04 +0000
+++ core/Engine.cpp 1970-01-01 00:00:00 +0000
@@ -1,13 +0,0 @@
-#include<yade/core/Engine.hpp>
-
-CREATE_LOGGER(Engine);
-
-void Engine::action(){
- LOG_FATAL("Engine "<<getClassName()<<" calling virtual method Engine::action(). Please submit bug report at http://bugs.launchpad.net/yade.";);
- throw std::logic_error("Engine::action() called.");
-}
-
-void Engine::explicitAction(){
- scene=Omega::instance().getScene().get(); action();
-}
-
=== modified file 'core/Engine.hpp'
--- core/Engine.hpp 2014-07-03 17:20:40 +0000
+++ core/Engine.hpp 2014-07-19 19:52:41 +0000
@@ -18,6 +18,8 @@
class Body;
class Scene;
+CREATE_LOGGER(Engine);
+
class Engine: public Serializable{
public:
// pointer to the simulation, set at every step by Scene::moveToNextTimeStep
@@ -29,14 +31,17 @@
virtual ~Engine() {};
virtual bool isActivated() { return true; };
- virtual void action();
+ virtual void action() {
+ LOG_FATAL("Engine "<<getClassName()<<" calling virtual method Engine::action(). Please submit bug report at http://bugs.launchpad.net/yade.";);
+ throw std::logic_error("Engine::action() called.");
+ }
private:
// py access funcs
TimingInfo::delta timingInfo_nsec_get(){return timingInfo.nsec;};
void timingInfo_nsec_set(TimingInfo::delta d){ timingInfo.nsec=d;}
long timingInfo_nExec_get(){return timingInfo.nExec;};
void timingInfo_nExec_set(long d){ timingInfo.nExec=d;}
- void explicitAction();
+ void explicitAction() {scene=Omega::instance().getScene().get(); action();};
DECLARE_LOGGER;
=== removed file 'core/FrontEnd.cpp'
--- core/FrontEnd.cpp 2007-03-31 06:18:17 +0000
+++ core/FrontEnd.cpp 1970-01-01 00:00:00 +0000
@@ -1,21 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Olivier Galizzi *
-* olivier.galizzi@xxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#include "FrontEnd.hpp"
-
-
-FrontEnd::FrontEnd()
-{
-}
-
-
-FrontEnd::~FrontEnd()
-{
-
-}
-
=== modified file 'core/FrontEnd.hpp'
--- core/FrontEnd.hpp 2010-11-07 11:46:20 +0000
+++ core/FrontEnd.hpp 2014-07-19 19:52:41 +0000
@@ -15,8 +15,8 @@
class FrontEnd : public Factorable
{
public :
- FrontEnd ();
- virtual ~FrontEnd ();
+ FrontEnd () {};
+ virtual ~FrontEnd () {};
virtual int run(int , char * []) { return -1;};
// called before actually invoking it
=== modified file 'core/Functor.hpp'
--- core/Functor.hpp 2014-07-02 16:11:24 +0000
+++ core/Functor.hpp 2014-07-19 19:52:41 +0000
@@ -23,7 +23,7 @@
shared_ptr<TimingDeltas> timingDeltas;
//! updated before every dispatch loop by the dispatcher; DO NOT ABUSE access to scene, except for getting global variables like scene->dt.
Scene* scene;
- virtual ~Functor(); // defined in Dispatcher.cpp
+ virtual ~Functor() {}; // defined in Dispatcher.cpp
YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(Functor,Serializable,"Function-like object that is called by Dispatcher, if types of arguments match those the Functor declares to accept.",
((string,label,,,"Textual label for this object; must be a valid python identifier, you can refer to it directly from python.")),
/*ctor*/
=== modified file 'core/Material.cpp'
--- core/Material.cpp 2014-07-03 17:20:40 +0000
+++ core/Material.cpp 2014-07-19 19:52:41 +0000
@@ -2,8 +2,6 @@
#include<yade/core/Material.hpp>
#include<yade/core/Scene.hpp>
-Material::~Material(){}
-
const shared_ptr<Material> Material::byId(int id, Scene* w_){
Scene* w=w_?w_:Omega::instance().getScene().get();
assert(id>=0 && (size_t)id<w->materials.size());
=== modified file 'core/Material.hpp'
--- core/Material.hpp 2014-07-03 17:20:40 +0000
+++ core/Material.hpp 2014-07-19 19:52:41 +0000
@@ -14,7 +14,7 @@
*/
class Material: public Serializable, public Indexable{
public:
- virtual ~Material();
+ virtual ~Material() {};
//! Function to return empty default-initialized instance of State that
// is supposed to go along with this Material. Don't override unless you need
=== removed file 'core/PartialEngine.cpp'
--- core/PartialEngine.cpp 2010-03-20 12:40:44 +0000
+++ core/PartialEngine.cpp 1970-01-01 00:00:00 +0000
@@ -1,10 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Olivier Galizzi *
-* olivier.galizzi@xxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#include "PartialEngine.hpp"
-
=== removed file 'core/Shape.cpp'
--- core/Shape.cpp 2010-02-06 16:01:07 +0000
+++ core/Shape.cpp 1970-01-01 00:00:00 +0000
@@ -1,2 +0,0 @@
-#include<yade/core/Shape.hpp>
-Shape::~Shape(){}
=== modified file 'core/Shape.hpp'
--- core/Shape.hpp 2010-11-07 11:46:20 +0000
+++ core/Shape.hpp 2014-07-19 19:52:41 +0000
@@ -18,7 +18,7 @@
class Shape: public Serializable, public Indexable {
public:
- ~Shape(); // vtable
+ ~Shape() {}; // vtable
#ifdef BV_FUNCTOR_CACHE
shared_ptr<BoundFunctor> boundFunctor;
#endif
=== removed file 'core/TimeStepper.cpp'
--- core/TimeStepper.cpp 2010-08-16 21:31:08 +0000
+++ core/TimeStepper.cpp 1970-01-01 00:00:00 +0000
@@ -1,29 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Olivier Galizzi *
-* olivier.galizzi@xxxxxxx *
-* Copyright (C) 2004 by Janek Kozicki *
-* cosurgi@xxxxxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#include<yade/core/TimeStepper.hpp>
-#include<yade/core/GlobalEngine.hpp>
-#include<yade/core/Scene.hpp>
-
-bool TimeStepper::isActivated()
-{
- return (active && (scene->iter % timeStepUpdateInterval == 0));
-}
-
-
-
-void TimeStepper::setActive(bool a, int nb)
-{
- active = a;
- if (nb>0)
- timeStepUpdateInterval = (unsigned int)nb;
-}
-
-
=== modified file 'core/TimeStepper.hpp'
--- core/TimeStepper.hpp 2010-08-24 12:54:14 +0000
+++ core/TimeStepper.hpp 2014-07-19 19:52:41 +0000
@@ -12,15 +12,16 @@
#include <vector>
#include "Interaction.hpp"
#include "GlobalEngine.hpp"
+#include "Scene.hpp"
class Body;
class TimeStepper: public GlobalEngine{
public:
virtual void computeTimeStep(Scene* ) { throw; };
- virtual bool isActivated();
+ virtual bool isActivated() {return (active && (scene->iter % timeStepUpdateInterval == 0));};
virtual void action() { computeTimeStep(scene);} ;
- void setActive(bool a, int nb=-1);
+ void setActive(bool a, int nb=-1) {active = a; if (nb>0) {timeStepUpdateInterval = (unsigned int)nb;}}
YADE_CLASS_BASE_DOC_ATTRS(
TimeStepper,GlobalEngine,"Engine defining time-step (fundamental class)",
=== modified file 'core/main/main.py.in'
--- core/main/main.py.in 2014-05-27 11:24:16 +0000
+++ core/main/main.py.in 2014-07-22 17:23:17 +0000
@@ -48,7 +48,7 @@
)
par.add_argument('-v','--version',help='Print version and exit.',dest='version',action='store_true')
par.add_argument('-j','--threads',help='Number of OpenMP threads to run; defaults to 1. Equivalent to setting OMP_NUM_THREADS environment variable.',dest='threads',type=int)
-par.add_argument('--cores',help='Set number of OpenMP threads (as \-\-threads) and in addition set affinity of threads to the cores given.',dest='cores',type=str)
+par.add_argument('--cores',help='Set number of OpenMP threads (as \-\-threads) and in addition set affinity of threads to the cores given. Please provide a string with comma-separated core-ids.',dest='cores',type=str)
par.add_argument('--update',help='Update deprecated class names in given script(s) using text search & replace. Changed files will be backed up with ~ suffix. Exit when done without running any simulation.',dest='updateScripts',nargs='+')
par.add_argument('--nice',help='Increase nice level (i.e. decrease priority) by given number.',dest='nice',type=int)
par.add_argument('-x',help='Exit when the script finishes',dest='exitAfter',action='store_true')
@@ -107,7 +107,7 @@
print 'WARNING: compiled without OpenMP, -j/--threads/--cores have no effect.'
# OpenMP env variables must be set before loading yade libs ("import yade" below)
-# changes have no effeect after libgomp initializes
+# changes have no effect after libgomp initializes
if opts.cores:
if opts.threads: print 'WARNING: --threads ignored, since --cores specified.'
try:
@@ -115,13 +115,13 @@
except ValueError:
raise ValueError('Invalid --cores specification %s, should be a comma-separated list of non-negative integers'%opts.cores)
opts.nthreads=len(cores)
- os.environ['GOMP_CPU_AFFINITY']=' '.join([str(cores[0])]+[str(c) for c in cores])
+ os.environ['GOMP_CPU_AFFINITY']=' '.join([str(c) for c in cores])
os.environ['OMP_NUM_THREADS']=str(len(cores))
elif opts.threads: os.environ['OMP_NUM_THREADS']=str(opts.threads)
else: os.environ['OMP_NUM_THREADS']='1'
if __name__ == "__main__": # do not print this while importing yade in other python application
- sys.stderr.write('Welcome to Yade '+version+debugbuild+'\n')
+ sys.stdout.write('Welcome to Yade '+version+debugbuild+'\n')
# initialization and c++ plugins import
import yade
=== modified file 'doc/references.bib'
--- doc/references.bib 2014-05-08 05:57:04 +0000
+++ doc/references.bib 2014-07-28 06:24:36 +0000
@@ -34,6 +34,16 @@
year = "2005"
}
+@Article{ Chan2011,
+ title = "Film drainage and coalescence between deformable drops and bubbles.",
+ author = "D. Chan and E. Klaseboer and R. Manica",
+ journal = "Soft Matter",
+ pages = "2235-2264",
+ volume = "7",
+ number = "6",
+ year = "2011"
+}
+
@Article{ Chareyre2002a,
title = "Theoretical versus experimental modeling of the anchorage capacity of geotextiles in trenches.",
author = "B. Chareyre and L. Briancon and P. Villard",
=== modified file 'doc/sphinx/conf.py'
--- doc/sphinx/conf.py 2014-06-07 08:15:18 +0000
+++ doc/sphinx/conf.py 2014-07-23 09:07:00 +0000
@@ -695,7 +695,7 @@
\text{\sffamily\bfseries\LARGE Authors}\\
\\
\text{\sffamily\bfseries\Large V\'{a}clav \v{S}milauer}\\
-\text{\sffamily\Large University of Innsbruck}\\
+\text{\sffamily\Large Freelance consultant (http://woodem.eu)}\\
\\
\text{\sffamily\bfseries\Large Emanuele Catalano}\\
\text{\sffamily\Large Grenoble INP, UJF, CNRS, lab. 3SR}\\
=== modified file 'doc/sphinx/installation.rst'
--- doc/sphinx/installation.rst 2014-06-25 18:18:39 +0000
+++ doc/sphinx/installation.rst 2014-07-25 08:01:21 +0000
@@ -179,12 +179,12 @@
^^^^^^^^^^^
You should create a separate build-place-folder, where Yade will be configured
-and where the source code will be compiled. Here is an example for a folderstructure:
+and where the source code will be compiled. Here is an example for a folder structure::
- myYade/ ## base directory
- trunk/ ## folder for sourcecode in which you use github
- build/ ## folder in which sources will be compiled; build-directory; use cmake here
- install/ ## installfolder
+ myYade/ ## base directory
+ trunk/ ## folder for sourcecode in which you use github
+ build/ ## folder in which sources will be compiled; build-directory; use cmake here
+ install/ ## install folder
Then inside this build-directory you should start cmake to configure the compilation process::
@@ -192,7 +192,7 @@
For the folder structure given above call the following command in folder "build":
- cmake -DINSTALL_PREFIX=../install ../trunk
+ cmake -DINSTALL_PREFIX=../install ../trunk
Additional options can be configured in the same line with the following
syntax::
@@ -227,19 +227,22 @@
make
-Installing performs with the following command::
-
- make install
-
-The "install" command will in fact also recompile if source files have been modified.
-Hence there is no absolute need to type the two commands separately.
-
The compilation process can take a long time, be patient. An additional
parameter on many cores systems ``-j`` can be added to decrease compilation time
and split the compilation on many cores. For example, on 4-core machines
it would be reasonable to set the parameter ``-j4``. Note, the Yade requires
approximately 2GB/core for compilation, otherwise the swap-file will be used
-and a compilation time dramatically increases. After compilation finished successfully
+and a compilation time dramatically increases.
+
+Installing performs with the following command::
+
+ make install
+
+The "install" command will in fact also recompile if source files have been modified.
+Hence there is no absolute need to type the two commands separately. You may receive make errors if you don't permission to write into the target folder.
+These errors are not critical but without writing permissions Yade won't be installed in /usr/local/bin/.
+
+After compilation finished successfully
the new built can be started by navigating to /path/to/installfolder/bin and calling yade via (based on version yade-2014-02-20.git-a7048f4)::
cd /path/to/installfolder/bin
@@ -247,7 +250,7 @@
For building the documentation you should at first execute the command "make install"
and then "make doc" to build it. The generated files will be stored in your current
-build directory/doc/sphinx/_build.
+build directory/doc/sphinx/_build. Once again writing permissions are necessary for installing into /usr/local/share/doc/.
"make manpage" command generates and moves manpages in a standard place.
"make check" command executes standard test to check the functionality of compiled program.
=== modified file 'doc/yade-conferences.bib'
--- doc/yade-conferences.bib 2014-04-02 15:19:41 +0000
+++ doc/yade-conferences.bib 2014-07-18 18:16:01 +0000
@@ -611,3 +611,23 @@
document_type={Conference Paper},
source={Scopus},
}
+
+@InProceedings{Yuan2014,
+ author={Chao Yuan and Bruno Chareyre and Félix Darve},
+ title={Pore-Scale Simulations of Drainage for Two-phase Flow in Dense Sphere Packings},
+ journal={XX . International Conference on Computational Methods in Water Resources},
+ year={2014},
+ pages={232},
+ url={http://cmwr14.de/images/bookofabstracts/CMWR14BookofAbstracts.pdf}
+}
+
+@InProceedings{Sweijen2014,
+ author={Thomas Sweijen and Majid Hassanizadeh and Bruno Chareyre},
+ title={Pore-Scale modeling of swelling porous media; application to super absorbent polymers},
+ journal={XX . International Conference on Computational Methods in Water Resources},
+ year={2014},
+ pages={227},
+ url={http://cmwr14.de/images/bookofabstracts/CMWR14BookofAbstracts.pdf}
+}
+
+
=== added file 'examples/bouncingbubble.py'
--- examples/bouncingbubble.py 1970-01-01 00:00:00 +0000
+++ examples/bouncingbubble.py 2014-07-29 19:02:55 +0000
@@ -0,0 +1,20 @@
+
+#Simulates a 5mm Diameter bubble in water rising and colliding with another bubble of the same diameter
+
+rad = 2.5e-3
+#O.materials.append(FrictMat(young=1e3,density=1000))
+O.materials.append(BubbleMat(density=1000,surfaceTension=71.97e-3))
+O.bodies.append([
+ utils.sphere(center=(0,0,0),radius=rad,fixed=True),
+ utils.sphere((0,0,-2*rad*1.1),rad)
+])
+O.dt = 1e-7
+
+O.engines=[
+ ForceResetter(),
+ InsertionSortCollider([Bo1_Sphere_Aabb()]),
+ InteractionLoop([Ig2_Sphere_Sphere_ScGeom()],[Ip2_BubbleMat_BubbleMat_BubblePhys()],[Law2_ScGeom_BubblePhys_Bubble()]),
+# InteractionLoop([Ig2_Sphere_Sphere_ScGeom()],[Ip2_FrictMat_FrictMat_FrictPhys()],[Law2_ScGeom_FrictPhys_CundallStrack()]),
+ NewtonIntegrator(damping=0.1,gravity=(0,0,9.81))
+]
+O.saveTmp()
=== modified file 'examples/gts-horse/gts-horse.py'
--- examples/gts-horse/gts-horse.py 2013-03-28 11:04:00 +0000
+++ examples/gts-horse/gts-horse.py 2014-07-29 17:26:54 +0000
@@ -40,9 +40,9 @@
ForceResetter(),
InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()],label='collider'),
InteractionLoop(
- [Ig2_Sphere_Sphere_L3Geom(),Ig2_Facet_Sphere_L3Geom()],
+ [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom()],
[Ip2_FrictMat_FrictMat_FrictPhys()],
- [Law2_L3Geom_FrictPhys_ElPerfPl()],
+ [Law2_ScGeom_FrictPhys_CundallStrack()],
),
NewtonIntegrator(damping=.1,gravity=[0,0,-5000]),
PyRunner(iterPeriod=1000,command='timing.stats(); O.pause();'),
=== modified file 'examples/stl-gts/gts-stl.py'
--- examples/stl-gts/gts-stl.py 2013-12-03 07:47:09 +0000
+++ examples/stl-gts/gts-stl.py 2014-07-29 17:26:54 +0000
@@ -20,9 +20,9 @@
ForceResetter(),
InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()],label='collider'),
InteractionLoop(
- [Ig2_Sphere_Sphere_L3Geom(),Ig2_Facet_Sphere_L3Geom()],
+ [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom()],
[Ip2_FrictMat_FrictMat_FrictPhys()],
- [Law2_L3Geom_FrictPhys_ElPerfPl()],
+ [Law2_ScGeom_FrictPhys_CundallStrack()],
),
NewtonIntegrator(damping=.1,gravity=[0,0,-500.0]),
PyRunner(iterPeriod=1000,command='timing.stats(); O.pause();'),
=== modified file 'examples/triax-tutorial/script-session1.py'
--- examples/triax-tutorial/script-session1.py 2014-01-03 20:43:21 +0000
+++ examples/triax-tutorial/script-session1.py 2014-08-01 15:37:51 +0000
@@ -68,6 +68,7 @@
## generate positions and input them in the simulation
sp.makeClumpCloud(mn,mx,[c1],periodic=False)
sp.toSimulation(material='spheres')
+ O.bodies.updateClumpProperties()#get more accurate clump masses/volumes/inertia
else:
sp.makeCloud(mn,mx,-1,0.3333,num_spheres,False, 0.95,seed=1) #"seed" make the "random" generation always the same
O.bodies.append([sphere(center,rad,material='spheres') for center,rad in sp])
=== modified file 'lib/base/Math.cpp'
--- lib/base/Math.cpp 2014-06-17 10:18:00 +0000
+++ lib/base/Math.cpp 2014-08-01 16:00:13 +0000
@@ -20,8 +20,8 @@
mask_t operator&(int i, const mask_t& g) { return g & i; }
mask_t operator|(const mask_t& g, int i) { return g | mask_t(i); }
mask_t operator|(int i, const mask_t& g) { return g | i; }
-bool operator||(const mask_t& g, bool b) { return (g == 0) || b; }
+bool operator||(const mask_t& g, bool b) { return (g != 0) || b; }
bool operator||(bool b, const mask_t& g) { return g || b; }
-bool operator&&(const mask_t& g, bool b) { return (g == 0) && b; }
+bool operator&&(const mask_t& g, bool b) { return (g != 0) && b; }
bool operator&&(bool b, const mask_t& g) { return g && b; }
#endif
=== removed file 'lib/multimethods/Indexable.cpp'
--- lib/multimethods/Indexable.cpp 2007-03-31 06:18:17 +0000
+++ lib/multimethods/Indexable.cpp 1970-01-01 00:00:00 +0000
@@ -1,34 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Janek Kozicki *
-* cosurgi@xxxxxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#include "Indexable.hpp"
-
-
-Indexable::Indexable ()
-{
-}
-
-
-Indexable::~Indexable ()
-{
-
-}
-
-
-void Indexable::createIndex ()
-{
- int& index = getClassIndex();
- if(index == -1) // assign new index
- {
- index = getMaxCurrentlyUsedClassIndex()+1;
- // so that other dispatchers will not fall in conflict with this index
- incrementMaxCurrentlyUsedClassIndex();
- }
-
-}
-
=== modified file 'lib/multimethods/Indexable.hpp'
--- lib/multimethods/Indexable.hpp 2012-02-16 16:05:15 +0000
+++ lib/multimethods/Indexable.hpp 2014-07-21 18:17:38 +0000
@@ -24,11 +24,19 @@
class Indexable
{
protected :
- void createIndex();
+ void createIndex() {
+ int& index = getClassIndex();
+ if(index == -1) // assign new index
+ {
+ index = getMaxCurrentlyUsedClassIndex()+1;
+ // so that other dispatchers will not fall in conflict with this index
+ incrementMaxCurrentlyUsedClassIndex();
+ }
+ }
public :
- Indexable ();
- virtual ~Indexable ();
+ Indexable () {};
+ virtual ~Indexable () {};
/// Returns the id of the current class. This id is set by a multimethod manager
virtual int& getClassIndex() { _THROW_NOT_OVERRIDDEN;};
=== removed file 'pkg/common/Aabb.cpp'
--- pkg/common/Aabb.cpp 2010-10-13 16:23:08 +0000
+++ pkg/common/Aabb.cpp 1970-01-01 00:00:00 +0000
@@ -1,13 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Olivier Galizzi *
-* olivier.galizzi@xxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#include "Aabb.hpp"
-
-Aabb::~Aabb(){}
-YADE_PLUGIN((Aabb));
-
=== modified file 'pkg/common/Aabb.hpp'
--- pkg/common/Aabb.hpp 2010-10-13 16:23:08 +0000
+++ pkg/common/Aabb.hpp 2014-07-19 19:52:41 +0000
@@ -18,11 +18,9 @@
*/
class Aabb : public Bound{
public :
- virtual ~Aabb();
+ virtual ~Aabb() {};
YADE_CLASS_BASE_DOC_ATTRS_CTOR(Aabb,Bound,"Axis-aligned bounding box, for use with :yref:`InsertionSortCollider`. (This class is quasi-redundant since min,max are already contained in :yref:`Bound` itself. That might change at some point, though.)",/*attrs*/,/*ctor*/createIndex(););
REGISTER_CLASS_INDEX(Aabb,Bound);
};
REGISTER_SERIALIZABLE(Aabb);
-
-
=== added file 'pkg/common/Bo1_Aabb.cpp'
--- pkg/common/Bo1_Aabb.cpp 1970-01-01 00:00:00 +0000
+++ pkg/common/Bo1_Aabb.cpp 2014-07-19 19:52:41 +0000
@@ -0,0 +1,87 @@
+/*************************************************************************
+* Copyright (C) 2004 by Olivier Galizzi *
+* olivier.galizzi@xxxxxxx *
+* *
+* This program is free software; it is licensed under the terms of the *
+* GNU General Public License v2 or later. See file LICENSE for details. *
+*************************************************************************/
+
+#include <yade/pkg/common/Bo1_Aabb.hpp>
+
+YADE_PLUGIN((Bo1_Sphere_Aabb)(Bo1_Facet_Aabb)(Bo1_Box_Aabb));
+
+void Bo1_Sphere_Aabb::go(const shared_ptr<Shape>& cm, shared_ptr<Bound>& bv, const Se3r& se3, const Body* b){
+ Sphere* sphere = static_cast<Sphere*>(cm.get());
+ if(!bv){ bv=shared_ptr<Bound>(new Aabb); }
+ Aabb* aabb=static_cast<Aabb*>(bv.get());
+ Vector3r halfSize = (aabbEnlargeFactor>0?aabbEnlargeFactor:1.)*Vector3r(sphere->radius,sphere->radius,sphere->radius);
+ if(!scene->isPeriodic){
+ aabb->min=se3.position-halfSize; aabb->max=se3.position+halfSize;
+ return;
+ }
+ // adjust box size along axes so that sphere doesn't stick out of the box even if sheared (i.e. parallelepiped)
+ if(scene->cell->hasShear()) {
+ Vector3r refHalfSize(halfSize);
+ const Vector3r& cos=scene->cell->getCos();
+ for(int i=0; i<3; i++){
+ //cerr<<"cos["<<i<<"]"<<cos[i]<<" ";
+ int i1=(i+1)%3,i2=(i+2)%3;
+ halfSize[i1]+=.5*refHalfSize[i1]*(1/cos[i]-1);
+ halfSize[i2]+=.5*refHalfSize[i2]*(1/cos[i]-1);
+ }
+ }
+ //cerr<<" || "<<halfSize<<endl;
+ aabb->min = scene->cell->unshearPt(se3.position)-halfSize;
+ aabb->max = scene->cell->unshearPt(se3.position)+halfSize;
+}
+
+void Bo1_Facet_Aabb::go( const shared_ptr<Shape>& cm
+ , shared_ptr<Bound>& bv
+ , const Se3r& se3
+ , const Body* b)
+{
+ if(!bv){ bv=shared_ptr<Bound>(new Aabb); }
+ Aabb* aabb=static_cast<Aabb*>(bv.get());
+ Facet* facet = static_cast<Facet*>(cm.get());
+ const Vector3r& O = se3.position;
+ Matrix3r facetAxisT=se3.orientation.toRotationMatrix();
+ const vector<Vector3r>& vertices=facet->vertices;
+ if(!scene->isPeriodic){
+ aabb->min=aabb->max = O + facetAxisT * vertices[0];
+ for (int i=1;i<3;++i)
+ {
+ Vector3r v = O + facetAxisT * vertices[i];
+ aabb->min = aabb->min.cwiseMin(v);
+ aabb->max = aabb->max.cwiseMax(v);
+ }
+ } else {
+ Real inf=std::numeric_limits<Real>::infinity();
+ aabb->min=Vector3r(inf,inf,inf); aabb->max=Vector3r(-inf,-inf,-inf);
+ for(int i=0; i<3; i++){
+ Vector3r v=scene->cell->unshearPt(O+facetAxisT*vertices[i]);
+ aabb->min=aabb->min.cwiseMin(v);
+ aabb->max=aabb->max.cwiseMax(v);
+ }
+ }
+}
+
+void Bo1_Box_Aabb::go( const shared_ptr<Shape>& cm,
+ shared_ptr<Bound>& bv,
+ const Se3r& se3,
+ const Body* b)
+{
+ Box* box = static_cast<Box*>(cm.get());
+ if(!bv){ bv=shared_ptr<Bound>(new Aabb); }
+ Aabb* aabb=static_cast<Aabb*>(bv.get());
+
+ if(scene->isPeriodic && scene->cell->hasShear()) throw logic_error(__FILE__ "Boxes not (yet?) supported in sheared cell.");
+
+ Matrix3r r=se3.orientation.toRotationMatrix();
+ Vector3r halfSize(Vector3r::Zero());
+ for( int i=0; i<3; ++i )
+ for( int j=0; j<3; ++j )
+ halfSize[i] += std::abs( r(i,j) * box->extents[j] );
+
+ aabb->min = se3.position-halfSize;
+ aabb->max = se3.position+halfSize;
+}
=== added file 'pkg/common/Bo1_Aabb.hpp'
--- pkg/common/Bo1_Aabb.hpp 1970-01-01 00:00:00 +0000
+++ pkg/common/Bo1_Aabb.hpp 2014-07-19 19:52:41 +0000
@@ -0,0 +1,66 @@
+/*************************************************************************
+* Copyright (C) 2004 by Olivier Galizzi *
+* olivier.galizzi@xxxxxxx *
+* *
+* This program is free software; it is licensed under the terms of the *
+* GNU General Public License v2 or later. See file LICENSE for details. *
+*************************************************************************/
+
+#pragma once
+
+#include <yade/pkg/common/Dispatching.hpp>
+#include <yade/pkg/common/Aabb.hpp>
+#include <yade/pkg/common/Sphere.hpp>
+#include <yade/pkg/common/Facet.hpp>
+#include <yade/pkg/common/Box.hpp>
+
+class Bo1_Sphere_Aabb : public BoundFunctor
+{
+ public :
+ void go(const shared_ptr<Shape>& cm, shared_ptr<Bound>& bv, const Se3r&, const Body*);
+ FUNCTOR1D(Sphere);
+ YADE_CLASS_BASE_DOC_ATTRS(Bo1_Sphere_Aabb,BoundFunctor,"Functor creating :yref:`Aabb` from :yref:`Sphere`.",
+ ((Real,aabbEnlargeFactor,((void)"deactivated",-1),,"Relative enlargement of the bounding box; deactivated if negative.\n\n.. note::\n\tThis attribute is used to create distant interaction, but is only meaningful with an :yref:`IGeomFunctor` which will not simply discard such interactions: :yref:`Ig2_Sphere_Sphere_ScGeom::interactionDetectionFactor` should have the same value as :yref:`aabbEnlargeFactor<Bo1_Sphere_Aabb::aabbEnlargeFactor>`."))
+ );
+};
+
+REGISTER_SERIALIZABLE(Bo1_Sphere_Aabb);
+
+
+/*************************************************************************
+* Copyright (C) 2008 by Sergei Dorofeenko *
+* sega@xxxxxxxxxxxxxxxx *
+* *
+* This program is free software; it is licensed under the terms of the *
+* GNU General Public License v2 or later. See file LICENSE for details. *
+*************************************************************************/
+
+class Bo1_Facet_Aabb : public BoundFunctor{
+ public:
+ void go(const shared_ptr<Shape>& cm, shared_ptr<Bound>& bv, const Se3r& se3, const Body*);
+ FUNCTOR1D(Facet);
+ YADE_CLASS_BASE_DOC(Bo1_Facet_Aabb,BoundFunctor,"Creates/updates an :yref:`Aabb` of a :yref:`Facet`.");
+};
+REGISTER_SERIALIZABLE(Bo1_Facet_Aabb);
+
+
+/*************************************************************************
+* Copyright (C) 2004 by Janek Kozicki *
+* cosurgi@xxxxxxxxxx *
+* *
+* This program is free software; it is licensed under the terms of the *
+* GNU General Public License v2 or later. See file LICENSE for details. *
+*************************************************************************/
+
+class Box;
+class Bo1_Box_Aabb : public BoundFunctor{
+ public:
+ void go(const shared_ptr<Shape>& cm, shared_ptr<Bound>& bv, const Se3r& se3, const Body*);
+ FUNCTOR1D(Box);
+ YADE_CLASS_BASE_DOC(Bo1_Box_Aabb,BoundFunctor,"Create/update an :yref:`Aabb` of a :yref:`Box`.");
+};
+
+REGISTER_SERIALIZABLE(Bo1_Box_Aabb);
+
+
+
=== removed file 'pkg/common/Bo1_Box_Aabb.cpp'
--- pkg/common/Bo1_Box_Aabb.cpp 2014-07-03 07:16:58 +0000
+++ pkg/common/Bo1_Box_Aabb.cpp 1970-01-01 00:00:00 +0000
@@ -1,35 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Olivier Galizzi *
-* olivier.galizzi@xxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#include "Bo1_Box_Aabb.hpp"
-#include<yade/pkg/common/Box.hpp>
-#include<yade/pkg/common/Aabb.hpp>
-
-
-void Bo1_Box_Aabb::go( const shared_ptr<Shape>& cm,
- shared_ptr<Bound>& bv,
- const Se3r& se3,
- const Body* b)
-{
- Box* box = static_cast<Box*>(cm.get());
- if(!bv){ bv=shared_ptr<Bound>(new Aabb); }
- Aabb* aabb=static_cast<Aabb*>(bv.get());
-
- if(scene->isPeriodic && scene->cell->hasShear()) throw logic_error(__FILE__ "Boxes not (yet?) supported in sheared cell.");
-
- Matrix3r r=se3.orientation.toRotationMatrix();
- Vector3r halfSize(Vector3r::Zero());
- for( int i=0; i<3; ++i )
- for( int j=0; j<3; ++j )
- halfSize[i] += std::abs( r(i,j) * box->extents[j] );
-
- aabb->min = se3.position-halfSize;
- aabb->max = se3.position+halfSize;
-}
-
-YADE_PLUGIN((Bo1_Box_Aabb));
=== removed file 'pkg/common/Bo1_Box_Aabb.hpp'
--- pkg/common/Bo1_Box_Aabb.hpp 2011-01-09 16:34:50 +0000
+++ pkg/common/Bo1_Box_Aabb.hpp 1970-01-01 00:00:00 +0000
@@ -1,24 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Janek Kozicki *
-* cosurgi@xxxxxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#pragma once
-
-
-#include<yade/pkg/common/Dispatching.hpp>
-
-class Box;
-class Bo1_Box_Aabb : public BoundFunctor{
- public:
- void go(const shared_ptr<Shape>& cm, shared_ptr<Bound>& bv, const Se3r& se3, const Body*);
- FUNCTOR1D(Box);
- YADE_CLASS_BASE_DOC(Bo1_Box_Aabb,BoundFunctor,"Create/update an :yref:`Aabb` of a :yref:`Box`.");
-};
-
-REGISTER_SERIALIZABLE(Bo1_Box_Aabb);
-
-
=== removed file 'pkg/common/Bo1_Facet_Aabb.cpp'
--- pkg/common/Bo1_Facet_Aabb.cpp 2013-03-07 18:28:01 +0000
+++ pkg/common/Bo1_Facet_Aabb.cpp 1970-01-01 00:00:00 +0000
@@ -1,43 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2008 by Sergei Dorofeenko *
-* sega@xxxxxxxxxxxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#include <yade/pkg/common/Facet.hpp>
-#include <yade/pkg/common/Bo1_Facet_Aabb.hpp>
-#include <yade/pkg/common/Aabb.hpp>
-
-void Bo1_Facet_Aabb::go( const shared_ptr<Shape>& cm
- , shared_ptr<Bound>& bv
- , const Se3r& se3
- , const Body* b)
-{
- if(!bv){ bv=shared_ptr<Bound>(new Aabb); }
- Aabb* aabb=static_cast<Aabb*>(bv.get());
- Facet* facet = static_cast<Facet*>(cm.get());
- const Vector3r& O = se3.position;
- Matrix3r facetAxisT=se3.orientation.toRotationMatrix();
- const vector<Vector3r>& vertices=facet->vertices;
- if(!scene->isPeriodic){
- aabb->min=aabb->max = O + facetAxisT * vertices[0];
- for (int i=1;i<3;++i)
- {
- Vector3r v = O + facetAxisT * vertices[i];
- aabb->min = aabb->min.cwiseMin(v);
- aabb->max = aabb->max.cwiseMax(v);
- }
- } else {
- Real inf=std::numeric_limits<Real>::infinity();
- aabb->min=Vector3r(inf,inf,inf); aabb->max=Vector3r(-inf,-inf,-inf);
- for(int i=0; i<3; i++){
- Vector3r v=scene->cell->unshearPt(O+facetAxisT*vertices[i]);
- aabb->min=aabb->min.cwiseMin(v);
- aabb->max=aabb->max.cwiseMax(v);
- }
- }
-}
-
-YADE_PLUGIN((Bo1_Facet_Aabb));
=== removed file 'pkg/common/Bo1_Facet_Aabb.hpp'
--- pkg/common/Bo1_Facet_Aabb.hpp 2010-11-12 08:03:16 +0000
+++ pkg/common/Bo1_Facet_Aabb.hpp 1970-01-01 00:00:00 +0000
@@ -1,21 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2008 by Sergei Dorofeenko *
-* sega@xxxxxxxxxxxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#pragma once
-
-#include <yade/pkg/common/Dispatching.hpp>
-
-class Bo1_Facet_Aabb : public BoundFunctor{
- public:
- void go(const shared_ptr<Shape>& cm, shared_ptr<Bound>& bv, const Se3r& se3, const Body*);
- FUNCTOR1D(Facet);
- YADE_CLASS_BASE_DOC(Bo1_Facet_Aabb,BoundFunctor,"Creates/updates an :yref:`Aabb` of a :yref:`Facet`.");
-};
-REGISTER_SERIALIZABLE(Bo1_Facet_Aabb);
-
-
=== removed file 'pkg/common/Bo1_Sphere_Aabb.cpp'
--- pkg/common/Bo1_Sphere_Aabb.cpp 2010-11-07 11:46:20 +0000
+++ pkg/common/Bo1_Sphere_Aabb.cpp 1970-01-01 00:00:00 +0000
@@ -1,38 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Olivier Galizzi *
-* olivier.galizzi@xxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#include "Bo1_Sphere_Aabb.hpp"
-#include<yade/pkg/common/Sphere.hpp>
-#include<yade/pkg/common/Aabb.hpp>
-
-void Bo1_Sphere_Aabb::go(const shared_ptr<Shape>& cm, shared_ptr<Bound>& bv, const Se3r& se3, const Body* b){
- Sphere* sphere = static_cast<Sphere*>(cm.get());
- if(!bv){ bv=shared_ptr<Bound>(new Aabb); }
- Aabb* aabb=static_cast<Aabb*>(bv.get());
- Vector3r halfSize = (aabbEnlargeFactor>0?aabbEnlargeFactor:1.)*Vector3r(sphere->radius,sphere->radius,sphere->radius);
- if(!scene->isPeriodic){
- aabb->min=se3.position-halfSize; aabb->max=se3.position+halfSize;
- return;
- }
- // adjust box size along axes so that sphere doesn't stick out of the box even if sheared (i.e. parallelepiped)
- if(scene->cell->hasShear()) {
- Vector3r refHalfSize(halfSize);
- const Vector3r& cos=scene->cell->getCos();
- for(int i=0; i<3; i++){
- //cerr<<"cos["<<i<<"]"<<cos[i]<<" ";
- int i1=(i+1)%3,i2=(i+2)%3;
- halfSize[i1]+=.5*refHalfSize[i1]*(1/cos[i]-1);
- halfSize[i2]+=.5*refHalfSize[i2]*(1/cos[i]-1);
- }
- }
- //cerr<<" || "<<halfSize<<endl;
- aabb->min = scene->cell->unshearPt(se3.position)-halfSize;
- aabb->max = scene->cell->unshearPt(se3.position)+halfSize;
-}
-
-YADE_PLUGIN((Bo1_Sphere_Aabb));
=== removed file 'pkg/common/Bo1_Sphere_Aabb.hpp'
--- pkg/common/Bo1_Sphere_Aabb.hpp 2012-06-28 19:34:29 +0000
+++ pkg/common/Bo1_Sphere_Aabb.hpp 1970-01-01 00:00:00 +0000
@@ -1,26 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Olivier Galizzi *
-* olivier.galizzi@xxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#pragma once
-
-#include<yade/pkg/common/Dispatching.hpp>
-#include<yade/pkg/common/Sphere.hpp>
-
-class Bo1_Sphere_Aabb : public BoundFunctor
-{
- public :
- void go(const shared_ptr<Shape>& cm, shared_ptr<Bound>& bv, const Se3r&, const Body*);
- FUNCTOR1D(Sphere);
- YADE_CLASS_BASE_DOC_ATTRS(Bo1_Sphere_Aabb,BoundFunctor,"Functor creating :yref:`Aabb` from :yref:`Sphere`.",
- ((Real,aabbEnlargeFactor,((void)"deactivated",-1),,"Relative enlargement of the bounding box; deactivated if negative.\n\n.. note::\n\tThis attribute is used to create distant interaction, but is only meaningful with an :yref:`IGeomFunctor` which will not simply discard such interactions: :yref:`Ig2_Sphere_Sphere_ScGeom::interactionDetectionFactor` should have the same value as :yref:`aabbEnlargeFactor<Bo1_Sphere_Aabb::aabbEnlargeFactor>`."))
- );
-};
-
-REGISTER_SERIALIZABLE(Bo1_Sphere_Aabb);
-
-
=== removed file 'pkg/common/BoundaryController.cpp'
--- pkg/common/BoundaryController.cpp 2014-07-03 17:20:40 +0000
+++ pkg/common/BoundaryController.cpp 1970-01-01 00:00:00 +0000
@@ -1,3 +0,0 @@
-#include<yade/pkg/common/BoundaryController.hpp>
-void BoundaryController::action(){ throw std::runtime_error("BoundaryController must not be used in simulations directly (BoundaryController::action called)."); }
-YADE_PLUGIN((BoundaryController));
=== modified file 'pkg/common/BoundaryController.hpp'
--- pkg/common/BoundaryController.hpp 2010-10-13 16:23:08 +0000
+++ pkg/common/BoundaryController.hpp 2014-07-21 09:09:14 +0000
@@ -1,7 +1,9 @@
#pragma once
#include<yade/core/GlobalEngine.hpp>
class BoundaryController: public GlobalEngine{
- virtual void action();
+ virtual void action() {
+ { throw std::runtime_error("BoundaryController must not be used in simulations directly (BoundaryController::action called)."); }
+ }
YADE_CLASS_BASE_DOC(BoundaryController,GlobalEngine,"Base for engines controlling boundary conditions of simulations. Not to be used directly.");
};
REGISTER_SERIALIZABLE(BoundaryController);
=== removed file 'pkg/common/Box.cpp'
--- pkg/common/Box.cpp 2010-10-13 16:23:08 +0000
+++ pkg/common/Box.cpp 1970-01-01 00:00:00 +0000
@@ -1,15 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Olivier Galizzi *
-* olivier.galizzi@xxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#include "Box.hpp"
-
-Box::~Box ()
-{
-}
-
-YADE_PLUGIN((Box));
=== modified file 'pkg/common/Box.hpp'
--- pkg/common/Box.hpp 2010-10-13 16:23:08 +0000
+++ pkg/common/Box.hpp 2014-07-21 09:09:14 +0000
@@ -14,7 +14,7 @@
class Box: public Shape{
public:
Box(const Vector3r& _extents): extents(_extents){}
- virtual ~Box ();
+ virtual ~Box () {};
YADE_CLASS_BASE_DOC_ATTRS_CTOR(Box,Shape,"Box (cuboid) particle geometry. (Avoid using in new code, prefer :yref:`Facet` instead.",
((Vector3r,extents,,,"Half-size of the cuboid")),
/* ctor */ createIndex();
=== removed file 'pkg/common/Callbacks.cpp'
--- pkg/common/Callbacks.cpp 2011-01-09 16:34:50 +0000
+++ pkg/common/Callbacks.cpp 1970-01-01 00:00:00 +0000
@@ -1,13 +0,0 @@
-#include<yade/pkg/common/Callbacks.hpp>
-
-#ifdef YADE_BODY_CALLBACK
- BodyCallback::~BodyCallback(){};
-#endif
-
-IntrCallback::~IntrCallback(){};
-
-YADE_PLUGIN((IntrCallback)
- #ifdef YADE_BODY_CALLBACK
- (BodyCallback)
- #endif
-);
=== modified file 'pkg/common/Callbacks.hpp'
--- pkg/common/Callbacks.hpp 2014-07-03 17:20:40 +0000
+++ pkg/common/Callbacks.hpp 2014-07-21 09:09:14 +0000
@@ -9,7 +9,7 @@
class IntrCallback: public Serializable{
public:
- virtual ~IntrCallback(); // vtable
+ virtual ~IntrCallback() {}; // vtable
typedef void(*FuncPtr)(IntrCallback*,Interaction*);
// should be set at every step by InteractionLoop
Scene* scene;
@@ -26,7 +26,7 @@
#ifdef YADE_BODY_CALLBACKS
class BodyCallback: public Serializable{
public:
- virtual ~BodyCallback(); // vtable
+ virtual ~BodyCallback() {}; // vtable
typedef void(*FuncPtr)(BodyCallback*,Body*);
// set at every step, before stepInit() is called
Scene* scene;
=== removed file 'pkg/common/CylScGeom6D.cpp'
--- pkg/common/CylScGeom6D.cpp 2012-04-10 19:19:27 +0000
+++ pkg/common/CylScGeom6D.cpp 1970-01-01 00:00:00 +0000
@@ -1,24 +0,0 @@
-#include "CylScGeom6D.hpp"
-#include<yade/core/Omega.hpp>
-#include<yade/core/Scene.hpp>
-
-CylScGeom::~CylScGeom() {}
-CylScGeom6D::~CylScGeom6D() {}
-
-
-
-void CylScGeom6D::precomputeRotations(const State& rbp1, const State& rbp2, bool isNew, bool creep) {
- initRotations(rbp1,rbp2);
-}
-
-void CylScGeom6D::initRotations(const State& state1, const State& state2)
-{
- initialOrientation1 = state1.ori;
- initialOrientation2 = state2.ori;
- twist=0;
- bending=Vector3r::Zero();
- twistCreep=Quaternionr(1.0,0.0,0.0,0.0);
-}
-
-YADE_PLUGIN((CylScGeom6D)(CylScGeom));
-
=== modified file 'pkg/common/CylScGeom6D.hpp'
--- pkg/common/CylScGeom6D.hpp 2014-05-06 15:32:52 +0000
+++ pkg/common/CylScGeom6D.hpp 2014-07-21 09:09:14 +0000
@@ -10,7 +10,7 @@
State fictiousState;
// shared_ptr<Interaction> duplicate;
- virtual ~CylScGeom ();
+ virtual ~CylScGeom () {};
YADE_CLASS_BASE_DOC_ATTRS_CTOR(CylScGeom,ScGeom,"Geometry of a cylinder-sphere contact.",
((bool,onNode,false,,"contact on node?"))
((int,isDuplicate,0,,"this flag is turned true (1) automatically if the contact is shared between two chained cylinders. A duplicated interaction will be skipped once by the constitutive law, so that only one contact at a time is effective. If isDuplicate=2, it means one of the two duplicates has no longer geometric interaction, and should be erased by the constitutive laws."))
@@ -28,9 +28,17 @@
class CylScGeom6D: public ScGeom6D {
public:
- virtual ~CylScGeom6D();
- void precomputeRotations(const State& rbp1, const State& rbp2, bool isNew, bool creep=false);
- void initRotations(const State& rbp1, const State& rbp2);
+ virtual ~CylScGeom6D() {};
+ void precomputeRotations(const State& rbp1, const State& rbp2, bool isNew, bool creep=false) {
+ initRotations(rbp1,rbp2);
+ }
+ void initRotations(const State& rbp1, const State& rbp2){
+ initialOrientation1 = rbp1.ori;
+ initialOrientation2 = rbp2.ori;
+ twist=0;
+ bending=Vector3r::Zero();
+ twistCreep=Quaternionr(1.0,0.0,0.0,0.0);
+ }
State fictiousState;
YADE_CLASS_BASE_DOC_ATTRS_INIT_CTOR_PY(CylScGeom6D,ScGeom6D,"Class representing :yref:`geometry<IGeom>` of two :yref:`bodies<Body>` in contact. The contact has 6 DOFs (normal, 2×shear, twist, 2xbending) and uses :yref:`ScGeom` incremental algorithm for updating shear.",
((bool,onNode,false,,"contact on node?"))
=== modified file 'pkg/common/Cylinder.cpp'
--- pkg/common/Cylinder.cpp 2014-07-03 07:16:58 +0000
+++ pkg/common/Cylinder.cpp 2014-07-18 18:18:50 +0000
@@ -658,7 +658,7 @@
}
}
-void Law2_CylScGeom_FrictPhys_CundallStrack::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
+bool Law2_CylScGeom_FrictPhys_CundallStrack::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
int id1 = contact->getId1(), id2 = contact->getId2();
CylScGeom* geom= static_cast<CylScGeom*>(ig.get());
@@ -667,13 +667,11 @@
if (neverErase) {
phys->shearForce = Vector3r::Zero();
phys->normalForce = Vector3r::Zero();}
- else scene->interactions->requestErase(contact);
- return;}
+ else return false;}
if (geom->isDuplicate) {
if (id2!=geom->trueInt) {
//cerr<<"skip duplicate "<<id1<<" "<<id2<<endl;
- if (geom->isDuplicate==2) {/*cerr<<"erase duplicate "<<id1<<" "<<id2<<endl;*/scene->interactions->requestErase(contact);}
- return;}
+ if (geom->isDuplicate==2) return false;}
}
Real& un=geom->penetrationDepth;
phys->normalForce=phys->kn*std::max(un,(Real) 0)*geom->normal;
@@ -721,10 +719,11 @@
scene->forces.addTorque(id1,(geom->radius1-0.5*geom->penetrationDepth)* geom->normal.cross(force));
scene->forces.addTorque(id2,(geom->radius2-0.5*geom->penetrationDepth)* geom->normal.cross(force));
}
+ return true;
}
-void Law2_CylScGeom6D_CohFrictPhys_CohesionMoment::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact) {
+bool Law2_CylScGeom6D_CohFrictPhys_CohesionMoment::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact) {
int id1 = contact->getId1(), id2 = contact->getId2();
@@ -738,19 +737,17 @@
if (geom->isDuplicate) {
if (id2!=geom->trueInt) {
//cerr<<"skip duplicate "<<id1<<" "<<id2<<endl;
- if (geom->isDuplicate==2) {/*cerr<<"erase duplicate coh "<<id1<<" "<<id2<<endl;*/scene->interactions->requestErase(contact);}
- return;}
+ if (geom->isDuplicate==2) return false;}
}
if (currentContactPhysics->fragile && (-Fn)> currentContactPhysics->normalAdhesion) {
// BREAK due to tension
- scene->interactions->requestErase(contact);
+ return false;
} else {
if ((-Fn)> currentContactPhysics->normalAdhesion) {//normal plasticity
Fn=-currentContactPhysics->normalAdhesion;
currentContactPhysics->unp = un+currentContactPhysics->normalAdhesion/currentContactPhysics->kn;
- if (currentContactPhysics->unpMax && currentContactPhysics->unp<currentContactPhysics->unpMax)
- scene->interactions->requestErase(contact);
+ if (currentContactPhysics->unpMax && currentContactPhysics->unp<currentContactPhysics->unpMax) return false;
}
currentContactPhysics->normalForce = Fn*geom->normal;
Vector3r& shearForce = geom->rotate(currentContactPhysics->shearForce);
@@ -794,11 +791,12 @@
}
//applyForceAtContactPoint(-currentContactPhysics->normalForce-shearForce, geom->contactPoint, id1, de1->se3.position, id2, de2->se3.position);
}
+ return true;
}
-void Law2_ChCylGeom6D_CohFrictPhys_CohesionMoment::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
+bool Law2_ChCylGeom6D_CohFrictPhys_CohesionMoment::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
int id1 = contact->getId1(), id2 = contact->getId2();
ChCylGeom6D* geom= YADE_CAST<ChCylGeom6D*>(ig.get());
@@ -816,24 +814,15 @@
if (contact->isFresh(scene)) shearForce = Vector3r::Zero(); //contact nouveau => force tengentielle = 0,0,0
Real un = geom->penetrationDepth; //un : interpenetration
Real Fn = currentContactPhysics->kn*(un-currentContactPhysics->unp); //Fn : force normale
- /*if (geom->isDuplicate) {
- if (id2!=geom->trueInt) {
- //cerr<<"skip duplicate "<<id1<<" "<<id2<<endl;
- if (geom->isDuplicate==2) {cerr<<"erase duplicate coh "<<id1<<" "<<id2<<endl;scene->interactions->requestErase(contact);}
- return;}
- }
-*/
- if (currentContactPhysics->fragile && (-Fn)> currentContactPhysics->normalAdhesion) {
- // BREAK due to tension
- scene->interactions->requestErase(contact);
- } else {
+ if (currentContactPhysics->fragile && (-Fn)> currentContactPhysics->normalAdhesion) return false; // BREAK due to tension
+ else {
if ((-Fn)> currentContactPhysics->normalAdhesion) {//normal plasticity
Fn=-currentContactPhysics->normalAdhesion;
currentContactPhysics->unp = un+currentContactPhysics->normalAdhesion/currentContactPhysics->kn;
if (currentContactPhysics->unpMax && currentContactPhysics->unp<currentContactPhysics->unpMax)
- scene->interactions->requestErase(contact);
- }
+ return false;
+ }
currentContactPhysics->normalForce = Fn*geom->normal;
@@ -885,6 +874,7 @@
scene->forces.addTorque(id2,(geom->radius2-0.5*geom->penetrationDepth)* geom->normal.cross(force));
}
}
+ return true;
}
=== modified file 'pkg/common/Cylinder.hpp'
--- pkg/common/Cylinder.hpp 2014-07-14 20:08:33 +0000
+++ pkg/common/Cylinder.hpp 2014-07-18 18:18:50 +0000
@@ -218,7 +218,7 @@
class Law2_CylScGeom_FrictPhys_CundallStrack: public LawFunctor{
public:
//OpenMPAccumulator<Real> plasticDissipation;
- virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
//Real elasticEnergy ();
//Real getPlasticDissipation();
//void initPlasticDissipation(Real initVal=0);
@@ -239,7 +239,7 @@
class Law2_CylScGeom6D_CohFrictPhys_CohesionMoment: public LawFunctor {
public:
//OpenMPAccumulator<Real> plasticDissipation;
- virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
//Real elasticEnergy ();
//Real getPlasticDissipation();
//void initPlasticDissipation(Real initVal=0);
@@ -265,7 +265,7 @@
class Law2_ChCylGeom6D_CohFrictPhys_CohesionMoment: public LawFunctor {
public:
//OpenMPAccumulator<Real> plasticDissipation;
- virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(Law2_ChCylGeom6D_CohFrictPhys_CohesionMoment,LawFunctor,"Law for linear compression, and Mohr-Coulomb plasticity surface without cohesion.\nThis law implements the classical linear elastic-plastic law from [CundallStrack1979]_ (see also [Pfc3dManual30]_). The normal force is (with the convention of positive tensile forces) $F_n=\\min(k_n u_n, 0)$. The shear force is $F_s=k_s u_s$, the plasticity condition defines the maximum value of the shear force : $F_s^{\\max}=F_n\\tan(\\phi)$, with $\\phi$ the friction angle.\n\n.. note::\n This law is well tested in the context of triaxial simulation, and has been used for a number of published results (see e.g. [Scholtes2009b]_ and other papers from the same authors). It is generalised by :yref:`Law2_ScGeom6D_CohFrictPhys_CohesionMoment`, which adds cohesion and moments at contact.",
((bool,neverErase,false,,"Keep interactions even if particles go away from each other (only in case another constitutive law is in the scene, e.g. :yref:`Law2_ScGeom_CapillaryPhys_Capillarity`)"))
((bool,traceEnergy,false,Attr::hidden,"Define the total energy dissipated in plastic slips at all contacts."))
=== modified file 'pkg/common/Dispatching.hpp'
--- pkg/common/Dispatching.hpp 2014-02-22 04:56:33 +0000
+++ pkg/common/Dispatching.hpp 2014-07-18 18:18:50 +0000
@@ -52,7 +52,7 @@
class LawFunctor: public Functor2D<
/*dispatch types*/ IGeom,IPhys,
- /*return type*/ void,
+ /*return type*/ bool,
/*argument types*/ TYPELIST_3(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*)
>{
public: virtual ~LawFunctor();
=== removed file 'pkg/common/ElastMat.cpp'
--- pkg/common/ElastMat.cpp 2010-11-07 11:46:20 +0000
+++ pkg/common/ElastMat.cpp 1970-01-01 00:00:00 +0000
@@ -1,5 +0,0 @@
-// 2009 © Václav Šmilauer <eudoxos@xxxxxxxx>
-#include<yade/pkg/common/ElastMat.hpp>
-YADE_PLUGIN((ElastMat)(FrictMat));
-ElastMat::~ElastMat(){}
-FrictMat::~FrictMat(){}
=== modified file 'pkg/common/ElastMat.hpp'
--- pkg/common/ElastMat.hpp 2014-07-03 17:20:40 +0000
+++ pkg/common/ElastMat.hpp 2014-07-19 19:52:41 +0000
@@ -5,7 +5,7 @@
/*! Elastic material */
class ElastMat: public Material{
public:
- virtual ~ElastMat();
+ virtual ~ElastMat() {};
YADE_CLASS_BASE_DOC_ATTRS_CTOR(ElastMat,Material,"Purely elastic material. The material parameters may have different meanings depending on the :yref:`IPhysFunctor` used : true Young and Poisson in :yref:`Ip2_FrictMat_FrictMat_MindlinPhys`, or contact stiffnesses in :yref:`Ip2_FrictMat_FrictMat_FrictPhys`.",
((Real,young,1e9,,"elastic modulus [Pa]. It has different meanings depending on the Ip functor."))
((Real,poisson,.25,,"Poisson's ratio or the ratio between shear and normal stiffness [-]. It has different meanings depending on the Ip functor. ")),
@@ -18,7 +18,7 @@
/*! Granular material */
class FrictMat: public ElastMat{
public:
- virtual ~FrictMat();
+ virtual ~FrictMat() {};
YADE_CLASS_BASE_DOC_ATTRS_CTOR(FrictMat,ElastMat,"Elastic material with contact friction. See also :yref:`ElastMat`.",
((Real,frictionAngle,.5,,"Contact friction angle (in radians). Hint : use 'radians(degreesValue)' in python scripts.")),
createIndex();
=== removed file 'pkg/common/FieldApplier.cpp'
--- pkg/common/FieldApplier.cpp 2010-11-07 11:46:20 +0000
+++ pkg/common/FieldApplier.cpp 1970-01-01 00:00:00 +0000
@@ -1,5 +0,0 @@
-#include<yade/pkg/common/FieldApplier.hpp>
-#include<stdexcept>
-void FieldApplier::action(){ throw std::runtime_error("FieldApplier must not be used in simulations directly (FieldApplier::action called)."); }
-YADE_PLUGIN((FieldApplier));
-
=== modified file 'pkg/common/FieldApplier.hpp'
--- pkg/common/FieldApplier.hpp 2010-11-13 21:11:39 +0000
+++ pkg/common/FieldApplier.hpp 2014-07-21 09:09:14 +0000
@@ -1,7 +1,11 @@
#pragma once
-#include<yade/core/GlobalEngine.hpp>
+#include <yade/core/GlobalEngine.hpp>
+#include <stdexcept>
+
class FieldApplier: public GlobalEngine{
- virtual void action();
+ virtual void action() {
+ throw std::runtime_error("FieldApplier must not be used in simulations directly (FieldApplier::action called).");
+ }
YADE_CLASS_BASE_DOC_ATTRS(FieldApplier,GlobalEngine,"Base for engines applying force files on particles. Not to be used directly.",
((int,fieldWorkIx,-1,(Attr::hidden|Attr::noSave),"Index for the work done by this field, if tracking energies."))
);
=== removed file 'pkg/common/ForceResetter.cpp'
--- pkg/common/ForceResetter.cpp 2013-03-06 17:30:45 +0000
+++ pkg/common/ForceResetter.cpp 1970-01-01 00:00:00 +0000
@@ -1,10 +0,0 @@
-#include<yade/pkg/common/ForceResetter.hpp>
-#include<yade/core/Scene.hpp>
-
-YADE_PLUGIN((ForceResetter));
-
-void ForceResetter::action(){
- scene->forces.reset(scene->iter);
- if(scene->trackEnergy) scene->energy->resetResettables();
-}
-
=== modified file 'pkg/common/ForceResetter.hpp'
--- pkg/common/ForceResetter.hpp 2010-10-29 10:12:44 +0000
+++ pkg/common/ForceResetter.hpp 2014-07-21 09:09:14 +0000
@@ -1,11 +1,15 @@
#pragma once
#include<yade/core/GlobalEngine.hpp>
+#include<yade/core/Scene.hpp>
class Scene;
class ForceResetter: public GlobalEngine{
public:
- virtual void action();
+ virtual void action() {
+ scene->forces.reset(scene->iter);
+ if(scene->trackEnergy) scene->energy->resetResettables();
+ }
YADE_CLASS_BASE_DOC(ForceResetter,GlobalEngine,"Reset all forces stored in Scene::forces (``O.forces`` in python). Typically, this is the first engine to be run at every step. In addition, reset those energies that should be reset, if energy tracing is enabled.");
};
REGISTER_SERIALIZABLE(ForceResetter);
=== removed file 'pkg/common/GLDrawFunctors.cpp'
--- pkg/common/GLDrawFunctors.cpp 2010-11-19 12:30:08 +0000
+++ pkg/common/GLDrawFunctors.cpp 1970-01-01 00:00:00 +0000
@@ -1,7 +0,0 @@
-#include<yade/pkg/common/GLDrawFunctors.hpp>
-#ifdef YADE_OPENGL
- YADE_PLUGIN(
- (GlBoundFunctor)(GlShapeFunctor)(GlIGeomFunctor)(GlIPhysFunctor)(GlStateFunctor)
- (GlBoundDispatcher)(GlShapeDispatcher)(GlIGeomDispatcher)(GlIPhysDispatcher)(GlStateDispatcher)
- );
-#endif
=== removed file 'pkg/common/Gl1_Aabb.cpp'
--- pkg/common/Gl1_Aabb.cpp 2010-11-19 12:30:08 +0000
+++ pkg/common/Gl1_Aabb.cpp 1970-01-01 00:00:00 +0000
@@ -1,33 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Olivier Galizzi *
-* olivier.galizzi@xxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#ifdef YADE_OPENGL
-
-#include"Gl1_Aabb.hpp"
-#include<yade/lib/opengl/OpenGLWrapper.hpp>
-#include<yade/pkg/common/Aabb.hpp>
-#include<yade/core/Scene.hpp>
-
-void Gl1_Aabb::go(const shared_ptr<Bound>& bv, Scene* scene){
- Aabb* aabb = static_cast<Aabb*>(bv.get());
- glColor3v(bv->color);
- // glDisable(GL_LIGHTING);
- if(!scene->isPeriodic){
- glTranslatev(Vector3r(.5*(aabb->min+aabb->max)));
- glScalev(Vector3r(aabb->max-aabb->min));
- } else {
- glTranslatev(Vector3r(scene->cell->shearPt(scene->cell->wrapPt(.5*(aabb->min+aabb->max)))));
- glMultMatrixd(scene->cell->getGlShearTrsfMatrix());
- glScalev(Vector3r(aabb->max-aabb->min));
- }
- glutWireCube(1);
-}
-
-YADE_PLUGIN((Gl1_Aabb));
-
-#endif /* YADE_OPENGL */
=== removed file 'pkg/common/Gl1_Aabb.hpp'
--- pkg/common/Gl1_Aabb.hpp 2011-01-10 08:08:23 +0000
+++ pkg/common/Gl1_Aabb.hpp 1970-01-01 00:00:00 +0000
@@ -1,22 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Olivier Galizzi *
-* olivier.galizzi@xxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#pragma once
-
-#include<yade/pkg/common/GLDrawFunctors.hpp>
-#include<yade/pkg/common/Aabb.hpp>
-
-class Gl1_Aabb: public GlBoundFunctor{
- public:
- virtual void go(const shared_ptr<Bound>&, Scene*);
- RENDERS(Aabb);
- YADE_CLASS_BASE_DOC(Gl1_Aabb,GlBoundFunctor,"Render Axis-aligned bounding box (:yref:`Aabb`).");
-};
-REGISTER_SERIALIZABLE(Gl1_Aabb);
-
-
=== removed file 'pkg/common/Gl1_Box.cpp'
--- pkg/common/Gl1_Box.cpp 2010-11-19 12:30:08 +0000
+++ pkg/common/Gl1_Box.cpp 1970-01-01 00:00:00 +0000
@@ -1,29 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Olivier Galizzi *
-* olivier.galizzi@xxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#ifdef YADE_OPENGL
-
-#include "Gl1_Box.hpp"
-#include<yade/pkg/common/Box.hpp>
-#include<yade/lib/opengl/OpenGLWrapper.hpp>
-
-void Gl1_Box::go(const shared_ptr<Shape>& cg, const shared_ptr<State>&,bool wire,const GLViewInfo&)
-{
- glColor3v(cg->color);
-
- Vector3r &extents = (static_cast<Box*>(cg.get()))->extents;
-
- glScalef(2*extents[0],2*extents[1],2*extents[2]);
-
- if (wire) glutWireCube(1);
- else glutSolidCube(1);
-}
-
-YADE_PLUGIN((Gl1_Box));
-
-#endif
=== removed file 'pkg/common/Gl1_Box.hpp'
--- pkg/common/Gl1_Box.hpp 2011-01-10 08:08:23 +0000
+++ pkg/common/Gl1_Box.hpp 1970-01-01 00:00:00 +0000
@@ -1,23 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2004 by Olivier Galizzi *
-* olivier.galizzi@xxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#pragma once
-
-#include<yade/pkg/common/GLDrawFunctors.hpp>
-#include<yade/pkg/common/Box.hpp>
-
-class Gl1_Box : public GlShapeFunctor{
- public :
- virtual void go(const shared_ptr<Shape>&, const shared_ptr<State>&,bool,const GLViewInfo&);
- RENDERS(Box);
- YADE_CLASS_BASE_DOC(Gl1_Box,GlShapeFunctor,"Renders :yref:`Box` object");
-};
-
-REGISTER_SERIALIZABLE(Gl1_Box);
-
-
=== removed file 'pkg/common/Gl1_Facet.cpp'
--- pkg/common/Gl1_Facet.cpp 2011-02-14 08:05:09 +0000
+++ pkg/common/Gl1_Facet.cpp 1970-01-01 00:00:00 +0000
@@ -1,61 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2008 by Sergei Dorofeenko *
-* sega@xxxxxxxxxxxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#ifdef YADE_OPENGL
-
-#include<yade/pkg/common/Gl1_Facet.hpp>
-#include<yade/pkg/common/Facet.hpp>
-#include<yade/lib/opengl/OpenGLWrapper.hpp>
-
-bool Gl1_Facet::normals=false;
-
-void Gl1_Facet::go(const shared_ptr<Shape>& cm, const shared_ptr<State>& ,bool wire,const GLViewInfo&)
-{
- Facet* facet = static_cast<Facet*>(cm.get());
- const vector<Vector3r>& vertices = facet->vertices;
- const Vector3r* ne = facet->ne;
- const Real& icr = facet->icr;
-
- if(cm->wire || wire){
- // facet
- glBegin(GL_LINE_LOOP);
- glColor3v(normals ? Vector3r(1,0,0): cm->color);
- glVertex3v(vertices[0]);
- glVertex3v(vertices[1]);
- glVertex3v(vertices[2]);
- glEnd();
- if(!normals) return;
- // facet's normal
- glBegin(GL_LINES);
- glColor3(0.0,0.0,1.0);
- glVertex3(0.0,0.0,0.0);
- glVertex3v(facet->normal);
- glEnd();
- // normal of edges
- glColor3(0.0,0.0,1.0);
- glBegin(GL_LINES);
- glVertex3(0.0,0.0,0.0); glVertex3v(Vector3r(icr*ne[0]));
- glVertex3(0.0,0.0,0.0); glVertex3v(Vector3r(icr*ne[1]));
- glVertex3(0.0,0.0,0.0); glVertex3v(Vector3r(icr*ne[2]));
- glEnd();
- } else {
- glDisable(GL_CULL_FACE);
- Vector3r normal=(facet->vertices[1]-facet->vertices[0]).cross(facet->vertices[2]-facet->vertices[1]); normal.normalize();
- glColor3v(cm->color);
- glBegin(GL_TRIANGLES);
- glNormal3v(normal); // this makes every triangle different WRT the light direction; important!
- glVertex3v(facet->vertices[0]);
- glVertex3v(facet->vertices[1]);
- glVertex3v(facet->vertices[2]);
- glEnd();
- }
-}
-
-YADE_PLUGIN((Gl1_Facet));
-
-#endif /* YADE_OPENGL */
=== removed file 'pkg/common/Gl1_Facet.hpp'
--- pkg/common/Gl1_Facet.hpp 2011-01-09 16:34:50 +0000
+++ pkg/common/Gl1_Facet.hpp 1970-01-01 00:00:00 +0000
@@ -1,26 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2008 by Sergei Dorofeenko *
-* sega@xxxxxxxxxxxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#pragma once
-
-#include<yade/pkg/common/GLDrawFunctors.hpp>
-#include<yade/pkg/common/Facet.hpp>
-
-class Gl1_Facet : public GlShapeFunctor
-{
- public:
- virtual void go(const shared_ptr<Shape>&, const shared_ptr<State>&,bool,const GLViewInfo&);
- RENDERS(Facet);
- YADE_CLASS_BASE_DOC_STATICATTRS(Gl1_Facet,GlShapeFunctor,"Renders :yref:`Facet` object",
- ((bool,normals,false,,"In wire mode, render normals of facets and edges; facet's :yref:`colors<Shape::color>` are disregarded in that case."))
- );
-};
-
-REGISTER_SERIALIZABLE(Gl1_Facet);
-
-
=== added file 'pkg/common/Gl1_Primitives.cpp'
--- pkg/common/Gl1_Primitives.cpp 1970-01-01 00:00:00 +0000
+++ pkg/common/Gl1_Primitives.cpp 2014-07-19 21:36:19 +0000
@@ -0,0 +1,225 @@
+/*************************************************************************
+* Copyright (C) 2004 by Olivier Galizzi *
+* olivier.galizzi@xxxxxxx *
+* *
+* Copyright (C) 2008 by Sergei Dorofeenko *
+* sega@xxxxxxxxxxxxxxxx *
+* *
+* © 2008 Václav Šmilauer <eudoxos@xxxxxxxx> *
+* *
+* © 2008 Bruno Chareyre <bruno.chareyre@xxxxxxxxxxx> *
+* *
+* This program is free software; it is licensed under the terms of the *
+* GNU General Public License v2 or later. See file LICENSE for details. *
+*************************************************************************/
+
+#ifdef YADE_OPENGL
+
+#include "Gl1_Primitives.hpp"
+#include <yade/lib/opengl/OpenGLWrapper.hpp>
+#include <yade/core/Scene.hpp>
+
+YADE_PLUGIN((Gl1_Aabb)(Gl1_Box)(Gl1_Facet));
+YADE_PLUGIN((GlBoundFunctor)(GlShapeFunctor)(GlIGeomFunctor)(GlIPhysFunctor)(GlStateFunctor)
+ (GlBoundDispatcher)(GlShapeDispatcher)(GlIGeomDispatcher)(GlIPhysDispatcher)
+ (GlStateDispatcher));
+
+
+void Gl1_Aabb::go(const shared_ptr<Bound>& bv, Scene* scene){
+ Aabb* aabb = static_cast<Aabb*>(bv.get());
+ glColor3v(bv->color);
+ if(!scene->isPeriodic){
+ glTranslatev(Vector3r(.5*(aabb->min+aabb->max)));
+ glScalev(Vector3r(aabb->max-aabb->min));
+ } else {
+ glTranslatev(Vector3r(scene->cell->shearPt(scene->cell->wrapPt(.5*(aabb->min+aabb->max)))));
+ glMultMatrixd(scene->cell->getGlShearTrsfMatrix());
+ glScalev(Vector3r(aabb->max-aabb->min));
+ }
+ glutWireCube(1);
+}
+
+void Gl1_Box::go(const shared_ptr<Shape>& cg, const shared_ptr<State>&,bool wire,const GLViewInfo&)
+{
+ glColor3v(cg->color);
+ Vector3r &extents = (static_cast<Box*>(cg.get()))->extents;
+ glScalef(2*extents[0],2*extents[1],2*extents[2]);
+ if (wire) glutWireCube(1);
+ else glutSolidCube(1);
+}
+
+
+bool Gl1_Facet::normals=false;
+
+void Gl1_Facet::go(const shared_ptr<Shape>& cm, const shared_ptr<State>& ,bool wire,const GLViewInfo&)
+{
+ Facet* facet = static_cast<Facet*>(cm.get());
+ const vector<Vector3r>& vertices = facet->vertices;
+ const Vector3r* ne = facet->ne;
+ const Real& icr = facet->icr;
+
+ if(cm->wire || wire){
+ // facet
+ glBegin(GL_LINE_LOOP);
+ glColor3v(normals ? Vector3r(1,0,0): cm->color);
+ glVertex3v(vertices[0]);
+ glVertex3v(vertices[1]);
+ glVertex3v(vertices[2]);
+ glEnd();
+ if(!normals) return;
+ // facet's normal
+ glBegin(GL_LINES);
+ glColor3(0.0,0.0,1.0);
+ glVertex3(0.0,0.0,0.0);
+ glVertex3v(facet->normal);
+ glEnd();
+ // normal of edges
+ glColor3(0.0,0.0,1.0);
+ glBegin(GL_LINES);
+ glVertex3(0.0,0.0,0.0); glVertex3v(Vector3r(icr*ne[0]));
+ glVertex3(0.0,0.0,0.0); glVertex3v(Vector3r(icr*ne[1]));
+ glVertex3(0.0,0.0,0.0); glVertex3v(Vector3r(icr*ne[2]));
+ glEnd();
+ } else {
+ glDisable(GL_CULL_FACE);
+ Vector3r normal=(facet->vertices[1]-facet->vertices[0]).cross(facet->vertices[2]-facet->vertices[1]); normal.normalize();
+ glColor3v(cm->color);
+ glBegin(GL_TRIANGLES);
+ glNormal3v(normal); // this makes every triangle different WRT the light direction; important!
+ glVertex3v(facet->vertices[0]);
+ glVertex3v(facet->vertices[1]);
+ glVertex3v(facet->vertices[2]);
+ glEnd();
+ }
+}
+
+// Spheres==========================================
+
+bool Gl1_Sphere::wire;
+bool Gl1_Sphere::stripes;
+int Gl1_Sphere::glutSlices;
+int Gl1_Sphere::glutStacks;
+Real Gl1_Sphere::quality;
+bool Gl1_Sphere::localSpecView;
+vector<Vector3r> Gl1_Sphere::vertices, Gl1_Sphere::faces;
+int Gl1_Sphere::glStripedSphereList=-1;
+int Gl1_Sphere::glGlutSphereList=-1;
+Real Gl1_Sphere::prevQuality=0;
+
+void Gl1_Sphere::go(const shared_ptr<Shape>& cm, const shared_ptr<State>& ,bool wire2, const GLViewInfo&)
+{
+ glClearDepth(1.0f);
+ glEnable(GL_NORMALIZE);
+
+ Real r=(static_cast<Sphere*>(cm.get()))->radius;
+ glColor3v(cm->color);
+ if (wire || wire2) glutWireSphere(r,quality*glutSlices,quality*glutStacks);
+ else {
+ //Check if quality has been modified or if previous lists are invalidated (e.g. by creating a new qt view), then regenerate lists
+ bool somethingChanged = (std::abs(quality-prevQuality)>0.001 || glIsList(glStripedSphereList)!=GL_TRUE);
+ if (somethingChanged) {initStripedGlList(); initGlutGlList(); prevQuality=quality;}
+ glScalef(r,r,r);
+ if(stripes) glCallList(glStripedSphereList);
+ else glCallList(glGlutSphereList);
+ }
+ return;
+}
+YADE_PLUGIN((Gl1_Sphere));
+
+void Gl1_Sphere::subdivideTriangle(Vector3r& v1,Vector3r& v2,Vector3r& v3, int depth){
+ Vector3r v;
+ //Change color only at the appropriate level, i.e. 8 times in total, since we draw 8 mono-color sectors one after another
+ if (depth==int(quality) || quality<=0){
+ v = (v1+v2+v3)/3.0;
+ GLfloat matEmit[4];
+ if (v[1]*v[0]*v[2]>0){
+ matEmit[0] = 0.3;
+ matEmit[1] = 0.3;
+ matEmit[2] = 0.3;
+ matEmit[3] = 1.f;
+ }else{
+ matEmit[0] = 0.15;
+ matEmit[1] = 0.15;
+ matEmit[2] = 0.15;
+ matEmit[3] = 0.2;
+ }
+ glMaterialfv(GL_FRONT, GL_EMISSION, matEmit);
+ }
+ if (depth==1){//Then display 4 triangles
+ Vector3r v12 = v1+v2;
+ Vector3r v23 = v2+v3;
+ Vector3r v31 = v3+v1;
+ v12.normalize();
+ v23.normalize();
+ v31.normalize();
+ //Use TRIANGLE_STRIP for faster display of adjacent facets
+ glBegin(GL_TRIANGLE_STRIP);
+ glNormal3v(v1); glVertex3v(v1);
+ glNormal3v(v31); glVertex3v(v31);
+ glNormal3v(v12); glVertex3v(v12);
+ glNormal3v(v23); glVertex3v(v23);
+ glNormal3v(v2); glVertex3v(v2);
+ glEnd();
+ //terminate with this triangle left behind
+ glBegin(GL_TRIANGLES);
+ glNormal3v(v3); glVertex3v(v3);
+ glNormal3v(v23); glVertex3v(v23);
+ glNormal3v(v31); glVertex3v(v31);
+ glEnd();
+ return;
+ }
+ Vector3r v12 = v1+v2;
+ Vector3r v23 = v2+v3;
+ Vector3r v31 = v3+v1;
+ v12.normalize();
+ v23.normalize();
+ v31.normalize();
+ subdivideTriangle(v1,v12,v31,depth-1);
+ subdivideTriangle(v2,v23,v12,depth-1);
+ subdivideTriangle(v3,v31,v23,depth-1);
+ subdivideTriangle(v12,v23,v31,depth-1);
+}
+
+void Gl1_Sphere::initStripedGlList() {
+ if (!vertices.size()){//Fill vectors with vertices and facets
+ //Define 6 points for +/- coordinates
+ vertices.push_back(Vector3r(-1,0,0));//0
+ vertices.push_back(Vector3r(1,0,0));//1
+ vertices.push_back(Vector3r(0,-1,0));//2
+ vertices.push_back(Vector3r(0,1,0));//3
+ vertices.push_back(Vector3r(0,0,-1));//4
+ vertices.push_back(Vector3r(0,0,1));//5
+ //Define 8 sectors of the sphere
+ faces.push_back(Vector3r(3,4,1));
+ faces.push_back(Vector3r(3,0,4));
+ faces.push_back(Vector3r(3,5,0));
+ faces.push_back(Vector3r(3,1,5));
+ faces.push_back(Vector3r(2,1,4));
+ faces.push_back(Vector3r(2,4,0));
+ faces.push_back(Vector3r(2,0,5));
+ faces.push_back(Vector3r(2,5,1));
+ }
+ //Generate the list. Only once for each qtView, or more if quality is modified.
+ glDeleteLists(glStripedSphereList,1);
+ glStripedSphereList = glGenLists(1);
+ glNewList(glStripedSphereList,GL_COMPILE);
+ glEnable(GL_LIGHTING);
+ glShadeModel(GL_SMOOTH);
+ // render the sphere now
+ for (int i=0;i<8;i++)
+ subdivideTriangle(vertices[(unsigned int)faces[i][0]],vertices[(unsigned int)faces[i][1]],vertices[(unsigned int)faces[i][2]],1+ (int) quality);
+ glEndList();
+
+}
+
+void Gl1_Sphere::initGlutGlList(){
+ //Generate the "no-stripes" display list, each time quality is modified
+ glDeleteLists(glGlutSphereList,1);
+ glGlutSphereList = glGenLists(1);
+ glNewList(glGlutSphereList,GL_COMPILE);
+ glEnable(GL_LIGHTING);
+ glShadeModel(GL_SMOOTH);
+ glutSolidSphere(1.0,max(quality*glutSlices,(Real)2.),max(quality*glutStacks,(Real)3.));
+ glEndList();
+}
+#endif /* YADE_OPENGL */
=== added file 'pkg/common/Gl1_Primitives.hpp'
--- pkg/common/Gl1_Primitives.hpp 1970-01-01 00:00:00 +0000
+++ pkg/common/Gl1_Primitives.hpp 2014-07-19 21:36:19 +0000
@@ -0,0 +1,81 @@
+/*************************************************************************
+* Copyright (C) 2004 by Olivier Galizzi *
+* olivier.galizzi@xxxxxxx *
+* *
+* Copyright (C) 2008 by Sergei Dorofeenko *
+* sega@xxxxxxxxxxxxxxxx *
+* *
+* © 2008 Václav Šmilauer <eudoxos@xxxxxxxx> *
+* *
+* © 2008 Bruno Chareyre <bruno.chareyre@xxxxxxxxxxx> *
+* *
+* This program is free software; it is licensed under the terms of the *
+* GNU General Public License v2 or later. See file LICENSE for details. *
+*************************************************************************/
+
+#pragma once
+
+#include<yade/pkg/common/GLDrawFunctors.hpp>
+#include<yade/pkg/common/Aabb.hpp>
+#include<yade/pkg/common/Box.hpp>
+#include<yade/pkg/common/Facet.hpp>
+#include<yade/pkg/common/Sphere.hpp>
+
+class Gl1_Aabb: public GlBoundFunctor{
+ public:
+ virtual void go(const shared_ptr<Bound>&, Scene*);
+ RENDERS(Aabb);
+ YADE_CLASS_BASE_DOC(Gl1_Aabb,GlBoundFunctor,"Render Axis-aligned bounding box (:yref:`Aabb`).");
+};
+REGISTER_SERIALIZABLE(Gl1_Aabb);
+
+
+class Gl1_Box : public GlShapeFunctor{
+ public :
+ virtual void go(const shared_ptr<Shape>&, const shared_ptr<State>&,bool,const GLViewInfo&);
+ RENDERS(Box);
+ YADE_CLASS_BASE_DOC(Gl1_Box,GlShapeFunctor,"Renders :yref:`Box` object");
+};
+
+REGISTER_SERIALIZABLE(Gl1_Box);
+
+class Gl1_Facet : public GlShapeFunctor
+{
+ public:
+ virtual void go(const shared_ptr<Shape>&, const shared_ptr<State>&,bool,const GLViewInfo&);
+ RENDERS(Facet);
+ YADE_CLASS_BASE_DOC_STATICATTRS(Gl1_Facet,GlShapeFunctor,"Renders :yref:`Facet` object",
+ ((bool,normals,false,,"In wire mode, render normals of facets and edges; facet's :yref:`colors<Shape::color>` are disregarded in that case."))
+ );
+};
+
+REGISTER_SERIALIZABLE(Gl1_Facet);
+
+
+class Gl1_Sphere : public GlShapeFunctor{
+ private:
+ // for stripes
+ static vector<Vector3r> vertices, faces;
+ static int glStripedSphereList;
+ static int glGlutSphereList;
+ void subdivideTriangle(Vector3r& v1,Vector3r& v2,Vector3r& v3, int depth);
+ //Generate GlList for GLUT sphere
+ void initGlutGlList();
+ //Generate GlList for sliced spheres
+ void initStripedGlList();
+ //for regenerating glutSphere list if needed
+ static Real prevQuality;
+ public:
+ virtual void go(const shared_ptr<Shape>&, const shared_ptr<State>&,bool,const GLViewInfo&);
+ YADE_CLASS_BASE_DOC_STATICATTRS(Gl1_Sphere,GlShapeFunctor,"Renders :yref:`Sphere` object",
+ ((Real,quality,1.0,,"Change discretization level of spheres. quality>1 for better image quality, at the price of more cpu/gpu usage, 0<quality<1 for faster rendering. If mono-color spheres are displayed (:yref:`Gl1_Sphere::stripes` = False), quality mutiplies :yref:`Gl1_Sphere::glutSlices` and :yref:`Gl1_Sphere::glutStacks`. If striped spheres are displayed (:yref:`Gl1_Sphere::stripes` = True), only integer increments are meaningfull : quality=1 and quality=1.9 will give the same result, quality=2 will give finer result."))
+ ((bool,wire,false,,"Only show wireframe (controlled by ``glutSlices`` and ``glutStacks``."))
+ ((bool,stripes,false,,"In non-wire rendering, show stripes clearly showing particle rotation."))
+ ((bool,localSpecView,true,,"Compute specular light in local eye coordinate system."))
+ ((int,glutSlices,12,(Attr::noSave | Attr::readonly),"Base number of sphere slices, multiplied by :yref:`Gl1_Sphere::quality` before use); not used with ``stripes`` (see `glut{Solid,Wire}Sphere reference <http://www.opengl.org/documentation/specs/glut/spec3/node81.html>`_)"))
+ ((int,glutStacks,6,(Attr::noSave | Attr::readonly),"Base number of sphere stacks, multiplied by :yref:`Gl1_Sphere::quality` before use; not used with ``stripes`` (see `glut{Solid,Wire}Sphere reference <http://www.opengl.org/documentation/specs/glut/spec3/node81.html>`_)"))
+ );
+ RENDERS(Sphere);
+};
+
+REGISTER_SERIALIZABLE(Gl1_Sphere);
=== removed file 'pkg/common/Gl1_Sphere.cpp'
--- pkg/common/Gl1_Sphere.cpp 2014-07-03 07:16:58 +0000
+++ pkg/common/Gl1_Sphere.cpp 1970-01-01 00:00:00 +0000
@@ -1,182 +0,0 @@
-/*************************************************************************
-* © 2004 Olivier Galizzi <olivier.galizzi@xxxxxxx> *
-* © 2008 Václav Šmilauer <eudoxos@xxxxxxxx> *
-* © 2008 Bruno Chareyre <bruno.chareyre@xxxxxxxxxxx> *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#ifdef YADE_OPENGL
-
-#include "Gl1_Sphere.hpp"
-#include<yade/pkg/common/Sphere.hpp>
-#include<yade/lib/opengl/OpenGLWrapper.hpp>
-
-bool Gl1_Sphere::wire;
-bool Gl1_Sphere::stripes;
-int Gl1_Sphere::glutSlices;
-int Gl1_Sphere::glutStacks;
-Real Gl1_Sphere::quality;
-bool Gl1_Sphere::localSpecView;
-vector<Vector3r> Gl1_Sphere::vertices, Gl1_Sphere::faces;
-int Gl1_Sphere::glStripedSphereList=-1;
-int Gl1_Sphere::glGlutSphereList=-1;
-Real Gl1_Sphere::prevQuality=0;
-
-void Gl1_Sphere::go(const shared_ptr<Shape>& cm, const shared_ptr<State>& ,bool wire2, const GLViewInfo&)
-{
- glClearDepth(1.0f);
- glEnable(GL_NORMALIZE);
-
- Real r=(static_cast<Sphere*>(cm.get()))->radius;
- glColor3v(cm->color);
- if (wire || wire2) glutWireSphere(r,quality*glutSlices,quality*glutStacks);
- else {
- //Check if quality has been modified or if previous lists are invalidated (e.g. by creating a new qt view), then regenerate lists
- bool somethingChanged = (std::abs(quality-prevQuality)>0.001 || glIsList(glStripedSphereList)!=GL_TRUE);
- if (somethingChanged) {initStripedGlList(); initGlutGlList(); prevQuality=quality;}
- glScalef(r,r,r);
- if(stripes) glCallList(glStripedSphereList);
- else glCallList(glGlutSphereList);
- }
- return;
-}
-YADE_PLUGIN((Gl1_Sphere));
-
-void Gl1_Sphere::subdivideTriangle(Vector3r& v1,Vector3r& v2,Vector3r& v3, int depth){
- Vector3r v;
- //Change color only at the appropriate level, i.e. 8 times in total, since we draw 8 mono-color sectors one after another
- if (depth==int(quality) || quality<=0){
- v = (v1+v2+v3)/3.0;
- GLfloat matEmit[4];
- if (v[1]*v[0]*v[2]>0){
- matEmit[0] = 0.3;
- matEmit[1] = 0.3;
- matEmit[2] = 0.3;
- matEmit[3] = 1.f;
- }else{
- matEmit[0] = 0.15;
- matEmit[1] = 0.15;
- matEmit[2] = 0.15;
- matEmit[3] = 0.2;
- }
- glMaterialfv(GL_FRONT, GL_EMISSION, matEmit);
- }
- if (depth==1){//Then display 4 triangles
- Vector3r v12 = v1+v2;
- Vector3r v23 = v2+v3;
- Vector3r v31 = v3+v1;
- v12.normalize();
- v23.normalize();
- v31.normalize();
- //Use TRIANGLE_STRIP for faster display of adjacent facets
- glBegin(GL_TRIANGLE_STRIP);
- glNormal3v(v1); glVertex3v(v1);
- glNormal3v(v31); glVertex3v(v31);
- glNormal3v(v12); glVertex3v(v12);
- glNormal3v(v23); glVertex3v(v23);
- glNormal3v(v2); glVertex3v(v2);
- glEnd();
- //terminate with this triangle left behind
- glBegin(GL_TRIANGLES);
- glNormal3v(v3); glVertex3v(v3);
- glNormal3v(v23); glVertex3v(v23);
- glNormal3v(v31); glVertex3v(v31);
- glEnd();
- return;
- }
- Vector3r v12 = v1+v2;
- Vector3r v23 = v2+v3;
- Vector3r v31 = v3+v1;
- v12.normalize();
- v23.normalize();
- v31.normalize();
- subdivideTriangle(v1,v12,v31,depth-1);
- subdivideTriangle(v2,v23,v12,depth-1);
- subdivideTriangle(v3,v31,v23,depth-1);
- subdivideTriangle(v12,v23,v31,depth-1);
-}
-
-void Gl1_Sphere::initStripedGlList() {
- if (!vertices.size()){//Fill vectors with vertices and facets
- //Define 6 points for +/- coordinates
- vertices.push_back(Vector3r(-1,0,0));//0
- vertices.push_back(Vector3r(1,0,0));//1
- vertices.push_back(Vector3r(0,-1,0));//2
- vertices.push_back(Vector3r(0,1,0));//3
- vertices.push_back(Vector3r(0,0,-1));//4
- vertices.push_back(Vector3r(0,0,1));//5
- //Define 8 sectors of the sphere
- faces.push_back(Vector3r(3,4,1));
- faces.push_back(Vector3r(3,0,4));
- faces.push_back(Vector3r(3,5,0));
- faces.push_back(Vector3r(3,1,5));
- faces.push_back(Vector3r(2,1,4));
- faces.push_back(Vector3r(2,4,0));
- faces.push_back(Vector3r(2,0,5));
- faces.push_back(Vector3r(2,5,1));
- }
- //Generate the list. Only once for each qtView, or more if quality is modified.
- glDeleteLists(glStripedSphereList,1);
- glStripedSphereList = glGenLists(1);
- glNewList(glStripedSphereList,GL_COMPILE);
- glEnable(GL_LIGHTING);
- glShadeModel(GL_SMOOTH);
- // render the sphere now
- for (int i=0;i<8;i++)
- subdivideTriangle(vertices[(unsigned int)faces[i][0]],vertices[(unsigned int)faces[i][1]],vertices[(unsigned int)faces[i][2]],1+ (int) quality);
- glEndList();
-
-}
-
-void Gl1_Sphere::initGlutGlList(){
- //Generate the "no-stripes" display list, each time quality is modified
- glDeleteLists(glGlutSphereList,1);
- glGlutSphereList = glGenLists(1);
- glNewList(glGlutSphereList,GL_COMPILE);
- glEnable(GL_LIGHTING);
- glShadeModel(GL_SMOOTH);
- glutSolidSphere(1.0,max(quality*glutSlices,(Real)2.),max(quality*glutStacks,(Real)3.));
- glEndList();
-}
-
-#endif /* YADE_OPENGL */
-
-// ///The old Galizzi's lists
-// if (!vertices.size()) {
-// Real X = 0.525731112119133606;
-// Real Z = 0.850650808352039932;
-// vertices.push_back(Vector3r(-X,0,Z));//0
-// vertices.push_back(Vector3r(X,0,Z));//1
-// vertices.push_back(Vector3r(-X,0,-Z));//2
-// vertices.push_back(Vector3r(X,0,-Z));//3
-// vertices.push_back(Vector3r(0,Z,X));//4
-// vertices.push_back(Vector3r(0,Z,-X));//5
-// vertices.push_back(Vector3r(0,-Z,X));//6
-// vertices.push_back(Vector3r(0,-Z,-X));//7
-// vertices.push_back(Vector3r(Z,X,0));//8
-// vertices.push_back(Vector3r(-Z,X,0));//9
-// vertices.push_back(Vector3r(Z,-X,0));//10
-// vertices.push_back(Vector3r(-Z,-X,0));//11
-// faces.push_back(Vector3r(0,4,1));
-// faces.push_back(Vector3r(0,9,4));
-// faces.push_back(Vector3r(9,5,4));
-// faces.push_back(Vector3r(4,5,8));
-// faces.push_back(Vector3r(4,8,1));
-// faces.push_back(Vector3r(8,10,1));
-// faces.push_back(Vector3r(8,3,10));
-// faces.push_back(Vector3r(5,3,8));
-// faces.push_back(Vector3r(5,2,3));
-// faces.push_back(Vector3r(2,7,3));
-// faces.push_back(Vector3r(7,10,3));
-// faces.push_back(Vector3r(7,6,10));
-// faces.push_back(Vector3r(7,11,6));
-// faces.push_back(Vector3r(11,0,6));
-// faces.push_back(Vector3r(0,1,6));
-// faces.push_back(Vector3r(6,1,10));
-// faces.push_back(Vector3r(9,0,11));
-// faces.push_back(Vector3r(9,11,2));
-// faces.push_back(Vector3r(9,2,5));
-// faces.push_back(Vector3r(7,2,11));}
-
=== removed file 'pkg/common/Gl1_Sphere.hpp'
--- pkg/common/Gl1_Sphere.hpp 2012-11-22 14:32:53 +0000
+++ pkg/common/Gl1_Sphere.hpp 1970-01-01 00:00:00 +0000
@@ -1,44 +0,0 @@
-/*************************************************************************
-* © 2004 Olivier Galizzi <olivier.galizzi@xxxxxxx> *
-* © 2008 Václav Šmilauer <eudoxos@xxxxxxxx> *
-* © 2008 Bruno Chareyre <bruno.chareyre@xxxxxxxxxxx> *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-*************************************************************************/
-
-#pragma once
-
-#include<yade/pkg/common/GLDrawFunctors.hpp>
-#include<yade/pkg/common/Sphere.hpp>
-
-
-class Gl1_Sphere : public GlShapeFunctor{
- private:
- // for stripes
- static vector<Vector3r> vertices, faces;
- static int glStripedSphereList;
- static int glGlutSphereList;
- void subdivideTriangle(Vector3r& v1,Vector3r& v2,Vector3r& v3, int depth);
-// void drawSphere(const Vector3r& color);
- //Generate GlList for GLUT sphere
- void initGlutGlList();
- //Generate GlList for sliced spheres
- void initStripedGlList();
- //for regenerating glutSphere list if needed
- static Real prevQuality;
- public:
- virtual void go(const shared_ptr<Shape>&, const shared_ptr<State>&,bool,const GLViewInfo&);
- YADE_CLASS_BASE_DOC_STATICATTRS(Gl1_Sphere,GlShapeFunctor,"Renders :yref:`Sphere` object",
- ((Real,quality,1.0,,"Change discretization level of spheres. quality>1 for better image quality, at the price of more cpu/gpu usage, 0<quality<1 for faster rendering. If mono-color spheres are displayed (:yref:`Gl1_Sphere::stripes` = False), quality mutiplies :yref:`Gl1_Sphere::glutSlices` and :yref:`Gl1_Sphere::glutStacks`. If striped spheres are displayed (:yref:`Gl1_Sphere::stripes` = True), only integer increments are meaningfull : quality=1 and quality=1.9 will give the same result, quality=2 will give finer result."))
- ((bool,wire,false,,"Only show wireframe (controlled by ``glutSlices`` and ``glutStacks``."))
- ((bool,stripes,false,,"In non-wire rendering, show stripes clearly showing particle rotation."))
- ((bool,localSpecView,true,,"Compute specular light in local eye coordinate system."))
- ((int,glutSlices,12,(Attr::noSave | Attr::readonly),"Base number of sphere slices, multiplied by :yref:`Gl1_Sphere::quality` before use); not used with ``stripes`` (see `glut{Solid,Wire}Sphere reference <http://www.opengl.org/documentation/specs/glut/spec3/node81.html>`_)"))
- ((int,glutStacks,6,(Attr::noSave | Attr::readonly),"Base number of sphere stacks, multiplied by :yref:`Gl1_Sphere::quality` before use; not used with ``stripes`` (see `glut{Solid,Wire}Sphere reference <http://www.opengl.org/documentation/specs/glut/spec3/node81.html>`_)"))
- );
- RENDERS(Sphere);
-};
-REGISTER_SERIALIZABLE(Gl1_Sphere);
-
-
=== modified file 'pkg/common/Grid.cpp'
--- pkg/common/Grid.cpp 2014-07-03 11:31:29 +0000
+++ pkg/common/Grid.cpp 2014-07-25 16:12:46 +0000
@@ -242,7 +242,7 @@
const shared_ptr<Interaction> intr = scene->interactions->find(c->id1,gridNo1->ConnList[i]->getId());
if(intr && intr->isReal()){
shared_ptr<ScGridCoGeom> intrGeom=YADE_PTR_CAST<ScGridCoGeom>(intr->geom);
- if(!intrGeom->isDuplicate==1){ //skip contact.
+ if(!(intrGeom->isDuplicate==1)){ //skip contact.
if (isNew) {return false;}
else {scm->isDuplicate=1;}/*cout<<"Declare "<<c->id1<<"-"<<c->id2<<" as duplicated."<<endl;*/
}
@@ -277,7 +277,7 @@
const shared_ptr<Interaction> intr = scene->interactions->find(c->id1,gridNo2->ConnList[i]->getId());
if(intr && intr->isReal()){
shared_ptr<ScGridCoGeom> intrGeom=YADE_PTR_CAST<ScGridCoGeom>(intr->geom);
- if(!intrGeom->isDuplicate==1){
+ if(!(intrGeom->isDuplicate==1)){
if (isNew) return false;
else scm->isDuplicate=1;/*cout<<"Declare "<<c->id1<<"-"<<c->id2<<" as duplicated."<<endl;*/
}
@@ -397,7 +397,7 @@
//!################## Laws #####################
//! O/
-void Law2_ScGridCoGeom_FrictPhys_CundallStrack::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
+bool Law2_ScGridCoGeom_FrictPhys_CundallStrack::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
int id1 = contact->getId1(), id2 = contact->getId2();
ScGridCoGeom* geom= static_cast<ScGridCoGeom*>(ig.get());
FrictPhys* phys = static_cast<FrictPhys*>(ip.get());
@@ -405,16 +405,11 @@
if (neverErase) {
phys->shearForce = Vector3r::Zero();
phys->normalForce = Vector3r::Zero();}
- else scene->interactions->requestErase(contact);
- return;}
+ else return false;}
if (geom->isDuplicate) {
if (id2!=geom->trueInt) {
//cerr<<"skip duplicate "<<id1<<" "<<id2<<endl;
- if (geom->isDuplicate==2) {
- //cerr<<"erase duplicate "<<id1<<" "<<id2<<endl;
- scene->interactions->requestErase(contact);
- }
- return;
+ if (geom->isDuplicate==2) return false;
}
}
Real& un=geom->penetrationDepth;
@@ -451,11 +446,12 @@
scene->forces.addTorque(geom->id3,(1-geom->relPos)*twist);
scene->forces.addForce(geom->id4,(-geom->relPos)*force);
scene->forces.addTorque(geom->id4,geom->relPos*twist);
+ return true;
}
YADE_PLUGIN((Law2_ScGridCoGeom_FrictPhys_CundallStrack));
-void Law2_ScGridCoGeom_CohFrictPhys_CundallStrack::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
+bool Law2_ScGridCoGeom_CohFrictPhys_CundallStrack::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
const int &id1 = contact->getId1();
const int &id2 = contact->getId2();
ScGridCoGeom* geom = YADE_CAST<ScGridCoGeom*> (ig.get());
@@ -464,11 +460,7 @@
if (geom->isDuplicate) {
if (id2!=geom->trueInt) {
//cerr<<"skip duplicate "<<id1<<" "<<id2<<endl;
- if (geom->isDuplicate==2) {
- //cerr<<"erase duplicate "<<id1<<" "<<id2<<endl;
- scene->interactions->requestErase(contact);
- }
- return;
+ if (geom->isDuplicate==2) return false;
}
}
@@ -479,13 +471,13 @@
if (phys->fragile && (-Fn)> phys->normalAdhesion) {
// BREAK due to tension
- scene->interactions->requestErase(contact); return;
+ return false;
} else {
if ((-Fn)> phys->normalAdhesion) {//normal plasticity
Fn=-phys->normalAdhesion;
phys->unp = un+phys->normalAdhesion/phys->kn;
if (phys->unpMax && phys->unp<phys->unpMax)
- scene->interactions->requestErase(contact); return;
+ return false;
}
phys->normalForce = Fn*geom->normal;
Vector3r& shearForce = geom->rotate(phys->shearForce);
@@ -520,11 +512,12 @@
scene->forces.addTorque(geom->id3,(1-geom->relPos)*twist);
scene->forces.addForce(geom->id4,(-geom->relPos)*force);
scene->forces.addTorque(geom->id4,geom->relPos*twist);
+ return true;
}
}
YADE_PLUGIN((Law2_ScGridCoGeom_CohFrictPhys_CundallStrack));
-void Law2_GridCoGridCoGeom_FrictPhys_CundallStrack::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
+bool Law2_GridCoGridCoGeom_FrictPhys_CundallStrack::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
int id1 = contact->getId1(), id2 = contact->getId2();
id_t id11 = (static_cast<GridConnection*>((&Body::byId(id1)->shape)->get()))->node1->getId();
id_t id12 = (static_cast<GridConnection*>((&Body::byId(id1)->shape)->get()))->node2->getId();
@@ -536,8 +529,7 @@
if (neverErase) {
phys->shearForce = Vector3r::Zero();
phys->normalForce = Vector3r::Zero();}
- else scene->interactions->requestErase(contact);
- return;}
+ else return false;}
Real& un=geom->penetrationDepth;
phys->normalForce=phys->kn*std::max(un,(Real) 0)*geom->normal;
@@ -579,6 +571,7 @@
scene->forces.addTorque(id12,geom->relPos1*torque1);
scene->forces.addTorque(id21,(1-geom->relPos2)*torque2);
scene->forces.addTorque(id22,geom->relPos2*torque2);
+ return true;
}
YADE_PLUGIN((Law2_GridCoGridCoGeom_FrictPhys_CundallStrack));
//!################## Bounds #####################
=== modified file 'pkg/common/Grid.hpp'
--- pkg/common/Grid.hpp 2014-07-14 20:08:33 +0000
+++ pkg/common/Grid.hpp 2014-07-18 18:18:50 +0000
@@ -168,7 +168,7 @@
//! O/
class Law2_ScGridCoGeom_FrictPhys_CundallStrack: public LawFunctor{
public:
- virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(Law2_ScGridCoGeom_FrictPhys_CundallStrack,LawFunctor,"Law between a frictional :yref:`GridConnection` and a frictional :yref:`Sphere`. Almost the same than :yref:`Law2_ScGeom_FrictPhys_CundallStrack`, but the force is divided and applied on the two :yref:`GridNodes<GridNode>` only.",
((bool,neverErase,false,,"Keep interactions even if particles go away from each other (only in case another constitutive law is in the scene, e.g. :yref:`Law2_ScGeom_CapillaryPhys_Capillarity`)"))
((int,plastDissipIx,-1,(Attr::hidden|Attr::noSave),"Index for plastic dissipation (with O.trackEnergy)"))
@@ -181,7 +181,7 @@
class Law2_ScGridCoGeom_CohFrictPhys_CundallStrack: public LawFunctor{
public:
- virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(Law2_ScGridCoGeom_CohFrictPhys_CundallStrack,LawFunctor,"Law between a cohesive frictional :yref:`GridConnection` and a cohesive frictional :yref:`Sphere`. Almost the same than :yref:`Law2_ScGeom6D_CohFrictPhys_CohesionMoment`, but THE ROTATIONAL MOMENTS ARE NOT COMPUTED.",
((bool,neverErase,false,,"Keep interactions even if particles go away from each other (only in case another constitutive law is in the scene, e.g. :yref:`Law2_ScGeom_CapillaryPhys_Capillarity`)"))
((int,plastDissipIx,-1,(Attr::hidden|Attr::noSave),"Index for plastic dissipation (with O.trackEnergy)"))
@@ -194,7 +194,7 @@
//! -/-
class Law2_GridCoGridCoGeom_FrictPhys_CundallStrack: public Law2_ScGeom_FrictPhys_CundallStrack{
public:
- virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
YADE_CLASS_BASE_DOC_ATTRS(Law2_GridCoGridCoGeom_FrictPhys_CundallStrack,Law2_ScGeom_FrictPhys_CundallStrack,"Frictional elastic contact law between two :yref:`gridConnection` . See :yref:`Law2_ScGeom_FrictPhys_CundallStrack` for more details.",
/*ATTRS*/
);
=== modified file 'pkg/common/InteractionLoop.cpp'
--- pkg/common/InteractionLoop.cpp 2014-06-06 19:26:47 +0000
+++ pkg/common/InteractionLoop.cpp 2014-07-18 18:18:50 +0000
@@ -137,7 +137,8 @@
assert(!swap); // reverse call would make no sense, as the arguments are of different types
}
assert(I->functorCache.constLaw);
- I->functorCache.constLaw->go(I->geom,I->phys,I.get());
+ //If the functor return false, the interaction is reset
+ if (!I->functorCache.constLaw->go(I->geom,I->phys,I.get())) scene->interactions->requestErase(I);
// process callbacks for this interaction
if(!I->isReal()) continue; // it is possible that Law2_ functor called requestErase, hence this check
=== modified file 'pkg/common/MatchMaker.cpp'
--- pkg/common/MatchMaker.cpp 2014-07-03 17:20:40 +0000
+++ pkg/common/MatchMaker.cpp 2014-07-23 09:11:35 +0000
@@ -1,9 +1,8 @@
// 2010 © Václav Šmilauer <eudoxos@xxxxxxxx>
-#include<yade/pkg/common/MatchMaker.hpp>
+#include <yade/pkg/common/MatchMaker.hpp>
YADE_PLUGIN((MatchMaker));
-MatchMaker::~MatchMaker(){}
Real MatchMaker::operator()(int id1, int id2, Real val1, Real val2) const {
FOREACH(const Vector3r& m, matches){
=== modified file 'pkg/common/MatchMaker.hpp'
--- pkg/common/MatchMaker.hpp 2014-07-03 17:20:40 +0000
+++ pkg/common/MatchMaker.hpp 2014-07-21 09:09:14 +0000
@@ -24,7 +24,7 @@
bool fbNeedsValues;
#endif
public:
- virtual ~MatchMaker();
+ virtual ~MatchMaker() {};
MatchMaker(std::string _algo): algo(_algo){ postLoad(*this); }
MatchMaker(Real _val): algo("val"), val(_val){ postLoad(*this); }
Real computeFallback(Real val1, Real val2) const ;
=== removed file 'pkg/common/NormShearPhys.cpp'
--- pkg/common/NormShearPhys.cpp 2010-10-13 16:23:08 +0000
+++ pkg/common/NormShearPhys.cpp 1970-01-01 00:00:00 +0000
@@ -1,8 +0,0 @@
-#include"NormShearPhys.hpp"
-YADE_PLUGIN((NormPhys)(NormShearPhys));
-/* At least one virtual function must be in the shared object; let's put empty desctructors here
- * Otherwise downcasting via dynamic_cast will not work (no vtable in the shared lib?)
- */
-NormPhys::~NormPhys(){};
-NormShearPhys::~NormShearPhys(){};
-
=== modified file 'pkg/common/NormShearPhys.hpp'
--- pkg/common/NormShearPhys.hpp 2012-02-14 16:51:38 +0000
+++ pkg/common/NormShearPhys.hpp 2014-07-21 09:09:14 +0000
@@ -6,7 +6,7 @@
class NormPhys:public IPhys {
public:
- virtual ~NormPhys();
+ virtual ~NormPhys() {};
YADE_CLASS_BASE_DOC_ATTRS_CTOR(NormPhys,IPhys,"Abstract class for interactions that have normal stiffness.",
((Real,kn,0,,"Normal stiffness"))
((Vector3r,normalForce,Vector3r::Zero(),,"Normal force after previous step (in global coordinates).")),
@@ -18,7 +18,7 @@
class NormShearPhys: public NormPhys{
public:
- virtual ~NormShearPhys();
+ virtual ~NormShearPhys() {};
YADE_CLASS_BASE_DOC_ATTRS_CTOR(NormShearPhys,NormPhys,
"Abstract class for interactions that have shear stiffnesses, in addition to normal stiffness. This class is used in the PFC3d-style stiffness timestepper.",
((Real,ks,0,,"Shear stiffness"))
=== modified file 'pkg/common/OpenGLRenderer.cpp'
--- pkg/common/OpenGLRenderer.cpp 2014-05-23 13:05:19 +0000
+++ pkg/common/OpenGLRenderer.cpp 2014-07-31 16:09:05 +0000
@@ -9,6 +9,7 @@
#include<yade/core/Timing.hpp>
#include<yade/core/Scene.hpp>
#include<yade/pkg/common/Aabb.hpp>
+#include<yade/lib/pyutil/gil.hpp>
#ifdef __APPLE__
# include <OpenGL/glu.h>
@@ -140,6 +141,7 @@
void OpenGLRenderer::render(const shared_ptr<Scene>& _scene,Body::id_t selection){
+ gilLock lockgil;
if(!initDone) init();
assert(initDone);
selId = selection;
=== modified file 'pkg/common/ParallelEngine.cpp'
--- pkg/common/ParallelEngine.cpp 2014-07-03 17:20:40 +0000
+++ pkg/common/ParallelEngine.cpp 2014-07-23 09:11:35 +0000
@@ -1,11 +1,10 @@
-#include"ParallelEngine.hpp"
+#include <yade/pkg/common/ParallelEngine.hpp>
+YADE_PLUGIN((ParallelEngine));
#ifdef YADE_OPENMP
#include<omp.h>
#endif
-YADE_PLUGIN((ParallelEngine));
-
//! ParallelEngine's pseudo-ctor (factory), taking nested lists of slave engines (might be moved to real ctor perhaps)
shared_ptr<ParallelEngine> ParallelEngine_ctor_list(const boost::python::list& slaves){ shared_ptr<ParallelEngine> instance(new ParallelEngine); instance->slaves_set(slaves); return instance; }
=== removed file 'pkg/common/PeriodicEngines.cpp'
--- pkg/common/PeriodicEngines.cpp 2010-11-07 11:46:20 +0000
+++ pkg/common/PeriodicEngines.cpp 1970-01-01 00:00:00 +0000
@@ -1,3 +0,0 @@
-#include<yade/pkg/common/PeriodicEngines.hpp>
-YADE_PLUGIN((PeriodicEngine));
-PeriodicEngine::~PeriodicEngine(){};
=== modified file 'pkg/common/PeriodicEngines.hpp'
--- pkg/common/PeriodicEngines.hpp 2013-03-26 12:03:38 +0000
+++ pkg/common/PeriodicEngines.hpp 2014-07-19 19:52:41 +0000
@@ -9,7 +9,7 @@
class PeriodicEngine: public GlobalEngine {
public:
static Real getClock(){ timeval tp; gettimeofday(&tp,NULL); return tp.tv_sec+tp.tv_usec/1e6; }
- virtual ~PeriodicEngine(); // vtable
+ virtual ~PeriodicEngine() {}; // vtable
virtual bool isActivated(){
const Real& virtNow=scene->time;
Real realNow=getClock();
@@ -67,44 +67,3 @@
);
};
REGISTER_SERIALIZABLE(PeriodicEngine);
-
-#if 0
- /*!
- PeriodicEngine but with constraint that may be stretched by a given stretchFactor (default 2).
- Limits for each periodicity criterion may be set and the mayStretch bool says whether the period
- can be stretched (default: doubled) without active criteria getting off limits.
-
- stretchFactor must be positive; if >1, period is stretched, for <1, it is shrunk.
-
- Limit consistency (whether actual period is not over/below the limit) is checked: period is set to the
- limit value if we are off. If the limit is zero, however, and the period is non-zero, the limit is set
- to the period value (therefore, if you initialize only iterPeriod, you will get what you expect: engine
- running at iterPeriod).
-
- Note: the logic here is probably too complicated to be practical, although it works. Chances are that
- if unused, it will be removed from the codebase.
- */
- class StretchPeriodicEngine: public PeriodicEngine{
- public:
- StretchPeriodicEngine(): PeriodicEngine(), realLim(0.), virtLim(0.), iterLim(0), stretchFactor(2.), mayStretch(false){}
- Real realLim, virtLim; long iterLim;
- Real stretchFactor;
- bool mayStretch;
- virtual bool isActivated(){
- assert(stretchFactor>0);
- if(iterLim==0 && iterPeriod!=0){iterLim=iterPeriod;} else if(iterLim!=0 && iterPeriod==0){iterPeriod=iterLim;}
- if(realLim==0 && realPeriod!=0){realLim=realPeriod;} else if(realLim!=0 && realPeriod==0){realPeriod=realLim;}
- if(virtLim==0 && virtPeriod!=0){virtLim=virtPeriod;} else if(virtLim!=0 && virtPeriod==0){virtPeriod=virtLim;}
- if(stretchFactor>1){iterPeriod=min(iterPeriod,iterLim); realPeriod=min(realPeriod,realLim); virtPeriod=min(virtPeriod,virtLim);}
- else {iterPeriod=max(iterPeriod,iterLim); realPeriod=max(realPeriod,realLim); virtPeriod=max(virtPeriod,virtLim);}
- mayStretch=((virtPeriod<0 || (stretchFactor>1 ? stretchFactor*virtPeriod<=virtLim : stretchFactor*virtPeriod>=virtLim))
- && (realPeriod<0 || (stretchFactor>1 ? stretchFactor*realPeriod<=realLim : stretchFactor*realPeriod>=realLim))
- && (iterPeriod<0 || (stretchFactor>1 ? stretchFactor*iterPeriod<=iterLim : stretchFactor*iterPeriod>=iterLim)));
- return PeriodicEngine::isActivated();
- }
- REGISTER_ATTRIBUTES(PeriodicEngine,(realLim)(virtLim)(iterLim)(mayStretch)(stretchFactor));
- REGISTER_CLASS_NAME(StretchPeriodicEngine);
- REGISTER_BASE_CLASS_NAME(PeriodicEngine);
- };
- REGISTER_SERIALIZABLE(StretchPeriodicEngine);
-#endif
=== removed file 'pkg/common/PyRunner.cpp'
--- pkg/common/PyRunner.cpp 2010-11-07 11:46:20 +0000
+++ pkg/common/PyRunner.cpp 1970-01-01 00:00:00 +0000
@@ -1,2 +0,0 @@
-#include<yade/pkg/common/PyRunner.hpp>
-YADE_PLUGIN((PyRunner));
=== modified file 'pkg/common/PyRunner.hpp'
--- pkg/common/PyRunner.hpp 2010-11-07 11:46:20 +0000
+++ pkg/common/PyRunner.hpp 2014-07-19 19:52:41 +0000
@@ -15,4 +15,3 @@
);
};
REGISTER_SERIALIZABLE(PyRunner);
-
=== removed file 'pkg/common/Recorder.cpp'
--- pkg/common/Recorder.cpp 2014-05-23 13:05:19 +0000
+++ pkg/common/Recorder.cpp 1970-01-01 00:00:00 +0000
@@ -1,15 +0,0 @@
-// 2009 © Václav Šmilauer <eudoxos@xxxxxxxx>
-#include<yade/pkg/common/Recorder.hpp>
-YADE_PLUGIN((Recorder));
-
-Recorder::~Recorder(){}
-void Recorder::openAndCheck(){
- assert(!out.is_open());
-
- std::string fileTemp = file;
- if (addIterNum) fileTemp+="-" + boost::lexical_cast<string>(scene->iter);
-
- if(fileTemp.empty()) throw ios_base::failure(__FILE__ ": Empty filename.");
- out.open(fileTemp.c_str(), truncate ? fstream::trunc : fstream::app);
- if(!out.good()) throw ios_base::failure(__FILE__ ": I/O error opening file `"+fileTemp+"'.");
-}
=== modified file 'pkg/common/Recorder.hpp'
--- pkg/common/Recorder.hpp 2014-07-03 17:20:40 +0000
+++ pkg/common/Recorder.hpp 2014-07-21 09:09:14 +0000
@@ -2,12 +2,21 @@
#pragma once
#include<yade/pkg/common/PeriodicEngines.hpp>
class Recorder: public PeriodicEngine{
- void openAndCheck();
+ void openAndCheck() {
+ assert(!out.is_open());
+
+ std::string fileTemp = file;
+ if (addIterNum) fileTemp+="-" + boost::lexical_cast<string>(scene->iter);
+
+ if(fileTemp.empty()) throw ios_base::failure(__FILE__ ": Empty filename.");
+ out.open(fileTemp.c_str(), truncate ? fstream::trunc : fstream::app);
+ if(!out.good()) throw ios_base::failure(__FILE__ ": I/O error opening file `"+fileTemp+"'.");
+ }
protected:
//! stream object that derived engines should write to
std::ofstream out;
public:
- virtual ~Recorder(); // vtable
+ virtual ~Recorder() {};
virtual bool isActivated(){
if(PeriodicEngine::isActivated()){
if(!out.is_open()) openAndCheck();
=== removed file 'pkg/common/Sphere.cpp'
--- pkg/common/Sphere.cpp 2010-10-13 16:23:08 +0000
+++ pkg/common/Sphere.cpp 1970-01-01 00:00:00 +0000
@@ -1,3 +0,0 @@
-#include "Sphere.hpp"
-Sphere::~Sphere(){}
-YADE_PLUGIN((Sphere));
=== modified file 'pkg/common/Sphere.hpp'
--- pkg/common/Sphere.hpp 2010-10-13 16:23:08 +0000
+++ pkg/common/Sphere.hpp 2014-07-19 19:52:41 +0000
@@ -8,7 +8,7 @@
class Sphere: public Shape{
public:
Sphere(Real _radius): radius(_radius){ createIndex(); }
- virtual ~Sphere ();
+ virtual ~Sphere () {};
YADE_CLASS_BASE_DOC_ATTRS_CTOR(Sphere,Shape,"Geometry of spherical particle.",
((Real,radius,NaN,,"Radius [m]")),
createIndex(); /*ctor*/
=== removed file 'pkg/common/StepDisplacer.cpp'
--- pkg/common/StepDisplacer.cpp 2010-12-21 12:19:50 +0000
+++ pkg/common/StepDisplacer.cpp 1970-01-01 00:00:00 +0000
@@ -1,24 +0,0 @@
-// 2008 © Václav Šmilauer <eudoxos@xxxxxxxx>
-#include"StepDisplacer.hpp"
-#include<yade/core/State.hpp>
-#include<yade/core/Scene.hpp>
-CREATE_LOGGER(StepDisplacer);
-YADE_PLUGIN((StepDisplacer));
-
-void StepDisplacer::action(){
- FOREACH(Body::id_t id, ids){
- const shared_ptr<Body>& b=Body::byId(id,scene);
- if(setVelocities){
- const Real& dt=scene->dt;
- b->state->vel=mov/dt;
- AngleAxisr aa(rot); aa.axis().normalize();
- b->state->angVel=aa.axis()*aa.angle()/dt;
- LOG_DEBUG("Angular velocity set to "<<aa.axis()*aa.angle()/dt<<". Axis="<<aa.axis()<<", angle="<<aa.angle());
- }
- if(!setVelocities){
- b->state->pos+=mov;
- b->state->ori=rot*b->state->ori;
- }
- }
-}
-
=== modified file 'pkg/common/StepDisplacer.hpp'
--- pkg/common/StepDisplacer.hpp 2011-02-21 14:45:36 +0000
+++ pkg/common/StepDisplacer.hpp 2014-07-21 09:09:14 +0000
@@ -1,11 +1,27 @@
// 2008 © Václav Šmilauer <eudoxos@xxxxxxxx>
#pragma once
-
#include<yade/core/PartialEngine.hpp>
+#include<yade/core/State.hpp>
+#include<yade/core/Scene.hpp>
class StepDisplacer: public PartialEngine {
public:
- virtual void action();
+ virtual void action() {
+ FOREACH(Body::id_t id, ids){
+ const shared_ptr<Body>& b=Body::byId(id,scene);
+ if(setVelocities){
+ const Real& dt=scene->dt;
+ b->state->vel=mov/dt;
+ AngleAxisr aa(rot); aa.axis().normalize();
+ b->state->angVel=aa.axis()*aa.angle()/dt;
+ LOG_DEBUG("Angular velocity set to "<<aa.axis()*aa.angle()/dt<<". Axis="<<aa.axis()<<", angle="<<aa.angle());
+ }
+ if(!setVelocities){
+ b->state->pos+=mov;
+ b->state->ori=rot*b->state->ori;
+ }
+ }
+ }
YADE_CLASS_BASE_DOC_ATTRS(StepDisplacer,PartialEngine,"Apply generalized displacement (displacement or rotation) stepwise on subscribed bodies. Could be used for purposes of contact law tests (by moving one sphere compared to another), but in this case, see rather :yref:`LawTester`",
((Vector3r,mov,Vector3r::Zero(),,"Linear displacement step to be applied per iteration, by addition to :yref:`State.pos`."))
((Quaternionr,rot,Quaternionr::Identity(),,"Rotation step to be applied per iteration (via rotation composition with :yref:`State.ori`)."))
@@ -14,5 +30,3 @@
DECLARE_LOGGER;
};
REGISTER_SERIALIZABLE(StepDisplacer);
-
-
=== removed file 'pkg/common/TorqueEngine.cpp'
--- pkg/common/TorqueEngine.cpp 2010-10-13 16:23:08 +0000
+++ pkg/common/TorqueEngine.cpp 1970-01-01 00:00:00 +0000
@@ -1,12 +0,0 @@
-#include"TorqueEngine.hpp"
-#include<yade/core/Scene.hpp>
-
-void TorqueEngine::action(){
- FOREACH(const Body::id_t id, ids){
- // check that body really exists?
- scene->forces.addTorque(id,moment);
- }
-}
-
-YADE_PLUGIN((TorqueEngine));
-
=== modified file 'pkg/common/TorqueEngine.hpp'
--- pkg/common/TorqueEngine.hpp 2010-10-13 16:23:08 +0000
+++ pkg/common/TorqueEngine.hpp 2014-07-21 09:09:14 +0000
@@ -8,11 +8,17 @@
#pragma once
-#include<yade/core/PartialEngine.hpp>
+#include <yade/core/PartialEngine.hpp>
+#include <yade/core/Scene.hpp>
class TorqueEngine: public PartialEngine{
public:
- virtual void action();
+ virtual void action() {
+ FOREACH(const Body::id_t id, ids){
+ // check that body really exists?
+ scene->forces.addTorque(id,moment);
+ }
+ }
YADE_CLASS_BASE_DOC_ATTRS(TorqueEngine,PartialEngine,"Apply given torque (momentum) value at every subscribed particle, at every step.",
((Vector3r,moment,Vector3r::Zero(),,"Torque value to be applied."))
);
=== added file 'pkg/common/common.cpp'
--- pkg/common/common.cpp 1970-01-01 00:00:00 +0000
+++ pkg/common/common.cpp 2014-07-23 09:11:35 +0000
@@ -0,0 +1,31 @@
+// =======================================================
+// Some plugins from removed CPP-fiels
+
+#include <yade/pkg/dem/DemXDofGeom.hpp>
+#include <yade/pkg/common/TorqueEngine.hpp>
+#include <yade/pkg/common/ForceResetter.hpp>
+#include <yade/pkg/common/FieldApplier.hpp>
+#include <yade/pkg/common/Callbacks.hpp>
+#include <yade/pkg/common/BoundaryController.hpp>
+#include <yade/pkg/common/NormShearPhys.hpp>
+#include <yade/pkg/common/Recorder.hpp>
+#include <yade/pkg/common/CylScGeom6D.hpp>
+#include <yade/pkg/common/Box.hpp>
+#include <yade/pkg/common/StepDisplacer.hpp>
+#include <yade/pkg/common/PeriodicEngines.hpp>
+#include <yade/pkg/common/ElastMat.hpp>
+#include <yade/pkg/common/PyRunner.hpp>
+#include <yade/pkg/common/Sphere.hpp>
+#include <yade/pkg/common/Aabb.hpp>
+
+YADE_PLUGIN((IntrCallback)
+ #ifdef YADE_BODY_CALLBACK
+ (BodyCallback)
+ #endif
+);
+
+YADE_PLUGIN((ForceResetter)(TorqueEngine)(FieldApplier)(BoundaryController)
+ (NormPhys)(NormShearPhys)(Recorder)(CylScGeom6D)(CylScGeom)(Box)
+ (StepDisplacer)(GenericSpheresContact)
+ (PeriodicEngine)(Sphere)(Aabb)(ElastMat)(FrictMat)(PyRunner)
+ );
=== modified file 'pkg/dem/BubbleMat.cpp'
--- pkg/dem/BubbleMat.cpp 2014-07-03 07:16:58 +0000
+++ pkg/dem/BubbleMat.cpp 2014-07-28 06:41:18 +0000
@@ -13,54 +13,70 @@
shared_ptr<BubblePhys> phys(new BubblePhys());
interaction->phys = phys;
- BubbleMat* mat1 = YADE_CAST<BubbleMat*>(m1.get());
- BubbleMat* mat2 = YADE_CAST<BubbleMat*>(m2.get());
-
- // averaging over both materials
- phys->surfaceTension = .5*(mat1->surfaceTension + mat2->surfaceTension);
}
-
/********************** BubblePhys ****************************/
CREATE_LOGGER(BubblePhys);
-Real BubblePhys::computeForce(Real penetrationDepth, Real surfaceTension, Real rAvg, int newtonIter, Real newtonTol) {
- if (penetrationDepth <= 0.0) { return 0.0; }
-
- Real f,df,ll,retOld,residual;
- Real c1 = 1./(4*Mathr::PI*surfaceTension);
- Real c2 = 1./(8*Mathr::PI*surfaceTension*rAvg);
- Real ret=1./c2;
- int i = 0;
- do {
- retOld = ret;
- ll = log( ret*c2 );
- f = penetrationDepth - ret*c1*ll;
- df = -c1*(ll + 1);
- ret -= f/df;
- residual = std::abs(ret - retOld)/retOld;
- if (i++ > newtonIter) {
- throw runtime_error("BubblePhys::computeForce: no convergence\n");
- }
- } while (residual > newtonTol);
- return ret;
-}
-
-
-
+void BubblePhys::computeCoeffs(Real pctMaxForce, Real surfaceTension, Real c1)
+{
+ Real Fmax = pctMaxForce*c1*rAvg;
+ Real logPct = log(pctMaxForce/4);
+ Dmax = pctMaxForce*rAvg*logPct;
+ Real dfdd = c1/(logPct+1);
+ coeffB = dfdd/Fmax;
+ coeffA = Fmax/exp(coeffB*Dmax);
+ fN = 0.1*Fmax;
+}
+
+Real BubblePhys::computeForce(Real separation, Real surfaceTension, Real rAvg, int newtonIter, Real newtonTol, Real c1, Real fN, BubblePhys* phys) {
+
+ if (separation >= phys->Dmax) {
+ Real f,df,g,retOld,residual;
+ Real c2 = 1./(4*c1*rAvg);
+ Real ret = fN;
+ int i = 0;
+ do { //[Chan2011], Loop solves modified form of equation (25) using Newton-Raphson method
+ retOld = ret;
+ g = log(ret*c2);
+ f = separation*c1 - ret*g;
+ df = g+1;
+ ret += f/df;
+ if(ret <= 0.0){ //Need to make sure ret > 0, otherwise the next iteration will try to evaluate the natural logarithm of a negative number, which results in NaN
+ ret = 0.9*fabs(ret);
+ residual = newtonTol*2; //Also, if, by chance, retOld = 0.9*ret it would cause the loop to exit based on the boolean residual > newtonTol, so we force the next iteration by setting residual = newtonTol*2
+ }
+ else {residual = fabs(ret - retOld)/retOld;}
+ if (i++ > newtonIter) {
+ throw runtime_error("BubblePhys::computeForce: no convergence\n");
+ }
+ } while (residual > newtonTol);
+ return ret;
+ }
+ else { //Artificial Extension of [Chan2011] equation 25 to approximiate behaviour outside the valid regime (large penetration cases)
+ return phys->coeffA*exp(phys->coeffB*separation);
+ }
+}
/********************** Law2_ScGeom_BubblePhys_Bubble ****************************/
CREATE_LOGGER(Law2_ScGeom_BubblePhys_Bubble);
-void Law2_ScGeom_BubblePhys_Bubble::go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I){
+bool Law2_ScGeom_BubblePhys_Bubble::go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I){
ScGeom* geom=static_cast<ScGeom*>(_geom.get());
BubblePhys* phys=static_cast<BubblePhys*>(_phys.get());
-
+
+ if(geom->penetrationDepth <= 0.0) {
+ return false;
+ }
+
if (I->isFresh(scene)) {
+ c1 = 2*Mathr::PI*surfaceTension;
phys->rAvg = .5*(geom->refR1 + geom->refR2);
+ phys->computeCoeffs(pctMaxForce,surfaceTension,c1);
}
Real &fN = phys->fN;
- fN = phys->computeForce(geom->penetrationDepth, phys->surfaceTension, phys->rAvg, phys->newtonIter, phys->newtonTol);
+ fN = phys->computeForce(-geom->penetrationDepth, surfaceTension, phys->rAvg, phys->newtonIter, phys->newtonTol, c1, fN, phys);
+ phys->fN = fN;
Vector3r &normalForce = phys->normalForce;
normalForce = fN*geom->normal;
@@ -75,4 +91,5 @@
scene->forces.addTorque(id1,(geom->radius1-0.5*geom->penetrationDepth)* geom->normal.cross(normalForce));
scene->forces.addTorque(id2,(geom->radius2-0.5*geom->penetrationDepth)* geom->normal.cross(normalForce));
}
+ return true;
}
=== modified file 'pkg/dem/BubbleMat.hpp'
--- pkg/dem/BubbleMat.hpp 2013-07-05 15:54:47 +0000
+++ pkg/dem/BubbleMat.hpp 2014-07-28 06:24:36 +0000
@@ -15,7 +15,7 @@
((Real,surfaceTension,0.07197,,"The surface tension in the fluid surrounding the bubbles. The default value is that of water at 25 degrees Celcius."))
,
createIndex();
- density=1; // TODO density default value
+ density=1000;
);
REGISTER_CLASS_INDEX(BubbleMat,Material);
};
@@ -24,9 +24,14 @@
/********************** BubblePhys ****************************/
class BubblePhys : public IPhys {
+ private:
+
+ Real coeffA,coeffB; //Coefficents for artificial curve
+
public:
- static Real computeForce(Real penetrationDepth, Real surfaceTension, Real rAvg, int newtonIter, Real newtonTol);
+ void computeCoeffs(Real pctMaxForce,Real surfaceTension, Real c1);
+ static Real computeForce(Real separation, Real surfaceTension, Real rAvg, int newtonIter, Real newtonTol, Real c1, Real fN, BubblePhys* phys);
virtual ~BubblePhys(){};
YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(BubblePhys,IPhys,"Physics of bubble-bubble interactions, for use with BubbleMat",
@@ -34,6 +39,7 @@
((Real,surfaceTension,NaN,,"Surface tension of the surrounding liquid"))
((Real,fN,NaN,,"Contact normal force"))
((Real,rAvg,NaN,,"Average radius of the two interacting bubbles"))
+ ((Real,Dmax,NaN,,"Maximum penetrationDepth of the bubbles before the force displacement curve changes to an artificial exponential curve. Setting this value will have no effect. See Law2_ScGeom_BubblePhys_Bubble::pctMaxForce for more information"))
((Real,newtonIter,50,,"Maximum number of force iterations allowed"))
((Real,newtonTol,1e-6,,"Convergence criteria for force iterations"))
,
@@ -63,10 +69,16 @@
/********************** Law2_ScGeom_BubblePhys_Bubble ****************************/
class Law2_ScGeom_BubblePhys_Bubble : public LawFunctor{
+ private:
+
+ Real c1; //Coeff used for many contacts
+
public:
- void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* interaction);
+ bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* interaction);
FUNCTOR2D(GenericSpheresContact,BubblePhys);
YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(Law2_ScGeom_BubblePhys_Bubble,LawFunctor,"Constitutive law for Bubble model.",
+ ((Real,pctMaxForce,0.1,,"Chan[2011] states the contact law is valid only for small interferences; therefore an exponential force-displacement curve models the contact stiffness outside that regime (large penetration). This artificial stiffening ensures that bubbles will not pass through eachother or completely overlap during the simulation. The maximum force is Fmax = (2*pi*surfaceTension*rAvg). pctMaxForce is the percentage of the maximum force dictates the separation threshold, Dmax, for each contact. Penetrations less than Dmax calculate the reaction force from the derived contact law, while penetrations equal to or greater than Dmax calculate the reaction force from the artificial exponential curve."))
+ ((Real,surfaceTension,0.07197,,"The surface tension in the liquid surrounding the bubbles. The default value is that of water at 25 degrees Celcius."))
,
/*ctor*/,
);
=== modified file 'pkg/dem/CapillaryTriaxialTest.cpp'
--- pkg/dem/CapillaryTriaxialTest.cpp 2014-07-14 20:08:34 +0000
+++ pkg/dem/CapillaryTriaxialTest.cpp 2014-07-19 19:52:41 +0000
@@ -25,8 +25,7 @@
#include<yade/pkg/common/InsertionSortCollider.hpp>
#include<yade/core/Interaction.hpp>
#include<yade/pkg/common/Dispatching.hpp>
-#include<yade/pkg/common/Bo1_Sphere_Aabb.hpp>
-#include<yade/pkg/common/Bo1_Box_Aabb.hpp>
+#include<yade/pkg/common/Bo1_Aabb.hpp>
#include<yade/pkg/common/GravityEngines.hpp>
#include<yade/pkg/dem/NewtonIntegrator.hpp>
=== modified file 'pkg/dem/CohesiveFrictionalContactLaw.cpp'
--- pkg/dem/CohesiveFrictionalContactLaw.cpp 2014-07-14 20:08:33 +0000
+++ pkg/dem/CohesiveFrictionalContactLaw.cpp 2014-07-18 18:18:50 +0000
@@ -99,7 +99,7 @@
functor->go(I->geom, I->phys, I.get());}
}
-void Law2_ScGeom6D_CohFrictPhys_CohesionMoment::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact)
+bool Law2_ScGeom6D_CohFrictPhys_CohesionMoment::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact)
{
const Real& dt = scene->dt;
const int &id1 = contact->getId1();
@@ -114,13 +114,13 @@
if (phys->fragile && (-Fn)> phys->normalAdhesion) {
// BREAK due to tension
- scene->interactions->requestErase(contact); return;
+ return false;
} else {
if ((-Fn)> phys->normalAdhesion) {//normal plasticity
Fn=-phys->normalAdhesion;
phys->unp = un+phys->normalAdhesion/phys->kn;
if (phys->unpMax && phys->unp<phys->unpMax)
- scene->interactions->requestErase(contact); return;
+ return false;
}
phys->normalForce = Fn*geom->normal;
State* de1 = Body::byId(id1,scene)->state.get();
@@ -223,6 +223,7 @@
}
/// Moment law END ///
}
+ return true;
}
=== modified file 'pkg/dem/CohesiveFrictionalContactLaw.hpp'
--- pkg/dem/CohesiveFrictionalContactLaw.hpp 2014-07-14 20:08:33 +0000
+++ pkg/dem/CohesiveFrictionalContactLaw.hpp 2014-07-18 18:18:50 +0000
@@ -86,7 +86,7 @@
Real totalElastEnergy();
Real getPlasticDissipation();
void initPlasticDissipation(Real initVal=0);
- virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(Law2_ScGeom6D_CohFrictPhys_CohesionMoment,LawFunctor,"Law for linear traction-compression-bending-twisting, with cohesion+friction and Mohr-Coulomb plasticity surface. This law adds adhesion and moments to :yref:`Law2_ScGeom_FrictPhys_CundallStrack`.\n\nThe normal force is (with the convention of positive tensile forces) $F_n=min(k_n*u_n, a_n)$, with $a_n$ the normal adhesion. The shear force is $F_s=k_s*u_s$, the plasticity condition defines the maximum value of the shear force, by default $F_s^{max}=F_n*tan(\\phi)+a_s$, with $\\phi$ the friction angle and $a_s$ the shear adhesion. If :yref:`CohFrictPhys::cohesionDisableFriction` is True, friction is ignored as long as adhesion is active, and the maximum shear force is only $F_s^{max}=a_s$.\n\nIf the maximum tensile or maximum shear force is reached and :yref:`CohFrictPhys::fragile` =True (default), the cohesive link is broken, and $a_n, a_s$ are set back to zero. If a tensile force is present, the contact is lost, else the shear strength is $F_s^{max}=F_n*tan(\\phi)$. If :yref:`CohFrictPhys::fragile` =False (in course of implementation), the behaviour is perfectly plastic, and the shear strength is kept constant.\n\nIf :yref:`Law2_ScGeom6D_CohFrictPhys_CohesionMoment::momentRotationLaw` =True, bending and twisting moments are computed using a linear law with moduli respectively $k_t$ and $k_r$ (the two values are the same currently), so that the moments are : $M_b=k_b*\\Theta_b$ and $M_t=k_t*\\Theta_t$, with $\\Theta_{b,t}$ the relative rotations between interacting bodies. The maximum value of moments can be defined and takes the form of rolling friction. Cohesive -type moment may also be included in the future.\n\nCreep at contact is implemented in this law, as defined in [Hassan2010]_. If activated, there is a viscous behaviour of the shear and twisting components, and the evolution of the elastic parts of shear displacement and relative twist is given by $du_{s,e}/dt=-F_s/\\nu_s$ and $d\\Theta_{t,e}/dt=-M_t/\\nu_t$.",
((bool,neverErase,false,,"Keep interactions even if particles go away from each other (only in case another constitutive law is in the scene, e.g. :yref:`Law2_ScGeom_CapillaryPhys_Capillarity`)"))
((bool,always_use_moment_law,false,,"If true, use bending/twisting moments at all contacts. If false, compute moments only for cohesive contacts."))
=== modified file 'pkg/dem/CohesiveTriaxialTest.cpp'
--- pkg/dem/CohesiveTriaxialTest.cpp 2014-07-14 20:08:34 +0000
+++ pkg/dem/CohesiveTriaxialTest.cpp 2014-07-19 19:52:41 +0000
@@ -29,8 +29,7 @@
#include<yade/core/Body.hpp>
#include<yade/pkg/common/Box.hpp>
#include<yade/pkg/common/Sphere.hpp>
-#include<yade/pkg/common/Bo1_Sphere_Aabb.hpp>
-#include<yade/pkg/common/Bo1_Box_Aabb.hpp>
+#include<yade/pkg/common/Bo1_Aabb.hpp>
#include<yade/pkg/common/ForceResetter.hpp>
=== modified file 'pkg/dem/ConcretePM.cpp'
--- pkg/dem/ConcretePM.cpp 2014-07-03 07:16:58 +0000
+++ pkg/dem/ConcretePM.cpp 2014-07-18 18:18:50 +0000
@@ -285,7 +285,7 @@
#define NNAN(a) YADE_VERIFY(!isnan(a));
#define NNANV(v) YADE_VERIFY(!isnan(v[0])); assert(!isnan(v[1])); assert(!isnan(v[2]));
-void Law2_ScGeom_CpmPhys_Cpm::go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I){
+bool Law2_ScGeom_CpmPhys_Cpm::go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I){
TIMING_DELTAS_START();
ScGeom* geom=static_cast<ScGeom*>(_geom.get());
CpmPhys* phys=static_cast<CpmPhys*>(_phys.get());
@@ -413,8 +413,7 @@
{ boost::mutex::scoped_lock lock(st1->updateMutex); st1->numBrokenCohesive += 1; /* st1->epsPlBroken += epsPlSum; */ }
{ boost::mutex::scoped_lock lock(st2->updateMutex); st2->numBrokenCohesive += 1; /* st2->epsPlBroken += epsPlSum; */ }
/* } */
- scene->interactions->requestErase(I);
- return;
+ return false;
}
Fn = sigmaN*crossSection; phys->normalForce = -Fn*geom->normal;
@@ -438,6 +437,7 @@
scene->forces.addTorque(id2,(geom->radius2+.5*(phys->refPD-geom->penetrationDepth))*geom->normal.cross(f));
}
TIMING_DELTAS_CHECKPOINT("rest");
+ return true;
}
=== modified file 'pkg/dem/ConcretePM.hpp'
--- pkg/dem/ConcretePM.hpp 2014-07-14 20:08:34 +0000
+++ pkg/dem/ConcretePM.hpp 2014-07-18 18:18:50 +0000
@@ -255,7 +255,7 @@
class Law2_ScGeom_CpmPhys_Cpm: public LawFunctor{
public:
- void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
Real elasticEnergy();
Real yieldSigmaTMagnitude(Real sigmaN, Real omega, Real undamagedCohesion, Real tanFrictionAngle) {
=== removed file 'pkg/dem/DemXDofGeom.cpp'
--- pkg/dem/DemXDofGeom.cpp 2013-08-29 10:30:31 +0000
+++ pkg/dem/DemXDofGeom.cpp 1970-01-01 00:00:00 +0000
@@ -1,3 +0,0 @@
-#include"DemXDofGeom.hpp"
-YADE_PLUGIN((GenericSpheresContact));
-GenericSpheresContact::~GenericSpheresContact(){}
=== modified file 'pkg/dem/DemXDofGeom.hpp'
--- pkg/dem/DemXDofGeom.hpp 2013-08-29 10:30:31 +0000
+++ pkg/dem/DemXDofGeom.hpp 2014-07-21 09:09:14 +0000
@@ -16,6 +16,6 @@
);
REGISTER_CLASS_INDEX(GenericSpheresContact,IGeom);
- virtual ~GenericSpheresContact(); // vtable
+ virtual ~GenericSpheresContact() {};
};
REGISTER_SERIALIZABLE(GenericSpheresContact);
=== modified file 'pkg/dem/DomainLimiter.cpp'
--- pkg/dem/DomainLimiter.cpp 2014-07-03 07:16:58 +0000
+++ pkg/dem/DomainLimiter.cpp 2014-07-23 09:11:35 +0000
@@ -1,4 +1,5 @@
#include<yade/pkg/dem/DomainLimiter.hpp>
+#include<yade/pkg/dem/DemXDofGeom.hpp>
#include<yade/pkg/dem/Shop.hpp>
YADE_PLUGIN((DomainLimiter)(LawTester)
=== modified file 'pkg/dem/ElasticContactLaw.cpp'
--- pkg/dem/ElasticContactLaw.cpp 2013-03-06 17:30:45 +0000
+++ pkg/dem/ElasticContactLaw.cpp 2014-07-19 14:26:26 +0000
@@ -46,7 +46,7 @@
}
CREATE_LOGGER(Law2_ScGeom_FrictPhys_CundallStrack);
-void Law2_ScGeom_FrictPhys_CundallStrack::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
+bool Law2_ScGeom_FrictPhys_CundallStrack::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
int id1 = contact->getId1(), id2 = contact->getId2();
ScGeom* geom= static_cast<ScGeom*>(ig.get());
@@ -55,8 +55,8 @@
if (neverErase) {
phys->shearForce = Vector3r::Zero();
phys->normalForce = Vector3r::Zero();}
- else scene->interactions->requestErase(contact);
- return;}
+ else return false;
+ }
Real& un=geom->penetrationDepth;
phys->normalForce=phys->kn*std::max(un,(Real) 0)*geom->normal;
@@ -96,9 +96,10 @@
scene->forces.addTorque(id1,(geom->radius1-0.5*geom->penetrationDepth)* geom->normal.cross(force));
scene->forces.addTorque(id2,(geom->radius2-0.5*geom->penetrationDepth)* geom->normal.cross(force));
}
+ return true;
}
-void Law2_ScGeom_ViscoFrictPhys_CundallStrack::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
+bool Law2_ScGeom_ViscoFrictPhys_CundallStrack::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
ScGeom* geom= static_cast<ScGeom*>(ig.get());
ViscoFrictPhys* phys = static_cast<ViscoFrictPhys*>(ip.get());
if (shearCreep) {
@@ -106,5 +107,5 @@
geom->rotate(phys->creepedShear);
phys->creepedShear+= creepStiffness*phys->ks*(phys->shearForce-phys->creepedShear)*dt/viscosity;
phys->shearForce -= phys->ks*((phys->shearForce-phys->creepedShear)*dt/viscosity);}
- Law2_ScGeom_FrictPhys_CundallStrack::go(ig,ip,contact);
+ return Law2_ScGeom_FrictPhys_CundallStrack::go(ig,ip,contact);
}
=== modified file 'pkg/dem/ElasticContactLaw.hpp'
--- pkg/dem/ElasticContactLaw.hpp 2014-07-14 20:08:34 +0000
+++ pkg/dem/ElasticContactLaw.hpp 2014-07-18 18:18:50 +0000
@@ -16,7 +16,7 @@
class Law2_ScGeom_FrictPhys_CundallStrack: public LawFunctor{
public:
OpenMPAccumulator<Real> plasticDissipation;
- virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
Real elasticEnergy ();
Real getPlasticDissipation();
void initPlasticDissipation(Real initVal=0);
@@ -38,7 +38,7 @@
class Law2_ScGeom_ViscoFrictPhys_CundallStrack: public Law2_ScGeom_FrictPhys_CundallStrack{
public:
- virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(Law2_ScGeom_ViscoFrictPhys_CundallStrack,Law2_ScGeom_FrictPhys_CundallStrack,"Law similar to :yref:`Law2_ScGeom_FrictPhys_CundallStrack` with the addition of shear creep at contacts.",
((bool,shearCreep,false,," "))
((Real,viscosity,1,," "))
=== modified file 'pkg/dem/FrictViscoPM.cpp'
--- pkg/dem/FrictViscoPM.cpp 2014-01-24 21:59:41 +0000
+++ pkg/dem/FrictViscoPM.cpp 2014-07-18 18:18:50 +0000
@@ -154,7 +154,7 @@
CREATE_LOGGER(Law2_ScGeom_FrictViscoPhys_CundallStrackVisco);
-void Law2_ScGeom_FrictViscoPhys_CundallStrackVisco::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact)
+bool Law2_ScGeom_FrictViscoPhys_CundallStrackVisco::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact)
{
LOG_TRACE( "Law2_ScGeom_FrictViscoPhys_CundallStrackVisco::go - create interaction physics" );
@@ -166,8 +166,7 @@
if (neverErase) {
phys->shearForce = Vector3r::Zero();
phys->normalForce = Vector3r::Zero();}
- else scene->interactions->requestErase(contact);
- return;
+ else return false;
}
Real& un=geom->penetrationDepth;
phys->normalForce = phys->kn*std::max(un,(Real) 0) * geom->normal;
@@ -224,5 +223,6 @@
scene->forces.addTorque(id1,(geom->radius1-0.5*geom->penetrationDepth)* geom->normal.cross(force));
scene->forces.addTorque(id2,(geom->radius2-0.5*geom->penetrationDepth)* geom->normal.cross(force));
}
+ return true;
}
=== modified file 'pkg/dem/FrictViscoPM.hpp'
--- pkg/dem/FrictViscoPM.hpp 2014-07-14 20:08:34 +0000
+++ pkg/dem/FrictViscoPM.hpp 2014-07-18 18:18:50 +0000
@@ -95,7 +95,7 @@
class Law2_ScGeom_FrictViscoPhys_CundallStrackVisco: public LawFunctor{
public:
OpenMPAccumulator<Real> plasticDissipation;
- virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
Real elasticEnergy ();
Real getPlasticDissipation();
void initPlasticDissipation(Real initVal=0);
=== modified file 'pkg/dem/HertzMindlin.cpp'
--- pkg/dem/HertzMindlin.cpp 2014-07-14 20:08:33 +0000
+++ pkg/dem/HertzMindlin.cpp 2014-07-25 16:12:46 +0000
@@ -146,20 +146,20 @@
return adhesionEnergy;
}
-void Law2_ScGeom_MindlinPhys_MindlinDeresiewitz::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
+bool Law2_ScGeom_MindlinPhys_MindlinDeresiewitz::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
Body::id_t id1(contact->getId1()), id2(contact->getId2());
ScGeom* geom = static_cast<ScGeom*>(ig.get());
MindlinPhys* phys=static_cast<MindlinPhys*>(ip.get());
const Real uN=geom->penetrationDepth;
if (uN<0) {
- if (neverErase) {phys->shearForce = phys->normalForce = Vector3r::Zero(); phys->kn=phys->ks=0; return;}
- else {scene->interactions->requestErase(contact); return;}
+ if (neverErase) {phys->shearForce = phys->normalForce = Vector3r::Zero(); phys->kn=phys->ks=0; return true;}
+ else {return false;}
}
// normal force
Real Fn=phys->kno*pow(uN,3/2.);
phys->normalForce=Fn*geom->normal;
// exactly zero would not work with the shear formulation, and would give zero shear force anyway
- if(Fn==0) return;
+ if(Fn==0) return true;
//phys->kn=3./2.*phys->kno*std::pow(uN,0.5); // update stiffness, not needed
// contact radius
@@ -186,16 +186,17 @@
scene->forces.addForce(id2,-f);
scene->forces.addTorque(id1,(geom->radius1-.5*geom->penetrationDepth)*geom->normal.cross(f));
scene->forces.addTorque(id2,(geom->radius2-.5*geom->penetrationDepth)*geom->normal.cross(f));
+ return true;
}
-void Law2_ScGeom_MindlinPhys_HertzWithLinearShear::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
+bool Law2_ScGeom_MindlinPhys_HertzWithLinearShear::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
Body::id_t id1(contact->getId1()), id2(contact->getId2());
ScGeom* geom = static_cast<ScGeom*>(ig.get());
MindlinPhys* phys=static_cast<MindlinPhys*>(ip.get());
const Real uN=geom->penetrationDepth;
if (uN<0) {
- if (neverErase) {phys->shearForce = phys->normalForce = Vector3r::Zero(); phys->kn=phys->ks=0; return;}
- else {scene->interactions->requestErase(contact); return;}
+ if (neverErase) {phys->shearForce = phys->normalForce = Vector3r::Zero(); phys->kn=phys->ks=0; return true;}
+ else return false;
}
// normal force
Real Fn=phys->kno*pow(uN,3/2.);
@@ -228,13 +229,14 @@
scene->forces.addForce(id2,-f);
scene->forces.addTorque(id1,(geom->radius1-.5*geom->penetrationDepth)*geom->normal.cross(f));
scene->forces.addTorque(id2,(geom->radius2-.5*geom->penetrationDepth)*geom->normal.cross(f));
+ return true;
}
/******************** Law2_ScGeom_MindlinPhys_Mindlin *********/
CREATE_LOGGER(Law2_ScGeom_MindlinPhys_Mindlin);
-void Law2_ScGeom_MindlinPhys_Mindlin::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
+bool Law2_ScGeom_MindlinPhys_Mindlin::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
const Real& dt = scene->dt; // get time step
Body::id_t id1 = contact->getId1(); // get id body 1
@@ -262,8 +264,8 @@
Real uN = scg->penetrationDepth; // get overlapping
if (uN<0) {
- if (neverErase) {phys->shearForce = phys->normalForce = Vector3r::Zero(); phys->kn=phys->ks=0; return;}
- else {scene->interactions->requestErase(contact); return;}
+ if (neverErase) {phys->shearForce = phys->normalForce = Vector3r::Zero(); phys->kn=phys->ks=0; return true;}
+ else return false;
}
/* Hertz-Mindlin's formulation (PFC)
Note that the normal stiffness here is a secant value (so as it is cannot be used in the GSTS)
@@ -518,7 +520,7 @@
scene->forces.addTorque(id1,-moment);
scene->forces.addTorque(id2,moment);
}
-
+ return true;
// update variables
//phys->prevNormal = scg->normal;
}
=== modified file 'pkg/dem/HertzMindlin.hpp'
--- pkg/dem/HertzMindlin.hpp 2014-07-14 20:08:34 +0000
+++ pkg/dem/HertzMindlin.hpp 2014-07-23 09:07:11 +0000
@@ -1,7 +1,7 @@
// 2010 © Chiara Modenese <c.modenese@xxxxxxxxx>
//
/*
-=== HIGH LEVEL OVERVIEW OF Mindlin ===
+=== HIGH LEVEL OVERVIEW OF MINDLIN ===
Mindlin is a set of classes to include the Hertz-Mindlin formulation for the contact stiffnesses.
The DMT formulation is also considered (for adhesive particles, rigid and small bodies).
@@ -51,11 +51,11 @@
((bool,isAdhesive,false,,"bool to identify if the contact is adhesive, that is to say if the contact force is attractive"))
((bool,isSliding,false,,"check if the contact is sliding (useful to calculate the ratio of sliding contacts)"))
- // Contact damping coefficients as for linear elastic contact law
- ((Real,betan,0.0,,"Fraction of the viscous damping coefficient (normal direction) equal to $\\frac{c_{n}}{C_{n,crit}}$."))
- ((Real,betas,0.0,,"Fraction of the viscous damping coefficient (shear direction) equal to $\\frac{c_{s}}{C_{s,crit}}$."))
+ // Contact damping ratio as for linear elastic contact law
+ ((Real,betan,0.0,,"Normal Damping Ratio. Fraction of the viscous damping coefficient (normal direction) equal to $\\frac{c_{n}}{C_{n,crit}}$."))
+ ((Real,betas,0.0,,"Shear Damping Ratio. Fraction of the viscous damping coefficient (shear direction) equal to $\\frac{c_{s}}{C_{s,crit}}$."))
- // Contact damping coefficient for non-linear elastic contact law (of Hertz-Mindlin type)
+ // Contact damping ratio for non-linear elastic contact law (of Hertz-Mindlin type)
((Real,alpha,0.0,,"Constant coefficient to define contact viscous damping for non-linear elastic force-displacement relationship."))
// temporary
@@ -75,15 +75,28 @@
public :
virtual void go(const shared_ptr<Material>& b1, const shared_ptr<Material>& b2, const shared_ptr<Interaction>& interaction);
FUNCTOR2D(FrictMat,FrictMat);
- YADE_CLASS_BASE_DOC_ATTRS(Ip2_FrictMat_FrictMat_MindlinPhys,IPhysFunctor,"Calculate some physical parameters needed to obtain the normal and shear stiffnesses according to the Hertz-Mindlin's formulation (as implemented in PFC).\n\nViscous parameters can be specified either using coefficients of restitution ($e_n$, $e_s$) or viscous damping coefficient ($\\beta_n$, $\\beta_s$). The following rules apply:\n#. If the $\\beta_n$ ($\\beta_s$) coefficient is given, it is assigned to :yref:`MindlinPhys.betan` (:yref:`MindlinPhys.betas`) directly.\n#. If $e_n$ is given, :yref:`MindlinPhys.betan` is computed using $\\beta_n=-(\\log e_n)/\\sqrt{\\pi^2+(\\log e_n)^2}$. The same applies to $e_s$, :yref:`MindlinPhys.betas`.\n#. It is an error (exception) to specify both $e_n$ and $\\beta_n$ ($e_s$ and $\\beta_s$).\n#. If neither $e_n$ nor $\\beta_n$ is given, zero value for :yref:`MindlinPhys.betan` is used; there will be no viscous effects.\n#.If neither $e_s$ nor $\\beta_s$ is given, the value of :yref:`MindlinPhys.betan` is used for :yref:`MindlinPhys.betas` as well.\n\nThe $e_n$, $\\beta_n$, $e_s$, $\\beta_s$ are :yref:`MatchMaker` objects; they can be constructed from float values to always return constant value.\n\nSee :ysrc:`scripts/test/shots.py` for an example of specifying $e_n$ based on combination of parameters.",
+ YADE_CLASS_BASE_DOC_ATTRS(
+ Ip2_FrictMat_FrictMat_MindlinPhys,IPhysFunctor,"Calculate some physical parameters needed to obtain \
+the normal and shear stiffnesses according to the Hertz-Mindlin formulation (as implemented in PFC).\n\n\
+Viscous parameters can be specified either using coefficients of restitution ($e_n$, $e_s$) or viscous \
+damping ratio ($\\beta_n$, $\\beta_s$). The following rules apply:\n#. If the $\\beta_n$ ($\\beta_s$) \
+ratio is given, it is assigned to :yref:`MindlinPhys.betan` (:yref:`MindlinPhys.betas`) directly.\n#. \
+If $e_n$ is given, :yref:`MindlinPhys.betan` is computed using $\\beta_n=-(\\log e_n)/\\sqrt{\\pi^2+(\\log e_n)^2}$. \
+The same applies to $e_s$, :yref:`MindlinPhys.betas`.\n#. It is an error (exception) to specify both $e_n$ \
+and $\\beta_n$ ($e_s$ and $\\beta_s$).\n#. If neither $e_n$ nor $\\beta_n$ is given, zero value \
+for :yref:`MindlinPhys.betan` is used; there will be no viscous effects.\n#.If neither $e_s$ nor $\\beta_s$ \
+is given, the value of :yref:`MindlinPhys.betan` is used for :yref:`MindlinPhys.betas` as well.\n\nThe \
+$e_n$, $\\beta_n$, $e_s$, $\\beta_s$ are :yref:`MatchMaker` objects; they can be constructed from float \
+values to always return constant value.\n\nSee :ysrc:`scripts/test/shots.py` for an example of specifying \
+$e_n$ based on combination of parameters.",
((Real,gamma,0.0,,"Surface energy parameter [J/m^2] per each unit contact surface, to derive DMT formulation from HM"))
((Real,eta,0.0,,"Coefficient to determine the plastic bending moment"))
((Real,krot,0.0,,"Rotational stiffness for moment contact law"))
((Real,ktwist,0.0,,"Torsional stiffness for moment contact law"))
((shared_ptr<MatchMaker>,en,,,"Normal coefficient of restitution $e_n$."))
((shared_ptr<MatchMaker>,es,,,"Shear coefficient of restitution $e_s$."))
- ((shared_ptr<MatchMaker>,betan,,,"Normal viscous damping coefficient $\\beta_n$."))
- ((shared_ptr<MatchMaker>,betas,,,"Shear viscous damping coefficient $\\beta_s$."))
+ ((shared_ptr<MatchMaker>,betan,,,"Normal viscous damping ratio $\\beta_n$."))
+ ((shared_ptr<MatchMaker>,betas,,,"Shear viscous damping ratio $\\beta_s$."))
((shared_ptr<MatchMaker>,frictAngle,,,"Instance of :yref:`MatchMaker` determining how to compute the friction angle of an interaction. If ``None``, minimum value is used."))
);
DECLARE_LOGGER;
@@ -93,7 +106,7 @@
class Law2_ScGeom_MindlinPhys_MindlinDeresiewitz: public LawFunctor{
public:
- virtual void go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
+ virtual bool go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
FUNCTOR2D(ScGeom,MindlinPhys);
YADE_CLASS_BASE_DOC_ATTRS(Law2_ScGeom_MindlinPhys_MindlinDeresiewitz,LawFunctor,
"Hertz-Mindlin contact law with partial slip solution, as described in [Thornton1991]_.",
@@ -104,7 +117,7 @@
class Law2_ScGeom_MindlinPhys_HertzWithLinearShear: public LawFunctor{
public:
- virtual void go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
+ virtual bool go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
FUNCTOR2D(ScGeom,MindlinPhys);
YADE_CLASS_BASE_DOC_ATTRS(Law2_ScGeom_MindlinPhys_HertzWithLinearShear,LawFunctor,
"Constitutive law for the Hertz formulation (using :yref:`MindlinPhys.kno`) and linear beahvior in shear (using :yref:`MindlinPhys.kso` for stiffness and :yref:`FrictPhys.tangensOfFrictionAngle`). \n\n.. note:: No viscosity or damping. If you need those, look at :yref:`Law2_ScGeom_MindlinPhys_Mindlin`, which also includes non-linear Mindlin shear.",
@@ -119,7 +132,7 @@
class Law2_ScGeom_MindlinPhys_Mindlin: public LawFunctor{
public:
- virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
Real normElastEnergy();
Real adhesionEnergy();
@@ -200,8 +213,8 @@
((Real,ktwist,0.0,,"Torsional stiffness for moment contact law"))
((shared_ptr<MatchMaker>,en,,,"Normal coefficient of restitution $e_n$."))
((shared_ptr<MatchMaker>,es,,,"Shear coefficient of restitution $e_s$."))
- ((shared_ptr<MatchMaker>,betan,,,"Normal viscous damping coefficient $\\beta_n$."))
- ((shared_ptr<MatchMaker>,betas,,,"Shear viscous damping coefficient $\\beta_s$."))
+ ((shared_ptr<MatchMaker>,betan,,,"Normal viscous damping ratio $\\beta_n$."))
+ ((shared_ptr<MatchMaker>,betas,,,"Shear viscous damping ratio $\\beta_s$."))
);
DECLARE_LOGGER;
};
=== modified file 'pkg/dem/InelastCohFrictPM.cpp'
--- pkg/dem/InelastCohFrictPM.cpp 2014-07-14 20:08:34 +0000
+++ pkg/dem/InelastCohFrictPM.cpp 2014-07-18 18:18:50 +0000
@@ -91,7 +91,7 @@
}
-void Law2_ScGeom6D_InelastCohFrictPhys_CohesionMoment::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact)
+bool Law2_ScGeom6D_InelastCohFrictPhys_CohesionMoment::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact)
{
//FIXME : non cohesive contact are not implemented, it would be useful to use setCohesionNow, setCohesionOnNewContacts etc ...
const int &id1 = contact->getId1();
@@ -112,8 +112,7 @@
if(-un>phys->maxExten || phys->isBroken){//plastic failure.
phys->isBroken=1;
phys->normalForce=phys->shearForce=phys->moment_twist=phys->moment_bending=Vector3r(0,0,0);
- scene->interactions->requestErase(contact);
- return;
+ return false;
}
Fn=phys->knT*un; //elasticity
if(-Fn>phys->maxElT || phys->onPlastT){ //so we are on plastic deformation.
@@ -144,9 +143,9 @@
phys->isBroken=1;
phys->normalForce=phys->shearForce=phys->moment_twist=phys->moment_bending=Vector3r(0,0,0);
if(geom->penetrationDepth<=0){ //do not erase the contact while penetrationDepth<0 because it would be recreated at next timestep.
- scene->interactions->requestErase(contact);
+ return false;
}
- return;
+ return true;
}
Fn=phys->knC*un;
if(Fn>phys->maxElC || phys->onPlastC){
@@ -268,4 +267,5 @@
applyForceAtContactPoint(phys->normalForce+phys->shearForce, geom->contactPoint, id1, de1->se3.position, id2, de2->se3.position + (scene->isPeriodic ? scene->cell->intrShiftPos(contact->cellDist): Vector3r::Zero()));
scene->forces.addTorque(id1,-phys->moment_bending-phys->moment_twist);
scene->forces.addTorque(id2,phys->moment_bending+phys->moment_twist);
+ return true;
}
=== modified file 'pkg/dem/InelastCohFrictPM.hpp'
--- pkg/dem/InelastCohFrictPM.hpp 2014-07-14 20:08:34 +0000
+++ pkg/dem/InelastCohFrictPM.hpp 2014-07-18 18:18:50 +0000
@@ -137,7 +137,7 @@
public:
Real normElastEnergy();
Real shearElastEnergy();
- virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(Law2_ScGeom6D_InelastCohFrictPhys_CohesionMoment,LawFunctor,"This law is currently under developpement. Final version and documentation will come before the end of 2014.",
,,
.def("normElastEnergy",&Law2_ScGeom6D_InelastCohFrictPhys_CohesionMoment::normElastEnergy,"Compute normal elastic energy.")
=== modified file 'pkg/dem/JointedCohesiveFrictionalPM.cpp'
--- pkg/dem/JointedCohesiveFrictionalPM.cpp 2014-07-03 07:16:58 +0000
+++ pkg/dem/JointedCohesiveFrictionalPM.cpp 2014-07-18 18:18:50 +0000
@@ -11,7 +11,7 @@
/********************** Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM ****************************/
CREATE_LOGGER(Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM);
-void Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
+bool Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
const int &id1 = contact->getId1();
const int &id2 = contact->getId2();
@@ -42,7 +42,7 @@
if ( smoothJoint && phys->isOnJoint ) {
if ( phys->more || ( phys->jointCumulativeSliding > (2*min(geom->radius1,geom->radius2)) ) ) {
- scene->interactions->requestErase(contact); return;
+ return false;
} else {
D = phys->initD - std::abs((b1->state->pos - b2->state->pos).dot(phys->jointNormal));
}
@@ -54,7 +54,7 @@
if (D < 0) { //spheres do not "touch" !
if ( !phys->isCohesive )
{
- scene->interactions->requestErase(contact); return;
+ return false;
}
if ( phys->isCohesive && (phys->FnMax>0) && (std::abs(D)>Dtensile) ) {
@@ -78,7 +78,7 @@
cracksFileExist=true;
// delete contact
- scene->interactions->requestErase(contact); return;
+ return false;
}
}
@@ -143,7 +143,7 @@
// delete contact if in tension, set the contact properties to friction if in compression
if ( D < 0 )
{
- scene->interactions->requestErase(contact); return;
+ return false;
}
else
{
@@ -165,11 +165,12 @@
scene->forces.addForce (id2, f);
// simple solution to avoid torque computation for particles interacting on a smooth joint
- if ( (phys->isOnJoint)&&(smoothJoint) ) return;
+ if ( (phys->isOnJoint)&&(smoothJoint) ) return true;
/// those lines are needed if rootBody->forces.addForce and rootBody->forces.addMoment are used instead of applyForceAtContactPoint -> NOTE need to check for accuracy!!!
scene->forces.addTorque(id1,(geom->radius1-0.5*geom->penetrationDepth)* geom->normal.cross(-f));
scene->forces.addTorque(id2,(geom->radius2-0.5*geom->penetrationDepth)* geom->normal.cross(-f));
+ return true;
}
=== modified file 'pkg/dem/JointedCohesiveFrictionalPM.hpp'
--- pkg/dem/JointedCohesiveFrictionalPM.hpp 2014-06-10 14:27:43 +0000
+++ pkg/dem/JointedCohesiveFrictionalPM.hpp 2014-07-24 11:35:36 +0000
@@ -59,7 +59,7 @@
((Real,initD,0,,"equilibrium distance for interacting particles. Computed as the interparticular distance at first contact detection."))
((bool,isCohesive,false,,"If false, particles interact in a frictional way. If true, particles are bonded regarding the given :yref:`cohesion<JCFpmMat.cohesion>` and :yref:`tensile strength<JCFpmMat.tensileStrength>` (or their jointed variants)."))
((bool,more,false,,"specifies if the interaction is crossed by more than 3 joints. If true, interaction is deleted (temporary solution)."))
- ((bool,isOnJoint,false,,"defined as true when both interacting particles are :yref:`on joint<JCFpmState.onJoint>` and are in opposite sides of the joint surface. In this case, mechanical parameters of the interaction are derived from the ''joint...'' material properties of the particles, and normal of the interaction is re-oriented (see also :yref:`Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM`)."))
+ ((bool,isOnJoint,false,,"defined as true when both interacting particles are :yref:`on joint<JCFpmState.onJoint>` and are in opposite sides of the joint surface. In this case, mechanical parameters of the interaction are derived from the ''joint...'' material properties of the particles. Furthermore, the normal of the interaction may be re-oriented (see :yref:`Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM.smoothJoint`)."))
((Real,tanFrictionAngle,0,,"tangent of Coulomb friction angle for this interaction (auto. computed). [-]"))
((Real,crossSection,0,,"crossSection=pi*Rmin^2. [m2]"))
((Real,FnMax,0,,"positiv value computed from :yref:`tensile strength<JCFpmMat.tensileStrength>` (or joint variant) to define the maximum admissible normal force in traction: Fn >= -FnMax. [N]"))
@@ -95,7 +95,7 @@
/** 2d functor creating the interaction law (Law2) based on SphereContactGeometry (ScGeom) and JCFpmPhys of 2 bodies, returning type JointedCohesiveFrictionalPM */
class Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM: public LawFunctor{
public:
- virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
FUNCTOR2D(ScGeom,JCFpmPhys);
YADE_CLASS_BASE_DOC_ATTRS(Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM,LawFunctor,"Interaction law for cohesive frictional material, e.g. rock, possibly presenting joint surfaces, that can be mechanically described with a smooth contact logic [Ivars2011]_ (implemented in Yade in [Scholtes2012]_). See examples/jointedCohesiveFrictionalPM for script examples. Joint surface definitions (through stl meshes or direct definition with gts module) are illustrated there.",
=== modified file 'pkg/dem/L3Geom.cpp'
--- pkg/dem/L3Geom.cpp 2014-07-03 17:20:40 +0000
+++ pkg/dem/L3Geom.cpp 2014-07-25 16:12:46 +0000
@@ -129,7 +129,7 @@
Vector3r normTwistVec=avgNormal*scene->dt*.5*avgNormal.dot(state1.angVel+state2.angVel);
// compute relative velocity
// noRatch: take radius or current distance as the branch vector; see discussion in ScGeom::precompute (avoidGranularRatcheting)
- Vector3r c1x=((noRatch && !r1>0) ? ( r1*normal).eval() : (contPt-state1.pos).eval()); // used only for sphere-sphere
+ Vector3r c1x=((noRatch && !(r1>0)) ? ( r1*normal).eval() : (contPt-state1.pos).eval()); // used only for sphere-sphere
Vector3r c2x=(noRatch ? (-r2*normal).eval() : (contPt-state2.pos+shift2).eval());
//Vector3r state2velCorrected=state2.vel+(scene->isPeriodic?scene->cell->intrShiftVel(I->cellDist):Vector3r::Zero()); // velocity of the second particle, corrected with meanfield velocity if necessary
//cerr<<"correction "<<(scene->isPeriodic?scene->cell->intrShiftVel(I->cellDist):Vector3r::Zero())<<endl;
@@ -284,14 +284,14 @@
return true;
}
-void Law2_L3Geom_FrictPhys_ElPerfPl::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* I){
+bool Law2_L3Geom_FrictPhys_ElPerfPl::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* I){
L3Geom* geom=static_cast<L3Geom*>(ig.get()); FrictPhys* phys=static_cast<FrictPhys*>(ip.get());
// compute force
Vector3r& localF(geom->F);
localF=geom->relU().cwiseProduct(Vector3r(phys->kn,phys->ks,phys->ks));
// break if necessary
- if(localF[0]>0 && !noBreak){ scene->interactions->requestErase(I); return; }
+ if(localF[0]>0 && !noBreak) return false;
if(!noSlip){
// plastic slip, if necessary; non-zero elastic limit only for compression
@@ -308,10 +308,11 @@
if(scene->trackEnergy) { scene->energy->add(0.5*(pow(geom->relU()[0],2)*phys->kn+(pow(geom->relU()[1],2)+pow(geom->relU()[2],2))*phys->ks),"elastPotential",elastPotentialIx,/*reset at every timestep*/true); }
// apply force: this converts the force to global space, updates NormShearPhys::{normal,shear}Force, applies to particles
geom->applyLocalForce(localF,I,scene,phys);
+ return true;
}
-void Law2_L6Geom_FrictPhys_Linear::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* I){
+bool Law2_L6Geom_FrictPhys_Linear::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* I){
L6Geom& geom=ig->cast<L6Geom>(); FrictPhys& phys=ip->cast<FrictPhys>();
// simple linear relationships
@@ -319,6 +320,7 @@
Vector3r localT=charLen*(geom.relPhi().cwiseProduct(Vector3r(phys.kn,phys.ks,phys.ks)));
geom.applyLocalForceTorque(localF,localT,I,scene,static_cast<NormShearPhys*>(ip.get()));
+ return true;
}
#ifdef YADE_OPENGL
=== modified file 'pkg/dem/L3Geom.hpp'
--- pkg/dem/L3Geom.hpp 2014-01-20 16:14:58 +0000
+++ pkg/dem/L3Geom.hpp 2014-07-18 18:18:50 +0000
@@ -174,7 +174,7 @@
struct Law2_L3Geom_FrictPhys_ElPerfPl: public LawFunctor{
- virtual void go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
+ virtual bool go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
FUNCTOR2D(L3Geom,FrictPhys);
YADE_CLASS_BASE_DOC_ATTRS(Law2_L3Geom_FrictPhys_ElPerfPl,LawFunctor,"Basic law for testing :yref:`L3Geom`; it bears no cohesion (unless *noBreak* is ``True``), and plastic slip obeys the Mohr-Coulomb criterion (unless *noSlip* is ``True``).",
((bool,noBreak,false,,"Do not break contacts when particles separate."))
@@ -186,7 +186,7 @@
REGISTER_SERIALIZABLE(Law2_L3Geom_FrictPhys_ElPerfPl);
struct Law2_L6Geom_FrictPhys_Linear: public Law2_L3Geom_FrictPhys_ElPerfPl{
- virtual void go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
+ virtual bool go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
FUNCTOR2D(L6Geom,FrictPhys);
YADE_CLASS_BASE_DOC_ATTRS(Law2_L6Geom_FrictPhys_Linear,Law2_L3Geom_FrictPhys_ElPerfPl,"Basic law for testing :yref:`L6Geom` -- linear in both normal and shear sense, without slip or breakage.",
((Real,charLen,1,,"Characteristic length with the meaning of the stiffness ratios bending/shear and torsion/normal."))
=== modified file 'pkg/dem/LudingPM.cpp'
--- pkg/dem/LudingPM.cpp 2013-10-11 19:22:48 +0000
+++ pkg/dem/LudingPM.cpp 2014-07-18 18:18:50 +0000
@@ -71,7 +71,7 @@
return 2.0*a;
}
-void Law2_ScGeom_LudingPhys_Basic::go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I){
+bool Law2_ScGeom_LudingPhys_Basic::go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I){
const ScGeom& geom=*static_cast<ScGeom*>(_geom.get());
LudingPhys& phys=*static_cast<LudingPhys*>(_phys.get());
@@ -82,8 +82,7 @@
const Real Delt = geom.penetrationDepth;
if (Delt<0) {
- scene->interactions->requestErase(I);
- return;
+ return false;
};
const BodyContainer& bodies = *scene->bodies;
@@ -195,6 +194,6 @@
addTorque(id1,-c1x.cross(f),scene);
addTorque(id2, c2x.cross(f),scene);
}
-
+ return true;
}
=== modified file 'pkg/dem/LudingPM.hpp'
--- pkg/dem/LudingPM.hpp 2013-10-11 14:42:51 +0000
+++ pkg/dem/LudingPM.hpp 2014-07-18 18:18:50 +0000
@@ -59,7 +59,7 @@
class Law2_ScGeom_LudingPhys_Basic: public LawFunctor {
public :
- virtual void go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
+ virtual bool go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
private:
Real calculateCapillarForce(const ScGeom& geom, LudingPhys& phys);
FUNCTOR2D(ScGeom,LudingPhys);
=== modified file 'pkg/dem/NormalInelasticPM.cpp'
--- pkg/dem/NormalInelasticPM.cpp 2014-07-14 20:08:34 +0000
+++ pkg/dem/NormalInelasticPM.cpp 2014-07-18 18:18:50 +0000
@@ -14,7 +14,7 @@
YADE_PLUGIN((NormalInelasticMat)(NormalInelasticityPhys)(Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity)(Ip2_2xNormalInelasticMat_NormalInelasticityPhys));
-void Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity::go(shared_ptr<IGeom>& iG, shared_ptr<IPhys>& iP, Interaction* contact)
+bool Law2_ScGeom6D_NormalInelasticityPhys_NormalInelasticity::go(shared_ptr<IGeom>& iG, shared_ptr<IPhys>& iP, Interaction* contact)
{
int id1 = contact->getId1();
int id2 = contact->getId2();
@@ -39,8 +39,7 @@
// Check if there is a real overlap or not. The Ig2... seems to let exist interactions with negative un (= no overlap). Such interactions seem then to have to be deleted here.
if ( un < 0 )
{
- scene->interactions->requestErase(contact);// this, among other things, resets the interaction : geometry and physics variables (as forces, ...) are reset to defaut values
- return;
+ return false;// this, among other things, resets the interaction : geometry and physics variables (as forces, ...) are reset to defaut values
}
@@ -87,10 +86,7 @@
// ******** Tangential force ******* //
if ( un < 0 )
{ // BREAK due to tension
- scene->interactions->requestErase(contact); return;
- // probably not useful anymore
- currentContactPhysics->normalForce = Vector3r::Zero();
- currentContactPhysics->shearForce = Vector3r::Zero();
+ return false;
}
else
{
@@ -145,6 +141,7 @@
}
// ******** Moment law END ******* //
}
+ return true;
}
=== modified file 'pkg/dem/NormalInelasticPM.hpp'
--- pkg/dem/NormalInelasticPM.hpp 2014-07-14 20:08:34 +0000
+++ pkg/dem/NormalInelasticPM.hpp 2014-07-18 18:18:50 +0000
@@ -66,7 +66,7 @@
,maxFs; // maximum value of shear force according to Coulomb-like criterion
Real un; // value of interpenetration in the interaction
public :
- virtual void go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
+ virtual bool go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
FUNCTOR2D(ScGeom,NormalInelasticityPhys);
=== modified file 'pkg/dem/PeriIsoCompressor.hpp'
--- pkg/dem/PeriIsoCompressor.hpp 2014-05-22 15:26:55 +0000
+++ pkg/dem/PeriIsoCompressor.hpp 2014-07-21 09:09:14 +0000
@@ -76,11 +76,8 @@
class Peri3dController: public BoundaryController{
public:
- Vector6r stressOld;//, stressGoal, strainGoal;
- //Vector6i pe, ps;
+ Vector6r stressOld;
Matrix3r sigma, epsilon, epsilonRate, rot, nonrot;
- //int pathSizes[6], pathsCounter[6], lenPe, lenPs;
- //vector<Vector2r>* paths[6];
virtual void action();
YADE_CLASS_BASE_DOC_ATTRS_CTOR_PY(Peri3dController,BoundaryController,"Experimental controller of full strain/stress tensors on periodic cell. Detailed documentation is in py/_extraDocs.py.",
@@ -114,15 +111,10 @@
((Vector6i,pathsCounter,Vector6i::Zero(),Attr::readonly,"Peri3dController internal variable"))
((int,lenPe,NaN,Attr::readonly,"Peri3dController internal variable"))
((int,lenPs,NaN,Attr::readonly,"Peri3dController internal variable"))
- // not yet used
- //((Real,currUnbalanced,NaN,,"current unbalanced force |yupdate|"))
- //((Real,maxUnbalanced,1e-4,,"Maximum unbalanced force"))
,
/*ctor*/
,
/*py*/
-
);
- DECLARE_LOGGER;
};
REGISTER_SERIALIZABLE(Peri3dController);
=== modified file 'pkg/dem/Polyhedra.cpp'
--- pkg/dem/Polyhedra.cpp 2014-07-03 07:16:58 +0000
+++ pkg/dem/Polyhedra.cpp 2014-07-18 18:18:50 +0000
@@ -413,11 +413,11 @@
//**************************************************************************************
// Apply forces on polyhedrons in collision based on geometric configuration
-void PolyhedraVolumetricLaw::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* I){
+bool PolyhedraVolumetricLaw::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* I){
- if (!I->geom) {return;}
+ if (!I->geom) {return true;}
const shared_ptr<PolyhedraGeom>& contactGeom(YADE_PTR_DYN_CAST<PolyhedraGeom>(I->geom));
- if(!contactGeom) {return;}
+ if(!contactGeom) {return true;}
const Body::id_t idA=I->getId1(), idB=I->getId2();
const shared_ptr<Body>& A=Body::byId(idA), B=Body::byId(idB);
@@ -425,13 +425,11 @@
//erase the interaction when aAbB shows separation, otherwise keep it to be able to store previous separating plane for fast detection of separation
if (A->bound->min[0] >= B->bound->max[0] || B->bound->min[0] >= A->bound->max[0] || A->bound->min[1] >= B->bound->max[1] || B->bound->min[1] >= A->bound->max[1] || A->bound->min[2] >= B->bound->max[2] || B->bound->min[2] >= A->bound->max[2]) {
- scene->interactions->requestErase(I);
- phys->normalForce = Vector3r(0.,0.,0.); phys->shearForce = Vector3r(0.,0.,0.);
- return;
+ return false;
}
//zero penetration depth means no interaction force
- if(!(contactGeom->equivalentPenetrationDepth > 1E-18) || !(contactGeom->penetrationVolume > 0)) {phys->normalForce = Vector3r(0.,0.,0.); phys->shearForce = Vector3r(0.,0.,0.); return;}
+ if(!(contactGeom->equivalentPenetrationDepth > 1E-18) || !(contactGeom->penetrationVolume > 0)) {phys->normalForce = Vector3r(0.,0.,0.); phys->shearForce = Vector3r(0.,0.,0.); return true;}
Vector3r normalForce=contactGeom->normal*contactGeom->penetrationVolume*phys->kn;
//shear force: in case the polyhdras are separated and come to contact again, one should not use the previous shear force
@@ -491,6 +489,7 @@
//needed to be able to acces interaction forces in other parts of yade
phys->normalForce = normalForce;
phys->shearForce = shearForce;
+ return true;
}
#endif // YADE_CGAL
=== modified file 'pkg/dem/Polyhedra.hpp'
--- pkg/dem/Polyhedra.hpp 2014-07-03 17:20:40 +0000
+++ pkg/dem/Polyhedra.hpp 2014-07-18 18:18:50 +0000
@@ -242,7 +242,7 @@
class PolyhedraVolumetricLaw: public LawFunctor{
OpenMPAccumulator<Real> plasticDissipation;
- virtual void go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
+ virtual bool go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
Real elasticEnergy ();
Real getPlasticDissipation();
void initPlasticDissipation(Real initVal=0);
=== modified file 'pkg/dem/Shop.hpp'
--- pkg/dem/Shop.hpp 2014-07-03 17:20:40 +0000
+++ pkg/dem/Shop.hpp 2014-08-01 17:25:38 +0000
@@ -139,7 +139,8 @@
static void setContactFriction(Real angleRad);
//! Homothetic change of sizes of spheres and clumps
- static void growParticles(Real multiplier, bool updateMass, bool dynamicOnly, unsigned int discretization);
+ static void growParticles(Real multiplier, bool updateMass, bool dynamicOnly);
//! Change of size of a single sphere or a clump
+ // DEPREC, update wrt growParticles()
static void growParticle(Body::id_t bodyID, Real multiplier, bool updateMass);
};
=== modified file 'pkg/dem/Shop_01.cpp'
--- pkg/dem/Shop_01.cpp 2014-07-14 20:08:34 +0000
+++ pkg/dem/Shop_01.cpp 2014-07-19 19:52:41 +0000
@@ -16,8 +16,7 @@
#include"yade/pkg/dem/ViscoelasticPM.hpp"
#include"yade/pkg/dem/CapillaryPhys.hpp"
-#include"yade/pkg/common/Bo1_Sphere_Aabb.hpp"
-#include"yade/pkg/common/Bo1_Box_Aabb.hpp"
+#include"yade/pkg/common/Bo1_Aabb.hpp"
#include"yade/pkg/dem/NewtonIntegrator.hpp"
#include"yade/pkg/dem/Ig2_Sphere_Sphere_ScGeom.hpp"
#include"yade/pkg/dem/Ig2_Box_Sphere_ScGeom.hpp"
=== modified file 'pkg/dem/Shop_02.cpp'
--- pkg/dem/Shop_02.cpp 2014-07-14 20:08:34 +0000
+++ pkg/dem/Shop_02.cpp 2014-08-01 17:25:38 +0000
@@ -15,8 +15,7 @@
#include<yade/pkg/dem/ViscoelasticPM.hpp>
#include<yade/pkg/dem/CapillaryPhys.hpp>
-#include<yade/pkg/common/Bo1_Sphere_Aabb.hpp>
-#include<yade/pkg/common/Bo1_Box_Aabb.hpp>
+#include<yade/pkg/common/Bo1_Aabb.hpp>
#include<yade/pkg/dem/NewtonIntegrator.hpp>
#include<yade/pkg/dem/Ig2_Sphere_Sphere_ScGeom.hpp>
#include<yade/pkg/dem/Ig2_Box_Sphere_ScGeom.hpp>
@@ -335,7 +334,7 @@
// *** Compute total fabric tensor from the two tensors above ***/
Matrix3r fabricTot(Matrix3r::Zero());
int q(0);
- if(!count==0){ // compute only if there are some interactions
+ if(!count){ // compute only if there are some interactions
q=nStrong*1./count;
fabricTot=(1-q)*fabricWeak+q*fabricStrong;
}
@@ -457,32 +456,31 @@
}
}
-void Shop::growParticles(Real multiplier, bool updateMass, bool dynamicOnly, unsigned int discretization)
+void Shop::growParticles(Real multiplier, bool updateMass, bool dynamicOnly)
{
Scene* scene = Omega::instance().getScene().get();
+ const int sphereIdx = Sphere::getClassIndexStatic();
FOREACH(const shared_ptr<Body>& b,*scene->bodies){
if (dynamicOnly && !b->isDynamic()) continue;
- int ci=b->shape->getClassIndex();
- if(b->isClump() || ci==GridNode::getClassIndexStatic() || ci==GridConnection::getClassIndexStatic()) continue;
- if (updateMass) {b->state->mass*=pow(multiplier,3); b->state->inertia*=pow(multiplier,5);}
- (YADE_CAST<Sphere*> (b->shape.get()))->radius *= multiplier;
- // Clump volume variation with homothetic displacement from its center
- if (b->isClumpMember()) b->state->pos += (multiplier-1) * (b->state->pos - Body::byId(b->clumpId, scene)->state->pos);
- }
- FOREACH(const shared_ptr<Body>& b,*scene->bodies){
- if(b->isClump()){
- Clump* clumpSt = YADE_CAST<Clump*>(b->shape.get());
- clumpSt->updateProperties(b, discretization);
+ //We grow only spheres and clumps
+ if(b->isClump() or sphereIdx == b->shape->getClassIndex()){
+ if (updateMass) {b->state->mass*=pow(multiplier,3); b->state->inertia*=pow(multiplier,5);}
+ // for clumps we updated inertia, nothing else to do
+ if (b->isClump()) continue;
+ // for spheres, we update radius
+ (YADE_CAST<Sphere*> (b->shape.get()))->radius *= multiplier;
+ // and if they are clump members,clump volume variation with homothetic displacement of all members
+ if (b->isClumpMember()) b->state->pos += (multiplier-1) * (b->state->pos - Body::byId(b->clumpId, scene)->state->pos);
}
}
FOREACH(const shared_ptr<Interaction>& ii, *scene->interactions){
- int ci=(*(scene->bodies))[ii->getId1()]->shape->getClassIndex();
- if(ci==GridNode::getClassIndexStatic() || ci==GridConnection::getClassIndexStatic()) continue;
+ int ci1=(*(scene->bodies))[ii->getId1()]->shape->getClassIndex();
+ int ci2=(*(scene->bodies))[ii->getId2()]->shape->getClassIndex();
if (ii->isReal()) {
GenericSpheresContact* contact = YADE_CAST<GenericSpheresContact*>(ii->geom.get());
- if (!dynamicOnly || (*(scene->bodies))[ii->getId1()]->isDynamic())
+ if ((!dynamicOnly || (*(scene->bodies))[ii->getId1()]->isDynamic()) && ci1==sphereIdx)
contact->refR1 = YADE_CAST<Sphere*>((* (scene->bodies))[ii->getId1()]->shape.get())->radius;
- if (!dynamicOnly || (*(scene->bodies))[ii->getId2()]->isDynamic())
+ if ((!dynamicOnly || (*(scene->bodies))[ii->getId2()]->isDynamic()) && ci2==sphereIdx)
contact->refR2 = YADE_CAST<Sphere*>((* (scene->bodies))[ii->getId2()]->shape.get())->radius;
const shared_ptr<FrictPhys>& contactPhysics = YADE_PTR_CAST<FrictPhys>(ii->phys);
contactPhysics->kn*=multiplier; contactPhysics->ks*=multiplier;
=== modified file 'pkg/dem/SimpleShear.cpp'
--- pkg/dem/SimpleShear.cpp 2014-07-14 20:08:34 +0000
+++ pkg/dem/SimpleShear.cpp 2014-07-19 19:52:41 +0000
@@ -29,8 +29,8 @@
#include<yade/pkg/dem/Ig2_Sphere_Sphere_ScGeom.hpp>
#include<yade/pkg/dem/Ig2_Box_Sphere_ScGeom.hpp>
-#include<yade/pkg/common/Bo1_Sphere_Aabb.hpp>
-#include<yade/pkg/common/Bo1_Box_Aabb.hpp>
+#include<yade/pkg/common/Bo1_Aabb.hpp>
+#include<yade/pkg/common/Bo1_Aabb.hpp>
#include<yade/core/Body.hpp>
#include<yade/pkg/common/Box.hpp>
=== modified file 'pkg/dem/Tetra.cpp'
--- pkg/dem/Tetra.cpp 2014-07-03 07:16:58 +0000
+++ pkg/dem/Tetra.cpp 2014-07-18 18:18:50 +0000
@@ -1208,13 +1208,12 @@
#ifdef YADE_CGAL
-void Law2_TTetraSimpleGeom_NormPhys_Simple::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
+bool Law2_TTetraSimpleGeom_NormPhys_Simple::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
int id1 = contact->getId1(), id2 = contact->getId2();
TTetraSimpleGeom* geom= static_cast<TTetraSimpleGeom*>(ig.get());
NormPhys* phys = static_cast<NormPhys*>(ip.get());
if ( geom->flag == 0 || geom->penetrationVolume <= 0. ) {
- scene->interactions->requestErase(contact);
- return;
+ return false;
}
Real& un=geom->penetrationVolume;
phys->normalForce=phys->kn*std::max(un,(Real) 0)*geom->normal;
@@ -1223,6 +1222,7 @@
State* de2 = Body::byId(id2,scene)->state.get();
applyForceAtContactPoint(-phys->normalForce, geom->contactPoint, id1, de1->se3.position, id2, de2->se3.position);
// TODO periodic
+ return true;
}
#endif
=== modified file 'pkg/dem/Tetra.hpp'
--- pkg/dem/Tetra.hpp 2014-07-03 17:20:40 +0000
+++ pkg/dem/Tetra.hpp 2014-07-18 18:18:50 +0000
@@ -169,7 +169,7 @@
class Law2_TTetraSimpleGeom_NormPhys_Simple: public LawFunctor{
public:
- virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
YADE_CLASS_BASE_DOC(Law2_TTetraSimpleGeom_NormPhys_Simple,LawFunctor,"EXPERIMENTAL. TODO");
FUNCTOR2D(TTetraSimpleGeom,NormPhys);
DECLARE_LOGGER;
=== modified file 'pkg/dem/TriaxialStressController.cpp'
--- pkg/dem/TriaxialStressController.cpp 2013-07-16 09:57:52 +0000
+++ pkg/dem/TriaxialStressController.cpp 2014-08-01 17:27:01 +0000
@@ -15,6 +15,7 @@
#include<assert.h>
#include<yade/core/Scene.hpp>
#include<yade/pkg/dem/Shop.hpp>
+#include<yade/core/Clump.hpp>
CREATE_LOGGER(TriaxialStressController);
YADE_PLUGIN((TriaxialStressController));
@@ -105,14 +106,18 @@
if (first) {
BodyContainer::iterator bi = scene->bodies->begin();
BodyContainer::iterator biEnd = scene->bodies->end();
- spheresVolume = 0;
- for ( ; bi!=biEnd; ++bi )
- {
- if((*bi)->isClump()) continue;
+
+ particlesVolume = 0;
+ for ( ; bi!=biEnd; ++bi ) {
const shared_ptr<Body>& b = *bi;
- if ( b->isDynamic() || b->isClumpMember() ) {
+ if (b->isClump()) {
+ const shared_ptr<Clump>& clump = YADE_PTR_CAST<Clump>(b->shape);
+ const shared_ptr<Body>& member = Body::byId(clump->members.begin()->first,scene);
+ particlesVolume += b->state->mass / member->material->density;
+ }
+ else if (b->isDynamic() && !b->isClumpMember()) {
const shared_ptr<Sphere>& sphere = YADE_PTR_CAST<Sphere> ( b->shape );
- spheresVolume += 1.3333333*Mathr::PI*pow ( sphere->radius, 3 );
+ particlesVolume += 1.3333333*Mathr::PI*pow ( sphere->radius, 3 );
}
}
first = false;
@@ -121,7 +126,7 @@
max_vel2=3 * height /(height+width+depth)*max_vel;
max_vel3 =3 * depth /(height+width+depth)*max_vel;
- porosity = ( boxVolume - spheresVolume ) /boxVolume;
+ porosity = ( boxVolume - particlesVolume ) /boxVolume;
position_top = p_top->se3.position.y();
position_bottom = p_bottom->se3.position.y();
position_right = p_right->se3.position.x();
@@ -229,33 +234,8 @@
void TriaxialStressController::controlInternalStress ( Real multiplier )
{
- spheresVolume *= pow ( multiplier,3 );
- BodyContainer::iterator bi = scene->bodies->begin();
- BodyContainer::iterator biEnd = scene->bodies->end();
- for ( ; bi!=biEnd ; ++bi )
- {
- if ( ( *bi )->isDynamic() )
- {
- ( static_cast<Sphere*> ( ( *bi )->shape.get() ) )->radius *= multiplier;
- (*bi)->state->mass*=pow(multiplier,3);
- (*bi)->state->inertia*=pow(multiplier,5);
-
- }
- }
- InteractionContainer::iterator ii = scene->interactions->begin();
- InteractionContainer::iterator iiEnd = scene->interactions->end();
- for (; ii!=iiEnd ; ++ii)
- {
- if ((*ii)->isReal()) {
- ScGeom* contact = static_cast<ScGeom*>((*ii)->geom.get());
- if ((*(scene->bodies))[(*ii)->getId1()]->isDynamic())
- contact->radius1 = static_cast<Sphere*>((* (scene->bodies))[(*ii)->getId1()]->shape.get())->radius;
- if ((* (scene->bodies))[(*ii)->getId2()]->isDynamic())
- contact->radius2 = static_cast<Sphere*>((* (scene->bodies))[(*ii)->getId2()]->shape.get())->radius;
- const shared_ptr<FrictPhys>& contactPhysics = YADE_PTR_CAST<FrictPhys>((*ii)->phys);
- contactPhysics->kn*=multiplier; contactPhysics->ks*=multiplier;
- }
- }
+ particlesVolume *= pow ( multiplier,3 );
+ Shop::growParticles(multiplier,true,true);
}
/*!
=== modified file 'pkg/dem/TriaxialStressController.hpp'
--- pkg/dem/TriaxialStressController.hpp 2014-05-23 13:05:19 +0000
+++ pkg/dem/TriaxialStressController.hpp 2014-08-01 15:37:51 +0000
@@ -40,8 +40,8 @@
//! The values of stresses
Vector3r stress [6];
Vector3r force [6];
- //! Value of spheres volume (solid volume)
- Real spheresVolume;
+ //! Value of particles volume (solid volume of clumps and spheres)
+ Real particlesVolume;
//! Value of box volume
Real boxVolume;
//! Sample porosity
@@ -142,7 +142,8 @@
.def_readonly("strainRate",&TriaxialStressController::getStrainRate,"Current strain rate in a vector d/dt(exx,eyy,ezz).")
.def_readonly("porosity",&TriaxialStressController::porosity,"Porosity of the packing.")
.def_readonly("boxVolume",&TriaxialStressController::boxVolume,"Total packing volume.")
- .def_readonly("spheresVolume",&TriaxialStressController::spheresVolume,"Total volume pf spheres.")
+ .def_readonly("particlesVolume",&TriaxialStressController::particlesVolume,"Total volume of particles (clumps and spheres).")
+ .def_readonly("spheresVolume",&TriaxialStressController::particlesVolume,"Shorthand for :yref:`TriaxialStressController::particlesVolume`")
.def_readonly("max_vel1",&TriaxialStressController::max_vel1,"see :yref:`TriaxialStressController::max_vel` |ycomp|")
.def_readonly("max_vel2",&TriaxialStressController::max_vel2,"see :yref:`TriaxialStressController::max_vel` |ycomp|")
.def_readonly("max_vel3",&TriaxialStressController::max_vel3,"see :yref:`TriaxialStressController::max_vel` |ycomp|")
=== modified file 'pkg/dem/TriaxialTest.cpp'
--- pkg/dem/TriaxialTest.cpp 2014-07-14 20:08:34 +0000
+++ pkg/dem/TriaxialTest.cpp 2014-07-19 19:52:41 +0000
@@ -34,9 +34,7 @@
#include<yade/pkg/dem/Ig2_Sphere_Sphere_ScGeom.hpp>
#include<yade/pkg/dem/Ig2_Box_Sphere_ScGeom.hpp>
#include<yade/pkg/dem/Ig2_Facet_Sphere_ScGeom.hpp>
-#include<yade/pkg/common/Bo1_Sphere_Aabb.hpp>
-#include<yade/pkg/common/Bo1_Box_Aabb.hpp>
-#include<yade/pkg/common/Bo1_Facet_Aabb.hpp>
+#include<yade/pkg/common/Bo1_Aabb.hpp>
#include<yade/pkg/common/Wall.hpp>
#include <boost/numeric/conversion/bounds.hpp>
=== modified file 'pkg/dem/VTKRecorder.cpp'
--- pkg/dem/VTKRecorder.cpp 2014-07-17 06:53:17 +0000
+++ pkg/dem/VTKRecorder.cpp 2014-07-17 08:39:25 +0000
@@ -449,7 +449,7 @@
radii->InsertNextValue(sphere->radius);
if (recActive[REC_BSTRESS]) {
- const Matrix3r& bStress = - bStresses[b->getId()]; // compressive states are negativ for getStressLWForEachBody; I want them as positiv
+ const Matrix3r& bStress = bStresses[b->getId()];
Eigen::SelfAdjointEigenSolver<Matrix3r> solver(bStress); // bStress is probably not symmetric (= self-adjoint for real matrices), but the solver hopefully works (considering only one half of bStress). And, moreover, existence of (real) eigenvalues is not sure for not symmetric bStress..
Matrix3r dirAll = solver.eigenvectors();
Vector3r eigenVal = solver.eigenvalues(); // cf http://eigen.tuxfamily.org/dox/classEigen_1_1SelfAdjointEigenSolver.html#a30caf3c3884a7f4a46b8ec94efd23c5e to be sure that eigenVal[i] * dirAll.col(i) = bStress * dirAll.col(i)
=== modified file 'pkg/dem/VTKRecorder.hpp'
--- pkg/dem/VTKRecorder.hpp 2014-07-16 08:34:55 +0000
+++ pkg/dem/VTKRecorder.hpp 2014-07-17 08:39:25 +0000
@@ -23,7 +23,7 @@
((bool,multiblock,false,,"Use multi-block (``.vtm``) files to store data, rather than separate ``.vtu`` files."))
#endif
((string,fileName,"",,"Base file name; it will be appended with {spheres,intrs,facets}-243100.vtu (unless *multiblock* is ``True``) depending on active recorders and step number (243100 in this case). It can contain slashes, but the directory must exist already."))
- ((vector<string>,recorders,vector<string>(1,string("all")),,"List of active recorders (as strings). ``all`` (the default value) enables all base and generic recorders.\n\n.. admonition:: Base recorders\n\n\tBase recorders save the geometry (unstructured grids) on which other data is defined. They are implicitly activated by many of the other recorders. Each of them creates a new file (or a block, if :yref:`multiblock <VTKRecorder.multiblock>` is set).\n\n\t``spheres``\n\t\tSaves positions and radii (``radii``) of :yref:`spherical<Sphere>` particles.\n\t``facets``\n\t\tSave :yref:`facets<Facet>` positions (vertices).\n\t``boxes``\n\t\tSave :yref:`boxes<Box>` positions (edges).\n\t``intr``\n\t\tStore interactions as lines between nodes at respective particles positions. Additionally stores magnitude of normal (``forceN``) and shear (``absForceT``) forces on interactions (the :yref:`geom<Interaction.geom> must be of type :yref:`NormShearPhys`). \n\n.. admonition:: Generic recorders\n\n\tGeneric recorders do not depend on specific model being used and save commonly useful data.\n\n\t``id``\n\t\tSaves id's (field ``id``) of spheres; active only if ``spheres`` is active.\n\t``mass``\n\t\tSaves masses (field ``mass``) of spheres; active only if ``spheres`` is active.\n\t``clumpId``\n\t\tSaves id's of clumps to which each sphere belongs (field ``clumpId``); active only if ``spheres`` is active.\n\t``colors``\n\t\tSaves colors of :yref:`spheres<Sphere>` and of :yref:`facets<Facet>` (field ``color``); only active if ``spheres`` or ``facets`` are activated.\n\t``mask``\n\t\tSaves groupMasks of :yref:`spheres<Sphere>` and of :yref:`facets<Facet>` (field ``mask``); only active if ``spheres`` or ``facets`` are activated.\n\t``materialId``\n\t\tSaves materialID of :yref:`spheres<Sphere>` and of :yref:`facets<Facet>`; only active if ``spheres`` or ``facets`` are activated.\n\t``velocity``\n\t\tSaves linear and angular velocities of spherical particles as Vector3 and length(fields ``linVelVec``, ``linVelLen`` and ``angVelVec``, ``angVelLen`` respectively``); only effective with ``spheres``.\n\t``stress``\n\t\tSaves stresses of :yref:`spheres<Sphere>` and of :yref:`facets<Facet>` as Vector3 and length; only active if ``spheres`` or ``facets`` are activated.\n\t``pericell``\n\t\tSaves the shape of the cell (simulation has to be periodic).\n\t``bstresses``\n\t\tSaves per-particle principal stresses (sigI >= sigII >= sigIII) and associated principal directions (dirI/II/III). Per-particle stress tensor is computed from :yref:`bodyStressTensors<yade.utils.bodyStressTensors>` considering positiv values for compressive states.\n\n.. admonition:: Specific recorders\n\n\tThe following should only be activated in appropriate cases, otherwise crashes can occur due to violation of type presuppositions.\n\n\t``cpm``\n\t\tSaves data pertaining to the :yref:`concrete model<Law2_ScGeom_CpmPhys_Cpm>`: ``cpmDamage`` (normalized residual strength averaged on particle), ``cpmStress`` (stress on particle); ``intr`` is activated automatically by ``cpm``\n\t``wpm``\n\t\tSaves data pertaining to the :yref:`wire particle model<Law2_ScGeom_WirePhys_WirePM>`: ``wpmForceNFactor`` shows the loading factor for the wire, e.g. normal force divided by threshold normal force.\n\t``jcfpm``\n\t\tSaves data pertaining to the :yref:`rock (smooth)-jointed model<Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM>`: ``damage`` is defined by :yref:`JCFpmState.tensBreak` + :yref:`JCFpmState.shearBreak`; ``intr`` is activated automatically by ``jcfpm``, and :yref:`on joint<JCFpmPhys.isOnJoint>` or :yref:`cohesive<JCFpmPhys.isCohesive>` interactions can be vizualized.\n\t``cracks``\n\t\tSaves other data pertaining to the :yref:`rock model<Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM>`: ``cracks`` shows locations where cohesive bonds failed during the simulation, with their types (0/1 for tensile/shear breakages), their sizes (0.5*(R1+R2)), and their normal directions. The :yref:`corresponding attribute<Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM.recordCracks>` has to be activated, and Key attributes have to be consistent.\n\n"))
+ ((vector<string>,recorders,vector<string>(1,string("all")),,"List of active recorders (as strings). ``all`` (the default value) enables all base and generic recorders.\n\n.. admonition:: Base recorders\n\n\tBase recorders save the geometry (unstructured grids) on which other data is defined. They are implicitly activated by many of the other recorders. Each of them creates a new file (or a block, if :yref:`multiblock <VTKRecorder.multiblock>` is set).\n\n\t``spheres``\n\t\tSaves positions and radii (``radii``) of :yref:`spherical<Sphere>` particles.\n\t``facets``\n\t\tSave :yref:`facets<Facet>` positions (vertices).\n\t``boxes``\n\t\tSave :yref:`boxes<Box>` positions (edges).\n\t``intr``\n\t\tStore interactions as lines between nodes at respective particles positions. Additionally stores magnitude of normal (``forceN``) and shear (``absForceT``) forces on interactions (the :yref:`geom<Interaction.geom> must be of type :yref:`NormShearPhys`). \n\n.. admonition:: Generic recorders\n\n\tGeneric recorders do not depend on specific model being used and save commonly useful data.\n\n\t``id``\n\t\tSaves id's (field ``id``) of spheres; active only if ``spheres`` is active.\n\t``mass``\n\t\tSaves masses (field ``mass``) of spheres; active only if ``spheres`` is active.\n\t``clumpId``\n\t\tSaves id's of clumps to which each sphere belongs (field ``clumpId``); active only if ``spheres`` is active.\n\t``colors``\n\t\tSaves colors of :yref:`spheres<Sphere>` and of :yref:`facets<Facet>` (field ``color``); only active if ``spheres`` or ``facets`` are activated.\n\t``mask``\n\t\tSaves groupMasks of :yref:`spheres<Sphere>` and of :yref:`facets<Facet>` (field ``mask``); only active if ``spheres`` or ``facets`` are activated.\n\t``materialId``\n\t\tSaves materialID of :yref:`spheres<Sphere>` and of :yref:`facets<Facet>`; only active if ``spheres`` or ``facets`` are activated.\n\t``velocity``\n\t\tSaves linear and angular velocities of spherical particles as Vector3 and length(fields ``linVelVec``, ``linVelLen`` and ``angVelVec``, ``angVelLen`` respectively``); only effective with ``spheres``.\n\t``stress``\n\t\tSaves stresses of :yref:`spheres<Sphere>` and of :yref:`facets<Facet>` as Vector3 and length; only active if ``spheres`` or ``facets`` are activated.\n\t``pericell``\n\t\tSaves the shape of the cell (simulation has to be periodic).\n\t``bstresses``\n\t\tSaves per-particle principal stresses (sigI >= sigII >= sigIII) and associated principal directions (dirI/II/III). Per-particle stress tensors are given by :yref:`bodyStressTensors<yade.utils.bodyStressTensors>` (positive values for tensile states).\n\n.. admonition:: Specific recorders\n\n\tThe following should only be activated in appropriate cases, otherwise crashes can occur due to violation of type presuppositions.\n\n\t``cpm``\n\t\tSaves data pertaining to the :yref:`concrete model<Law2_ScGeom_CpmPhys_Cpm>`: ``cpmDamage`` (normalized residual strength averaged on particle), ``cpmStress`` (stress on particle); ``intr`` is activated automatically by ``cpm``\n\t``wpm``\n\t\tSaves data pertaining to the :yref:`wire particle model<Law2_ScGeom_WirePhys_WirePM>`: ``wpmForceNFactor`` shows the loading factor for the wire, e.g. normal force divided by threshold normal force.\n\t``jcfpm``\n\t\tSaves data pertaining to the :yref:`rock (smooth)-jointed model<Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM>`: ``damage`` is defined by :yref:`JCFpmState.tensBreak` + :yref:`JCFpmState.shearBreak`; ``intr`` is activated automatically by ``jcfpm``, and :yref:`on joint<JCFpmPhys.isOnJoint>` or :yref:`cohesive<JCFpmPhys.isCohesive>` interactions can be vizualized.\n\t``cracks``\n\t\tSaves other data pertaining to the :yref:`rock model<Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM>`: ``cracks`` shows locations where cohesive bonds failed during the simulation, with their types (0/1 for tensile/shear breakages), their sizes (0.5*(R1+R2)), and their normal directions. The :yref:`corresponding attribute<Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM.recordCracks>` has to be activated, and Key attributes have to be consistent.\n\n"))
((string,Key,"",,"Necessary if :yref:`recorders<VTKRecorder.recorders>` contains 'cracks'. A string specifying the name of file 'cracks___.txt' that is considered in this case (see :yref:`corresponding attribute<Law2_ScGeom_JCFpmPhys_JointedCohesiveFrictionalPM.Key>`)."))
((int,mask,0,,"If mask defined, only bodies with corresponding groupMask will be exported. If 0, all bodies will be exported.")),
/*ctor*/
=== modified file 'pkg/dem/ViscoelasticCapillarPM.cpp'
--- pkg/dem/ViscoelasticCapillarPM.cpp 2014-06-18 07:35:57 +0000
+++ pkg/dem/ViscoelasticCapillarPM.cpp 2014-07-21 09:09:14 +0000
@@ -62,11 +62,11 @@
}
/* Law2_ScGeom_ViscElCapPhys_Basic */
-void Law2_ScGeom_ViscElCapPhys_Basic::go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I) {
+bool Law2_ScGeom_ViscElCapPhys_Basic::go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I) {
Vector3r force = Vector3r::Zero();
- const int id1 = I->getId1();
- const int id2 = I->getId2();
+ const id_t id1 = I->getId1();
+ const id_t id2 = I->getId2();
const ScGeom& geom=*static_cast<ScGeom*>(_geom.get());
Scene* scene=Omega::instance().getScene().get();
@@ -111,7 +111,7 @@
addForce (id1,-phys.normalForce,scene);
addForce (id2, phys.normalForce,scene);
};
- return;
+ return true;
} else {
if (phys.liqBridgeActive) {
VLiqBridg -= phys.Vb;
@@ -127,8 +127,7 @@
scene->delIntrs.push_back(B1);
scene->delIntrs.push_back(B2);
#endif
- scene->interactions->requestErase(I);
- return;
+ return false;
};
};
@@ -150,6 +149,7 @@
addTorque(id1, torque1,scene);
addTorque(id2, torque2,scene);
}
+ return true;
}
Real Law2_ScGeom_ViscElCapPhys_Basic::critDist(const Real& Vb, const Real& R, const Real& Theta) {
=== modified file 'pkg/dem/ViscoelasticCapillarPM.hpp'
--- pkg/dem/ViscoelasticCapillarPM.hpp 2014-06-25 14:46:14 +0000
+++ pkg/dem/ViscoelasticCapillarPM.hpp 2014-07-18 18:18:50 +0000
@@ -61,7 +61,7 @@
/// Constitutive law
class Law2_ScGeom_ViscElCapPhys_Basic: public LawFunctor {
public :
- virtual void go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
+ virtual bool go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
static Real Willett_numeric_f (const ScGeom& geom, ViscElCapPhys& phys);
static Real Willett_analytic_f (const ScGeom& geom, ViscElCapPhys& phys);
static Real Weigert_f (const ScGeom& geom, ViscElCapPhys& phys);
=== modified file 'pkg/dem/ViscoelasticPM.cpp'
--- pkg/dem/ViscoelasticPM.cpp 2014-07-03 07:16:58 +0000
+++ pkg/dem/ViscoelasticPM.cpp 2014-07-18 18:18:50 +0000
@@ -29,7 +29,7 @@
}
/* Law2_ScGeom_ViscElPhys_Basic */
-void Law2_ScGeom_ViscElPhys_Basic::go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I) {
+bool Law2_ScGeom_ViscElPhys_Basic::go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I) {
Vector3r force = Vector3r::Zero();
Vector3r torque1 = Vector3r::Zero();
Vector3r torque2 = Vector3r::Zero();
@@ -41,11 +41,8 @@
addForce (id2, force,scene);
addTorque(id1, torque1,scene);
addTorque(id2, torque2,scene);
- return;
- } else {
- scene->interactions->requestErase(I);
- return;
- }
+ return true;
+ } else return false;
}
bool computeForceTorqueViscEl(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I, Vector3r & force, Vector3r & torque1, Vector3r & torque2) {
=== modified file 'pkg/dem/ViscoelasticPM.hpp'
--- pkg/dem/ViscoelasticPM.hpp 2014-06-03 08:53:31 +0000
+++ pkg/dem/ViscoelasticPM.hpp 2014-07-18 18:18:50 +0000
@@ -90,7 +90,7 @@
/// This class provides linear viscoelastic contact model
class Law2_ScGeom_ViscElPhys_Basic: public LawFunctor {
public :
- virtual void go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
+ virtual bool go(shared_ptr<IGeom>&, shared_ptr<IPhys>&, Interaction*);
public :
FUNCTOR2D(ScGeom,ViscElPhys);
YADE_CLASS_BASE_DOC(Law2_ScGeom_ViscElPhys_Basic,LawFunctor,"Linear viscoelastic model operating on :yref:`ScGeom` and :yref:`ViscElPhys`. The model is mostly based on the paper for For details see Pournin [Pournin2001]_ .");
=== modified file 'pkg/dem/WirePM.cpp'
--- pkg/dem/WirePM.cpp 2014-07-03 07:16:58 +0000
+++ pkg/dem/WirePM.cpp 2014-07-18 18:18:50 +0000
@@ -80,7 +80,7 @@
/********************** Law2_ScGeom_WirePhys_WirePM ****************************/
CREATE_LOGGER(Law2_ScGeom_WirePhys_WirePM);
-void Law2_ScGeom_WirePhys_WirePM::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
+bool Law2_ScGeom_WirePhys_WirePM::go(shared_ptr<IGeom>& ig, shared_ptr<IPhys>& ip, Interaction* contact){
LOG_TRACE( "Law2_ScGeom_WirePhys_WirePM::go - contact law" );
@@ -102,8 +102,7 @@
/* check whether the particles are linked or not */
if ( !phys->isLinked ) { // destroy the interaction before calculation
- scene->interactions->requestErase(contact);
- return;
+ return false;
}
if ( (phys->isLinked) && (D < DFValues.back()(0)) ) { // spheres are linked but failure because of reaching maximal admissible displacement
phys->isLinked=false;
@@ -112,8 +111,7 @@
WireState* st2=dynamic_cast<WireState*>(b2->state.get());
st1->numBrokenLinks+=1;
st2->numBrokenLinks+=1;
- scene->interactions->requestErase(contact);
- return;
+ return false;
}
/* compute normal force Fn */
@@ -163,6 +161,7 @@
/* set shear force to zero */
phys->shearForce = Vector3r::Zero();
+ return true;
}
=== modified file 'pkg/dem/WirePM.hpp'
--- pkg/dem/WirePM.hpp 2013-08-14 08:01:55 +0000
+++ pkg/dem/WirePM.hpp 2014-07-18 18:18:50 +0000
@@ -116,7 +116,7 @@
/** 2d functor creating the interaction law (Law2) based on SphereContactGeometry (ScGeom) and WirePhys of 2 bodies, returning type WirePM */
class Law2_ScGeom_WirePhys_WirePM: public LawFunctor{
public:
- virtual void go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
+ virtual bool go(shared_ptr<IGeom>& _geom, shared_ptr<IPhys>& _phys, Interaction* I);
FUNCTOR2D(ScGeom,WirePhys);
=== removed file 'pkg/lbm/LBMbody.cpp'
--- pkg/lbm/LBMbody.cpp 2014-03-31 16:01:57 +0000
+++ pkg/lbm/LBMbody.cpp 1970-01-01 00:00:00 +0000
@@ -1,16 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2009-2012 by Franck Lominé *
-* franck.lomine@xxxxxxxxxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-* *
-*************************************************************************/
-#ifdef LBM_ENGINE
-
-#include "LBMbody.hpp"
-
-YADE_PLUGIN((LBMbody));
-LBMbody::~LBMbody(){};
-
-#endif //LBM_ENGINE
=== modified file 'pkg/lbm/LBMbody.hpp'
--- pkg/lbm/LBMbody.hpp 2014-03-31 16:01:57 +0000
+++ pkg/lbm/LBMbody.hpp 2014-07-21 18:14:48 +0000
@@ -14,7 +14,7 @@
class LBMbody: public Serializable{
public:
- virtual ~LBMbody();
+ virtual ~LBMbody() {};
//Real radius(){return ext[0];}
bool isBox(){if(type==1)return true; else return false;}
bool isPtc(){if(type==2)return true; else return false;}
=== removed file 'pkg/lbm/LBMlink.cpp'
--- pkg/lbm/LBMlink.cpp 2014-03-31 16:01:57 +0000
+++ pkg/lbm/LBMlink.cpp 1970-01-01 00:00:00 +0000
@@ -1,17 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2009-2012 by Franck Lominé *
-* franck.lomine@xxxxxxxxxxxxxx *
-* *
-* This program is free software; it is licensed under the terms of the *
-* GNU General Public License v2 or later. See file LICENSE for details. *
-* *
-*************************************************************************/
-#ifdef LBM_ENGINE
-
-#include "LBMlink.hpp"
-
-YADE_PLUGIN((LBMlink));
-LBMlink::~LBMlink(){};
-void LBMlink::ReinitDynamicalProperties(){sid=-1;fid=-1;idx_sigma_i=-1;isBd=false;VbMid=Vector3r::Zero();DistMid=Vector3r::Zero();ct=0.;return;}
-
-#endif //LBM_ENGINE
=== modified file 'pkg/lbm/LBMlink.hpp'
--- pkg/lbm/LBMlink.hpp 2014-03-31 16:01:57 +0000
+++ pkg/lbm/LBMlink.hpp 2014-07-21 18:14:48 +0000
@@ -14,8 +14,15 @@
class LBMlink: public Serializable{
public:
- void ReinitDynamicalProperties();
- virtual ~LBMlink();
+ void ReinitDynamicalProperties() {
+ sid=-1; fid=-1;
+ idx_sigma_i=-1; isBd=false;
+ VbMid=Vector3r::Zero();
+ DistMid=Vector3r::Zero();
+ ct=0.;
+ return;
+ };
+ virtual ~LBMlink() {};
YADE_CLASS_BASE_DOC_ATTRS_CTOR(LBMlink,Serializable,
"Link class for Lattice Boltzmann Method ",
=== modified file 'pkg/lbm/LBMnode.cpp'
--- pkg/lbm/LBMnode.cpp 2014-03-31 16:01:57 +0000
+++ pkg/lbm/LBMnode.cpp 2014-07-21 18:14:48 +0000
@@ -9,10 +9,14 @@
#ifdef LBM_ENGINE
#include "LBMnode.hpp"
+#include "LBMlink.hpp"
+#include "LBMbody.hpp"
-YADE_PLUGIN((LBMnode));
+YADE_PLUGIN((LBMnode)(LBMlink)(LBMbody));
LBMnode::~LBMnode(){};
+
+
void LBMnode::MixteBC(string lbmodel,Real density, Vector3r U, string where){
Real rhoVx=density*U.x();
Real rhoVy=density*U.y();
=== modified file 'pkg/pfv/FlowEngine.hpp.in'
--- pkg/pfv/FlowEngine.hpp.in 2014-06-20 15:28:21 +0000
+++ pkg/pfv/FlowEngine.hpp.in 2014-08-01 11:16:21 +0000
@@ -40,6 +40,7 @@
#include<yade/lib/triangulation/FlowBoundingSphere.hpp>
#include<yade/lib/triangulation/PeriodicFlow.hpp>
#include<yade/pkg/pfv/FlowEngine.hpp>
+#include<yade/core/Clump.hpp>
template<class _CellInfo, class _VertexInfo, class _Tesselation=CGT::_Tesselation<CGT::TriangulationTypes<_VertexInfo,_CellInfo> >, class solverT=DEFAULTSOLVER >
class TemplateFlowEngine_@TEMPLATE_FLOW_NAME@ : public PartialEngine
@@ -62,7 +63,7 @@
shared_ptr<FlowSolver> backgroundSolver;
volatile bool backgroundCompleted;
Cell cachedCell;
- struct posData {Body::id_t id; Vector3r pos; Real radius; bool isSphere; bool exists; posData(){exists=0;}};
+ struct posData {Body::id_t id; Vector3r pos; Real radius; bool isSphere; bool isClump; bool exists; posData(){exists=0;}};
vector<posData> positionBufferCurrent;//reflect last known positions before we start computations
vector<posData> positionBufferParallel;//keep the positions from a given step for multithread factorization
//copy positions in a buffer for faster and/or parallel access
=== modified file 'pkg/pfv/FlowEngine.ipp.in'
--- pkg/pfv/FlowEngine.ipp.in 2014-07-03 07:16:58 +0000
+++ pkg/pfv/FlowEngine.ipp.in 2014-08-01 11:16:21 +0000
@@ -275,14 +275,20 @@
buffer.clear();
buffer.resize(scene->bodies->size());
shared_ptr<Sphere> sph ( new Sphere );
- const int Sph_Index = sph->getClassIndexStatic();
+ const int Sph_Index = sph->getClassIndexStatic();
FOREACH ( const shared_ptr<Body>& b, *scene->bodies ) {
- if (!b || ignoredBody==b->getId()) continue;
- posData& dat = buffer[b->getId()];
+ if (!b || ignoredBody==b->getId() || b->isClumpMember()) continue;
+ posData& dat = buffer[b->getId()];
dat.id=b->getId();
dat.pos=b->state->pos;
dat.isSphere= (b->shape->getClassIndex() == Sph_Index);
+ dat.isClump = b->isClump();
if (dat.isSphere) dat.radius = YADE_CAST<Sphere*>(b->shape.get())->radius;
+ if (dat.isClump) {
+ const shared_ptr<Clump>& clump = YADE_PTR_CAST<Clump>(b->shape);
+ const shared_ptr<Body>& member = Body::byId(clump->members.begin()->first,scene);
+ dat.radius = pow( (3*b->state->mass)/(4*Mathr::PI*member->material->density) , 1.0/3.0); //use equivalent radius of clump (just valid for nearly spherical clumps)
+ }
dat.exists=true;
}
}
@@ -293,7 +299,7 @@
solver->xMin = Mathr::MAX_REAL, solver->xMax = -Mathr::MAX_REAL, solver->yMin = Mathr::MAX_REAL, solver->yMax = -Mathr::MAX_REAL, solver->zMin = Mathr::MAX_REAL, solver->zMax = -Mathr::MAX_REAL;
FOREACH ( const posData& b, buffer ) {
if ( !b.exists ) continue;
- if ( b.isSphere ) {
+ if ( b.isSphere || b.isClump ) {
const Real& rad = b.radius;
const Real& x = b.pos[0];
const Real& y = b.pos[1];
@@ -358,11 +364,9 @@
///Using one-by-one insertion
vector<posData>& buffer = multithread ? positionBufferParallel : positionBufferCurrent;
FOREACH ( const posData& b, buffer ) {
- if ( !b.exists ) continue;
- if ( b.isSphere ) {
- if (b.id==ignoredBody) continue;
- flow.tesselation().insert ( b.pos[0], b.pos[1], b.pos[2], b.radius, b.id );}
- }
+ if ( !b.exists || b.id==ignoredBody ) continue;
+ if ( b.isSphere || b.isClump ) flow.tesselation().insert ( b.pos[0], b.pos[1], b.pos[2], b.radius, b.id );
+ }
flow.tesselation().redirected=true;//By inserting one-by-one, we already redirected
flow.shearLubricationForces.resize ( flow.tesselation().maxId+1 );
flow.shearLubricationTorques.resize ( flow.tesselation().maxId+1 );
=== modified file 'py/_utils.cpp'
--- py/_utils.cpp 2014-07-03 17:20:40 +0000
+++ py/_utils.cpp 2014-08-01 17:25:38 +0000
@@ -495,7 +495,7 @@
py::def("calm",Shop__calm,(py::arg("mask")=-1),"Set translational and rotational velocities of bodies to zero. Applied to all bodies by default. To calm only some bodies, use mask parameter, it will calm only bodies with groupMask compatible to given value");
py::def("setNewVerticesOfFacet",setNewVerticesOfFacet,(py::arg("b"),py::arg("v1"),py::arg("v2"),py::arg("v3")),"Sets new vertices (in global coordinates) to given facet.");
py::def("setContactFriction",Shop::setContactFriction,py::arg("angleRad"),"Modify the friction angle (in radians) inside the material classes and existing contacts. The friction for non-dynamic bodies is not modified.");
- py::def("growParticles",Shop::growParticles,(py::args("multiplier"), py::args("updateMass")=true, py::args("dynamicOnly")=true, py::args("discretization")=0), "Change the size of spheres and sphere clumps by the multiplier. If updateMass=True, then the mass is updated. dynamicOnly=True is mandatory in many cases since in current implementation the function would crash on the non-spherical and non-dynamic bodies (e.g. facets, boxes, etc.). For clumps the masses and inertias are adapted automatically when discretization>0 (for details of inertia tensor integration scheme see :yref:`clump()<BodyContainer.clump>`).");
+ py::def("growParticles",Shop::growParticles,(py::args("multiplier"), py::args("updateMass")=true, py::args("dynamicOnly")=true), "Change the size of spheres and sphere clumps by the multiplier. If updateMass=True, then the mass and inertia are updated. dynamicOnly=True will select dynamic bodies.");
py::def("growParticle",Shop::growParticle,(py::args("bodyID"),py::args("multiplier"), py::args("updateMass")=true), "Change the size of a single sphere (to be implemented: single clump). If updateMass=True, then the mass is updated.");
py::def("intrsOfEachBody",intrsOfEachBody,"returns list of lists of interactions of each body");
py::def("numIntrsOfEachBody",numIntrsOfEachBody,"returns list of number of interactions of each body");
=== modified file 'py/tests/__init__.py'
--- py/tests/__init__.py 2013-08-29 10:30:31 +0000
+++ py/tests/__init__.py 2014-07-22 17:23:17 +0000
@@ -1,7 +1,7 @@
# encoding: utf-8
# 2009 © Václav Šmilauer <eudoxos@xxxxxxxx>
"""All defined functionality tests for yade."""
-import unittest,inspect
+import unittest,inspect,sys
# add any new test suites to the list here, so that they are picked up by testAll
allTests=['wrapper','core','pbc','clump','cohesive-chain']
@@ -25,7 +25,7 @@
@param module: fully-qualified module name, e.g. yade.tests.wrapper
"""
suite=unittest.defaultTestLoader().loadTestsFromName(module)
- return unittest.TextTestRunner(verbosity=2).run(suite)
+ return unittest.TextTestRunner(stream=sys.stdout,verbosity=2).run(suite)
def testAll():
"""Run all tests defined in all yade.tests.* modules and return
@@ -34,7 +34,7 @@
import doctest
for mod in allModules:
suite.addTest(doctest.DocTestSuite(mod))
- return unittest.TextTestRunner(verbosity=2).run(suite)
+ return unittest.TextTestRunner(stream=sys.stdout,verbosity=2).run(suite)
=== modified file 'py/utils.py'
--- py/utils.py 2014-07-01 05:36:28 +0000
+++ py/utils.py 2014-07-22 17:23:17 +0000
@@ -37,12 +37,8 @@
(1, 2, 3)
those variables will be save in the .xml file, when the simulation itself is saved. To recover those variables once the .xml is loaded again, use
-
- >>> loadVars('something')
-
- and they will be defined in the yade.params.\ *mark* module. The *loadNow* parameter calls :yref:`yade.utils.loadVars` after saving automatically.
-
- If 'something' already exists, given variables will be inserted.
+ ``loadVars('something')``and they will be defined in the yade.params.\ *mark* module. The *loadNow* parameter calls :yref:`yade.utils.loadVars`
+ after saving automatically. If 'something' already exists, given variables will be inserted.
"""
import cPickle
try:
=== modified file 'scripts/checks-and-tests/checks/DEM-PFV-check.py'
--- scripts/checks-and-tests/checks/DEM-PFV-check.py 2014-05-27 11:29:38 +0000
+++ scripts/checks-and-tests/checks/DEM-PFV-check.py 2014-07-28 06:50:17 +0000
@@ -3,7 +3,7 @@
# the test is based on examples/FluidCouplingPFV/oedometer.py, only slightly simplified and using less particles
-if ('PFVflow' in features):
+if ('PFVFLOW' in features):
errors=0
tolerance=0.01
@@ -21,7 +21,6 @@
sp=pack.SpherePack()
sp.makeCloud(mn,mx,-1,0.3333,num_spheres,False, 0.95,seed=0) #"seed" is not enough for portable determinism it seems, let us use a data file
- checksPath="/home/3S-LAB/bchareyre/yade/yade-git/trunk/scripts/checks-and-tests/checks"
sp.load(checksPath+'/data/100spheres')
sp.toSimulation(material='spheres')
@@ -136,8 +135,8 @@
flow.forceMetis=True
O.run(201,1)
if not flow.metisUsed():
- print "DEM-PFV: Metis is not used during cholmod's reordering although explicitely enabled, something wrong with libraries"
- errors+=1
+ print "DEM-PFV: Metis is not used during cholmod's reordering although explicitly enabled, something wrong with libraries"
+ #errors+=1
if (errors):
resultStatus +=1 #Test is failed
=== modified file 'scripts/checks-and-tests/checks/checkList.py'
--- scripts/checks-and-tests/checks/checkList.py 2013-09-09 07:37:01 +0000
+++ scripts/checks-and-tests/checks/checkList.py 2014-07-28 06:50:17 +0000
@@ -6,7 +6,7 @@
resultStatus = 0
nFailed=0
-skipScripts = ['checkList.py', 'DEM-PFV-check.py']
+skipScripts = ['checkList.py']
for script in scriptsToRun:
if (script[len(script)-3:]==".py" and script not in skipScripts):
