yade-users team mailing list archive

[liucheng83 <question247705@xxxxxxxxxxxxxxxxxxxxx>]

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

    Status: Answered => Open

liucheng83 is still having a problem:
Hi Bruno,
Thank you!

>>>(1) The default viscosity of fluid is 1, not 0.001Pa.s for water, why?
>>Some fluids are not water.

So if the viscosity of fluid is 0.001Pa.s , I can set the viscosity
value to 0.001, is it true?

And I have test the problems similar to the sedimentation in "4 Immersed
granular deposition" [Catalano2014a] Catalano, E., Chareyre, B.,
Barthélémy, E. (2014), Pore-scale modeling of fluid-particles
interaction and emerging poromechanical effects.

It works for some time, but it breaks down 35s , could you please give some advice for the following scripts,  thank you!
## ---------------------------------start

# basic simulation showing sphere falling ball gravity,
# bouncing against another sphere representing the support
from yade import pack
# DATA COMPONENTS
num_spheres=1000# number of spheres
young=1e6
compFricDegree = 3 # initial contact friction during the confining phase
finalFricDegree = 30 # contact friction during the deviatoric loading
mn,mx=Vector3(0,0,0),Vector3(2,1.5,0.15) # 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,0.03,0.333,num_spheres,False, 0.61,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=0, # If true the confining pressure is generated by growing particles
)

newton=NewtonIntegrator(damping=0.1,gravity=(0,-9.8,0))

# simulation loop -- see presentation for the explanation
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.wall_bottom_activated=False
triax.wall_top_activated=False
triax.wall_left_activated=False
triax.wall_right_activated=False
triax.wall_back_activated=False
triax.wall_front_activated=False

# Blocked certain degress of freedom to make 2D-Model in plane-XY
##for k in O.bodies:
## if isinstance(k.shape, Sphere): k.state.blockedDOFs='zXY'

# set timestep to a fraction of the critical timestep
# the fraction is very small, so that the simulation is not too fast
# and the motion can be observed
O.dt=.5e-4*utils.PWaveTimeStep()

O.run(1,1)
print "###------  ###"

# save the simulation, so that it can be reloaded later, for experimentation
O.saveTmp()

###---------------------------------------------------------------------------------------------------add
by lc

####A. Check bulk modulus of the dry material from load/unload cycles
###triax.stressMask=2
###triax.goal1=triax.goal3=0

###triax.internalCompaction=False
###triax.wall_bottom_activated=False
####load
###triax.goal2=11000; O.run(2000,1)
####unload
###triax.goal2=10000; O.run(2000,1)
####load
###triax.goal2=11000; O.run(2000,1)
###e22=triax.strain[1]
####unload
###triax.goal2=10000; O.run(2000,1)

###e22=e22-triax.strain[1]
###modulus = 1000./abs(e22)

#B. Activate flow engine and set boundary conditions in order to get permeability
flow.dead=0
flow.defTolerance=0.3
flow.meshUpdateInterval=200
flow.useSolver=3
flow.permeabilityFactor=1
flow.viscosity=100
flow.bndCondIsPressure=[0,0,0,1,0,0]
flow.bndCondValue=[0,0,0,0,0,0]
flow.boundaryUseMaxMin=[0,0,0,0,0,0]
O.dt=0.1e-3
O.dynDt=False

O.run(1,1)

Qin = flow.getBoundaryFlux(2)
Qout = flow.getBoundaryFlux(3)
permeability = abs(Qin)/O.dt
print "Qin=",Qin," Qout=",Qout," permeability=",permeability

#newton.damping=0

zeroTime=O.time


from yade import plot

## a function saving variables
def history():
	plot.addData(P1 = flow.averageSlicePressure(0.000),P2 = flow.averageSlicePressure(0.3),P3 = flow.averageSlicePressure(0.6),P4 = flow.averageSlicePressure(0.9),P5 = flow.averageSlicePressure(1.2),P6 = flow.averageSlicePressure(1.5),t=O.time)
  	#plot.addData(Qin = flow.getBoundaryFlux(5),Qout = flow.getBoundaryFlux(4),t=O.time,permeability = abs(flow.getBoundaryFlux(5))/2/0.75/O.dt)
  	#plot.addData(e22=triax.strain[1],t=O.time,s22=-triax.stress(2)[1],p=flow.MeasurePorePressure((0.5,0.5,0.5)))

O.engines=O.engines+[PyRunner(iterPeriod=200,command='history()',label='recorder')]
##make nice animations:
#O.engines=O.engines+[PyRunner(iterPeriod=200,command='flow.saveVtk()')]

from yade import plot
plot.plots={'t':('P1','P2','P3','P4','P5','P6')}
plot.plot()
O.saveTmp()
O.timingEnabled=1
from yade import timing
print "starting oedometer simulation"
O.run(200,1)
timing.stats()

## Make more steps to see the convergence to the stationnary solution

##----------------------------------------------------------------------------------------------------------------------------------------------end

-- 
You received this question notification because you are a member of
yade-users, which is an answer contact for Yade.