2D Granular Chute Flow

Don’t respond unless its going to the entire group. I have no personal qualm with you. That being said, your method for calculating a diffusion constant still doesn’t make sense. The slope of the MSD (measured on ‘long’ time scales) only holds for a single diffusion coefficient that does not depend on space. So if you’re trying to get at the diffusion coefficient it makes more sense to do canonical flow problems like simple shear. You can report what you want, but it won’t be a diffusion coefficient. That IS why they don’t report a diffusion coefficient.

I reiterate that no one is going to track such a specific problem of matching with literature down - that’s a matter of contacting authors directly.

I personally trust LAMMPS granular models and have used them a great deal, and tested the models down to matching on the level of single collisions.

It’s not really a diffusion coefficient in the strictest sense of the word in the way that a molecular solute has a diffusion coefficient. But in case it’d be of interest, I was only exploring ‘diffusion’ in granular flows due to this 2014 paper. Right now, I am putting aside trying to obtain the diffusion coefficient until I’ve read enough papers on previous studies related to it.

My question this time is about chute flows in LAMMPS. As I’ve mentioned previously, when I run my simulations I always get flagged by an error which reads ‘Lost atoms: original x current y.’ The documentation page for this specific error message attributes the loss of atoms to bad dynamics and I cannot quite figure out why since my script is rather simple – it’s just a flow of 9900 0.5mm beads (discs) over a bed of discs with twice the diameter:

units si
dimension 2
atom_style sphere
atom_modify id yes
boundary p fs p
newton off
comm_modify vel yes

variable diam equal 0.0005
lattice sq ${diam}

region hill block 0 100 0 100 0 1 units lattice
create_box 1 hill

However, I wrote another script in lj units while retaining g=9.8 and the 2:1 diameter ratio of the roughened surface to the flowing assembly (contained in data.2dposgen data file which containing the positions of the atoms) and the problem goes away.

(sorry if you’re seeing this msg twice. the first one wasn’t sent. here it is again in its entirety)

It’s not really a diffusion coefficient in the strictest sense of the word in the way that a molecular solute has a diffusion coefficient. But in case it’d be of interest, I was only exploring ‘diffusion’ in granular flows due to this 2014 paper. Right now, I am putting aside trying to obtain the diffusion coefficient until I’ve read enough papers on previous studies related to it.

My question this time is about chute flows in LAMMPS. As I’ve mentioned previously, when I run my simulations I always get flagged by an error which reads ‘Lost atoms: original x current y.’ The documentation page for this specific error message attributes the loss of atoms to bad dynamics and I cannot quite figure out why since my script is rather simple – it’s just a flow of 9900 0.5mm beads (discs) over a bed of discs with twice the diameter:

units si
dimension 2
atom_style sphere
atom_modify id yes
boundary p fs p
newton off
comm_modify vel yes

variable diam equal 0.0005
lattice sq ${diam}
region hill block 0 100 0 100 0 1 units lattice
create_box 1 hill

neighbor 0.001 bin
neigh_modify delay 0

variable kn equal 200000
variable dcnormal equal 33.5

variable dctangent equal 0
variable syc equal 0.5
variable tdf equal 1

pair_style gran/hertz/history {kn} NULL {dcnormal} {dctangent} {syc} ${tdf}
pair_coeff * *
timestep 0.0001

create_atoms 1 region hill
group gm id > 100
group bed id <> 1 100
set group bed diameter 0.001 density/disc 2450
set group gm diameter ${diam} density/disc 2450

fix 1 gm nve/sphere disc
fix 2 bed freeze
fix 4 all wall/gran hertz/history {kn} NULL {dcnormal} NULL {syc} {tdf} yplane 0 NULL
reset_timestep 0
fix 7 all gravity 9.8 chute 18
thermo 100

thermo_style custom step atoms ke pe etotal

thermo_modify lost warn norm no
dump dc1 all custom 100 diff.incflo.grandiff* id x y
fix 5 all enforce2d
run 5000

However, I wrote another script in lj units while (retaining g=9.8) and the 2:1 diameter ratio of the roughened surface to the flowing assembly (positions contained in data.2dposgen ) and the problem goes away.

units lj
dimension 2
atom_style sphere
boundary p fs p

newton off
comm_modify vel yes
read_data data.2dposgen
group gm id > 100
group bed id <> 1 100
neighbor 0.1 bin
neigh_modify delay 0
variable kn equal 200000
variable dcnormal equal 33.5

variable dctangent equal 0
variable syc equal 0.5
variable tdf equal 0
pair_style gran/hooke/history {kn} NULL {dcnormal} {dctangent} {syc} ${tdf}
pair_coeff * *
timestep 0.0001
fix 1 gm nve/sphere disc
fix 2 bed freeze
fix 4 all gravity 9.8 chute 23.0
thermo 100
thermo_style custom step atoms ke pe etotal
thermo_modify lost warn norm no
fix 5 all enforce2d
run 110000

What could be the possible reason/s for the ‘loss of atoms’ in si and the non-appearance of this error in lj units?

Lastly, the documentation page for fix id gravity g chute theta has g=1 as a default value and is also g value in the sample chute input script that come with LAMMPS. Since g is in units of force/mass and force in lj units is f*/m = fsigma/(mass*epsilon) with mass=epsilon=sigma=1, would it not make sense to set the value of g to 9.8 instead of 1?

Thank you for your patience and guidance.

Hi Eric. My messages got trimmed into sections (i don’t know why). Expand the thread to view it whole. It picks up from our previous exchange last 4 Feb 2019.

So this may seem like a LAMMPS problem on the surface, but its really quite fundamental.

I highly suggest non-dimensionalizing your equations of motion, and looking at all the inputs - what changes between setting units to SI or LJ? What is not being scaled consistently? This should clear up your issue. You’ll see that it doesn’t really matter what the units are, or the particular values you give for the size of a particle or gravity. So long as everything is scaled correctly there will be no difference in the dynamics. This is useful for example when MD folks set sigma = 1. In SI units its really on the scale of 10^(-10) m, but using such small numbers is not really advantageous from a numerical perspective - but you are still able to transform between units to say get the values in SI.

You are in luck that there are actually very few nondimensional parameters in this set-up. (k_springD_particle/massgravity) is one of them.

So this may seem like a LAMMPS problem on the surface, but its really quite fundamental.

I highly suggest non-dimensionalizing your equations of motion, and looking at all the inputs - what changes between setting units to SI or LJ? What is not being scaled consistently? This should clear up your issue. You'll see that it doesn't really matter what the units are, or the particular values you give for the size of a particle or gravity. So long as everything is scaled correctly there will be no difference in the dynamics. This is useful for example when MD folks set sigma = 1. In SI units its really on the scale of 10^(-10) m, but using such small numbers is not really advantageous from a numerical perspective - but you are still able to transform between units to say get the values in SI.

FWIW, with floating point math, the order of magnitude, doesn't make
much of a difference for as long as numbers *can* be represented.
however, what is often a big problem when switching between reduced
units and non-reduced units, is the unit of time. as that is not
reduced by itself, but derived from length, mass, and energy (which
makes sense, since e.g. increasing the mass slows down motion, which
is equivalent to a shorter time unit).

axel.

Great points.