Merge pull request #4370 from jtclemm/small-patches

Fixing potential segfaults in granular/nstencil classes and other minor clean ups
This commit is contained in:
Axel Kohlmeyer
2024-11-08 06:05:49 -05:00
committed by GitHub
14 changed files with 21055 additions and 8342 deletions

View File

@ -5,7 +5,11 @@ The BPM package implements bonded particle models which can be used to
simulate mesoscale solids. Solids are constructed as a collection of
particles, which each represent a coarse-grained region of space much
larger than the atomistic scale. Particles within a solid region are
then connected by a network of bonds to provide solid elasticity.
then connected by a network of bonds to model solid elasticity.
There are many names for methods that are based on similar (or
equivalent) capabilities to those in this package, including, but not
limited to, cohesive beam models, bonded DEMs, lattice spring models,
mass spring models, and lattice particle methods.
Unlike traditional bonds in molecular dynamics, the equilibrium bond
length can vary between bonds. Bonds store the reference state. This

File diff suppressed because it is too large Load Diff

View File

@ -11,7 +11,7 @@ mass 1 1.0
mass 2 1.0
velocity all create 2.0 87287 loop geom
pair_style hybrid/overlay lj/cut 2.5 tracker
pair_style hybrid/overlay lj/cut 2.5 tracker trackfix 1000 id1 id2 time/created time/broken time/total time/min 0.5
pair_coeff * * lj/cut 1.0 1.0 2.5
pair_coeff * * tracker 2.5
@ -19,10 +19,9 @@ neighbor 0.3 bin
neigh_modify every 20 delay 0 check no
fix 1 all nve
fix 2 all pair/tracker 1000 id1 id2 time/created time/broken time/total time/min 0.5
fix 3 all ave/histo 1000 1 1000 0 6 30 f_2[5] mode vector kind local file lifetime_hist.dat
fix 2 all ave/histo 1000 1 1000 0 6 30 f_trackfix[5] mode vector kind local file lifetime_hist.dat
dump 1 all local 1000 contact_history.dat f_2[1] f_2[2] f_2[3] f_2[4]
dump 1 all local 1000 contact_history.dat f_trackfix[1] f_trackfix[2] f_trackfix[3] f_trackfix[4]
dump_modify 1 header no
# compute 1 all property/local patom1 patom2

View File

@ -1,4 +1,4 @@
# Histogrammed data for fix 3
# Histogrammed data for fix 2
# TimeStep Number-of-bins Total-counts Missing-counts Min-value Max-value
# Bin Coord Count Count/Total
0 30 0 0 1e+20 -1e+20
@ -32,34 +32,34 @@
28 5.5 0 0
29 5.7 0 0
30 5.9 0 0
1000 30 8122 0 0.5 5
1000 30 528 20285 1 1000
1 0.1 0 0
2 0.3 0 0
3 0.5 910 0.112041
4 0.7 1253 0.154272
5 0.9 953 0.117336
6 1.1 747 0.0919724
7 1.3 559 0.0688254
8 1.5 501 0.0616843
9 1.7 421 0.0518345
10 1.9 356 0.0438316
11 2.1 300 0.0369367
12 2.3 281 0.0345974
13 2.5 242 0.0297956
14 2.7 226 0.0278257
15 2.9 175 0.0215464
16 3.1 168 0.0206846
17 3.3 162 0.0199458
18 3.5 129 0.0158828
19 3.7 151 0.0185915
20 3.9 137 0.0168678
21 4.1 98 0.012066
22 4.3 104 0.0128047
23 4.5 83 0.0102192
24 4.7 77 0.00948042
25 4.9 86 0.0105885
26 5.1 3 0.000369367
3 0.5 0 0
4 0.7 0 0
5 0.9 0 0
6 1.1 17 0.032197
7 1.3 0 0
8 1.5 0 0
9 1.7 0 0
10 1.9 0 0
11 2.1 51 0.0965909
12 2.3 0 0
13 2.5 0 0
14 2.7 0 0
15 2.9 0 0
16 3.1 80 0.151515
17 3.3 0 0
18 3.5 0 0
19 3.7 0 0
20 3.9 0 0
21 4.1 91 0.172348
22 4.3 0 0
23 4.5 0 0
24 4.7 0 0
25 4.9 0 0
26 5.1 117 0.221591
27 5.3 0 0
28 5.5 0 0
29 5.7 0 0
30 5.9 0 0
30 5.9 172 0.325758

