From cdbec96e0703dab98f16e61d0c6d115baea91844 Mon Sep 17 00:00:00 2001 From: Axel Kohlmeyer Date: Wed, 21 Aug 2024 09:00:11 -0400 Subject: [PATCH] enable and apply clang-format --- src/group.cpp | 727 +++++++++++++++++++++++++------------------------- 1 file changed, 368 insertions(+), 359 deletions(-) diff --git a/src/group.cpp b/src/group.cpp index baf9f377ff..909e741c6b 100644 --- a/src/group.cpp +++ b/src/group.cpp @@ -1,4 +1,3 @@ -// clang-format off /* ---------------------------------------------------------------------- LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator https://www.lammps.org/, Sandia National Laboratories @@ -20,18 +19,18 @@ #include "domain.h" #include "dump.h" #include "error.h" +#include "exceptions.h" #include "fix.h" #include "force.h" #include "input.h" -#include "math_extra.h" #include "math_eigen.h" +#include "math_extra.h" #include "memory.h" #include "modify.h" #include "output.h" #include "region.h" #include "tokenizer.h" #include "variable.h" -#include "exceptions.h" #include #include @@ -43,8 +42,8 @@ using namespace LAMMPS_NS; static constexpr int MAX_GROUP = 32; static constexpr double EPSILON = 1.0e-6; -enum{NONE,TYPE,MOLECULE,ID}; -enum{LT,LE,GT,GE,EQ,NEQ,BETWEEN}; +enum { NONE, TYPE, MOLECULE, ID }; +enum { LT, LE, GT, GE, EQ, NEQ, BETWEEN }; static constexpr double BIG = 1.0e20; @@ -54,9 +53,9 @@ static constexpr double BIG = 1.0e20; Group::Group(LAMMPS *lmp) : Pointers(lmp) { - MPI_Comm_rank(world,&me); + MPI_Comm_rank(world, &me); - names = new char*[MAX_GROUP]; + names = new char *[MAX_GROUP]; bitmask = new int[MAX_GROUP]; inversemask = new int[MAX_GROUP]; dynamic = new int[MAX_GROUP]; @@ -78,11 +77,11 @@ Group::Group(LAMMPS *lmp) : Pointers(lmp) Group::~Group() { - for (int i = 0; i < MAX_GROUP; i++) delete [] names[i]; - delete [] names; - delete [] bitmask; - delete [] inversemask; - delete [] dynamic; + for (int i = 0; i < MAX_GROUP; i++) delete[] names[i]; + delete[] names; + delete[] bitmask; + delete[] inversemask; + delete[] dynamic; } /* ---------------------------------------------------------------------- @@ -93,39 +92,37 @@ void Group::assign(int narg, char **arg) { int i; - if (domain->box_exist == 0) - error->all(FLERR,"Group command before simulation box is defined"); + if (domain->box_exist == 0) error->all(FLERR, "Group command before simulation box is defined"); if (narg < 2) utils::missing_cmd_args(FLERR, "group", error); // delete the group if not being used elsewhere // clear mask of each atom assigned to this group - if (strcmp(arg[1],"delete") == 0) { - if (narg != 2) error->all(FLERR,"Illegal group delete command: too many arguments"); + if (strcmp(arg[1], "delete") == 0) { + if (narg != 2) error->all(FLERR, "Illegal group delete command: too many arguments"); int igroup = find(arg[0]); - if (igroup == -1) error->all(FLERR,"Could not find group delete group ID {}",arg[0]); - if (igroup == 0) error->all(FLERR,"Cannot delete group all"); + if (igroup == -1) error->all(FLERR, "Could not find group delete group ID {}", arg[0]); + if (igroup == 0) error->all(FLERR, "Cannot delete group all"); for (const auto &i : modify->get_fix_list()) if (i->igroup == igroup) - error->all(FLERR,"Cannot delete group {} currently used by fix ID {}",arg[0],i->id); + error->all(FLERR, "Cannot delete group {} currently used by fix ID {}", arg[0], i->id); for (const auto &i : modify->get_compute_list()) if (i->igroup == igroup) - error->all(FLERR,"Cannot delete group {} currently used by compute ID {}",arg[0],i->id); + error->all(FLERR, "Cannot delete group {} currently used by compute ID {}", arg[0], i->id); for (const auto &i : output->get_dump_list()) if (i->igroup == igroup) - error->all(FLERR,"Cannot delete group {} currently used by dump ID {}",arg[0],i->id); - if (atom->firstgroupname && strcmp(arg[0],atom->firstgroupname) == 0) - error->all(FLERR,"Cannot delete group {} currently used by atom_modify first",arg[0]); + error->all(FLERR, "Cannot delete group {} currently used by dump ID {}", arg[0], i->id); + if (atom->firstgroupname && strcmp(arg[0], atom->firstgroupname) == 0) + error->all(FLERR, "Cannot delete group {} currently used by atom_modify first", arg[0]); int *mask = atom->mask; int nlocal = atom->nlocal; int bits = inversemask[igroup]; for (i = 0; i < nlocal; i++) mask[i] &= bits; - if (dynamic[igroup]) - modify->delete_fix(std::string("GROUP_") + names[igroup]); + if (dynamic[igroup]) modify->delete_fix(std::string("GROUP_") + names[igroup]); - delete [] names[igroup]; + delete[] names[igroup]; names[igroup] = nullptr; dynamic[igroup] = 0; ngroup--; @@ -135,10 +132,10 @@ void Group::assign(int narg, char **arg) // clear the group - if (strcmp(arg[1],"clear") == 0) { + if (strcmp(arg[1], "clear") == 0) { int igroup = find(arg[0]); - if (igroup == -1) error->all (FLERR,"Could not find group clear group ID {}",arg[0]); - if (igroup == 0) error->all (FLERR,"Cannot clear group all"); + if (igroup == -1) error->all(FLERR, "Could not find group clear group ID {}", arg[0]); + if (igroup == 0) error->all(FLERR, "Cannot clear group all"); int *mask = atom->mask; int nlocal = atom->nlocal; @@ -155,7 +152,7 @@ void Group::assign(int narg, char **arg) bool created = false; if (igroup == -1) { - if (ngroup == MAX_GROUP) error->all(FLERR,"Too many groups (max {})", MAX_GROUP); + if (ngroup == MAX_GROUP) error->all(FLERR, "Too many groups (max {})", MAX_GROUP); igroup = find_unused(); names[igroup] = utils::strdup(arg[0]); ngroup++; @@ -171,81 +168,96 @@ void Group::assign(int narg, char **arg) // style = region // add to group if atom is in region - if (strcmp(arg[1],"region") == 0) { + if (strcmp(arg[1], "region") == 0) { - if (narg != 3) error->all(FLERR,"Illegal group region command"); + if (narg != 3) error->all(FLERR, "Illegal group region command"); auto region = domain->get_region_by_id(arg[2]); - if (!region) error->all(FLERR,"Region {} for group region does not exist", arg[2]); + if (!region) error->all(FLERR, "Region {} for group region does not exist", arg[2]); region->init(); region->prematch(); for (i = 0; i < nlocal; i++) - if (region->match(x[i][0],x[i][1],x[i][2])) - mask[i] |= bit; + if (region->match(x[i][0], x[i][1], x[i][2])) mask[i] |= bit; - // create an empty group + // create an empty group - } else if (strcmp(arg[1],"empty") == 0) { + } else if (strcmp(arg[1], "empty") == 0) { - if (narg != 2) error->all(FLERR,"Illegal group empty command"); + if (narg != 2) error->all(FLERR, "Illegal group empty command"); // nothing else to do here - // style = type, molecule, id - // add to group if atom matches type/molecule/id or condition + // style = type, molecule, id + // add to group if atom matches type/molecule/id or condition - } else if (strcmp(arg[1],"type") == 0 || strcmp(arg[1],"molecule") == 0 || - strcmp(arg[1],"id") == 0) { + } else if (strcmp(arg[1], "type") == 0 || strcmp(arg[1], "molecule") == 0 || + strcmp(arg[1], "id") == 0) { - if (narg < 3) utils::missing_cmd_args(FLERR, std::string("group ")+arg[1], error); + if (narg < 3) utils::missing_cmd_args(FLERR, std::string("group ") + arg[1], error); - int category=NONE; - if (strcmp(arg[1],"type") == 0) category = TYPE; - else if (strcmp(arg[1],"molecule") == 0) category = MOLECULE; - else if (strcmp(arg[1],"id") == 0) category = ID; + int category = NONE; + if (strcmp(arg[1], "type") == 0) + category = TYPE; + else if (strcmp(arg[1], "molecule") == 0) + category = MOLECULE; + else if (strcmp(arg[1], "id") == 0) + category = ID; if ((category == MOLECULE) && (!atom->molecule_flag)) - error->all(FLERR,"Group molecule command requires atom attribute molecule"); + error->all(FLERR, "Group molecule command requires atom attribute molecule"); if ((category == ID) && (!atom->tag_enable)) - error->all(FLERR,"Group id command requires atom IDs"); + error->all(FLERR, "Group id command requires atom IDs"); // args = logical condition if (narg > 3 && - (strcmp(arg[2],"<") == 0 || strcmp(arg[2],">") == 0 || - strcmp(arg[2],"<=") == 0 || strcmp(arg[2],">=") == 0 || - strcmp(arg[2],"==") == 0 || strcmp(arg[2],"!=") == 0 || - strcmp(arg[2],"<>") == 0)) { + (strcmp(arg[2], "<") == 0 || strcmp(arg[2], ">") == 0 || strcmp(arg[2], "<=") == 0 || + strcmp(arg[2], ">=") == 0 || strcmp(arg[2], "==") == 0 || strcmp(arg[2], "!=") == 0 || + strcmp(arg[2], "<>") == 0)) { int condition = -1; - if (strcmp(arg[2],"<") == 0) condition = LT; - else if (strcmp(arg[2],"<=") == 0) condition = LE; - else if (strcmp(arg[2],">") == 0) condition = GT; - else if (strcmp(arg[2],">=") == 0) condition = GE; - else if (strcmp(arg[2],"==") == 0) condition = EQ; - else if (strcmp(arg[2],"!=") == 0) condition = NEQ; - else if (strcmp(arg[2],"<>") == 0) condition = BETWEEN; - else error->all(FLERR,"Illegal group command"); + if (strcmp(arg[2], "<") == 0) + condition = LT; + else if (strcmp(arg[2], "<=") == 0) + condition = LE; + else if (strcmp(arg[2], ">") == 0) + condition = GT; + else if (strcmp(arg[2], ">=") == 0) + condition = GE; + else if (strcmp(arg[2], "==") == 0) + condition = EQ; + else if (strcmp(arg[2], "!=") == 0) + condition = NEQ; + else if (strcmp(arg[2], "<>") == 0) + condition = BETWEEN; + else + error->all(FLERR, "Illegal group command"); - tagint bound1,bound2; + tagint bound1, bound2; if (category == TYPE) bound1 = (tagint) utils::expand_type_int(FLERR, arg[3], Atom::ATOM, lmp); - else bound1 = utils::tnumeric(FLERR, arg[3], false, lmp); + else + bound1 = utils::tnumeric(FLERR, arg[3], false, lmp); bound2 = -1; if (condition == BETWEEN) { - if (narg != 5) error->all(FLERR,"Illegal group command"); + if (narg != 5) error->all(FLERR, "Illegal group command"); if (category == TYPE) bound2 = (tagint) utils::expand_type_int(FLERR, arg[4], Atom::ATOM, lmp); - else bound2 = utils::tnumeric(FLERR, arg[4], false, lmp); - } else if (narg != 4) error->all(FLERR,"Illegal group command"); + else + bound2 = utils::tnumeric(FLERR, arg[4], false, lmp); + } else if (narg != 4) + error->all(FLERR, "Illegal group command"); int *attribute = nullptr; tagint *tattribute = nullptr; - if (category == TYPE) attribute = atom->type; - else if (category == MOLECULE) tattribute = atom->molecule; - else if (category == ID) tattribute = atom->tag; + if (category == TYPE) + attribute = atom->type; + else if (category == MOLECULE) + tattribute = atom->molecule; + else if (category == ID) + tattribute = atom->tag; // add to group if meets condition @@ -270,8 +282,7 @@ void Group::assign(int narg, char **arg) if (attribute[i] != bound1) mask[i] |= bit; } else if (condition == BETWEEN) { for (i = 0; i < nlocal; i++) - if (attribute[i] >= bound1 && attribute[i] <= bound2) - mask[i] |= bit; + if (attribute[i] >= bound1 && attribute[i] <= bound2) mask[i] |= bit; } } else { if (condition == LT) { @@ -294,22 +305,24 @@ void Group::assign(int narg, char **arg) if (tattribute[i] != bound1) mask[i] |= bit; } else if (condition == BETWEEN) { for (i = 0; i < nlocal; i++) - if (tattribute[i] >= bound1 && tattribute[i] <= bound2) - mask[i] |= bit; + if (tattribute[i] >= bound1 && tattribute[i] <= bound2) mask[i] |= bit; } } - // args = list of values + // args = list of values } else { int *attribute = nullptr; tagint *tattribute = nullptr; - if (category == TYPE) attribute = atom->type; - else if (category == MOLECULE) tattribute = atom->molecule; - else if (category == ID) tattribute = atom->tag; + if (category == TYPE) + attribute = atom->type; + else if (category == MOLECULE) + tattribute = atom->molecule; + else if (category == ID) + tattribute = atom->tag; char *typestr = nullptr; - tagint start,stop,delta; + tagint start, stop, delta; for (int iarg = 2; iarg < narg; iarg++) { delta = 1; @@ -320,21 +333,21 @@ void Group::assign(int narg, char **arg) } if (typestr == nullptr) { try { - ValueTokenizer values(arg[iarg],":"); + ValueTokenizer values(arg[iarg], ":"); start = values.next_tagint(); - if (utils::strmatch(arg[iarg],"^-?\\d+$")) { + if (utils::strmatch(arg[iarg], "^-?\\d+$")) { stop = start; - } else if (utils::strmatch(arg[iarg],"^-?\\d+:-?\\d+$")) { + } else if (utils::strmatch(arg[iarg], "^-?\\d+:-?\\d+$")) { stop = values.next_tagint(); - } else if (utils::strmatch(arg[iarg],"^-?\\d+:-?\\d+:\\d+$")) { + } else if (utils::strmatch(arg[iarg], "^-?\\d+:-?\\d+:\\d+$")) { stop = values.next_tagint(); delta = values.next_tagint(); - } else throw TokenizerException("Syntax error",""); + } else + throw TokenizerException("Syntax error", ""); } catch (TokenizerException &e) { - error->all(FLERR,"Incorrect range string '{}': {}",arg[iarg],e.what()); + error->all(FLERR, "Incorrect range string '{}': {}", arg[iarg], e.what()); } - if (delta < 1) - error->all(FLERR,"Illegal range increment value"); + if (delta < 1) error->all(FLERR, "Illegal range increment value"); } // add to group if attribute matches value or sequence @@ -342,33 +355,35 @@ void Group::assign(int narg, char **arg) if (attribute) { for (i = 0; i < nlocal; i++) if (attribute[i] >= start && attribute[i] <= stop && - (attribute[i]-start) % delta == 0) mask[i] |= bit; + (attribute[i] - start) % delta == 0) + mask[i] |= bit; } else { for (i = 0; i < nlocal; i++) if (tattribute[i] >= start && tattribute[i] <= stop && - (tattribute[i]-start) % delta == 0) mask[i] |= bit; + (tattribute[i] - start) % delta == 0) + mask[i] |= bit; } } delete[] typestr; } - // style = variable - // add to group if atom-atyle variable is non-zero + // style = variable + // add to group if atom-atyle variable is non-zero - } else if (strcmp(arg[1],"variable") == 0) { + } else if (strcmp(arg[1], "variable") == 0) { int ivar = input->variable->find(arg[2]); - if (ivar < 0) error->all(FLERR,"Variable name {} for group does not exist",arg[2]); + if (ivar < 0) error->all(FLERR, "Variable name {} for group does not exist", arg[2]); if (!input->variable->atomstyle(ivar)) - error->all(FLERR,"Variable {} for group is invalid style",arg[2]); + error->all(FLERR, "Variable {} for group is invalid style", arg[2]); double *aflag; // aflag = evaluation of per-atom variable - memory->create(aflag,nlocal,"group:aflag"); - input->variable->compute_atom(ivar,0,aflag,1,0); + memory->create(aflag, nlocal, "group:aflag"); + input->variable->compute_atom(ivar, 0, aflag, 1, 0); // add to group if per-atom variable evaluated to non-zero @@ -377,34 +392,35 @@ void Group::assign(int narg, char **arg) memory->destroy(aflag); - // style = include + // style = include - } else if (strcmp(arg[1],"include") == 0) { + } else if (strcmp(arg[1], "include") == 0) { - if (narg != 3) error->all(FLERR,"Illegal group include command"); - if (strcmp(arg[2],"molecule") == 0) { + if (narg != 3) error->all(FLERR, "Illegal group include command"); + if (strcmp(arg[2], "molecule") == 0) { if (!atom->molecule_flag) - error->all(FLERR,"Group include molecule command requires atom attribute molecule"); + error->all(FLERR, "Group include molecule command requires atom attribute molecule"); - add_molecules(igroup,bit); + add_molecules(igroup, bit); - } else error->all(FLERR,"Unknown group include keyword {}",arg[2]); + } else + error->all(FLERR, "Unknown group include keyword {}", arg[2]); - // style = subtract + // style = subtract - } else if (strcmp(arg[1],"subtract") == 0) { + } else if (strcmp(arg[1], "subtract") == 0) { if (narg < 4) utils::missing_cmd_args(FLERR, "group subtract", error); - int length = narg-2; + int length = narg - 2; std::vector list(length); int jgroup; for (int iarg = 2; iarg < narg; iarg++) { jgroup = find(arg[iarg]); - if (jgroup == -1) error->all(FLERR,"Group ID {} does not exist", arg[iarg]); - if (dynamic[jgroup]) error->all(FLERR,"Cannot subtract dynamic groups"); - list[iarg-2] = jgroup; + if (jgroup == -1) error->all(FLERR, "Group ID {} does not exist", arg[iarg]); + if (dynamic[jgroup]) error->all(FLERR, "Cannot subtract dynamic groups"); + list[iarg - 2] = jgroup; } // add to group if in 1st group in list @@ -425,21 +441,21 @@ void Group::assign(int narg, char **arg) if (mask[i] & otherbit) mask[i] &= inverse; } - // style = union + // style = union - } else if (strcmp(arg[1],"union") == 0) { + } else if (strcmp(arg[1], "union") == 0) { if (narg < 3) utils::missing_cmd_args(FLERR, "group union", error); - int length = narg-2; + int length = narg - 2; std::vector list(length); int jgroup; for (int iarg = 2; iarg < narg; iarg++) { jgroup = find(arg[iarg]); - if (jgroup == -1) error->all(FLERR,"Group ID {} does not exist",arg[iarg]); - if (dynamic[jgroup]) error->all(FLERR,"Cannot union groups from a dynamic group"); - list[iarg-2] = jgroup; + if (jgroup == -1) error->all(FLERR, "Group ID {} does not exist", arg[iarg]); + if (dynamic[jgroup]) error->all(FLERR, "Cannot union groups from a dynamic group"); + list[iarg - 2] = jgroup; } // add to group if in any other group in list @@ -452,27 +468,26 @@ void Group::assign(int narg, char **arg) if (mask[i] & otherbit) mask[i] |= bit; } - // style = intersect + // style = intersect - } else if (strcmp(arg[1],"intersect") == 0) { + } else if (strcmp(arg[1], "intersect") == 0) { if (narg < 4) utils::missing_cmd_args(FLERR, "group intersect", error); - int length = narg-2; + int length = narg - 2; std::vector list(length); int jgroup; for (int iarg = 2; iarg < narg; iarg++) { jgroup = find(arg[iarg]); - if (jgroup == -1) error->all(FLERR,"Group ID {} does not exist",arg[iarg]); - if (dynamic[jgroup]) - error->all(FLERR,"Cannot intersect groups using a dynamic group"); - list[iarg-2] = jgroup; + if (jgroup == -1) error->all(FLERR, "Group ID {} does not exist", arg[iarg]); + if (dynamic[jgroup]) error->all(FLERR, "Cannot intersect groups using a dynamic group"); + list[iarg - 2] = jgroup; } // add to group if in all groups in list - int otherbit,ok,ilist; + int otherbit, ok, ilist; for (i = 0; i < nlocal; i++) { ok = 1; @@ -483,50 +498,48 @@ void Group::assign(int narg, char **arg) if (ok) mask[i] |= bit; } - // style = dynamic - // create a new FixGroup to dynamically determine atoms in group + // style = dynamic + // create a new FixGroup to dynamically determine atoms in group - } else if (strcmp(arg[1],"dynamic") == 0) { + } else if (strcmp(arg[1], "dynamic") == 0) { - if (narg < 4) error->all(FLERR,"Illegal group command"); - if (strcmp(arg[0],arg[2]) == 0) - error->all(FLERR,"Group dynamic cannot reference itself"); + if (narg < 4) error->all(FLERR, "Illegal group command"); + if (strcmp(arg[0], arg[2]) == 0) error->all(FLERR, "Group dynamic cannot reference itself"); if (find(arg[2]) < 0) - error->all(FLERR,"Group dynamic parent group {} does not exist", arg[2]); - if (igroup == 0) error->all(FLERR,"Group all cannot be made dynamic"); + error->all(FLERR, "Group dynamic parent group {} does not exist", arg[2]); + if (igroup == 0) error->all(FLERR, "Group all cannot be made dynamic"); // if group is already dynamic, delete existing FixGroup - if (dynamic[igroup]) - modify->delete_fix(std::string("GROUP_") + names[igroup]); + if (dynamic[igroup]) modify->delete_fix(std::string("GROUP_") + names[igroup]); dynamic[igroup] = 1; std::string fixcmd = "GROUP_"; - fixcmd += fmt::format("{} {} GROUP",names[igroup],arg[2]); + fixcmd += fmt::format("{} {} GROUP", names[igroup], arg[2]); for (i = 3; i < narg; i++) fixcmd += std::string(" ") + arg[i]; modify->add_fix(fixcmd); - // style = static - // remove dynamic FixGroup if necessary + // style = static + // remove dynamic FixGroup if necessary - } else if (strcmp(arg[1],"static") == 0) { + } else if (strcmp(arg[1], "static") == 0) { - if (narg != 2) error->all(FLERR,"Illegal group static command"); + if (narg != 2) error->all(FLERR, "Illegal group static command"); - if (dynamic[igroup]) - modify->delete_fix(std::string("GROUP_") + names[igroup]); + if (dynamic[igroup]) modify->delete_fix(std::string("GROUP_") + names[igroup]); dynamic[igroup] = 0; - // not a valid group style + // not a valid group style - } else error->all(FLERR,"Unknown group command keyword: {}",arg[1]); + } else + error->all(FLERR, "Unknown group command keyword: {}", arg[1]); - } catch (LAMMPSException & e) { + } catch (LAMMPSException &e) { // undo created group in case of an error if (created) { - delete [] names[igroup]; + delete[] names[igroup]; names[igroup] = nullptr; ngroup--; } @@ -537,17 +550,18 @@ void Group::assign(int narg, char **arg) int n; n = 0; - for (i = 0; i < nlocal; i++) if (mask[i] & bit) n++; + for (i = 0; i < nlocal; i++) + if (mask[i] & bit) n++; double rlocal = n; double all; - MPI_Allreduce(&rlocal,&all,1,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(&rlocal, &all, 1, MPI_DOUBLE, MPI_SUM, world); if (me == 0) { if (dynamic[igroup]) - utils::logmesg(lmp,"dynamic group {} defined\n",names[igroup]); + utils::logmesg(lmp, "dynamic group {} defined\n", names[igroup]); else - utils::logmesg(lmp,"{:.15g} atoms in group {}\n",all,names[igroup]); + utils::logmesg(lmp, "{:.15g} atoms in group {}\n", all, names[igroup]); } } @@ -558,12 +572,10 @@ void Group::assign(int narg, char **arg) void Group::assign(const std::string &groupcmd) { auto args = utils::split_words(groupcmd); - std::vector newarg(args.size()); - int i=0; - for (const auto &arg : args) { - newarg[i++] = (char *)arg.c_str(); - } - assign(args.size(),newarg.data()); + std::vector newarg(args.size()); + int i = 0; + for (const auto &arg : args) { newarg[i++] = (char *) arg.c_str(); } + assign(args.size(), newarg.data()); } /* ---------------------------------------------------------------------- @@ -580,7 +592,7 @@ void Group::create(const std::string &name, int *flag) int igroup = find(name); if (igroup == -1) { - if (ngroup == MAX_GROUP) error->all(FLERR,"Too many groups (max {})", MAX_GROUP); + if (ngroup == MAX_GROUP) error->all(FLERR, "Too many groups (max {})", MAX_GROUP); igroup = find_unused(); names[igroup] = utils::strdup(name); ngroup++; @@ -617,7 +629,7 @@ int Group::find_or_create(const char *name) int igroup = find(name); if (igroup >= 0) return igroup; - if (ngroup == MAX_GROUP) error->all(FLERR,"Too many groups (max {})", MAX_GROUP); + if (ngroup == MAX_GROUP) error->all(FLERR, "Too many groups (max {})", MAX_GROUP); igroup = find_unused(); names[igroup] = utils::strdup(name); ngroup++; @@ -646,7 +658,7 @@ void Group::add_molecules(int /*igroup*/, int bit) { // hash = unique molecule IDs of atoms already in group - hash = new std::map(); + hash = new std::map(); tagint *molecule = atom->molecule; int *mask = atom->mask; @@ -663,14 +675,14 @@ void Group::add_molecules(int /*igroup*/, int bit) int n = hash->size(); tagint *list; - memory->create(list,n,"group:list"); + memory->create(list, n, "group:list"); n = 0; - std::map::iterator pos; + std::map::iterator pos; for (pos = hash->begin(); pos != hash->end(); ++pos) list[n++] = pos->first; molbit = bit; - comm->ring(n,sizeof(tagint),list,1,molring,nullptr,(void *)this); + comm->ring(n, sizeof(tagint), list, 1, molring, nullptr, (void *) this); delete hash; memory->destroy(list); @@ -687,7 +699,7 @@ void Group::molring(int n, char *cbuf, void *ptr) { auto gptr = (Group *) ptr; auto list = (tagint *) cbuf; - std::map *hash = gptr->hash; + std::map *hash = gptr->hash; int nlocal = gptr->atom->nlocal; tagint *molecule = gptr->atom->molecule; int *mask = gptr->atom->mask; @@ -707,7 +719,7 @@ void Group::molring(int n, char *cbuf, void *ptr) void Group::write_restart(FILE *fp) { - fwrite(&ngroup,sizeof(int),1,fp); + fwrite(&ngroup, sizeof(int), 1, fp); // use count to not change restart format with deleted groups // remove this on next major release @@ -715,11 +727,13 @@ void Group::write_restart(FILE *fp) int n; int count = 0; for (int i = 0; i < MAX_GROUP; i++) { - if (names[i]) n = strlen(names[i]) + 1; - else n = 0; - fwrite(&n,sizeof(int),1,fp); + if (names[i]) + n = strlen(names[i]) + 1; + else + n = 0; + fwrite(&n, sizeof(int), 1, fp); if (n) { - fwrite(names[i],sizeof(char),n,fp); + fwrite(names[i], sizeof(char), n, fp); count++; } if (count == ngroup) break; @@ -733,15 +747,15 @@ void Group::write_restart(FILE *fp) void Group::read_restart(FILE *fp) { - int i,n; + int i, n; // delete existing group names // atom masks will be overwritten by reading of restart file - for (i = 0; i < MAX_GROUP; i++) delete [] names[i]; + for (i = 0; i < MAX_GROUP; i++) delete[] names[i]; - if (me == 0) utils::sfread(FLERR,&ngroup,sizeof(int),1,fp,nullptr,error); - MPI_Bcast(&ngroup,1,MPI_INT,0,world); + if (me == 0) utils::sfread(FLERR, &ngroup, sizeof(int), 1, fp, nullptr, error); + MPI_Bcast(&ngroup, 1, MPI_INT, 0, world); // use count to not change restart format with deleted groups // remove this on next major release @@ -752,14 +766,15 @@ void Group::read_restart(FILE *fp) names[i] = nullptr; continue; } - if (me == 0) utils::sfread(FLERR,&n,sizeof(int),1,fp,nullptr,error); - MPI_Bcast(&n,1,MPI_INT,0,world); + if (me == 0) utils::sfread(FLERR, &n, sizeof(int), 1, fp, nullptr, error); + MPI_Bcast(&n, 1, MPI_INT, 0, world); if (n) { names[i] = new char[n]; - if (me == 0) utils::sfread(FLERR,names[i],sizeof(char),n,fp,nullptr,error); - MPI_Bcast(names[i],n,MPI_CHAR,0,world); + if (me == 0) utils::sfread(FLERR, names[i], sizeof(char), n, fp, nullptr, error); + MPI_Bcast(names[i], n, MPI_CHAR, 0, world); count++; - } else names[i] = nullptr; + } else + names[i] = nullptr; } } @@ -775,7 +790,7 @@ bigint Group::count_all() { bigint nme = atom->nlocal; bigint nall; - MPI_Allreduce(&nme,&nall,1,MPI_LMP_BIGINT,MPI_SUM,world); + MPI_Allreduce(&nme, &nall, 1, MPI_LMP_BIGINT, MPI_SUM, world); return nall; } @@ -796,7 +811,7 @@ bigint Group::count(int igroup) bigint nsingle = n; bigint nall; - MPI_Allreduce(&nsingle,&nall,1,MPI_LMP_BIGINT,MPI_SUM,world); + MPI_Allreduce(&nsingle, &nall, 1, MPI_LMP_BIGINT, MPI_SUM, world); return nall; } @@ -815,11 +830,11 @@ bigint Group::count(int igroup, Region *region) int n = 0; for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) n++; + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) n++; bigint nsingle = n; bigint nall; - MPI_Allreduce(&nsingle,&nall,1,MPI_LMP_BIGINT,MPI_SUM,world); + MPI_Allreduce(&nsingle, &nall, 1, MPI_LMP_BIGINT, MPI_SUM, world); return nall; } @@ -849,7 +864,7 @@ double Group::mass(int igroup) } double all; - MPI_Allreduce(&one,&all,1,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(&one, &all, 1, MPI_DOUBLE, MPI_SUM, world); return all; } @@ -874,16 +889,14 @@ double Group::mass(int igroup, Region *region) if (rmass) { for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) - one += rmass[i]; + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) one += rmass[i]; } else { for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) - one += mass[type[i]]; + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) one += mass[type[i]]; } double all; - MPI_Allreduce(&one,&all,1,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(&one, &all, 1, MPI_DOUBLE, MPI_SUM, world); return all; } @@ -904,7 +917,7 @@ double Group::charge(int igroup) if (mask[i] & groupbit) qone += q[i]; double qall; - MPI_Allreduce(&qone,&qall,1,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(&qone, &qall, 1, MPI_DOUBLE, MPI_SUM, world); return qall; } @@ -924,11 +937,10 @@ double Group::charge(int igroup, Region *region) double qone = 0.0; for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) - qone += q[i]; + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) qone += q[i]; double qall; - MPI_Allreduce(&qone,&qall,1,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(&qone, &qall, 1, MPI_DOUBLE, MPI_SUM, world); return qall; } @@ -951,12 +963,12 @@ void Group::bounds(int igroup, double *minmax) for (int i = 0; i < nlocal; i++) { if (mask[i] & groupbit) { - extent[0] = MIN(extent[0],x[i][0]); - extent[1] = MAX(extent[1],x[i][0]); - extent[2] = MIN(extent[2],x[i][1]); - extent[3] = MAX(extent[3],x[i][1]); - extent[4] = MIN(extent[4],x[i][2]); - extent[5] = MAX(extent[5],x[i][2]); + extent[0] = MIN(extent[0], x[i][0]); + extent[1] = MAX(extent[1], x[i][0]); + extent[2] = MIN(extent[2], x[i][1]); + extent[3] = MAX(extent[3], x[i][1]); + extent[4] = MIN(extent[4], x[i][2]); + extent[5] = MAX(extent[5], x[i][2]); } } @@ -968,7 +980,7 @@ void Group::bounds(int igroup, double *minmax) extent[2] = -extent[2]; extent[4] = -extent[4]; - MPI_Allreduce(extent,minmax,6,MPI_DOUBLE,MPI_MAX,world); + MPI_Allreduce(extent, minmax, 6, MPI_DOUBLE, MPI_MAX, world); minmax[0] = -minmax[0]; minmax[2] = -minmax[2]; @@ -994,13 +1006,13 @@ void Group::bounds(int igroup, double *minmax, Region *region) int nlocal = atom->nlocal; for (int i = 0; i < nlocal; i++) { - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) { - extent[0] = MIN(extent[0],x[i][0]); - extent[1] = MAX(extent[1],x[i][0]); - extent[2] = MIN(extent[2],x[i][1]); - extent[3] = MAX(extent[3],x[i][1]); - extent[4] = MIN(extent[4],x[i][2]); - extent[5] = MAX(extent[5],x[i][2]); + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) { + extent[0] = MIN(extent[0], x[i][0]); + extent[1] = MAX(extent[1], x[i][0]); + extent[2] = MIN(extent[2], x[i][1]); + extent[3] = MAX(extent[3], x[i][1]); + extent[4] = MIN(extent[4], x[i][2]); + extent[5] = MAX(extent[5], x[i][2]); } } @@ -1012,7 +1024,7 @@ void Group::bounds(int igroup, double *minmax, Region *region) extent[2] = -extent[2]; extent[4] = -extent[4]; - MPI_Allreduce(extent,minmax,6,MPI_DOUBLE,MPI_MAX,world); + MPI_Allreduce(extent, minmax, 6, MPI_DOUBLE, MPI_MAX, world); minmax[0] = -minmax[0]; minmax[2] = -minmax[2]; @@ -1048,7 +1060,7 @@ void Group::xcm(int igroup, double masstotal, double *cm) for (int i = 0; i < nlocal; i++) if (mask[i] & groupbit) { massone = rmass[i]; - domain->unmap(x[i],image[i],unwrap); + domain->unmap(x[i], image[i], unwrap); cmone[0] += unwrap[0] * massone; cmone[1] += unwrap[1] * massone; cmone[2] += unwrap[2] * massone; @@ -1057,14 +1069,14 @@ void Group::xcm(int igroup, double masstotal, double *cm) for (int i = 0; i < nlocal; i++) if (mask[i] & groupbit) { massone = mass[type[i]]; - domain->unmap(x[i],image[i],unwrap); + domain->unmap(x[i], image[i], unwrap); cmone[0] += unwrap[0] * massone; cmone[1] += unwrap[1] * massone; cmone[2] += unwrap[2] * massone; } } - MPI_Allreduce(cmone,cm,3,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(cmone, cm, 3, MPI_DOUBLE, MPI_SUM, world); if (masstotal > 0.0) { cm[0] /= masstotal; cm[1] /= masstotal; @@ -1100,25 +1112,25 @@ void Group::xcm(int igroup, double masstotal, double *cm, Region *region) if (rmass) { for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) { + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) { massone = rmass[i]; - domain->unmap(x[i],image[i],unwrap); + domain->unmap(x[i], image[i], unwrap); cmone[0] += unwrap[0] * massone; cmone[1] += unwrap[1] * massone; cmone[2] += unwrap[2] * massone; } } else { for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) { + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) { massone = mass[type[i]]; - domain->unmap(x[i],image[i],unwrap); + domain->unmap(x[i], image[i], unwrap); cmone[0] += unwrap[0] * massone; cmone[1] += unwrap[1] * massone; cmone[2] += unwrap[2] * massone; } } - MPI_Allreduce(cmone,cm,3,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(cmone, cm, 3, MPI_DOUBLE, MPI_SUM, world); if (masstotal > 0.0) { cm[0] /= masstotal; cm[1] /= masstotal; @@ -1143,28 +1155,28 @@ void Group::vcm(int igroup, double masstotal, double *cm) double *rmass = atom->rmass; int nlocal = atom->nlocal; - double p[3],massone; + double p[3], massone; p[0] = p[1] = p[2] = 0.0; if (rmass) { for (int i = 0; i < nlocal; i++) if (mask[i] & groupbit) { massone = rmass[i]; - p[0] += v[i][0]*massone; - p[1] += v[i][1]*massone; - p[2] += v[i][2]*massone; + p[0] += v[i][0] * massone; + p[1] += v[i][1] * massone; + p[2] += v[i][2] * massone; } } else { for (int i = 0; i < nlocal; i++) if (mask[i] & groupbit) { massone = mass[type[i]]; - p[0] += v[i][0]*massone; - p[1] += v[i][1]*massone; - p[2] += v[i][2]*massone; + p[0] += v[i][0] * massone; + p[1] += v[i][1] * massone; + p[2] += v[i][2] * massone; } } - MPI_Allreduce(p,cm,3,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(p, cm, 3, MPI_DOUBLE, MPI_SUM, world); if (masstotal > 0.0) { cm[0] /= masstotal; cm[1] /= masstotal; @@ -1191,28 +1203,28 @@ void Group::vcm(int igroup, double masstotal, double *cm, Region *region) double *rmass = atom->rmass; int nlocal = atom->nlocal; - double p[3],massone; + double p[3], massone; p[0] = p[1] = p[2] = 0.0; if (rmass) { for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) { + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) { massone = rmass[i]; - p[0] += v[i][0]*massone; - p[1] += v[i][1]*massone; - p[2] += v[i][2]*massone; + p[0] += v[i][0] * massone; + p[1] += v[i][1] * massone; + p[2] += v[i][2] * massone; } } else { for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) { + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) { massone = mass[type[i]]; - p[0] += v[i][0]*massone; - p[1] += v[i][1]*massone; - p[2] += v[i][2]*massone; + p[0] += v[i][0] * massone; + p[1] += v[i][1] * massone; + p[2] += v[i][2] * massone; } } - MPI_Allreduce(p,cm,3,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(p, cm, 3, MPI_DOUBLE, MPI_SUM, world); if (masstotal > 0.0) { cm[0] /= masstotal; cm[1] /= masstotal; @@ -1242,7 +1254,7 @@ void Group::fcm(int igroup, double *cm) flocal[2] += f[i][2]; } - MPI_Allreduce(flocal,cm,3,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(flocal, cm, 3, MPI_DOUBLE, MPI_SUM, world); } /* ---------------------------------------------------------------------- @@ -1263,13 +1275,13 @@ void Group::fcm(int igroup, double *cm, Region *region) flocal[0] = flocal[1] = flocal[2] = 0.0; for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) { + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) { flocal[0] += f[i][0]; flocal[1] += f[i][1]; flocal[2] += f[i][2]; } - MPI_Allreduce(flocal,cm,3,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(flocal, cm, 3, MPI_DOUBLE, MPI_SUM, world); } /* ---------------------------------------------------------------------- @@ -1292,17 +1304,15 @@ double Group::ke(int igroup) if (rmass) { for (int i = 0; i < nlocal; i++) if (mask[i] & groupbit) - one += (v[i][0]*v[i][0] + v[i][1]*v[i][1] + v[i][2]*v[i][2]) * - rmass[i]; + one += (v[i][0] * v[i][0] + v[i][1] * v[i][1] + v[i][2] * v[i][2]) * rmass[i]; } else { for (int i = 0; i < nlocal; i++) if (mask[i] & groupbit) - one += (v[i][0]*v[i][0] + v[i][1]*v[i][1] + v[i][2]*v[i][2]) * - mass[type[i]]; + one += (v[i][0] * v[i][0] + v[i][1] * v[i][1] + v[i][2] * v[i][2]) * mass[type[i]]; } double all; - MPI_Allreduce(&one,&all,1,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(&one, &all, 1, MPI_DOUBLE, MPI_SUM, world); all *= 0.5 * force->mvv2e; return all; } @@ -1328,18 +1338,16 @@ double Group::ke(int igroup, Region *region) if (rmass) { for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) - one += (v[i][0]*v[i][0] + v[i][1]*v[i][1] + v[i][2]*v[i][2]) * - rmass[i]; + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) + one += (v[i][0] * v[i][0] + v[i][1] * v[i][1] + v[i][2] * v[i][2]) * rmass[i]; } else { for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) - one += (v[i][0]*v[i][0] + v[i][1]*v[i][1] + v[i][2]*v[i][2]) * - mass[type[i]]; + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) + one += (v[i][0] * v[i][0] + v[i][1] * v[i][1] + v[i][2] * v[i][2]) * mass[type[i]]; } double all; - MPI_Allreduce(&one,&all,1,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(&one, &all, 1, MPI_DOUBLE, MPI_SUM, world); all *= 0.5 * force->mvv2e; return all; } @@ -1362,24 +1370,26 @@ double Group::gyration(int igroup, double masstotal, double *cm) double *rmass = atom->rmass; int nlocal = atom->nlocal; - double dx,dy,dz,massone; + double dx, dy, dz, massone; double unwrap[3]; double rg = 0.0; for (int i = 0; i < nlocal; i++) if (mask[i] & groupbit) { - domain->unmap(x[i],image[i],unwrap); + domain->unmap(x[i], image[i], unwrap); dx = unwrap[0] - cm[0]; dy = unwrap[1] - cm[1]; dz = unwrap[2] - cm[2]; - if (rmass) massone = rmass[i]; - else massone = mass[type[i]]; - rg += (dx*dx + dy*dy + dz*dz) * massone; + if (rmass) + massone = rmass[i]; + else + massone = mass[type[i]]; + rg += (dx * dx + dy * dy + dz * dz) * massone; } double rg_all; - MPI_Allreduce(&rg,&rg_all,1,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(&rg, &rg_all, 1, MPI_DOUBLE, MPI_SUM, world); - if (masstotal > 0.0) return sqrt(rg_all/masstotal); + if (masstotal > 0.0) return sqrt(rg_all / masstotal); return 0.0; } @@ -1402,24 +1412,26 @@ double Group::gyration(int igroup, double masstotal, double *cm, Region *region) double *rmass = atom->rmass; int nlocal = atom->nlocal; - double dx,dy,dz,massone; + double dx, dy, dz, massone; double unwrap[3]; double rg = 0.0; for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) { - domain->unmap(x[i],image[i],unwrap); + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) { + domain->unmap(x[i], image[i], unwrap); dx = unwrap[0] - cm[0]; dy = unwrap[1] - cm[1]; dz = unwrap[2] - cm[2]; - if (rmass) massone = rmass[i]; - else massone = mass[type[i]]; - rg += (dx*dx + dy*dy + dz*dz) * massone; + if (rmass) + massone = rmass[i]; + else + massone = mass[type[i]]; + rg += (dx * dx + dy * dy + dz * dz) * massone; } double rg_all; - MPI_Allreduce(&rg,&rg_all,1,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(&rg, &rg_all, 1, MPI_DOUBLE, MPI_SUM, world); - if (masstotal > 0.0) return sqrt(rg_all/masstotal); + if (masstotal > 0.0) return sqrt(rg_all / masstotal); return 0.0; } @@ -1442,7 +1454,7 @@ void Group::angmom(int igroup, double *cm, double *lmom) double *rmass = atom->rmass; int nlocal = atom->nlocal; - double dx,dy,dz,massone; + double dx, dy, dz, massone; double unwrap[3]; double p[3]; @@ -1450,18 +1462,20 @@ void Group::angmom(int igroup, double *cm, double *lmom) for (int i = 0; i < nlocal; i++) if (mask[i] & groupbit) { - domain->unmap(x[i],image[i],unwrap); + domain->unmap(x[i], image[i], unwrap); dx = unwrap[0] - cm[0]; dy = unwrap[1] - cm[1]; dz = unwrap[2] - cm[2]; - if (rmass) massone = rmass[i]; - else massone = mass[type[i]]; - p[0] += massone * (dy*v[i][2] - dz*v[i][1]); - p[1] += massone * (dz*v[i][0] - dx*v[i][2]); - p[2] += massone * (dx*v[i][1] - dy*v[i][0]); + if (rmass) + massone = rmass[i]; + else + massone = mass[type[i]]; + p[0] += massone * (dy * v[i][2] - dz * v[i][1]); + p[1] += massone * (dz * v[i][0] - dx * v[i][2]); + p[2] += massone * (dx * v[i][1] - dy * v[i][0]); } - MPI_Allreduce(p,lmom,3,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(p, lmom, 3, MPI_DOUBLE, MPI_SUM, world); } /* ---------------------------------------------------------------------- @@ -1484,26 +1498,28 @@ void Group::angmom(int igroup, double *cm, double *lmom, Region *region) double *rmass = atom->rmass; int nlocal = atom->nlocal; - double dx,dy,dz,massone; + double dx, dy, dz, massone; double unwrap[3]; double p[3]; p[0] = p[1] = p[2] = 0.0; for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) { - domain->unmap(x[i],image[i],unwrap); + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) { + domain->unmap(x[i], image[i], unwrap); dx = unwrap[0] - cm[0]; dy = unwrap[1] - cm[1]; dz = unwrap[2] - cm[2]; - if (rmass) massone = rmass[i]; - else massone = mass[type[i]]; - p[0] += massone * (dy*v[i][2] - dz*v[i][1]); - p[1] += massone * (dz*v[i][0] - dx*v[i][2]); - p[2] += massone * (dx*v[i][1] - dy*v[i][0]); + if (rmass) + massone = rmass[i]; + else + massone = mass[type[i]]; + p[0] += massone * (dy * v[i][2] - dz * v[i][1]); + p[1] += massone * (dz * v[i][0] - dx * v[i][2]); + p[2] += massone * (dx * v[i][1] - dy * v[i][0]); } - MPI_Allreduce(p,lmom,3,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(p, lmom, 3, MPI_DOUBLE, MPI_SUM, world); } /* ---------------------------------------------------------------------- @@ -1522,7 +1538,7 @@ void Group::torque(int igroup, double *cm, double *tq) imageint *image = atom->image; int nlocal = atom->nlocal; - double dx,dy,dz; + double dx, dy, dz; double unwrap[3]; double tlocal[3]; @@ -1530,16 +1546,16 @@ void Group::torque(int igroup, double *cm, double *tq) for (int i = 0; i < nlocal; i++) if (mask[i] & groupbit) { - domain->unmap(x[i],image[i],unwrap); + domain->unmap(x[i], image[i], unwrap); dx = unwrap[0] - cm[0]; dy = unwrap[1] - cm[1]; dz = unwrap[2] - cm[2]; - tlocal[0] += dy*f[i][2] - dz*f[i][1]; - tlocal[1] += dz*f[i][0] - dx*f[i][2]; - tlocal[2] += dx*f[i][1] - dy*f[i][0]; + tlocal[0] += dy * f[i][2] - dz * f[i][1]; + tlocal[1] += dz * f[i][0] - dx * f[i][2]; + tlocal[2] += dx * f[i][1] - dy * f[i][0]; } - MPI_Allreduce(tlocal,tq,3,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(tlocal, tq, 3, MPI_DOUBLE, MPI_SUM, world); } /* ---------------------------------------------------------------------- @@ -1559,24 +1575,24 @@ void Group::torque(int igroup, double *cm, double *tq, Region *region) imageint *image = atom->image; int nlocal = atom->nlocal; - double dx,dy,dz; + double dx, dy, dz; double unwrap[3]; double tlocal[3]; tlocal[0] = tlocal[1] = tlocal[2] = 0.0; for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) { - domain->unmap(x[i],image[i],unwrap); + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) { + domain->unmap(x[i], image[i], unwrap); dx = unwrap[0] - cm[0]; dy = unwrap[1] - cm[1]; dz = unwrap[2] - cm[2]; - tlocal[0] += dy*f[i][2] - dz*f[i][1]; - tlocal[1] += dz*f[i][0] - dx*f[i][2]; - tlocal[2] += dx*f[i][1] - dy*f[i][0]; + tlocal[0] += dy * f[i][2] - dz * f[i][1]; + tlocal[1] += dz * f[i][0] - dx * f[i][2]; + tlocal[2] += dx * f[i][1] - dy * f[i][0]; } - MPI_Allreduce(tlocal,tq,3,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(tlocal, tq, 3, MPI_DOUBLE, MPI_SUM, world); } /* ---------------------------------------------------------------------- @@ -1586,7 +1602,7 @@ void Group::torque(int igroup, double *cm, double *tq, Region *region) void Group::inertia(int igroup, double *cm, double itensor[3][3]) { - int i,j; + int i, j; int groupbit = bitmask[igroup]; @@ -1598,34 +1614,35 @@ void Group::inertia(int igroup, double *cm, double itensor[3][3]) double *rmass = atom->rmass; int nlocal = atom->nlocal; - double dx,dy,dz,massone; + double dx, dy, dz, massone; double unwrap[3]; double ione[3][3]; for (i = 0; i < 3; i++) - for (j = 0; j < 3; j++) - ione[i][j] = 0.0; + for (j = 0; j < 3; j++) ione[i][j] = 0.0; for (i = 0; i < nlocal; i++) if (mask[i] & groupbit) { - domain->unmap(x[i],image[i],unwrap); + domain->unmap(x[i], image[i], unwrap); dx = unwrap[0] - cm[0]; dy = unwrap[1] - cm[1]; dz = unwrap[2] - cm[2]; - if (rmass) massone = rmass[i]; - else massone = mass[type[i]]; - ione[0][0] += massone * (dy*dy + dz*dz); - ione[1][1] += massone * (dx*dx + dz*dz); - ione[2][2] += massone * (dx*dx + dy*dy); - ione[0][1] -= massone * dx*dy; - ione[1][2] -= massone * dy*dz; - ione[0][2] -= massone * dx*dz; + if (rmass) + massone = rmass[i]; + else + massone = mass[type[i]]; + ione[0][0] += massone * (dy * dy + dz * dz); + ione[1][1] += massone * (dx * dx + dz * dz); + ione[2][2] += massone * (dx * dx + dy * dy); + ione[0][1] -= massone * dx * dy; + ione[1][2] -= massone * dy * dz; + ione[0][2] -= massone * dx * dz; } ione[1][0] = ione[0][1]; ione[2][1] = ione[1][2]; ione[2][0] = ione[0][2]; - MPI_Allreduce(&ione[0][0],&itensor[0][0],9,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(&ione[0][0], &itensor[0][0], 9, MPI_DOUBLE, MPI_SUM, world); } /* ---------------------------------------------------------------------- @@ -1635,7 +1652,7 @@ void Group::inertia(int igroup, double *cm, double itensor[3][3]) void Group::inertia(int igroup, double *cm, double itensor[3][3], Region *region) { - int i,j; + int i, j; int groupbit = bitmask[igroup]; region->prematch(); @@ -1648,34 +1665,35 @@ void Group::inertia(int igroup, double *cm, double itensor[3][3], Region *region double *rmass = atom->rmass; int nlocal = atom->nlocal; - double dx,dy,dz,massone; + double dx, dy, dz, massone; double unwrap[3]; double ione[3][3]; for (i = 0; i < 3; i++) - for (j = 0; j < 3; j++) - ione[i][j] = 0.0; + for (j = 0; j < 3; j++) ione[i][j] = 0.0; for (i = 0; i < nlocal; i++) - if (mask[i] & groupbit && region->match(x[i][0],x[i][1],x[i][2])) { - domain->unmap(x[i],image[i],unwrap); + if (mask[i] & groupbit && region->match(x[i][0], x[i][1], x[i][2])) { + domain->unmap(x[i], image[i], unwrap); dx = unwrap[0] - cm[0]; dy = unwrap[1] - cm[1]; dz = unwrap[2] - cm[2]; - if (rmass) massone = rmass[i]; - else massone = mass[type[i]]; - ione[0][0] += massone * (dy*dy + dz*dz); - ione[1][1] += massone * (dx*dx + dz*dz); - ione[2][2] += massone * (dx*dx + dy*dy); - ione[0][1] -= massone * dx*dy; - ione[1][2] -= massone * dy*dz; - ione[0][2] -= massone * dx*dz; + if (rmass) + massone = rmass[i]; + else + massone = mass[type[i]]; + ione[0][0] += massone * (dy * dy + dz * dz); + ione[1][1] += massone * (dx * dx + dz * dz); + ione[2][2] += massone * (dx * dx + dy * dy); + ione[0][1] -= massone * dx * dy; + ione[1][2] -= massone * dy * dz; + ione[0][2] -= massone * dx * dz; } ione[1][0] = ione[0][1]; ione[2][1] = ione[1][2]; ione[2][0] = ione[0][2]; - MPI_Allreduce(&ione[0][0],&itensor[0][0],9,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(&ione[0][0], &itensor[0][0], 9, MPI_DOUBLE, MPI_SUM, world); } /* ---------------------------------------------------------------------- @@ -1684,17 +1702,16 @@ void Group::inertia(int igroup, double *cm, double itensor[3][3], Region *region void Group::omega(double *angmom, double inertia[3][3], double *w) { - double idiag[3],ex[3],ey[3],ez[3],cross[3]; - double evectors[3][3],inverse[3][3]; + double idiag[3], ex[3], ey[3], ez[3], cross[3]; + double evectors[3][3], inverse[3][3]; // determinant = triple product of rows of inertia matrix - double determinant = inertia[0][0]*inertia[1][1]*inertia[2][2] + - inertia[0][1]*inertia[1][2]*inertia[2][0] + - inertia[0][2]*inertia[1][0]*inertia[2][1] - - inertia[0][0]*inertia[1][2]*inertia[2][1] - - inertia[0][1]*inertia[1][0]*inertia[2][2] - - inertia[2][0]*inertia[1][1]*inertia[0][2]; + double determinant = inertia[0][0] * inertia[1][1] * inertia[2][2] + + inertia[0][1] * inertia[1][2] * inertia[2][0] + + inertia[0][2] * inertia[1][0] * inertia[2][1] - + inertia[0][0] * inertia[1][2] * inertia[2][1] - + inertia[0][1] * inertia[1][0] * inertia[2][2] - inertia[2][0] * inertia[1][1] * inertia[0][2]; // non-singular I matrix // use L = Iw, inverting I to solve for w @@ -1702,37 +1719,29 @@ void Group::omega(double *angmom, double inertia[3][3], double *w) if (determinant > EPSILON) { - inverse[0][0] = inertia[1][1]*inertia[2][2] - inertia[1][2]*inertia[2][1]; - inverse[0][1] = -(inertia[0][1]*inertia[2][2] - - inertia[0][2]*inertia[2][1]); - inverse[0][2] = inertia[0][1]*inertia[1][2] - inertia[0][2]*inertia[1][1]; + inverse[0][0] = inertia[1][1] * inertia[2][2] - inertia[1][2] * inertia[2][1]; + inverse[0][1] = -(inertia[0][1] * inertia[2][2] - inertia[0][2] * inertia[2][1]); + inverse[0][2] = inertia[0][1] * inertia[1][2] - inertia[0][2] * inertia[1][1]; - inverse[1][0] = -(inertia[1][0]*inertia[2][2] - - inertia[1][2]*inertia[2][0]); - inverse[1][1] = inertia[0][0]*inertia[2][2] - inertia[0][2]*inertia[2][0]; - inverse[1][2] = -(inertia[0][0]*inertia[1][2] - - inertia[0][2]*inertia[1][0]); + inverse[1][0] = -(inertia[1][0] * inertia[2][2] - inertia[1][2] * inertia[2][0]); + inverse[1][1] = inertia[0][0] * inertia[2][2] - inertia[0][2] * inertia[2][0]; + inverse[1][2] = -(inertia[0][0] * inertia[1][2] - inertia[0][2] * inertia[1][0]); - inverse[2][0] = inertia[1][0]*inertia[2][1] - inertia[1][1]*inertia[2][0]; - inverse[2][1] = -(inertia[0][0]*inertia[2][1] - - inertia[0][1]*inertia[2][0]); - inverse[2][2] = inertia[0][0]*inertia[1][1] - inertia[0][1]*inertia[1][0]; + inverse[2][0] = inertia[1][0] * inertia[2][1] - inertia[1][1] * inertia[2][0]; + inverse[2][1] = -(inertia[0][0] * inertia[2][1] - inertia[0][1] * inertia[2][0]); + inverse[2][2] = inertia[0][0] * inertia[1][1] - inertia[0][1] * inertia[1][0]; for (int i = 0; i < 3; i++) - for (int j = 0; j < 3; j++) - inverse[i][j] /= determinant; + for (int j = 0; j < 3; j++) inverse[i][j] /= determinant; - w[0] = inverse[0][0]*angmom[0] + inverse[0][1]*angmom[1] + - inverse[0][2]*angmom[2]; - w[1] = inverse[1][0]*angmom[0] + inverse[1][1]*angmom[1] + - inverse[1][2]*angmom[2]; - w[2] = inverse[2][0]*angmom[0] + inverse[2][1]*angmom[1] + - inverse[2][2]*angmom[2]; + w[0] = inverse[0][0] * angmom[0] + inverse[0][1] * angmom[1] + inverse[0][2] * angmom[2]; + w[1] = inverse[1][0] * angmom[0] + inverse[1][1] * angmom[1] + inverse[1][2] * angmom[2]; + w[2] = inverse[2][0] * angmom[0] + inverse[2][1] * angmom[1] + inverse[2][2] * angmom[2]; - // handle (nearly) singular I matrix - // typically due to 2-atom group or linear molecule - // use jacobi3() and angmom_to_omega() to calculate valid omega - // less exact answer than matrix inversion, due to iterative Jacobi method + // handle (nearly) singular I matrix + // typically due to 2-atom group or linear molecule + // use jacobi3() and angmom_to_omega() to calculate valid omega + // less exact answer than matrix inversion, due to iterative Jacobi method } else { int ierror = MathEigen::jacobi3(inertia, idiag, evectors); @@ -1751,21 +1760,21 @@ void Group::omega(double *angmom, double inertia[3][3], double *w) // enforce 3 evectors as a right-handed coordinate system // flip 3rd vector if needed - MathExtra::cross3(ex,ey,cross); - if (MathExtra::dot3(cross,ez) < 0.0) MathExtra::negate3(ez); + MathExtra::cross3(ex, ey, cross); + if (MathExtra::dot3(cross, ez) < 0.0) MathExtra::negate3(ez); // if any principal moment < scaled EPSILON, set to 0.0 double max; - max = MAX(idiag[0],idiag[1]); - max = MAX(max,idiag[2]); + max = MAX(idiag[0], idiag[1]); + max = MAX(max, idiag[2]); - if (idiag[0] < EPSILON*max) idiag[0] = 0.0; - if (idiag[1] < EPSILON*max) idiag[1] = 0.0; - if (idiag[2] < EPSILON*max) idiag[2] = 0.0; + if (idiag[0] < EPSILON * max) idiag[0] = 0.0; + if (idiag[1] < EPSILON * max) idiag[1] = 0.0; + if (idiag[2] < EPSILON * max) idiag[2] = 0.0; // calculate omega using diagonalized inertia matrix - MathExtra::angmom_to_omega(angmom,ex,ey,ez,idiag,w); + MathExtra::angmom_to_omega(angmom, ex, ey, ez, idiag, w); } }