rename "set spin" to "set spin/atom" and "set spin/random" to "set spin/atom/random"

this also adds unit tests and improves error checking and messages
This commit is contained in:
Axel Kohlmeyer
2022-08-28 18:57:19 -04:00
parent bb6ea422b7
commit 74bbbb9e17
3 changed files with 178 additions and 21 deletions

View File

@ -15,7 +15,7 @@ Syntax
* one or more keyword/value pairs may be appended * one or more keyword/value pairs may be appended
* keyword = *type* or *type/fraction* or *type/ratio* or *type/subset* * keyword = *type* or *type/fraction* or *type/ratio* or *type/subset*
or *mol* or *x* or *y* or *z* or *vx* or *vy* or *vz* or *charge* or or *mol* or *x* or *y* or *z* or *vx* or *vy* or *vz* or *charge* or
*dipole* or *dipole/random* or *quat* or *spin* or *spin/random* or *dipole* or *dipole/random* or *quat* or *spin/atom* or *spin/atom/random* or
*quat* or *quat/random* or *diameter* or *shape* or *quat* or *quat/random* or *diameter* or *shape* or
*length* or *tri* or *theta* or *theta/random* or *angmom* or *length* or *tri* or *theta* or *theta/random* or *angmom* or
*omega* or *mass* or *density* or *density/disc* or *omega* or *mass* or *density* or *density/disc* or
@ -55,11 +55,11 @@ Syntax
*dipole/random* value = seed Dlen *dipole/random* value = seed Dlen
seed = random # seed (positive integer) for dipole moment orientations seed = random # seed (positive integer) for dipole moment orientations
Dlen = magnitude of dipole moment (dipole units) Dlen = magnitude of dipole moment (dipole units)
*spin* values = g x y z *spin/atom* values = g x y z
g = magnitude of magnetic spin vector (in Bohr magneton's unit) g = magnitude of magnetic spin vector (in Bohr magneton's unit)
x,y,z = orientation of magnetic spin vector x,y,z = orientation of magnetic spin vector
any of x,y,z can be an atom-style variable (see below) any of x,y,z can be an atom-style variable (see below)
*spin/random* value = seed Dlen *spin/atom/random* value = seed Dlen
seed = random # seed (positive integer) for magnetic spin orientations seed = random # seed (positive integer) for magnetic spin orientations
Dlen = magnitude of magnetic spin vector (in Bohr magneton's unit) Dlen = magnitude of magnetic spin vector (in Bohr magneton's unit)
*quat* values = a b c theta *quat* values = a b c theta
@ -277,14 +277,18 @@ the orientation of a particular atom is the same, regardless of how
many processors are being used. This keyword does not allow use of an many processors are being used. This keyword does not allow use of an
atom-style variable. atom-style variable.
Keyword *spin* uses the specified g value to set the magnitude of the .. versionchanged:: TBD
Keyword *spin/atom* uses the specified g value to set the magnitude of the
magnetic spin vectors, and the x,y,z values as components of a vector magnetic spin vectors, and the x,y,z values as components of a vector
to set as the orientation of the magnetic spin vectors of the selected to set as the orientation of the magnetic spin vectors of the selected
atoms. atoms. This keyword was previously called *spin*.
Keyword *spin/random* randomizes the orientation of the magnetic spin .. versionchanged:: TBD
Keyword *spin/atom/random* randomizes the orientation of the magnetic spin
vectors for the selected atoms and sets the magnitude of each to the vectors for the selected atoms and sets the magnitude of each to the
specified *Dlen* value. specified *Dlen* value. This keyword was previously called *spin/random*.
Keyword *quat* uses the specified values to create a quaternion Keyword *quat* uses the specified values to create a quaternion
(4-vector) that represents the orientation of the selected atoms. The (4-vector) that represents the orientation of the selected atoms. The

View File

