Welcome to ask.woodem.org. You may post when you login through your GitHub account.

General Questions about modelling rigid bodys with woo

0 votes
84 views
asked Jul 5, 2016 by Clemens1992 (290 points)
reopened Jul 13, 2016 by Clemens1992

Hi everybody,

i'm trying to use woo for the simulation of a rigid body contact problem.

Basically there are five elements in my simulation:

  • 2 passive elements, whose DoF are blocked and they represent the ground
  • 2 active elements, whose DoF are imposed
  • One main element (which is kind of a "bar" - rectangular solid) which should interact with the other elements 

My MWE looks like this:

import woo
from woo.core import *
from woo.dem import *
import woo.gl
import math
from woo.fem import *
from math import *
import numpy as np
from minieigen import *
import csv

#Set the materials
matSteel = FrictMat(density = 7850.,young = 210e6,ktDivKn=0.1, tanPhi=math.tan(20./180*pi))
matAlu   = FrictMat(density = 2700.,young = 70e6,ktDivKn=0.1, tanPhi=math.tan(20./180*pi))
matBar   = matAlu

#Master-Scene 
S=woo.master.scene=Scene(fields=[DemField(gravity=(0,0,-9.81))],dtSafety=.6,dt=1e-5,stopAtTime = 0.2)
S.dt = 1e-5  

#Import STL's, define centralNode for InterpolatedMotion() and addClumped()
ZL  = woo.utils.importSTL('ZL.stl',mat=matSteel, shift = [0.32814+0.005*0.3,-0.035,-0.42562], thickness=.002,scale=1,ori=Quaternion((1,0,0),0),color=-.25)
ZR  = woo.utils.importSTL('ZL.stl',mat=matSteel, shift = [0.32814+0.005*0.3,+0.035,-0.42562],thickness=.002,scale=1,ori=Quaternion((1,0,0),0),color=-.25)
S.lab.ZL = Node(pos=(0,0,0),dem=ClumpData(blocked='xyzXYZ'))
S.lab.ZR = Node(pos=(0,0,0),dem=ClumpData(blocked='xyzXYZ'))

S.dem.par.addClumped(ZL, centralNode  = S.lab.ZL)
S.dem.par.addClumped(ZL,centralNode  = S.lab.ZR)

#Create Bar als triangulatedBox and define mass and inertia within addClumped()
rDim        = Vector3(0.03,0.1,0.015)
rCent       = Vector3(0.31314,0,-0.41312+0.002/2+0.5/1000*1.2)
halfThick   = 0.5/1000. #0.5mm radius an vertices
nodesBar = woo.triangulated.box(dim=rDim,center=rCent,halfThick=halfThick,mat=matBar,fixed=False,wire=False,mask=DemField.defaultMovableMask)
density     = 1e3 
m           = rDim[0]*rDim[1]*rDim[2]*density
I           = m*(1/12.)*Vector3(rDim[1]**2+rDim[2]**2,rDim[0]**2+rDim[2]**2,rDim[1]**2+rDim[2]**2)

bar         = S.dem.par.addClumped(nodesBar,centralNode=Node(pos=rCent,dem=ClumpData(mass=m,inertia=I,blocked='')))

#All Nodes colllected here, so later initialized Nodes will not move
S.dem.collectNodes()

#Import 'Ground' - passive element, should not move, but active elements should collide with it
Ground_Left  = woo.utils.importSTL("GL.STL",mat=matSteel,scale = 1, shift = [0.2015,+0.02,-0.42812900],thickness=.002)
Ground_Right = woo.utils.importSTL("GL.STL",mat=matSteel,scale = 1, shift = [0.2015,-0.02,-0.42812900],thickness=.002)
S.dem.par.addClumped(Ground_Left,centralNode = Node(pos = Vector3(0,0,0), dem = ClumpData(mass = 1, inertia = I, blocked='xyzXYZ')))
S.dem.par.addClumped(Ground_Right,centralNode = Node(pos = Vector3(0,0,0), dem = ClumpData(mass = 1, inertia = I, blocked='xyzXYZ')))
#S.dem.par.add(Ground_Left)
#S.dem.par.add(Ground_Right)

#------------------------------------------------------------------------------
#                       Import Motion for imposed elements
#------------------------------------------------------------------------------
Motion      = np.loadtxt('motion.txt',delimiter=',')
Time        = Motion[:,0]
xyz_Left    = []
xyz_Right   = []
Motion_ori  = [Quaternion.Identity]

for i in range(0,len(Motion)):
    xyz_Left.append(tuple([Motion[i,1]-0.022+0.005/1.2,-0.035,Motion[i,3]+0.0239505+0.005])) 
    xyz_Right.append(tuple([Motion[i,1]-0.022+0.005/1.2,+0.035,Motion[i,3]+0.0229505+0.005])) 
    Motion_ori.append(tuple([tuple([1,0,0]),0]))
   
