diff --git a/src/MC/fix_neighbor_swap.cpp b/src/MC/fix_neighbor_swap.cpp index 703207783a..95e74bee22 100755 --- a/src/MC/fix_neighbor_swap.cpp +++ b/src/MC/fix_neighbor_swap.cpp @@ -1,4 +1,3 @@ -// clang-format off /* ---------------------------------------------------------------------- LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator https://www.lammps.org/, Sandia National Laboratories @@ -23,6 +22,7 @@ #include "bond.h" #include "comm.h" #include "compute.h" +#include "compute_voronoi_atom.h" #include "dihedral.h" #include "domain.h" #include "error.h" @@ -39,37 +39,36 @@ #include "random_park.h" #include "region.h" #include "update.h" -#include "compute_voronoi_atom.h" -#include #include #include +#include #include using namespace LAMMPS_NS; using namespace FixConst; static const char cite_fix_neighbor_swap_c[] = - "fix neighbor/swap command: doi:10.1016/j.commatsci.2022.111929\n\n" - "@Article{Tavenner2023111929,\n" - " author = {Jacob P. Tavenner and Mikhail I. Mendelev and John W. Lawson},\n" - " title = {Molecular dynamics based kinetic Monte Carlo simulation for accelerated diffusion},\n" - " journal = {Computational Materials Science},\n" - " year = {2023},\n" - " volume = {218},\n" - " pages = {111929}\n" - " url = {https://www.sciencedirect.com/science/article/pii/S0927025622006401}\n" - "}\n\n"; + "fix neighbor/swap command: doi:10.1016/j.commatsci.2022.111929\n\n" + "@Article{Tavenner2023111929,\n" + " author = {Jacob P. Tavenner and Mikhail I. Mendelev and John W. Lawson},\n" + " title = {Molecular dynamics based kinetic Monte Carlo simulation for accelerated " + "diffusion},\n" + " journal = {Computational Materials Science},\n" + " year = {2023},\n" + " volume = {218},\n" + " pages = {111929}\n" + " url = {https://www.sciencedirect.com/science/article/pii/S0927025622006401}\n" + "}\n\n"; /* ---------------------------------------------------------------------- */ FixNeighborSwap::FixNeighborSwap(LAMMPS *lmp, int narg, char **arg) : - Fix(lmp, narg, arg), region(nullptr), idregion(nullptr), type_list(nullptr), - qtype(nullptr), c_voro(nullptr), voro_neighbor_list(nullptr), - sqrt_mass_ratio(nullptr), local_swap_iatom_list(nullptr), - random_equal(nullptr), random_unequal(nullptr), c_pe(nullptr) + Fix(lmp, narg, arg), region(nullptr), idregion(nullptr), type_list(nullptr), qtype(nullptr), + c_voro(nullptr), voro_neighbor_list(nullptr), sqrt_mass_ratio(nullptr), + local_swap_iatom_list(nullptr), random_equal(nullptr), random_unequal(nullptr), c_pe(nullptr) { - if (narg < 10) utils::missing_cmd_args(FLERR,"fix neighbor/swap", error); + if (narg < 10) utils::missing_cmd_args(FLERR, "fix neighbor/swap", error); dynamic_group_allow = 1; @@ -82,34 +81,34 @@ FixNeighborSwap::FixNeighborSwap(LAMMPS *lmp, int narg, char **arg) : // required args - nevery = utils::inumeric(FLERR,arg[3],false,lmp); - ncycles = utils::inumeric(FLERR,arg[4],false,lmp); - seed = utils::inumeric(FLERR,arg[5],false,lmp); - double temperature = utils::numeric(FLERR,arg[6],false,lmp); + nevery = utils::inumeric(FLERR, arg[3], false, lmp); + ncycles = utils::inumeric(FLERR, arg[4], false, lmp); + seed = utils::inumeric(FLERR, arg[5], false, lmp); + double temperature = utils::numeric(FLERR, arg[6], false, lmp); - if (nevery <= 0) error->all(FLERR,"Illegal fix neighbor/swap command nevery value"); - if (ncycles < 0) error->all(FLERR,"Illegal fix neighbor/swap command ncycles value"); - if (seed <= 0) error->all(FLERR,"Illegal fix neighbor/swap command seed value"); - if (temperature <= 0.0) error->all(FLERR,"Illegal fix neighbor/swap command temperature value"); + if (nevery <= 0) error->all(FLERR, "Illegal fix neighbor/swap command nevery value"); + if (ncycles < 0) error->all(FLERR, "Illegal fix neighbor/swap command ncycles value"); + if (seed <= 0) error->all(FLERR, "Illegal fix neighbor/swap command seed value"); + if (temperature <= 0.0) error->all(FLERR, "Illegal fix neighbor/swap command temperature value"); beta = 1.0 / (force->boltz * temperature); - memory->create(type_list,atom->ntypes,"neighbor/swap:type_list"); - memory->create(rate_list,atom->ntypes,"neighbor/swap:rate_list"); + memory->create(type_list, atom->ntypes, "neighbor/swap:type_list"); + memory->create(rate_list, atom->ntypes, "neighbor/swap:rate_list"); // read options from end of input line - options(narg-7,&arg[7]); + options(narg - 7, &arg[7]); - if (voro_flag != 1) error->all(FLERR,"Voronoi compute required for fix neighbor/swap command"); + if (voro_flag != 1) error->all(FLERR, "Voronoi compute required for fix neighbor/swap command"); // random number generator, same for all procs - random_equal = new RanPark(lmp,seed); + random_equal = new RanPark(lmp, seed); // random number generator, not the same for all procs - random_unequal = new RanPark(lmp,seed); + random_unequal = new RanPark(lmp, seed); // set up reneighboring @@ -130,9 +129,10 @@ FixNeighborSwap::FixNeighborSwap(LAMMPS *lmp, int narg, char **arg) : // set comm size needed by this Fix - if (atom->q_flag) comm_forward = 2; - else comm_forward = 1; - + if (atom->q_flag) + comm_forward = 2; + else + comm_forward = 1; } /* ---------------------------------------------------------------------- */ @@ -158,7 +158,7 @@ FixNeighborSwap::~FixNeighborSwap() void FixNeighborSwap::options(int narg, char **arg) { - if (narg < 0) error->all(FLERR,"Illegal fix neighbor/swap command\n"); + if (narg < 0) error->all(FLERR, "Illegal fix neighbor/swap command\n"); ke_flag = 1; diff_flag = 0; @@ -168,64 +168,64 @@ void FixNeighborSwap::options(int narg, char **arg) int iarg = 0; while (iarg < narg) { - if (strcmp(arg[iarg],"region") == 0) { - if (iarg+2 > narg) error->all(FLERR,"Illegal fix neighbor/swap command"); - region = domain->get_region_by_id(arg[iarg+1]); - if (!region) error->all(FLERR,"Region ID for fix neighbor/swap does not exist"); + if (strcmp(arg[iarg], "region") == 0) { + if (iarg + 2 > narg) error->all(FLERR, "Illegal fix neighbor/swap command"); + region = domain->get_region_by_id(arg[iarg + 1]); + if (!region) error->all(FLERR, "Region ID for fix neighbor/swap does not exist"); idregion = utils::strdup(arg[iarg + 1]); iarg += 2; - } else if (strcmp(arg[iarg],"ke") == 0) { - if (iarg + 2 > narg) error->all(FLERR,"Illegal fix neighbor/swap command"); - ke_flag = utils::logical(FLERR,arg[iarg+1],false,lmp); + } else if (strcmp(arg[iarg], "ke") == 0) { + if (iarg + 2 > narg) error->all(FLERR, "Illegal fix neighbor/swap command"); + ke_flag = utils::logical(FLERR, arg[iarg + 1], false, lmp); iarg += 2; - } else if (strcmp(arg[iarg],"types") == 0) { - if (iarg + 3 > narg) error->all(FLERR,"Illegal fix neighbor/swap command"); - if (diff_flag != 0) error->all(FLERR,"Illegal fix neighbor/swap command"); + } else if (strcmp(arg[iarg], "types") == 0) { + if (iarg + 3 > narg) error->all(FLERR, "Illegal fix neighbor/swap command"); + if (diff_flag != 0) error->all(FLERR, "Illegal fix neighbor/swap command"); iarg++; nswaptypes = 0; while (iarg < narg) { if (isalpha(arg[iarg][0])) break; - if (nswaptypes >= atom->ntypes) error->all(FLERR,"Illegal fix neighbor/swap command"); - type_list[nswaptypes] = utils::numeric(FLERR,arg[iarg],false,lmp); + if (nswaptypes >= atom->ntypes) error->all(FLERR, "Illegal fix neighbor/swap command"); + type_list[nswaptypes] = utils::numeric(FLERR, arg[iarg], false, lmp); nswaptypes++; iarg++; } - } else if (strcmp(arg[iarg],"voro") == 0) { - if (iarg + 2 > narg) error->all(FLERR,"Illegal fix neighbor/swap command"); + } else if (strcmp(arg[iarg], "voro") == 0) { + if (iarg + 2 > narg) error->all(FLERR, "Illegal fix neighbor/swap command"); - int icompute = modify->find_compute(utils::strdup(arg[iarg+1])); + int icompute = modify->find_compute(utils::strdup(arg[iarg + 1])); - if (icompute < 0) error->all(FLERR,"Could not find neighbor compute ID"); + if (icompute < 0) error->all(FLERR, "Could not find neighbor compute ID"); c_voro = modify->compute[icompute]; if (c_voro->local_flag == 0) - error->all(FLERR,"Neighbor compute does not compute local info"); + error->all(FLERR, "Neighbor compute does not compute local info"); if (c_voro->size_local_cols != 3) - error->all(FLERR,"Neighbor compute does not give i, j, size as expected"); + error->all(FLERR, "Neighbor compute does not give i, j, size as expected"); voro_flag = 1; iarg += 2; - } else if (strcmp(arg[iarg],"diff") == 0) { - if (iarg + 2 > narg) error->all(FLERR,"Illegal fix neighbor/swap command"); - if (nswaptypes != 0) error->all(FLERR,"Illegal fix neighbor/swap command"); - type_list[nswaptypes] = utils::numeric(FLERR,arg[iarg+1],false,lmp); + } else if (strcmp(arg[iarg], "diff") == 0) { + if (iarg + 2 > narg) error->all(FLERR, "Illegal fix neighbor/swap command"); + if (nswaptypes != 0) error->all(FLERR, "Illegal fix neighbor/swap command"); + type_list[nswaptypes] = utils::numeric(FLERR, arg[iarg + 1], false, lmp); diff_flag = 1; nswaptypes++; iarg += 2; - } else if (strcmp(arg[iarg],"rates") == 0) { - if (iarg+atom->ntypes >= narg) error->all(FLERR,"Illegal fix neighbor/swap command"); + } else if (strcmp(arg[iarg], "rates") == 0) { + if (iarg + atom->ntypes >= narg) error->all(FLERR, "Illegal fix neighbor/swap command"); iarg++; int i = 0; while (iarg < narg) { if (isalpha(arg[iarg][0])) break; - if (i >= atom->ntypes) error->all(FLERR,"Illegal fix neighbor/swap command"); - rate_list[i] = utils::numeric(FLERR,arg[iarg],false,lmp); + if (i >= atom->ntypes) error->all(FLERR, "Illegal fix neighbor/swap command"); + rate_list[i] = utils::numeric(FLERR, arg[iarg], false, lmp); i++; iarg++; } rates_flag = 1; - if (i != atom->ntypes) error->all(FLERR,"Illegal fix neighbor/swap command"); - } - else error->all(FLERR,"Illegal fix neighbor/swap command"); + if (i != atom->ntypes) error->all(FLERR, "Illegal fix neighbor/swap command"); + } else + error->all(FLERR, "Illegal fix neighbor/swap command"); } } @@ -246,7 +246,8 @@ void FixNeighborSwap::init() int *type = atom->type; - if (nswaptypes < 2 && !diff_flag) error->all(FLERR,"Must specify at least 2 types in fix neighbor/swap command"); + if (nswaptypes < 2 && !diff_flag) + error->all(FLERR, "Must specify at least 2 types in fix neighbor/swap command"); // set index and check validity of region @@ -257,12 +258,12 @@ void FixNeighborSwap::init() for (int iswaptype = 0; iswaptype < nswaptypes; iswaptype++) if (type_list[iswaptype] <= 0 || type_list[iswaptype] > atom->ntypes) - error->all(FLERR,"Invalid atom type in fix neighbor/swap command"); + error->all(FLERR, "Invalid atom type in fix neighbor/swap command"); if (atom->q_flag) { - double qmax,qmin; - int firstall,first; - memory->create(qtype,nswaptypes,"neighbor/swap:qtype"); + double qmax, qmin; + int firstall, first; + memory->create(qtype, nswaptypes, "neighbor/swap:qtype"); for (int iswaptype = 0; iswaptype < nswaptypes; iswaptype++) { first = 1; for (int i = 0; i < atom->nlocal; i++) { @@ -272,24 +273,27 @@ void FixNeighborSwap::init() qtype[iswaptype] = atom->q[i]; first = 0; } else if (qtype[iswaptype] != atom->q[i]) - error->one(FLERR,"All atoms of a swapped type must have the same charge."); + error->one(FLERR, "All atoms of a swapped type must have the same charge."); } } } - MPI_Allreduce(&first,&firstall,1,MPI_INT,MPI_MIN,world); - if (firstall) error->all(FLERR,"At least one atom of each swapped type must be present to define charges."); + MPI_Allreduce(&first, &firstall, 1, MPI_INT, MPI_MIN, world); + if (firstall) + error->all(FLERR, + "At least one atom of each swapped type must be present to define charges."); if (first) qtype[iswaptype] = -DBL_MAX; - MPI_Allreduce(&qtype[iswaptype],&qmax,1,MPI_DOUBLE,MPI_MAX,world); + MPI_Allreduce(&qtype[iswaptype], &qmax, 1, MPI_DOUBLE, MPI_MAX, world); if (first) qtype[iswaptype] = DBL_MAX; - MPI_Allreduce(&qtype[iswaptype],&qmin,1,MPI_DOUBLE,MPI_MIN,world); - if (qmax != qmin) error->all(FLERR,"All atoms of a swapped type must have same charge."); + MPI_Allreduce(&qtype[iswaptype], &qmin, 1, MPI_DOUBLE, MPI_MIN, world); + if (qmax != qmin) error->all(FLERR, "All atoms of a swapped type must have same charge."); } } - memory->create(sqrt_mass_ratio,atom->ntypes+1,atom->ntypes+1,"neighbor/swap:sqrt_mass_ratio"); + memory->create(sqrt_mass_ratio, atom->ntypes + 1, atom->ntypes + 1, + "neighbor/swap:sqrt_mass_ratio"); for (int itype = 1; itype <= atom->ntypes; itype++) for (int jtype = 1; jtype <= atom->ntypes; jtype++) - sqrt_mass_ratio[itype][jtype] = sqrt(atom->mass[itype]/atom->mass[jtype]); + sqrt_mass_ratio[itype][jtype] = sqrt(atom->mass[itype] / atom->mass[jtype]); // check to see if itype and jtype cutoffs are the same // if not, reneighboring will be needed between swaps @@ -314,9 +318,9 @@ void FixNeighborSwap::init() if ((mask[i] == groupbit) && (mask[i] && firstgroupbit)) flag = 1; int flagall; - MPI_Allreduce(&flag,&flagall,1,MPI_INT,MPI_SUM,world); + MPI_Allreduce(&flag, &flagall, 1, MPI_INT, MPI_SUM, world); - if (flagall) error->all(FLERR,"Cannot do neighbor/swap on atoms in atom_modify first group"); + if (flagall) error->all(FLERR, "Cannot do neighbor/swap on atoms in atom_modify first group"); } } @@ -336,7 +340,7 @@ void FixNeighborSwap::pre_exchange() domain->pbc(); comm->exchange(); comm->borders(); - if (domain->triclinic) domain->lamda2x(atom->nlocal+atom->nghost); + if (domain->triclinic) domain->lamda2x(atom->nlocal + atom->nghost); if (modify->n_pre_neighbor) modify->pre_neighbor(); neighbor->build(1); @@ -369,7 +373,7 @@ int FixNeighborSwap::attempt_swap() // int nlocal = atom->nlocal; tagint *id = atom->tag; - if ( niswap == 0 ) return 0; + if (niswap == 0) return 0; // pre-swap energy @@ -385,7 +389,7 @@ int FixNeighborSwap::attempt_swap() // build nearest-neighbor list based on atom i build_i_neighbor_list(i); - if ( njswap <= 0 ) return 0; + if (njswap <= 0) return 0; // pick a neighbor atom j based on i neighbor list jtype_selected = -1; @@ -395,16 +399,14 @@ int FixNeighborSwap::attempt_swap() int jtype = jtype_selected; // Accept swap if types are equal, no change to system - if ( itype == jtype ) { - return 1; - } + if (itype == jtype) { return 1; } // swap their properties - if ( i >= 0 ) { + if (i >= 0) { atom->type[i] = jtype; if (atom->q_flag) atom->q[i] = qtype[jtype_selected]; } - if ( j >= 0 ) { + if (j >= 0) { atom->type[j] = itype; if (atom->q_flag) atom->q[j] = qtype[0]; } @@ -418,7 +420,7 @@ int FixNeighborSwap::attempt_swap() domain->pbc(); comm->exchange(); comm->borders(); - if (domain->triclinic) domain->lamda2x(atom->nlocal+atom->nghost); + if (domain->triclinic) domain->lamda2x(atom->nlocal + atom->nghost); if (modify->n_pre_neighbor) modify->pre_neighbor(); neighbor->build(1); } else { @@ -432,16 +434,15 @@ int FixNeighborSwap::attempt_swap() // if swap accepted, return 1 // if ke_flag, rescale atom velocities - if (random_equal->uniform() < - exp(beta*(energy_before - energy_after))) { + if (random_equal->uniform() < exp(beta * (energy_before - energy_after))) { update_iswap_atoms_list(); if (ke_flag) { - if ( i >= 0 ) { + if (i >= 0) { atom->v[i][0] *= sqrt_mass_ratio[itype][jtype]; atom->v[i][1] *= sqrt_mass_ratio[itype][jtype]; atom->v[i][2] *= sqrt_mass_ratio[itype][jtype]; } - if ( j >= 0 ) { + if (j >= 0) { atom->v[j][0] *= sqrt_mass_ratio[jtype][itype]; atom->v[j][1] *= sqrt_mass_ratio[jtype][itype]; atom->v[j][2] *= sqrt_mass_ratio[jtype][itype]; @@ -456,11 +457,11 @@ int FixNeighborSwap::attempt_swap() // do not need to re-call comm->borders() and rebuild neighbor list // since will be done on next cycle or in Verlet when this fix finishes - if ( i >= 0 ) { + if (i >= 0) { atom->type[i] = itype; if (atom->q_flag) atom->q[i] = qtype[0]; } - if ( j >= 0 ) { + if (j >= 0) { atom->type[j] = jtype; if (atom->q_flag) atom->q[j] = qtype[jtype_selected]; } @@ -479,16 +480,16 @@ double FixNeighborSwap::energy_full() if (modify->n_pre_force) modify->pre_force(vflag); - if (force->pair) force->pair->compute(eflag,vflag); + if (force->pair) force->pair->compute(eflag, vflag); if (atom->molecular != Atom::ATOMIC) { - if (force->bond) force->bond->compute(eflag,vflag); - if (force->angle) force->angle->compute(eflag,vflag); - if (force->dihedral) force->dihedral->compute(eflag,vflag); - if (force->improper) force->improper->compute(eflag,vflag); + if (force->bond) force->bond->compute(eflag, vflag); + if (force->angle) force->angle->compute(eflag, vflag); + if (force->dihedral) force->dihedral->compute(eflag, vflag); + if (force->improper) force->improper->compute(eflag, vflag); } - if (force->kspace) force->kspace->compute(eflag,vflag); + if (force->kspace) force->kspace->compute(eflag, vflag); if (modify->n_post_force_any) modify->post_force(vflag); @@ -507,15 +508,14 @@ int FixNeighborSwap::pick_i_swap_atom() int id_center_local = -1; int i = -1; - int iwhichglobal = static_cast (niswap*random_equal->uniform()); - if ((iwhichglobal >= niswap_before) && - (iwhichglobal < niswap_before + niswap_local)) { + int iwhichglobal = static_cast(niswap * random_equal->uniform()); + if ((iwhichglobal >= niswap_before) && (iwhichglobal < niswap_before + niswap_local)) { int iwhichlocal = iwhichglobal - niswap_before; i = local_swap_iatom_list[iwhichlocal]; id_center_local = id[i]; - MPI_Allreduce(&id[i],&id_center,1,MPI_INT,MPI_MAX,world); + MPI_Allreduce(&id[i], &id_center, 1, MPI_INT, MPI_MAX, world); } else { - MPI_Allreduce(&id[i],&id_center,1,MPI_INT,MPI_MAX,world); + MPI_Allreduce(&id[i], &id_center, 1, MPI_INT, MPI_MAX, world); } return i; @@ -530,7 +530,7 @@ int FixNeighborSwap::pick_j_swap_neighbor(int i) int jtype_selected_local = -1; // Generate random double from 0 to maximum global probability - double selected_prob = static_cast (global_probability*random_equal->uniform()); + double selected_prob = static_cast(global_probability * random_equal->uniform()); // Find which local swap atom corresponds to probability if ((selected_prob >= prev_probability) && @@ -539,28 +539,26 @@ int FixNeighborSwap::pick_j_swap_neighbor(int i) for (int n = 0; n < njswap_local; n++) { if (search_prob > local_swap_probability[n]) { search_prob -= local_swap_probability[n]; - } - else{ + } else { j = local_swap_neighbor_list[n]; jtype_selected_local = local_swap_type_list[n]; - MPI_Allreduce(&jtype_selected_local,&jtype_selected,1,MPI_INT,MPI_MAX,world); + MPI_Allreduce(&jtype_selected_local, &jtype_selected, 1, MPI_INT, MPI_MAX, world); return j; } } - error->all(FLERR,"Did not select local neighbor swap atom"); + error->all(FLERR, "Did not select local neighbor swap atom"); } - MPI_Allreduce(&jtype_selected_local,&jtype_selected,1,MPI_INT,MPI_MAX,world); + MPI_Allreduce(&jtype_selected_local, &jtype_selected, 1, MPI_INT, MPI_MAX, world); return j; } /* ---------------------------------------------------------------------- ------------------------------------------------------------------------- */ -double FixNeighborSwap::get_distance(double* i, double* j) +double FixNeighborSwap::get_distance(double *i, double *j) { - double r = sqrt(MathSpecial::square((i[0] - j[0])) + - MathSpecial::square((i[1] - j[1])) + + double r = sqrt(MathSpecial::square((i[0] - j[0])) + MathSpecial::square((i[1] - j[1])) + MathSpecial::square((i[2] - j[2]))); return r; } @@ -580,16 +578,16 @@ void FixNeighborSwap::build_i_neighbor_list(int i_center) memory->sfree(local_swap_neighbor_list); atom_swap_nmax = atom->nmax; - local_swap_neighbor_list = (int *) memory->smalloc(atom_swap_nmax*sizeof(int), - "MCSWAP:local_swap_neighbor_list"); + local_swap_neighbor_list = + (int *) memory->smalloc(atom_swap_nmax * sizeof(int), "MCSWAP:local_swap_neighbor_list"); memory->sfree(local_swap_probability); - local_swap_probability = (double *) memory->smalloc(atom_swap_nmax*sizeof(double), - "MCSWAP:local_swap_probability_list"); + local_swap_probability = (double *) memory->smalloc(atom_swap_nmax * sizeof(double), + "MCSWAP:local_swap_probability_list"); memory->sfree(local_swap_type_list); - local_swap_type_list = (int *) memory->smalloc(atom_swap_nmax*sizeof(int), - "MCSWAP:local_swap_type_list"); + local_swap_type_list = + (int *) memory->smalloc(atom_swap_nmax * sizeof(int), "MCSWAP:local_swap_type_list"); // Compute voronoi and access neighbor list @@ -605,11 +603,10 @@ void FixNeighborSwap::build_i_neighbor_list(int i_center) int temp_j = -1; // Find local voronoi entry with selected central atom - if ( (int)voro_neighbor_list[n][0] == id_center ) { + if ((int) voro_neighbor_list[n][0] == id_center) { temp_j_id = voro_neighbor_list[n][1]; temp_j = -1; - } else if ( ((int)voro_neighbor_list[n][1] == id_center) && - ( i_center < 0 ) ) { + } else if (((int) voro_neighbor_list[n][1] == id_center) && (i_center < 0)) { temp_j_id = voro_neighbor_list[n][0]; temp_j = -1; } else { @@ -618,39 +615,37 @@ void FixNeighborSwap::build_i_neighbor_list(int i_center) // Find which local atom corresponds to neighbor for (int j = 0; j < nlocal; j++) { - if ( temp_j_id == id[j] ) { + if (temp_j_id == id[j]) { temp_j = j; break; } } // If temp_j not on this processor, skip - if ( temp_j < 0 ) continue; + if (temp_j < 0) continue; if (region) { - if (region->match(x[temp_j][0],x[temp_j][1],x[temp_j][2]) == 1) { + if (region->match(x[temp_j][0], x[temp_j][1], x[temp_j][2]) == 1) { if (atom->mask[temp_j] & groupbit) { if (diff_flag) { // Calculate distance from i to each j, adjust probability of selection // Get distance if own centr atom double r = INFINITY; - if ( i_center >= 0 ) { - double r = get_distance(x[temp_j], x[i_center]); - } + if (i_center >= 0) { double r = get_distance(x[temp_j], x[i_center]); } // Get local id of ghost center atom when ghost - for (int i=nlocal; i < nlocal+nghost; i++) { - if ( (id[i] == id_center) && - (get_distance(x[temp_j], x[i]) < r) ) { - r = get_distance(x[temp_j], x[i]); + for (int i = nlocal; i < nlocal + nghost; i++) { + if ((id[i] == id_center) && (get_distance(x[temp_j], x[i]) < r)) { + r = get_distance(x[temp_j], x[i]); } } if (rates_flag) { - local_swap_probability[njswap_local] = rate_list[type[temp_j] - 1]*exp(-MathSpecial::square(r/3.0)); + local_swap_probability[njswap_local] = + rate_list[type[temp_j] - 1] * exp(-MathSpecial::square(r / 3.0)); } else { - local_swap_probability[njswap_local] = exp(-MathSpecial::square(r/3.0)); + local_swap_probability[njswap_local] = exp(-MathSpecial::square(r / 3.0)); } local_probability += local_swap_probability[njswap_local]; local_swap_type_list[njswap_local] = type[temp_j]; @@ -662,22 +657,20 @@ void FixNeighborSwap::build_i_neighbor_list(int i_center) // Calculate distance from i to each j, adjust probability of selection // Get distance if own center atom double r = INFINITY; - if ( i_center >= 0 ) { - double r = get_distance(x[temp_j], x[i_center]); - } + if (i_center >= 0) { double r = get_distance(x[temp_j], x[i_center]); } // Get local id of ghost center atom when ghost - for (int i=nlocal; i < nlocal+nghost; i++) { - if ( (id[i] == id_center) && - (get_distance(x[temp_j], x[i]) < r) ) { - r = get_distance(x[temp_j], x[i]); + for (int i = nlocal; i < nlocal + nghost; i++) { + if ((id[i] == id_center) && (get_distance(x[temp_j], x[i]) < r)) { + r = get_distance(x[temp_j], x[i]); } } if (rates_flag) { - local_swap_probability[njswap_local] = rate_list[type[temp_j] - 1]*exp(-MathSpecial::square(r/3.0)); + local_swap_probability[njswap_local] = + rate_list[type[temp_j] - 1] * exp(-MathSpecial::square(r / 3.0)); } else { - local_swap_probability[njswap_local] = exp(-MathSpecial::square(r/3.0)); + local_swap_probability[njswap_local] = exp(-MathSpecial::square(r / 3.0)); } local_probability += local_swap_probability[njswap_local]; @@ -695,20 +688,19 @@ void FixNeighborSwap::build_i_neighbor_list(int i_center) // Calculate distance from i to each j, adjust probability of selection // Get distance if own center atom double r = INFINITY; - if ( i_center >= 0 ) { - r = get_distance(x[temp_j], x[i_center]); - } + if (i_center >= 0) { r = get_distance(x[temp_j], x[i_center]); } // Get local id of ghost center atoms - for (int i=nlocal; i < nlocal+nghost; i++) { - if ( (id[i] == id_center) && (get_distance(x[temp_j], x[i]) < r) ) + for (int i = nlocal; i < nlocal + nghost; i++) { + if ((id[i] == id_center) && (get_distance(x[temp_j], x[i]) < r)) r = get_distance(x[temp_j], x[i]); } if (rates_flag) { - local_swap_probability[njswap_local] = rate_list[type[temp_j] - 1]*exp(-MathSpecial::square(r/3.0)); - } else{ - local_swap_probability[njswap_local] = exp(-MathSpecial::square(r/3.0)); + local_swap_probability[njswap_local] = + rate_list[type[temp_j] - 1] * exp(-MathSpecial::square(r / 3.0)); + } else { + local_swap_probability[njswap_local] = exp(-MathSpecial::square(r / 3.0)); } local_probability += local_swap_probability[njswap_local]; @@ -721,22 +713,20 @@ void FixNeighborSwap::build_i_neighbor_list(int i_center) // Calculate distance from i to each j, adjust probability of selection // Get distance if own center atom double r = INFINITY; - if ( i_center >= 0 ) { - double r = get_distance(x[temp_j], x[i_center]); - } + if (i_center >= 0) { double r = get_distance(x[temp_j], x[i_center]); } // Get local id of ghost center atom when ghost - for (int i=nlocal; i < nlocal+nghost; i++) { - if ( (id[i] == id_center) && - (get_distance(x[temp_j], x[i]) < r)) { - r = get_distance(x[temp_j], x[i]); + for (int i = nlocal; i < nlocal + nghost; i++) { + if ((id[i] == id_center) && (get_distance(x[temp_j], x[i]) < r)) { + r = get_distance(x[temp_j], x[i]); } } if (rates_flag) { - local_swap_probability[njswap_local] = rate_list[type[temp_j] - 1]*exp(-MathSpecial::square(r/3.0)); + local_swap_probability[njswap_local] = + rate_list[type[temp_j] - 1] * exp(-MathSpecial::square(r / 3.0)); } else { - local_swap_probability[njswap_local] = exp(-MathSpecial::square(r/3.0)); + local_swap_probability[njswap_local] = exp(-MathSpecial::square(r / 3.0)); } local_probability += local_swap_probability[njswap_local]; @@ -750,12 +740,12 @@ void FixNeighborSwap::build_i_neighbor_list(int i_center) } } - MPI_Allreduce(&njswap_local,&njswap,1,MPI_INT,MPI_SUM,world); - MPI_Scan(&njswap_local,&njswap_before,1,MPI_INT,MPI_SUM,world); + MPI_Allreduce(&njswap_local, &njswap, 1, MPI_INT, MPI_SUM, world); + MPI_Scan(&njswap_local, &njswap_before, 1, MPI_INT, MPI_SUM, world); njswap_before -= njswap_local; - MPI_Allreduce(&local_probability,&global_probability,1,MPI_DOUBLE,MPI_SUM,world); - MPI_Scan(&local_probability,&prev_probability,1,MPI_DOUBLE,MPI_SUM,world); + MPI_Allreduce(&local_probability, &global_probability, 1, MPI_DOUBLE, MPI_SUM, world); + MPI_Scan(&local_probability, &prev_probability, 1, MPI_DOUBLE, MPI_SUM, world); prev_probability -= local_probability; } @@ -772,8 +762,8 @@ void FixNeighborSwap::update_iswap_atoms_list() if (atom->nmax > atom_swap_nmax) { memory->sfree(local_swap_iatom_list); atom_swap_nmax = atom->nmax; - local_swap_iatom_list = (int *) memory->smalloc(atom_swap_nmax*sizeof(int), - "MCSWAP:local_swap_iatom_list"); + local_swap_iatom_list = + (int *) memory->smalloc(atom_swap_nmax * sizeof(int), "MCSWAP:local_swap_iatom_list"); } niswap_local = 0; @@ -781,9 +771,9 @@ void FixNeighborSwap::update_iswap_atoms_list() if (region) { for (int i = 0; i < nlocal; i++) { - if (region->match(x[i][0],x[i][1],x[i][2]) == 1) { + if (region->match(x[i][0], x[i][1], x[i][2]) == 1) { if (atom->mask[i] & groupbit) { - if (type[i] == type_list[0]) { + if (type[i] == type_list[0]) { local_swap_iatom_list[niswap_local] = i; niswap_local++; } @@ -794,7 +784,7 @@ void FixNeighborSwap::update_iswap_atoms_list() } else { for (int i = 0; i < nlocal; i++) { if (atom->mask[i] & groupbit) { - if (type[i] == type_list[0]) { + if (type[i] == type_list[0]) { local_swap_iatom_list[niswap_local] = i; niswap_local++; } @@ -802,16 +792,17 @@ void FixNeighborSwap::update_iswap_atoms_list() } } - MPI_Allreduce(&niswap_local,&niswap,1,MPI_INT,MPI_SUM,world); - MPI_Scan(&niswap_local,&niswap_before,1,MPI_INT,MPI_SUM,world); + MPI_Allreduce(&niswap_local, &niswap, 1, MPI_INT, MPI_SUM, world); + MPI_Scan(&niswap_local, &niswap_before, 1, MPI_INT, MPI_SUM, world); niswap_before -= niswap_local; } /* ---------------------------------------------------------------------- */ -int FixNeighborSwap::pack_forward_comm(int n, int *list, double *buf, int /*pbc_flag*/, int * /*pbc*/) +int FixNeighborSwap::pack_forward_comm(int n, int *list, double *buf, int /*pbc_flag*/, + int * /*pbc*/) { - int i,j,m; + int i, j, m; int *type = atom->type; double *q = atom->q; @@ -838,7 +829,7 @@ int FixNeighborSwap::pack_forward_comm(int n, int *list, double *buf, int /*pbc_ void FixNeighborSwap::unpack_forward_comm(int n, int first, double *buf) { - int i,m,last; + int i, m, last; int *type = atom->type; double *q = atom->q; @@ -848,12 +839,11 @@ void FixNeighborSwap::unpack_forward_comm(int n, int first, double *buf) if (atom->q_flag) { for (i = first; i < last; i++) { - type[i] = static_cast (buf[m++]); + type[i] = static_cast(buf[m++]); q[i] = buf[m++]; } } else { - for (i = first; i < last; i++) - type[i] = static_cast (buf[m++]); + for (i = first; i < last; i++) type[i] = static_cast(buf[m++]); } } @@ -874,7 +864,7 @@ double FixNeighborSwap::compute_vector(int n) double FixNeighborSwap::memory_usage() { - double bytes = (double)atom_swap_nmax * sizeof(int); + double bytes = (double) atom_swap_nmax * sizeof(int); return bytes; } @@ -895,8 +885,8 @@ void FixNeighborSwap::write_restart(FILE *fp) if (comm->me == 0) { int size = n * sizeof(double); - fwrite(&size,sizeof(int),1,fp); - fwrite(list,sizeof(double),n,fp); + fwrite(&size, sizeof(int), 1, fp); + fwrite(list, sizeof(double), n, fp); } } @@ -909,10 +899,10 @@ void FixNeighborSwap::restart(char *buf) int n = 0; double *list = (double *) buf; - seed = static_cast (list[n++]); + seed = static_cast(list[n++]); random_equal->reset(seed); - seed = static_cast (list[n++]); + seed = static_cast(list[n++]); random_unequal->reset(seed); next_reneighbor = (bigint) ubuf(list[n++]).i; @@ -922,5 +912,5 @@ void FixNeighborSwap::restart(char *buf) bigint ntimestep_restart = (bigint) ubuf(list[n++]).i; if (ntimestep_restart != update->ntimestep) - error->all(FLERR,"Must not reset timestep when restarting fix neighbor/swap"); + error->all(FLERR, "Must not reset timestep when restarting fix neighbor/swap"); } diff --git a/src/MC/fix_neighbor_swap.h b/src/MC/fix_neighbor_swap.h index df75505a6d..826ab96a65 100755 --- a/src/MC/fix_neighbor_swap.h +++ b/src/MC/fix_neighbor_swap.h @@ -40,17 +40,17 @@ class FixNeighborSwap : public Fix { private: int nevery, seed; - int ke_flag; // yes = conserve ke, no = do not conserve ke - int diff_flag; // yes = simulate diffusion of central atom, no = swap only to certain types - int rates_flag; // yes = use modified type rates, no = swap rates are equivilent across types - int voro_flag; // yes = use given voronoi calculation, no = use internal voronoi calculation + int ke_flag; // yes = conserve ke, no = do not conserve ke + int diff_flag; // yes = simulate diffusion of central atom, no = swap only to certain types + int rates_flag; // yes = use modified type rates, no = swap rates are equivilent across types + int voro_flag; // yes = use given voronoi calculation, no = use internal voronoi calculation int ncycles; int niswap, njswap; // # of i,j swap atoms on all procs int niswap_local, njswap_local; // # of swap atoms on this proc int niswap_before, njswap_before; // # of swap atoms on procs < this proc // int global_i_ID; // global id of selected i atom - class Region *region; // swap region - char *idregion; // swap region id + class Region *region; // swap region + char *idregion; // swap region id int nswaptypes; int jtype_selected; @@ -68,19 +68,18 @@ class FixNeighborSwap : public Fix { int atom_swap_nmax; double beta; - double local_probability; // Total swap probability stored on this proc - double global_probability; // Total swap probability across all proc - double prev_probability; // Swap probability on proc < this proc + double local_probability; // Total swap probability stored on this proc + double global_probability; // Total swap probability across all proc + double prev_probability; // Swap probability on proc < this proc double *qtype; double energy_stored; double **sqrt_mass_ratio; double **voro_neighbor_list; int *local_swap_iatom_list; int *local_swap_neighbor_list; - int *local_swap_type_list; // Type list index of atoms stored on this proc + int *local_swap_type_list; // Type list index of atoms stored on this proc double *local_swap_probability; - class RanPark *random_equal; class RanPark *random_unequal;