@ -46,9 +46,9 @@ enum{ATOM_SELECT,MOL_SELECT,TYPE_SELECT,GROUP_SELECT,REGION_SELECT};
enum{TYPE,TYPE_FRACTION,TYPE_RATIO,TYPE_SUBSET, enum{TYPE,TYPE_FRACTION,TYPE_RATIO,TYPE_SUBSET,
MOLECULE,X,Y,Z,VX,VY,VZ,CHARGE,MASS,SHAPE,LENGTH,TRI, MOLECULE,X,Y,Z,VX,VY,VZ,CHARGE,MASS,SHAPE,LENGTH,TRI,
DIPOLE,DIPOLE_RANDOM,SPIN,SPIN_RANDOM,QUAT,QUAT_RANDOM, DIPOLE,DIPOLE_RANDOM,SPIN_ATOM,SPIN_RANDOM,SPIN_ELECTRON,RADIUS_ELECTRON,
THETA,THETA_RANDOM,ANGMOM,OMEGA, QUAT,QUAT_RANDOM,THETA,THETA_RANDOM,ANGMOM,OMEGA,
DIAMETER,DENSITY,VOLUME,IMAGE,BOND,ANGLE,DIHEDRAL,IMPROPER, DIAMETER,RADIUS_ATOM,DENSITY,VOLUME,IMAGE,BOND,ANGLE,DIHEDRAL,IMPROPER,
SPH_E,SPH_CV,SPH_RHO,EDPD_TEMP,EDPD_CV,CC,SMD_MASS_DENSITY, SPH_E,SPH_CV,SPH_RHO,EDPD_TEMP,EDPD_CV,CC,SMD_MASS_DENSITY,
SMD_CONTACT_RADIUS,DPDTHETA,EPSILON,IVEC,DVEC,IARRAY,DARRAY}; SMD_CONTACT_RADIUS,DPDTHETA,EPSILON,IVEC,DVEC,IARRAY,DARRAY};
@ -265,8 +265,10 @@ void Set::command(int narg, char **arg)
setrandom(DIPOLE_RANDOM); setrandom(DIPOLE_RANDOM);
iarg += 3; iarg += 3;
} else if (strcmp(arg[iarg],"spin") == 0) { } else if ((strcmp(arg[iarg],"spin") == 0) || (strcmp(arg[iarg],"spin/atom") == 0)) {
if (iarg+4 > narg) error->all(FLERR,"Illegal set command"); if ((strcmp(arg[iarg],"spin") == 0) && (comm->me == 0))
error->warning(FLERR, "Set attribute spin is deprecated. Please use spin/atom instead.");
if (iarg+5 > narg) utils::missing_cmd_args(FLERR, "set spin/atom", error);
if (utils::strmatch(arg[iarg+1],"^v_")) varparse(arg[iarg+1],1); if (utils::strmatch(arg[iarg+1],"^v_")) varparse(arg[iarg+1],1);
else dvalue = utils::numeric(FLERR,arg[iarg+1],false,lmp); else dvalue = utils::numeric(FLERR,arg[iarg+1],false,lmp);
if (utils::strmatch(arg[iarg+2],"^v_")) varparse(arg[iarg+2],2); if (utils::strmatch(arg[iarg+2],"^v_")) varparse(arg[iarg+2],2);
@ -275,24 +277,50 @@ void Set::command(int narg, char **arg)
else yvalue = utils::numeric(FLERR,arg[iarg+3],false,lmp); else yvalue = utils::numeric(FLERR,arg[iarg+3],false,lmp);
if (utils::strmatch(arg[iarg+4],"^v_")) varparse(arg[iarg+4],4); if (utils::strmatch(arg[iarg+4],"^v_")) varparse(arg[iarg+4],4);
else zvalue = utils::numeric(FLERR,arg[iarg+4],false,lmp); else zvalue = utils::numeric(FLERR,arg[iarg+4],false,lmp);
if ((xvalue == 0.0) && (yvalue == 0.0) && (zvalue == 0.0))
error->all(FLERR,"At least one spin vector component must be non-zero");
if (!atom->sp_flag) if (!atom->sp_flag)
error->all(FLERR,"Cannot set this attribute for this atom style"); error->all(FLERR,"Cannot set attribute {} for atom style {}", arg[iarg], atom->get_style());
set(SPIN); if (dvalue <= 0.0)
error->all(FLERR,"Invalid spin magnitude {} in set {} command", dvalue, arg[iarg]);
set(SPIN_ATOM);
iarg += 5; iarg += 5;
} else if (strcmp(arg[iarg],"spin/random") == 0) { } else if ((strcmp(arg[iarg],"spin/random") == 0) ||
if (iarg+3 > narg) error->all(FLERR,"Illegal set command"); (strcmp(arg[iarg],"spin/atom/random") == 0)) {
if (iarg+3 > narg) utils::missing_cmd_args(FLERR, "set spin/atom/random", error);
ivalue = utils::inumeric(FLERR,arg[iarg+1],false,lmp); ivalue = utils::inumeric(FLERR,arg[iarg+1],false,lmp);
dvalue = utils::numeric(FLERR,arg[iarg+2],false,lmp); dvalue = utils::numeric(FLERR,arg[iarg+2],false,lmp);
if ((strcmp(arg[iarg],"spin/random") == 0) && (comm->me == 0))
error->warning(FLERR, "Set attribute spin/random is deprecated. "
"Please use spin/atom/random instead.");
if (!atom->sp_flag) if (!atom->sp_flag)
error->all(FLERR,"Cannot set this attribute for this atom style"); error->all(FLERR,"Cannot set attribute {} for atom style {}", arg[iarg], atom->get_style());
if (ivalue <= 0) if (ivalue <= 0)
error->all(FLERR,"Invalid random number seed in set command"); error->all(FLERR,"Invalid random number seed {} in set {} command", ivalue, arg[iarg]);
if (dvalue <= 0.0) if (dvalue <= 0.0)
error->all(FLERR,"Invalid dipole length in set command"); error->all(FLERR,"Invalid spin magnitude {} in set {} command", dvalue, arg[iarg]);
setrandom(SPIN_RANDOM); setrandom(SPIN_RANDOM);
iarg += 3; iarg += 3;
} else if (strcmp(arg[iarg],"radius/electron") == 0) {
if (iarg+2 > narg) utils::missing_cmd_args(FLERR, "set radius/electron", error);
if (utils::strmatch(arg[iarg+1],"^v_")) varparse(arg[iarg+1],1);
else dvalue = utils::numeric(FLERR,arg[iarg+1],false,lmp);
if (!atom->eradius_flag)
error->all(FLERR,"Cannot set attribute {} for atom style {}", arg[iarg], atom->get_style());
set(RADIUS_ELECTRON);
iarg += 2;
} else if (strcmp(arg[iarg],"spin/electron") == 0) {
if (iarg+2 > narg) utils::missing_cmd_args(FLERR, "set spin/electron", error);
if (utils::strmatch(arg[iarg+1],"^v_")) varparse(arg[iarg+1],1);
else dvalue = utils::numeric(FLERR,arg[iarg+1],false,lmp);
if (!atom->spin_flag)
error->all(FLERR,"Cannot set attribute {} for atom style {}", arg[iarg], atom->get_style());
set(SPIN_ELECTRON);
iarg += 2;
} else if (strcmp(arg[iarg],"quat") == 0) { } else if (strcmp(arg[iarg],"quat") == 0) {
if (iarg+5 > narg) utils::missing_cmd_args(FLERR, "set quat", error); if (iarg+5 > narg) utils::missing_cmd_args(FLERR, "set quat", error);
if (utils::strmatch(arg[iarg+1],"^v_")) varparse(arg[iarg+1],1); if (utils::strmatch(arg[iarg+1],"^v_")) varparse(arg[iarg+1],1);
@ -714,7 +742,7 @@ void Set::selection(int n)
} else if (style == GROUP_SELECT) { } else if (style == GROUP_SELECT) {
int igroup = group->find(id); int igroup = group->find(id);
if (igroup == -1) error->all(FLERR,"Could not find set group ID"); if (igroup == -1) error->all(FLERR,"Could not find set group ID {}", id);
int groupbit = group->bitmask[igroup]; int groupbit = group->bitmask[igroup];
int *mask = atom->mask; int *mask = atom->mask;
@ -947,7 +975,7 @@ void Set::set(int keyword)
// set magnetic moments // set magnetic moments
else if (keyword == SPIN) { else if (keyword == SPIN_ATOM) {
double **sp = atom->sp; double **sp = atom->sp;
double inorm = 1.0/sqrt(xvalue*xvalue+yvalue*yvalue+zvalue*zvalue); double inorm = 1.0/sqrt(xvalue*xvalue+yvalue*yvalue+zvalue*zvalue);
sp[i][0] = inorm*xvalue; sp[i][0] = inorm*xvalue;

View File

@ -87,6 +87,11 @@ TEST_F(SetTest, NoBoxNoAtoms)
TEST_FAILURE(".*ERROR: Unknown set command style: xxx.*", command("set xxx 1 x 0.0");); TEST_FAILURE(".*ERROR: Unknown set command style: xxx.*", command("set xxx 1 x 0.0"););
TEST_FAILURE(".*ERROR: Set keyword or custom property yyy does not exist.*", TEST_FAILURE(".*ERROR: Set keyword or custom property yyy does not exist.*",
command("set type 1 yyy 0.0");); command("set type 1 yyy 0.0"););
TEST_FAILURE(".*ERROR: Cannot set attribute spin/atom for atom style atomic.*",
command("set atom * spin/atom 1.0 1.0 0.0 0.0"););
TEST_FAILURE(".*ERROR: Cannot set attribute spin/atom/random for atom style atomic.*",
command("set atom * spin/atom/random 436273456 1.0"););
} }
TEST_F(SetTest, StylesTypes) TEST_F(SetTest, StylesTypes)
@ -161,6 +166,126 @@ TEST_F(SetTest, StylesTypes)
for (int i = 0; i < 8; ++i) for (int i = 0; i < 8; ++i)
sum += (atom->type[i] == 2) ? 1 : 0; sum += (atom->type[i] == 2) ? 1 : 0;
ASSERT_EQ(sum, 4); ASSERT_EQ(sum, 4);
TEST_FAILURE(".*ERROR: Numeric index 9 is out of bounds .1-8.*", command("set type 9 x 0.0"););
TEST_FAILURE(".*ERROR: Invalid range string: 3:10.*", command("set type 3:10 x 0.0"););
TEST_FAILURE(".*ERROR: Could not find set group ID nope.*", command("set group nope x 0.0"););
}
TEST_F(SetTest, PosVelCharge)
{
atomic_system("charge");
ASSERT_EQ(atom->natoms, 8);
BEGIN_HIDE_OUTPUT();
command("set group top charge 1.0");
command("set atom 5*8 charge -1.0");
END_HIDE_OUTPUT();
EXPECT_EQ(atom->q[0], 1);
EXPECT_EQ(atom->q[1], 1);
EXPECT_EQ(atom->q[2], 0);
EXPECT_EQ(atom->q[3], 0);
EXPECT_EQ(atom->q[4], -1);
EXPECT_EQ(atom->q[5], -1);
EXPECT_EQ(atom->q[6], -1);
EXPECT_EQ(atom->q[7], -1);
BEGIN_HIDE_OUTPUT();
command("variable xpos atom 0.5-x");
command("variable ypos atom y*0.5");
command("set atom * x v_xpos y v_ypos z 0.5");
command("set group all vx v_xpos vy v_ypos vz 0.5");
END_HIDE_OUTPUT();
EXPECT_EQ(atom->x[0][0], 0.375);
EXPECT_EQ(atom->x[0][1], 0.0625);
EXPECT_EQ(atom->x[0][2], 0.5);
EXPECT_EQ(atom->x[1][0], -0.625);
EXPECT_EQ(atom->x[1][1], 0.0625);
EXPECT_EQ(atom->x[1][2], 0.5);
EXPECT_EQ(atom->x[2][0], 0.375);
EXPECT_EQ(atom->x[2][1], 0.5625);
EXPECT_EQ(atom->x[2][2], 0.5);
EXPECT_EQ(atom->x[3][0], -0.625);
EXPECT_EQ(atom->x[3][1], 0.5625);
EXPECT_EQ(atom->x[3][2], 0.5);
EXPECT_EQ(atom->x[4][0], 0.375);
EXPECT_EQ(atom->x[4][1], 0.0625);
EXPECT_EQ(atom->x[4][2], 0.5);
EXPECT_EQ(atom->x[5][0], -0.625);
EXPECT_EQ(atom->x[5][1], 0.0625);
EXPECT_EQ(atom->x[5][2], 0.5);
EXPECT_EQ(atom->x[6][0], 0.375);
EXPECT_EQ(atom->x[6][1], 0.5625);
EXPECT_EQ(atom->x[6][2], 0.5);
EXPECT_EQ(atom->x[7][0], -0.625);
EXPECT_EQ(atom->x[7][1], 0.5625);
EXPECT_EQ(atom->x[7][2], 0.5);
EXPECT_EQ(atom->v[0][0], 0.125);
EXPECT_EQ(atom->v[0][1], 0.03125);
EXPECT_EQ(atom->v[0][2], 0.5);
EXPECT_EQ(atom->v[1][0], 1.125);
EXPECT_EQ(atom->v[1][1], 0.03125);
EXPECT_EQ(atom->v[1][2], 0.5);
EXPECT_EQ(atom->v[2][0], 0.125);
EXPECT_EQ(atom->v[2][1], 0.28125);
EXPECT_EQ(atom->v[2][2], 0.5);
EXPECT_EQ(atom->v[3][0], 1.125);
EXPECT_EQ(atom->v[3][1], 0.28125);
EXPECT_EQ(atom->v[3][2], 0.5);
EXPECT_EQ(atom->v[4][0], 0.125);
EXPECT_EQ(atom->v[4][1], 0.03125);
EXPECT_EQ(atom->v[4][2], 0.5);
EXPECT_EQ(atom->v[5][0], 1.125);
EXPECT_EQ(atom->v[5][1], 0.03125);
EXPECT_EQ(atom->v[5][2], 0.5);
EXPECT_EQ(atom->v[6][0], 0.125);
EXPECT_EQ(atom->v[6][1], 0.28125);
EXPECT_EQ(atom->v[6][2], 0.5);
EXPECT_EQ(atom->v[7][0], 1.125);
EXPECT_EQ(atom->v[7][1], 0.28125);
EXPECT_EQ(atom->v[7][2], 0.5);
}
TEST_F(SetTest, SpinPackage)
{
if (!Info::has_package("SPIN")) GTEST_SKIP();
atomic_system("spin");
ASSERT_EQ(atom->natoms, 8);
BEGIN_HIDE_OUTPUT();
command("set atom 1*2 spin/atom 0.5 0.1 0.5 -0.1");
command("set atom 8 spin/atom/random 23974 0.25");
END_HIDE_OUTPUT();
constexpr double vx = 0.1;
constexpr double vy = 0.5;
constexpr double vz = -0.1;
constexpr double norm = 1.0 / sqrt(vx * vx + vy * vy + vz * vz);
EXPECT_EQ(atom->sp[0][0], vx * norm);
EXPECT_EQ(atom->sp[0][1], vy * norm);
EXPECT_EQ(atom->sp[0][2], vz * norm);
EXPECT_EQ(atom->sp[0][3], 0.5);
EXPECT_EQ(atom->sp[1][0], vx * norm);
EXPECT_EQ(atom->sp[1][1], vy * norm);
EXPECT_EQ(atom->sp[1][2], vz * norm);
EXPECT_EQ(atom->sp[1][3], 0.5);
EXPECT_NE(atom->sp[7][0], 0.0);
EXPECT_NE(atom->sp[7][1], 0.0);
EXPECT_NE(atom->sp[7][2], 0.0);
EXPECT_EQ(atom->sp[7][3], 0.25);
for (int i = 2; i < 7; ++i)
for (int j = 0; j < 4; ++j)
EXPECT_EQ(atom->sp[i][j], 0);
TEST_FAILURE(".*ERROR: Invalid spin magnitude -0.1 in set spin/atom command.*",
command("set atom * spin/atom -0.1 1.0 0.0 0.0"););
TEST_FAILURE(".*ERROR: At least one spin vector component must be non-zero.*",
command("set atom * spin/atom 1.0 0.0 0.0 0.0"););
TEST_FAILURE(".*ERROR: Invalid spin magnitude -0.2 in set spin/atom/random command.*",
command("set atom * spin/atom/random 436273456 -0.2"););
TEST_FAILURE(".*ERROR: Invalid random number seed 0 in set spin/atom/random command.*",
command("set atom * spin/atom/random 0 1.0"););
} }
} // namespace LAMMPS_NS } // namespace LAMMPS_NS