del Motion_ori[-1] 

#Set motion with the help of InterpolatedMotion() - centralNode is needed
S.lab.ZL.dem.impose = InterpolatedMotion(t0 = 0, poss = xyz_Left, oris=Motion_ori,times=time)
S.lab.ZR.dem.impose = InterpolatedMotion(t0 = 0, poss = xyz_Right,oris=Motion_ori,times=time)

#Init Hertz Model for Contact
modelHertz=woo.models.ContactModelSelector(
    name='Hertz',
    damping=0,
    poisson = 0.3,
    restitution = 0.8)
    
#Simulation-Engine
S.engines=DemField.minimalEngines(model = modelHertz)+[IntraForce([In2_Facet()])]+ [PyRunner(10,'S.plot.autoData()')]
S.saveTmp()

My problem is, that I don't understand some central points of this modelling technique and I hope someone can help me out :)

1.) I created the active element Bar as a triangulated box and added it to the simulation as clumpedData. But I don't understand what this command does internal. The shape of the bar is triangulated, so does addClumped() initialize small spheres on each node of the triangles? Analogously the question remains for the imported STL-files.

2.) I created a central node for the Bar. It is possible to get its position by calling bar.pos, but I get an error if I try to call bar.vel for its velocity. How do I get the velocity (or acceleration) of the bar?

3.) I have some troubles regarding the right init positions of the bodys. So e.g. I would like to put the bar on the ground at the beginning of the simulation. Right now it doesn't work well, because I don't understand in what relation I define the central nodes of the imported STL-files. In my CAD-programme I defined a point of origin of my model, but is this also the imported frame base? In addition to that I think I have to correct the shift-positions of the bodys, because of the thickness I set to the Clumps. Is this consideration right?

4.) My last question relates to the general modelling approach. Is this the right way to simulate rigid body interaction? In most papers it is suggested to build a rigid body out of multiple spheres, which are connected by spring(-damper)-elements. Is this possible with woo? Especially this is an important question, because I would like to quantify the contact-forces.

I hope I could make my problems clear and I hope someone can help me with these questions.

Thanks in advance!

Cheers, Clemens

2 Answers

+1 vote
answered Jul 7, 2016 by eudoxos (44,890 points)
selected Jul 13, 2016 by Clemens1992
 
Best answer

Hi Clements,

thank you for your well-formulated questions. I was able to run your MWE after some adjustments and using STL files you sent by e-mail, though there is some coordinate mismatch so it looks a bit funny, but works:

screenshot animation

1.) I created the active element Bar as a triangulated box and added it to the simulation as clumpedData. But I don't understand what this command does internal. The shape of the bar is triangulated, so does addClumped() initialize small spheres on each node of the triangles? Analogously the question remains for the imported STL-files.

Clumping means that the body behaves as a rigid thing. This is useful if you don't want to do solid deformation calculations (of the internals; they cost lot of time plus have problem on their own) but still have the whole body moving as response to forces (contacts, gravity etc). addClumped tells all the particles and nodes that they a part of this clump and will move with it. Interally, clumping has an effect on motion integration: before computing acceleration, forces are collected from all clumped particles to the central node, its movement is computed, and then all clumped nodes move with the central node (so that the rigid shape is retained).

Clumping is not useful (=will have no effect) for particles with don't move at all - such as the passive ground. If you just block all DOFs, that will do.

2.) I created a central node for the Bar. It is possible to get its position by calling bar.pos, but I get an error if I try to call bar.vel for its velocity. How do I get the velocity (or acceleration) of the bar?

The node is a Node instance and if you check the docs, node has only position and orientation. There is more in DemData attached to it which you access as node.dem. If your node is "n", you get position as n.pos and velocity as n.dem.vel and so on. Acceleration is not stored anywhere, it is computed from force at every step, used to update vel/pos and discarded.

3.) I have some troubles regarding the right init positions of the bodys. So e.g. I would like to put the bar on the ground at the beginning of the simulation. Right now it doesn't work well, because I don't understand in what relation I define the central nodes of the imported STL-files. In my CAD-programme I defined a point of origin of my model, but is this also the imported frame base? In addition to that I think I have to correct the shift-positions of the bodys, because of the thickness I set to the Clumps. Is this consideration right?

Central node (some people like to call it "handle" for the clump) is positioned in the simulation space (global coordinates). woo.utils.importSTL uses global coordinates, but "scale, shift and ori are applied in this order before nodal coordinates are computed". In short, importSTL can transform coords used in the STL file to the simulation space; and central node is positioned in that simulation space.

If you add thickness to the surface (because only facets with non-zero thickness can have meaningful contacts), you need to adjust positions, because thickness will be there in excess to the size of the zero-thickenss surface in the STL file. This is not specific to clumps, thickness is a property of facets (which you use as components of the clump), internally stored in Facet.halfThick.

