diff --git a/src/atom_vec_sphere_bpm.cpp b/src/atom_vec_sphere_bpm.cpp deleted file mode 100644 index 3fe372e340..0000000000 --- a/src/atom_vec_sphere_bpm.cpp +++ /dev/null @@ -1,254 +0,0 @@ -// clang-format off -/* ---------------------------------------------------------------------- - LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator - https://www.lammps.org/, Sandia National Laboratories - Steve Plimpton, sjplimp@sandia.gov - - Copyright (2003) Sandia Corporation. Under the terms of Contract - DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains - certain rights in this software. This software is distributed under - the GNU General Public License. - - See the README file in the top-level LAMMPS directory. -------------------------------------------------------------------------- */ - -#include "atom_vec_sphere_bpm.h" - -#include "atom.h" -#include "error.h" -#include "fix.h" -#include "fix_adapt.h" -#include "math_const.h" -#include "modify.h" -#include "utils.h" - -#include - -using namespace LAMMPS_NS; -using namespace MathConst; - -/* ---------------------------------------------------------------------- */ - -AtomVecSphereBPM::AtomVecSphereBPM(LAMMPS *lmp) : AtomVec(lmp) -{ - mass_type = PER_ATOM; - molecular = Atom::MOLECULAR; - bonds_allow = 1; - - atom->molecule_flag = 1; - atom->sphere_flag = 1; - atom->radius_flag = atom->rmass_flag = atom->omega_flag = - atom->torque_flag = atom->quat_flag = 1; - - // strings with peratom variables to include in each AtomVec method - // strings cannot contain fields in corresponding AtomVec default strings - // order of fields in a string does not matter - // except: fields_data_atom & fields_data_vel must match data file - - fields_grow = (char *) - "molecule num_bond bond_type bond_atom nspecial special radius rmass omega torque quat"; - fields_copy = (char *) - "molecule num_bond bond_type bond_atom nspecial special radius rmass omega quat"; - fields_comm = (char *) ""; - fields_comm_vel = (char *) "omega quat"; - fields_reverse = (char *) "torque"; - fields_border = (char *) "molecule radius rmass"; - fields_border_vel = (char *) "molecule radius rmass omega quat"; - fields_exchange = (char *) - "molecule num_bond bond_type bond_atom nspecial special radius rmass omega quat"; - fields_restart = (char *) - "molecule num_bond bond_type bond_atom radius rmass omega quat"; - fields_create = (char *) - "molecule num_bond nspecial radius rmass omega quat"; - fields_data_atom = (char *) "id molecule type radius rmass x"; - fields_data_vel = (char *) "id v omega"; - - bond_per_atom = 0; - bond_negative = NULL; -} - -/* ---------------------------------------------------------------------- - process sub-style args - optional arg = 0/1 for static/dynamic particle radii -------------------------------------------------------------------------- */ - -void AtomVecSphereBPM::process_args(int narg, char **arg) -{ - if (narg != 0 && narg != 1) - error->all(FLERR,"Illegal atom_style sphere/bpm command"); - - radvary = 0; - if (narg == 1) { - radvary = utils::numeric(FLERR,arg[0],true,lmp); - if (radvary < 0 || radvary > 1) - error->all(FLERR,"Illegal atom_style sphere/bpm command"); - } - - // dynamic particle radius and mass must be communicated every step - - if (radvary) { - fields_comm = (char *) "radius rmass"; - fields_comm_vel = (char *) "radius rmass omega"; - } - - // delay setting up of fields until now - - setup_fields(); -} - -/* ---------------------------------------------------------------------- */ - -void AtomVecSphereBPM::init() -{ - AtomVec::init(); - - // check if optional radvary setting should have been set to 1 - - for (int i = 0; i < modify->nfix; i++) - if (strcmp(modify->fix[i]->style,"adapt") == 0) { - FixAdapt *fix = (FixAdapt *) modify->fix[i]; - if (fix->diamflag && radvary == 0) - error->all(FLERR,"Fix adapt changes particle radii " - "but atom_style sphere is not dynamic"); - } -} - -/* ---------------------------------------------------------------------- - set local copies of all grow ptrs used by this class, except defaults - needed in replicate when 2 atom classes exist and it calls pack_restart() -------------------------------------------------------------------------- */ - -void AtomVecSphereBPM::grow_pointers() -{ - radius = atom->radius; - rmass = atom->rmass; - omega = atom->omega; - quat = atom->quat; - - num_bond = atom->num_bond; - bond_type = atom->bond_type; - nspecial = atom->nspecial; -} - - -/* ---------------------------------------------------------------------- - initialize non-zero atom quantities -------------------------------------------------------------------------- */ - -void AtomVecSphereBPM::create_atom_post(int ilocal) -{ - radius[ilocal] = 0.5; - rmass[ilocal] = 4.0*MY_PI/3.0 * 0.5*0.5*0.5; - - quat[ilocal][0] = 1.0; - quat[ilocal][1] = 0.0; - quat[ilocal][2] = 0.0; - quat[ilocal][3] = 0.0; - -} - -/* ---------------------------------------------------------------------- - modify values for AtomVec::pack_restart() to pack -------------------------------------------------------------------------- */ - -void AtomVecSphereBPM::pack_restart_pre(int ilocal) -{ - // insure bond_negative vector is needed length - - if (bond_per_atom < atom->bond_per_atom) { - delete [] bond_negative; - bond_per_atom = atom->bond_per_atom; - bond_negative = new int[bond_per_atom]; - } - - // flip any negative types to positive and flag which ones - - any_bond_negative = 0; - for (int m = 0; m < num_bond[ilocal]; m++) { - if (bond_type[ilocal][m] < 0) { - bond_negative[m] = 1; - bond_type[ilocal][m] = -bond_type[ilocal][m]; - any_bond_negative = 1; - } else bond_negative[m] = 0; - } -} - -/* ---------------------------------------------------------------------- - unmodify values packed by AtomVec::pack_restart() -------------------------------------------------------------------------- */ - -void AtomVecSphereBPM::pack_restart_post(int ilocal) -{ - // restore the flagged types to their negative values - - if (any_bond_negative) { - for (int m = 0; m < num_bond[ilocal]; m++) - if (bond_negative[m]) bond_type[ilocal][m] = -bond_type[ilocal][m]; - } -} - -/* ---------------------------------------------------------------------- - initialize other atom quantities after AtomVec::unpack_restart() -------------------------------------------------------------------------- */ - -void AtomVecSphereBPM::unpack_restart_init(int ilocal) -{ - nspecial[ilocal][0] = 0; - nspecial[ilocal][1] = 0; - nspecial[ilocal][2] = 0; -} - -/* ---------------------------------------------------------------------- - modify what AtomVec::data_atom() just unpacked - or initialize other atom quantities -------------------------------------------------------------------------- */ - -void AtomVecSphereBPM::data_atom_post(int ilocal) -{ - radius_one = 0.5 * atom->radius[ilocal]; - radius[ilocal] = radius_one; - if (radius_one > 0.0) - rmass[ilocal] *= 4.0*MY_PI/3.0 * radius_one*radius_one*radius_one; - - if (rmass[ilocal] <= 0.0) - error->one(FLERR,"Invalid density in Atoms section of data file"); - - omega[ilocal][0] = 0.0; - omega[ilocal][1] = 0.0; - omega[ilocal][2] = 0.0; - - quat[ilocal][0] = 1.0; - quat[ilocal][1] = 0.0; - quat[ilocal][2] = 0.0; - quat[ilocal][3] = 0.0; - - num_bond[ilocal] = 0; - nspecial[ilocal][0] = 0; - nspecial[ilocal][1] = 0; - nspecial[ilocal][2] = 0; -} - -/* ---------------------------------------------------------------------- - modify values for AtomVec::pack_data() to pack -------------------------------------------------------------------------- */ - -void AtomVecSphereBPM::pack_data_pre(int ilocal) -{ - radius_one = radius[ilocal]; - rmass_one = rmass[ilocal]; - - radius[ilocal] *= 2.0; - if (radius_one!= 0.0) - rmass[ilocal] = - rmass_one / (4.0*MY_PI/3.0 * radius_one*radius_one*radius_one); -} - -/* ---------------------------------------------------------------------- - unmodify values packed by AtomVec::pack_data() -------------------------------------------------------------------------- */ - -void AtomVecSphereBPM::pack_data_post(int ilocal) -{ - radius[ilocal] = radius_one; - rmass[ilocal] = rmass_one; -} diff --git a/src/atom_vec_sphere_bpm.h b/src/atom_vec_sphere_bpm.h deleted file mode 100644 index 2a8aa36d32..0000000000 --- a/src/atom_vec_sphere_bpm.h +++ /dev/null @@ -1,83 +0,0 @@ -/* -*- c++ -*- ---------------------------------------------------------- - LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator - https://www.lammps.org/, Sandia National Laboratories - Steve Plimpton, sjplimp@sandia.gov - - Copyright (2003) Sandia Corporation. Under the terms of Contract - DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains - certain rights in this software. This software is distributed under - the GNU General Public License. - - See the README file in the top-level LAMMPS directory. -------------------------------------------------------------------------- */ - -#ifdef ATOM_CLASS -// clang-format off -AtomStyle(sphere/bpm,AtomVecSphereBPM) -// clang-format on -#else - -#ifndef LMP_ATOM_VEC_SPHERE_BPM_H -#define LMP_ATOM_VEC_SPHERE_BPM_H - -#include "atom_vec.h" - -namespace LAMMPS_NS { - -class AtomVecSphereBPM : public AtomVec { - public: - AtomVecSphereBPM(class LAMMPS *); - void process_args(int, char **) override; - void init() override; - - void grow_pointers() override; - void create_atom_post(int) override; - void pack_restart_pre(int) override; - void pack_restart_post(int) override; - void unpack_restart_init(int) override; - void data_atom_post(int) override; - void pack_data_pre(int) override; - void pack_data_post(int) override; - - - private: - int *num_bond; - int **bond_type; - int **nspecial; - - double *radius,*rmass; - double **omega, **torque, **quat; - - int any_bond_negative; - int bond_per_atom; - int *bond_negative; - - int radvary; - double radius_one,rmass_one; -}; - -} // namespace LAMMPS_NS - -#endif -#endif - -/* ERROR/WARNING messages: - -E: Per-processor system is too big - -The number of owned atoms plus ghost atoms on a single -processor must fit in 32-bit integer. - -E: Invalid atom type in Atoms section of data file - -Atom types must range from 1 to specified # of types. - -E: Invalid radius in Atoms section of data file - -Radius must be >= 0.0. - -E: Invalid density in Atoms section of data file - -Density value cannot be <= 0.0. - -*/ diff --git a/src/fix_nve_sphere_bpm.cpp b/src/fix_nve_sphere_bpm.cpp deleted file mode 100644 index 66801cc0b6..0000000000 --- a/src/fix_nve_sphere_bpm.cpp +++ /dev/null @@ -1,161 +0,0 @@ -// clang-format off -/* ---------------------------------------------------------------------- - LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator - https://www.lammps.org/, Sandia National Laboratories - Steve Plimpton, sjplimp@sandia.gov - - Copyright (2003) Sandia Corporation. Under the terms of Contract - DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains - certain rights in this software. This software is distributed under - the GNU General Public License. - - See the README file in the top-level LAMMPS directory. -------------------------------------------------------------------------- */ - -#include "fix_nve_sphere_bpm.h" - -#include "atom.h" -#include "atom_vec.h" -#include "domain.h" -#include "error.h" -#include "force.h" -#include "math_extra.h" - -#include -#include - -using namespace LAMMPS_NS; -using namespace FixConst; -using namespace MathExtra; - - -/* ---------------------------------------------------------------------- */ - -FixNVESphereBPM::FixNVESphereBPM(LAMMPS *lmp, int narg, char **arg) : - FixNVE(lmp, narg, arg) -{ - if (narg < 3) error->all(FLERR,"Illegal fix nve/sphere/bpm command"); - - time_integrate = 1; - - // process extra keywords - // inertia = moment of inertia prefactor for sphere or disc - - inertia = 0.4; - - int iarg = 3; - while (iarg < narg) { - if (strcmp(arg[iarg],"disc")==0) { - inertia = 0.5; - if (domain->dimension != 2) - error->all(FLERR,"Fix nve/sphere/bpm disc requires 2d simulation"); - iarg++; - } - else error->all(FLERR,"Illegal fix nve/sphere/bpm command"); - } - - inv_inertia = 1.0/inertia; - - // error checks - - if (!atom->quat_flag || !atom->sphere_flag) - error->all(FLERR,"Fix nve/sphere/bpm requires atom style sphere/bpm"); -} - -/* ---------------------------------------------------------------------- */ - -void FixNVESphereBPM::init() -{ - FixNVE::init(); - - // check that all particles are finite-size spheres - // no point particles allowed - - double *radius = atom->radius; - int *mask = atom->mask; - int nlocal = atom->nlocal; - - for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit) - if (radius[i] == 0.0) - error->one(FLERR,"Fix nve/sphere/bpm requires extended particles"); -} - -/* ---------------------------------------------------------------------- */ - -void FixNVESphereBPM::initial_integrate(int /*vflag*/) -{ - double dtq,dtfm,dtirotate,particle_inertia; - - double **x = atom->x; - double **v = atom->v; - double **f = atom->f; - double **omega = atom->omega; - double **torque = atom->torque; - double **quat = atom->quat; - double *radius = atom->radius; - double *rmass = atom->rmass; - int *mask = atom->mask; - int nlocal = atom->nlocal; - if (igroup == atom->firstgroup) nlocal = atom->nfirst; - - // set timestep here since dt may have changed or come via rRESPA - dtq = 0.5 * dtv; - - // update v,x,omega,quat for all particles - // d_omega/dt = torque / inertia - - for (int i = 0; i < nlocal; i++) { - if (mask[i] & groupbit) { - dtfm = dtf / rmass[i]; - v[i][0] += dtfm * f[i][0]; - v[i][1] += dtfm * f[i][1]; - v[i][2] += dtfm * f[i][2]; - x[i][0] += dtv * v[i][0]; - x[i][1] += dtv * v[i][1]; - x[i][2] += dtv * v[i][2]; - - particle_inertia = inertia*(radius[i]*radius[i]*rmass[i]); - dtirotate = dtf / particle_inertia; - omega[i][0] += dtirotate * torque[i][0]; - omega[i][1] += dtirotate * torque[i][1]; - omega[i][2] += dtirotate * torque[i][2]; - - MathExtra::richardson_sphere(quat[i],omega[i],dtq); - } - } -} - -/* ---------------------------------------------------------------------- */ - -void FixNVESphereBPM::final_integrate() -{ - double dtfm,dtirotate,particle_inertia; - - double **v = atom->v; - double **f = atom->f; - double **omega = atom->omega; - double **torque = atom->torque; - double *rmass = atom->rmass; - double *radius = atom->radius; - int *mask = atom->mask; - int nlocal = atom->nlocal; - if (igroup == atom->firstgroup) nlocal = atom->nfirst; - - // update v,omega for all particles - // d_omega/dt = torque / inertia - - for (int i = 0; i < nlocal; i++) - if (mask[i] & groupbit) { - dtfm = dtf / rmass[i]; - v[i][0] += dtfm * f[i][0]; - v[i][1] += dtfm * f[i][1]; - v[i][2] += dtfm * f[i][2]; - - particle_inertia = inertia*(radius[i]*radius[i]*rmass[i]); - dtirotate = dtf / particle_inertia; - omega[i][0] += dtirotate * torque[i][0]; - omega[i][1] += dtirotate * torque[i][1]; - omega[i][2] += dtirotate * torque[i][2]; - } -} diff --git a/src/fix_nve_sphere_bpm.h b/src/fix_nve_sphere_bpm.h deleted file mode 100644 index ee6f1b452c..0000000000 --- a/src/fix_nve_sphere_bpm.h +++ /dev/null @@ -1,71 +0,0 @@ -/* -*- c++ -*- ---------------------------------------------------------- - LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator - https://www.lammps.org/, Sandia National Laboratories - Steve Plimpton, sjplimp@sandia.gov - - Copyright (2003) Sandia Corporation. Under the terms of Contract - DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains - certain rights in this software. This software is distributed under - the GNU General Public License. - - See the README file in the top-level LAMMPS directory. -------------------------------------------------------------------------- */ - -#ifdef FIX_CLASS -// clang-format off -FixStyle(nve/sphere/bpm,FixNVESphereBPM) -// clang-format on -#else - -#ifndef LMP_FIX_NVE_SPHERE_BPM_H -#define LMP_FIX_NVE_SPHERE_BPM_H - -#include "fix_nve.h" - -namespace LAMMPS_NS { - -class FixNVESphereBPM : public FixNVE { - public: - FixNVESphereBPM(class LAMMPS *, int, char **); - - void init() override; - void initial_integrate(int) override; - void final_integrate() override; - - protected: - double inertia, inv_inertia; - int extra; - int dlm; -}; - -} // namespace LAMMPS_NS - -#endif -#endif - -/* ERROR/WARNING messages: - -E: Illegal ... command - -Self-explanatory. Check the input script syntax and compare to the -documentation for the command. You can use -echo screen as a -command-line option when running LAMMPS to see the offending line. - -E: Fix nve/sphere/bpm disc requires 2d simulation - -UNDOCUMENTED - -E: Fix nve/sphere/bpm requires atom style sphere/bpm - -Self-explanatory. - -E: Fix nve/sphere/bpm update dipole requires atom attribute mu - -An atom style with this attribute is needed. - -E: Fix nve/sphere/bpm requires extended particles - -This fix can only be used for particles of a finite size. - - -*/