← Back to team overview

yade-dev team mailing list archive

[Branch ~yade-pkg/yade/git-trunk] Rev 3762: Export normal and viscous part of the visco-elastic contact

 

------------------------------------------------------------
revno: 3762
committer: Anton Gladky <gladky.anton@xxxxxxxxx>
timestamp: Wed 2015-12-16 07:39:37 +0100
message:
  Export normal and viscous part of the visco-elastic contact
added:
  examples/simple-scene/2SpheresNormVisc.py
modified:
  pkg/dem/ViscoelasticPM.cpp
  pkg/dem/ViscoelasticPM.hpp


--
lp:yade
https://code.launchpad.net/~yade-pkg/yade/git-trunk

Your team Yade developers is subscribed to branch lp:yade.
To unsubscribe from this branch go to https://code.launchpad.net/~yade-pkg/yade/git-trunk/+edit-subscription
=== added file 'examples/simple-scene/2SpheresNormVisc.py'
--- examples/simple-scene/2SpheresNormVisc.py	1970-01-01 00:00:00 +0000
+++ examples/simple-scene/2SpheresNormVisc.py	2015-12-16 06:39:37 +0000
@@ -0,0 +1,53 @@
+#!/usr/bin/env python
+# encoding: utf-8
+#
+# This script plots the viscous and normal forces of the contact
+# of 2 spheres
+
+from yade import utils, plot, qt
+o = Omega()
+
+# Physical parameters
+fr = 0.5;rho=2000
+tc = 0.0001; en = 0.7; et = 0.7;
+o.dt = 0.000002*tc
+Rad = 2.0e-3
+
+# Add material
+mat1 = O.materials.append(ViscElMat(frictionAngle=fr,tc=tc,en=en,et=et))
+
+# Add spheres
+id1 = O.bodies.append(sphere(center=[0,0,0],radius=Rad,material=mat1,fixed=True))
+id2 = O.bodies.append(sphere(center=[0,0,(2.0*Rad)],radius=Rad,material=mat1,fixed=False))
+
+O.bodies[id2].state.vel[2] = -1.0
+# Add engines
+o.engines = [
+  ForceResetter(),
+  InsertionSortCollider([Bo1_Sphere_Aabb()]),
+  InteractionLoop(
+    [Ig2_Sphere_Sphere_ScGeom()],
+    [Ip2_ViscElMat_ViscElMat_ViscElPhys()],
+    [Law2_ScGeom_ViscElPhys_Basic()],
+  ),
+  NewtonIntegrator(damping=0,gravity=[0,0,-9.81]),
+  PyRunner(command='addPlotData()',iterPeriod=100),
+]
+
+# Function to add data to plot
+def addPlotData(): 
+  try:
+    delta = (O.bodies[id2].state.pos[2]-O.bodies[id1].state.pos[2])-(2*Rad)
+    plot.addData(delta=delta, time1=O.time, time2=O.time, time3=O.time, time4=O.time,
+              Fn = O.interactions[0,1].phys.Fn,
+              Fv = O.interactions[0,1].phys.Fv,
+              deltaDot = O.bodies[id2].state.vel[2] - O.bodies[id1].state.vel[2])
+  except:
+    pass
+  
+plot.plots={'time1':('delta'), 'time2':('deltaDot'), 'time3':('Fn'), 'time4':('Fv')}; plot.plot()
+
+O.run(1)
+qt.View()
+
+#O.wait() ; plot.saveGnuplot('sim-data_Sphere')

=== modified file 'pkg/dem/ViscoelasticPM.cpp'
--- pkg/dem/ViscoelasticPM.cpp	2015-06-02 17:09:57 +0000
+++ pkg/dem/ViscoelasticPM.cpp	2015-12-16 06:39:37 +0000
@@ -101,7 +101,9 @@
 		// Prevent appearing of attraction forces due to a viscous component
 		// [Radjai2011], page 3, equation [1.7]
 		// [Schwager2007]
-		const Real normForceReal = phys.kn * geom.penetrationDepth + phys.cn * normalVelocity;
+		phys.Fn = phys.kn * geom.penetrationDepth;
+		phys.Fv = phys.cn * normalVelocity;
+		const Real normForceReal = phys.Fn + phys.Fv;
 		if (normForceReal < 0) {
 			phys.normalForce = Vector3r::Zero();
 		} else {

=== modified file 'pkg/dem/ViscoelasticPM.hpp'
--- pkg/dem/ViscoelasticPM.hpp	2015-06-04 21:02:26 +0000
+++ pkg/dem/ViscoelasticPM.hpp	2015-12-16 06:39:37 +0000
@@ -55,6 +55,8 @@
 	YADE_CLASS_BASE_DOC_ATTRS_CTOR(ViscElPhys,FrictPhys,"IPhys created from :yref:`ViscElMat`, for use with :yref:`Law2_ScGeom_ViscElPhys_Basic`.",
 		((Real,cn,NaN,,"Normal viscous constant"))
 		((Real,cs,NaN,,"Shear viscous constant"))
+		((Real,Fn,0.0,,"Normal force of the contact"))
+		((Real,Fv,0.0,,"Viscous force of the contact"))
 		((Real,mR,0.0,,"Rolling resistance, see [Zhou1999536]_."))
 #ifdef YADE_SPH
 		((bool,SPHmode,false,,"True, if SPH-mode is enabled."))