View File

@ -0,0 +1,108 @@
LAMMPS (29 Aug 2024 - Development - patch_29Aug2024-731-g5ce635757f)
# 3d Lennard-Jones melt with tracking
units lj
atom_style atomic
lattice fcc 0.8442
Lattice spacing in x,y,z = 1.6795962 1.6795962 1.6795962
region box block 0 4 0 4 0 4
create_box 2 box
Created orthogonal box = (0 0 0) to (6.7183848 6.7183848 6.7183848)
1 by 1 by 1 MPI processor grid
create_atoms 2 box
Created 256 atoms
using lattice units in orthogonal box = (0 0 0) to (6.7183848 6.7183848 6.7183848)
create_atoms CPU = 0.000 seconds
mass 1 1.0
mass 2 1.0
velocity all create 2.0 87287 loop geom
pair_style hybrid/overlay lj/cut 2.5 tracker trackfix 1000 id1 id2 time/created time/broken time/total time/min 0.5
pair_coeff * * lj/cut 1.0 1.0 2.5
pair_coeff * * tracker 2.5
neighbor 0.3 bin
neigh_modify every 20 delay 0 check no
fix 1 all nve
fix 2 all ave/histo 1000 1 1000 0 6 30 f_trackfix[5] mode vector kind local file lifetime_hist.dat
dump 1 all local 1000 contact_history.dat f_trackfix[1] f_trackfix[2] f_trackfix[3] f_trackfix[4]
dump_modify 1 header no
# compute 1 all property/local patom1 patom2
# dump 2 all local 1 pairs.dat c_1[1] c_1[2]
# dump 3 all atom 1 atoms.dat
thermo 50
run 1000 #0
Generated 0 of 1 mixed pair_coeff terms from geometric mixing rule
Neighbor list info ...
update: every = 20 steps, delay = 0 steps, check = no
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 2.8
ghost atom cutoff = 2.8
binsize = 1.4, bins = 5 5 5
2 neighbor lists, perpetual/occasional/extra = 2 0 0
(1) pair lj/cut, perpetual
attributes: half, newton on
pair build: half/bin/atomonly/newton
stencil: half/bin/3d
bin: standard
(2) pair tracker, perpetual, copy from (1)
attributes: half, newton on
pair build: copy
stencil: none
bin: none
Per MPI rank memory allocation (min/avg/max) = 8.679 | 8.679 | 8.679 Mbytes
Step Temp E_pair E_mol TotEng Press
0 2 -6.7733681 0 -3.7850868 -4.5535126
50 1.0540428 -5.3670859 0 -3.7921977 2.3386375
100 1.0402713 -5.3493439 0 -3.7950323 2.4553748
150 1.0570745 -5.3738961 0 -3.7944781 2.3767396
200 1.0431846 -5.3518647 0 -3.7932002 2.5010135
250 1.070121 -5.3902744 0 -3.791363 2.4908658
300 1.0667723 -5.3866302 0 -3.7927224 2.3589344
350 1.000601 -5.2859643 0 -3.7909257 2.9065274
400 0.99256113 -5.2738812 0 -3.7908553 2.8595867
450 1.0482542 -5.357452 0 -3.7912128 2.4707397
500 1.0196176 -5.3123538 0 -3.7889017 2.7230338
550 0.98274535 -5.2586303 0 -3.7902706 2.9156947
600 1.0683914 -5.3863229 0 -3.789996 2.3002719
650 1.0130779 -5.303917 0 -3.7902362 2.8726423
700 1.0583333 -5.3737358 0 -3.792437 2.5770307
750 0.98274506 -5.2612464 0 -3.792887 2.9447027
800 1.0294191 -5.332001 0 -3.7939042 2.5293193
850 0.99240027 -5.2735754 0 -3.7907899 2.7672711
900 1.0293488 -5.3306241 0 -3.7926323 2.6054041
950 0.97137182 -5.2424403 0 -3.7910742 3.129989
1000 1.0009431 -5.2864286 0 -3.7908788 2.7536598
Loop time of 0.286151 on 1 procs for 1000 steps with 256 atoms
Performance: 1509690.552 tau/day, 3494.654 timesteps/s, 894.631 katom-step/s
99.4% CPU use with 1 MPI tasks x no OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 0.21578 | 0.21578 | 0.21578 | 0.0 | 75.41
Neigh | 0.020262 | 0.020262 | 0.020262 | 0.0 | 7.08
Comm | 0.0051679 | 0.0051679 | 0.0051679 | 0.0 | 1.81
Output | 0.021204 | 0.021204 | 0.021204 | 0.0 | 7.41
Modify | 0.02307 | 0.02307 | 0.02307 | 0.0 | 8.06
Other | | 0.0006654 | | | 0.23
Nlocal: 256 ave 256 max 256 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost: 1327 ave 1327 max 1327 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs: 9612 ave 9612 max 9612 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Total # of neighbors = 9612
Ave neighs/atom = 37.546875
Neighbor list builds = 50
Dangerous builds not checked
Total wall time: 0:00:00

