Yade

O.load but still in the first code

Asked by jamespaul

Hi,

In my first script,only if maxUnbalanced<3e-3 ,the process will enter the next function and save the .gz file.But when i run the second script,the unbalancedForce is largger than 3e-3.Then until unbalancedForce is smaller than 3e-3 will start the compression.

1.Why the unbalancedForce become largger?

2.I've rewritten the triaxFinished function. Why still in the first compress until the unbalancedForce is small enough?

Thanks

James
###############
first script:

# encoding: utf-8
from yade import pack, qt, plot

sigmaIso=-25000

O.periodic=True

frictionAngle1=0.0286
frictionAngle11=14.03624347

sphere=O.materials.append(FrictMat(young=64e9,poisson=0.12,density=2650,frictionAngle=radians(frictionAngle1),label='sphere'))

multiple=.001
length=3
psdSizes=[.18*multiple,.2*multiple,.21*multiple,.25*multiple,.27*multiple,.30*multiple,.31*multiple,.34*multiple,.38*multiple]
psdCumm=[0.11,0.27,0.45,0.73,0.88,0.95,0.98,0.99,1]

sp=pack.SpherePack()
sp.makeCloud((0,0,0),(length*multiple,length*multiple,length*multiple),psdSizes=psdSizes,psdCumm=psdCumm,distributeMass=True,periodic=True)
n=sp.toSimulation()

O.engines=[
 ForceResetter(),
 InsertionSortCollider([Bo1_Sphere_Aabb()]),
 InteractionLoop(
    [Ig2_Sphere_Sphere_ScGeom()],
    [Ip2_FrictMat_FrictMat_FrictPhys()],
    [Law2_ScGeom_FrictPhys_CundallStrack()]
  ),
  PeriTriaxController(label='triax',
    goal=(sigmaIso,sigmaIso,sigmaIso),stressMask=7,
    dynCell=True,maxStrainRate=(1000,1000,1000),
    maxUnbalanced=3e-3,relStressTol=1e-3,
    doneHook='changeMaterial()'
  ),
  NewtonIntegrator(damping=.1),
 PyRunner(command='num()',initRun=1,nDo=1),
]
O.dt=1*PWaveTimeStep()

def num():
    print(O.bodies[-1].id)

def changeMaterial():
  O.cell.trsf=Matrix3.Identity
  for i in n:
          O.bodies[i].material.frictionAngle=radians(frictionAngle11)
 triax.doneHook='triaxFinished()'

def triaxFinished():
 O.cell.trsf=Matrix3.Identity
     print 'Finished'
     O.pause()

qt.Controller()
O.run()
O.wait()
O.save('four.yade.gz')
O.saveTmp()

#################
scecond script:
#!/usr/bin/python
from yade import pack, qt, plot

O.load('four.yade.gz')
O.run()
O.dt=1e-13

O.engines += [PyRunner(command='addPlotData2()',iterPeriod=1000)];

def triaxFinished():
 O.cell.trsf=Matrix3.Identity
 triax.goal=(0,0,-4000000)
 triax.stressMask=4
 triax.maxStrainRate=(0,0,30)
 triax.maxUnbalanced=1e-8
 triax.doneHook='functionToRunWhenFinished()'

def functionToRunWhenFinished():
    print 'Finished'
    O.pause()

def addPlotData2():
  plot.addData(unbalanced=unbalancedForce(),i=O.iter,
    sxx=-triax.stress[0],syy=-triax.stress[1],szz=-triax.stress[2],
    exx=-triax.strain[0],eyy=-triax.strain[1],ezz=-triax.strain[2],
    q=-triax.stress[2]+0.5*(triax.stress[0]+triax.stress[1]),
    p=-(triax.stress[0]+triax.stress[1]+triax.stress[2])/3,
  k=(triax.stress[0]*0.5+triax.stress[1]*0.5)/triax.stress[2]
  )
 print('sxx=',-triax.stress[0],'syy=',-triax.stress[1],'szz=',-triax.stress[2])
 print(porosity())

plot.plots={'i':('unbalanced',),'i ':('sxx','syy','szz'),' i':('exx','eyy','ezz'),' p ':('q'),}
plot.plot()

Question information

Language:
English Edit question
Status:
Solved
For:
Yade Edit question
Assignee:
No assignee Edit question
Solved by:
Jan Stránský
Solved:
Last query:
Last reply:
Revision history for this message
Jan Stránský (honzik) said : 		#1

Hi James,

> unbalancedForce is largger than 3e-3

please provide exact value

> 1.Why the unbalancedForce become largger?

from my experience, the unbalanced force has oscillating value. So once you stop the simulation when it reaches a small value, it is likely that after re-run, it exceeds this limit again

> 2.I've rewritten the triaxFinished function. Why still in the first compress until the unbalancedForce is small enough?
>
> def triaxFinished():
> ...
> triax.doneHook='functionToRunWhenFinished()'
>
> def functionToRunWhenFinished():
> print 'Finished'
> O.pause()

I would say, that in triaxFinished, the new doneHook is set, but it is not executed (there is no reason) and the simulation continues.. Why not to put "print 'Finished'" and "O.pause()" directly in triaxFinished?

One more note:

