diff --git a/src/dump.cpp b/src/dump.cpp index 639dc2fc07..43fdb775b5 100644 --- a/src/dump.cpp +++ b/src/dump.cpp @@ -918,9 +918,22 @@ void Dump::modify_params(int narg, char **arg) else delay_flag = 0; iarg += 2; - } else if (strcmp(arg[iarg],"header") == 0) { + } else if (strcmp(arg[iarg],"delta") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command"); - header_flag = utils::logical(FLERR,arg[iarg+1],false,lmp); + int idump; + for (idump = 0; idump < output->ndump; idump++) + if (strcmp(id,output->dump[idump]->id) == 0) break; + double delta; + if (strstr(arg[iarg+1],"v_") == arg[iarg+1]) { + delete [] output->var_dump[idump]; + output->var_dump[idump] = utils::strdup(&arg[iarg+1][2]); + delta = 0.0; + } else { + delta = utils::numeric(FLERR,arg[iarg+1],false,lmp); + if (delta <= 0.0) error->all(FLERR,"Illegal dump_modify command"); + } + output->mode_dump[idump] = 1; + output->delta_dump[idump] = delta; iarg += 2; } else if (strcmp(arg[iarg],"every") == 0) { @@ -937,28 +950,10 @@ void Dump::modify_params(int narg, char **arg) n = utils::inumeric(FLERR,arg[iarg+1],false,lmp); if (n <= 0) error->all(FLERR,"Illegal dump_modify command"); } - output->time_dump[idump] = 0; + output->mode_dump[idump] = 0; output->every_dump[idump] = n; iarg += 2; - } else if (strcmp(arg[iarg],"delta") == 0) { - if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command"); - int idump; - for (idump = 0; idump < output->ndump; idump++) - if (strcmp(id,output->dump[idump]->id) == 0) break; - double delta; - if (strstr(arg[iarg+1],"v_") == arg[iarg+1]) { - delete [] output->var_dump[idump]; - output->var_dump[idump] = utils::strdup(&arg[iarg+1][2]); - delta = 0.0; - } else { - delta = utils::numeric(FLERR,arg[iarg+1],false,lmp); - if (delta <= 0.0) error->all(FLERR,"Illegal dump_modify command"); - } - output->time_dump[idump] = 1; - output->delta_dump[idump] = delta; - iarg += 2; - } else if (strcmp(arg[iarg],"fileper") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command"); if (!multiproc) @@ -1026,6 +1021,11 @@ void Dump::modify_params(int narg, char **arg) iarg += n; } + } else if (strcmp(arg[iarg],"header") == 0) { + if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command"); + header_flag = utils::logical(FLERR,arg[iarg+1],false,lmp); + iarg += 2; + } else if (strcmp(arg[iarg],"maxfiles") == 0) { if (iarg+2 > narg) error->all(FLERR,"Illegal dump_modify command"); if (!multifile) diff --git a/src/output.cpp b/src/output.cpp index 0c283efa55..a47a4f8ec1 100644 --- a/src/output.cpp +++ b/src/output.cpp @@ -34,6 +34,7 @@ using namespace LAMMPS_NS; #define DELTA 1 +#define EPSDT 1.0e-6 /* ---------------------------------------------------------------------- initialize all output @@ -59,10 +60,11 @@ Output::Output(LAMMPS *lmp) : Pointers(lmp) ndump = 0; max_dump = 0; - time_dump = nullptr; + mode_dump = nullptr; every_dump = nullptr; delta_dump = nullptr; next_dump = nullptr; + next_time_dump = nullptr; last_dump = nullptr; var_dump = nullptr; ivar_dump = nullptr; @@ -94,10 +96,11 @@ Output::~Output() if (thermo) delete thermo; delete [] var_thermo; - memory->destroy(time_dump); + memory->destroy(mode_dump); memory->destroy(every_dump); memory->destroy(delta_dump); memory->destroy(next_dump); + memory->destroy(next_time_dump); memory->destroy(last_dump); for (int i = 0; i < ndump; i++) delete [] var_dump[i]; memory->sfree(var_dump); @@ -130,8 +133,8 @@ void Output::init() for (int i = 0; i < ndump; i++) dump[i]->init(); for (int i = 0; i < ndump; i++) - if ((time_dump[i] == 0 && every_dump[i] == 0) || - (time_dump[i] == 1 && delta_dump[i] == 0.0)} { + if ((mode_dump[i] == 0 && every_dump[i] == 0) || + (mode_dump[i] == 1 && delta_dump[i] == 0.0)) { ivar_dump[i] = input->variable->find(var_dump[i]); if (ivar_dump[i] < 0) error->all(FLERR,"Variable name for dump every or delta does not exist"); @@ -170,9 +173,9 @@ void Output::setup(int memflag) // current timestep is multiple of every and last dump not >= this step // this is first run after dump created and firstflag is set // note that variable freq will not write unless triggered by firstflag - // set next_dump to multiple of every or variable value + // set next_dump, and also next_time_dump for mode_dump = 1 // set next_dump_any to smallest next_dump - // wrap dumps that invoke computes and variable eval with clear/add + // wrap dumps that invoke computes or do variable eval with clear/add // if dump not written now, use addstep_compute_all() since don't know // what computes the dump write would invoke // if no dumps, set next_dump_any to last+1 so will not influence next @@ -183,9 +186,15 @@ void Output::setup(int memflag) for (int idump = 0; idump < ndump; idump++) { if (dump[idump]->clearstep || var_dump[idump]) modify->clearstep_compute(); + writeflag = 0; - if (every_dump[idump] && ntimestep % every_dump[idump] == 0 && - last_dump[idump] != ntimestep) writeflag = 1; + if (mode_dump[idump] == 0) { + if (every_dump[idump] && (ntimestep % every_dump[idump] == 0) && + last_dump[idump] != ntimestep) writeflag = 1; + } else { + if (delta_dump[idump] >= 0.0 && last_dump[idump] != ntimestep) + writeflag = 1; + } if (last_dump[idump] < 0 && dump[idump]->first_flag == 1) writeflag = 1; if (writeflag) { @@ -193,21 +202,15 @@ void Output::setup(int memflag) last_dump[idump] = ntimestep; } - if (every_dump[idump]) - next_dump[idump] = - (ntimestep/every_dump[idump])*every_dump[idump] + every_dump[idump]; - else { - bigint nextdump = static_cast - (input->variable->compute_equal(ivar_dump[idump])); - if (nextdump <= ntimestep) - error->all(FLERR,"Dump every variable returned a bad timestep"); - next_dump[idump] = nextdump; - } + // set next_dump and next_time_dump, 0 arg for setup() + + calculate_next_dump(0,idump,ntimestep); if (dump[idump]->clearstep || var_dump[idump]) { if (writeflag) modify->addstep_compute(next_dump[idump]); else modify->addstep_compute_all(next_dump[idump]); } + if (idump) next_dump_any = MIN(next_dump_any,next_dump[idump]); else next_dump_any = next_dump[0]; } @@ -288,43 +291,40 @@ void Output::setup(int memflag) void Output::write(bigint ntimestep) { - // next_dump does not force output on last step of run - // wrap dumps that invoke computes or eval of variable with clear/add + // perform dump if its next_dump = current ntimestep + // but not if it was already written on this step + // set next_dump and also next_time_dump for mode_dump = 1 + // set next_dump_any to smallest next_dump + // wrap dumps that invoke computes or do variable eval with clear/add + // if dump not written now, use addstep_compute_all() since don't know + // what computes the dump write would invoke + + int writeflag; if (next_dump_any == ntimestep) { for (int idump = 0; idump < ndump; idump++) { if (next_dump[idump] == ntimestep) { if (dump[idump]->clearstep || var_dump[idump]) modify->clearstep_compute(); - if (last_dump[idump] != ntimestep) { + + writeflag = 0; + if (last_dump[idump] != ntimestep) writeflag = 1; + + if (writeflag) { dump[idump]->write(); last_dump[idump] = ntimestep; } - if (time_dump[idump] == 0) { - if (every_dump[idump]) next_dump[idump] += every_dump[idump]; - else { - bigint nextdump = static_cast - (input->variable->compute_equal(ivar_dump[idump])); - if (nextdump <= ntimestep) - error->all(FLERR,"Dump every variable returned a bad timestep"); - next_dump[idump] = nextdump; - } - } else { - if (delta_dump[idump] > 0.0) { - bigint nextstep = - next_dump[idump] += every_dump[idump]; - } else { - double nexttime = input->variable->compute_equal(ivar_dump[idump]); - if (nexttime <= current_time) // NOTE: what is current time - error->all(FLERR,"Dump delta variable returned a bad time"); - next_dump[idump] = nextdump; - } - } + // set next_dump and next_time_dump, 1 arg for write() + + calculate_next_dump(1,idump,ntimestep); - if (dump[idump]->clearstep || var_dump[idump]) - modify->addstep_compute(next_dump[idump]); + if (dump[idump]->clearstep || var_dump[idump]) { + if (writeflag) modify->addstep_compute(next_dump[idump]); + else modify->addstep_compute_all(next_dump[idump]); + } } + if (idump) next_dump_any = MIN(next_dump_any,next_dump[idump]); else next_dump_any = next_dump[0]; } @@ -332,7 +332,7 @@ void Output::write(bigint ntimestep) // next_restart does not force output on last step of run // for toggle = 0, replace "*" with current timestep in restart filename - // eval of variable may invoke computes so wrap with clear/add + // next restart variable may invoke computes so wrap with clear/add if (next_restart == ntimestep) { if (next_restart_single == ntimestep) { @@ -417,6 +417,88 @@ void Output::write_dump(bigint ntimestep) } } +/* ---------------------------------------------------------------------- + calculate when next dump occurs for Dump instance idump + operates in one of two modes, based on mode_dump flag + for timestep mode, set next_dump + for simulation time mode, set next_time_dump and next_dump +------------------------------------------------------------------------- */ + + void Output::calculate_next_dump(int which, int idump, bigint ntimestep) +{ + // dump mode is by timestep + // just set next_dump + + if (mode_dump[idump] == 0) { + + // for setup, make next_dump a multiple of every_dump + + if (every_dump[idump]) { + if (which == 1) next_dump[idump] += every_dump[idump]; + else + next_dump[idump] = + (ntimestep/every_dump[idump])*every_dump[idump] + every_dump[idump]; + } else { + next_dump[idump] = static_cast + (input->variable->compute_equal(ivar_dump[idump])); + if (next_dump[idump] <= ntimestep) + error->all(FLERR,"Dump every variable returned a bad timestep"); + } + + // dump mode is by simulation time + // set next_time_dump and next_dump + + } else { + bigint nextdump; + double nexttime; + double tcurrent = update->atime + + (ntimestep - update->atimestep) * update->dt; + + // for setup, make nexttime a multiple of delta_dump + + if (delta_dump[idump] > 0.0) { + if (which == 1) nexttime = next_time_dump[idump] + delta_dump[idump]; + else + nexttime = static_cast (tcurrent/delta_dump[idump]) * + delta_dump[idump] + delta_dump[idump]; + + nextdump = ntimestep + + static_cast ((nexttime - tcurrent + EPSDT*update->dt) / + update->dt); + + // 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) / delta_dump[idump]); + nexttime = nexttime + (multiple+1)*delta_dump[idump]; + nextdump = ntimestep + + static_cast ((nexttime - tcurrent + EPSDT*update->dt) / + update->dt); + } + + } else { + nexttime = input->variable->compute_equal(ivar_dump[idump]); + if (nexttime <= tcurrent) + error->all(FLERR,"Dump delta variable returned a bad time"); + nextdump = ntimestep + + static_cast ((nexttime - tcurrent + EPSDT*update->dt) / + update->dt); + if (nextdump <= ntimestep) + error->all(FLERR,"Dump delta variable too small for next timestep"); + } + + next_time_dump[idump] = nexttime; + next_dump[idump] = nextdump; + + //printf("END time %20.16g step %ld ratio %g\n", + // next_time_dump[idump],next_dump[idump], + // next_time_dump[idump]/update->dt/(next_dump[idump]+1)); + } +} + /* ---------------------------------------------------------------------- force restart file(s) to be written called from PRD and TAD @@ -568,10 +650,11 @@ void Output::add_dump(int narg, char **arg) max_dump += DELTA; dump = (Dump **) memory->srealloc(dump,max_dump*sizeof(Dump *),"output:dump"); - memory->grow(time_dump,max_dump,"output:time_dump"); + memory->grow(mode_dump,max_dump,"output:mode_dump"); memory->grow(every_dump,max_dump,"output:every_dump"); memory->grow(delta_dump,max_dump,"output:delta_dump"); memory->grow(next_dump,max_dump,"output:next_dump"); + memory->grow(next_time_dump,max_dump,"output:next_time_dump"); memory->grow(last_dump,max_dump,"output:last_dump"); var_dump = (char **) memory->srealloc(var_dump,max_dump*sizeof(char *),"output:var_dump"); @@ -587,7 +670,7 @@ void Output::add_dump(int narg, char **arg) // initialize per-dump data to suitable default values - time_dump[ndump] = 0; + mode_dump[ndump] = 0; every_dump[ndump] = utils::inumeric(FLERR,arg[3],false,lmp); if (every_dump[ndump] <= 0) error->all(FLERR,"Illegal dump command"); delta_dump[ndump] = 0.0; @@ -646,10 +729,11 @@ void Output::delete_dump(char *id) for (int i = idump+1; i < ndump; i++) { dump[i-1] = dump[i]; - time_dump[i-1] = time_dump[i]; + mode_dump[i-1] = mode_dump[i]; every_dump[i-1] = every_dump[i]; delta_dump[i-1] = delta_dump[i]; next_dump[i-1] = next_dump[i]; + next_time_dump[i-1] = next_time_dump[i]; last_dump[i-1] = last_dump[i]; var_dump[i-1] = var_dump[i]; ivar_dump[i-1] = ivar_dump[i]; diff --git a/src/output.h b/src/output.h index b47ec84dea..ed5bfc162b 100644 --- a/src/output.h +++ b/src/output.h @@ -36,13 +36,14 @@ class Output : protected Pointers { int ndump; // # of Dumps defined int max_dump; // max size of Dump list bigint next_dump_any; // next timestep for any Dump - int *time_dump; // 0/1 if write every N timesteps or Delta in sim time - int *every_dump; // write every this many timesteps, 0 if var - double *delta_dump; // write every this delta sim time, 0.0 if var - bigint *next_dump; // next timestep to do each Dump + int *mode_dump; // 0/1 if write every N timesteps or Delta in sim time + int *every_dump; // dump every this many timesteps, 0 if variable + double *delta_dump; // dump every this delta sim time, 0.0 if variable + bigint *next_dump; // next timestep to perform dump + double *next_time_dump; // next simulation time to perform dump (mode = 1) bigint *last_dump; // last timestep each snapshot was output - char **var_dump; // variable name for dump frequency (steps or sim time) - int *ivar_dump; // variable index for dump frequency + char **var_dump; // variable name for next dump (steps or sim time) + int *ivar_dump; // variable index of var_dump name Dump **dump; // list of defined Dumps int restart_flag; // 1 if any restart files are written @@ -89,6 +90,7 @@ class Output : protected Pointers { private: template static Dump *dump_creator(LAMMPS *, int, char **); + void calculate_next_dump(int, int, bigint); }; } // namespace LAMMPS_NS