This issue was already happening with previous versions of NUFEB, could it be a problem of using a Lennard-Jones with atoms of type coccus?

atom.in

NUFEB Simulation 

       0 atoms
       3 atom types

       0   2e-05  xlo xhi
       0   2e-05  ylo yhi
       0   2e-05  zlo zhi

A small 8 atoms sphere in the middle of the simulation puts the pressure at 2.9760175e+17.
```
# NUFEB simulation with methanogens and GAC floor
```

units si # using si units atom_style coccus # using nufeb atom style

map array: find atoms using indices sort 1000 5.0e-6: sort every 1000 steps with 5.0e-6 binsize

atom_modify map array sort 1000 1e-6

periodic boundaries in x and y fixed boundary in z

boundary ff ff ff

forces between local and ghost atoms are computed in each processor without communication

newton off

processors * # processor grid

communicate velocities for ghost atoms

comm_modify vel yes

guarantee that enough atoms are communicated to correctly compute

comm_modify cutoff 2e-6

read_data atom.in

Shift the lattice grid 0.5 unit so that atom can be created in the center of the grid

lattice sc 1e-6 origin 0.5 0.5 0.5 region reg block 0 20 0 20 1 20

create_atoms 1 random 80 31324 reg

region gac_sphere1 sphere 10 10 10 1 create_atoms 3 region gac_sphere1

set type 1 density 150 set type 1 diameter 1.3e-6 set type 1 mass 1.725e-16 set type 1 outer_diameter 1.3e-6 set type 1 outer_density 30

define attributes for type 3, set gac diameter to 1e-6, make it consistent with lattice size

set type 3 density 150 set type 3 mass 1.725e-16 set type 3 diameter 1e-6 set type 3 outer_diameter 1e-6 set type 3 outer_density 30

set the groups

group HET type 1 group EPS type 2 group GAC type 3

setting neighbour skin distance and style

neighbor 7e-7 bin

rebuild neighbour list if any atom had moved more than half the skin distance

neigh_modify check yes

defining grid sytle, substrate names, and grid size

grid_style nufeb/chemostat 4 sub o2 no2 no3 4e-6

set diffusion boundary conditions and initial concentrations (liquid:kg/m3)

grid_modify set sub pp pp nd 1e-4
grid_modify set o2 pp pp nd 1e-4 grid_modify set no2 pp pp nd 1e-4 grid_modify set no3 pp pp nd 1e-4

define pair styles

pair_style hybrid gran/hooke/history 1e-4 NULL 1e-5 NULL 0.0 1 lj/cut 2.5e-6 pair_coeff 2 2 gran/hooke/history pair_coeff 3 3 lj/cut 1.0 1.0e-6 2.5e-6 pair_coeff *2 3 lj/cut 1.0e-6 1.0e-6 2.5e-6

NVE integration with maximum distance limit -> only update position for non gac atoms

fix nve1 HET nve/limit 1e-8 fix nve2 EPS nve/limit 1e-8

heterotrophs growth

fix growth_het HET nufeb/growth/het sub 3.5e-5 o2 0 no2 0 no3 0 & growth 0.00028 yield 0.61 decay 0.0 epsyield 0.18 anoxic 0.0 epsdens 30

heterotrophs division, division diameter: 1.36e-6

fix div HET nufeb/division/coccus 1.36e-6 1234

EPS secretion from heterotrophs

fix eps_ext HET nufeb/eps_extract 2 EPS 1.3 30 2345

diffusion reaction fixes

fix diff_sub all nufeb/diffusion_reaction sub 1.6e-9

fix div HET nufeb/division/coccus 1.36e-6 1234 fix eps_ext HET nufeb/eps_extract 2 EPS 1.3 30 5678

mechanical model fixes

fix wall all wall/gran hooke/history 1e-3 NULL 1e-4 NULL 0 0 zplane 0.0 8e-5 fix eps_adh all nufeb/adhesion/eps EPS 1e-6 fix vis all viscous 1e-5

pressure computation

compute vol all nufeb/volume compute ke all ke variable one equal 1.0 compute press all pressure NULL pair vol v_one variable press equal "(c_ke + c_press) / (3.0 * c_vol)" variable mass equal "mass(all)"

