diff --git a/src/output.cpp b/src/output.cpp index 1036ff26e5..3fa0d87cf3 100644 --- a/src/output.cpp +++ b/src/output.cpp @@ -545,11 +545,9 @@ void Output::calculate_next_dump(int which, int idump, bigint ntimestep) // if delta is too small to reach next timestep, use multiple of delta if (nextdump == ntimestep) { - double tnext = update->atime + - (ntimestep+1 - update->atimestep) * update->dt; - int multiple = static_cast - ((tnext - nexttime) / every_time_dump[idump]); - nexttime = nexttime + (multiple+1)*every_time_dump[idump]; + double tnext = update->atime + (ntimestep + 1 - update->atimestep) * update->dt; + int multiple = static_cast((tnext - nexttime) / every_time_dump[idump]); + nexttime = nexttime + (multiple + 1) * every_time_dump[idump]; nextdump = ntimestep + static_cast ((nexttime - tcurrent - EPSDT*update->dt) / update->dt) + 1; } diff --git a/src/rerun.cpp b/src/rerun.cpp index 50088a9c2c..2edf6eba86 100644 --- a/src/rerun.cpp +++ b/src/rerun.cpp @@ -17,17 +17,20 @@ #include "domain.h" #include "error.h" #include "finish.h" +#include "input.h" #include "integrate.h" #include "modify.h" #include "output.h" #include "read_dump.h" #include "timer.h" #include "update.h" +#include "variable.h" #include using namespace LAMMPS_NS; +#define EPSDT 1.0e-6 /* ---------------------------------------------------------------------- */ Rerun::Rerun(LAMMPS *lmp) : Command(lmp) {} @@ -165,6 +168,78 @@ void Rerun::command(int narg, char **arg) modify->init(); update->integrate->setup_minimal(1); modify->end_of_step(); + + // fix up the "next_dump" settings for dumps if the sequence differs from what is read in + + for (int idump = 0; idump < output->ndump; ++idump) { + // dumps triggered by timestep + if (output->mode_dump[idump] == 0) { + // rerun has advanced the timestep faster than the dump expected so we need to catch it up + if (output->next_dump[idump] < ntimestep) { + // equidistant dumps + if (output->every_dump[idump]) { + // if current step compatible with dump frequency, adjust next_dump setting for dump + if (ntimestep % output->every_dump[idump] == 0) output->next_dump[idump] = ntimestep; + } else { + // next dump is determined by variable. + // advance next timestep computation from variable until it is equal or larger + // than the current dump timestep; trigger dump only if equal. + bigint savedstep = update->ntimestep; + update->ntimestep = output->next_dump[idump]; + bigint nextdump; + do { + nextdump = (bigint) input->variable->compute_equal(output->ivar_dump[idump]); + update->ntimestep = nextdump; + } while (nextdump < ntimestep); + output->next_dump[idump] = nextdump; + update->ntimestep = savedstep; + } + } + } else { + // dumps triggered by time + double tcurrent = update->atime + (ntimestep - update->atimestep) * update->dt; + // rerun time has moved beyond expected time for dump so we need to catch it up + if (output->next_time_dump[idump] < tcurrent) { + // equidistant dumps in time + if (output->every_time_dump[idump] > 0.0) { + // trigger dump if current time is within +/- half a timestep of the every interval + double every = output->every_time_dump[idump]; + double rest = fabs(tcurrent/every - round(tcurrent/every)) * every/update->dt; + if (rest < 0.5) { + output->next_dump[idump] = ntimestep; + output->next_time_dump[idump] = tcurrent; + } else { + double nexttime = (floor(tcurrent/every) + 1.0) * every; + output->next_dump[idump] = update->ntimestep + + (bigint) ((nexttime - (update->atime + (update->ntimestep - update->atimestep) * + update->dt) - EPSDT*update->dt) / update->dt) + 1; + output->next_time_dump[idump] = nexttime; + } + } else { + // next dump time is determined by variable. + // advance next time computation from variable until is is equal or larger + // than the current time/timestep + bigint savedstep = update->ntimestep; + update->ntimestep = output->next_dump[idump]; + double nexttime = output->next_time_dump[idump]; + bigint nextstep = output->next_dump[idump]; + while (nextstep < ntimestep) { + nexttime = input->variable->compute_equal(output->ivar_dump[idump]); + nextstep = update->ntimestep + + (bigint) ((nexttime - (update->atime + (update->ntimestep - update->atimestep) * + update->dt) - EPSDT*update->dt) / update->dt) + 1; + update->ntimestep = nextstep; + }; + if (ntimestep > 0) { + output->next_time_dump[idump] = nexttime; + output->next_dump[idump] = nextstep; + } + update->ntimestep = savedstep; + } + } + } + } + output->next_dump_any = ntimestep; if (firstflag) output->setup(); else if (output->next) output->write(ntimestep); @@ -172,7 +247,7 @@ void Rerun::command(int narg, char **arg) firstflag = 0; ntimestep = rd->next(ntimestep,last,nevery,nskip); if (stopflag && ntimestep > stop) - error->all(FLERR,"Read rerun dump file timestep > specified stop"); + error->all(FLERR,"Read rerun dump file timestep {} > specified stop {}", ntimestep, stop); if (ntimestep < 0) break; }