← Back to team overview

yade-users team mailing list archive

Re: [Question #679320]: squashing modelling

 

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

    Status: Open => Answered

Chareyre proposed the following answer:
For 2/
You write force component in some variable at some point, then you record
that variable. Of course it's constant.
Like saying that speed of wind at time t_0 is independent of time t, which
is true and expected.
Bruno

Le mar. 19 mars. 2019 17:08, Rioual <question679320@xxxxxxxxxxxxxxxxxxxxx>
a écrit :

> Question #679320 on Yade changed:
> https://answers.launchpad.net/yade/+question/679320
>
>     Status: Answered => Open
>
> Rioual is still having a problem:
> Hello Bruno
>
> Thank's for your input.
>
> 1- I tried the simulation with different smaller time scales 1e-4, 1e-5,
> 1e-6 and the error messages appear at fairly the same simulation time
> apparently:
> "
> CHOLMOD warning: matrix not positive definite
> something went wrong in Cholesky factorization, use LDLt as fallback this
> time
> 65 : Vh==NULL!! id=65 Point=-0.803797 -0.683743 1.61915 rad=0.0666365
> 147 : Vh==NULL!! id=147 Point=-0.927542 -0.647778 1.92058 rad=0.0638237
> ...etc"
>
> 2- Before the error messages, the compression is processed as can be
> observed on the graphical output (motion of the boundaries) but nothing
> changes when I print the position of the upper plate and the force on
> this plate. See a minimal working example code of the compression
> process below that shows both problems 1 and 2:
>
>
> ****************************************************************************
>
> ## ______________   First section, similar to
> triax-tutorial/script-session1.py  _________________
> from yade import pack
>
> num_spheres=1000# number of spheres
> young=1e6
> C=1e6
> compFricDegree = 3 # initial contact friction during the confining phase
> finalFricDegree = 11 # contact friction during the deviatoric loading (for
> ice!!)
> mn,mx=Vector3(0,0,0),Vector3(1,1,1) # corners of the initial packing
>
>
> O.materials.append(FrictMat(young=young,poisson=0.5,frictionAngle=radians(compFricDegree),density=2600,label='spheres'))
>
> O.materials.append(FrictMat(young=young,poisson=0.5,frictionAngle=0,density=0,label='walls'))
> walls=aabbWalls([mn,mx],thickness=0,material='walls')
> wallIds=O.bodies.append(walls)
>
> sp=pack.SpherePack()
> sp.makeCloud(mn,mx,-1,0.3333,num_spheres,False, 0.95,seed=1) #"seed" make
> the "random" generation always the same
> sp.toSimulation(material='spheres')
>
> triax=TriaxialStressController(
>         maxMultiplier=1.+2e4/young, # spheres growing factor (fast growth)
>         finalMaxMultiplier=1.+2e3/young, # spheres growing factor (slow
> growth)
>         thickness = 0,
>         stressMask = 7,
>         max_vel = 0.005,
>         internalCompaction=True, # If true the confining pressure is
> generated by growing particles
> )
>
> newton=NewtonIntegrator(damping=0.2)
>
> 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()],label="iloop"
>         ),
>         FlowEngine(dead=1,label="flow"),#introduced as a dead engine for
> the moment, see 2nd section
>
> GlobalStiffnessTimeStepper(active=1,timeStepUpdateInterval=100,timestepSafetyCoefficient=0.8),
>         triax,
>         newton
> ]
>
> triax.goal1=triax.goal2=triax.goal3=-10000
>
> while 1:
>   print '***'
>   O.run(1000, True)
>   unb=unbalancedForce()
>   if unb<0.001 and abs(-10000-triax.meanStress)/10000<0.001:
>     break
>
> setContactFriction(radians(finalFricDegree))
>
> ##Second section: Activate flow engine and set boundary conditions for
> the compression (squashing) test, impermeable on each face:
>
> flow.dead=0
> flow.defTolerance=0.3
> flow.meshUpdateInterval=200
> flow.useSolver=3
> flow.permeabilityFactor=1
> flow.viscosity=10
> flow.boundaryUseMaxMin=[0,0,0,0,0,0]
> O.dt=1e-4
> O.dynDt=False
>
> flow.bndCondIsPressure=[False,False,False,False,False,False]
> #flow.bndCondValue=[0,0,0,0,0,0]
>
> newton.damping=0
>
> #aditional compressibility of the fluid
> flow.fluidBulkModulus=1e6
>
> triax.stressMask=5
> triax.internalCompaction=False
> triax.goal2=-0.1
> triax.goal1= 0
> triax.goal3= 0
> triax.wall_left_activated= True
> triax.wall_right_activated= True
> triax.wall_front_activated= True
> triax.wall_back_activated= True
> triax.wall_top_activated= True
>
> global plate
> plate = O.bodies[triax.wall_top_id]
>
> Ftot0 = O.forces.f(plate.id)[0]
> Ftot1 = O.forces.f(plate.id)[1]
> Ftot2 = O.forces.f(plate.id)[2]
>
> xtot = plate.state.pos[0]
> ytot = plate.state.pos[1]
> ztot = plate.state.pos[2]
>
>
> ## a function printing variables
> def history():
>         print 't=', O.time, 'Force on the plate (Ftot)=', Ftot0, Ftot1,
> Ftot2
>         print 't=', O.time, 'Position of the plate (xtot,ytot,ztot)=',
> xtot, ytot, ztot
>
>
> O.engines=O.engines+[PyRunner(iterPeriod=200,command='history()',label='recorder')]
>
> from yade import timing
> print "starting squashing simulation"
> O.run(200,1)
> timing.stats()
>
>
> ************************************************************************************************
>
>
> Best wishes,
>
> Fr.
>
> --
> 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     : yade-users@xxxxxxxxxxxxxxxxxxx
> Unsubscribe : https://launchpad.net/~yade-users
> More help   : https://help.launchpad.net/ListHelp
>
>
>

-- 
You received this question notification because your team yade-users is
an answer contact for Yade.