variable nmet equal "count(HET)" variable neps equal "count(EPS)" compute mycon all nufeb/ave_conc

file output

shell mkdir vtk_gac_rigids dump 1 all vtk 10 vtk_gac_rigids/dump.vtu id type diameter dump 2 all grid/vtk 10 vtk_gacrigids/dump%_.vti con

thermo output

thermo_style custom step atoms v_press v_mass v_nmet v_neps c_mycon[*] thermo 1 thermo_modify lost ignore flush yes

issue run command

run_style nufeb diffdt 1e-4 difftol 1e-12 pairdt 2 pairtol 0 pairmax 500 diffmax 50000 timestep 1000 run 300

2. The weirdest one. Simply creating random atoms of type 3 puts the pressure at ```-5.2855696e+15```.

NUFEB simulation with methanogens and GAC floor

units si # using si units atom_style coccus # using nufeb atom style

map array: find atoms using indices sort 1000 5.0e-6: sort every 1000 steps with 5.0e-6 binsize

atom_modify map array sort 1000 1e-6

periodic boundaries in x and y fixed boundary in z

boundary ff ff ff

forces between local and ghost atoms are computed in each processor without communication

newton off

processors * # processor grid

communicate velocities for ghost atoms

comm_modify vel yes

guarantee that enough atoms are communicated to correctly compute

comm_modify cutoff 2e-6

read_data atom.in

Shift the lattice grid 0.5 unit so that atom can be created in the center of the grid

lattice sc 1e-6 origin 0.5 0.5 0.5 region reg block 0 20 0 20 1 20

create_atoms 1 random 80 31324 reg

create_atoms 3 random 80 31325 reg

set type 1 density 150 set type 1 diameter 1.3e-6 set type 1 mass 1.725e-16 set type 1 outer_diameter 1.3e-6 set type 1 outer_density 30

define attributes for type 3, set gac diameter to 1e-6, make it consistent with lattice size

set type 3 density 150 set type 3 mass 1.725e-16 set type 3 diameter 1e-6 set type 3 outer_diameter 1e-6 set type 3 outer_density 30

set the groups

group HET type 1 group EPS type 2 group GAC type 3

setting neighbour skin distance and style

neighbor 7e-7 bin

rebuild neighbour list if any atom had moved more than half the skin distance

neigh_modify check yes

defining grid sytle, substrate names, and grid size

grid_style nufeb/chemostat 4 sub o2 no2 no3 4e-6

set diffusion boundary conditions and initial concentrations (liquid:kg/m3)

grid_modify set sub pp pp nd 1e-4
grid_modify set o2 pp pp nd 1e-4 grid_modify set no2 pp pp nd 1e-4 grid_modify set no3 pp pp nd 1e-4

define pair styles

pair_style hybrid gran/hooke/history 1e-4 NULL 1e-5 NULL 0.0 1 lj/cut 2.5e-6 pair_coeff 2 2 gran/hooke/history pair_coeff 3 3 lj/cut 1.0 1.0e-6 2.5e-6 pair_coeff *2 3 lj/cut 1.0e-6 1.0e-6 2.5e-6

NVE integration with maximum distance limit -> only update position for non gac atoms

fix nve1 HET nve/limit 1e-8 fix nve2 EPS nve/limit 1e-8

heterotrophs growth

fix growth_het HET nufeb/growth/het sub 3.5e-5 o2 0 no2 0 no3 0 & growth 0.00028 yield 0.61 decay 0.0 epsyield 0.18 anoxic 0.0 epsdens 30

heterotrophs division, division diameter: 1.36e-6

fix div HET nufeb/division/coccus 1.36e-6 1234

EPS secretion from heterotrophs

fix eps_ext HET nufeb/eps_extract 2 EPS 1.3 30 2345

diffusion reaction fixes

fix diff_sub all nufeb/diffusion_reaction sub 1.6e-9

fix div HET nufeb/division/coccus 1.36e-6 1234 fix eps_ext HET nufeb/eps_extract 2 EPS 1.3 30 5678

mechanical model fixes

fix wall all wall/gran hooke/history 1e-3 NULL 1e-4 NULL 0 0 zplane 0.0 8e-5 fix eps_adh all nufeb/adhesion/eps EPS 1e-6 fix vis all viscous 1e-5