View File

@ -1,107 +0,0 @@
LAMMPS (8 Apr 2021)
# 3d Lennard-Jones melt with tracking
units lj
atom_style atomic
lattice fcc 0.8442
Lattice spacing in x,y,z = 1.6795962 1.6795962 1.6795962
region box block 0 4 0 4 0 4
create_box 2 box
Created orthogonal box = (0.0000000 0.0000000 0.0000000) to (6.7183848 6.7183848 6.7183848)
1 by 1 by 1 MPI processor grid
create_atoms 2 box
Created 256 atoms
create_atoms CPU = 0.001 seconds
mass 1 1.0
mass 2 1.0
velocity all create 2.0 87287 loop geom
pair_style hybrid/overlay lj/cut 2.5 tracker
pair_coeff * * lj/cut 1.0 1.0 2.5
pair_coeff * * tracker 2.5
neighbor 0.3 bin
neigh_modify every 20 delay 0 check no
fix 1 all nve
fix 2 all pair/tracker 1000 id1 id2 time/created time/broken time/total time/min 0.5
fix 3 all ave/histo 1000 1 1000 0 6 30 f_2[5] mode vector kind local file lifetime_hist.dat
dump 1 all local 1000 contact_history.dat f_2[1] f_2[2] f_2[3] f_2[4]
dump_modify 1 header no
# compute 1 all property/local patom1 patom2
# dump 2 all local 1 pairs.dat c_1[1] c_1[2]
# dump 3 all atom 1 atoms.dat
thermo 50
run 1000 #0
Neighbor list info ...
update every 20 steps, delay 0 steps, check no
max neighbors/atom: 2000, page size: 100000
master list distance cutoff = 2.8
ghost atom cutoff = 2.8
binsize = 1.4, bins = 5 5 5
2 neighbor lists, perpetual/occasional/extra = 2 0 0
(1) pair lj/cut, perpetual
attributes: half, newton on
pair build: half/bin/atomonly/newton
stencil: half/bin/3d
bin: standard
(2) pair tracker, perpetual
attributes: half, newton on, history
pair build: half/bin/atomonly/newton
stencil: half/bin/3d
bin: standard
Per MPI rank memory allocation (min/avg/max) = 10.46 | 10.46 | 10.46 Mbytes
Step Temp E_pair E_mol TotEng Press
0 2 -6.7733681 0 -3.7850868 -4.5535126
50 1.0540428 -5.3670859 0 -3.7921977 2.3386375
100 1.0402713 -5.3493439 0 -3.7950323 2.4553748
150 1.0570745 -5.3738961 0 -3.7944781 2.3767396
200 1.0431846 -5.3518647 0 -3.7932002 2.5010135
250 1.070121 -5.3902744 0 -3.791363 2.4908658
300 1.0667723 -5.3866302 0 -3.7927224 2.3589344
350 1.000601 -5.2859643 0 -3.7909257 2.9065274
400 0.99256113 -5.2738812 0 -3.7908553 2.8595867
450 1.0482542 -5.357452 0 -3.7912128 2.4707397
500 1.0196176 -5.3123538 0 -3.7889017 2.7230338
550 0.98274535 -5.2586303 0 -3.7902706 2.9156947
600 1.0683914 -5.3863229 0 -3.789996 2.3002719
650 1.0130779 -5.303917 0 -3.7902362 2.8726423
700 1.0583333 -5.3737358 0 -3.792437 2.5770307
750 0.98274506 -5.2612464 0 -3.792887 2.9447027
800 1.0294191 -5.332001 0 -3.7939042 2.5293193
850 0.99240027 -5.2735754 0 -3.7907899 2.7672711
900 1.0293488 -5.3306241 0 -3.7926323 2.6054041
950 0.97137182 -5.2424403 0 -3.7910742 3.129989
1000 1.0009431 -5.2864286 0 -3.7908788 2.7536598
Loop time of 0.310363 on 1 procs for 1000 steps with 256 atoms
Performance: 1391918.252 tau/day, 3222.033 timesteps/s
100.7% CPU use with 1 MPI tasks x no OpenMP threads
MPI task timing breakdown:
Section | min time | avg time | max time |%varavg| %total
---------------------------------------------------------------
Pair | 0.22132 | 0.22132 | 0.22132 | 0.0 | 71.31
Neigh | 0.03458 | 0.03458 | 0.03458 | 0.0 | 11.14
Comm | 0.0087938 | 0.0087938 | 0.0087938 | 0.0 | 2.83
Output | 0.014075 | 0.014075 | 0.014075 | 0.0 | 4.54
Modify | 0.02818 | 0.02818 | 0.02818 | 0.0 | 9.08
Other | | 0.003411 | | | 1.10
Nlocal: 256.000 ave 256 max 256 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Nghost: 1327.00 ave 1327 max 1327 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Neighs: 9612.00 ave 9612 max 9612 min
Histogram: 1 0 0 0 0 0 0 0 0 0
Total # of neighbors = 9612
Ave neighs/atom = 37.546875
Neighbor list builds = 50
Dangerous builds not checked
Total wall time: 0:00:00

