# How to detect the Contact Geometry and Contact Force?

185 views
asked Oct 13, 2016

Hi Vaclav.

As your request, I poste the question here.

So I repeat what I asked in the previous questions in order to let everyone have a better view.

Here I provide this Rough Picture.

My question is this:

1. How can I determine the particle whether it contacts the membrane?

2. After determining the number of particles contacting the membrane, how can I make a sum of the conctact force which comes from the contact between particle and membrane?

3.

sp=woo.pack.SpherePack()
sp.makeCloud((0,0,0.1),(3.4,3.4,1.1),rMean=.3*xmax/xdiv,rRelFuzz=.3,periodic=False,num=10)
sp.toSimulation(S,mat=FrictMat(young=200e6,density=2650,ktDivKn=0.25,tanPhi=0.3))

Here I want to know if the total mass of the spheres = volume*num*density ?

Is there any method I can set the total mass of the particles?

Thanks

Best Regards

Feixiang Xuan

commented Oct 15, 2016 by (770 points)

Add 1 more question:

I can access the node pos. However, I think the Number of the node is random instead of fixed.

For example, the four corner has fixed the translational D.o.F, which means the Pos[2] should be zero. If I type S.dem.nodes[0].pos[2], the displacement is not zero meaning this node is not what we consider at the corner. Besides, I checked the 'Inspector' for the dem.nodes and I found every time I run the simulation, the nodes position is re-arranged again.

So my question is that is there any command which can make the coordinate of the node fixed? Or could you tell me how does the code arrange the coordinate of the node?

Best Regards

Feixiang Xuan

commented Oct 17, 2016 by (770 points)

could you tell me the right way to access the node force? I want the reaction force data.

I tried S.dem.node.force[2] it doesn't work.

Thanks

Feixiang Xuan

answered Oct 17, 2016 by (48,170 points)

Hi, please post unrelated questions separately, otherwise it gets quite messy. Using your question numbers:

1. to access all contacts with the memberane, you can either loop through all membrane elements and then over all its contacts, or over all contacts and filter out only those with membrane:

# option one
for c in S.dem.con:
# skip contacts where neither particles is a Membrane
if not isinstance(c.pA.shape,woo.fem.Membrane) and not isinstance(c.pB.shape,woo.dem.Membrane): continue
# do something, such as sum forces etc.

# option two
for p in S.dem.par:
if not isinstance(p.shape,woo.fem.Membrane): continue
for c in p.con:
# do something with c, which is a contact between Membrane and something

2. I'd guess every time you run the script, the numbering will be the same -- deterministic but arbitrary (maybe not if python sequences come into play). But node numbering is an implementation detail and is not guaranteed to be anything. You should loop over all nodes and find the corner nodes by comparing coordinates -- something like this (make sure you do it just once and store the result e.g. in S.lab, as it traverses all nodes every time called):

from minieigen import *
def findNodesAround(S.dem,coord,eps=1e-3):
# construct box around the coordinate
box=AlignedBox3(Vector3(coord)-eps*Vector3.Ones,Vector3(coord)+eps*Vector3.Ones)
# return all nodes within the box, the caller should get the item (s)he wants
return [n for n in S.dem.nodes if n.pos in box]
S.lab.cornerNodes=[findNodesAround(S.dem,coord)[0] for coord in (0,0,0),(0,1,0),(1,0,0),(0,0,1)] # in this example, corners of unit square in xy-plane

3. To set the total mass, find volume and compute density from required mass and volume (to compute the volume, call SpherePack.solidVolume, obviously after calling makeCloud and before calling toSimulation). If you use more sophisticated generators, you can set the mass you require directly (see tutorial on inlets).

4. Nodal force: Node.dem.force (DemData.force, read the docs, really) so it will be S.dem.nodes[i].dem.force if you start with Scene, i being the node number.

Cheers, v.

commented Oct 19, 2016 by (770 points)

Hi Vaclav
But here for the node numbering (2nd question), I tried your method here it gives me the exact node <node@address coord> instead of the node id (number). For instance S.dem.nodes[i].pos[2], I need the i.
Maybe I didn't either understand your method or you misunderstand my meaning. Anyhow, I solve this question by myself through the following code:

    # Find the interesting nodes id through coordinates
i=0
for i in range(len(S.dem.nodes)):
if (S.dem.nodes[i].pos[0]==0) and (S.dem.nodes[i].pos[1]==0):
print('c_1',i)
c_1=i
elif (S.dem.nodes[i].pos[0]==0) and (S.dem.nodes[i].pos[1]==2.8):
print('c_2',i)
c_2=i
elif (S.dem.nodes[i].pos[0]==2.8) and (S.dem.nodes[i].pos[1]==0):
print('c_3',i)
c_3=i
elif (S.dem.nodes[i].pos[0]==2.8) and (S.dem.nodes[i].pos[1]==2.8):
print('c_4',i)
c_4=i
else: i=i+1


Now I'm going to do the contact loop, before doing it, could you give me some details explanation about the 'PA 'and 'PB'?

Thanks very much
Best Regards

Feixiang Xuan

commented Oct 20, 2016 by (48,170 points)

Hi, that's correct, you are getting the Node instance, which is what you want anyway. You can store it and re-use it later (if you want to keep it across saves, put it in S.lab) as a variable. Nodes don't have IDs, they are numbered, but the number can change anytime (if nodes are deleted, for example), it is just an implementation detail you should not be relying upon; why do you need it anyway? (Particles, on the other hand, do have guaranteed-stable IDs).

For pA and pB see documentation, it is the way to access particles in the contact (they are always 2).

commented Oct 20, 2016 by (770 points)

I need these nodes numbering because I want to get the results of the nodal displacement and nodal reaction forces varying with time instead of a final answer.

Frankly speaking, I need to plot these results in a diagram. So there are several node position that I'm interested in.

Thanks

commented Oct 20, 2016 by (48,170 points)

You don't need to know number of nodes, they have no meaning, it is just and internal detail. You need to get those corner nodes, which I what showed you before. So when you get the nodes, just ask them for displacement. Provided you have S.lab.cornerNodes as I showed you above, you can do something like:

def addPlotData(S):
# do whatever you need
cn=S.lab.cornerNodes # we are lazy
# plot z-positions and z-forces of corner nodes