pressure computation

variable nmet equal "count(HET)" variable neps equal "count(EPS)" compute mycon all nufeb/ave_conc

file output

shell mkdir vtk_gac_rigids dump 1 all vtk 10 vtk_gac_rigids/dump.vtu id type diameter dump 2 all grid/vtk 10 vtk_gacrigids/dump%_.vti con

thermo output

thermo_style custom step atoms v_press v_mass v_nmet v_neps c_mycon[*] thermo 1 thermo_modify lost ignore flush yes

issue run command

run_style nufeb diffdt 1e-4 difftol 1e-12 pairdt 2 pairtol 0 pairmax 500 diffmax 50000 timestep 1000 run 300

3. Mountains of immovable atoms of type 3 put the pressure at ```-nan```.

NUFEB simulation with methanogens and GAC floor

units si # using si units atom_style coccus # using nufeb atom style

map array: find atoms using indices sort 1000 5.0e-6: sort every 1000 steps with 5.0e-6 binsize

atom_modify map array sort 1000 1e-6

periodic boundaries in x and y fixed boundary in z

boundary ff ff ff

forces between local and ghost atoms are computed in each processor without communication

newton off

processors * # processor grid

communicate velocities for ghost atoms

comm_modify vel yes

guarantee that enough atoms are communicated to correctly compute

comm_modify cutoff 2e-6

read_data atom.in

Shift the lattice grid 0.5 unit so that atom can be created in the center of the grid

lattice sc 1e-6 origin 0.5 0.5 0.5 region reg block 0 20 0 20 1 20

create_atoms 1 random 80 31324 reg

region gac_floor block 0 20 0 20 0 1 create_atoms 3 region gac_floor region gac_limits block 0 20 0 20 0 20 variable xx internal 0.0 variable yy internal 0.0 variable zz internal 0.0 variable v equal "v_zz < sin(sin(5e5v_xx))sin(cos(5e5*v_yy))/10e4" create_atoms 3 region gac_limits var v set x xx set y yy set z zz

set type 1 density 150 set type 1 diameter 1.3e-6 set type 1 mass 1.725e-16 set type 1 outer_diameter 1.3e-6 set type 1 outer_density 30

define attributes for type 3, set gac diameter to 1e-6, make it consistent with lattice size

set type 3 density 150 set type 3 mass 1.725e-16 set type 3 diameter 1e-6 set type 3 outer_diameter 1e-6 set type 3 outer_density 30

set the groups

group HET type 1 group EPS type 2 group GAC type 3

setting neighbour skin distance and style

neighbor 7e-7 bin

rebuild neighbour list if any atom had moved more than half the skin distance

neigh_modify check yes

defining grid sytle, substrate names, and grid size

grid_style nufeb/chemostat 4 sub o2 no2 no3 4e-6

set diffusion boundary conditions and initial concentrations (liquid:kg/m3)

grid_modify set sub pp pp nd 1e-4
grid_modify set o2 pp pp nd 1e-4 grid_modify set no2 pp pp nd 1e-4 grid_modify set no3 pp pp nd 1e-4

define pair styles

pair_style hybrid gran/hooke/history 1e-4 NULL 1e-5 NULL 0.0 1 lj/cut 2.5e-6 pair_coeff 2 2 gran/hooke/history pair_coeff 3 3 lj/cut 1.0 1.0e-6 2.5e-6 pair_coeff *2 3 lj/cut 1.0e-6 1.0e-6 2.5e-6

NVE integration with maximum distance limit -> only update position for non gac atoms

fix nve1 HET nve/limit 1e-8 fix nve2 EPS nve/limit 1e-8

heterotrophs growth

fix growth_het HET nufeb/growth/het sub 3.5e-5 o2 0 no2 0 no3 0 & growth 0.00028 yield 0.61 decay 0.0 epsyield 0.18 anoxic 0.0 epsdens 30

heterotrophs division, division diameter: 1.36e-6

fix div HET nufeb/division/coccus 1.36e-6 1234

EPS secretion from heterotrophs

fix eps_ext HET nufeb/eps_extract 2 EPS 1.3 30 2345

diffusion reaction fixes

fix diff_sub all nufeb/diffusion_reaction sub 1.6e-9

