errors in plot, print, save file, micro data, time step, 3D show

Hello,

if this is really the code you are using, I got:
- a few SyntaxError because of missing commas.
- Spherepack AttributeError (should be SpherePack)
- NewtonIntegrato NameError (should be NewtonIntegrator)
Please next time try to post exactly the code you are using and make sure it is working.

> 1- When I run my code, no plot will show up however I coded plotting commands

Without O.run(), I got the plot. What is the behavior on your side without O.run?

> 2 ... when I try to modify my code a little

please provide actual code. Little modifications and have significant consequences.

> 3- The timp steps are so low and run will take a long time (more than 1 hour)

GlobalStiffnessTimeStepper should be reliable, so dt is close to the critical time step. Without changing the physics of the system (increase density or decrease stiffness), this is just how DEM works..

> 4- I cannot save the dada as a text file that I coded its command at the end of my code

please be more specific what is the problem. The file is not save? is saved but is empty? Do you get some error? If yes, which one? ....
I got the saved file OK..

> 5- I do not know how can I get micro information such as contact forces through my code

e.g. (depending on your specific needs)
for i in O.interactions:
   id1,id2 = i.id1,i.id2
   f = i.phys.normalForce + i.phys.shearForce

Basically from Python you have access to any data in the simulation.

> 6- I do not see on the screen the printed test text I used in my code.

Could not reproduce it. I got 4 warnings concerning SpherePack and then ****** TEST TEXT ****** line

> 7- please let me know if I have any other problem with my code

O.run() should not be directly followed by other commands. O.run() just starts the simulation and Python continues to evaluate further commands, which are evaluated "at random time" w.r.t Yade running.
Use O.run(N,True) # running N time steps and continue after it is finished
or
O.run(); O.wait()
This start running, but before evaluating other commands, Python waits for Yade to stop (e.g. some O.pause() inside PyRunner)

cheers
Jan

Thank you so much Jan. I appreciate your time and help. mostly you solved my problem. Just few things related to my current code in below:

1- Yes, I get an empty file which only contains the first line of title as e11, e22 .... , but no data is in the file (I want to make a clear and organized text file for post proccesing)
2- When I add "PyRunner(command='checkUnbalanced()',realPeriod=2)" to code in O.engines or those 2 print lines at the end of the code, the plot window doesn't show up
3- still I don't get any output for contact forces (let say data in saved file or plot)
4- in the plot command, the 'i' vs 's11,s22,s33' don't show up
5- I don't know the direction of 1,2,3 axises. is it like 3 for vertical axis, 2 for the horizontal axis and 1 for the third one?

and this is the output I get on Terminal:
Engine::action(): This engine is deprecated, please switch to TriaxialStressController if you expect long term support.
yade.plot: creating new line for kinetic
/home/ehsan/yade/install/lib/x86_64-linux-gnu/yade-2019-08-08.git-775ae74/py/yade/plot.py:505: MatplotlibDeprecationWarning:
The 'verts' kwarg was deprecated in Matplotlib 3.0 and will be removed in 3.2. Use 'marker' instead.
  scatter=ax.scatter(scatterPt[0],scatterPt[1],s=60,color=line.get_color(),**scatterMarkerKw)
yade.plot: creating new line for gravWork
yade.plot: creating new line for elastPotential
yade.plot: creating new line for plastDissip
yade.plot: creating new line for nonviscDamp
In [1]: The constructor with a shareWidget is deprecated, use the regular contructor instead.

this is the code:

from yade import pack, plot, qt

sp=pack.SpherePack()
mn, mx=Vector3(0,0,0), Vector3(10,10,10)
sp.makeCloud(minCorner=mn,maxCorner=mx,rRelFuzz=0.2,num=2000,porosity=0.4)

O.materials.append(FrictMat(young=15e6,poisson=0.4,frictionAngle=radians(30),density=2600,label='spheres'))
O.materials.append(FrictMat(young=15e6,poisson=0.4,frictionAngle=0,density=0,label='frictionless'))

O.bodies.append([sphere(center,rad,material='spheres') for center, rad in sp])

walls=aabbWalls(thickness=1e-10,material='frictionless')
wallIds=O.bodies.append(walls)

triax=TriaxialCompressionEngine(
 wall_bottom_id=wallIds[2],
 wall_top_id=wallIds[3],
 wall_left_id=wallIds[0],
 wall_right_id=wallIds[1],
 wall_back_id=wallIds[4],
 wall_front_id=wallIds[5],
 internalCompaction=False,
 sigmaIsoCompaction=-50e3,
 sigmaLateralConfinement=-50e3,
 max_vel=10,
 strainRate=0.01,
 label="triax",
)

O.engines=[
 ForceResetter(),
 InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Box_Aabb()]),
 InteractionLoop(
 [Ig2_Sphere_Sphere_ScGeom(),Ig2_Box_Sphere_ScGeom()],
 [Ip2_FrictMat_FrictMat_FrictPhys()],
 [Law2_ScGeom_FrictPhys_CundallStrack()]
 ),
 GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.8),
 triax,
 NewtonIntegrator(damping=0.4)
 #PyRunner(command='checkUnbalanced()',realPeriod=2),
]

O.engines=O.engines[0:5]+[PyRunner(iterPeriod=20,command='history()',label='recorder')]+O.engines[5:7]
O.dt=.5*PWaveTimeStep()

def history():
 plot.addData(unbalanced=unbalancedForce(),i=O.iter,**O.energy,
  e11=-triax.strain[0], e22=-triax.strain[1], e33=-triax.strain[2],
  s11=-triax.stress(0)[0], s22=-triax.stress(2)[1], s33=-triax.stress(4)[2],
  Etot=O.energy.total())

