# Soil Settlement

106 views

Thanks Valclav.

Actually I really don't want to bother you because I know you also need to work and are busy as everyone instead of  Code babysitter.

Every question I asked here I tried to solve by myself with many times. The Plotting question I asked before is not so crucial for me now.

If you have time, please have a look at this code about the Soil Settelment to a RC plate, I think there's some problem during the engine creation.

from woo.core import *
from woo.dem import *
from woo.fem import *
import woo
import woo.gl
import math
from math import pi
from minieigen import *
import woo.dem
import woo.pack, woo.utils, numpy
woo.master.usesApi=10103
# create the scene
S=woo.master.scene=Scene(fields=[DemField(gravity=(0,0,-9.81))])
# graphic library
S.gl.demField.nodes=True
S.gl.node.wd=4
S.gl.membrane.node=False
S.gl.node.len=.05
S.gl.membrane.phiScale=0.
S.gl.demField.glyph=woo.gl.Gl1_DemField.glyphForce
S.gl.demField.nodes=False
S.gl.membrane.uScale=0.
S.gl.membrane.relPhi=0.
S.dem.gravity=(0,0,-30)

S.gl.demField.colorBy='disp'
# linear contact model
model=woo.models.ContactModelSelector(name='linear')
# create the box
wall=woo.dem.Wall.makeBox(
((0,0,0.1),(3.4,3.4,1.1)),
(1,1,0,1,1,0)
)
# create the mesh
mat_membrane=FrictMat(
young=30e9,
tanPhi=.3,
ktDivKn=.25,
density=2500
)
xmax,ymax=3.4,3.4
xdiv,ydiv=35,35
ff=woo.pack.gtsSurface2Facets(woo.pack.sweptPolylines2gtsSurface([[(x,y,0) for x in numpy.linspace(0,xmax,num=xdiv)] for y in numpy.linspace(0,ymax,num=ydiv)]),flex=True,halfThick=0.1)
# set boundary conditions
for n in S.dem.nodes:
n.dem.blocked=''
if (n.pos[0]>.25 and n.pos[0]<.35) and (n.pos[1]>.25 and n.pos[1]<.35): n.dem.blocked='xyz'
if (n.pos[0]>.25 and n.pos[0]<.35) and (n.pos[1]>3.05 and n.pos[1]<3.15): n.dem.blocked='xyz'
if (n.pos[0]>3.05 and n.pos[0]<3.15) and (n.pos[1]>.25 and n.pos[1]<.35): n.dem.blocked='xyz'
if (n.pos[0]>3.05 and n.pos[0]<3.15) and (n.pos[1]>3.05 and n.pos[1]<3.15): n.dem.blocked='xyz'

#S.engines=[Leapfrog(reset=True,damping=.05),IntraForce([In2_Membrane_ElastMat(thickness=.2,bending=True,bendThickness=.2)])]
# engines we want
S.engines=DemField.minimalEngines(model=model)+[
Leapfrog(reset=True,damping=.1),
InsertionSortCollider([Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()],verletDist=0.01),
ContactLoop([Cg2_Sphere_Sphere_L6Geom(),Cg2_Facet_Sphere_L6Geom(),Cg2_Wall_Sphere_L6Geom(),Cg2_Wall_Facet_L6Geom()],[Cp2_FrictMat_FrictPhys()],[Law2_L6Geom_FrictPhys_IdealElPl()],applyForces=False),
# forces are applied in IntraForce
IntraForce([In2_Membrane_FrictMat(thickness=0.2,bending=True,bendThickness=0.2)]),
# VtkExport(out='/tmp/membrane',stepPeriod=100,what=VtkExport.spheres|VtkExport.mesh)
]

# soil grains
sp=woo.pack.SpherePack()
sp.makeCloud((0,0,0.1),(3.4,3.4,1.1),rMean=.3*xmax/xdiv,rRelFuzz=.3,periodic=False,num=1000)
sp.toSimulation(S,mat=FrictMat(young=200e6,density=2650,ktDivKn=0.25,tanPhi=0.3))
# set fake inertia and mass to the node
for n in S.dem.nodes: DemData.setOriMassInertia(n)

S.dt=min(1e-4,.7*woo.utils.pWaveDt(S))

S.saveTmp()

import woo.qt
woo.qt.Controller()
woo.qt.View()

I need the soil grains settle to the RC membrane by the gravity, but the RC membrane's deflection is so large with a wrong value.

If I just make the RC membrane without soil grains (sphere.pack) and create the engine only with

S.engines=DemField.minimalEngines(damping=.01,verletDist=-0.01)+[IntraForce([In2_Membrane_ElastMat(thickness=0.2,bending=True,bendThickness=0.2)])]

Then the deflection is right with around 0.31 mm.

Could you help me solve this problem?

It's really grateful and a little bit sorry to bother you again.

Best Regards

Feixiang Xuan

+1 vote
answered Oct 7, 2016 by (46,490 points)
selected Oct 21, 2016 by eudoxos

Hi Feixiang/Xuan (sorry, what is your first name? some places in Asia write family name -- first name, some the other way),

sorry for the belated answer. You can use the "Inspector" button in the controller to see what engines you have.

DemField.minimalEngines(...) returns a pre-made set of engines which are usually needed for DEM simulations (integrator, collider, contact loop, dt control) with some pre-set parameters according to the contact model you choose and other params.

If you say

S.engines=DemField.minimalEngines(model=model)+[
Leapfrog(reset=True,damping=.1),
InsertionSortCollider([     Bo1_Sphere_Aabb(),Bo1_Facet_Aabb(),Bo1_Wall_Aabb()],verletDist=0.01),
ContactLoop([Cg2_Sphere_Sphere_L6Geom(),Cg2_Facet_Sphere_L6Geom(),Cg2_Wall_Sphere_L6Geom(),Cg2_Wall_Facet_L6Geom()],[Cp2_FrictMat_FrictPhys()],[Law2_L6Geom_FrictPhys_IdealElPl()],applyForces=False), # forces are applied in IntraForce IntraForce([In2_Membrane_FrictMat(thickness=0.2,bending=True,bendThickness=0.2)]), # VtkExport(out='/tmp/membrane',stepPeriod=100,what=VtkExport.spheres|VtkExport.mesh) ]

then if you look into Inspector, you will see there are two Leapfrog instances, two colliders, two contact loops. This will make a whole mess, as contact forces will be applied twice within one timestep and so on.

So either specify all engines by hand (without using minimalEngines), which needs some understanding of the internals (which is useful for you, anyway), or use minimalEngines but only add IntraForce([In2_Membrane_ElastMat(...)]). ContactLoop applies forces (by default, with applyForces=True) on uninodal particles (spheres, in your case), while IntraForce (with the default setting contacts=True) will take care of applying contact forces on the membrane.

HTH, Vaclav

commented Oct 13, 2016 by (46,490 points)

Hi, could you please open a separate question for this? It is a good question which deserves to be better visible. Cheers, v.

commented Oct 20, 2016 by (46,490 points)
edited Oct 20, 2016 by eudoxos

Hi, there is some stress implementation in the latest commit https://github.com/woodem/woo/commit/70e4969f4cb36aac91bd4b94bf5a408671560fd2 , if you pull from git, you will get it. Please read thoroughly https://woodem.org/theory/membrane-element.html#stresses and look at examples/membrane-stress.py as well.

I can't guarantee it is correct, you better test it yourself and also look at the implementation. Don't have time for more. Please really next time open separate questions for unrelated topics, I had hard time finding this question on stress. Cheers, v.

commented Oct 20, 2016 by (770 points)

Hi Vaclav