How to add RadialForceEngine

Asked by liucheng83

Hi, all
I want to add RadialForceEngine to one particle with ids=pmove (id=3) to overcome the gravity effect, but when I define the RadialForceEngine with ids=pmove, the particle is gone away, and form inspection, its velocity is nan, why? Here is the code,
anybody can give me some advice, thank you!
-----------------------------------------------------------------
# -*- coding: utf-8 -*-
# encoding: utf-8
from yade import utils, ymport, qt, plot

# basic simulation showing sphere falling ball gravity,
# bouncing against another sphere representing the support

interactionRadius=1.5
iterPeriod_checkInteraction=1000

# DATA COMPONENTS

O.materials.append(FrictMat(young=8e6,poisson=0.3,frictionAngle=radians(20),density=2700,label='spheres'))
# add 2 particles to the simulation
# they the default material (utils.defaultMat)
O.bodies.append([
 # fixed: particle's position in space will not change (support)
 utils.sphere(center=(0,0,0),radius=.5,fixed=True),
 utils.sphere(center=(0,3,0),radius=.5,fixed=True),
 # this particles is free, subject to dynamics
 utils.sphere((0,0,5),.5),
 #utils.sphere((0,3,20),.5)
])

pmove=O.bodies.append([utils.sphere((0,3,5),.5)
])

# FUNCTIONAL COMPONENTS

# simulation loop -- see presentation for the explanation
O.engines=[
 ForceResetter(),
 InsertionSortCollider([Bo1_Sphere_Aabb(aabbEnlargeFactor=interactionRadius,label='aabb')]),
 InteractionLoop(
  [Ig2_Sphere_Sphere_ScGeom(interactionDetectionFactor=interactionRadius,label='Ig2ssGeom')],
  [Ip2_FrictMat_FrictMat_FrictPhys()],
  [Law2_ScGeom_FrictPhys_CundallStrack()]
 ),
 # apply gravity force to particles
 #GravityEngine is deprecated
 #GravityEngine(gravity=(0,0,-9.81)),
 RadialForceEngine(ids=pmove,axisDir=(0,0,1),axisPt=(0,3,5),fNorm=O.bodies[3].state.mass*9.81*0.5,label="radialForce"),
 # damping: numerical dissipation of energy
 NewtonIntegrator(damping=0.1,gravity=(0,0,-9.81)),
 PyRunner(command='checkInteraction()',iterPeriod=iterPeriod_checkInteraction,label='checker')
]

# set timestep to a fraction of the critical timestep
# the fraction is very small, so that the simulation is not too fast
# and the motion can be observed
O.dt=.5e-4*utils.PWaveTimeStep()

#### to see it
v=qt.Controller()
v=qt.View()
rr=qt.Renderer()
rr.intrAllWire=True

# save the simulation, so that it can be reloaded later, for experimentation

def checkInteraction():
 print O.iter
# if O.interactions[3,1].isReal:
# i=O.interactions[3,1]
# print O.iter
# else:
# print O.iter

#run one single step
O.step()
#### initializes now the interaction detection factor
#
#aabb.aabbEnlargeFactor=-1.
#Ig2ssGeom.interactionDetectionFactor=-1.
O.saveTmp()

Question information

Language:
English Edit question
Status:
Solved
For:
Yade Edit question
Assignee:
No assignee Edit question
Solved by:
Anton Gladky
Solved:
Last query:
Last reply:
Revision history for this message
Anton Gladky (gladky-anton) said :
#1

> I want to add RadialForceEngine to one particle with ids=pmove (id=3) to overcome the gravity effect

Hi,

if removing of gravity for 1 particle is the only goal for you, you do
not need to use this engine. Use GravityEngine instead of gravity
parameter for NewtonIntegrator, but apply mask-parameter.

It can disable the gravity for "unmasked" particles. Nans in
RadialForceEngine should be investigated, it can be a bug.

Cheers,

Anton

Revision history for this message
Best Anton Gladky (gladky-anton) said :
#3

Why do you use RadialForceEngine instead of ForceEngine?

https://github.com/yade/trunk/blob/master/pkg/common/ForceEngine.cpp#L36

   Vector3r radial=(pos - (axisPt+axisDir * /* t */
((pos-axisPt).dot(axisDir)))).normalized();

This is, probably, the "blame-line", but it is difficult for me
quickly to understand, how it
works, sorry.

Anton

Revision history for this message
liucheng83 (lcheng83) said :
#4

Hi, Anto . I will try the ForceEngine, a good advice to me, thank you!

Revision history for this message
liucheng83 (lcheng83) said :
#5

Thanks Anton Gladky, that solved my question.