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 #include "write_data.h"
16
17 #include "angle.h"
18 #include "atom.h"
19 #include "atom_vec.h"
20 #include "bond.h"
21 #include "comm.h"
22 #include "dihedral.h"
23 #include "domain.h"
24 #include "error.h"
25 #include "fix.h"
26 #include "force.h"
27 #include "improper.h"
28 #include "memory.h"
29 #include "modify.h"
30 #include "output.h"
31 #include "pair.h"
32 #include "thermo.h"
33 #include "update.h"
34
35 #include <cstring>
36
37 using namespace LAMMPS_NS;
38
39 enum{II,IJ};
40 enum{ELLIPSOID,LINE,TRIANGLE,BODY}; // also in AtomVecHybrid
41
42 /* ---------------------------------------------------------------------- */
43
WriteData(LAMMPS * lmp)44 WriteData::WriteData(LAMMPS *lmp) : Command(lmp)
45 {
46 MPI_Comm_rank(world,&me);
47 MPI_Comm_size(world,&nprocs);
48 }
49
50 /* ----------------------------------------------------------------------
51 called as write_data command in input script
52 ------------------------------------------------------------------------- */
53
command(int narg,char ** arg)54 void WriteData::command(int narg, char **arg)
55 {
56 if (domain->box_exist == 0)
57 error->all(FLERR,"Write_data command before simulation box is defined");
58
59 if (narg < 1) error->all(FLERR,"Illegal write_data command");
60
61 // if filename contains a "*", replace with current timestep
62
63 std::string file = arg[0];
64 std::size_t found = file.find('*');
65 if (found != std::string::npos)
66 file.replace(found,1,fmt::format("{}",update->ntimestep));
67
68 // read optional args
69 // noinit is a hidden arg, only used by -r command-line switch
70
71 pairflag = II;
72 coeffflag = 1;
73 fixflag = 1;
74 int noinit = 0;
75
76 int iarg = 1;
77 while (iarg < narg) {
78 if (strcmp(arg[iarg],"pair") == 0) {
79 if (iarg+2 > narg) error->all(FLERR,"Illegal write_data command");
80 if (strcmp(arg[iarg+1],"ii") == 0) pairflag = II;
81 else if (strcmp(arg[iarg+1],"ij") == 0) pairflag = IJ;
82 else error->all(FLERR,"Illegal write_data command");
83 iarg += 2;
84 } else if (strcmp(arg[iarg],"noinit") == 0) {
85 noinit = 1;
86 iarg++;
87 } else if (strcmp(arg[iarg],"nocoeff") == 0) {
88 coeffflag = 0;
89 iarg++;
90 } else if (strcmp(arg[iarg],"nofix") == 0) {
91 fixflag = 0;
92 iarg++;
93 } else error->all(FLERR,"Illegal write_data command");
94 }
95
96 // init entire system since comm->exchange is done
97 // comm::init needs neighbor::init needs pair::init needs kspace::init, etc
98 // exception is when called by -r command-line switch
99 // then write_data immediately follows reading of restart file
100 // assume that read_restart initialized necessary values
101 // if don't make exception:
102 // pair->init() can fail due to various unset values:
103 // e.g. pair hybrid coeffs, dpd ghost-atom velocity setting
104
105 if (noinit == 0) {
106 if (comm->me == 0) utils::logmesg(lmp,"System init for write_data ...\n");
107 lmp->init();
108
109 // move atoms to new processors before writing file
110 // do setup_pre_exchange to force update of per-atom info if needed
111 // enforce PBC in case atoms are outside box
112 // call borders() to rebuild atom map since exchange() destroys map
113
114 modify->setup_pre_exchange();
115 if (domain->triclinic) domain->x2lamda(atom->nlocal);
116 domain->pbc();
117 domain->reset_box();
118 comm->setup();
119 comm->exchange();
120 comm->borders();
121 if (domain->triclinic) domain->lamda2x(atom->nlocal+atom->nghost);
122 }
123
124 write(file);
125 }
126
127 /* ----------------------------------------------------------------------
128 called from command()
129 might later let it be directly called within run/minimize loop
130 ------------------------------------------------------------------------- */
131
write(const std::string & file)132 void WriteData::write(const std::string &file)
133 {
134 // special case where reneighboring is not done in integrator
135 // on timestep data file is written (due to build_once being set)
136 // if box is changing, must be reset, else data file will have
137 // wrong box size and atoms will be lost when data file is read
138 // other calls to pbc and domain and comm are not made,
139 // b/c they only make sense if reneighboring is actually performed
140
141 //if (neighbor->build_once) domain->reset_box();
142
143 // natoms = sum of nlocal = value to write into data file
144 // if unequal and thermo lostflag is "error", don't write data file
145
146 bigint nblocal = atom->nlocal;
147 bigint natoms;
148 MPI_Allreduce(&nblocal,&natoms,1,MPI_LMP_BIGINT,MPI_SUM,world);
149 if (natoms != atom->natoms && output->thermo->lostflag == Thermo::ERROR)
150 error->all(FLERR,"Atom count is inconsistent, cannot write data file");
151
152 // sum up bond,angle,dihedral,improper counts
153 // may be different than atom->nbonds,nangles, etc. if broken/turned-off
154
155 if (atom->molecular == Atom::MOLECULAR && (atom->nbonds || atom->nbondtypes)) {
156 nbonds_local = atom->avec->pack_bond(nullptr);
157 MPI_Allreduce(&nbonds_local,&nbonds,1,MPI_LMP_BIGINT,MPI_SUM,world);
158 }
159 if (atom->molecular == Atom::MOLECULAR && (atom->nangles || atom->nangletypes)) {
160 nangles_local = atom->avec->pack_angle(nullptr);
161 MPI_Allreduce(&nangles_local,&nangles,1,MPI_LMP_BIGINT,MPI_SUM,world);
162 }
163
164 if (atom->molecular == Atom::MOLECULAR && (atom->ndihedrals || atom->ndihedraltypes)) {
165 ndihedrals_local = atom->avec->pack_dihedral(nullptr);
166 MPI_Allreduce(&ndihedrals_local,&ndihedrals,1,MPI_LMP_BIGINT,MPI_SUM,world);
167 }
168
169 if (atom->molecular == Atom::MOLECULAR && (atom->nimpropers || atom->nimpropertypes)) {
170 nimpropers_local = atom->avec->pack_improper(nullptr);
171 MPI_Allreduce(&nimpropers_local,&nimpropers,1,MPI_LMP_BIGINT,MPI_SUM,world);
172 }
173
174 // open data file
175
176 if (me == 0) {
177 fp = fopen(file.c_str(),"w");
178 if (fp == nullptr)
179 error->one(FLERR,"Cannot open data file {}: {}",
180 file, utils::getsyserror());
181 }
182
183 // proc 0 writes header, ntype-length arrays, force fields
184
185 if (me == 0) {
186 header();
187 type_arrays();
188 if (coeffflag) force_fields();
189 }
190
191 // per atom info in Atoms and Velocities sections
192
193 if (natoms) atoms();
194 if (natoms) velocities();
195
196 // molecular topology info if defined
197 // do not write molecular topology for atom_style template
198
199 if (atom->molecular == Atom::MOLECULAR) {
200 if (atom->nbonds && nbonds) bonds();
201 if (atom->nangles && nangles) angles();
202 if (atom->ndihedrals) dihedrals();
203 if (atom->nimpropers) impropers();
204 }
205
206 // bonus info if defined
207
208 if (natoms && atom->ellipsoid_flag) bonus(ELLIPSOID);
209 if (natoms && atom->line_flag) bonus(LINE);
210 if (natoms && atom->tri_flag) bonus(TRIANGLE);
211 if (natoms && atom->body_flag) bonus(BODY);
212
213 // extra sections managed by fixes
214
215 if (fixflag)
216 for (int i = 0; i < modify->nfix; i++)
217 if (modify->fix[i]->wd_section)
218 for (int m = 0; m < modify->fix[i]->wd_section; m++) fix(i,m);
219
220 // close data file
221
222 if (me == 0) fclose(fp);
223 }
224
225 /* ----------------------------------------------------------------------
226 proc 0 writes out data file header
227 ------------------------------------------------------------------------- */
228
header()229 void WriteData::header()
230 {
231 fmt::print(fp,"LAMMPS data file via write_data, version {}, "
232 "timestep = {}\n\n",lmp->version,update->ntimestep);
233
234 fmt::print(fp,"{} atoms\n{} atom types\n",atom->natoms,atom->ntypes);
235
236 // only write out number of types for atom style template
237
238 if (atom->molecular == Atom::MOLECULAR) {
239 if (atom->nbonds || atom->nbondtypes)
240 fmt::print(fp,"{} bonds\n{} bond types\n",
241 nbonds,atom->nbondtypes);
242 if (atom->nangles || atom->nangletypes)
243 fmt::print(fp,"{} angles\n{} angle types\n",
244 nangles,atom->nangletypes);
245 if (atom->ndihedrals || atom->ndihedraltypes)
246 fmt::print(fp,"{} dihedrals\n{} dihedral types\n",
247 ndihedrals,atom->ndihedraltypes);
248 if (atom->nimpropers || atom->nimpropertypes)
249 fmt::print(fp,"{} impropers\n{} improper types\n",
250 nimpropers,atom->nimpropertypes);
251 }
252
253 if (atom->molecular == Atom::TEMPLATE) {
254 if (atom->nbondtypes) fmt::print(fp,"{} bond types\n",atom->nbondtypes);
255 if (atom->nangletypes) fmt::print(fp,"{} angle types\n",atom->nangletypes);
256 if (atom->ndihedraltypes) fmt::print(fp,"{} dihedral types\n",atom->ndihedraltypes);
257 if (atom->nimpropertypes) fmt::print(fp,"{} improper types\n",atom->nimpropertypes);
258 }
259
260 // bonus info
261
262 if (atom->ellipsoid_flag) fmt::print(fp,"{} ellipsoids\n",atom->nellipsoids);
263 if (atom->line_flag) fmt::print(fp,"{} lines\n",atom->nlines);
264 if (atom->tri_flag) fmt::print(fp,"{} triangles\n",atom->ntris);
265 if (atom->body_flag) fmt::print(fp,"{} bodies\n",atom->nbodies);
266
267 // fix info
268
269 if (fixflag)
270 for (int i = 0; i < modify->nfix; i++)
271 if (modify->fix[i]->wd_header)
272 for (int m = 0; m < modify->fix[i]->wd_header; m++)
273 modify->fix[i]->write_data_header(fp,m);
274
275 // box info
276
277 auto box = fmt::format("\n{} {} xlo xhi"
278 "\n{} {} ylo yhi"
279 "\n{} {} zlo zhi\n",
280 domain->boxlo[0],domain->boxhi[0],
281 domain->boxlo[1],domain->boxhi[1],
282 domain->boxlo[2],domain->boxhi[2]);
283 if (domain->triclinic)
284 box += fmt::format("{} {} {} xy xz yz\n",
285 domain->xy,domain->xz,domain->yz);
286 fputs(box.c_str(),fp);
287 }
288
289 /* ----------------------------------------------------------------------
290 proc 0 writes out any type-based arrays that are defined
291 ------------------------------------------------------------------------- */
292
type_arrays()293 void WriteData::type_arrays()
294 {
295 if (atom->mass) {
296 double *mass = atom->mass;
297 fputs("\nMasses\n\n",fp);
298 for (int i = 1; i <= atom->ntypes; i++)
299 fmt::print(fp,"{} {:.16g}\n",i,mass[i]);
300 }
301 }
302
303 /* ----------------------------------------------------------------------
304 proc 0 writes out force field info
305 ------------------------------------------------------------------------- */
306
force_fields()307 void WriteData::force_fields()
308 {
309 if (force->pair && force->pair->writedata) {
310 if (pairflag == II) {
311 if ((comm->me == 0) && (force->pair->mixed_flag == 0))
312 error->warning(FLERR,"Not all mixed pair coeffs generated from mixing. "
313 "Use write_data with 'pair ij' option to store all pair coeffs.");
314 fmt::print(fp,"\nPair Coeffs # {}\n\n", force->pair_style);
315 force->pair->write_data(fp);
316 } else if (pairflag == IJ) {
317 fmt::print(fp,"\nPairIJ Coeffs # {}\n\n", force->pair_style);
318 force->pair->write_data_all(fp);
319 }
320 }
321 if (force->bond && force->bond->writedata && atom->nbondtypes) {
322 fmt::print(fp,"\nBond Coeffs # {}\n\n", force->bond_style);
323 force->bond->write_data(fp);
324 }
325 if (force->angle && force->angle->writedata && atom->nangletypes) {
326 fmt::print(fp,"\nAngle Coeffs # {}\n\n", force->angle_style);
327 force->angle->write_data(fp);
328 }
329 if (force->dihedral && force->dihedral->writedata && atom->ndihedraltypes) {
330 fmt::print(fp,"\nDihedral Coeffs # {}\n\n", force->dihedral_style);
331 force->dihedral->write_data(fp);
332 }
333 if (force->improper && force->improper->writedata && atom->nimpropertypes) {
334 fmt::print(fp,"\nImproper Coeffs # {}\n\n", force->improper_style);
335 force->improper->write_data(fp);
336 }
337 }
338
339 /* ----------------------------------------------------------------------
340 write out Atoms section of data file
341 ------------------------------------------------------------------------- */
342
atoms()343 void WriteData::atoms()
344 {
345 // communication buffer for all my Atom info
346 // maxrow X ncol = largest buffer needed by any proc
347
348 int ncol = atom->avec->size_data_atom + 3;
349 int sendrow = atom->nlocal;
350 int maxrow;
351 MPI_Allreduce(&sendrow,&maxrow,1,MPI_INT,MPI_MAX,world);
352
353 double **buf;
354 if (me == 0) memory->create(buf,MAX(1,maxrow),ncol,"write_data:buf");
355 else memory->create(buf,MAX(1,sendrow),ncol,"write_data:buf");
356
357 // pack my atom data into buf
358
359 atom->avec->pack_data(buf);
360
361 // write one chunk of atoms per proc to file
362 // proc 0 pings each proc, receives its chunk, writes to file
363 // all other procs wait for ping, send their chunk to proc 0
364
365 int tmp,recvrow;
366
367 if (me == 0) {
368 MPI_Status status;
369 MPI_Request request;
370
371 fmt::print(fp,"\nAtoms # {}\n\n",atom->atom_style);
372 for (int iproc = 0; iproc < nprocs; iproc++) {
373 if (iproc) {
374 MPI_Irecv(&buf[0][0],maxrow*ncol,MPI_DOUBLE,iproc,0,world,&request);
375 MPI_Send(&tmp,0,MPI_INT,iproc,0,world);
376 MPI_Wait(&request,&status);
377 MPI_Get_count(&status,MPI_DOUBLE,&recvrow);
378 recvrow /= ncol;
379 } else recvrow = sendrow;
380
381 atom->avec->write_data(fp,recvrow,buf);
382 }
383
384 } else {
385 MPI_Recv(&tmp,0,MPI_INT,0,0,world,MPI_STATUS_IGNORE);
386 MPI_Rsend(&buf[0][0],sendrow*ncol,MPI_DOUBLE,0,0,world);
387 }
388
389 memory->destroy(buf);
390 }
391
392 /* ----------------------------------------------------------------------
393 write out Velocities section of data file
394 ------------------------------------------------------------------------- */
395
velocities()396 void WriteData::velocities()
397 {
398 // communication buffer for all my Atom info
399 // maxrow X ncol = largest buffer needed by any proc
400
401 int ncol = atom->avec->size_velocity + 1;
402 int sendrow = atom->nlocal;
403 int maxrow;
404 MPI_Allreduce(&sendrow,&maxrow,1,MPI_INT,MPI_MAX,world);
405
406 double **buf;
407 if (me == 0) memory->create(buf,MAX(1,maxrow),ncol,"write_data:buf");
408 else memory->create(buf,MAX(1,sendrow),ncol,"write_data:buf");
409
410 // pack my velocity data into buf
411
412 atom->avec->pack_vel(buf);
413
414 // write one chunk of velocities per proc to file
415 // proc 0 pings each proc, receives its chunk, writes to file
416 // all other procs wait for ping, send their chunk to proc 0
417
418 int tmp,recvrow;
419
420 if (me == 0) {
421 MPI_Status status;
422 MPI_Request request;
423
424 fputs("\nVelocities\n\n",fp);
425 for (int iproc = 0; iproc < nprocs; iproc++) {
426 if (iproc) {
427 MPI_Irecv(&buf[0][0],maxrow*ncol,MPI_DOUBLE,iproc,0,world,&request);
428 MPI_Send(&tmp,0,MPI_INT,iproc,0,world);
429 MPI_Wait(&request,&status);
430 MPI_Get_count(&status,MPI_DOUBLE,&recvrow);
431 recvrow /= ncol;
432 } else recvrow = sendrow;
433
434 atom->avec->write_vel(fp,recvrow,buf);
435 }
436
437 } else {
438 MPI_Recv(&tmp,0,MPI_INT,0,0,world,MPI_STATUS_IGNORE);
439 MPI_Rsend(&buf[0][0],sendrow*ncol,MPI_DOUBLE,0,0,world);
440 }
441
442 memory->destroy(buf);
443 }
444
445 /* ----------------------------------------------------------------------
446 write out Bonds section of data file
447 ------------------------------------------------------------------------- */
448
bonds()449 void WriteData::bonds()
450 {
451 // communication buffer for all my Bond info
452 // maxrow X ncol = largest buffer needed by any proc
453
454 int ncol = 3;
455 int sendrow = static_cast<int> (nbonds_local);
456 int maxrow;
457 MPI_Allreduce(&sendrow,&maxrow,1,MPI_INT,MPI_MAX,world);
458
459 tagint **buf;
460 if (me == 0) memory->create(buf,MAX(1,maxrow),ncol,"write_data:buf");
461 else memory->create(buf,MAX(1,sendrow),ncol,"write_data:buf");
462
463 // pack my bond data into buf
464
465 atom->avec->pack_bond(buf);
466
467 // write one chunk of info per proc to file
468 // proc 0 pings each proc, receives its chunk, writes to file
469 // all other procs wait for ping, send their chunk to proc 0
470
471 int tmp,recvrow;
472
473 int index = 1;
474 if (me == 0) {
475 MPI_Status status;
476 MPI_Request request;
477
478 fputs("\nBonds\n\n",fp);
479 for (int iproc = 0; iproc < nprocs; iproc++) {
480 if (iproc) {
481 MPI_Irecv(&buf[0][0],maxrow*ncol,MPI_LMP_TAGINT,iproc,0,world,&request);
482 MPI_Send(&tmp,0,MPI_INT,iproc,0,world);
483 MPI_Wait(&request,&status);
484 MPI_Get_count(&status,MPI_LMP_TAGINT,&recvrow);
485 recvrow /= ncol;
486 } else recvrow = sendrow;
487
488 atom->avec->write_bond(fp,recvrow,buf,index);
489 index += recvrow;
490 }
491
492 } else {
493 MPI_Recv(&tmp,0,MPI_INT,0,0,world,MPI_STATUS_IGNORE);
494 MPI_Rsend(&buf[0][0],sendrow*ncol,MPI_LMP_TAGINT,0,0,world);
495 }
496
497 memory->destroy(buf);
498 }
499
500 /* ----------------------------------------------------------------------
501 write out Angles section of data file
502 ------------------------------------------------------------------------- */
503
angles()504 void WriteData::angles()
505 {
506 // communication buffer for all my Angle info
507 // maxrow X ncol = largest buffer needed by any proc
508
509 int ncol = 4;
510 int sendrow = static_cast<int> (nangles_local);
511 int maxrow;
512 MPI_Allreduce(&sendrow,&maxrow,1,MPI_INT,MPI_MAX,world);
513
514 tagint **buf;
515 if (me == 0) memory->create(buf,MAX(1,maxrow),ncol,"write_data:buf");
516 else memory->create(buf,MAX(1,sendrow),ncol,"write_data:buf");
517
518 // pack my angle data into buf
519
520 atom->avec->pack_angle(buf);
521
522 // write one chunk of info per proc to file
523 // proc 0 pings each proc, receives its chunk, writes to file
524 // all other procs wait for ping, send their chunk to proc 0
525
526 int tmp,recvrow;
527
528 int index = 1;
529 if (me == 0) {
530 MPI_Status status;
531 MPI_Request request;
532
533 fputs("\nAngles\n\n",fp);
534 for (int iproc = 0; iproc < nprocs; iproc++) {
535 if (iproc) {
536 MPI_Irecv(&buf[0][0],maxrow*ncol,MPI_LMP_TAGINT,iproc,0,world,&request);
537 MPI_Send(&tmp,0,MPI_INT,iproc,0,world);
538 MPI_Wait(&request,&status);
539 MPI_Get_count(&status,MPI_LMP_TAGINT,&recvrow);
540 recvrow /= ncol;
541 } else recvrow = sendrow;
542
543 atom->avec->write_angle(fp,recvrow,buf,index);
544 index += recvrow;
545 }
546
547 } else {
548 MPI_Recv(&tmp,0,MPI_INT,0,0,world,MPI_STATUS_IGNORE);
549 MPI_Rsend(&buf[0][0],sendrow*ncol,MPI_LMP_TAGINT,0,0,world);
550 }
551
552 memory->destroy(buf);
553 }
554
555 /* ----------------------------------------------------------------------
556 write out Dihedrals section of data file
557 ------------------------------------------------------------------------- */
558
dihedrals()559 void WriteData::dihedrals()
560 {
561 // communication buffer for all my Dihedral info
562 // maxrow X ncol = largest buffer needed by any proc
563
564 int ncol = 5;
565 int sendrow = static_cast<int> (ndihedrals_local);
566 int maxrow;
567 MPI_Allreduce(&sendrow,&maxrow,1,MPI_INT,MPI_MAX,world);
568
569 tagint **buf;
570 if (me == 0) memory->create(buf,MAX(1,maxrow),ncol,"write_data:buf");
571 else memory->create(buf,MAX(1,sendrow),ncol,"write_data:buf");
572
573 // pack my dihedral data into buf
574
575 atom->avec->pack_dihedral(buf);
576
577 // write one chunk of info per proc to file
578 // proc 0 pings each proc, receives its chunk, writes to file
579 // all other procs wait for ping, send their chunk to proc 0
580
581 int tmp,recvrow;
582
583 int index = 1;
584 if (me == 0) {
585 MPI_Status status;
586 MPI_Request request;
587
588 fputs("\nDihedrals\n\n",fp);
589 for (int iproc = 0; iproc < nprocs; iproc++) {
590 if (iproc) {
591 MPI_Irecv(&buf[0][0],maxrow*ncol,MPI_LMP_TAGINT,iproc,0,world,&request);
592 MPI_Send(&tmp,0,MPI_INT,iproc,0,world);
593 MPI_Wait(&request,&status);
594 MPI_Get_count(&status,MPI_LMP_TAGINT,&recvrow);
595 recvrow /= ncol;
596 } else recvrow = sendrow;
597
598 atom->avec->write_dihedral(fp,recvrow,buf,index);
599 index += recvrow;
600 }
601
602 } else {
603 MPI_Recv(&tmp,0,MPI_INT,0,0,world,MPI_STATUS_IGNORE);
604 MPI_Rsend(&buf[0][0],sendrow*ncol,MPI_LMP_TAGINT,0,0,world);
605 }
606
607 memory->destroy(buf);
608 }
609
610 /* ----------------------------------------------------------------------
611 write out Impropers section of data file
612 ------------------------------------------------------------------------- */
613
impropers()614 void WriteData::impropers()
615 {
616 // communication buffer for all my Improper info
617 // maxrow X ncol = largest buffer needed by any proc
618
619 int ncol = 5;
620 int sendrow = static_cast<int> (nimpropers_local);
621 int maxrow;
622 MPI_Allreduce(&sendrow,&maxrow,1,MPI_INT,MPI_MAX,world);
623
624 tagint **buf;
625 if (me == 0) memory->create(buf,MAX(1,maxrow),ncol,"write_data:buf");
626 else memory->create(buf,MAX(1,sendrow),ncol,"write_data:buf");
627
628 // pack my improper data into buf
629
630 atom->avec->pack_improper(buf);
631
632 // write one chunk of info per proc to file
633 // proc 0 pings each proc, receives its chunk, writes to file
634 // all other procs wait for ping, send their chunk to proc 0
635
636 int tmp,recvrow;
637
638 int index = 1;
639 if (me == 0) {
640 MPI_Status status;
641 MPI_Request request;
642
643 fputs("\nImpropers\n\n",fp);
644 for (int iproc = 0; iproc < nprocs; iproc++) {
645 if (iproc) {
646 MPI_Irecv(&buf[0][0],maxrow*ncol,MPI_LMP_TAGINT,iproc,0,world,&request);
647 MPI_Send(&tmp,0,MPI_INT,iproc,0,world);
648 MPI_Wait(&request,&status);
649 MPI_Get_count(&status,MPI_LMP_TAGINT,&recvrow);
650 recvrow /= ncol;
651 } else recvrow = sendrow;
652
653 atom->avec->write_improper(fp,recvrow,buf,index);
654 index += recvrow;
655 }
656
657 } else {
658 MPI_Recv(&tmp,0,MPI_INT,0,0,world,MPI_STATUS_IGNORE);
659 MPI_Rsend(&buf[0][0],sendrow*ncol,MPI_LMP_TAGINT,0,0,world);
660 }
661
662 memory->destroy(buf);
663 }
664
665 /* ----------------------------------------------------------------------
666 write out Bonus sections of data file
667 flag indicates which bonus section it is
668 ------------------------------------------------------------------------- */
669
bonus(int flag)670 void WriteData::bonus(int flag)
671 {
672 // communication buffer for all my Bonus info
673 // maxvalues = largest buffer needed by any proc
674
675 int nvalues = atom->avec->pack_data_bonus(nullptr,flag);
676 int maxvalues;
677 MPI_Allreduce(&nvalues,&maxvalues,1,MPI_INT,MPI_MAX,world);
678
679 double *buf = nullptr;
680 if (me == 0) memory->create(buf,MAX(1,maxvalues),"write_data:buf");
681 else memory->create(buf,MAX(1,nvalues),"write_data:buf");
682
683 // pack my bonus data into buf
684
685 atom->avec->pack_data_bonus(buf,flag);
686
687 // write one chunk of info per proc to file
688 // proc 0 pings each proc, receives its chunk, writes to file
689 // all other procs wait for ping, send their chunk to proc 0
690
691 int tmp;
692
693 if (me == 0) {
694 MPI_Status status;
695 MPI_Request request;
696
697 if (flag == ELLIPSOID) fputs("\nEllipsoids\n\n",fp);
698 if (flag == LINE) fputs("\nLines\n\n",fp);
699 if (flag == TRIANGLE) fputs("\nTriangles\n\n",fp);
700 if (flag == BODY) fputs("\nBodies\n\n",fp);
701
702 for (int iproc = 0; iproc < nprocs; iproc++) {
703 if (iproc) {
704 MPI_Irecv(buf,maxvalues,MPI_DOUBLE,iproc,0,world,&request);
705 MPI_Send(&tmp,0,MPI_INT,iproc,0,world);
706 MPI_Wait(&request,&status);
707 MPI_Get_count(&status,MPI_DOUBLE,&nvalues);
708 }
709
710 atom->avec->write_data_bonus(fp,nvalues,buf,flag);
711 }
712
713 } else {
714 MPI_Recv(&tmp,0,MPI_INT,0,0,world,MPI_STATUS_IGNORE);
715 MPI_Rsend(buf,nvalues,MPI_DOUBLE,0,0,world);
716 }
717
718 memory->destroy(buf);
719 }
720
721 /* ----------------------------------------------------------------------
722 write out Mth section of data file owned by Fix ifix
723 ------------------------------------------------------------------------- */
724
fix(int ifix,int mth)725 void WriteData::fix(int ifix, int mth)
726 {
727 // communication buffer for Fix info
728 // maxrow X ncol = largest buffer needed by any proc
729
730 int sendrow,ncol;
731 modify->fix[ifix]->write_data_section_size(mth,sendrow,ncol);
732 int maxrow;
733 MPI_Allreduce(&sendrow,&maxrow,1,MPI_INT,MPI_MAX,world);
734
735 double **buf;
736 if (me == 0) memory->create(buf,MAX(1,maxrow),ncol,"write_data:buf");
737 else memory->create(buf,MAX(1,sendrow),ncol,"write_data:buf");
738
739 // pack my fix data into buf
740
741 modify->fix[ifix]->write_data_section_pack(mth,buf);
742
743 // write one chunk of info per proc to file
744 // proc 0 pings each proc, receives its chunk, writes to file
745 // all other procs wait for ping, send their chunk to proc 0
746
747 int tmp,recvrow;
748
749 int index = 1;
750 if (me == 0) {
751 MPI_Status status;
752 MPI_Request request;
753
754 modify->fix[ifix]->write_data_section_keyword(mth,fp);
755 for (int iproc = 0; iproc < nprocs; iproc++) {
756 if (iproc) {
757 MPI_Irecv(&buf[0][0],maxrow*ncol,MPI_DOUBLE,iproc,0,world,&request);
758 MPI_Send(&tmp,0,MPI_INT,iproc,0,world);
759 MPI_Wait(&request,&status);
760 MPI_Get_count(&status,MPI_DOUBLE,&recvrow);
761 recvrow /= ncol;
762 } else recvrow = sendrow;
763
764 modify->fix[ifix]->write_data_section(mth,fp,recvrow,buf,index);
765 index += recvrow;
766 }
767
768 } else {
769 MPI_Recv(&tmp,0,MPI_INT,0,0,world,MPI_STATUS_IGNORE);
770 MPI_Rsend(&buf[0][0],sendrow*ncol,MPI_DOUBLE,0,0,world);
771 }
772
773 memory->destroy(buf);
774 }
775