diff --git a/src/MANYBODY/pair_eam_cd.cpp b/src/MANYBODY/pair_eam_cd.cpp index 46d439b879..dc1ec360c6 100644 --- a/src/MANYBODY/pair_eam_cd.cpp +++ b/src/MANYBODY/pair_eam_cd.cpp @@ -499,10 +499,13 @@ void PairEAMCD::read_h_coeff(char *filename) // Seek to end of file, read last part into a buffer and // then skip over lines in buffer until reaching the end. - platform::fseek(fptr, platform::END_OF_FILE); - platform::fseek(fptr, platform::ftell(fptr) - MAXLINE); + if ( (platform::fseek(fptr, platform::END_OF_FILE) < 0) + || (platform::fseek(fptr, platform::ftell(fptr) - MAXLINE) < 0)) + error->one(FLERR,"Failure to seek to end-of-file for reading h(x) coeffs: {}", + utils::getsyserror()); + char *buf = new char[MAXLINE+1]; - fread(buf, 1, MAXLINE, fptr); + utils::sfread(FLERR, buf, 1, MAXLINE, fptr, filename, error); buf[MAXLINE] = '\0'; // must 0-terminate buffer for string processing Tokenizer lines(buf, "\n"); delete[] buf; diff --git a/src/REAXFF/fix_reaxff_species.cpp b/src/REAXFF/fix_reaxff_species.cpp index 4f44cc7c64..bed7fabb20 100644 --- a/src/REAXFF/fix_reaxff_species.cpp +++ b/src/REAXFF/fix_reaxff_species.cpp @@ -246,8 +246,10 @@ FixReaxFFSpecies::~FixReaxFFSpecies() if (posflag && multipos_opened) fclose(pos); } - modify->delete_compute(fmt::format("SPECATOM_{}",id)); - modify->delete_fix(fmt::format("SPECBOND_{}",id)); + try { + modify->delete_compute(fmt::format("SPECATOM_{}",id)); + modify->delete_fix(fmt::format("SPECBOND_{}",id)); + } catch (std::exception &) {} } /* ---------------------------------------------------------------------- */ diff --git a/src/fix_move.cpp b/src/fix_move.cpp index 72fd9b75d2..f7bc4d3640 100644 --- a/src/fix_move.cpp +++ b/src/fix_move.cpp @@ -522,7 +522,7 @@ void FixMove::initial_integrate(int /*vflag*/) } } - // for wiggle: X = X0 + A sin(w*dt) + // for wiggle: X = X0 + A sin(w*dt) } else if (mstyle == WIGGLE) { double arg = omega_rotate * delta; @@ -578,19 +578,19 @@ void FixMove::initial_integrate(int /*vflag*/) } } - // for rotate by right-hand rule around omega: - // P = point = vector = point of rotation - // R = vector = axis of rotation - // w = omega of rotation (from period) - // X0 = xoriginal = initial coord of atom - // R0 = runit = unit vector for R - // D = X0 - P = vector from P to X0 - // C = (D dot R0) R0 = projection of atom coord onto R line - // A = D - C = vector from R line to X0 - // B = R0 cross A = vector perp to A in plane of rotation - // A,B define plane of circular rotation around R line - // X = P + C + A cos(w*dt) + B sin(w*dt) - // V = w R0 cross (A cos(w*dt) + B sin(w*dt)) + // for rotate by right-hand rule around omega: + // P = point = vector = point of rotation + // R = vector = axis of rotation + // w = omega of rotation (from period) + // X0 = xoriginal = initial coord of atom + // R0 = runit = unit vector for R + // D = X0 - P = vector from P to X0 + // C = (D dot R0) R0 = projection of atom coord onto R line + // A = D - C = vector from R line to X0 + // B = R0 cross A = vector perp to A in plane of rotation + // A,B define plane of circular rotation around R line + // X = P + C + A cos(w*dt) + B sin(w*dt) + // V = w R0 cross (A cos(w*dt) + B sin(w*dt)) } else if (mstyle == ROTATE) { double arg = omega_rotate * delta; @@ -707,10 +707,10 @@ void FixMove::initial_integrate(int /*vflag*/) } } - // for variable: compute x,v from variables - // NOTE: also allow for changes to extra attributes? - // omega, angmom, theta, quat - // only necessary if prescribed motion involves rotation + // for variable: compute x,v from variables + // NOTE: also allow for changes to extra attributes? + // omega, angmom, theta, quat + // only necessary if prescribed motion involves rotation } else if (mstyle == VARIABLE) { @@ -778,21 +778,16 @@ void FixMove::initial_integrate(int /*vflag*/) } else if (vxvarstr) { if (vxvarstyle == EQUAL) v[i][0] = vx; else v[i][0] = velocity[i][0]; - if (rmass) { - x[i][0] += dtv * v[i][0]; - } else { - x[i][0] += dtv * v[i][0]; - } + x[i][0] += dtv * v[i][0]; } else { if (rmass) { dtfm = dtf / rmass[i]; v[i][0] += dtfm * f[i][0]; - x[i][0] += dtv * v[i][0]; } else { dtfm = dtf / mass[type[i]]; v[i][0] += dtfm * f[i][0]; - x[i][0] += dtv * v[i][0]; } + x[i][0] += dtv * v[i][0]; } if (yvarstr && vyvarstr) { @@ -806,21 +801,16 @@ void FixMove::initial_integrate(int /*vflag*/) } else if (vyvarstr) { if (vyvarstyle == EQUAL) v[i][1] = vy; else v[i][1] = velocity[i][1]; - if (rmass) { - x[i][1] += dtv * v[i][1]; - } else { - x[i][1] += dtv * v[i][1]; - } + x[i][1] += dtv * v[i][1]; } else { if (rmass) { dtfm = dtf / rmass[i]; v[i][1] += dtfm * f[i][1]; - x[i][1] += dtv * v[i][1]; } else { dtfm = dtf / mass[type[i]]; v[i][1] += dtfm * f[i][1]; - x[i][1] += dtv * v[i][1]; } + x[i][1] += dtv * v[i][1]; } if (zvarstr && vzvarstr) { @@ -834,21 +824,16 @@ void FixMove::initial_integrate(int /*vflag*/) } else if (vzvarstr) { if (vzvarstyle == EQUAL) v[i][2] = vz; else v[i][2] = velocity[i][2]; - if (rmass) { - x[i][2] += dtv * v[i][2]; - } else { - x[i][2] += dtv * v[i][2]; - } + x[i][2] += dtv * v[i][2]; } else { if (rmass) { dtfm = dtf / rmass[i]; v[i][2] += dtfm * f[i][2]; - x[i][2] += dtv * v[i][2]; } else { dtfm = dtf / mass[type[i]]; v[i][2] += dtfm * f[i][2]; - x[i][2] += dtv * v[i][2]; } + x[i][2] += dtv * v[i][2]; } domain->remap_near(x[i],xold); diff --git a/src/output.cpp b/src/output.cpp index b2b8011842..449077d1fb 100644 --- a/src/output.cpp +++ b/src/output.cpp @@ -41,6 +41,8 @@ using namespace LAMMPS_NS; #define DELTA 1 #define EPSDT 1.0e-6 +enum {SETUP, WRITE, RESET_DT}; + /* ---------------------------------------------------------------------- initialize all output ------------------------------------------------------------------------- */ @@ -232,7 +234,7 @@ void Output::setup(int memflag) // only do this if dump written or dump has not been written yet if (writeflag || last_dump[idump] < 0) - calculate_next_dump(0,idump,ntimestep); + calculate_next_dump(SETUP,idump,ntimestep); // if dump not written now, use addstep_compute_all() // since don't know what computes the dump will invoke @@ -361,7 +363,7 @@ void Output::setup(int memflag) dump[idump]->write(); last_dump[idump] = ntimestep; - calculate_next_dump(1,idump,ntimestep); + calculate_next_dump(WRITE,idump,ntimestep); if (mode_dump[idump] == 0 && (dump[idump]->clearstep || var_dump[idump])) @@ -468,12 +470,12 @@ void Output::setup(int memflag) for timestep mode, set next_dump for simulation time mode, set next_time_dump and next_dump which flag depends on caller - 0 = from setup() at start of run - 1 = from write() during run each time a dump file is written - 2 = from reset_dt() called from fix dt/reset when it changes timestep size + SETUP = from setup() at start of run + WRITE = from write() during run each time a dump file is written + RESET_DT = from reset_dt() called from fix dt/reset when it changes timestep size ------------------------------------------------------------------------- */ - void Output::calculate_next_dump(int which, int idump, bigint ntimestep) +void Output::calculate_next_dump(int which, int idump, bigint ntimestep) { // dump mode is by timestep // just set next_dump @@ -482,13 +484,13 @@ void Output::setup(int memflag) if (every_dump[idump]) { - // which = 0: nextdump = next multiple of every_dump - // which = 1: increment nextdump by every_dump + // which = SETUP: nextdump = next multiple of every_dump + // which = WRITE: increment nextdump by every_dump - if (which == 0) + if (which == SETUP) next_dump[idump] = (ntimestep/every_dump[idump])*every_dump[idump] + every_dump[idump]; - else if (which == 1) + else if (which == WRITE) next_dump[idump] += every_dump[idump]; } else { @@ -510,17 +512,28 @@ void Output::setup(int memflag) if (every_time_dump[idump] > 0.0) { - // which = 0: nexttime = next multiple of every_time_dump - // which = 1: increment nexttime by every_time_dump - // which = 2: no change to previous nexttime (only timestep has changed) + // which = SETUP: nexttime = next multiple of every_time_dump + // which = WRITE: increment nexttime by every_time_dump + // which = RESET_DT: no change to previous nexttime (only timestep has changed) - if (which == 0) + switch (which) { + case SETUP: nexttime = static_cast (tcurrent/every_time_dump[idump]) * every_time_dump[idump] + every_time_dump[idump]; - else if (which == 1) + break; + + case WRITE: nexttime = next_time_dump[idump] + every_time_dump[idump]; - else if (which == 2) + break; + + case RESET_DT: nexttime = next_time_dump[idump]; + break; + + default: + nexttime = 0; + error->all(FLERR,"Unexpected argument to calculate_next_dump"); + } nextdump = ntimestep + static_cast ((nexttime - tcurrent - EPSDT*update->dt) / @@ -541,10 +554,10 @@ void Output::setup(int memflag) } else { - // do not re-evaulate variable for which = 2, leave nexttime as-is + // do not re-evaulate variable for which = RESET_DT, leave nexttime as-is // unless next_time_dump < 0.0, which means variable never yet evaluated - if (which < 2 || next_time_dump[idump] < 0.0) { + if (which < RESET_DT || next_time_dump[idump] < 0.0) { nexttime = input->variable->compute_equal(ivar_dump[idump]); } else nexttime = next_time_dump[idump]; @@ -703,7 +716,7 @@ void Output::reset_dt() // since timestep change affects next step if (next_dump[idump] != ntimestep) - calculate_next_dump(2,idump,update->ntimestep); + calculate_next_dump(RESET_DT,idump,update->ntimestep); if (dump[idump]->clearstep || var_dump[idump]) next_time_dump_any = MIN(next_time_dump_any,next_dump[idump]);