> O.load('four.yade.gz')
> O.run()
> # some settings

I would first set everything properly and then run O.run afterwards, otherwise you cannot be sure when the settings are applied (might be before the first iteration, but it is possible that the settings are applied after unknown number of time steps)

cheers
Jan

Revision history for this message
jamespaul (jamespauljames) said : 		#2

Thanks Jan,thank you very much,

>please provide exact value

0.016

> Why not to put "print 'Finished'" and "O.pause()" directly in triaxFinished?

I correct it.

def triaxFinished():
 O.cell.trsf=Matrix3.Identity
...
 triax.maxStrainRate=(0,0,.2)
 triax.doneHook='print "Hydrostatic load reached."; O.pause()'

>Though i have rewritten triaxFinished() in the second script,but the program still execute PeriTriaxController() in the first code.I print O.engines in the second script,

[<ForceResetter instance at 0x5615b0329cd0>, <InsertionSortCollider instance at 0x5615b1600660>, <InteractionLoop instance at 0x5615b13334f0>, <PeriTriaxController instance at 0x5615af815100>, <NewtonIntegrator instance at 0x5615b10d2880>, <PyRunner instance at 0x5615b129b1c0>, <PyRunner instance at 0x5615b1299ef0>]

How can i kill the PeriTriaxController()?

Revision history for this message
Jan Stránský (honzik) said : 		#3

> >please provide exact value
> 0.016

Do you have an idea what is the evolution of unbalanced force before? like how much it is oscillating etc.?

> Though i have rewritten triaxFinished() in the second script,but the program still execute PeriTriaxController() in the first code.
> How can i kill the PeriTriaxController()?

triax.dead = True

but then triaxFinished will never be called in the second script..

Jan

Revision history for this message
jamespaul (jamespauljames) said : 		#4

> Do you have an idea what is the evolution of unbalanced force before? like how much it is oscillating etc.?

Yes!I know about servo systems.

> The reason i need to separate the PeriTriaxController() into 2 parts is the stiffness of particle is large.First isotropic compression, then one-dimensional compression.If O.dt is big,the second part will distort.If O.dt is small,it will waste time in the first part.So,is there any method to save simulation time,such as change O.dt to a the specified value in the second part.

Revision history for this message
Best Jan Stránský (honzik) said : 		#5

>> Do you have an idea what is the evolution of unbalanced force before? like how much it is oscillating etc.?
> Yes!I know about servo systems.

Sorry for not being specific enough.. I meant this specific case. Maybe the unbalanced force oscillates like between 0.005 to 0.2 and then the behavior what you described is ok..
If it monotoncally converges to 3e-3 and then is 0.16, then it is another case

Sorry, I did not get the second part at all..

> The reason i need to separate the PeriTriaxController() into 2 parts is the stiffness of particle is large. First isotropic compression, then one-dimensional compression.

"large" is too vague term. The value "young=64e9" is quite normal I would say.
I see no reason how stiffness implies the need of "separate the PeriTriaxController() into 2 parts"..
Maybe the change of loading is the problem, see the last paragraph.

> If O.dt is big,the second part will distort.

yes, but then also the first part should distort..

> If O.dt is small,it will waste time in the first part. So,is there any method to save simulation time,such as change O.dt to a the specified value in the second part.

simulation time (O.time) is of course saved.
you can set O.dt to whatever value..

Do you know how PeriTriaxController behaves when changing stressMask and goal?
I don't use it and have no idea..
If this would be a problem, you can replace the controller by
###
triax = PeriTriaxController(...) # completely new instance with no history
O.engines = O.engines[:3] + [triax] + O.engines[4:]
###
It is a good idea to test both approaches (creating a new controller and modifying the existing one) to see if there is a significant difference.

Jan

Revision history for this message
Chareyre (bruno-chareyre-9) said : 		#6

My two cents: this is a typical situation where making the whole thing
happen in a unique script would help clarify. Inserting save/load in
between phase transitions only makes things more confused.

Le mar. 26 févr. 2019 14:08, jamespaul <email address hidden>
a écrit :

> Question #678799 on Yade changed:
> https://answers.launchpad.net/yade/+question/678799
>
> jamespaul posted a new comment:
> > Do you have an idea what is the evolution of unbalanced force before?
> like how much it is oscillating etc.?
>
> Yes!I know about servo systems.
>
> > The reason i need to separate the PeriTriaxController() into 2 parts
> is the stiffness of particle is large.First isotropic compression, then
> one-dimensional compression.If O.dt is big,the second part will
> distort.If O.dt is small,it will waste time in the first part.So,is
> there any method to save simulation time,such as change O.dt to a the
> specified value in the second part.
>
> --
> You received this question notification because your team yade-users is
> an answer contact for Yade.
>
> _______________________________________________
> Mailing list: https://launchpad.net/~yade-users
> Post to : <email address hidden>
> Unsubscribe : https://launchpad.net/~yade-users
> More help : https://help.launchpad.net/ListHelp
>
>
>

Revision history for this message
jamespaul (jamespauljames) said : 		#7

Thanks Jan Stránský, that solved my question.

Revision history for this message
jamespaul (jamespauljames) said : 		#8

Thanks Jan and Chareyre!