ICODE-ADC/lammps
4
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

add sanity check on radscale value, set radius also for quasi-random mode

Browse Source
This commit is contained in:
Axel Kohlmeyer
2022-05-17 15:04:42 -04:00
parent 173e80a970
commit e7d072c593
1 changed files with 6 additions and 2 deletions
Download Patch File Download Diff File Expand all files Collapse all files

8
src/create_atoms.cpp
View File

@ -276,11 +276,12 @@ void CreateAtoms::command(int narg, char **arg)
iarg += 3; iarg += 3;
} else if (strcmp(arg[iarg], "radscale") == 0) { } else if (strcmp(arg[iarg], "radscale") == 0) {
if (iarg + 2 > narg) utils::missing_cmd_args(FLERR, "create_atoms radscale", error); if (iarg + 2 > narg) utils::missing_cmd_args(FLERR, "create_atoms radscale", error);
if ((style != MESH) && (mesh_style != BISECTION)) if (style != MESH)
error->all(FLERR, "Create_atoms radscale can only be used with mesh style in bisect mode"); error->all(FLERR, "Create_atoms radscale can only be used with mesh style");
if (!atom->radius_flag) if (!atom->radius_flag)
error->all(FLERR, "Must have atom attribute radius to set radscale factor"); error->all(FLERR, "Must have atom attribute radius to set radscale factor");
radscale = utils::numeric(FLERR, arg[iarg + 1], false, lmp); radscale = utils::numeric(FLERR, arg[iarg + 1], false, lmp);
if (radscale <= 0.0) error->all(FLERR, "Illegal create_atoms radscale value: {}", radscale);
iarg += 2; iarg += 2;
} else } else
error->all(FLERR, "Illegal create_atoms command option {}", arg[iarg]); error->all(FLERR, "Illegal create_atoms command option {}", arg[iarg]);
@ -975,6 +976,8 @@ int CreateAtoms::add_quasirandom(const double vert[3][3], tagint molid)
MathExtra::cross3(ab, ac, temp); MathExtra::cross3(ab, ac, temp);
area = 0.5 * MathExtra::len3(temp); area = 0.5 * MathExtra::len3(temp);
int nparticles = ceil(mesh_density * area); int nparticles = ceil(mesh_density * area);
// estimate radius from number of particles and area
double rad = sqrt(area/MY_PI/nparticles);
for (int i = 0; i < nparticles; i++) { for (int i = 0; i < nparticles; i++) {
// Define point in unit square // Define point in unit square
@ -1009,6 +1012,7 @@ int CreateAtoms::add_quasirandom(const double vert[3][3], tagint molid)
atom->avec->create_atom(ntype, point); atom->avec->create_atom(ntype, point);
int idx = atom->nlocal - 1; int idx = atom->nlocal - 1;
if (atom->molecule_flag) atom->molecule[idx] = molid; if (atom->molecule_flag) atom->molecule[idx] = molid;
if (atom->radius_flag) atom->radius[idx] = rad * radscale;
} }
} }