View File

@ -146,11 +146,17 @@ void BondBPM::init_style()
"weights = 1,1,1");
if (id_fix_dummy_special && break_flag) {
id_fix_update_special_bonds = utils::strdup("BPM_UPDATE_SPECIAL_BONDS");
auto newfix = modify->replace_fix(
// check if an update fix already exists, if so use it
auto fixes = modify->get_fix_by_style("UPDATE_SPECIAL_BONDS");
if (fixes.size() > 0 ) {
fix_update_special_bonds = dynamic_cast<FixUpdateSpecialBonds *>(fixes[0]);
} else {
id_fix_update_special_bonds = utils::strdup("BPM_UPDATE_SPECIAL_BONDS");
auto newfix = modify->replace_fix(
id_fix_dummy_special,
fmt::format("{} all UPDATE_SPECIAL_BONDS", id_fix_update_special_bonds), 1);
fix_update_special_bonds = dynamic_cast<FixUpdateSpecialBonds *>(newfix);
fix_update_special_bonds = dynamic_cast<FixUpdateSpecialBonds *>(newfix);
}
delete[] id_fix_dummy_special;
id_fix_dummy_special = nullptr;
}
@ -162,10 +168,12 @@ void BondBPM::init_style()
error->all(FLERR, "Bond style bpm requires 1-3 and 1-4 special weights of 1.0");
}
if (force->angle || force->dihedral || force->improper)
error->all(FLERR, "Bond style bpm cannot be used with 3,4-body interactions");
if (atom->molecular == 2)
error->all(FLERR, "Bond style bpm cannot be used with atom style template");
if (break_flag) {
if (force->angle || force->dihedral || force->improper)
error->all(FLERR, "Bond style bpm cannot break with 3,4-body interactions");
if (atom->molecular == 2)
error->all(FLERR, "Bond style bpm cannot break with atom style template");
}
// find all instances of bond history to delete/shift data
// (bond hybrid may create multiple)

