remove ecouple variable from Fix

This commit is contained in:
Plimpton
2021-01-21 11:56:33 -07:00
parent 182eb35f1a
commit 6e3b9307a4
8 changed files with 39 additions and 36 deletions

View File

@ -47,7 +47,8 @@ class Fix : protected Pointers {
int energy_atom_flag; // 1 if contributes to peratom eng
int virial_global_flag; // 1 if contributes to global virial
int virial_peratom_flag; // 1 if contributes to peratom virial
int ecouple_flag; // 1 if thermostat which accumulates eng to ecouple
int ecouple_flag; // 1 if thermostat fix outputs cumulative
// reservoir energy via compute_scalar()
int time_integrate; // 1 if performs time integration, 0 if no
int rigid_flag; // 1 if integrates rigid bodies, 0 if not
int no_change_box; // 1 if cannot swap ortho <-> triclinic
@ -102,9 +103,8 @@ class Fix : protected Pointers {
int comm_reverse; // size of reverse communication (0 if none)
int comm_border; // size of border communication (0 if none)
double ecouple; // cumulative energy added to reservoir by thermostatting
double virial[6]; // virial for this timestep
double *eatom,**vatom; // per-atom energy/virial for this timestep
double virial[6]; // virial for this timestep
double *eatom,**vatom; // per-atom energy/virial for this timestep
int centroidstressflag; // centroid stress compared to two-body stress
// CENTROID_SAME = same as two-body stress

View File

@ -53,7 +53,8 @@ enum{CONSTANT,EQUAL,ATOM};
FixLangevin::FixLangevin(LAMMPS *lmp, int narg, char **arg) :
Fix(lmp, narg, arg),
gjfflag(0), gfactor1(nullptr), gfactor2(nullptr), ratio(nullptr), tstr(nullptr),
flangevin(nullptr), tforce(nullptr), franprev(nullptr), lv(nullptr), id_temp(nullptr), random(nullptr)
flangevin(nullptr), tforce(nullptr), franprev(nullptr),
lv(nullptr), id_temp(nullptr), random(nullptr)
{
if (narg < 7) error->all(FLERR,"Illegal fix langevin command");
@ -159,7 +160,7 @@ FixLangevin::FixLangevin(LAMMPS *lmp, int narg, char **arg) :
id_temp = nullptr;
temperature = nullptr;
ecouple = 0.0;
energy = 0.0;
// flangevin is unallocated until first call to setup()
// compute_scalar checks for this and returns 0.0
@ -990,7 +991,7 @@ void FixLangevin::end_of_step()
}
}
ecouple += energy_onestep*update->dt;
energy += energy_onestep*update->dt;
}
/* ---------------------------------------------------------------------- */
@ -1068,7 +1069,7 @@ double FixLangevin::compute_scalar()
if (mask[i] & groupbit)
energy_onestep += flangevin[i][0]*v[i][0] + flangevin[i][1]*v[i][1] +
flangevin[i][2]*v[i][2];
ecouple = 0.5*energy_onestep*update->dt;
energy = 0.5*energy_onestep*update->dt;
} else {
for (int i = 0; i < nlocal; i++)
if (mask[i] & groupbit) {
@ -1079,13 +1080,13 @@ double FixLangevin::compute_scalar()
if (tbiasflag)
temperature->restore_bias(i, lv[i]);
}
ecouple = -0.5*energy_onestep*update->dt;
energy = -0.5*energy_onestep*update->dt;
}
}
// convert midstep energy back to previous fullstep energy
double energy_me = ecouple - 0.5*energy_onestep*update->dt;
double energy_me = energy - 0.5*energy_onestep*update->dt;
double energy_all;
MPI_Allreduce(&energy_me,&energy_all,1,MPI_DOUBLE,MPI_SUM,world);

View File