fix div HET nufeb/division/coccus 1.36e-6 1234 fix eps_ext HET nufeb/eps_extract 2 EPS 1.3 30 5678

mechanical model fixes

fix wall all wall/gran hooke/history 1e-3 NULL 1e-4 NULL 0 0 zplane 0.0 8e-5 fix eps_adh all nufeb/adhesion/eps EPS 1e-6 fix vis all viscous 1e-5

pressure computation

variable nmet equal "count(HET)" variable neps equal "count(EPS)" compute mycon all nufeb/ave_conc

file output

shell mkdir vtk_gac_rigids dump 1 all vtk 10 vtk_gac_rigids/dump.vtu id type diameter dump 2 all grid/vtk 10 vtk_gacrigids/dump%_.vti con

thermo output

thermo_style custom step atoms v_press v_mass v_nmet v_neps c_mycon[*] thermo 1 thermo_modify lost ignore flush yes

issue run command

run_style nufeb diffdt 1e-4 difftol 1e-12 pairdt 2 pairtol 0 pairmax 500 diffmax 50000 timestep 1000 run 300

nufeb / NUFEB-2

Pressure at -nan or very high values (negative and positive) #14

map array: find atoms using indices sort 1000 5.0e-6: sort every 1000 steps with 5.0e-6 binsize

periodic boundaries in x and y fixed boundary in z

forces between local and ghost atoms are computed in each processor without communication

communicate velocities for ghost atoms

guarantee that enough atoms are communicated to correctly compute

Shift the lattice grid 0.5 unit so that atom can be created in the center of the grid

define attributes for type 3, set gac diameter to 1e-6, make it consistent with lattice size

set the groups

setting neighbour skin distance and style

rebuild neighbour list if any atom had moved more than half the skin distance

defining grid sytle, substrate names, and grid size

set diffusion boundary conditions and initial concentrations (liquid:kg/m3)

define pair styles

NVE integration with maximum distance limit -> only update position for non gac atoms

heterotrophs growth

heterotrophs division, division diameter: 1.36e-6

EPS secretion from heterotrophs

diffusion reaction fixes

mechanical model fixes

pressure computation

file output

thermo output

issue run command

NUFEB simulation with methanogens and GAC floor

map array: find atoms using indices sort 1000 5.0e-6: sort every 1000 steps with 5.0e-6 binsize

periodic boundaries in x and y fixed boundary in z

forces between local and ghost atoms are computed in each processor without communication

communicate velocities for ghost atoms

guarantee that enough atoms are communicated to correctly compute

Shift the lattice grid 0.5 unit so that atom can be created in the center of the grid

define attributes for type 3, set gac diameter to 1e-6, make it consistent with lattice size

set the groups

setting neighbour skin distance and style

rebuild neighbour list if any atom had moved more than half the skin distance

defining grid sytle, substrate names, and grid size

set diffusion boundary conditions and initial concentrations (liquid:kg/m3)

define pair styles

NVE integration with maximum distance limit -> only update position for non gac atoms

heterotrophs growth

heterotrophs division, division diameter: 1.36e-6

EPS secretion from heterotrophs

diffusion reaction fixes

mechanical model fixes

pressure computation

file output

thermo output

issue run command

NUFEB simulation with methanogens and GAC floor

map array: find atoms using indices sort 1000 5.0e-6: sort every 1000 steps with 5.0e-6 binsize

periodic boundaries in x and y fixed boundary in z

forces between local and ghost atoms are computed in each processor without communication

communicate velocities for ghost atoms

guarantee that enough atoms are communicated to correctly compute

Shift the lattice grid 0.5 unit so that atom can be created in the center of the grid

define attributes for type 3, set gac diameter to 1e-6, make it consistent with lattice size

set the groups

setting neighbour skin distance and style

rebuild neighbour list if any atom had moved more than half the skin distance

defining grid sytle, substrate names, and grid size

set diffusion boundary conditions and initial concentrations (liquid:kg/m3)

define pair styles

NVE integration with maximum distance limit -> only update position for non gac atoms

heterotrophs growth

heterotrophs division, division diameter: 1.36e-6

EPS secretion from heterotrophs

diffusion reaction fixes

mechanical model fixes

pressure computation

file output

thermo output

issue run command