View File

@ -454,14 +454,12 @@ void BondBPMRotational::damping_forces(int i1, int i2, int type, double *rhat, d
void BondBPMRotational::compute(int eflag, int vflag)
{
if (!fix_bond_history->stored_flag) {
fix_bond_history->stored_flag = true;
store_data();
}
if (hybrid_flag)
fix_bond_history->compress_history();
if (hybrid_flag) fix_bond_history->compress_history();
int i1, i2, itmp, n, type;
double r[3], r0[3], rhat[3];
@ -544,33 +542,36 @@ void BondBPMRotational::compute(int eflag, int vflag)
// Apply forces and torques to particles
// ------------------------------------------------------//
if (newton_bond || i1 < nlocal) {
f[i1][0] -= force1on2[0] * smooth;
f[i1][1] -= force1on2[1] * smooth;
f[i1][2] -= force1on2[2] * smooth;
MathExtra::scale3(smooth, force1on2);
torque[i1][0] += torque2on1[0] * smooth;
torque[i1][1] += torque2on1[1] * smooth;
torque[i1][2] += torque2on1[2] * smooth;
if (newton_bond || i1 < nlocal) {
f[i1][0] -= force1on2[0];
f[i1][1] -= force1on2[1];
f[i1][2] -= force1on2[2];
MathExtra::scale3(smooth, torque2on1);
torque[i1][0] += torque2on1[0];
torque[i1][1] += torque2on1[1];
torque[i1][2] += torque2on1[2];
}
if (newton_bond || i2 < nlocal) {
f[i2][0] += force1on2[0] * smooth;
f[i2][1] += force1on2[1] * smooth;
f[i2][2] += force1on2[2] * smooth;
f[i2][0] += force1on2[0];
f[i2][1] += force1on2[1];
f[i2][2] += force1on2[2];
torque[i2][0] += torque1on2[0] * smooth;
torque[i2][1] += torque1on2[1] * smooth;
torque[i2][2] += torque1on2[2] * smooth;
MathExtra::scale3(smooth, torque1on2);
torque[i2][0] += torque1on2[0];
torque[i2][1] += torque1on2[1];
torque[i2][2] += torque1on2[2];
}
if (evflag)
ev_tally_xyz(i1, i2, nlocal, newton_bond, 0.0, -force1on2[0] * smooth, -force1on2[1] * smooth,
-force1on2[2] * smooth, r[0], r[1], r[2]);
ev_tally_xyz(i1, i2, nlocal, newton_bond, 0.0, -force1on2[0], -force1on2[1],
-force1on2[2], r[0], r[1], r[2]);
}
if (hybrid_flag)
fix_bond_history->uncompress_history();
if (hybrid_flag) fix_bond_history->uncompress_history();
}
/* ---------------------------------------------------------------------- */
@ -817,21 +818,22 @@ double BondBPMRotational::single(int type, double rsq, int i, int j, double &ffo
// set single_extra quantities
MathExtra::scale3(smooth, force1on2);
svector[0] = r0_mag;
if (flipped) {
svector[1] = -r0[0];
svector[2] = -r0[1];
svector[3] = -r0[2];
svector[4] = force1on2[0] * smooth;
svector[5] = force1on2[1] * smooth;
svector[6] = force1on2[2] * smooth;
svector[4] = force1on2[0];
svector[5] = force1on2[1];
svector[6] = force1on2[2];
} else {
svector[1] = r0[0];
svector[2] = r0[1];
svector[3] = r0[2];
svector[4] = -force1on2[0] * smooth;
svector[5] = -force1on2[1] * smooth;
svector[6] = -force1on2[2] * smooth;
svector[4] = -force1on2[0];
svector[5] = -force1on2[1];
svector[6] = -force1on2[2];
}
return 0.0;

View File

@ -149,8 +149,6 @@ void BondBPMSpring::store_data()
void BondBPMSpring::compute(int eflag, int vflag)
{
if (hybrid_flag) fix_bond_history->compress_history();
int i, bond_change_flag;
double *vol0, *vol;
@ -195,6 +193,8 @@ void BondBPMSpring::compute(int eflag, int vflag)
}
}
if (hybrid_flag) fix_bond_history->compress_history();
int i1, i2, itmp, n, type;
double delx, dely, delz, delvx, delvy, delvz;
double e, rsq, r, r0, rinv, smooth, fbond, dot;

View File

@ -20,8 +20,6 @@
#include <cmath>
#include <cmath>
using namespace LAMMPS_NS;
using namespace Granular_NS;
using namespace MathConst;
@ -159,6 +157,8 @@ GranSubModDampingCoeffRestitution::GranSubModDampingCoeffRestitution(GranularMod
allow_cohesion = 0;
}
/* ---------------------------------------------------------------------- */
void GranSubModDampingCoeffRestitution::init()
{
// Calculate prefactor, assume Hertzian as default
@ -172,6 +172,8 @@ void GranSubModDampingCoeffRestitution::init()
}
}
/* ---------------------------------------------------------------------- */
double GranSubModDampingCoeffRestitution::calculate_forces()
{
// in case argument <= 0 due to precision issues

View File

@ -12,6 +12,7 @@
------------------------------------------------------------------------- */
#include "gran_sub_mod_normal.h"
#include "error.h"
#include "granular_model.h"
#include "math_const.h"

View File

@ -77,10 +77,11 @@ void GranSubModRollingSDS::calculate_forces()
Frcrit = mu * gm->normal_model->get_fncrit();
hist_temp[0] = gm->history[rhist0];
hist_temp[1] = gm->history[rhist1];
hist_temp[2] = gm->history[rhist2];
if (gm->history_update) {
hist_temp[0] = gm->history[rhist0];
hist_temp[1] = gm->history[rhist1];
hist_temp[2] = gm->history[rhist2];
rolldotn = dot3(hist_temp, gm->nx);
frameupdate = (fabs(rolldotn) * k) > (EPSILON * Frcrit);

View File

@ -231,13 +231,13 @@ void PairSPHTaitwater::coeff(int narg, char **arg)
allocate();
int ilo, ihi, jlo, jhi;
utils::bounds(FLERR,arg[0], 1, atom->ntypes, ilo, ihi, error);
utils::bounds(FLERR,arg[1], 1, atom->ntypes, jlo, jhi, error);
utils::bounds(FLERR, arg[0], 1, atom->ntypes, ilo, ihi, error);
utils::bounds(FLERR, arg[1], 1, atom->ntypes, jlo, jhi, error);
double rho0_one = utils::numeric(FLERR,arg[2],false,lmp);
double soundspeed_one = utils::numeric(FLERR,arg[3],false,lmp);
double viscosity_one = utils::numeric(FLERR,arg[4],false,lmp);
double cut_one = utils::numeric(FLERR,arg[5],false,lmp);
double rho0_one = utils::numeric(FLERR, arg[2], false, lmp);
double soundspeed_one = utils::numeric(FLERR, arg[3], false, lmp);
double viscosity_one = utils::numeric(FLERR, arg[4], false, lmp);
double cut_one = utils::numeric(FLERR, arg[5], false, lmp);
double B_one = soundspeed_one * soundspeed_one * rho0_one / 7.0;
int count = 0;
@ -254,7 +254,7 @@ void PairSPHTaitwater::coeff(int narg, char **arg)
}
if (count == 0)
error->all(FLERR,"Incorrect args for pair coefficients");
error->all(FLERR, "Incorrect args for pair coefficients");
}
/* ----------------------------------------------------------------------
@ -264,7 +264,7 @@ void PairSPHTaitwater::coeff(int narg, char **arg)
double PairSPHTaitwater::init_one(int i, int j)
{
if (setflag[i][j] == 0) {
error->all(FLERR,"All pair sph/taitwater coeffs are set");
error->all(FLERR, "All pair sph/taitwater coeffs are not set");
}
cut[j][i] = cut[i][j];

View File

@ -223,15 +223,15 @@ void NStencil::create_setup()
// smax = max possible size of entire 3d stencil
// stencil will be empty if cutneighmax = 0.0
sx = static_cast<int>(cutneighmax*bininvx);
if (sx*binsizex < cutneighmax) sx++;
sy = static_cast<int>(cutneighmax*bininvy);
if (sy*binsizey < cutneighmax) sy++;
sz = static_cast<int>(cutneighmax*bininvz);
if (sz*binsizez < cutneighmax) sz++;
sx = static_cast<int>(cutneighmax * bininvx);
if (sx * binsizex < cutneighmax) sx++;
sy = static_cast<int>(cutneighmax * bininvy);
if (sy * binsizey < cutneighmax) sy++;
sz = static_cast<int>(cutneighmax * bininvz);
if (sz * binsizez < cutneighmax) sz++;
if (dimension == 2) sz = 0;
int smax = (2*sx+1) * (2*sy+1) * (2*sz+1);
int smax = (2 * sx + 1) * (2 * sy + 1) * (2 * sz + 1);
// reallocate stencil structs if necessary
// for BIN and MULTI_OLD styles
@ -240,10 +240,10 @@ void NStencil::create_setup()
if (smax > maxstencil) {
maxstencil = smax;
memory->destroy(stencil);
memory->create(stencil,maxstencil,"neighstencil:stencil");
memory->create(stencil, maxstencil, "neighstencil:stencil");
if (xyzflag) {
memory->destroy(stencilxyz);
memory->create(stencilxyz,maxstencil,3,"neighstencil:stencilxyz");
memory->create(stencilxyz, maxstencil, 3, "neighstencil:stencilxyz");
}
}
@ -251,9 +251,9 @@ void NStencil::create_setup()
int i;
int n = atom->ntypes;
if (maxstencil_multi_old == 0) {
nstencil_multi_old = new int[n+1];
stencil_multi_old = new int*[n+1];
distsq_multi_old = new double*[n+1];
nstencil_multi_old = new int[n + 1];
stencil_multi_old = new int*[n + 1];
distsq_multi_old = new double*[n + 1];
for (i = 1; i <= n; i++) {
nstencil_multi_old[i] = 0;
stencil_multi_old[i] = nullptr;
@ -265,9 +265,9 @@ void NStencil::create_setup()
for (i = 1; i <= n; i++) {
memory->destroy(stencil_multi_old[i]);
memory->destroy(distsq_multi_old[i]);
memory->create(stencil_multi_old[i],maxstencil_multi_old,
memory->create(stencil_multi_old[i], maxstencil_multi_old,
"neighstencil:stencil_multi_old");
memory->create(distsq_multi_old[i],maxstencil_multi_old,
memory->create(distsq_multi_old[i], maxstencil_multi_old,
"neighstencil:distsq_multi_old");
}
}
@ -280,7 +280,7 @@ void NStencil::create_setup()
if (nb) copy_bin_info_multi();
// Deallocate arrays if previously allocated
if((n > maxcollections) && stencil_multi){
if ((n > maxcollections) && stencil_multi){
memory->destroy(nstencil_multi);
for (i = 0; i < maxcollections; i++) {
for (j = 0; j < maxcollections; j++)
@ -339,10 +339,10 @@ void NStencil::create_setup()
memory->create(maxstencil_multi, n, n, "neighstencil::maxstencil_multi");
memory->create(nstencil_multi, n, n, "neighstencil::nstencil_multi");
stencil_multi = new int**[n]();
for (i = 0; i < n; ++i) {
for (i = 0; i < n; i++) {
stencil_multi[i] = new int*[n]();
for (j = 0; j < n; ++j) {
maxstencil_multi[i][j] = 0;
for (j = 0; j < n; j++) {
maxstencil_multi[i][j] = 0;
nstencil_multi[i][j] = 0;
stencil_multi[i][j] = nullptr;
}
@ -378,23 +378,23 @@ void NStencil::create_setup()
stencil_range = sqrt(cutcollectionsq[i][j]);
sx = static_cast<int> (stencil_range*bininvx_multi[bin_collection]);
if (sx*binsizex_multi[bin_collection] < stencil_range) sx++;
sy = static_cast<int> (stencil_range*bininvy_multi[bin_collection]);
if (sy*binsizey_multi[bin_collection] < stencil_range) sy++;
sz = static_cast<int> (stencil_range*bininvz_multi[bin_collection]);
if (sz*binsizez_multi[bin_collection] < stencil_range) sz++;
sx = static_cast<int> (stencil_range * bininvx_multi[bin_collection]);
if (sx * binsizex_multi[bin_collection] < stencil_range) sx++;
sy = static_cast<int> (stencil_range * bininvy_multi[bin_collection]);
if (sy * binsizey_multi[bin_collection] < stencil_range) sy++;
sz = static_cast<int> (stencil_range * bininvz_multi[bin_collection]);
if (sz * binsizez_multi[bin_collection] < stencil_range) sz++;
if (dimension == 2) sz = 0;
stencil_sx_multi[i][j] = sx;
stencil_sy_multi[i][j] = sy;
stencil_sz_multi[i][j] = sz;
smax = ((2*sx+1) * (2*sy+1) * (2*sz+1));
smax = ((2 * sx + 1) * (2 * sy + 1) * (2 * sz + 1));
if (smax > maxstencil_multi[i][j]) {
maxstencil_multi[i][j] = smax;
if(stencil_multi[i][j])
if (stencil_multi[i][j])
memory->destroy(stencil_multi[i][j]);
memory->create(stencil_multi[i][j], smax,
"neighstencil::stencil_multi");
@ -457,16 +457,20 @@ double NStencil::memory_usage()
{
double bytes = 0;
if (neighstyle == Neighbor::BIN) {
bytes += memory->usage(stencil,maxstencil);
bytes += memory->usage(stencilxyz,maxstencil,3);
bytes += memory->usage(stencil, maxstencil);
bytes += memory->usage(stencilxyz, maxstencil, 3);
} else if (neighstyle == Neighbor::MULTI_OLD) {
bytes += (double)atom->ntypes*maxstencil_multi_old * sizeof(int);
bytes += (double)atom->ntypes*maxstencil_multi_old * sizeof(double);
bytes += (double)atom->ntypes * maxstencil_multi_old * sizeof(int);
bytes += (double)atom->ntypes * maxstencil_multi_old * sizeof(double);
} else if (neighstyle == Neighbor::MULTI) {
int n = ncollections;
for (int i = 0; i < n; i++) {
for (int j = 0; j < n; j++) {
bytes += (double)maxstencil_multi[i][j] * sizeof(int);
// May not yet be initialized for occasional lists
if (maxstencil_multi) {
for (int i = 0; i < n; i++) {
for (int j = 0; j < n; j++) {
bytes += (double)maxstencil_multi[i][j] * sizeof(int);
}
}
}
}