← Back to team overview

yade-users team mailing list archive

Re: [Question #248996]: Viscoelastic Force Model

 

Question #248996 on Yade changed:
https://answers.launchpad.net/yade/+question/248996

Anton Gladky proposed the following answer:
Sorry Andy, I do not get your idea. There is an updated script, which
lets you to analyze all parameters after the second impact e.g. between
interactions 14400-14800.

If you look at the blue line (normal force), it is always 0, when
there is no penetration.
And we never get an attractive force during penetration, what is
restricted by lines
105-106 in ViscoelasticPM.cpp. From my point of view, all is right.

> its posible to see that the force is actually acting over the particle, because the velocity
> also change before the particles make contact.

Sure, you have a gravity there.

Anton

================================================

#!/usr/bin/env python
# -*- coding: utf-8 -*-
from yade import pack, plot
import numpy as np

Corr='-test'

Arch1='energ'+Corr+'.out'
deltaT=5e-5
grav=-9.81
zi=10.0
sumado=Vector3(0.0,0.0,.0)
# "PARTÍCULAS"
idmatVisEl=O.materials.append(ViscElMat(kn=1e9,cn=300e4,ks=1.4e4,cs=200e4,frictionAngle=0.36,label='matVisEl'))

EsfMo = O.bodies.append(utils.sphere(center=(0,0,2.0),material='matVisEl',
radius=1.0,fixed=False))
EsFix = O.bodies.append(utils.sphere(center=(0,0,0),
material='matVisEl',radius=1.0,fixed=True))

O.bodies[EsfMo].state.vel=Vector3(0,0,-12.528367)

O.engines=[
   ForceResetter(),
   InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb()]),
   InteractionLoop(
      [Ig2_Sphere_Sphere_ScGeom(),Ig2_Facet_Sphere_ScGeom()],
      [Ip2_ViscElMat_ViscElMat_ViscElPhys()],
      [Law2_ScGeom_ViscElPhys_Basic()]
   ),
        NewtonIntegrator(gravity=(0,0,grav),damping=0.0),
        PyRunner(command='dissip(sumado)',iterPeriod=1,label='potDisip'),#Cuando
comienzo a calcular las interacciones
        PyRunner(command='addPlotData()',iterPeriod=1,dead=True,
label='PlotAdd'),
]

from yade import qt

qt.View()


# Function to add data to plot
def addPlotData():
  posZ = O.bodies[EsfMo].state.pos[2]
  velZ = O.bodies[EsfMo].state.vel[2]
  PDcalc = 2.0 - posZ
  PDgeom = 0.0

  try:
    PDgeom = O.interactions[0,1].geom.penetrationDepth
  except:
    pass

  Fn = 0.0
  try:
    Fn = O.interactions[0,1].phys.normalForce[2]
  except:
    pass

  FnCalculated = 5e8*PDgeom - 1.5e6*velZ

  plot.addData(it1=O.iter, PDgeom = PDgeom, velZ = velZ/100.0, FnYade
=  -Fn/1e8, FnCalculated = FnCalculated/1e8)


plot.plots={'it1':('PDgeom','velZ','FnYade','FnCalculated')}; plot.plot()

def dissip(sumado):
        g = open (Arch1,'a')
        v=O.bodies[0].state.vel[2]
        z=O.bodies[0].state.pos[2]
        m=O.bodies[0].state.mass
        Fn=0.0
        for i in O.interactions:
                PD=i.geom.penetrationDepth
                YFn=i.phys.normalForce
                Fn=YFn[2]

g.write(str(O.iter)+'\t'+str(z)+'\t'+str(v)+'\t'+'\t'+str(Fn)+'\t'+str(PD)+'\n')
        g.closed
O.dt=deltaT

O.run(14400, True)
O.wait
PlotAdd.dead=False
O.run(400, True)
O.wait
plot.saveGnuplot('sim-data_Sphere')

-- 
You received this question notification because you are a member of
yade-users, which is an answer contact for Yade.