From 7918919d303d058a15713e43af0c6659a5a4180e Mon Sep 17 00:00:00 2001 From: Steve Plimpton Date: Mon, 6 Jul 2020 17:30:45 -0600 Subject: [PATCH] enable write_data for atom styles with bonus data --- src/atom_vec.cpp | 10 ++ src/atom_vec.h | 4 + src/atom_vec_body.cpp | 39 ++++++++ src/atom_vec_body.h | 3 + src/atom_vec_ellipsoid.cpp | 46 +++++++++ src/atom_vec_ellipsoid.h | 3 + src/atom_vec_hybrid.cpp | 196 +++++++++++++++++++++++-------------- src/atom_vec_hybrid.h | 4 + src/atom_vec_line.cpp | 54 ++++++++++ src/atom_vec_line.h | 3 + src/atom_vec_tri.cpp | 61 ++++++++++++ src/atom_vec_tri.h | 3 + src/read_data.cpp | 2 +- src/write_data.cpp | 146 ++++++++++++++++----------- src/write_data.h | 1 + 15 files changed, 447 insertions(+), 128 deletions(-) diff --git a/src/atom_vec.cpp b/src/atom_vec.cpp index d3c336769b..1d2c99dd17 100644 --- a/src/atom_vec.cpp +++ b/src/atom_vec.cpp @@ -2304,6 +2304,16 @@ void AtomVec::write_improper(FILE *fp, int n, tagint **buf, int index) } } +/* ---------------------------------------------------------------------- + return size_data_bonus + only AtomVecHybrid overrides this, so it can select which bonus data via flag +------------------------------------------------------------------------- */ + +int AtomVec::size_data_bonus_query(int /*flag*/) +{ + return size_data_bonus; +} + /* ---------------------------------------------------------------------- return # of bytes of allocated memory ------------------------------------------------------------------------- */ diff --git a/src/atom_vec.h b/src/atom_vec.h index 8ccf922c4b..71057848a2 100644 --- a/src/atom_vec.h +++ b/src/atom_vec.h @@ -145,6 +145,10 @@ class AtomVec : protected Pointers { virtual int pack_improper(tagint **); virtual void write_improper(FILE *, int, tagint **, int); + virtual int size_data_bonus_query(int); + virtual int pack_data_bonus(double **, int) {return 0;} + virtual void write_data_bonus(FILE *, int, double **, int) {} + virtual int property_atom(char *) {return -1;} virtual void pack_property_atom(int, double *, int, int) {} diff --git a/src/atom_vec_body.cpp b/src/atom_vec_body.cpp index 89fe1681cb..18312482d1 100644 --- a/src/atom_vec_body.cpp +++ b/src/atom_vec_body.cpp @@ -579,6 +579,45 @@ void AtomVecBody::pack_data_pre(int ilocal) else body[ilocal] = 1; } +/* ---------------------------------------------------------------------- + pack bonus body info for writing to data file + if buf is NULL, just return count of bodies +------------------------------------------------------------------------- */ + +int AtomVecBody::pack_data_bonus(double **buf, int /*flag*/) +{ + int i,j; + + tagint *tag = atom->tag; + int nlocal = atom->nlocal; + + // NOTE: needs to call Body sub-class to fill buffer + + int m = 0; + for (i = 0; i < nlocal; i++) { + if (body[i] < 0) continue; + if (buf) { + buf[m][0] = ubuf(tag[i]).d; + j = body[i]; + } + m++; + } + + return m; +} + +/* ---------------------------------------------------------------------- + write bonus body info to data file +------------------------------------------------------------------------- */ + +void AtomVecBody::write_data_bonus(FILE *fp, int n, double **buf, int /*flag*/) +{ + // NOTE: needs to call Body sub-class to do the write + + for (int i = 0; i < n; i++) { + } +} + /* ---------------------------------------------------------------------- unmodify values packed by AtomVec::pack_data() ------------------------------------------------------------------------- */ diff --git a/src/atom_vec_body.h b/src/atom_vec_body.h index a47cfb3b54..04838074d9 100644 --- a/src/atom_vec_body.h +++ b/src/atom_vec_body.h @@ -63,6 +63,9 @@ class AtomVecBody : public AtomVec { void pack_data_pre(int); void pack_data_post(int); + int pack_data_bonus(double **, int); + void write_data_bonus(FILE *, int, double **, int); + // methods used by other classes to query/set body info double radius_body(int, int, int *, double *); diff --git a/src/atom_vec_ellipsoid.cpp b/src/atom_vec_ellipsoid.cpp index a0cded549f..4190f9ff6e 100644 --- a/src/atom_vec_ellipsoid.cpp +++ b/src/atom_vec_ellipsoid.cpp @@ -25,6 +25,7 @@ #include "memory.h" #include "error.h" #include "utils.h" +#include "fmt/format.h" using namespace LAMMPS_NS; using namespace MathConst; @@ -481,6 +482,51 @@ void AtomVecEllipsoid::pack_data_post(int ilocal) rmass[ilocal] = rmass_one; } +/* ---------------------------------------------------------------------- + pack bonus ellipsoid info for writing to data file + if buf is NULL, just return count of ellipsoids +------------------------------------------------------------------------- */ + +int AtomVecEllipsoid::pack_data_bonus(double **buf, int /*flag*/) +{ + int i,j; + + tagint *tag = atom->tag; + int nlocal = atom->nlocal; + + int m = 0; + for (i = 0; i < nlocal; i++) { + if (ellipsoid[i] < 0) continue; + if (buf) { + buf[m][0] = ubuf(tag[i]).d; + j = ellipsoid[i]; + buf[m][1] = bonus[j].shape[0]; + buf[m][2] = bonus[j].shape[1]; + buf[m][3] = bonus[j].shape[2]; + buf[m][4] = bonus[j].quat[0]; + buf[m][5] = bonus[j].quat[1]; + buf[m][6] = bonus[j].quat[2]; + buf[m][7] = bonus[j].quat[3]; + } + m++; + } + + return m; +} + +/* ---------------------------------------------------------------------- + write bonus ellipsoid info to data file +------------------------------------------------------------------------- */ + +void AtomVecEllipsoid::write_data_bonus(FILE *fp, int n, double **buf, int /*flag*/) +{ + for (int i = 0; i < n; i++) { + fmt::print(fp,"{} {} {} {} {} {} {} {}", + (tagint) ubuf(buf[i][0]).i,buf[i][1],buf[i][2],buf[i][3], + buf[i][4],buf[i][5],buf[i][6],buf[i][7]); + } +} + /* ---------------------------------------------------------------------- set shape values in bonus data for particle I oriented aligned with xyz axes diff --git a/src/atom_vec_ellipsoid.h b/src/atom_vec_ellipsoid.h index 4f63112ff7..91191d5f48 100644 --- a/src/atom_vec_ellipsoid.h +++ b/src/atom_vec_ellipsoid.h @@ -56,6 +56,9 @@ class AtomVecEllipsoid : public AtomVec { void pack_data_pre(int); void pack_data_post(int); + int pack_data_bonus(double **, int); + void write_data_bonus(FILE *, int, double **, int); + // unique to AtomVecEllipsoid void set_shape(int, double, double, double); diff --git a/src/atom_vec_hybrid.cpp b/src/atom_vec_hybrid.cpp index 4149a362e6..2548287f07 100644 --- a/src/atom_vec_hybrid.cpp +++ b/src/atom_vec_hybrid.cpp @@ -23,6 +23,7 @@ using namespace LAMMPS_NS; #define NFIELDSTRINGS 12 // # of field strings +enum{ELLIPSOID,LINE,TRIANGLE,BODY}; // also in WriteData /* ---------------------------------------------------------------------- */ @@ -272,77 +273,6 @@ void AtomVecHybrid::force_clear(int n, size_t nbytes) if (styles[k]->forceclearflag) styles[k]->force_clear(n,nbytes); } -/* ---------------------------------------------------------------------- - modify values for AtomVec::pack_restart() to pack -------------------------------------------------------------------------- */ - -void AtomVecHybrid::pack_restart_pre(int ilocal) -{ - for (int k = 0; k < nstyles; k++) - styles[k]->pack_restart_pre(ilocal); -} - -/* ---------------------------------------------------------------------- - unmodify values packed by AtomVec::pack_restart() -------------------------------------------------------------------------- */ - -void AtomVecHybrid::pack_restart_post(int ilocal) -{ - for (int k = 0; k < nstyles; k++) - styles[k]->pack_restart_post(ilocal); -} - -/* ---------------------------------------------------------------------- - initialize other atom quantities after AtomVec::unpack_restart() -------------------------------------------------------------------------- */ - -void AtomVecHybrid::unpack_restart_init(int ilocal) -{ - for (int k = 0; k < nstyles; k++) - styles[k]->unpack_restart_init(ilocal); -} - -/* ---------------------------------------------------------------------- - initialize non-zero atom quantities -------------------------------------------------------------------------- */ - -void AtomVecHybrid::create_atom_post(int ilocal) -{ - for (int k = 0; k < nstyles; k++) - styles[k]->create_atom_post(ilocal); -} - -/* ---------------------------------------------------------------------- - modify what AtomVec::data_atom() just unpacked - or initialize other atom quantities -------------------------------------------------------------------------- */ - -void AtomVecHybrid::data_atom_post(int ilocal) -{ - for (int k = 0; k < nstyles; k++) - styles[k]->data_atom_post(ilocal); -} - -/* ---------------------------------------------------------------------- - modify values for AtomVec::pack_data() to pack -------------------------------------------------------------------------- */ - -void AtomVecHybrid::pack_data_pre(int ilocal) -{ - for (int k = 0; k < nstyles; k++) - styles[k]->pack_data_pre(ilocal); -} - -/* ---------------------------------------------------------------------- - unmodify values packed by AtomVec::pack_data() -------------------------------------------------------------------------- */ - -void AtomVecHybrid::pack_data_post(int ilocal) -{ - for (int k = 0; k < nstyles; k++) - styles[k]->pack_data_post(ilocal); -} - /* ---------------------------------------------------------------------- */ void AtomVecHybrid::copy_bonus(int i, int j, int delflag) @@ -457,6 +387,130 @@ bigint AtomVecHybrid::memory_usage_bonus() return bytes; } +/* ---------------------------------------------------------------------- + modify values for AtomVec::pack_restart() to pack +------------------------------------------------------------------------- */ + +void AtomVecHybrid::pack_restart_pre(int ilocal) +{ + for (int k = 0; k < nstyles; k++) + styles[k]->pack_restart_pre(ilocal); +} + +/* ---------------------------------------------------------------------- + unmodify values packed by AtomVec::pack_restart() +------------------------------------------------------------------------- */ + +void AtomVecHybrid::pack_restart_post(int ilocal) +{ + for (int k = 0; k < nstyles; k++) + styles[k]->pack_restart_post(ilocal); +} + +/* ---------------------------------------------------------------------- + initialize other atom quantities after AtomVec::unpack_restart() +------------------------------------------------------------------------- */ + +void AtomVecHybrid::unpack_restart_init(int ilocal) +{ + for (int k = 0; k < nstyles; k++) + styles[k]->unpack_restart_init(ilocal); +} + +/* ---------------------------------------------------------------------- + initialize non-zero atom quantities +------------------------------------------------------------------------- */ + +void AtomVecHybrid::create_atom_post(int ilocal) +{ + for (int k = 0; k < nstyles; k++) + styles[k]->create_atom_post(ilocal); +} + +/* ---------------------------------------------------------------------- + modify what AtomVec::data_atom() just unpacked + or initialize other atom quantities +------------------------------------------------------------------------- */ + +void AtomVecHybrid::data_atom_post(int ilocal) +{ + for (int k = 0; k < nstyles; k++) + styles[k]->data_atom_post(ilocal); +} + +/* ---------------------------------------------------------------------- + modify values for AtomVec::pack_data() to pack +------------------------------------------------------------------------- */ + +void AtomVecHybrid::pack_data_pre(int ilocal) +{ + for (int k = 0; k < nstyles; k++) + styles[k]->pack_data_pre(ilocal); +} + +/* ---------------------------------------------------------------------- + unmodify values packed by AtomVec::pack_data() +------------------------------------------------------------------------- */ + +void AtomVecHybrid::pack_data_post(int ilocal) +{ + for (int k = 0; k < nstyles; k++) + styles[k]->pack_data_post(ilocal); +} + +/* ---------------------------------------------------------------------- + return size_data_bonus + match flag to sub-style +------------------------------------------------------------------------- */ + +int AtomVecHybrid::size_data_bonus_query(int flag) +{ + for (int k = 0; k < nstyles; k++) { + if (flag == ELLIPSOID && strcmp(keywords[k],"ellipsoid") == 0) + return styles[k]->size_data_bonus; + if (flag == LINE && strcmp(keywords[k],"line") == 0) + return styles[k]->size_data_bonus; + if (flag == TRIANGLE && strcmp(keywords[k],"tri") == 0) + return styles[k]->size_data_bonus; + // this will not work, body style does not set size_data_bonus + // if (flag == BODY && strcmp(keywords[k],"body") == 0) + // return styles[k]->size_data_bonus; + } +} + +/* ---------------------------------------------------------------------- + pack bonus info for writing to data file + match flag to sub-style +------------------------------------------------------------------------- */ + +int AtomVecHybrid::pack_data_bonus(double **buf, int flag) +{ + for (int k = 0; k < nstyles; k++) { + if (flag == ELLIPSOID && strcmp(keywords[k],"ellipsoid") == 0) + return styles[k]->pack_data_bonus(buf,flag); + if (flag == LINE && strcmp(keywords[k],"line") == 0) + return styles[k]->pack_data_bonus(buf,flag); + if (flag == TRIANGLE && strcmp(keywords[k],"tri") == 0) + return styles[k]->pack_data_bonus(buf,flag); + } +} + +/* ---------------------------------------------------------------------- + write bonus info to data file +------------------------------------------------------------------------- */ + +void AtomVecHybrid::write_data_bonus(FILE *fp, int n, double **buf, int flag) +{ + for (int k = 0; k < nstyles; k++) { + if (flag == ELLIPSOID && strcmp(keywords[k],"ellipsoid") == 0) + styles[k]->write_data_bonus(fp,n,buf,flag); + if (flag == LINE && strcmp(keywords[k],"line") == 0) + styles[k]->write_data_bonus(fp,n,buf,flag); + if (flag == TRIANGLE && strcmp(keywords[k],"tri") == 0) + styles[k]->write_data_bonus(fp,n,buf,flag); + } +} + /* ---------------------------------------------------------------------- assign an index to named atom property and return index returned value encodes which sub-style and index returned by sub-style diff --git a/src/atom_vec_hybrid.h b/src/atom_vec_hybrid.h index 69b0aab7b7..ff4f168214 100644 --- a/src/atom_vec_hybrid.h +++ b/src/atom_vec_hybrid.h @@ -58,6 +58,10 @@ class AtomVecHybrid : public AtomVec { void pack_data_pre(int); void pack_data_post(int); + int size_data_bonus_query(int); + int pack_data_bonus(double **, int); + void write_data_bonus(FILE *, int, double **, int); + int property_atom(char *); void pack_property_atom(int, double *, int, int); diff --git a/src/atom_vec_line.cpp b/src/atom_vec_line.cpp index 88d02161d5..6a368fe2da 100644 --- a/src/atom_vec_line.cpp +++ b/src/atom_vec_line.cpp @@ -22,6 +22,7 @@ #include "memory.h" #include "error.h" #include "utils.h" +#include "fmt/format.h" using namespace LAMMPS_NS; using namespace MathConst; @@ -454,6 +455,59 @@ void AtomVecLine::pack_data_post(int ilocal) rmass[ilocal] = rmass_one; } +/* ---------------------------------------------------------------------- + pack bonus line info for writing to data file + if buf is NULL, just return count of lines +------------------------------------------------------------------------- */ + +int AtomVecLine::pack_data_bonus(double **buf, int /*flag*/) +{ + int i,j; + double length,theta; + double xc,yc,x1,x2,y1,y2; + + double **x = atom->x; + tagint *tag = atom->tag; + int nlocal = atom->nlocal; + + int m = 0; + for (i = 0; i < nlocal; i++) { + if (line[i] < 0) continue; + if (buf) { + buf[m][0] = ubuf(tag[i]).d; + j = line[i]; + length = bonus[j].length; + theta = bonus[j].theta; + xc = x[i][0]; + yc = x[i][1]; + x1 = xc - 0.5*cos(theta)*length; + y1 = yc - 0.5*sin(theta)*length; + x2 = xc + 0.5*cos(theta)*length; + y2 = yc + 0.5*sin(theta)*length; + buf[m][1] = x1; + buf[m][2] = y1; + buf[m][3] = x2; + buf[m][4] = y2; + } + m++; + } + + return m; +} + +/* ---------------------------------------------------------------------- + write bonus line info to data file +------------------------------------------------------------------------- */ + +void AtomVecLine::write_data_bonus(FILE *fp, int n, double **buf, int /*flag*/) +{ + for (int i = 0; i < n; i++) { + fmt::print(fp,"{} {} {} {} {}", + (tagint) ubuf(buf[i][0]).i, + buf[i][1],buf[i][2],buf[i][3],buf[i][4]); + } +} + /* ---------------------------------------------------------------------- set length value in bonus data for particle I oriented along x axis diff --git a/src/atom_vec_line.h b/src/atom_vec_line.h index 7bca58c64b..e5e32c06bf 100644 --- a/src/atom_vec_line.h +++ b/src/atom_vec_line.h @@ -56,6 +56,9 @@ class AtomVecLine : public AtomVec { void pack_data_pre(int); void pack_data_post(int); + int pack_data_bonus(double **, int); + void write_data_bonus(FILE *, int, double **, int); + // unique to AtomVecLine void set_length(int, double); diff --git a/src/atom_vec_tri.cpp b/src/atom_vec_tri.cpp index d5ca2342e8..16183995fe 100644 --- a/src/atom_vec_tri.cpp +++ b/src/atom_vec_tri.cpp @@ -23,6 +23,7 @@ #include "memory.h" #include "error.h" #include "utils.h" +#include "fmt/format.h" using namespace LAMMPS_NS; using namespace MathConst; @@ -685,6 +686,66 @@ void AtomVecTri::pack_data_post(int ilocal) rmass[ilocal] = rmass_one; } +/* ---------------------------------------------------------------------- + pack bonus tri info for writing to data file + if buf is NULL, just return count of lines +------------------------------------------------------------------------- */ + +int AtomVecTri::pack_data_bonus(double **buf, int /*flag*/) +{ + int i,j; + double xc,yc,zc; + double dc1[3],dc2[3],dc3[3]; + double p[3][3]; + + double **x = atom->x; + tagint *tag = atom->tag; + int nlocal = atom->nlocal; + + int m = 0; + for (i = 0; i < nlocal; i++) { + if (tri[i] < 0) continue; + if (buf) { + buf[m][0] = ubuf(tag[i]).d; + j = tri[i]; + MathExtra::quat_to_mat(bonus[j].quat,p); + MathExtra::matvec(p,bonus[j].c1,dc1); + MathExtra::matvec(p,bonus[j].c2,dc2); + MathExtra::matvec(p,bonus[j].c3,dc3); + xc = x[i][0]; + yc = x[i][1]; + zc = x[i][2]; + buf[m][1] = xc + dc1[0]; + buf[m][2] = yc + dc1[1]; + buf[m][3] = zc + dc1[2]; + buf[m][4] = xc + dc2[0]; + buf[m][5] = yc + dc2[1]; + buf[m][6] = zc + dc2[2]; + buf[m][7] = xc + dc3[0]; + buf[m][8] = yc + dc3[1]; + buf[m][9] = zc + dc3[2]; + } + m++; + } + + return m; +} + +/* ---------------------------------------------------------------------- + write bonus tri info to data file +------------------------------------------------------------------------- */ + +void AtomVecTri::write_data_bonus(FILE *fp, int n, double **buf, int /*flag*/) +{ + for (int i = 0; i < n; i++) { + fmt::print(fp,"{} {} {} {} {} {} {} {} {} {} {}", + (tagint) ubuf(buf[i][0]).i, + buf[i][1],buf[i][2],buf[i][3], + buf[i][4],buf[i][5],buf[i][6], + buf[i][7],buf[i][8],buf[i][9]); + } +} + /* ---------------------------------------------------------------------- set equilateral tri of size in bonus data for particle I oriented symmetrically in xy plane diff --git a/src/atom_vec_tri.h b/src/atom_vec_tri.h index 0b16285fe6..cbdbc73290 100644 --- a/src/atom_vec_tri.h +++ b/src/atom_vec_tri.h @@ -58,6 +58,9 @@ class AtomVecTri : public AtomVec { void pack_data_pre(int); void pack_data_post(int); + int pack_data_bonus(double **, int); + void write_data_bonus(FILE *, int, double **, int); + // unique to AtomVecTri void set_equilateral(int, double); diff --git a/src/read_data.cpp b/src/read_data.cpp index 244ce12625..02889c44b6 100644 --- a/src/read_data.cpp +++ b/src/read_data.cpp @@ -80,7 +80,7 @@ ReadData::ReadData(LAMMPS *lmp) : Pointers(lmp) fp = NULL; // customize for new sections - // pointers to atom styles that store extra info + // pointers to atom styles that store bonus info nellipsoids = 0; avec_ellipsoid = (AtomVecEllipsoid *) atom->style_match("ellipsoid"); diff --git a/src/write_data.cpp b/src/write_data.cpp index 5d67d4471c..8567325058 100644 --- a/src/write_data.cpp +++ b/src/write_data.cpp @@ -39,6 +39,7 @@ using namespace LAMMPS_NS; enum{II,IJ}; +enum{ELLIPSOID,LINE,TRIANGLE,BODY}; // also in AtomVecHybrid /* ---------------------------------------------------------------------- */ @@ -173,17 +174,6 @@ void WriteData::write(const std::string &file) MPI_Allreduce(&nimpropers_local,&nimpropers,1,MPI_LMP_BIGINT,MPI_SUM,world); } - // check for bonus data. - if (me == 0) { - if ((atom->nellipsoids > 0) - || (atom->nlines > 0) - || (atom->ntris > 0) - || (atom->nbodies > 0)) - error->warning(FLERR,"System has ellipsoids, lines, triangles, or bodies. " - "Those are not yet supported by write_data. The data file " - "will thus be incomplete."); - } - // open data file if (me == 0) { @@ -201,19 +191,30 @@ void WriteData::write(const std::string &file) if (coeffflag) force_fields(); } - // per atom info - // do not write molecular topology for atom_style template + // per atom info in Atoms and Velocities sections if (natoms) atoms(); if (natoms) velocities(); + + // molecular topology info if defined + // do not write molecular topology for atom_style template + if (atom->molecular == 1) { if (atom->nbonds && nbonds) bonds(); if (atom->nangles && nangles) angles(); if (atom->ndihedrals) dihedrals(); if (atom->nimpropers) impropers(); } + + // bonus info if defined + + if (natoms && atom->ellipsoid_flag) bonus(ELLIPSOID); + if (natoms && atom->line_flag) bonus(LINE); + if (natoms && atom->tri_flag) bonus(TRIANGLE); + if (natoms && atom->body_flag) bonus(BODY); // extra sections managed by fixes + if (fixflag) for (int i = 0; i < modify->nfix; i++) if (modify->fix[i]->wd_section) @@ -252,6 +253,15 @@ void WriteData::header() nimpropers,atom->nimpropertypes); } + // bonus info + + if (atom->ellipsoid_flag) fmt::print(fp,"{} ellipsoids\n",atom->nellipsoids); + if (atom->line_flag) fmt::print(fp,"{} lines\n",atom->nlines); + if (atom->tri_flag) fmt::print(fp,"{} triangles\n",atom->ntris); + if (atom->body_flag) fmt::print(fp,"{} bodies\n",atom->nbodies); + + // fix info + if (fixflag) for (int i = 0; i < modify->nfix; i++) if (modify->fix[i]->wd_header) @@ -322,10 +332,9 @@ void WriteData::force_fields() void WriteData::atoms() { // communication buffer for all my Atom info - // max_size = largest buffer needed by any proc + // maxrow X ncol = largest buffer needed by any proc int ncol = atom->avec->size_data_atom + 3; - int sendrow = atom->nlocal; int maxrow; MPI_Allreduce(&sendrow,&maxrow,1,MPI_INT,MPI_MAX,world); @@ -376,10 +385,9 @@ void WriteData::atoms() void WriteData::velocities() { // communication buffer for all my Atom info - // max_size = largest buffer needed by any proc + // maxrow X ncol = largest buffer needed by any proc int ncol = atom->avec->size_velocity + 1; - int sendrow = atom->nlocal; int maxrow; MPI_Allreduce(&sendrow,&maxrow,1,MPI_INT,MPI_MAX,world); @@ -430,6 +438,7 @@ void WriteData::velocities() void WriteData::bonds() { // communication buffer for all my Bond info + // maxrow X ncol = largest buffer needed by any proc int ncol = 3; int sendrow = static_cast (nbonds_local); @@ -484,6 +493,7 @@ void WriteData::bonds() void WriteData::angles() { // communication buffer for all my Angle info + // maxrow X ncol = largest buffer needed by any proc int ncol = 4; int sendrow = static_cast (nangles_local); @@ -538,27 +548,10 @@ void WriteData::angles() void WriteData::dihedrals() { // communication buffer for all my Dihedral info - // max_size = largest buffer needed by any proc + // maxrow X ncol = largest buffer needed by any proc int ncol = 5; - - tagint *tag = atom->tag; - int *num_dihedral = atom->num_dihedral; - tagint **dihedral_atom2 = atom->dihedral_atom2; - int nlocal = atom->nlocal; - int newton_bond = force->newton_bond; - - int i,j; - int sendrow = 0; - if (newton_bond) { - for (i = 0; i < nlocal; i++) - sendrow += num_dihedral[i]; - } else { - for (i = 0; i < nlocal; i++) - for (j = 0; j < num_dihedral[i]; j++) - if (tag[i] == dihedral_atom2[i][j]) sendrow++; - } - + int sendrow = static_cast (ndihedrals_local); int maxrow; MPI_Allreduce(&sendrow,&maxrow,1,MPI_INT,MPI_MAX,world); @@ -610,29 +603,11 @@ void WriteData::dihedrals() void WriteData::impropers() { // communication buffer for all my Improper info - // max_size = largest buffer needed by any proc + // maxrow X ncol = largest buffer needed by any proc int ncol = 5; - - tagint *tag = atom->tag; - int *num_improper = atom->num_improper; - tagint **improper_atom2 = atom->improper_atom2; - int nlocal = atom->nlocal; - int newton_bond = force->newton_bond; - - int i,j; - int sendrow = 0; - if (newton_bond) { - for (i = 0; i < nlocal; i++) - sendrow += num_improper[i]; - } else { - for (i = 0; i < nlocal; i++) - for (j = 0; j < num_improper[i]; j++) - if (tag[i] == improper_atom2[i][j]) sendrow++; - } - + int sendrow = static_cast (nimpropers_local); int maxrow; - MPI_Allreduce(&sendrow,&maxrow,1,MPI_INT,MPI_MAX,world); tagint **buf; if (me == 0) memory->create(buf,MAX(1,maxrow),ncol,"write_data:buf"); @@ -675,6 +650,64 @@ void WriteData::impropers() memory->destroy(buf); } +/* ---------------------------------------------------------------------- + write out Bonus sections of data file + flag indicates which bonus section it is +------------------------------------------------------------------------- */ + +void WriteData::bonus(int flag) +{ + // communication buffer for all my Bonus info + // maxrow X ncol = largest buffer needed by any proc + + int ncol = atom->avec->size_data_bonus_query(flag); + int sendrow = atom->avec->pack_data_bonus(NULL,flag); + int maxrow; + MPI_Allreduce(&sendrow,&maxrow,1,MPI_INT,MPI_MAX,world); + + double **buf; + if (me == 0) memory->create(buf,MAX(1,maxrow),ncol,"write_data:buf"); + else memory->create(buf,MAX(1,sendrow),ncol,"write_data:buf"); + + // pack my bonus data into buf + + atom->avec->pack_data_bonus(buf,flag); + + // write one chunk of info per proc to file + // proc 0 pings each proc, receives its chunk, writes to file + // all other procs wait for ping, send their chunk to proc 0 + + int tmp,recvrow; + + if (me == 0) { + MPI_Status status; + MPI_Request request; + + if (flag == ELLIPSOID) fprintf(fp,"\nEllipsoids\n\n"); + if (flag == LINE) fprintf(fp,"\nLines\n\n"); + if (flag == TRIANGLE) fprintf(fp,"\nTriangles\n\n"); + if (flag == BODY) fprintf(fp,"\nBodies\n\n"); + + for (int iproc = 0; iproc < nprocs; iproc++) { + if (iproc) { + MPI_Irecv(&buf[0][0],maxrow*ncol,MPI_DOUBLE,iproc,0,world,&request); + MPI_Send(&tmp,0,MPI_INT,iproc,0,world); + MPI_Wait(&request,&status); + MPI_Get_count(&status,MPI_DOUBLE,&recvrow); + recvrow /= ncol; + } else recvrow = sendrow; + + atom->avec->write_data_bonus(fp,recvrow,buf,flag); + } + + } else { + MPI_Recv(&tmp,0,MPI_INT,0,0,world,MPI_STATUS_IGNORE); + MPI_Rsend(&buf[0][0],sendrow*ncol,MPI_DOUBLE,0,0,world); + } + + memory->destroy(buf); +} + /* ---------------------------------------------------------------------- write out Mth section of data file owned by Fix ifix ------------------------------------------------------------------------- */ @@ -682,6 +715,7 @@ void WriteData::impropers() void WriteData::fix(int ifix, int mth) { // communication buffer for Fix info + // maxrow X ncol = largest buffer needed by any proc int sendrow,ncol; modify->fix[ifix]->write_data_section_size(mth,sendrow,ncol); diff --git a/src/write_data.h b/src/write_data.h index 5863b68292..e94e7f5595 100644 --- a/src/write_data.h +++ b/src/write_data.h @@ -51,6 +51,7 @@ class WriteData : protected Pointers { void angles(); void dihedrals(); void impropers(); + void bonus(int); void fix(int, int); };