@ -51,7 +51,8 @@ enum{ISO,ANISO,TRICLINIC};
FixNH::FixNH(LAMMPS *lmp, int narg, char **arg) :
Fix(lmp, narg, arg),
rfix(nullptr), id_dilate(nullptr), irregular(nullptr), id_temp(nullptr), id_press(nullptr),
rfix(nullptr), id_dilate(nullptr), irregular(nullptr),
id_temp(nullptr), id_press(nullptr),
eta(nullptr), eta_dot(nullptr), eta_dotdot(nullptr),
eta_mass(nullptr), etap(nullptr), etap_dot(nullptr), etap_dotdot(nullptr),
etap_mass(nullptr)
@ -784,7 +785,8 @@ void FixNH::setup(int /*vflag*/)
} else {
t0 = temperature->compute_scalar();
if (t0 < EPSILON)
error->all(FLERR, "Current temperature too close to zero, consider using ptemp setting");
error->all(FLERR,"Current temperature too close to zero, "
"consider using ptemp setting");
}
}
t_target = t0;

View File

@ -81,7 +81,7 @@ FixTempBerendsen::FixTempBerendsen(LAMMPS *lmp, int narg, char **arg) :
modify->add_compute(cmd);
tflag = 1;
ecouple = 0;
energy = 0;
}
/* ---------------------------------------------------------------------- */
@ -170,7 +170,7 @@ void FixTempBerendsen::end_of_step()
double lamda = sqrt(1.0 + update->dt/t_period*(t_target/t_current - 1.0));
double efactor = 0.5 * force->boltz * tdof;
ecouple += t_current * (1.0-lamda*lamda) * efactor;
energy += t_current * (1.0-lamda*lamda) * efactor;
double **v = atom->v;
int *mask = atom->mask;
@ -238,7 +238,7 @@ void FixTempBerendsen::reset_target(double t_new)
double FixTempBerendsen::compute_scalar()
{
return ecouple;
return energy;
}
/* ----------------------------------------------------------------------
@ -249,7 +249,7 @@ void FixTempBerendsen::write_restart(FILE *fp)
{
int n = 0;
double list[1];
list[n++] = ecouple;
list[n++] = energy;
if (comm->me == 0) {
int size = n * sizeof(double);
@ -266,7 +266,7 @@ void FixTempBerendsen::restart(char *buf)
{
double *list = (double *) buf;
ecouple = list[0];
energy = list[0];
}
/* ----------------------------------------------------------------------

View File

@ -92,7 +92,7 @@ FixTempCSLD::FixTempCSLD(LAMMPS *lmp, int narg, char **arg) :
vhold = nullptr;
nmax = -1;
ecouple = 0.0;
energy = 0.0;
}
/* ---------------------------------------------------------------------- */
@ -250,7 +250,7 @@ void FixTempCSLD::end_of_step()
// tally the kinetic energy transferred between heat bath and system
t_current = temperature->compute_scalar();
ecouple += ekin_old - t_current * 0.5 * temperature->dof * force->boltz;
energy += ekin_old - t_current * 0.5 * temperature->dof * force->boltz;
}
/* ---------------------------------------------------------------------- */
@ -294,7 +294,7 @@ void FixTempCSLD::reset_target(double t_new)
double FixTempCSLD::compute_scalar()
{
return ecouple;
return energy;
}
/* ----------------------------------------------------------------------
@ -304,11 +304,11 @@ double FixTempCSLD::compute_scalar()
void FixTempCSLD::write_restart(FILE *fp)
{
const int PRNGSIZE = 98+2+3;
int nsize = PRNGSIZE*comm->nprocs+2; // pRNG state per proc + nprocs + ecouple
int nsize = PRNGSIZE*comm->nprocs+2; // pRNG state per proc + nprocs + energy
double *list = nullptr;
if (comm->me == 0) {
list = new double[nsize];
list[0] = ecouple;
list[0] = energy;
list[1] = comm->nprocs;
}
double state[PRNGSIZE];
@ -331,7 +331,7 @@ void FixTempCSLD::restart(char *buf)
{
double *list = (double *) buf;
ecouple = list[0];
energy = list[0];
int nprocs = (int) list[1];
if (nprocs != comm->nprocs) {
if (comm->me == 0)

View File

@ -92,7 +92,7 @@ FixTempCSVR::FixTempCSVR(LAMMPS *lmp, int narg, char **arg) :
tflag = 1;
nmax = -1;
ecouple = 0.0;
energy = 0.0;
}
/* ---------------------------------------------------------------------- */
@ -206,7 +206,7 @@ void FixTempCSVR::end_of_step()
// tally the kinetic energy transferred between heat bath and system
ecouple += ekin_old * (1.0 - lamda*lamda);
energy += ekin_old * (1.0 - lamda*lamda);
}
/* ---------------------------------------------------------------------- */
@ -329,7 +329,7 @@ void FixTempCSVR::reset_target(double t_new)
double FixTempCSVR::compute_scalar()
{
return ecouple;
return energy;
}
/* ----------------------------------------------------------------------
@ -339,11 +339,11 @@ double FixTempCSVR::compute_scalar()
void FixTempCSVR::write_restart(FILE *fp)
{
const int PRNGSIZE = 98+2+3;
int nsize = PRNGSIZE*comm->nprocs+2; // pRNG state per proc + nprocs + ecouple
int nsize = PRNGSIZE*comm->nprocs+2; // pRNG state per proc + nprocs + energy
double *list = nullptr;
if (comm->me == 0) {
list = new double[nsize];
list[0] = ecouple;
list[0] = energy;
list[1] = comm->nprocs;
}
double state[PRNGSIZE];
@ -366,7 +366,7 @@ void FixTempCSVR::restart(char *buf)
{
double *list = (double *) buf;
ecouple = list[0];
energy = list[0];
int nprocs = (int) list[1];
if (nprocs != comm->nprocs) {
if (comm->me == 0)

View File

@ -78,7 +78,7 @@ FixTempRescale::FixTempRescale(LAMMPS *lmp, int narg, char **arg) :
modify->add_compute(cmd);
tflag = 1;
ecouple = 0.0;
energy = 0.0;
}
/* ---------------------------------------------------------------------- */
@ -171,7 +171,7 @@ void FixTempRescale::end_of_step()
int *mask = atom->mask;
int nlocal = atom->nlocal;
ecouple += (t_current-t_target) * efactor;
energy += (t_current-t_target) * efactor;
if (which == NOBIAS) {
for (int i = 0; i < nlocal; i++) {
@ -236,7 +236,7 @@ void FixTempRescale::reset_target(double t_new)
double FixTempRescale::compute_scalar()
{
return ecouple;
return energy;
}
/* ----------------------------------------------------------------------
@ -247,7 +247,7 @@ void FixTempRescale::write_restart(FILE *fp)
{
int n = 0;
double list[1];
list[n++] = ecouple;
list[n++] = energy;
if (comm->me == 0) {
int size = n * sizeof(double);
@ -265,7 +265,7 @@ void FixTempRescale::restart(char *buf)
int n = 0;
double *list = (double *) buf;
ecouple = list[n++];
energy = list[n++];
}
/* ----------------------------------------------------------------------

View File

@ -489,7 +489,7 @@ void Modify::end_of_step()
/* ----------------------------------------------------------------------
coupling energy call, only for relevant fixes
stored by each fix in ecouple variable
each thermostsat fix returns this via compute_scalar()
ecouple = cumulative energy added to reservoir by thermostatting
------------------------------------------------------------------------- */
@ -497,7 +497,7 @@ double Modify::energy_couple()
{
double energy = 0.0;
for (int i = 0; i < n_energy_couple; i++)
energy += fix[list_energy_couple[i]]->ecouple;
energy += fix[list_energy_couple[i]]->compute_scalar();
return energy;
}