1 // clang-format off
2 /* ----------------------------------------------------------------------
3    LAMMPS - Large-scale Atomic/Molecular Massively Parallel Simulator
4    https://www.lammps.org/, Sandia National Laboratories
5    Steve Plimpton, sjplimp@sandia.gov
6 
7    Copyright (2003) Sandia Corporation.  Under the terms of Contract
8    DE-AC04-94AL85000 with Sandia Corporation, the U.S. Government retains
9    certain rights in this software.  This software is distributed under
10    the GNU General Public License.
11 
12    See the README file in the top-level LAMMPS directory.
13 ------------------------------------------------------------------------- */
14 
15 /* ----------------------------------------------------------------------
16    Contributing author: Norbert Podhorszki (ORNL)
17 ------------------------------------------------------------------------- */
18 
19 #include "dump_atom_adios.h"
20 #include "atom.h"
21 #include "domain.h"
22 #include "error.h"
23 #include "group.h"
24 #include "memory.h"
25 #include "universe.h"
26 #include "update.h"
27 #include <cstring>
28 
29 #include "adios2.h"
30 
31 using namespace LAMMPS_NS;
32 
33 #define MAX_TEXT_HEADER_SIZE 4096
34 #define DUMP_BUF_CHUNK_SIZE 16384
35 #define DUMP_BUF_INCREMENT_SIZE 4096
36 
37 namespace LAMMPS_NS
38 {
39 class DumpAtomADIOSInternal
40 {
41 
42 public:
DumpAtomADIOSInternal()43     DumpAtomADIOSInternal() {};
44     ~DumpAtomADIOSInternal() = default;
45 
46     // name of adios group, referrable in adios2_config.xml
47     const std::string ioName = "atom";
48     adios2::ADIOS *ad = nullptr; // adios object
49     adios2::IO io;     // adios group of variables and attributes in this dump
50     adios2::Engine fh; // adios file/stream handle object
51     // one ADIOS output variable we need to change every step
52     adios2::Variable<double> varAtoms;
53 };
54 } // namespace LAMMPS_NS
55 
56 /* ---------------------------------------------------------------------- */
57 
DumpAtomADIOS(LAMMPS * lmp,int narg,char ** arg)58 DumpAtomADIOS::DumpAtomADIOS(LAMMPS *lmp, int narg, char **arg)
59 : DumpAtom(lmp, narg, arg)
60 {
61     internal = new DumpAtomADIOSInternal();
62     try {
63         internal->ad =
64             new adios2::ADIOS("adios2_config.xml", world, adios2::DebugON);
65     } catch (std::ios_base::failure &e) {
66         char str[256];
67         snprintf(str, sizeof(str), "ADIOS initialization failed with error: %s",
68                  e.what());
69         error->one(FLERR, str);
70     }
71 }
72 
73 /* ---------------------------------------------------------------------- */
74 
~DumpAtomADIOS()75 DumpAtomADIOS::~DumpAtomADIOS()
76 {
77     if (internal->fh) {
78         internal->fh.Close();
79     }
80     delete internal->ad;
81     delete internal;
82 }
83 
84 /* ---------------------------------------------------------------------- */
85 
openfile()86 void DumpAtomADIOS::openfile()
87 {
88     if (multifile) {
89         // if one file per timestep, replace '*' with current timestep
90         char *filestar = strdup(filename);
91         char *filecurrent = new char[strlen(filestar) + 16];
92         char *ptr = strchr(filestar, '*');
93         *ptr = '\0';
94         if (padflag == 0)
95             snprintf(filecurrent, sizeof(filecurrent), "%s" BIGINT_FORMAT "%s",
96                      filestar, update->ntimestep, ptr + 1);
97         else {
98             char bif[8], pad[16];
99             strcpy(bif, BIGINT_FORMAT);
100             snprintf(pad, sizeof(pad), "%%s%%0%d%s%%s", padflag, &bif[1]);
101             snprintf(filecurrent, sizeof(filecurrent), pad, filestar,
102                      update->ntimestep, ptr + 1);
103         }
104         internal->fh =
105             internal->io.Open(filecurrent, adios2::Mode::Write, world);
106         if (!internal->fh) {
107             char str[128];
108             snprintf(str, sizeof(str), "Cannot open dump file %s", filecurrent);
109             error->one(FLERR, str);
110         }
111         free(filestar);
112         delete[] filecurrent;
113     } else {
114         if (!singlefile_opened) {
115             internal->fh =
116                 internal->io.Open(filename, adios2::Mode::Write, world);
117             if (!internal->fh) {
118                 char str[128];
119                 snprintf(str, sizeof(str), "Cannot open dump file %s",
120                          filename);
121                 error->one(FLERR, str);
122             }
123             singlefile_opened = 1;
124         }
125     }
126 }
127 
128 /* ---------------------------------------------------------------------- */
129 
write()130 void DumpAtomADIOS::write()
131 {
132     if (domain->triclinic == 0) {
133         boxxlo = domain->boxlo[0];
134         boxxhi = domain->boxhi[0];
135         boxylo = domain->boxlo[1];
136         boxyhi = domain->boxhi[1];
137         boxzlo = domain->boxlo[2];
138         boxzhi = domain->boxhi[2];
139     } else {
140         boxxlo = domain->boxlo_bound[0];
141         boxxhi = domain->boxhi_bound[0];
142         boxylo = domain->boxlo_bound[1];
143         boxyhi = domain->boxhi_bound[1];
144         boxzlo = domain->boxlo_bound[2];
145         boxzhi = domain->boxhi_bound[2];
146         boxxy = domain->xy;
147         boxxz = domain->xz;
148         boxyz = domain->yz;
149     }
150 
151     // nme = # of dump lines this proc contributes to dump
152 
153     nme = count();
154 
155     // ntotal = total # of atoms in snapshot
156     // atomOffset = sum of # of atoms up to this proc (exclusive prefix sum)
157 
158     bigint bnme = nme;
159     MPI_Allreduce(&bnme, &ntotal, 1, MPI_LMP_BIGINT, MPI_SUM, world);
160 
161     bigint atomOffset; // sum of all atoms on processes 0..me-1
162     MPI_Scan(&bnme, &atomOffset, 1, MPI_LMP_BIGINT, MPI_SUM, world);
163     atomOffset -= nme; // exclusive prefix sum needed
164 
165     // Now we know the global size and the local subset size and offset
166     // of the atoms table
167     size_t nAtomsGlobal = static_cast<size_t>(ntotal);
168     size_t startRow = static_cast<size_t>(atomOffset);
169     size_t nAtomsLocal = static_cast<size_t>(nme);
170     size_t nColumns = static_cast<size_t>(size_one);
171     internal->varAtoms.SetShape({nAtomsGlobal, nColumns});
172     internal->varAtoms.SetSelection({{startRow, 0}, {nAtomsLocal, nColumns}});
173 
174     // insure buf is sized for packing
175     // adios does not limit per-process data size so nme*size_one is not
176     // constrained to int
177     // if sorting on IDs also request ID list from pack()
178     // sort buf as needed
179 
180     if (nme > maxbuf) {
181         maxbuf = nme;
182         memory->destroy(buf);
183         memory->create(buf, (maxbuf * size_one), "dump:buf");
184     }
185     if (sort_flag && sortcol == 0 && nme > maxids) {
186         maxids = nme;
187         memory->destroy(ids);
188         memory->create(ids, maxids, "dump:ids");
189     }
190 
191     if (sort_flag && sortcol == 0)
192         pack(ids);
193     else
194         pack(nullptr);
195     if (sort_flag)
196         sort();
197 
198     openfile();
199     internal->fh.BeginStep();
200     // write info on data as scalars (by me==0)
201     if (me == 0) {
202         internal->fh.Put<uint64_t>("ntimestep", update->ntimestep);
203         internal->fh.Put<int>("nprocs", nprocs);
204 
205         internal->fh.Put<double>("boxxlo", boxxlo);
206         internal->fh.Put<double>("boxxhi", boxxhi);
207         internal->fh.Put<double>("boxylo", boxylo);
208         internal->fh.Put<double>("boxyhi", boxyhi);
209         internal->fh.Put<double>("boxzlo", boxzlo);
210         internal->fh.Put<double>("boxzhi", boxzhi);
211 
212         if (domain->triclinic) {
213             internal->fh.Put<double>("boxxy", boxxy);
214             internal->fh.Put<double>("boxxz", boxxz);
215             internal->fh.Put<double>("boxyz", boxyz);
216         }
217     }
218     // Everyone needs to write scalar variables that are used as dimensions and
219     // offsets of arrays
220     internal->fh.Put<uint64_t>("natoms", ntotal);
221     internal->fh.Put<int>("ncolumns", size_one);
222     internal->fh.Put<uint64_t>("nme", bnme);
223     internal->fh.Put<uint64_t>("offset", atomOffset);
224     // now write the atoms
225     internal->fh.Put<double>(internal->varAtoms, buf);
226     internal->fh.EndStep(); // I/O will happen now...
227 
228     if (multifile) {
229         internal->fh.Close();
230     }
231 }
232 
233 /* ---------------------------------------------------------------------- */
234 
init_style()235 void DumpAtomADIOS::init_style()
236 {
237     if (image_flag == 0)
238         size_one = 5;
239     else
240         size_one = 8;
241 
242     // setup boundary string
243 
244     domain->boundary_string(boundstr);
245 
246     // remove % from filename since ADIOS always writes a global file with
247     // data/metadata
248     int len = strlen(filename);
249     char *ptr = strchr(filename, '%');
250     if (ptr) {
251         *ptr = '\0';
252         char *s = new char[len - 1];
253         snprintf(s, sizeof(s), "%s%s", filename, ptr + 1);
254         strncpy(filename, s, len);
255     }
256 
257     // setup column string
258 
259     std::vector<std::string> columnNames;
260 
261     if (scale_flag == 0 && image_flag == 0) {
262         columns = (char *)"id type x y z";
263         columnNames = {"id", "type", "x", "y", "z"};
264     } else if (scale_flag == 0 && image_flag == 1) {
265         columns = (char *)"id type x y z ix iy iz";
266         columnNames = {"id", "type", "x", "y", "z", "ix", "iy", "iz"};
267     } else if (scale_flag == 1 && image_flag == 0) {
268         columns = (char *)"id type xs ys zs";
269         columnNames = {"id", "type", "xs", "ys", "zs"};
270     } else if (scale_flag == 1 && image_flag == 1) {
271         columns = (char *)"id type xs ys zs ix iy iz";
272         columnNames = {"id", "type", "xs", "ys", "zs", "ix", "iy", "iz"};
273     }
274 
275     // setup function ptrs
276 
277     if (scale_flag == 1 && image_flag == 0 && domain->triclinic == 0)
278         pack_choice = &DumpAtomADIOS::pack_scale_noimage;
279     else if (scale_flag == 1 && image_flag == 1 && domain->triclinic == 0)
280         pack_choice = &DumpAtomADIOS::pack_scale_image;
281     else if (scale_flag == 1 && image_flag == 0 && domain->triclinic == 1)
282         pack_choice = &DumpAtomADIOS::pack_scale_noimage_triclinic;
283     else if (scale_flag == 1 && image_flag == 1 && domain->triclinic == 1)
284         pack_choice = &DumpAtomADIOS::pack_scale_image_triclinic;
285     else if (scale_flag == 0 && image_flag == 0)
286         pack_choice = &DumpAtomADIOS::pack_noscale_noimage;
287     else if (scale_flag == 0 && image_flag == 1)
288         pack_choice = &DumpAtomADIOS::pack_noscale_image;
289 
290     /* Define the group of variables for the atom style here since it's a fixed
291      * set */
292     internal->io = internal->ad->DeclareIO(internal->ioName);
293     if (!internal->io.InConfigFile()) {
294         // if not defined by user, we can change the default settings
295         // BPFile is the default writer
296         internal->io.SetEngine("BPFile");
297         int num_aggregators = multiproc;
298         if (num_aggregators == 0)
299             num_aggregators = 1;
300         char nstreams[128];
301         snprintf(nstreams, sizeof(nstreams), "%d", num_aggregators);
302         internal->io.SetParameters({{"substreams", nstreams}});
303         if (me == 0 && screen)
304             fprintf(
305                 screen,
306                 "ADIOS method for %s is n-to-m (aggregation with %s writers)\n",
307                 filename, nstreams);
308     }
309 
310     internal->io.DefineVariable<uint64_t>("ntimestep");
311     internal->io.DefineVariable<uint64_t>("natoms");
312 
313     internal->io.DefineVariable<int>("nprocs");
314     internal->io.DefineVariable<int>("ncolumns");
315 
316     internal->io.DefineVariable<double>("boxxlo");
317     internal->io.DefineVariable<double>("boxxhi");
318     internal->io.DefineVariable<double>("boxylo");
319     internal->io.DefineVariable<double>("boxyhi");
320     internal->io.DefineVariable<double>("boxzlo");
321     internal->io.DefineVariable<double>("boxzhi");
322 
323     internal->io.DefineVariable<double>("boxxy");
324     internal->io.DefineVariable<double>("boxxz");
325     internal->io.DefineVariable<double>("boxyz");
326 
327     internal->io.DefineAttribute<int>("triclinic", domain->triclinic);
328     internal->io.DefineAttribute<int>("scaled", scale_flag);
329     internal->io.DefineAttribute<int>("image", image_flag);
330 
331     int *boundaryptr = reinterpret_cast<int *>(domain->boundary);
332     internal->io.DefineAttribute<int>("boundary", boundaryptr, 6);
333 
334     size_t nColumns = static_cast<size_t>(size_one);
335     internal->io.DefineAttribute<std::string>("columns", columnNames.data(),
336                                               nColumns);
337     internal->io.DefineAttribute<std::string>("columnstr", columns);
338     internal->io.DefineAttribute<std::string>("boundarystr", boundstr);
339     internal->io.DefineAttribute<std::string>("LAMMPS/dump_style", "atom");
340     internal->io.DefineAttribute<std::string>("LAMMPS/version",
341                                               lmp->version);
342     internal->io.DefineAttribute<std::string>("LAMMPS/num_ver",
343                                               std::to_string(lmp->num_ver));
344 
345     internal->io.DefineVariable<uint64_t>(
346         "nme", {adios2::LocalValueDim}); // local dimension variable
347     internal->io.DefineVariable<uint64_t>(
348         "offset", {adios2::LocalValueDim}); // local dimension variable
349 
350     // atom table size is not known at the moment
351     // it will be correctly defined at the moment of write
352     size_t UnknownSizeYet = 1;
353     internal->varAtoms = internal->io.DefineVariable<double>(
354         "atoms", {UnknownSizeYet, nColumns}, {UnknownSizeYet, 0},
355         {UnknownSizeYet, nColumns});
356 }
357