4.) My last question relates to the general modelling approach. Is this the right way to simulate rigid body interaction? In most papers it is suggested to build a rigid body out of multiple spheres, which are connected by spring(-damper)-elements. Is this possible with woo? Especially this is an important question, because I would like to quantify the contact-forces.

There is not one way to simulate rigid body interactions.

The mass-spring models are one of the early methods, since it is very easy to compute (there is no 3d deformation, just springs). If you want to do that with Woo, go ahead with Sphere and Truss elements (not really tested, but I can make a small example if you need it), and also look at examples/circular-rod.py (that one is non-deforming). The disadvantage is that you won't get much control over type of contact between bodies such as friction etc.

Another method, which is what you use just now, is good when internal deformations of the particles are negligible in comparison to the rigid-body motion. Bodies are represented as surfaces only, and there are only contacts between surfaces.

The last one would be to include internal deformation in the computation; you can do that with Woo (with some drawbacks) if you also triangulate (tetrahedralize) the volume of your solids and use appropriate functors (elastic deformation only, so far at least), see examples/tube.py for a simple example.

For contacts of planar solids (which are geometrically conforming), using hertzian model is not necessary (that one is for non-conforming surfaces, which touch only in one point when undeformed), but it is not worse than any other - there are always inherent simplifications and the level of simplification you want is what guides your model selection.

Hope this helps,

Vaclav

commented Jul 7, 2016 by Clemens1992 (290 points)

Hey Vaclav,

thank you very much for your detailed answer! It's really helpful!

I will try to implement your hints and hopefully I can soon present a working model ;)

Clemens

0 votes
answered Jul 7, 2016 by Clemens1992 (290 points)
edited Jul 7, 2016 by Clemens1992

Hello Vaclav,

I tried some things, but I notice some curious things about the imposed Motion:

1.) At t = 0 my active element, whose motion is imposed, has the right orientation. Two timesteps further, the orientation changes obviusly, but S.lab.ZL.ori, still says Quaternion((1,0,0),0). I add two screenshots of this, to illustrate this. 

Orientation at t = 0

Init Orientation

Two timesteps away from t = 0

Imposed Motion

(The Position jumps, but thats still a problem of my initial shfit position).

Do you know, why the orientation changes, allthough I declared [1,0,0]),0] as  Quaternion to the imposed Motion?

2.) I tried to plot the velocity of the active element and the velocity of the bar. There two curious effects:

The velocity of the active element looks like this:

Theoretically the velocity should be a smooth function. So where do these steps come from? In imposedMotion() I only set the position - so are the discontinuities consequences of the numerical calculation of the derivation? 

Or do the collisions with the bar have an effect to the moved element? If this is true, how can I avoid this behavior?

The second curious thing is, that I can write "bar.dem.vel[0]" into the console and get a result. But if I try to plot it, I get an error message: 

"Traceback (most recent call last):   File "c:\src\woo\build\pyi.win32\woo.pyinstaller\out00-PYZ.pyz\woo.plot", line 190, in colDictUpdate  File "<string>", line 1, in <module>
NameError: name 'bar' is not defined
WARN: ignoring exception raised while evaluating auto-column ` - bar.dem.vel[0]' ($v_{Riegel}$ )."

What did I do wrong?

3.) My last question is less technically: Is there a way to change the perspective? Right now the 3D output displays everything with a vanishing point - is it possible to change the display to parallel projection?

Thanks in advance!

Clemens

commented Jul 12, 2016 by eudoxos (44,890 points)

Minor points (I will get to the other ones later):

  • Perspective: https://woodem.org/user/visualization.html#integrated-3d-view (short: press t to switch between perspective/orthographic camera)
  • NameError: the expression to plot is evaluated in different context than the source file (which is the same as the console), so the "bar" variable does not exist there; you have to access it through the Scene object (defined as S); so first put it in S.lab (as you do with S.lab.ZLinks and so on) and then say S.lab.bar.dem.vel[0].

HTH, v.

commented Jul 13, 2016 by Clemens1992 (290 points)

Thank you Vaclav!

Luckily I already solved the other described problems.

The changing orientation of my imposed element at the beginning of the simulation, seems to relate to my imposed motion. If the starting position isn't the same as the first imposed point, then the element jumps to this point and changes its orientation (thats my explanation).

I solved this (and the problem with the non-smooth velocity) by not imposing the motion with the help of a xyz-vector. Instead I used VariableVelocity() and it works perfectly.

Right now I'm struggling with the friciton law implemented in woo. Should I open a new thread or do you think it relates to this topic?

Cheers, Clemens

commented Jul 13, 2016 by eudoxos (44,890 points)

Hi Clemens, glad that you solved it. I don't think change in position can change orientation (those are orthogonal, and also handled separately in the code).

Open a new question for the frictional law, I am sure there will be a lots to write. And mark this question as answered. Thanks.

...