git-svn-id: svn://svn.icms.temple.edu/lammps-ro/trunk@6375 f3b2605a-c512-4ea7-a41b-209d697bcdaa
This commit is contained in:
sjplimp
37
src/USER-AWPMD/Install.sh
Normal file
37
src/USER-AWPMD/Install.sh
Normal file
@ -0,0 +1,37 @@
|
||||
# Install/unInstall package files in LAMMPS
|
||||
# edit Makefile.package to include/exclude ATC library
|
||||
|
||||
if (test $1 = 1) then
|
||||
|
||||
if (test -e ../Makefile.package) then
|
||||
sed -i -e 's/[^ \t]*awpmd[^ \t]* //g' ../Makefile.package
|
||||
sed -i -e 's|^PKG_INC =[ \t]*|&-I../../lib/awpmd/ivutils/include -I../../lib/awpmd/systems/interact |' ../Makefile.package
|
||||
sed -i -e 's|^PKG_PATH =[ \t]*|&-L../../lib/awpmd |' ../Makefile.package
|
||||
sed -i -e 's|^PKG_LIB =[ \t]*|&-lawpmd |' ../Makefile.package
|
||||
sed -i -e 's|^PKG_SYSPATH =[ \t]*|&$(user-awpmd_SYSPATH) |' ../Makefile.package
|
||||
sed -i -e 's|^PKG_SYSLIB =[ \t]*|&$(user-awpmd_SYSLIB) |' ../Makefile.package
|
||||
fi
|
||||
|
||||
cp atom_vec_wavepacket.cpp ..
|
||||
cp fix_nve_awpmd.cpp ..
|
||||
cp pair_awpmd_cut.cpp ..
|
||||
|
||||
cp atom_vec_wavepacket.h ..
|
||||
cp fix_nve_awpmd.h ..
|
||||
cp pair_awpmd_cut.h ..
|
||||
|
||||
elif (test $1 = 0) then
|
||||
|
||||
if (test -e ../Makefile.package) then
|
||||
sed -i -e 's/[^ \t]*awpmd[^ \t]* //g' ../Makefile.package
|
||||
fi
|
||||
|
||||
rm -f ../atom_vec_wavepacket.cpp
|
||||
rm -f ../fix_nve_awpmd.cpp
|
||||
rm -f ../pair_awpmd_cut.cpp
|
||||
|
||||
rm -f ../atom_vec_wavepacket.h
|
||||
rm -f ../fix_nve_awpmd.h
|
||||
rm -f ../pair_awpmd_cut.h
|
||||
|
||||
fi
|
||||
24
src/USER-AWPMD/README
Normal file
24
src/USER-AWPMD/README
Normal file
@ -0,0 +1,24 @@
|
||||
The files in this directory are a user-contributed package for LAMMPS.
|
||||
|
||||
The person who created these files is Ilya Valuev
|
||||
(valuev@physik.hu-berlin.de). Contact him directly if you have
|
||||
questions.
|
||||
|
||||
PACKAGE DESCRIPTION:
|
||||
|
||||
Contains a LAMMPS implementation of the Antisymmetrized Wave Packet
|
||||
Molecular Dynamics (AWPMD) method.
|
||||
|
||||
INSTALLATION:
|
||||
|
||||
- compile the awpmd library in lib/awpmd
|
||||
- follow a normal LAMMPS package installation: make yes-user-awpmd
|
||||
|
||||
OTHERS FILES INCLUDED:
|
||||
|
||||
User examples are under examples/USER/awpmd
|
||||
|
||||
ACKNOWLEDGMENTS:
|
||||
|
||||
Thanks to Steve Plimpton and Aidan Thompson for their help in
|
||||
adaptation of the AWPMD code to LAMMPS.
|
||||
1009
src/USER-AWPMD/atom_vec_wavepacket.cpp
Normal file
1009
src/USER-AWPMD/atom_vec_wavepacket.cpp
Normal file
File diff suppressed because it is too large
Load Diff
98
src/USER-AWPMD/atom_vec_wavepacket.h
Normal file
98
src/USER-AWPMD/atom_vec_wavepacket.h
Normal file
@ -0,0 +1,98 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, 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.
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Contributing author: Ilya Valuev (JIHT RAS)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
|
||||
#ifdef ATOM_CLASS
|
||||
|
||||
AtomStyle(wavepacket,AtomVecWavepacket)
|
||||
|
||||
#else
|
||||
|
||||
#ifndef LMP_ATOM_VEC_WAVEPACKET_H
|
||||
#define LMP_ATOM_VEC_WAVEPACKET_H
|
||||
|
||||
#include "atom_vec.h"
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
|
||||
class AtomVecWavepacket : public AtomVec {
|
||||
public:
|
||||
AtomVecWavepacket(class LAMMPS *, int, char **);
|
||||
~AtomVecWavepacket() {}
|
||||
void grow(int);
|
||||
void grow_reset();
|
||||
void copy(int, int, int);
|
||||
int pack_comm(int, int *, double *, int, int *);
|
||||
int pack_comm_vel(int, int *, double *, int, int *);
|
||||
int pack_comm_hybrid(int, int *, double *);
|
||||
void unpack_comm(int, int, double *);
|
||||
void unpack_comm_vel(int, int, double *);
|
||||
int unpack_comm_hybrid(int, int, double *);
|
||||
int pack_reverse(int, int, double *);
|
||||
int pack_reverse_hybrid(int, int, double *);
|
||||
void unpack_reverse(int, int *, double *);
|
||||
int unpack_reverse_hybrid(int, int *, double *);
|
||||
int pack_border(int, int *, double *, int, int *);
|
||||
int pack_border_vel(int, int *, double *, int, int *);
|
||||
int pack_border_hybrid(int, int *, double *);
|
||||
void unpack_border(int, int, double *);
|
||||
void unpack_border_vel(int, int, double *);
|
||||
int unpack_border_hybrid(int, int, double *);
|
||||
int pack_exchange(int, double *);
|
||||
int unpack_exchange(double *);
|
||||
int size_restart();
|
||||
int pack_restart(int, double *);
|
||||
int unpack_restart(double *);
|
||||
void create_atom(int, double *);
|
||||
void data_atom(double *, int, char **);
|
||||
int data_atom_hybrid(int, char **);
|
||||
void data_vel(int, char **);
|
||||
int data_vel_hybrid(int, char **);
|
||||
bigint memory_usage();
|
||||
|
||||
private:
|
||||
int *tag,*type,*mask,*image;
|
||||
double **x,**v,**f;
|
||||
|
||||
///\en spin: -1 or 1 for electron, 0 for ion (compatible with eff)
|
||||
int *spin;
|
||||
///\en charge: must be specified in the corresponding units (-1 for electron in real units, eff compatible)
|
||||
double *q;
|
||||
///\en width of the wavepacket (compatible with eff)
|
||||
double *eradius;
|
||||
///\en width velocity for the wavepacket (compatible with eff)
|
||||
double *ervel;
|
||||
///\en (generalized) force on width (compatible with eff)
|
||||
double *erforce;
|
||||
|
||||
// AWPMD- specific:
|
||||
///\en electron tag: must be the same for the WPs belonging to the same electron
|
||||
int *etag;
|
||||
///\en wavepacket split coeffcients: cre, cim, size is 2*N
|
||||
double *cs;
|
||||
///\en force on wavepacket split coeffcients: re, im, size is 2*N
|
||||
double *csforce;
|
||||
///\en (generalized) force on velocity, size is 3*N
|
||||
double *vforce;
|
||||
///\en (generalized) force on radius velocity, size is N
|
||||
double *ervelforce;
|
||||
};
|
||||
|
||||
}
|
||||
|
||||
#endif
|
||||
#endif
|
||||
155
src/USER-AWPMD/fix_nve_awpmd.cpp
Normal file
155
src/USER-AWPMD/fix_nve_awpmd.cpp
Normal file
@ -0,0 +1,155 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, 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.
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Contributing author: Ilya Valuev (JIHT RAS)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#include "math.h"
|
||||
#include "stdio.h"
|
||||
#include "string.h"
|
||||
#include "fix_nve_awpmd.h"
|
||||
#include "atom.h"
|
||||
#include "force.h"
|
||||
#include "update.h"
|
||||
#include "respa.h"
|
||||
#include "error.h"
|
||||
#include "math.h"
|
||||
|
||||
#include "TCP/wpmd_split.h"
|
||||
|
||||
|
||||
|
||||
|
||||
using namespace LAMMPS_NS;
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
FixNVEAwpmd::FixNVEAwpmd(LAMMPS *lmp, int narg, char **arg) :
|
||||
Fix(lmp, narg, arg)
|
||||
{
|
||||
if (!atom->wavepacket_flag)
|
||||
error->all("Fix nve/awpmd requires atom style wavepacket");
|
||||
//if (!atom->mass_type != 1)
|
||||
// error->all("Fix nve/awpmd requires per type mass");
|
||||
|
||||
time_integrate = 1;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
int FixNVEAwpmd::setmask()
|
||||
{
|
||||
int mask = 0;
|
||||
mask |= INITIAL_INTEGRATE;
|
||||
mask |= FINAL_INTEGRATE;
|
||||
mask |= INITIAL_INTEGRATE_RESPA;
|
||||
mask |= FINAL_INTEGRATE_RESPA;
|
||||
return mask;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixNVEAwpmd::init()
|
||||
{
|
||||
dtv = update->dt;
|
||||
dtf = 0.5 * update->dt * force->ftm2v;
|
||||
|
||||
if (strcmp(update->integrate_style,"respa") == 0)
|
||||
step_respa = ((Respa *) update->integrate)->step;
|
||||
|
||||
awpmd_pair=(PairAWPMDCut *)force->pair;
|
||||
awpmd_pair->wpmd->norm_needed=1;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
allow for only per-type mass
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void FixNVEAwpmd::initial_integrate(int vflag)
|
||||
{
|
||||
|
||||
|
||||
// update v,vr and x,radius of atoms in group
|
||||
|
||||
double **x = atom->x;
|
||||
double *eradius = atom->eradius;
|
||||
double **v = atom->v;
|
||||
double *ervel = atom->ervel;
|
||||
double **f = atom->f;
|
||||
double *erforce = atom->erforce;
|
||||
double *vforce=atom->vforce;
|
||||
double *ervelforce=atom->ervelforce;
|
||||
double *cs=atom->cs;
|
||||
double *csforce=atom->csforce;
|
||||
|
||||
|
||||
double *mass = atom->mass;
|
||||
int *spin = atom->spin;
|
||||
int *type = atom->type;
|
||||
int *mask = atom->mask;
|
||||
int nlocal = atom->nlocal;
|
||||
if (igroup == atom->firstgroup) nlocal = atom->nfirst;
|
||||
|
||||
// x + dt * [v + 0.5 * dt * (f / m)];
|
||||
|
||||
// simple Euler update
|
||||
for (int i = 0; i < nlocal; i++) {
|
||||
if (mask[i] & groupbit) {
|
||||
double dtfm = dtf / mass[type[i]];
|
||||
double dtfmr=dtfm;
|
||||
for(int j=0;j<3;j++){
|
||||
x[i][j] += dtv*vforce[3*i+j];
|
||||
v[i][j] += dtfm*f[i][j];
|
||||
}
|
||||
eradius[i]+= dtv*ervelforce[i];
|
||||
ervel[i] += dtfmr*erforce[i];
|
||||
}
|
||||
}
|
||||
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixNVEAwpmd::final_integrate(){}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixNVEAwpmd::initial_integrate_respa(int vflag, int ilevel, int iloop)
|
||||
{
|
||||
dtv = step_respa[ilevel];
|
||||
dtf = 0.5 * step_respa[ilevel] * force->ftm2v;
|
||||
|
||||
// innermost level - NVE update of v and x
|
||||
// all other levels - NVE update of v
|
||||
|
||||
if (ilevel == 0) initial_integrate(vflag);
|
||||
else final_integrate();
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixNVEAwpmd::final_integrate_respa(int ilevel, int iloop)
|
||||
{
|
||||
dtf = 0.5 * step_respa[ilevel] * force->ftm2v;
|
||||
final_integrate();
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void FixNVEAwpmd::reset_dt()
|
||||
{
|
||||
dtv = update->dt;
|
||||
dtf = 0.5 * update->dt * force->ftm2v;
|
||||
}
|
||||
|
||||
54
src/USER-AWPMD/fix_nve_awpmd.h
Normal file
54
src/USER-AWPMD/fix_nve_awpmd.h
Normal file
@ -0,0 +1,54 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, 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.
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Contributing author: Ilya Valuev (JIHT RAS)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#ifdef FIX_CLASS
|
||||
|
||||
FixStyle(nve/awpmd,FixNVEAwpmd)
|
||||
|
||||
#else
|
||||
|
||||
#ifndef LMP_FIX_NVE_awpmd_H
|
||||
#define LMP_FIX_NVE_awpmd_H
|
||||
|
||||
#include "fix.h"
|
||||
#include "pair_awpmd_cut.h"
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
|
||||
class FixNVEAwpmd : public Fix {
|
||||
public:
|
||||
FixNVEAwpmd(class LAMMPS *, int, char **);
|
||||
int setmask();
|
||||
virtual void init();
|
||||
virtual void initial_integrate(int);
|
||||
virtual void final_integrate();
|
||||
void initial_integrate_respa(int, int, int);
|
||||
void final_integrate_respa(int, int);
|
||||
void reset_dt();
|
||||
|
||||
protected:
|
||||
double dtv,dtf;
|
||||
double *step_respa;
|
||||
int mass_require;
|
||||
|
||||
PairAWPMDCut *awpmd_pair;
|
||||
};
|
||||
|
||||
}
|
||||
|
||||
#endif
|
||||
#endif
|
||||
757
src/USER-AWPMD/pair_awpmd_cut.cpp
Normal file
757
src/USER-AWPMD/pair_awpmd_cut.cpp
Normal file
@ -0,0 +1,757 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, 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.
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
Contributing author: Ilya Valuev
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
#include "math.h"
|
||||
#include "stdio.h"
|
||||
#include "stdlib.h"
|
||||
#include "string.h"
|
||||
#include "pair_awpmd_cut.h"
|
||||
#include "atom.h"
|
||||
#include "update.h"
|
||||
#include "min.h"
|
||||
#include "domain.h"
|
||||
#include "comm.h"
|
||||
#include "force.h"
|
||||
#include "neighbor.h"
|
||||
#include "neigh_list.h"
|
||||
#include "neigh_request.h"
|
||||
#include "memory.h"
|
||||
#include "error.h"
|
||||
|
||||
#include "TCP/wpmd_split.h"
|
||||
|
||||
|
||||
using namespace LAMMPS_NS;
|
||||
|
||||
#define MIN(a,b) ((a) < (b) ? (a) : (b))
|
||||
#define MAX(a,b) ((a) > (b) ? (a) : (b))
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
PairAWPMDCut::PairAWPMDCut(LAMMPS *lmp) : Pair(lmp)
|
||||
{
|
||||
single_enable = 0;
|
||||
|
||||
nmax = 0;
|
||||
min_var = NULL;
|
||||
min_varforce = NULL;
|
||||
nextra = 4;
|
||||
pvector = new double[nextra];
|
||||
|
||||
ermscale=1.;
|
||||
width_pbc=0.;
|
||||
wpmd= new AWPMD_split();
|
||||
|
||||
half_box_length=0;
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
PairAWPMDCut::~PairAWPMDCut()
|
||||
{
|
||||
delete [] pvector;
|
||||
memory->destroy(min_var);
|
||||
memory->destroy(min_varforce);
|
||||
|
||||
if (allocated) {
|
||||
memory->destroy(setflag);
|
||||
memory->destroy(cutsq);
|
||||
memory->destroy(cut);
|
||||
}
|
||||
|
||||
delete wpmd;
|
||||
}
|
||||
|
||||
|
||||
struct cmp_x{
|
||||
double tol;
|
||||
double **xx;
|
||||
cmp_x(double **xx_=NULL, double tol_=1e-12):xx(xx_),tol(tol_){}
|
||||
bool operator()(const pair<int,int> &left, const pair<int,int> &right) const {
|
||||
if(left.first==right.first){
|
||||
double d=xx[left.second][0]-xx[right.second][0];
|
||||
if(d<-tol)
|
||||
return true;
|
||||
else if(d>tol)
|
||||
return false;
|
||||
d=xx[left.second][1]-xx[right.second][1];
|
||||
if(d<-tol)
|
||||
return true;
|
||||
else if(d>tol)
|
||||
return false;
|
||||
d=xx[left.second][2]-xx[right.second][2];
|
||||
if(d<-tol)
|
||||
return true;
|
||||
else
|
||||
return false;
|
||||
}
|
||||
else
|
||||
return left.first<right.first;
|
||||
}
|
||||
};
|
||||
|
||||
/* ---------------------------------------------------------------------- */
|
||||
|
||||
void PairAWPMDCut::compute(int eflag, int vflag)
|
||||
{
|
||||
|
||||
// pvector = [KE, Pauli, ecoul, radial_restraint]
|
||||
for (int i=0; i<4; i++) pvector[i] = 0.0;
|
||||
|
||||
if (eflag || vflag)
|
||||
ev_setup(eflag,vflag);
|
||||
else
|
||||
evflag = vflag_fdotr = 0; //??
|
||||
|
||||
double **x = atom->x;
|
||||
double **f = atom->f;
|
||||
double *q = atom->q;
|
||||
double *erforce = atom->erforce;
|
||||
double *eradius = atom->eradius;
|
||||
int *spin = atom->spin;
|
||||
int *type = atom->type;
|
||||
int *etag = atom->etag;
|
||||
double **v = atom->v;
|
||||
|
||||
int nlocal = atom->nlocal;
|
||||
int nghost = atom->nghost;
|
||||
int ntot=nlocal+nghost;
|
||||
|
||||
int newton_pair = force->newton_pair;
|
||||
double qqrd2e = force->qqrd2e;
|
||||
|
||||
int inum = list->inum;
|
||||
int *ilist = list->ilist;
|
||||
int *numneigh = list->numneigh;
|
||||
int **firstneigh = list->firstneigh;
|
||||
|
||||
|
||||
|
||||
|
||||
// width pbc
|
||||
if(width_pbc<0)
|
||||
wpmd->Lextra=2*half_box_length;
|
||||
else
|
||||
wpmd->Lextra=width_pbc;
|
||||
|
||||
wpmd->newton_pair=newton_pair;
|
||||
|
||||
|
||||
|
||||
# if 1
|
||||
// mapping of the LAMMPS numbers to the AWPMC numbers
|
||||
vector<int> gmap(ntot,-1);
|
||||
|
||||
for (int ii = 0; ii < inum; ii++) {
|
||||
int i = ilist[ii];
|
||||
// local particles are all there
|
||||
gmap[i]=0;
|
||||
Vector_3 ri=Vector_3(x[i][0],x[i][1],x[i][2]);
|
||||
int itype = type[i];
|
||||
int *jlist = firstneigh[i];
|
||||
int jnum = numneigh[i];
|
||||
for (int jj = 0; jj < jnum; jj++) {
|
||||
int j = jlist[jj];
|
||||
j &= NEIGHMASK;
|
||||
if(j>=nlocal){ // this is a ghost
|
||||
Vector_3 rj=Vector_3(x[j][0],x[j][1],x[j][2]);
|
||||
int jtype = type[j];
|
||||
double rsq=(ri-rj).norm2();
|
||||
if (rsq < cutsq[itype][jtype])
|
||||
gmap[j]=0; //bingo, this ghost is really needed
|
||||
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
# else // old mapping
|
||||
// mapping of the LAMMPS numbers to the AWPMC numbers
|
||||
vector<int> gmap(ntot,-1);
|
||||
// map for filtering the clones out: [tag,image] -> id
|
||||
typedef map< pair<int,int>, int, cmp_x > map_t;
|
||||
cmp_x cmp(x);
|
||||
map_t idmap(cmp);
|
||||
for (int ii = 0; ii < inum; ii++) {
|
||||
int i = ilist[ii];
|
||||
// local particles are all there
|
||||
idmap[make_pair(atom->tag[i],i)]=i;
|
||||
bool i_local= i<nlocal ? true : false;
|
||||
if(i_local)
|
||||
gmap[i]=0;
|
||||
else if(gmap[i]==0) // this is a ghost which already has been tested
|
||||
continue;
|
||||
Vector_3 ri=Vector_3(x[i][0],x[i][1],x[i][2]);
|
||||
int itype = type[i];
|
||||
int *jlist = firstneigh[i];
|
||||
int jnum = numneigh[i];
|
||||
for (int jj = 0; jj < jnum; jj++) {
|
||||
int j = jlist[jj];
|
||||
j &= NEIGHMASK;
|
||||
|
||||
pair<map_t::iterator,bool> res=idmap.insert(make_pair(make_pair(atom->tag[j],j),j));
|
||||
bool have_it=!res.second;
|
||||
if(have_it){ // the clone of this particle is already listed
|
||||
if(res.first->second!=j) // check that was not the very same particle
|
||||
gmap[j]=-1; // filter out
|
||||
continue;
|
||||
}
|
||||
|
||||
bool j_local= j<nlocal ? true : false;
|
||||
if((i_local && !j_local) || (j_local && !i_local)){ // some of them is a ghost
|
||||
Vector_3 rj=Vector_3(x[j][0],x[j][1],x[j][2]);
|
||||
int jtype = type[j];
|
||||
double rsq=(ri-rj).norm2();
|
||||
if (rsq < cutsq[itype][jtype]){
|
||||
if(!i_local){
|
||||
gmap[i]=0; //bingo, this ghost is really needed
|
||||
break; // don't need to continue j loop
|
||||
}
|
||||
else
|
||||
gmap[j]=0; //bingo, this ghost is really needed
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
# endif
|
||||
// prepare the solver object
|
||||
wpmd->reset();
|
||||
|
||||
map<int,vector<int> > etmap;
|
||||
// add particles to the AWPMD solver object
|
||||
for (int i = 0; i < ntot; i++) {
|
||||
//int i = ilist[ii];
|
||||
if(gmap[i]<0) // this particle was filtered out
|
||||
continue;
|
||||
if(spin[i]==0) // this is an ion
|
||||
gmap[i]=wpmd->add_ion(q[i], Vector_3(x[i][0],x[i][1],x[i][2]),i<nlocal ? atom->tag[i] : -atom->tag[i]);
|
||||
else if(spin[i]==1 || spin[i]==-1){ // electron, sort them according to the tag
|
||||
etmap[etag[i]].push_back(i);
|
||||
}
|
||||
else
|
||||
error->all(fmt("Invalid spin value (%d) for particle %d !",spin[i],i));
|
||||
}
|
||||
// ion force vector
|
||||
Vector_3 *fi=NULL;
|
||||
if(wpmd->ni)
|
||||
fi= new Vector_3[wpmd->ni];
|
||||
|
||||
// adding electrons
|
||||
for(map<int,vector<int> >::iterator it=etmap.begin(); it!= etmap.end(); ++it){
|
||||
vector<int> &el=it->second;
|
||||
if(!el.size()) // should not happen
|
||||
continue;
|
||||
int s=spin[el[0]] >0 ? 0 : 1;
|
||||
wpmd->add_electron(s); // starts adding the spits
|
||||
for(size_t k=0;k<el.size();k++){
|
||||
int i=el[k];
|
||||
if(spin[el[0]]!=spin[i])
|
||||
error->all(fmt("WP splits for one electron should have the same spin (at particles %d, %d)!",el[0],i));
|
||||
double m= atom->mass ? atom->mass[type[i]] : force->e_mass;
|
||||
Vector_3 xx=Vector_3(x[i][0],x[i][1],x[i][2]);
|
||||
Vector_3 rv=m*Vector_3(v[i][0],v[i][1],v[i][2]);
|
||||
double pv=ermscale*m*atom->ervel[i];
|
||||
Vector_2 cc=Vector_2(atom->cs[2*i],atom->cs[2*i+1]);
|
||||
gmap[i]=wpmd->add_split(xx,rv,atom->eradius[i],pv,cc,1.,atom->q[i],i<nlocal ? atom->tag[i] : -atom->tag[i]);
|
||||
// resetting for the case constraints were applied
|
||||
v[i][0]=rv[0]/m;
|
||||
v[i][1]=rv[1]/m;
|
||||
v[i][2]=rv[2]/m;
|
||||
atom->ervel[i]=pv/(m*ermscale);
|
||||
}
|
||||
}
|
||||
wpmd->set_pbc(NULL); // not required for LAMMPS
|
||||
wpmd->interaction(0x1|0x4|0x10,fi);
|
||||
|
||||
// get forces from the AWPMD solver object
|
||||
for (int ii = 0; ii < inum; ii++) {
|
||||
int i = ilist[ii];
|
||||
if(gmap[i]<0) // this particle was filtered out
|
||||
continue;
|
||||
if(spin[i]==0){ // this is an ion, copying forces
|
||||
int ion=gmap[i];
|
||||
f[i][0]=fi[ion][0];
|
||||
f[i][0]=fi[ion][1];
|
||||
f[i][0]=fi[ion][2];
|
||||
}
|
||||
else { // electron
|
||||
int iel=gmap[i];
|
||||
int s=spin[i] >0 ? 0 : 1;
|
||||
wpmd->get_wp_force(s,iel,(Vector_3 *)f[i],(Vector_3 *)(atom->vforce+3*i),atom->erforce+i,atom->ervelforce+i,(Vector_2 *)(atom->csforce+2*i));
|
||||
}
|
||||
}
|
||||
|
||||
if(fi)
|
||||
delete [] fi;
|
||||
|
||||
// update LAMMPS energy
|
||||
if (eflag_either) {
|
||||
if (eflag_global){
|
||||
eng_coul+= wpmd->get_energy();
|
||||
// pvector = [KE, Pauli, ecoul, radial_restraint]
|
||||
pvector[0] = wpmd->Ee[0]+wpmd->Ee[1];
|
||||
pvector[2] = wpmd->Eii+wpmd->Eei[0]+wpmd->Eei[1]+wpmd->Eee;
|
||||
pvector[1] = pvector[0] + pvector[2] - wpmd->Edk - wpmd->Edc - wpmd->Eii; // All except diagonal terms
|
||||
pvector[3] = wpmd->Ew;
|
||||
}
|
||||
|
||||
if (eflag_atom) {
|
||||
// transfer per-atom energies here
|
||||
for (int i = 0; i < ntot; i++) {
|
||||
if(gmap[i]<0) // this particle was filtered out
|
||||
continue;
|
||||
if(spin[i]==0){
|
||||
eatom[i]=wpmd->Eiep[gmap[i]]+wpmd->Eiip[gmap[i]];
|
||||
}
|
||||
else {
|
||||
int s=spin[i] >0 ? 0 : 1;
|
||||
eatom[i]=wpmd->Eep[s][gmap[i]]+wpmd->Eeip[s][gmap[i]]+wpmd->Eeep[s][gmap[i]]+wpmd->Ewp[s][gmap[i]];
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
if (vflag_fdotr) {
|
||||
virial_fdotr_compute();
|
||||
if (flexible_pressure_flag)
|
||||
virial_eradius_compute();
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
electron width-specific contribution to global virial
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void PairAWPMDCut::virial_eradius_compute()
|
||||
{
|
||||
double *eradius = atom->eradius;
|
||||
double *erforce = atom->erforce;
|
||||
double e_virial;
|
||||
int *spin = atom->spin;
|
||||
|
||||
// sum over force on all particles including ghosts
|
||||
|
||||
if (neighbor->includegroup == 0) {
|
||||
int nall = atom->nlocal + atom->nghost;
|
||||
for (int i = 0; i < nall; i++) {
|
||||
if (spin[i]) {
|
||||
e_virial = erforce[i]*eradius[i]/3;
|
||||
virial[0] += e_virial;
|
||||
virial[1] += e_virial;
|
||||
virial[2] += e_virial;
|
||||
}
|
||||
}
|
||||
|
||||
// neighbor includegroup flag is set
|
||||
// sum over force on initial nfirst particles and ghosts
|
||||
|
||||
} else {
|
||||
int nall = atom->nfirst;
|
||||
for (int i = 0; i < nall; i++) {
|
||||
if (spin[i]) {
|
||||
e_virial = erforce[i]*eradius[i]/3;
|
||||
virial[0] += e_virial;
|
||||
virial[1] += e_virial;
|
||||
virial[2] += e_virial;
|
||||
}
|
||||
}
|
||||
|
||||
nall = atom->nlocal + atom->nghost;
|
||||
for (int i = atom->nlocal; i < nall; i++) {
|
||||
if (spin[i]) {
|
||||
e_virial = erforce[i]*eradius[i]/3;
|
||||
virial[0] += e_virial;
|
||||
virial[1] += e_virial;
|
||||
virial[2] += e_virial;
|
||||
}
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
allocate all arrays
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void PairAWPMDCut::allocate()
|
||||
{
|
||||
allocated = 1;
|
||||
int n = atom->ntypes;
|
||||
|
||||
memory->create(setflag,n+1,n+1,"pair:setflag");
|
||||
for (int i = 1; i <= n; i++)
|
||||
for (int j = i; j <= n; j++)
|
||||
setflag[i][j] = 0;
|
||||
|
||||
memory->create(cutsq,n+1,n+1,"pair:cutsq");
|
||||
memory->create(cut,n+1,n+1,"pair:cut");
|
||||
}
|
||||
|
||||
/* ---------------------------------------------------------------------
|
||||
global settings
|
||||
------------------------------------------------------------------------- */
|
||||
// the format is: pair_style awpmd/cut [<global_cutoff|-1> [command1] [command2] ...]
|
||||
// commands:
|
||||
// [hartree|dproduct|uhf] -- quantum approximation level (default is hartree)
|
||||
// [free|pbc <length|-1>|fix <w0|-1>|relax|harm <w0>] -- width restriction (default is free)
|
||||
// [ermscale <number>] -- scaling factor between electron mass and effective width mass (used for equations of motion only) (default is 1)
|
||||
// [flex_press] -- set flexible pressure flag
|
||||
// -1 for length means default setting (L/2 for cutoff and L for width PBC)
|
||||
void PairAWPMDCut::settings(int narg, char **arg){
|
||||
cut_global=-1.;
|
||||
ermscale=1.;
|
||||
width_pbc=0.;
|
||||
for(int i=0;i<narg;i++){
|
||||
// reading commands
|
||||
// first is cutoff value
|
||||
if(i==0){
|
||||
cut_global = force->numeric(arg[0]);
|
||||
continue;
|
||||
}
|
||||
if(!strcmp(arg[i],"hartree"))
|
||||
wpmd->approx=AWPMD::HARTREE;
|
||||
else if(!strcmp(arg[i],"dproduct"))
|
||||
wpmd->approx=AWPMD::DPRODUCT;
|
||||
else if(!strcmp(arg[i],"uhf"))
|
||||
wpmd->approx=AWPMD::UHF;
|
||||
else if(!strcmp(arg[i],"free"))
|
||||
wpmd->constraint=AWPMD::NONE;
|
||||
else if(!strcmp(arg[i],"fix")){
|
||||
wpmd->constraint=AWPMD::FIX;
|
||||
i++;
|
||||
if(i>=narg)
|
||||
error->all("Setting 'fix' should be followed by a number in awpmd/cut");
|
||||
wpmd->w0=force->numeric(arg[i]);
|
||||
}
|
||||
else if(!strcmp(arg[i],"harm")){
|
||||
wpmd->constraint=AWPMD::HARM;
|
||||
i++;
|
||||
if(i>=narg)
|
||||
error->all("Setting 'harm' should be followed by a number in awpmd/cut");
|
||||
wpmd->w0=force->numeric(arg[i]);
|
||||
wpmd->set_harm_constr(wpmd->w0);
|
||||
}
|
||||
else if(!strcmp(arg[i],"pbc")){
|
||||
i++;
|
||||
if(i>=narg)
|
||||
error->all("Setting 'pbc' should be followed by a number in awpmd/cut");
|
||||
width_pbc=force->numeric(arg[i]);
|
||||
}
|
||||
else if(!strcmp(arg[i],"relax"))
|
||||
wpmd->constraint=AWPMD::RELAX;
|
||||
else if(!strcmp(arg[i],"ermscale")){
|
||||
i++;
|
||||
if(i>=narg)
|
||||
error->all("Setting 'ermscale' should be followed by a number in awpmd/cut");
|
||||
ermscale=force->numeric(arg[i]);
|
||||
}
|
||||
else if(!strcmp(arg[i],"flex_press"))
|
||||
flexible_pressure_flag = 1;
|
||||
}
|
||||
|
||||
|
||||
// reset cutoffs that have been explicitly set
|
||||
/*
|
||||
if (allocated) {
|
||||
int i,j;
|
||||
for (i = 1; i <= atom->ntypes; i++)
|
||||
for (j = i+1; j <= atom->ntypes; j++)
|
||||
if (setflag[i][j]) cut[i][j] = cut_global;
|
||||
}*/
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
set coeffs for one or more type pairs
|
||||
------------------------------------------------------------------------- */
|
||||
// pair settings are as usual
|
||||
void PairAWPMDCut::coeff(int narg, char **arg)
|
||||
{
|
||||
if (narg < 2 || narg > 3) error->all("Incorrect args for pair coefficients");
|
||||
|
||||
/*if(domain->xperiodic == 1 || domain->yperiodic == 1 ||
|
||||
domain->zperiodic == 1) {*/
|
||||
double delx = domain->boxhi[0]-domain->boxlo[0];
|
||||
double dely = domain->boxhi[1]-domain->boxlo[1];
|
||||
double delz = domain->boxhi[2]-domain->boxlo[2];
|
||||
half_box_length = 0.5 * MIN(delx, MIN(dely, delz));
|
||||
//}
|
||||
if(cut_global<0)
|
||||
cut_global=half_box_length;
|
||||
|
||||
if (!allocated)
|
||||
allocate();
|
||||
else{
|
||||
int i,j;
|
||||
for (i = 1; i <= atom->ntypes; i++)
|
||||
for (j = i+1; j <= atom->ntypes; j++)
|
||||
if (setflag[i][j]) cut[i][j] = cut_global;
|
||||
}
|
||||
|
||||
int ilo,ihi,jlo,jhi;
|
||||
force->bounds(arg[0],atom->ntypes,ilo,ihi);
|
||||
force->bounds(arg[1],atom->ntypes,jlo,jhi);
|
||||
|
||||
double cut_one = cut_global;
|
||||
if (narg == 3) cut_one = atof(arg[2]);
|
||||
|
||||
int count = 0;
|
||||
for (int i = ilo; i <= ihi; i++) {
|
||||
for (int j = MAX(jlo,i); j <= jhi; j++) {
|
||||
cut[i][j] = cut_one;
|
||||
setflag[i][j] = 1;
|
||||
count++;
|
||||
}
|
||||
}
|
||||
|
||||
if (count == 0) error->all("Incorrect args for pair coefficients");
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
init specific to this pair style
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void PairAWPMDCut::init_style()
|
||||
{
|
||||
// error and warning checks
|
||||
|
||||
if (!atom->q_flag || !atom->spin_flag ||
|
||||
!atom->eradius_flag || !atom->erforce_flag ) // TO DO: adjust this to match approximation used
|
||||
error->all("Pair awpmd/cut requires atom attributes "
|
||||
"q, spin, eradius, erforce");
|
||||
|
||||
/*
|
||||
if(vflag_atom){ // can't compute virial per atom
|
||||
//warning->
|
||||
error->all("Pair style awpmd can't compute per atom virials");
|
||||
}*/
|
||||
|
||||
// add hook to minimizer for eradius and erforce
|
||||
|
||||
if (update->whichflag == 2)
|
||||
int ignore = update->minimize->request(this,1,0.01);
|
||||
|
||||
// make sure to use the appropriate timestep when using real units
|
||||
|
||||
/*if (update->whichflag == 1) {
|
||||
if (force->qqr2e == 332.06371 && update->dt == 1.0)
|
||||
error->all("You must lower the default real units timestep for pEFF ");
|
||||
}*/
|
||||
|
||||
// need a half neigh list and optionally a granular history neigh list
|
||||
|
||||
//int irequest = neighbor->request(this);
|
||||
|
||||
//if (atom->tag_enable == 0)
|
||||
// error->all("Pair style reax requires atom IDs");
|
||||
|
||||
//if (force->newton_pair == 0)
|
||||
//error->all("Pair style awpmd requires newton pair on");
|
||||
|
||||
//if (strcmp(update->unit_style,"real") != 0 && comm->me == 0)
|
||||
//error->warning("Not using real units with pair reax");
|
||||
|
||||
int irequest = neighbor->request(this);
|
||||
neighbor->requests[irequest]->newton = 2;
|
||||
|
||||
if(force->e_mass==0. || force->hhmrr2e==0. || force->mvh2r==0.)
|
||||
error->all("Pair style awpmd requires e_mass and conversions hhmrr2e, mvh2r to be properly set for unit system");
|
||||
|
||||
wpmd->me=force->e_mass;
|
||||
wpmd->h2_me=force->hhmrr2e/force->e_mass;
|
||||
wpmd->one_h=force->mvh2r;
|
||||
wpmd->coul_pref=force->qqrd2e;
|
||||
|
||||
wpmd->calc_ii=1;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
init for one type pair i,j and corresponding j,i
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
double PairAWPMDCut::init_one(int i, int j)
|
||||
{
|
||||
if (setflag[i][j] == 0)
|
||||
cut[i][j] = mix_distance(cut[i][i],cut[j][j]);
|
||||
|
||||
return cut[i][j];
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
proc 0 writes to restart file
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void PairAWPMDCut::write_restart(FILE *fp)
|
||||
{
|
||||
write_restart_settings(fp);
|
||||
|
||||
int i,j;
|
||||
for (i = 1; i <= atom->ntypes; i++)
|
||||
for (j = i; j <= atom->ntypes; j++) {
|
||||
fwrite(&setflag[i][j],sizeof(int),1,fp);
|
||||
if (setflag[i][j]) fwrite(&cut[i][j],sizeof(double),1,fp);
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
proc 0 reads from restart file, bcasts
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void PairAWPMDCut::read_restart(FILE *fp)
|
||||
{
|
||||
read_restart_settings(fp);
|
||||
allocate();
|
||||
|
||||
int i,j;
|
||||
int me = comm->me;
|
||||
for (i = 1; i <= atom->ntypes; i++)
|
||||
for (j = i; j <= atom->ntypes; j++) {
|
||||
if (me == 0) fread(&setflag[i][j],sizeof(int),1,fp);
|
||||
MPI_Bcast(&setflag[i][j],1,MPI_INT,0,world);
|
||||
if (setflag[i][j]) {
|
||||
if (me == 0) fread(&cut[i][j],sizeof(double),1,fp);
|
||||
MPI_Bcast(&cut[i][j],1,MPI_DOUBLE,0,world);
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
proc 0 writes to restart file
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void PairAWPMDCut::write_restart_settings(FILE *fp)
|
||||
{
|
||||
fwrite(&cut_global,sizeof(double),1,fp);
|
||||
fwrite(&offset_flag,sizeof(int),1,fp);
|
||||
fwrite(&mix_flag,sizeof(int),1,fp);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
proc 0 reads from restart file, bcasts
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void PairAWPMDCut::read_restart_settings(FILE *fp)
|
||||
{
|
||||
if (comm->me == 0) {
|
||||
fread(&cut_global,sizeof(double),1,fp);
|
||||
fread(&offset_flag,sizeof(int),1,fp);
|
||||
fread(&mix_flag,sizeof(int),1,fp);
|
||||
}
|
||||
MPI_Bcast(&cut_global,1,MPI_DOUBLE,0,world);
|
||||
MPI_Bcast(&offset_flag,1,MPI_INT,0,world);
|
||||
MPI_Bcast(&mix_flag,1,MPI_INT,0,world);
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
returns pointers to the log() of electron radius and corresponding force
|
||||
minimizer operates on log(radius) so radius never goes negative
|
||||
these arrays are stored locally by pair style
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void PairAWPMDCut::min_xf_pointers(int ignore, double **xextra, double **fextra)
|
||||
{
|
||||
// grow arrays if necessary
|
||||
// need to be atom->nmax in length
|
||||
int nvar=atom->nmax*(3+1+1+2); // w(1), vel(3), pw(1), cs(2)
|
||||
|
||||
if (nvar > nmax) {
|
||||
memory->destroy(min_var);
|
||||
memory->destroy(min_varforce);
|
||||
nmax = nvar;
|
||||
memory->create(min_var,nmax,"pair:min_var");
|
||||
memory->create(min_varforce,nmax,"pair:min_varforce");
|
||||
}
|
||||
|
||||
*xextra = min_var;
|
||||
*fextra = min_varforce;
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
minimizer requests the log() of electron radius and corresponding force
|
||||
calculate and store in min_eradius and min_erforce
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void PairAWPMDCut::min_xf_get(int ignore)
|
||||
{
|
||||
double *eradius = atom->eradius;
|
||||
double *erforce = atom->erforce;
|
||||
double **v=atom->v;
|
||||
double *vforce=atom->vforce;
|
||||
double *ervel=atom->ervel;
|
||||
double *ervelforce=atom->ervelforce;
|
||||
double *cs=atom->cs;
|
||||
double *csforce=atom->csforce;
|
||||
|
||||
int *spin = atom->spin;
|
||||
int nlocal = atom->nlocal;
|
||||
|
||||
for (int i = 0; i < nlocal; i++)
|
||||
if (spin[i]) {
|
||||
min_var[7*i] = log(eradius[i]);
|
||||
min_varforce[7*i] = eradius[i]*erforce[i];
|
||||
for(int j=0;j<3;j++){
|
||||
min_var[7*i+1+3*j] = v[i][j];
|
||||
min_varforce[7*i+1+3*j] = vforce[3*i+j];
|
||||
}
|
||||
min_var[7*i+4] = ervel[i];
|
||||
min_varforce[7*i+4] = ervelforce[i];
|
||||
min_var[7*i+5] = cs[2*i];
|
||||
min_varforce[7*i+5] = csforce[2*i];
|
||||
min_var[7*i+6] = cs[2*i+1];
|
||||
min_varforce[7*i+6] = csforce[2*i+1];
|
||||
|
||||
} else {
|
||||
for(int j=0;j<7;j++)
|
||||
min_var[7*i+j] = min_varforce[7*i+j] = 0.0;
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
propagate the minimizer values to the atom values
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
void PairAWPMDCut::min_x_set(int ignore)
|
||||
{
|
||||
double *eradius = atom->eradius;
|
||||
double **v=atom->v;
|
||||
double *ervel=atom->ervel;
|
||||
double *cs=atom->cs;
|
||||
|
||||
int *spin = atom->spin;
|
||||
int nlocal = atom->nlocal;
|
||||
|
||||
for (int i = 0; i < nlocal; i++) {
|
||||
if (spin[i]){
|
||||
eradius[i]=exp(min_var[7*i]);
|
||||
for(int j=0;j<3;j++)
|
||||
v[i][j]=min_var[7*i+1+3*j];
|
||||
ervel[i]=min_var[7*i+4];
|
||||
cs[2*i]=min_var[7*i+5];
|
||||
cs[2*i+1]=min_var[7*i+6];
|
||||
}
|
||||
}
|
||||
}
|
||||
|
||||
/* ----------------------------------------------------------------------
|
||||
memory usage of local atom-based arrays
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
double PairAWPMDCut::memory_usage()
|
||||
{
|
||||
double bytes = maxeatom * sizeof(double);
|
||||
bytes += maxvatom*6 * sizeof(double);
|
||||
bytes += 2 * nmax * sizeof(double);
|
||||
return bytes;
|
||||
}
|
||||
81
src/USER-AWPMD/pair_awpmd_cut.h
Normal file
81
src/USER-AWPMD/pair_awpmd_cut.h
Normal file
@ -0,0 +1,81 @@
|
||||
/* ----------------------------------------------------------------------
|
||||
LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
|
||||
http://lammps.sandia.gov, 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.
|
||||
------------------------------------------------------------------------- */
|
||||
/* ----------------------------------------------------------------------
|
||||
Contributing author: Ilya Valuev (JIHT RAS)
|
||||
------------------------------------------------------------------------- */
|
||||
|
||||
|
||||
#ifdef PAIR_CLASS
|
||||
|
||||
PairStyle(awpmd/cut,PairAWPMDCut)
|
||||
|
||||
#else
|
||||
|
||||
#ifndef LMP_PAIR_AWPMD_CUT_H
|
||||
#define LMP_PAIR_AWPMD_CUT_H
|
||||
|
||||
#include "pair.h"
|
||||
|
||||
|
||||
class AWPMD_split;
|
||||
|
||||
|
||||
namespace LAMMPS_NS {
|
||||
|
||||
class PairAWPMDCut : public Pair {
|
||||
friend class FixNVEAwpmd;
|
||||
public:
|
||||
PairAWPMDCut(class LAMMPS *);
|
||||
virtual ~PairAWPMDCut();
|
||||
virtual void compute(int, int);
|
||||
virtual void settings(int, char **);
|
||||
void coeff(int, char **);
|
||||
void init_style();
|
||||
void min_pointers(double **, double **);
|
||||
double init_one(int, int);
|
||||
void write_restart(FILE *);
|
||||
void read_restart(FILE *);
|
||||
virtual void write_restart_settings(FILE *);
|
||||
virtual void read_restart_settings(FILE *);
|
||||
|
||||
void min_xf_pointers(int, double **, double **);
|
||||
void min_xf_get(int);
|
||||
void min_x_set(int);
|
||||
double memory_usage();
|
||||
|
||||
private:
|
||||
|
||||
|
||||
int flexible_pressure_flag;
|
||||
double cut_global;
|
||||
double **cut;
|
||||
|
||||
|
||||
int nmax; // number of additional variables for minimizer
|
||||
double *min_var,*min_varforce; // additional variables for minimizer
|
||||
|
||||
void allocate();
|
||||
|
||||
void virial_eradius_compute();
|
||||
|
||||
|
||||
AWPMD_split *wpmd; // solver oblect
|
||||
double ermscale; // scale of width mass for motion
|
||||
double width_pbc; // setting for width pbc
|
||||
double half_box_length; // calculated by coeff function
|
||||
};
|
||||
|
||||
}
|
||||
|
||||
#endif
|
||||
#endif
|
||||
Reference in New Issue
Block a user