1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
|
# -*- encoding=utf-8 -*-
##
## SCRIPT TO TEST A NEW CONSTITUTIVE LAW (MINDLIN - nonlinear elastic model)
## list of engines
O.engines = [
ForceResetter(),
InsertionSortCollider([Bo1_Sphere_Aabb(), Bo1_Box_Aabb()]),
InteractionLoop([Ig2_Sphere_Sphere_ScGeom(), Ig2_Box_Sphere_ScGeom()], [Ip2_FrictMat_FrictMat_MindlinPhys()], [Law2_ScGeom_MindlinPhys_Mindlin()]),
NewtonIntegrator(damping=0.0, gravity=(10, 0, 0)),
###
### NOTE this extra engine:
###
### You want snapshot to be taken every 1 sec (realTimeLim) or every 50 iterations (iterLim),
### whichever comes soones. virtTimeLim attribute is unset, hence virtual time period is not taken into account.
PyRunner(iterPeriod=1, command='myAddPlotData()')
]
## define and append material
mat = FrictMat(young=600.0e6, poisson=0.6, density=2.60e3, frictionAngle=radians(26), label='Friction')
O.materials.append(mat)
## create two spheres (one will be fixed) and append them
s0 = sphere([0, 0, 0], 1, color=[0, 1, 0], fixed=False, wire=True, material='Friction')
s1 = sphere([2, 0, 0], 1, color=[0, 2, 0], fixed=True, wire=True, material='Friction')
O.bodies.append(s0)
O.bodies.append(s1)
## time step
O.dt = .2 * PWaveTimeStep()
O.saveTmp('Mindlin')
from yade import qt
qt.View()
qt.Controller()
############################################
##### now the part pertaining to plots #####
############################################
from yade import plot
## make one plot: step as function of fn
plot.plots = {'un': ('fn')}
## this function is called by plotDataCollector
## it should add data with the labels that we will plot
## if a datum is not specified (but exists), it will be NaN and will not be plotted
def myAddPlotData():
if O.interactions[0, 1].isReal:
i = O.interactions[0, 1]
## store some numbers under some labels
plot.addData(fn=i.phys.normalForce[0], step=O.iter, un=2 * s0.shape.radius - s1.state.pos[0] + s0.state.pos[0], kn=i.phys.kn)
O.run(100, True)
plot.plot(subPlots=False)
## We will have:
## 1) data in graphs (if you call plot.plot())
## 2) data in file (if you call plot.saveGnuplot('/tmp/a')
## 3) data in memory as plot.data['step'], plot.data['fn'], plot.data['un'], etc. under the labels they were saved
|
