← Back to team overview

yade-users team mailing list archive

Re: Coupling Yade and Comsol


Hello Andrea

Here are few more details on how to account for the fluid forces in Yade 
(I had to conclude quickly yesterday).
Currently, each type of simulation in Yade is defined in a PreProcessor 

If you look in e.g. TriaxialTest.cpp, you will see a function called 
TriaxialTest::createActors(shared_ptr<MetaBody>& rootBody);
"rootBody" is the simulation.
In this function you see a list of "rootBody->engines.push_back(XXX)". 
This list defines what has to be done at each timestep.
If you write the HydroForceEngine as described in my previous email, you 
will include it using rootBody->engines.push_back(HydroForceEngine). Of 
course you can write a test in the engine itself so that it do something 
every N timestep, rather than each timestep.

You can write a second engine - or even use the same one - to write 
(every N timestep?) the positions of the grains in the file that you 
will load in Comsol.

Yade will read the list of forces, run a number of steps, and write the 
new coordinates to a file. Next question is when and how to remove 
grains. You have to define the criterion that will define what a "free 
grain" is. Perhaps based on velocity or number of contacts? The easiest 
way to implement this will be to do nothing. In yade, let the removed 
grains be free, they will just fall appart. In comsol, just ignore them, 
so that they don't affect the flow  (this is just for preliminary tests 

I paste below some comments from Janek, they may be usefull to 
understand what is happening in ForceEngine.cpp.



The comments are refering to this line :

*static_cast<Force*>( ncb->actionParameters->find( b->getId() , actionParameterForce->getClassIndex() ).get() )->force += force;*

/yes, in fact this whole line is too long, and it should be less
complicated. I agree that it must be difficult to comprehend.

I should find some easier way to do the same thing, perhaps some macro..
I'm just not sure how to make it simpler. So until I'll find something -
this is the only way to access PhysicalAction (for example a class
Force, which derives from PhysicalAction)

[1]  static_cast<Force*>
[2]  (
[3]    ncb->actionParameters->find
[4]   ( 
[5]       b->getId() , 
[6]       actionParameterForce->getClassIndex() 
[7]   )
[8]   .get()
[9]  )
[10] ->force 
[11] += force;

[1] static_cast is used to switch between classes which are in the same
inheritance tree. Here the result of static_cast<Force*>(...) is a
pointer to class Force. for example we could write:

PhysicalAction* p =ncb->actionParameters->find(b->getId(),actionParameterForce->getClassIndex()).get();
Force* f = static_cast<Force*>(p);
f->force +=force;

[2] - [9] between the brackets is the single argument for the static_cast  :) 

[3] retreiving that single argument is little complicated: we find() in
actionParameters a class PhysicalAction that:

[5] corresponds to "the b->getID()" which is id number of a Body. so it
is a PhysicalAction deticated for that Body

[6] and it is also a Force. actionParameterForce is of type
shared_ptr<Force> and we use it ONLY to get the index number of Force.
For example Force has number 1, and Momentum has number 2.

Similar functionality would be provided in following way:

find( b->getId() , "Force" )

it means that you want a PhysicalAction "Force" for this Body. The point
here is that using a string "Force" is tad slower than using a single
number ( getClassIndex() ), which we know that it means "Force".

[4] - [7] arguments for find()

[8] find() actually returns type boost::shared_ptr<PhysicalAction>. That
shared_ptr is a so-called shared pointer. It is an improved version of a
pointer (like for example void* ). Yade is using shared pointers almost
everywhere bacause they are safer to deal with (no memory leaks  ;) 

so that .get() here is used to extract a normal pointer PhysicalAction*
from shared_ptr<PhysicalAction>

[9] so after all that we have a PhysicalAction*, about which we know -
is of type "Force" because we used getClassIndex() to tell that we
want a "Force" not "Momentum", and not "Stiffness", and now
static_cast<Force*> changes the type of PhysicalAction* to Force*

[10] so if now we have a pointer to class Force* - we can access its member -> force

[11] and do something with it/

Chareyre Bruno
Maître de conférence

Institut National Polytechnique de Grenoble
Laboratoire 3S (Soils Solids Structures) - bureau I08
BP 53 - 38041, Grenoble cedex 9 - France     
Tél :

Yade-users mailing list