tx258/lammps-gran-kokkos
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

add support for text output and restart output

Browse Source
This commit is contained in:
Steve Plimpton
2021-08-24 13:46:29 -06:00
parent 3d34ee5087
commit 588f821ac8
4 changed files with 163 additions and 42 deletions
Download Patch File Download Diff File Expand all files Collapse all files

69
src/EXTRA-FIX/fix_ttm.cpp
View File

@ -21,6 +21,7 @@
#include <cmath>
#include <cstring>
#include <cstdio>
#include "atom.h"
#include "force.h"
#include "update.h"
@ -52,7 +53,7 @@ FixTTM::FixTTM(LAMMPS *lmp, int narg, char **arg) :
T_electron(nullptr), T_electron_old(nullptr),
net_energy_transfer(nullptr), net_energy_transfer_all(nullptr)
{
if (narg != 14) error->all(FLERR,"Illegal fix ttm command");
if (narg < 13) error->all(FLERR,"Illegal fix ttm command");
vector_flag = 1;
size_vector = 2;
@ -73,9 +74,25 @@ FixTTM::FixTTM(LAMMPS *lmp, int narg, char **arg) :
nygrid = utils::inumeric(FLERR,arg[11],false,lmp);
nzgrid = utils::inumeric(FLERR,arg[12],false,lmp);
int n = strlen(arg[13]) + 1;
infile = new char[n];
strcpy(infile,arg[13]);
infile = outfile = NULL;
int iarg = 13;
while (iarg < narg) {
if (strcmp(arg[iarg],"infile") == 0) {
if (iarg+2 > narg) error->all(FLERR,"Illegal fix ttm command");
int n = strlen(arg[iarg+1]) + 1;
infile = new char[n];
strcpy(infile,arg[iarg+1]);
iarg += 2;
} else if (strcmp(arg[iarg],"outfile") == 0) {
if (iarg+3 > narg) error->all(FLERR,"Illegal fix ttm command");
outevery = utils::inumeric(FLERR,arg[iarg+1],false,lmp);
int n = strlen(arg[iarg+2]) + 1;
outfile = new char[n];
strcpy(outfile,arg[iarg+2]);
iarg += 3;
} else error->all(FLERR,"Illegal fix ttm command");
}
// error check
@ -162,6 +179,10 @@ void FixTTM::post_constructor()
allocate_grid();
// zero electron temperatures (default)
memset(&T_electron[0][0][0],0,ngridtotal*sizeof(double));
// zero net_energy_transfer
// in case compute_vector accesses it on timestep 0
@ -170,7 +191,7 @@ void FixTTM::post_constructor()
// set initial electron temperatures from user input file
read_electron_temperatures(infile);
if (infile) read_electron_temperatures(infile);
}
/* ---------------------------------------------------------------------- */
@ -413,18 +434,28 @@ void FixTTM::end_of_step()
(T_electron_old[iz][right_ynode][ix] +
T_electron_old[iz][left_ynode][ix] -
2*T_electron_old[iz][iy][ix])/dy/dy +
(T_electron_old[right_znode][iy][iz] +
(T_electron_old[right_znode][iy][ix] +
T_electron_old[left_znode][iy][ix] -
2*T_electron_old[iz][iy][ix])/dz/dz) -
(net_energy_transfer_all[iz][iy][ix])/del_vol);
}
}
// output of grid temperatures to file
if (outfile && (update->ntimestep % outevery == 0)) {
char *newfile = new char[strlen(outfile) + 16];
strcpy(newfile,outfile);
sprintf(newfile,"%s.%ld",outfile,update->ntimestep);
write_electron_temperatures((const char *) newfile);
}
}
/* ----------------------------------------------------------------------
read in initial electron temperatures from a user-specified file
only called by proc 0
only read by proc 0, grid values are Bcast to other procs
------------------------------------------------------------------------- */
void FixTTM::read_electron_temperatures(const char *filename)
@ -476,13 +507,35 @@ void FixTTM::read_electron_temperatures(const char *filename)
for (int iy = 0; iy < nygrid; iy++)
for (int ix = 0; ix < nxgrid; ix++)
if (T_initial_set[iz][iy][ix] == 0)
error->one(FLERR,"Initial temperatures not all set in fix ttm");
error->all(FLERR,"Fix ttm infile did not set all temperatures");
memory->destroy(T_initial_set);
MPI_Bcast(&T_electron[0][0][0],ngridtotal,MPI_DOUBLE,0,world);
}
/* ----------------------------------------------------------------------
write out current electron temperatures to user-specified file
only written by proc 0
------------------------------------------------------------------------- */
void FixTTM::write_electron_temperatures(const char *filename)
{
if (comm->me) return;
FILE *fp = fopen(filename,"w");
if (!fp) error->one(FLERR,"Fix ttm could not open output file");
int ix,iy,iz;
for (iz = 0; iz < nzgrid; iz++)
for (iy = 0; iy < nygrid; iy++)
for (ix = 0; ix < nxgrid; ix++)
fprintf(fp,"%d %d %d %20.16g\n",ix,iy,iz,T_electron[iz][iy][ix]);
fclose(fp);
}
/* ---------------------------------------------------------------------- */
void FixTTM::reset_dt()

