yade-users team mailing list archive

[内山康太郎 <question704090@xxxxxxxxxxxxxxxxxxxxx>]

New question #704090 on Yade:
https://answers.launchpad.net/yade/+question/704090

I am Japanese. Sorry for my rude English ...
and, I am new to yade .... I am also new to programming and DEM.

I want to apply a positive shear stress to stress σxz and then a negative shear force to stress σxz.
I want to .....
1st, isotropic compression

2nd, add positive shear
 (O.cell.velGrad = Matrix3(0, 0, 1, 0, 0, 0, 0, 0, 0, 0))

3rd. Proceed with analysis for positive the distorsion value (shear strain)
 ((O.cell.trsf[0, 2]) < .5)

4th. add negative shear
 (O.cell.velGrad = Matrix3(0, 0, -1, 0, 0, 0, 0, 0, 0, 0))

5th, Proceed with analysis for negative the distorsion value (shear strain)
 ((O.cell.trsf[0, 2]) >- .5)

6th, add positive shear
(O.cell.velGrad = Matrix3(0, 0, 1, 0, 0, 0, 0, 0, 0, 0, 0))

We want to loop from the second to the sixth.

but,I could not do this loop

When I run the following program, I get a loop where no negative shear strain is generated, only positive shear strain.

How can I generate negative shear strain?
and, Could you please give me some advice or program.
＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿＿
# encoding: utf-8

# script for periodic simple shear test, with periodic boundary
# first compresses to attain some isotropic stress (checkStress),
# then loads in shear (checkDistorsion)
#
# the initial packing is either regular (hexagonal), with empty bands along the boundary,
# or periodic random cloud of spheres
#
# material friction angle is initially set to zero, so that the resulting packing is dense
# (sphere rearrangement is easier if there is no friction)
#

# setup the periodic boundary
from __future__ import print_function

coords = [1+2*i for i in range(4)]
O.periodic = True
O.cell.hSize = Matrix3(2, 0, 0, 0, 2, 0, 0, 0, 2)

from yade import pack, plot

# the "if 0:" block will be never executed, therefore the "else:" block will be
# to use cloud instead of regular packing, change to "if 1:" or something similar
if 1:
    # create cloud of spheres and insert them into the simulation
    # we give corners, mean radius, radius variation
    sp = pack.SpherePack()
    sp.makeCloud((0, 0, 0), (.06, .06, .02), rMean=.00200, rRelFuzz=.2, periodic=True)
    # insert the packing into the simulation
    sp.toSimulation(color=(.5, .263, .0))  # pure blue
else:
    # in this case, add dense packing
    O.bodies.append(pack.regularHexa(pack.inAlignedBox((0, 0, 0), (2, 2, 2)), radius=.05, gap=0, color=(0, 0, 1)))

# create "dense" packing by setting friction to zero initially
O.materials[0].frictionAngle = 0

# simulation loop (will be run at every step)
O.engines += [
        ForceResetter(),
        InsertionSortCollider([Bo1_Sphere_Aabb()]),
        InteractionLoop(
                # interaction loop
                [Ig2_Sphere_Sphere_ScGeom()],
                [Ip2_FrictMat_FrictMat_FrictPhys()],
                [Law2_ScGeom_FrictPhys_CundallStrack()]
        ),
        NewtonIntegrator(damping=.4),
        # run checkStress function (defined below) every second
        # the label is arbitrary, and is used later to refer to this engine
        PyRunner(command='checkStress()', realPeriod=1, label='checker'),
        # record data for plotting every 100 steps; addData function is defined below
        PyRunner(command='addData()', iterPeriod=100)
]

# set the integration timestep to be 1/2 of the "critical" timestep
O.dt = .5 * PWaveTimeStep()

# prescribe isotropic normal deformation (constant strain rate)
# of the periodic cell
O.cell.velGrad = Matrix3(-10, 0, 0, 0, -10, 0, 0, 0, -10)

# when to stop the isotropic compression (used inside checkStress)



# called every second by the PyRunner engine
def checkStress():
    # stress tensor as the sum of normal and shear contributions
    # Matrix3.Zero is the intial value for sum(...
    if abs(getStress()[0,0]) < 1e4:
        O.step()
    # if mean stress is below (bigger in absolute value) limitMeanStress, start shearing
    # cyclic loading

    for cycle in range(6):
        # positive shear
            O.cell.velGrad = Matrix3(0,0, 1, 0,0,0, 0,0,0)
            while (O.cell.trsf[0, 2]) < .4:
               O.step()
            # negative shear
            O.cell.velGrad = Matrix3(0,0,-1, 0,0,0, 0,0,0)
            while (O.cell.trsf[0, 2]) > -.4:
                   O.step()

# called periodically to store data history
# called periodically to store data history
def addData():
    # get the stress tensor (as 3x3 matrix)
    stress = sum(normalShearStressTensors(), Matrix3.Zero)
    # give names to values we are interested in and save them
    plot.addData(exz=O.cell.trsf[0, 2], szz=stress[2, 2], sxz=stress[0, 2], tanPhi=(stress[0, 2] / stress[2, 2]) if stress[2, 2] != 0 else 0, i=O.iter)
    # color particles based on rotation amount
    for b in O.bodies:
        # rot() gives rotation vector between reference and current position
        b.shape.color = scalarOnColorScale(b.state.rot().norm(), 0, pi / 2.)

# define what to plot (3 plots in total)
## exz(i), [left y axis, separate by None:] szz(i), sxz(i)
## szz(exz), sxz(exz)
## tanPhi(i)
# note the space in 'i ' so that it does not overwrite the 'i' entry
plot.plots = {'i': ('exz', None), 'exz': ('szz', 'sxz'), 'i ': ('tanPhi',)}

# better show rotation of particles
Gl1_Sphere.stripes = True

# open the plot on the screen
plot.plot()

O.saveTmp()

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