yade-dev team mailing list archive
-
yade-dev team
-
Mailing list archive
-
Message #03604
[Branch ~yade-dev/yade/trunk] Rev 2082: - Suppression of CinemCisEngine, which even me almost never used.
------------------------------------------------------------
revno: 2082
committer: jduriez <jduriez@c1solimara-l>
branch nick: trunk
timestamp: Tue 2010-03-16 17:02:55 +0100
message:
- Suppression of CinemCisEngine, which even me almost never used.
- Add of KinemCTDEngine: performs constant tangential displacement (=oedometrical) compressions on the simple shear box
- move of KinemCND from pkg/common to pkg/dem. This Engine does not use really any DEM-specific concept, but I guess it will always be used in DEM simulations, so I thought
it will be more adequate in dem/
- modify of corresponding "include"...
- rewrite of KinemCNDEngine.hpp, with "good" macros
removed:
pkg/common/Engine/PartialEngine/CinemCisEngine.cpp
pkg/common/Engine/PartialEngine/CinemCisEngine.hpp
added:
pkg/dem/Engine/PartialEngine/KinemCTDEngine.cpp
pkg/dem/Engine/PartialEngine/KinemCTDEngine.hpp
renamed:
pkg/common/Engine/PartialEngine/KinemCNDEngine.cpp => pkg/dem/Engine/PartialEngine/KinemCNDEngine.cpp
pkg/common/Engine/PartialEngine/KinemCNDEngine.hpp => pkg/dem/Engine/PartialEngine/KinemCNDEngine.hpp
modified:
pkg/dem/Engine/PartialEngine/Disp2DPropLoadEngine.hpp
pkg/dem/PreProcessor/DirectShearCis.cpp
pkg/dem/PreProcessor/SimpleShear.cpp
py/system.py
pkg/dem/Engine/PartialEngine/KinemCNDEngine.hpp
--
lp:yade
https://code.launchpad.net/~yade-dev/yade/trunk
Your team Yade developers is subscribed to branch lp:yade.
To unsubscribe from this branch go to https://code.launchpad.net/~yade-dev/yade/trunk/+edit-subscription.
=== removed file 'pkg/common/Engine/PartialEngine/CinemCisEngine.cpp'
--- pkg/common/Engine/PartialEngine/CinemCisEngine.cpp 2010-02-09 20:22:04 +0000
+++ pkg/common/Engine/PartialEngine/CinemCisEngine.cpp 1970-01-01 00:00:00 +0000
@@ -1,104 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2008 by Jerome Duriez *
-* jerome.duriez@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. *
-*************************************************************************/
-
-#include<yade/pkg-common/RigidBodyParameters.hpp>
-#include<yade/core/Scene.hpp>
-#include<yade/lib-base/Math.hpp>
-
-#include"CinemCisEngine.hpp"
-
-
-CinemCisEngine::CinemCisEngine()
-{
- rotationAxis=Vector3r(0,0,1);
- shearSpeed=0;
- theta=0;
- thetalim=9.0/10.0*Mathr::PI/2.0;
- id_boxhaut=3;
-}
-
-
-
-
-void CinemCisEngine::applyCondition(Scene * body)
-{
- if(theta<=thetalim)
- {
- applyRotTranslation(body); // to let move the lateral walls
- applyTranslation(body); // to let move the upper wall
- }
-
-}
-
-
-
-void CinemCisEngine::applyRotTranslation(Scene * ncb)
-{
- shared_ptr<BodyContainer> bodies = ncb->bodies;
-
- Yplaqsup=((*bodies)[id_boxhaut]->physicalParameters.get())->se3.position.Y(); // the height of the sample, which may be different from the initial "height" defined in the Preprocessor
-
- subscribedBodies.clear();
- subscribedBodies.push_back(0); // 0 and 2 are in my case the ids of the lateral walls
- subscribedBodies.push_back(2);
- std::vector<int>::const_iterator ii = subscribedBodies.begin();
- std::vector<int>::const_iterator iiEnd = subscribedBodies.end();
-
- Real dt = Omega::instance().getTimeStep();
- Real dx = shearSpeed*dt;
- Real dtheta;
- dtheta=Mathr::ATan(dx/(Yplaqsup+dx*Mathr::Tan(theta)+Yplaqsup*Mathr::Pow(Mathr::Tan(theta),2)));
- Quaternionr q;
- q.FromAxisAngle(rotationAxis,-dtheta);
-
- for(;ii!=iiEnd;++ii)
- {
- RigidBodyParameters * rb = dynamic_cast<RigidBodyParameters*>((*bodies)[*ii]->physicalParameters.get());
-
-// La partie de rotation :
- rb->se3.orientation = q*rb->se3.orientation;
- rb->se3.orientation.Normalize();
-/* rb->angularVelocity = rotationAxis*angularVelocity;*/
- rb->velocity = Vector3r(0,0,0);
-// La partie de translation :
- Real YcentreW1; //the altitude of the center of both lateral walls
- YcentreW1=((*bodies)[0]->physicalParameters.get())->se3.position.Y();
- rb->se3.position+=dt*(shearSpeed*YcentreW1/Yplaqsup)*Vector3r(1,0,0);
- }
- theta+=dtheta;
-}
-
-
-void CinemCisEngine::applyTranslation(Scene * ncb)
-{
- shared_ptr<BodyContainer>& bodies = ncb->bodies;
-
- int id_boxhaut=3;
- Real dt = Omega::instance().getTimeStep();
-
- if(ParticleParameters* p = dynamic_cast<ParticleParameters*>((*bodies)[id_boxhaut]->physicalParameters.get()))
- {
- p->se3.position += dt*shearSpeed*Vector3r(1,0,0); // by def of shearSpeed
- p->velocity = shearSpeed*Vector3r(1,0,0);
- }
- else if(PhysicalParameters* b = dynamic_cast<PhysicalParameters*>((*bodies)[id_boxhaut]->physicalParameters.get()))
- { // NOT everyone has velocity !
- b->se3.position += dt*shearSpeed*Vector3r(1,0,0);
- }
- else
- {
- std::cerr << "TranslationEngine::applyCondition, WARNING! dynamic_cast failed! for id: " << id_boxhaut << std::endl;
- }
-
-}
-
-
-
-
-YADE_REQUIRE_FEATURE(PHYSPAR);
-
=== removed file 'pkg/common/Engine/PartialEngine/CinemCisEngine.hpp'
--- pkg/common/Engine/PartialEngine/CinemCisEngine.hpp 2010-02-09 20:22:04 +0000
+++ pkg/common/Engine/PartialEngine/CinemCisEngine.hpp 1970-01-01 00:00:00 +0000
@@ -1,45 +0,0 @@
-/*************************************************************************
-* Copyright (C) 2008 by Jerome Duriez *
-* jerome.duriez@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/core/Omega.hpp>
-#include<yade/core/PartialEngine.hpp>
-#include<yade/lib-base/Math.hpp>
-
-/*! \brief To apply a zero normal displacement shear for a parallelogram box
-
-This engine, used in simulations issued from "DirectShearCis" Preprocessor, allows to translate horizontally the upper plate while the lateral ones rotate so that they always keep contact with the lower and upper walls
-*/
-
-
-class CinemCisEngine : public PartialEngine
-{
- private :
- Real Yplaqsup// height of the upper wall, the engine cares to find its value
- ,theta; // the angle between a lateral plate and its original orienation : will increase from 0 to thetalim
-
- public :
- CinemCisEngine();
- void applyCondition(Scene * body);
-
- Real shearSpeed // to be defined in the PreProcessor
- ,thetalim // the maximum value of theta, at wich the displacement is stopped
- ;
- body_id_t id_boxhaut; // the id of the upper wall : defined in the constructor
- Vector3r rotationAxis; // defined in the constructor
- void applyRotTranslation(Scene *); // to let move (rotation combined with translation) the lateral walls
- void applyTranslation(Scene *); // to let move (translation) the upper wall
- REGISTER_ATTRIBUTES(PartialEngine,(shearSpeed)(rotationAxis)(theta)(thetalim)(id_boxhaut));
- REGISTER_CLASS_NAME(CinemCisEngine);
- REGISTER_BASE_CLASS_NAME(PartialEngine);
-};
-
-REGISTER_SERIALIZABLE(CinemCisEngine);
-
-
=== modified file 'pkg/dem/Engine/PartialEngine/Disp2DPropLoadEngine.hpp'
--- pkg/dem/Engine/PartialEngine/Disp2DPropLoadEngine.hpp 2010-03-16 12:36:54 +0000
+++ pkg/dem/Engine/PartialEngine/Disp2DPropLoadEngine.hpp 2010-03-16 16:02:55 +0000
@@ -6,8 +6,7 @@
* GNU General Public License v2 or later. See file LICENSE for details. *
*************************************************************************/
-#ifndef DIRECRESEARCHENGINE_HPP
-#define DIRECRESEARCHENGINE_HPP
+#pragma once
#include<yade/core/Omega.hpp>
#include<yade/core/PartialEngine.hpp>
@@ -79,5 +78,3 @@
REGISTER_SERIALIZABLE(Disp2DPropLoadEngine);
-#endif
-
=== renamed file 'pkg/common/Engine/PartialEngine/KinemCNDEngine.cpp' => 'pkg/dem/Engine/PartialEngine/KinemCNDEngine.cpp'
=== renamed file 'pkg/common/Engine/PartialEngine/KinemCNDEngine.hpp' => 'pkg/dem/Engine/PartialEngine/KinemCNDEngine.hpp'
--- pkg/common/Engine/PartialEngine/KinemCNDEngine.hpp 2010-03-10 11:37:28 +0000
+++ pkg/dem/Engine/PartialEngine/KinemCNDEngine.hpp 2010-03-16 16:02:55 +0000
@@ -52,6 +52,7 @@
string Key;
protected :
+// YADE_CLASS_BASE_DOCS_ATTRS(KinemCNDEngine,PartialEngine
REGISTER_ATTRIBUTES(PartialEngine,(shearSpeed)(gammalim)(gamma)(gamma_save)(temoin_save)(id_boxhaut)(id_boxleft)(id_boxright)(Key));
REGISTER_CLASS_NAME(KinemCNDEngine);
REGISTER_BASE_CLASS_NAME(PartialEngine);
=== added file 'pkg/dem/Engine/PartialEngine/KinemCTDEngine.cpp'
--- pkg/dem/Engine/PartialEngine/KinemCTDEngine.cpp 1970-01-01 00:00:00 +0000
+++ pkg/dem/Engine/PartialEngine/KinemCTDEngine.cpp 2010-03-16 16:02:55 +0000
@@ -0,0 +1,165 @@
+/*************************************************************************
+* Copyright (C) 2008 by Jerome Duriez *
+* jerome.duriez@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. *
+*************************************************************************/
+
+#include"KinemCTDEngine.hpp"
+#include<yade/core/State.hpp>
+#include<yade/pkg-common/Box.hpp>
+#include<yade/core/Scene.hpp>
+
+YADE_PLUGIN((KinemCTDEngine));
+
+KinemCTDEngine::~KinemCTDEngine()
+{
+}
+
+void KinemCTDEngine::applyCondition(Scene * ncb)
+{
+
+ leftbox = Body::byId(id_boxleft);
+ rightbox = Body::byId(id_boxright);
+ frontbox = Body::byId(id_boxfront);
+ backbox = Body::byId(id_boxback);
+ topbox = Body::byId(id_topbox);
+ boxbas = Body::byId(id_boxbas);
+
+// Computation of the current dimensions of the box : //
+ Real Xleft = leftbox->state->pos.X() + (YADE_CAST<Box*>(leftbox->shape.get()))->extents.X();
+
+ Real Xright = rightbox->state->pos.X() - (YADE_CAST<Box*>(rightbox->shape.get()))->extents.X();
+
+ Real Zfront = frontbox->state->pos.Z() - YADE_CAST<Box*>(frontbox->shape.get())->extents.Z();
+ Real Zback = backbox->state->pos.Z() + (YADE_CAST<Box*>(backbox->shape.get()))->extents.Z();
+
+ Real Scontact = (Xright-Xleft)*(Zfront-Zback); // that's so the value of section at the middle of the height of the box
+// End of computation of the current dimensions of the box //
+
+ ncb->forces.sync(); Vector3r F_sup=ncb->forces.getForce(id_topbox);
+ Real current_NormalForce=(ncb->forces.getForce(id_topbox)).Y();
+ current_sigma=current_NormalForce/(1000.0*Scontact); // so we have the current value of sigma, in kPa
+
+ if( ((compSpeed > 0) && (current_sigma < target_sigma)) || ((compSpeed < 0) && (current_sigma > target_sigma)) )
+ {
+ if(temoin!=0)
+ {
+// cout << "j'ai ici un temoin #0 visiblement. En effet temoin =" <<lexical_cast<string>(temoin) << endl;
+ temoin=0;
+// cout << "Maintenant (toujours dans le if temoin!=0), temoin =" <<lexical_cast<string>(temoin) << endl;
+ }
+ letMove(ncb);
+ }
+ else if (temoin==0)
+ {
+ stopMovement();
+// cout << "Mouvement stoppe, temoin = " << lexical_cast<string>(temoin) << endl;
+// cout << " Dans le if, temoin =" << lexical_cast<string>(temoin) << endl;
+ string f;
+ if (compSpeed > 0)
+ f="sigmax_";
+ else
+ f="Sigmin_";
+
+ Omega::instance().saveSimulation(Key + f +lexical_cast<string> (floor(target_sigma)) + "kPaReached.xml");
+ temoin=1;
+// cout << " Fin du if, temoin =" << lexical_cast<string>(temoin) << endl << endl;
+ }
+
+
+ for(unsigned int j=0;j<sigma_save.size();j++)
+ {
+ if( ( ( (compSpeed>0)&&(current_sigma > sigma_save[j]) ) || ((compSpeed<0)&&(current_sigma < sigma_save[j])) ) && (temoin_save[j]==0))
+ {
+ stopMovement();
+ Omega::instance().saveSimulation(Key + "SigInt_" +lexical_cast<string> (floor(current_sigma)) + "kPareached.xml");
+ temoin_save[j]=1;
+ }
+ }
+// cout << "Fin de ApplyCondi, temoin = " << lexical_cast<string>(temoin) << endl;
+
+}
+
+void KinemCTDEngine::letMove(Scene * ncb)
+{
+ computeAlpha();
+ Real dt = Omega::instance().getTimeStep()
+ ,dh=-compSpeed*dt
+ ,dgamma=0
+ ;
+
+ Real Ysup = topbox->state->pos.Y();
+ Real Ylat = leftbox->state->pos.Y();
+
+// Changes in vertical and horizontal position :
+
+ topbox->state->pos += Vector3r(dgamma,dh,0);
+
+ leftbox->state->pos += Vector3r(dgamma/2.0,dh/2.0,0);
+ rightbox->state->pos += Vector3r(dgamma/2.0,dh/2.0,0);
+
+ Real Ysup_mod = topbox->state->pos.Y();
+ Real Ylat_mod = leftbox->state->pos.Y();
+
+// with the corresponding velocities :
+ topbox->state->vel = Vector3r(dgamma/dt,dh/dt,0);
+ leftbox->state->vel = Vector3r((dgamma/2.0)/dt,dh/(2.0*dt),0);
+ rightbox->state->vel = Vector3r((dgamma/2.0)/dt,dh/(2.0*dt),0);
+
+// Then computation of the angle of the rotation to be done :
+ if (alpha == Mathr::PI/2.0) // Case of the very beginning
+ {
+ dalpha = - Mathr::ATan( dgamma / (Ysup_mod -Ylat_mod) );
+ }
+ else
+ {
+ Real A = (Ysup_mod - Ylat_mod) * 2.0*Mathr::Tan(alpha) / (2.0*(Ysup - Ylat) + dgamma*Mathr::Tan(alpha) );
+ dalpha = Mathr::ATan( (A - Mathr::Tan(alpha))/(1.0 + A * Mathr::Tan(alpha)));
+ }
+
+ Quaternionr qcorr;
+ qcorr.FromAxisAngle(Vector3r(0,0,1),dalpha);
+
+// On applique la rotation en changeant l'orientation des plaques et leurs vang
+ leftbox->state->ori = qcorr*leftbox->state->ori;
+ leftbox->state->angVel = Vector3r(0,0,1)*dalpha/dt;
+
+ rightbox->state->ori = qcorr*leftbox->state->ori;
+ rightbox->state->angVel = Vector3r(0,0,1)*dalpha/dt;
+
+}
+
+void KinemCTDEngine::computeAlpha()
+{
+ Quaternionr orientationLeftBox,orientationRightBox;
+ orientationLeftBox = leftbox->state->ori;
+ orientationRightBox = rightbox->state->ori;
+ if(orientationLeftBox!=orientationRightBox)
+ {
+ cout << "WARNING !!! your lateral boxes have not the same orientation, you're not in the case of a box imagined for creating these engines" << endl;
+ }
+ Vector3r axis;
+ Real angle;
+ orientationLeftBox.ToAxisAngle(axis,angle);
+ alpha=Mathr::PI/2.0-angle; // right if the initial orientation of the body (on the beginning of the simulation) is q =(1,0,0,0) = FromAxisAngle((0,0,1),0)
+}
+
+void KinemCTDEngine::stopMovement()
+{
+ // annulation de la vitesse de la plaque du haut
+ topbox->state->vel = Vector3r(0,0,0);
+
+ // de la plaque gauche
+ leftbox->state->vel = Vector3r(0,0,0);
+ leftbox->state->angVel = Vector3r(0,0,0);
+
+ // de la plaque droite
+ rightbox->state->vel = Vector3r(0,0,0);
+ rightbox->state->angVel = Vector3r(0,0,0);
+}
+
+
+
+
=== added file 'pkg/dem/Engine/PartialEngine/KinemCTDEngine.hpp'
--- pkg/dem/Engine/PartialEngine/KinemCTDEngine.hpp 1970-01-01 00:00:00 +0000
+++ pkg/dem/Engine/PartialEngine/KinemCTDEngine.hpp 2010-03-16 16:02:55 +0000
@@ -0,0 +1,64 @@
+/*************************************************************************
+* Copyright (C) 2008 by Jerome Duriez *
+* jerome.duriez@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/core/PartialEngine.hpp>
+#include<yade/core/Body.hpp>
+#include<Wm3Vector3.h>
+
+
+class KinemCTDEngine : public PartialEngine
+{
+ private :
+ Real alpha // angle from the lower plate to the left box (trigo wise), the Engine finds itself its value
+ ,dalpha // the increment over alpha to set, to cancel the offset due to vertical displacement of upper box
+ ,current_sigma // Computed in kPa
+ ;
+ void stopMovement(); // to cancel all the (rotation) velocities when Flim is reached
+
+ int temoin;
+
+ shared_ptr<Body> leftbox;
+ shared_ptr<Body> rightbox;
+ shared_ptr<Body> frontbox;
+ shared_ptr<Body> backbox;
+ shared_ptr<Body> boxbas;
+ shared_ptr<Body> topbox;
+
+ void letMove(Scene*);
+
+ public :
+ virtual ~KinemCTDEngine();
+
+ void applyCondition(Scene *)
+ ,computeAlpha()
+ ;
+
+ YADE_CLASS_BASE_DOC_ATTRS_CTOR(KinemCTDEngine,PartialEngine,
+ "To compress a simple shear sample by moving the upper box in a vertical way only, so that the tangential displacement (defined by the horizontal gap between the upper and lower boxes) remains constant.\n \t The lateral boxes move also to keep always contact. All that until this box is submitted to a given stress (=*target_sigma*). Moreover saves are executed at each value of stresses stored in the vector *sigma_save*, and at *target_sigma*",
+ ((Real,compSpeed,0.0,"(vertical) speed of the upper box : >0 for real compression, <0 for unloading [m/s]"))
+ ((std::vector<Real>,sigma_save,,"vector with the values of sigma at which a save of the simulation should be performed [kPa]"))
+ ((std::vector<Real>,temoin_save,,"vector (same length as 'sigma_save'), with 0 or 1 depending whether the save for the corresponding value of gamma has been done (1) or not (0)"))
+ ((Real,target_sigma,0.0,"the value of sigma at which the compression should stop [kPa]"))
+ ((body_id_t,id_topbox,3,"the id of the upper wall"))
+ ((body_id_t,id_boxbas,1,"the id of the lower wall"))
+ ((body_id_t,id_boxleft,0,"the id of the left wall"))
+ ((body_id_t,id_boxright,2,"the id of the right wall"))
+ ((body_id_t,id_boxfront,5,"the id of the wall in front of the sample"))
+ ((body_id_t,id_boxback,4,"the id of the wall at the back of the sample"))
+ ((string,Key,"","string to add at the names of the saved files")),
+ temoin=0;
+ )
+
+
+};
+
+REGISTER_SERIALIZABLE(KinemCTDEngine);
+
+
=== modified file 'pkg/dem/PreProcessor/DirectShearCis.cpp'
--- pkg/dem/PreProcessor/DirectShearCis.cpp 2010-03-10 11:37:28 +0000
+++ pkg/dem/PreProcessor/DirectShearCis.cpp 2010-03-16 16:02:55 +0000
@@ -42,7 +42,7 @@
#include<yade/pkg-common/CundallNonViscousDamping.hpp>
#include<yade/pkg-common/CundallNonViscousDamping.hpp>
#include<yade/pkg-common/GravityEngines.hpp>
-#include<yade/pkg-common/KinemCNDEngine.hpp>
+#include<yade/pkg-dem/KinemCNDEngine.hpp>
#include<yade/pkg-common/InteractionGeometryDispatcher.hpp>
=== modified file 'pkg/dem/PreProcessor/SimpleShear.cpp'
--- pkg/dem/PreProcessor/SimpleShear.cpp 2010-03-10 11:37:28 +0000
+++ pkg/dem/PreProcessor/SimpleShear.cpp 2010-03-16 16:02:55 +0000
@@ -37,7 +37,7 @@
#include<yade/pkg-dem/NewtonIntegrator.hpp>
#include<yade/pkg-common/GravityEngines.hpp>
-#include<yade/pkg-common/KinemCNDEngine.hpp>
+#include<yade/pkg-dem/KinemCNDEngine.hpp>
#include<yade/pkg-common/InteractionGeometryDispatcher.hpp>
=== modified file 'py/system.py'
--- py/system.py 2010-03-15 14:35:44 +0000
+++ py/system.py 2010-03-16 16:02:55 +0000
@@ -125,6 +125,7 @@
'CinemDNCEngine':'KinemCNDEngine', # Wed Mar 10 12:34:27 2010, jduriez@c1solimara-l
'CinemDTCEngine':'KinemCTDEngine', # Wed Mar 10 12:34:37 2010, jduriez@c1solimara-l
'Ip2_BMP_BMP_CSPhys':'Ip2_2xFrictMat_CSPhys', # Wed Mar 10 15:08:56 2010, eudoxos@frigo
+ 'CinemDTCEngine':'KinemCTDEngine', # Tue Mar 16 13:54:21 2010, jduriez@c1solimara-l
### END_RENAMED_CLASSES_LIST ### (do not delete this line; scripts/rename-class.py uses it
}