12
src/EXTRA-FIX/fix_ttm.h
View File

@ -54,22 +54,24 @@ class FixTTM : public Fix {
int nxgrid, nygrid, nzgrid; // size of global grid
int ngridtotal; // total size of global grid
int deallocate_flag;
int outflag;
int outflag,outevery;
double shift;
double e_energy,transfer_energy;
char *infile,*outfile;
class RanMars *random;
char *infile;
double *gfactor1, *gfactor2, *ratio, **flangevin;
double ***T_electron, ***T_electron_old;
double ***net_energy_transfer, ***net_energy_transfer_all;
double electronic_specific_heat, electronic_density;
double electronic_thermal_conductivity;
double gamma_p, gamma_s, v_0, v_0_sq;
double *gfactor1, *gfactor2, *ratio, **flangevin;
double ***T_electron, ***T_electron_old;
double ***net_energy_transfer, ***net_energy_transfer_all;
virtual void allocate_grid();
virtual void deallocate_grid();
virtual void read_electron_temperatures(const char *);
virtual void write_electron_temperatures(const char *);
};
} // namespace LAMMPS_NS

120
src/EXTRA-FIX/fix_ttm_grid.cpp
View File

@ -38,7 +38,9 @@ using namespace FixConst;
static constexpr int MAXLINE = 256;
static constexpr int CHUNK = 1024;
#define OFFSET 16384
//#define OFFSET 16384
#define OFFSET 0
/* ---------------------------------------------------------------------- */
@ -65,6 +67,10 @@ void FixTTMGrid::post_constructor()
allocate_grid();
// zero electron temperatures (default)
memset(&T_electron[nzlo_out][nylo_out][nxlo_out],0,ngridout*sizeof(double));
// zero net_energy_transfer
// in case compute_vector accesses it on timestep 0
@ -74,7 +80,7 @@ void FixTTMGrid::post_constructor()
// set initial electron temperatures from user input file
read_electron_temperatures(infile);
if (infile) read_electron_temperatures(infile);
}
/* ---------------------------------------------------------------------- */
@ -247,6 +253,16 @@ void FixTTMGrid::end_of_step()
gc->forward_comm(GridComm::FIX,this,1,sizeof(double),0,
gc_buf1,gc_buf2,MPI_DOUBLE);
// output of grid temperatures to file
if (outfile && (update->ntimestep % outevery == 0)) {
char *newfile = new char[strlen(outfile) + 16];
strcpy(newfile,outfile);
sprintf(newfile,"%s.%ld",outfile,update->ntimestep);
write_electron_temperatures((const char *) newfile);
}
}
/* ----------------------------------------------------------------------
@ -261,9 +277,10 @@ void FixTTMGrid::read_electron_temperatures(const char *filename)
int me = comm->me;
// initialize my own+ghost grid values to zero
memset(&T_electron[nzlo_out][nylo_out][nxlo_out],0,ngridout*sizeof(double));
int ***T_initial_set;
memory->create3d_offset(T_initial_set,nzlo_in,nzhi_in,nylo_in,nyhi_in,
nxlo_in,nxhi_in,"ttm/grid:T_initial_set");
memset(&T_initial_set[nzlo_in][nylo_in][nxlo_in],0,ngridmine*sizeof(int));
// proc 0 opens file
@ -320,9 +337,11 @@ void FixTTMGrid::read_electron_temperatures(const char *filename)
if (ix >= nxlo_in && ix <= nxhi_in &&
iy >= nylo_in && iy <= nyhi_in &&
iz >= nzlo_in && iz <= nzhi_in)
iz >= nzlo_in && iz <= nzhi_in) {
T_electron[iz][iy][ix] = utils::numeric(FLERR,values[3],true,lmp);
T_initial_set[iz][iy][ix] = 1;
}
buf = next + 1;
}
@ -338,20 +357,24 @@ void FixTTMGrid::read_electron_temperatures(const char *filename)
delete [] values;
delete [] buffer;
// check that all owned temperature values are > 0.0
// check completeness of input data
int flag = 0;
for (iz = nzlo_in; iz <= nzhi_in; iz++)
for (iy = nylo_in; iy <= nyhi_in; iy++)
for (ix = nxlo_in; ix <= nxhi_in; ix++)
if (T_electron[iz][iy][ix] <= 0.0) flag = 1;
if (T_initial_set[iz][iy][ix] == 0) {
flag = 1;
printf("UNSET %d %d %d\n",ix,iy,iz);
}
int flagall;
MPI_Allreduce(&flag,&flagall,1,MPI_INT,MPI_SUM,world);
if (flagall)
error->all(FLERR,
"Fix ttm infile did not set all temperatures or some are <= 0.0");
error->all(FLERR,"Fix ttm/grid infile did not set all temperatures");
memory->destroy3d_offset(T_initial_set,nzlo_in,nylo_in,nxlo_in);
// communicate new T_electron values to ghost grid points
@ -359,6 +382,23 @@ void FixTTMGrid::read_electron_temperatures(const char *filename)
gc_buf1,gc_buf2,MPI_DOUBLE);
}
/* ----------------------------------------------------------------------
write out current electron temperatures to user-specified file
only written by proc 0
------------------------------------------------------------------------- */
void FixTTMGrid::write_electron_temperatures(const char *filename)
{
if (comm->me == 0) {
FPout = fopen(filename,"w");
if (!FPout) error->one(FLERR,"Fix ttm/grid could not open output file");
}
gc->gather(GridComm::FIX,this,1,sizeof(double),1,NULL,MPI_DOUBLE);
if (comm->me == 0) fclose(FPout);
}
/* ----------------------------------------------------------------------
pack own values to buf to send to another proc
------------------------------------------------------------------------- */
@ -417,6 +457,14 @@ void FixTTMGrid::unpack_reverse_grid(int flag, void *vbuf, int nlist, int *list)
void FixTTMGrid::allocate_grid()
{
// partition global grid across procs
// n xyz lo/hi in = lower/upper bounds of global grid this proc owns
// indices range from 0 to N-1 inclusive in each dim
comm->partition_grid(nxgrid,nygrid,nzgrid,0.0,
nxlo_in,nxhi_in,nylo_in,nyhi_in,nzlo_in,nzhi_in);
/*
// global indices of grid range from 0 to N-1
// nlo_in,nhi_in = lower/upper limits of the 3d sub-brick of
// global grid that I own without ghost cells
@ -442,6 +490,7 @@ void FixTTMGrid::allocate_grid()
nzlo_in = static_cast<int> (comm->mysplit[1][0] * nzgrid);
nzhi_in = static_cast<int> (comm->mysplit[1][1] * nzgrid) - 1;
}
*/
// nlo,nhi = min/max index of global grid pt my owned atoms can be mapped to
// finite difference stencil requires extra grid pt around my owned grid pts
@ -480,6 +529,10 @@ void FixTTMGrid::allocate_grid()
error->one(FLERR,"Too many owned+ghost grid points in fix ttm");
ngridout = totalmine;
totalmine = (bigint) (nxhi_in-nxlo_in+1) * (nyhi_in-nylo_in+1) *
(nzhi_in-nzlo_in+1);
ngridmine = totalmine;
gc = new GridComm(lmp,world,nxgrid,nygrid,nzgrid,
nxlo_in,nxhi_in,nylo_in,nyhi_in,nzlo_in,nzhi_in,
nxlo_out,nxhi_out,nylo_out,nyhi_out,nzlo_out,nzhi_out);
@ -520,8 +573,9 @@ void FixTTMGrid::write_restart(FILE *fp)
{
int ix,iy,iz;
int rsize = nxgrid*nygrid*nzgrid + 4;
double *rlist;
memory->create(rlist,nxgrid*nygrid*nzgrid+4,"ttm/grid:rlist");
memory->create(rlist,rsize,"ttm/grid:rlist");
int n = 0;
rlist[n++] = nxgrid;
@ -533,15 +587,10 @@ void FixTTMGrid::write_restart(FILE *fp)
gc->gather(GridComm::FIX,this,1,sizeof(double),0,&rlist[4],MPI_DOUBLE);
for (iz = nzlo_in; iz <= nzhi_in; iz++)
for (iy = nylo_in; iy <= nyhi_in; iy++)
for (ix = nxlo_in; ix <= nxhi_in; ix++)
rlist[n++] = T_electron[iz][iy][ix];
if (comm->me == 0) {
int size = n * sizeof(double);
int size = rsize * sizeof(double);
fwrite(&size,sizeof(int),1,fp);
fwrite(rlist,sizeof(double),n,fp);
fwrite(rlist,sizeof(double),rsize,fp);
}
memory->destroy(rlist);
@ -588,6 +637,7 @@ void FixTTMGrid::restart(char *buf)
/* ----------------------------------------------------------------------
pack values from local grid into buf
used by which = 0 and 1
------------------------------------------------------------------------- */
void FixTTMGrid::pack_gather_grid(int which, void *vbuf)
@ -604,7 +654,8 @@ void FixTTMGrid::pack_gather_grid(int which, void *vbuf)
}
/* ----------------------------------------------------------------------
unpack values from buf into global gbuf based on their indices
which = 0: unpack values from buf into global gbuf based on their indices
which = 1: print values from buf to FPout file
------------------------------------------------------------------------- */
void FixTTMGrid::unpack_gather_grid(int which, void *vbuf, void *vgbuf,
@ -616,15 +667,28 @@ void FixTTMGrid::unpack_gather_grid(int which, void *vbuf, void *vgbuf,
double *buf = (double *) vbuf;
double *gbuf = (double *) vgbuf;
int iglobal;
int ilocal = 0;
if (which == 0) {
int iglobal;
int ilocal = 0;
for (iz = zlo; iz <= zhi; iz++)
for (iy = ylo; iy <= yhi; iy++)
for (ix = xlo; ix <= xhi; ix++) {
iglobal = nygrid*nxgrid*iz + nxgrid*iy + ix;
gbuf[iglobal] = buf[ilocal++];
}
for (iz = zlo; iz <= zhi; iz++)
for (iy = ylo; iy <= yhi; iy++)
for (ix = xlo; ix <= xhi; ix++) {
iglobal = nygrid*nxgrid*iz + nxgrid*iy + ix;
gbuf[iglobal] = buf[ilocal++];
}
} else if (which == 1) {
int ilocal = 0;
double value;
for (iz = zlo; iz <= zhi; iz++)
for (iy = ylo; iy <= yhi; iy++)
for (ix = xlo; ix <= xhi; ix++) {
value = buf[ilocal++];
fprintf(FPout,"%d %d %d %20.16g\n",ix,iy,iz,value);
}
}
}
/* ----------------------------------------------------------------------

4
src/EXTRA-FIX/fix_ttm_grid.h
View File

@ -48,11 +48,12 @@ class FixTTMGrid : public FixTTM {
double memory_usage();
private:
int ngridout;
int ngridmine,ngridout;
int nxlo_in,nxhi_in,nylo_in,nyhi_in,nzlo_in,nzhi_in;
int nxlo_out,nxhi_out,nylo_out,nyhi_out,nzlo_out,nzhi_out;
double delxinv,delyinv,delzinv;
double skin_original;
FILE *FPout;
class GridComm *gc;
int ngc_buf1,ngc_buf2;
@ -61,6 +62,7 @@ class FixTTMGrid : public FixTTM {
void allocate_grid();
void deallocate_grid();
void read_electron_temperatures(const char *);
void write_electron_temperatures(const char *);
};
} // namespace LAMMPS_NS