#O.trackEnergy=True

def checkUnbalanced():
 if unbalancedForce()<.05:
  O.pause()
  plot.saveDataTxt('bbb.txt.bz2')

for i in O.interactions:
   id1,id2 = i.id1,i.id2
   fn = i.phys.normalForce
   fs = i.phys.shearForce

#plot.plots={'i':('unbalanced'),'i':('s11','s22','s33'),'i':('e11','e22','e33'),'i':(O.energy.keys,None,'Etot')}
plot.plots={'i':('unbalanced'),'i':('s11','s22','s33'),'i':('e11','e22','e33'),'e22':('s22'),'e11':('s11')}

Gl1_Sphere.stripes=True

plot.plot()

O.saveTmp()

plot.saveDataTxt('results')

#print 'stress difference:', (triax.stress[2]-triax.stress[1])
#print 'mean stress:',(triax.stress[2]+triax.stress[1]+triax.stress[0])/3
#print 'porosity:', porosity()

Thanks again.

Hi,

> 1- Yes, I get an empty file which only contains the first line of title as e11, e22 .... , but no data is in the file (I want to make a clear and organized text file for post proccesing)

here you have no step/run in the script and at the end save the data, so the file is saved empty.
Try something like
O.run(100,True)
plot.saveDataTxt('results')

> 2- When I add "PyRunner(command='checkUnbalanced()',realPeriod=2)" to code in O.engines or those 2 print lines at the end of the code, the plot window doesn't show up

do you get any error (in the terminal)?

> 3- still I don't get any output for contact forces (let say data in saved file or plot)
> for i in O.interactions:
> id1,id2 = i.id1,i.id2
> fn = i.phys.normalForce
> fs = i.phys.shearForce

Because you did not output them in the script. I answered "how can I get micro information". By that code, you got fn and fs for all contacts and you can do whatever with it (save, plot, ...).
The correct output method really depends on what you want. Contact force history of a specific contact? All forces? At one time, multiple times? ...

> 4- in the plot command, the 'i' vs 's11,s22,s33' don't show up

because it is replaced by another 'i', use 'i', 'i ', 'i ', ... [1]

> 5- I don't know the direction of 1,2,3 axises. is it like 3 for vertical axis, 2 for the horizontal axis and 1 for the third one?

s11 is just a variable name, so 11 does not mean anything strict.
According to your history() function:
> e11=-triax.strain[0], e22=-triax.strain[1], e33=-triax.strain[2],
> s11=-triax.stress(0)[0], s22=-triax.stress(2)[1], s33=-triax.stress(4)[2]
11 = x, 22 = y, 33 = z

> triax=TriaxialCompressionEngine(
> Engine::action(): This engine is deprecated, please switch to TriaxialStressController if you expect long term support.

as the error says, consider using TriaxialStressController instead of TeiaxialCompressionEngine

cheers
Jan

[1] https://yade-dem.org/doc/user.html#multiple-figures

Hi Jan and thank you for your help.

1- problem solved

2- the output is:
(I add these lines to O.engines:
VTKRecorder(fileName='3d-vtk-',recorders=['all'],iterPeriod=1000),
qt.SnapshotEngine(fileBase='3d-',iterPeriod=200,label='snapshot'),
PyRunner(command='finish()',iterPeriod=20000)
PyRunner(command='checkUnbalanced()',realPeriod=2),)

ehsan@ehsan:~/Desktop/mycodes/test$ /home/ehsan/yade/install/bin/yade-2019-08-08.git-775ae74 test.py
Welcome to Yade 2019-08-08.git-775ae74
Using python version: 3.6.8 (default, Oct 7 2019, 12:59:55)
[GCC 8.3.0]
TCP python prompt on localhost:9000, auth cookie `skydue'
XMLRPC info provider on http://localhost:21000
Running script test.py
Traceback (most recent call last):
  File "/home/ehsan/yade/install/bin/yade-2019-08-08.git-775ae74", line 336, in runScript
    execfile(script,globals())
  File "/usr/lib/python3/dist-packages/past/builtins/misc.py", line 81, in execfile
    code = compile(source, filename, "exec")
  File "test.py", line 43
    VTKRecorder(fileName='3d-vtk-',recorders=['all'],iterPeriod=1000),
              ^
SyntaxError: invalid syntax
[[ ^L clears screen, ^U kills line. F12 controller, F11 3D view (press "h" in 3D view for help), F10 both, F9 generator, F8 plot. ]]

In [1]: <WARNING> InsertionSortCollider:261 virtual void InsertionSortCollider::action(): verletDist is set to 0 because no spheres were found. It will result in suboptimal performances, consider setting a positive verletDist in your script.
<WARNING> Shop:445 static Real Shop::PWaveTimeStep(boost::shared_ptr<Scene>): PWaveTimeStep has not found any suitable spherical body to calculate dt. dt is set to 1.0

3- lets say I want all forces . I want to plot the data and even save the data on a file.

4- problem is solved

5- problem solved

6- I will work on this "TriaxialStressController" and ask if I have any questions.

Thank you again for your valuable helps.

> I don't know the direction of 1,2,3 axises. is it like 3 for vertical axis, 2 for the horizontal axis and 1 for the third one?

For future reference: the vertical axis can be deduced from the gravity defined in users script. For instance if gravity=(10,0,0) you can be sure that 1 is the vertical direction.

More interestingly, 1,2,3 corresponds to x,y,z in the Qt window. Clicking those x,y,z buttons or checking the axis in 3D view (press "a") should help.

Thank you, that solved my problem. I open other topics for new questions, Thanks again