1 /* ----------------------------------------------------------------------
2 This is the
3
4 ██╗ ██╗ ██████╗ ██████╗ ██████╗ ██╗ ██╗████████╗███████╗
5 ██║ ██║██╔════╝ ██╔════╝ ██╔════╝ ██║ ██║╚══██╔══╝██╔════╝
6 ██║ ██║██║ ███╗██║ ███╗██║ ███╗███████║ ██║ ███████╗
7 ██║ ██║██║ ██║██║ ██║██║ ██║██╔══██║ ██║ ╚════██║
8 ███████╗██║╚██████╔╝╚██████╔╝╚██████╔╝██║ ██║ ██║ ███████║
9 ╚══════╝╚═╝ ╚═════╝ ╚═════╝ ╚═════╝ ╚═╝ ╚═╝ ╚═╝ ╚══════╝®
10
11 DEM simulation engine, released by
12 DCS Computing Gmbh, Linz, Austria
13 http://www.dcs-computing.com, office@dcs-computing.com
14
15 LIGGGHTS® is part of CFDEM®project:
16 http://www.liggghts.com | http://www.cfdem.com
17
18 Core developer and main author:
19 Christoph Kloss, christoph.kloss@dcs-computing.com
20
21 LIGGGHTS® is open-source, distributed under the terms of the GNU Public
22 License, version 2 or later. It is distributed in the hope that it will
23 be useful, but WITHOUT ANY WARRANTY; without even the implied warranty
24 of MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. You should have
25 received a copy of the GNU General Public License along with LIGGGHTS®.
26 If not, see http://www.gnu.org/licenses . See also top-level README
27 and LICENSE files.
28
29 LIGGGHTS® and CFDEM® are registered trade marks of DCS Computing GmbH,
30 the producer of the LIGGGHTS® software and the CFDEM®coupling software
31 See http://www.cfdem.com/terms-trademark-policy for details.
32
33 -------------------------------------------------------------------------
34 Contributing author and copyright for this file:
35 (if no contributing author is listed, this file has been contributed
36 by the core developer)
37
38 Arno Mayrhofer (DCS Computing GmbH, Linz)
39
40 Copyright 2016- DCS Computing GmbH, Linz
41 ------------------------------------------------------------------------- */
42
43 #ifdef LAMMPS_VTK
44
45 #include <cmath>
46 #include <stdlib.h>
47 #include <string.h>
48 #include "dump_local_gran.h"
49 #include "atom.h"
50 #include "force.h"
51 #include "domain.h"
52 #include "region.h"
53 #include "group.h"
54 #include "input.h"
55 #include "variable.h"
56 #include "update.h"
57 #include "modify.h"
58 #include "compute.h"
59 #include "compute_pair_gran_local.h"
60 #include "fix.h"
61 #include "memory.h"
62 #include "error.h"
63 #include <vector>
64 #include <sstream>
65 #include <vtkVersion.h>
66 #ifndef VTK_MAJOR_VERSION
67 #include <vtkConfigure.h>
68 #endif
69 #include <vtkPointData.h>
70 #include <vtkCellData.h>
71 #include <vtkLine.h>
72 #include <vtkDoubleArray.h>
73 #include <vtkIntArray.h>
74 #include <vtkStringArray.h>
75 #include <vtkPolyData.h>
76 #include <vtkRectilinearGrid.h>
77 #include <vtkHexahedron.h>
78 #include <vtkUnstructuredGrid.h>
79 #include <vtkInformation.h>
80
81 #ifdef vtkGenericDataArray_h
82 #define InsertNextTupleValue InsertNextTypedTuple
83 #endif
84
85 using namespace LAMMPS_NS;
86
87 enum{INT,DOUBLE,STRING}; // same as in DumpCFG
88 enum{X1,X2,CP,V1,V2,ID1,ID2,F,FN,FT,TORQUE,TORQUEN,TORQUET,AREA,DELTA,HEAT,MSID1,MSID2}; // dumps positions, force, normal and tangential forces, torque, normal and tangential torque
89
90 /* ---------------------------------------------------------------------- */
91
DumpLocalGran(LAMMPS * lmp,int _igroup,int _nclusterprocs,int _multiproc,int _nevery,int _filewriter,int _fileproc)92 DumpLocalGran::DumpLocalGran(LAMMPS *lmp, int _igroup, int _nclusterprocs, int _multiproc, int _nevery, int _filewriter, int _fileproc) :
93 Pointers(lmp),
94 nevery(_nevery),
95 nclusterprocs(_nclusterprocs),
96 multiproc(_multiproc),
97 filewriter(_filewriter),
98 fileproc(_fileproc),
99 iregion(-1),
100 idregion(NULL),
101 igroup(_igroup),
102 groupbit(group->bitmask[igroup]),
103 nchoose(0),
104 sortBuffer(NULL),
105 maxbuf(0),
106 buf(NULL),
107 size_one(0),
108 cpgl_(0),
109 n_calls_(0)
110 {
111 pack_choice.clear();
112 vtype.clear();
113 name.clear();
114 myarrays.clear();
115 }
116
117 /* ---------------------------------------------------------------------- */
118
~DumpLocalGran()119 DumpLocalGran::~DumpLocalGran()
120 {
121 delete [] idregion;
122 if (sortBuffer)
123 delete sortBuffer;
124 }
125
126 /* ---------------------------------------------------------------------- */
127
parse_parameters(const int narg,const char * const * const arg,bool optional_keyword,std::list<std::string> keyword_list)128 int DumpLocalGran::parse_parameters(const int narg, const char *const *const arg, bool optional_keyword, std::list<std::string> keyword_list)
129 {
130 int iarg = 0;
131
132 if (narg < 1)
133 error->all(FLERR, "dump local/gran is missing arguments");
134
135 if (strcmp(arg[iarg], "local_gran") != 0)
136 {
137 if (!optional_keyword)
138 error->all(FLERR, "Missing keyword \"local_gran\" in dump local/gran");
139 }
140 else
141 iarg++;
142
143 if (narg < 1+iarg)
144 error->all(FLERR, "dump local/gran is missing arguments");
145
146 Compute *comp = modify->find_compute_id(arg[iarg++]);
147 if(!comp || !dynamic_cast<ComputePairGranLocal*>(comp))
148 error->all(FLERR,"dump local/gran requires a valid ID of a compute pair/gran/local to be provided");
149
150 cpgl_ = static_cast<ComputePairGranLocal*>(comp);
151
152 // error if compute does not write pos
153 if(cpgl_->offset_x1() < 0 || cpgl_->offset_x2() < 0)
154 error->all(FLERR,"dump local/gran requires a valid ID of a compute pair/gran/local that writes the positions");
155
156 // do stuff which needs the cpgl_ ptr
157 size_one = cpgl_->get_nvalues();
158
159 // fill data into containers
160 define_properties();
161
162 if (filewriter) reset_vtk_data_containers();
163
164 if (narg > iarg)
165 {
166 std::list<std::string>::iterator it;
167 bool found_keyword = false;
168 for (it = keyword_list.begin(); it != keyword_list.end(); it++)
169 {
170 if (it->compare(arg[iarg]) == 0)
171 {
172 found_keyword = true;
173 break;
174 }
175 }
176 if (!found_keyword)
177 error->all(FLERR, "Could not find follow-up keyword in local/gran");
178 }
179
180 return iarg;
181 }
182
183 /* ---------------------------------------------------------------------- */
184
init_style()185 void DumpLocalGran::init_style()
186 {
187 // set index and check validity of region
188
189 if (iregion >= 0) {
190 iregion = domain->find_region(idregion);
191 if (iregion == -1)
192 error->all(FLERR,"Region ID for dump custom/vtk does not exist");
193 }
194 }
195
196 /* ---------------------------------------------------------------------- */
197
count()198 int DumpLocalGran::count()
199 {
200
201 n_calls_ = 0;
202
203 //TODO generalize
204 //also check if have same length
205 cpgl_->compute_local();
206 cpgl_->invoked_flag |= INVOKED_LOCAL;
207
208 nchoose = cpgl_->get_ncount();
209
210 return nchoose;
211 }
212
213 /* ---------------------------------------------------------------------- */
214
prepare_mbSet(vtkSmartPointer<vtkMultiBlockDataSet> mbSet,bool use_poly_data)215 void DumpLocalGran::prepare_mbSet(vtkSmartPointer<vtkMultiBlockDataSet> mbSet, bool use_poly_data)
216 {
217
218 // nme = # of dump lines this proc contributes to dump
219
220 int nme = count();
221
222 // ntotal = total # of dump lines in snapshot
223 // nmax = max # of dump lines on any proc
224
225 bigint bnme = nme;
226 bigint ntotal = nme;
227 MPI_Allreduce(&bnme,&ntotal,1,MPI_LMP_BIGINT,MPI_SUM,world);
228
229 int nmax;
230 if (multiproc != comm->nprocs) MPI_Allreduce(&nme,&nmax,1,MPI_INT,MPI_MAX,world);
231 else nmax = nme;
232
233 // ensure buf is sized for packing and communicating
234 // use nmax to ensure filewriter proc can receive info from others
235 // limit nmax*size_one to int since used as arg in MPI calls
236
237 if (nmax > maxbuf) {
238 if ((bigint) nmax * size_one > MAXSMALLINT)
239 error->all(FLERR,"Too much per-proc info for dump");
240 maxbuf = nmax;
241 memory->destroy(buf);
242 memory->create(buf,maxbuf*size_one,"dump:buf");
243
244 }
245
246 // ensure ids buffer is sized for sorting
247 if (sortBuffer)
248 sortBuffer->realloc_ids(nmax);
249
250 // pack my data into buf
251 // if sorting on IDs also request ID list from pack()
252 // sort buf as needed
253
254 if (sortBuffer)
255 pack(sortBuffer->get_ids());
256 else
257 pack(NULL);
258 if (sortBuffer)
259 sortBuffer->sort(buf, nme, maxbuf, size_one, ntotal);
260
261 // filewriter = 1 = this proc writes to file
262 // ping each proc in my cluster, receive its data, write data to file
263 // else wait for ping from fileproc, send my data to fileproc
264
265 int tmp,nlines;
266 MPI_Status status;
267 MPI_Request request;
268
269 // comm and output buf of doubles
270
271 if (filewriter) {
272 for (int iproc = 0; iproc < nclusterprocs; iproc++) {
273 if (iproc) {
274 MPI_Irecv(buf,maxbuf*size_one,MPI_DOUBLE,comm->me+iproc,0,world,&request);
275 MPI_Send(&tmp,0,MPI_INT,comm->me+iproc,0,world);
276 MPI_Wait(&request,&status);
277 MPI_Get_count(&status,MPI_DOUBLE,&nlines);
278 nlines /= size_one;
279 } else nlines = nme;
280
281 write_data(nlines, buf, mbSet, use_poly_data);
282 }
283 } else {
284 MPI_Recv(&tmp,0,MPI_INT,fileproc,0,world,&status);
285 MPI_Rsend(buf,nme*size_one,MPI_DOUBLE,fileproc,0,world);
286 }
287
288 }
289
290 /* ---------------------------------------------------------------------- */
291
pack(int * ids)292 void DumpLocalGran::pack(int *ids)
293 {
294 int n = 0;
295 for (std::map<int,FnPtrPack>::iterator it = pack_choice.begin(); it != pack_choice.end(); ++it)
296 {
297 (this->*(it->second))(n);
298
299 // increase n by length of data
300 if(vector_set.find(it->first) != vector_set.end())
301 n += 3;
302 else
303 n++;
304
305 }
306
307 // similar to dump local, IDs are not used here (since are atom IDS
308 }
309
310 /* ---------------------------------------------------------------------- */
311
buf2arrays(int n,double * mybuf)312 void DumpLocalGran::buf2arrays(int n, double *mybuf)
313 {
314
315 const bool have_cp = cpgl_->offset_contact_point() >= 0;
316
317 for (int idata=0; idata < n; ++idata) {
318
319 // stores the ID of newly added points
320 vtkIdType pid[2];
321
322 pid[0] = points->InsertNextPoint(mybuf[idata*size_one],mybuf[idata*size_one+1],mybuf[idata*size_one+2]);
323 pid[1] = points->InsertNextPoint(mybuf[idata*size_one+3],mybuf[idata*size_one+4],mybuf[idata*size_one+5]);
324
325 // define the line going from point pid[0] to pid[1]
326 vtkSmartPointer<vtkLine> line0 = vtkSmartPointer<vtkLine>::New();
327 line0->GetPointIds()->SetId(0,pid[0]);
328 line0->GetPointIds()->SetId(1,pid[1]);
329
330 lineCells->InsertNextCell(line0);
331 if(have_cp)
332 {
333 vtkIdType pidCP[1];
334 pidCP[0] = points->InsertNextPoint(mybuf[idata*size_one+6],mybuf[idata*size_one+7],mybuf[idata*size_one+8]);
335 line0->GetPointIds()->SetId(0,pid[0]);
336 line0->GetPointIds()->SetId(1,pidCP[0]);
337 lineCells->InsertNextCell(line0);
338 line0->GetPointIds()->SetId(0,pid[1]);
339 line0->GetPointIds()->SetId(1,pidCP[0]);
340 lineCells->InsertNextCell(line0);
341 }
342
343 int j = 6; // 0,1,2,3,4,5 = 2 * (x,y,z) handled just above
344 if(have_cp)
345 j += 3;
346 for (std::map<int, vtkSmartPointer<vtkAbstractArray> >::iterator it=myarrays.begin(); it!=myarrays.end(); ++it) {
347
348 vtkAbstractArray *paa = it->second;
349 if (it->second->GetNumberOfComponents() == 3) {
350 switch (vtype[it->first]) {
351 case INT:
352 {
353 int iv3[3] = { static_cast<int>(mybuf[idata*size_one+j ]),
354 static_cast<int>(mybuf[idata*size_one+j+1]),
355 static_cast<int>(mybuf[idata*size_one+j+2]) };
356 vtkIntArray *pia = static_cast<vtkIntArray*>(paa);
357 pia->InsertNextTupleValue(iv3);
358 if (have_cp)
359 {
360 pia->InsertNextTupleValue(iv3);
361 pia->InsertNextTupleValue(iv3);
362 }
363 break;
364 }
365 case DOUBLE:
366 {
367 vtkDoubleArray *pda = static_cast<vtkDoubleArray*>(paa);
368 pda->InsertNextTupleValue(&mybuf[idata*size_one+j]);
369 if (have_cp)
370 {
371 pda->InsertNextTupleValue(&mybuf[idata*size_one+j]);
372 pda->InsertNextTupleValue(&mybuf[idata*size_one+j]);
373 }
374 break;
375 }
376 }
377 j+=3;
378 } else {
379 switch (vtype[it->first]) {
380 case INT:
381 {
382 vtkIntArray *pia = static_cast<vtkIntArray*>(paa);
383 pia->InsertNextValue(mybuf[idata*size_one+j]);
384 if (have_cp)
385 {
386 pia->InsertNextValue(mybuf[idata*size_one+j]);
387 pia->InsertNextValue(mybuf[idata*size_one+j]);
388 }
389 break;
390 }
391 case DOUBLE:
392 {
393 vtkDoubleArray *pda = static_cast<vtkDoubleArray*>(paa);
394 pda->InsertNextValue(mybuf[idata*size_one+j]);
395 if (have_cp)
396 {
397 pda->InsertNextValue(mybuf[idata*size_one+j]);
398 pda->InsertNextValue(mybuf[idata*size_one+j]);
399 }
400 break;
401 }
402 }
403 ++j;
404 }
405 }
406 }
407 }
408
409 /* ---------------------------------------------------------------------- */
410
write_data(int n,double * mybuf,vtkSmartPointer<vtkMultiBlockDataSet> mbSet,bool use_poly_data)411 void DumpLocalGran::write_data(int n, double *mybuf, vtkSmartPointer<vtkMultiBlockDataSet> mbSet, bool use_poly_data)
412 {
413 ++n_calls_;
414 int cur_block = mbSet->GetNumberOfBlocks();
415
416 buf2arrays(n, mybuf);
417
418 if (n_calls_ < nclusterprocs)
419 return; // multiple processors but only proc 0 is a filewriter (-> nclusterprocs procs contribute to the filewriter's output data)
420
421 if (!use_poly_data)
422 {
423 vtkSmartPointer<vtkUnstructuredGrid> unstructuredGrid = vtkSmartPointer<vtkUnstructuredGrid>::New();
424 unstructuredGrid->SetPoints(points);
425 unstructuredGrid->SetCells(VTK_LINE, lineCells);
426
427 for (std::map<int, vtkSmartPointer<vtkAbstractArray> >::iterator it=myarrays.begin(); it!=myarrays.end(); ++it) {
428 unstructuredGrid->GetCellData()->AddArray(it->second);
429 }
430 mbSet->SetBlock(cur_block++, unstructuredGrid);
431
432 }
433 else
434 {
435 vtkSmartPointer<vtkPolyData> polyData = vtkSmartPointer<vtkPolyData>::New();
436 polyData->SetPoints(points);
437 polyData->SetLines(lineCells);
438
439 for (std::map<int, vtkSmartPointer<vtkAbstractArray> >::iterator it=myarrays.begin(); it!=myarrays.end(); ++it) {
440 polyData->GetCellData()->AddArray(it->second);
441 }
442 mbSet->SetBlock(cur_block++, polyData);
443 }
444 std::string name = "local_gran_";
445 name.append(cpgl_->id);
446 mbSet->GetMetaData(cur_block-1)->Set(mbSet->NAME(), name.c_str());
447
448 reset_vtk_data_containers();
449 }
450
451 /* ---------------------------------------------------------------------- */
452
reset_vtk_data_containers()453 void DumpLocalGran::reset_vtk_data_containers()
454 {
455 points = vtkSmartPointer<vtkPoints>::New();
456 lineCells = vtkSmartPointer<vtkCellArray>::New();
457
458 std::map<int,int>::iterator it=vtype.begin();
459
460 ++it;
461 ++it;
462 if(cpgl_->offset_contact_point() >= 0)
463 ++it;
464
465 for (; it!=vtype.end(); ++it) {
466
467 // part 1: add VTK array to myarrays
468 switch(vtype[it->first]) {
469 case INT:
470 myarrays[it->first] = vtkSmartPointer<vtkIntArray>::New();
471 break;
472 case DOUBLE:
473 myarrays[it->first] = vtkSmartPointer<vtkDoubleArray>::New();
474 break;
475 case STRING:
476 myarrays[it->first] = vtkSmartPointer<vtkStringArray>::New();
477 break;
478 }
479
480 // part 2: if vector, set length to 3; set name
481 if (vector_set.find(it->first) != vector_set.end()) {
482 myarrays[it->first]->SetNumberOfComponents(3);
483 myarrays[it->first]->SetName(name[it->first].c_str());
484 } else {
485 myarrays[it->first]->SetName(name[it->first].c_str());
486 }
487 }
488 }
489
490 /* ---------------------------------------------------------------------- */
491
492 // customize here
493
define_properties()494 void DumpLocalGran::define_properties()
495 {
496 pack_choice[X1] = &DumpLocalGran::pack_x1;
497 vtype[X1] = DOUBLE;
498 name[X1] = "pos1";
499 vector_set.insert(X1);
500
501 pack_choice[X2] = &DumpLocalGran::pack_x2;
502 vtype[X2] = DOUBLE;
503 name[X2] = "pos2";
504 vector_set.insert(X2);
505
506 if(cpgl_->offset_contact_point() >= 0)
507 {
508 pack_choice[CP] = &DumpLocalGran::pack_contact_point;
509 vtype[CP] = DOUBLE;
510 name[CP] = "contact_point";
511 vector_set.insert(CP);
512 }
513
514 if(cpgl_->offset_v1() >= 0)
515 {
516 pack_choice[V1] = &DumpLocalGran::pack_v1;
517 vtype[V1] = DOUBLE;
518 name[V1] = "vel1";
519 vector_set.insert(V1);
520 }
521
522 if(cpgl_->offset_v2() >= 0)
523 {
524 pack_choice[V2] = &DumpLocalGran::pack_v2;
525 vtype[V2] = DOUBLE;
526 name[V2] = "vel2";
527 vector_set.insert(V2);
528 }
529
530 if(cpgl_->offset_id1() >= 0)
531 {
532 pack_choice[ID1] = &DumpLocalGran::pack_id1;
533 vtype[ID1] = DOUBLE;
534 name[ID1] = "id1";
535 //scalar
536 }
537
538 if(cpgl_->offset_id2() >= 0)
539 {
540 pack_choice[ID2] = &DumpLocalGran::pack_id2;
541 vtype[ID2] = DOUBLE;
542 name[ID2] = "id2";
543 //scalar
544 }
545
546 if(cpgl_->offset_f() >= 0)
547 {
548 pack_choice[F] = &DumpLocalGran::pack_f;
549 vtype[F] = DOUBLE;
550 name[F] = "force";
551 vector_set.insert(F);
552 }
553
554 if(cpgl_->offset_fn() >= 0)
555 {
556 pack_choice[FN] = &DumpLocalGran::pack_fn;
557 vtype[FN] = DOUBLE;
558 name[FN] = "force_normal";
559 vector_set.insert(FN);
560 }
561
562 if(cpgl_->offset_ft() >= 0)
563 {
564 pack_choice[FT] = &DumpLocalGran::pack_ft;
565 vtype[FT] = DOUBLE;
566 name[FT] = "force_tangential";
567 vector_set.insert(FT);
568 }
569
570 if(cpgl_->offset_torque() >= 0)
571 {
572 pack_choice[TORQUE] = &DumpLocalGran::pack_torque;
573 vtype[TORQUE] = DOUBLE;
574 name[TORQUE] = "torque";
575 vector_set.insert(TORQUE);
576 }
577
578 if(cpgl_->offset_torquen() >= 0)
579 {
580 pack_choice[TORQUEN] = &DumpLocalGran::pack_torquen;
581 vtype[TORQUEN] = DOUBLE;
582 name[TORQUEN] = "torque_normal";
583 vector_set.insert(TORQUEN);
584 }
585
586 if(cpgl_->offset_torquet() >= 0)
587 {
588 pack_choice[TORQUET] = &DumpLocalGran::pack_torquet;
589 vtype[TORQUET] = DOUBLE;
590 name[TORQUET] = "torque_tangential";
591 vector_set.insert(TORQUET);
592 }
593
594 if(cpgl_->offset_area() >= 0)
595 {
596
597 pack_choice[AREA] = &DumpLocalGran::pack_area;
598 vtype[AREA] = DOUBLE;
599 name[AREA] = "contact_area";
600 //scalar
601 }
602
603 if(cpgl_->offset_delta() >= 0)
604 {
605 pack_choice[DELTA] = &DumpLocalGran::pack_delta;
606 vtype[DELTA] = DOUBLE;
607 name[DELTA] = "delta";
608 //scalar
609 }
610
611 if(cpgl_->offset_heat() >= 0)
612 {
613 pack_choice[HEAT] = &DumpLocalGran::pack_heat;
614 vtype[HEAT] = DOUBLE;
615 name[HEAT] = "heat_flux";
616 //scalar
617 }
618
619 if(cpgl_->offset_ms_id1() >= 0)
620 {
621 pack_choice[MSID1] = &DumpLocalGran::pack_ms_id1;
622 vtype[MSID1] = DOUBLE;
623 name[MSID1] = "ms_id1";
624 //scalar
625 }
626
627 if(cpgl_->offset_ms_id2() >= 0)
628 {
629 pack_choice[MSID2] = &DumpLocalGran::pack_ms_id2;
630 vtype[MSID2] = DOUBLE;
631 name[MSID2] = "ms_id2";
632 //scalar
633 }
634 }
635
636 /* ---------------------------------------------------------------------- */
637
modify_param(int narg,char ** arg)638 int DumpLocalGran::modify_param(int narg, char **arg)
639 {
640 if (strcmp(arg[0],"region") == 0) {
641 if (narg < 2) error->all(FLERR,"Illegal dump_modify command");
642 if (strcmp(arg[1],"none") == 0) iregion = -1;
643 else {
644 iregion = domain->find_region(arg[1]);
645 if (iregion == -1)
646 error->all(FLERR,"Dump_modify region ID does not exist");
647 delete [] idregion;
648 int n = strlen(arg[1]) + 1;
649 idregion = new char[n];
650 strcpy(idregion,arg[1]);
651 }
652 return 2;
653 }
654
655 if (strcmp(arg[0],"element") == 0) {
656 error->all(FLERR,"Dump local/gran does not support dump_modify 'element' ");
657 return 0;
658 }
659
660 if (strcmp(arg[0],"thresh") == 0) {
661 error->all(FLERR,"Dump local/gran does not support dump_modify 'thresh' ");
662 }
663
664 if (!sortBuffer)
665 sortBuffer = new SortBuffer(lmp, false);
666
667 return sortBuffer->modify_param(narg, arg);
668 }
669
670 /* ----------------------------------------------------------------------
671 return # of bytes of allocated memory in buf, choose, variable arrays
672 ------------------------------------------------------------------------- */
673
memory_usage()674 bigint DumpLocalGran::memory_usage()
675 {
676 bigint bytes = memory->usage(buf,maxbuf*size_one);
677 if (sortBuffer)
678 bytes += sortBuffer->memory_usage(size_one);
679 return bytes;
680 }
681
682 /* ----------------------------------------------------------------------
683 pack properties
684 TODO generalize function, and length of packing
685 ------------------------------------------------------------------------- */
686
687 // customize here
688
pack_x1(int n)689 void DumpLocalGran::pack_x1(int n)
690 {
691 int offset = cpgl_->offset_x1();
692
693 for (int i = 0; i < nchoose; i++) {
694 vectorCopy3D(&cpgl_->get_data()[i][offset],&buf[n]);
695 n += size_one;
696 }
697 }
698
pack_x2(int n)699 void DumpLocalGran::pack_x2(int n)
700 {
701 int offset = cpgl_->offset_x2();
702
703 for (int i = 0; i < nchoose; i++) {
704 vectorCopy3D(&cpgl_->get_data()[i][offset],&buf[n]);
705 n += size_one;
706 }
707 }
708
pack_v1(int n)709 void DumpLocalGran::pack_v1(int n)
710 {
711 int offset = cpgl_->offset_v1();
712
713 for (int i = 0; i < nchoose; i++) {
714 vectorCopy3D(&cpgl_->get_data()[i][offset],&buf[n]);
715 n += size_one;
716 }
717 }
718
pack_v2(int n)719 void DumpLocalGran::pack_v2(int n)
720 {
721 int offset = cpgl_->offset_v2();
722
723 for (int i = 0; i < nchoose; i++) {
724 vectorCopy3D(&cpgl_->get_data()[i][offset],&buf[n]);
725 n += size_one;
726 }
727 }
728
pack_id1(int n)729 void DumpLocalGran::pack_id1(int n)
730 {
731 int offset = cpgl_->offset_id1();
732
733 for (int i = 0; i < nchoose; i++) {
734 buf[n] = cpgl_->get_data()[i][offset];
735 n += size_one;
736 }
737 }
738
pack_id2(int n)739 void DumpLocalGran::pack_id2(int n)
740 {
741 int offset = cpgl_->offset_id2();
742
743 for (int i = 0; i < nchoose; i++) {
744 buf[n] = cpgl_->get_data()[i][offset];
745 n += size_one;
746 }
747 }
748
pack_f(int n)749 void DumpLocalGran::pack_f(int n)
750 {
751 int offset = cpgl_->offset_f();
752
753 for (int i = 0; i < nchoose; i++) {
754 vectorCopy3D(&cpgl_->get_data()[i][offset],&buf[n]);
755 n += size_one;
756 }
757 }
758
pack_fn(int n)759 void DumpLocalGran::pack_fn(int n)
760 {
761 int offset = cpgl_->offset_fn();
762
763 for (int i = 0; i < nchoose; i++) {
764 vectorCopy3D(&cpgl_->get_data()[i][offset],&buf[n]);
765 n += size_one;
766 }
767 }
768
pack_ft(int n)769 void DumpLocalGran::pack_ft(int n)
770 {
771 int offset = cpgl_->offset_ft();
772
773 for (int i = 0; i < nchoose; i++) {
774 vectorCopy3D(&cpgl_->get_data()[i][offset],&buf[n]);
775 n += size_one;
776 }
777 }
778
pack_torque(int n)779 void DumpLocalGran::pack_torque(int n)
780 {
781 int offset = cpgl_->offset_torque();
782
783 for (int i = 0; i < nchoose; i++) {
784 vectorCopy3D(&cpgl_->get_data()[i][offset],&buf[n]);
785 n += size_one;
786 }
787 }
788
pack_torquen(int n)789 void DumpLocalGran::pack_torquen(int n)
790 {
791 int offset = cpgl_->offset_torquen();
792
793 for (int i = 0; i < nchoose; i++) {
794 vectorCopy3D(&cpgl_->get_data()[i][offset],&buf[n]);
795 n += size_one;
796 }
797 }
798
pack_torquet(int n)799 void DumpLocalGran::pack_torquet(int n)
800 {
801 int offset = cpgl_->offset_torquet();
802
803 for (int i = 0; i < nchoose; i++) {
804 vectorCopy3D(&cpgl_->get_data()[i][offset],&buf[n]);
805 n += size_one;
806 }
807 }
808
pack_area(int n)809 void DumpLocalGran::pack_area(int n)
810 {
811 int offset = cpgl_->offset_area();
812
813 for (int i = 0; i < nchoose; i++) {
814 buf[n] = cpgl_->get_data()[i][offset];
815 n += size_one;
816 }
817 }
818
pack_delta(int n)819 void DumpLocalGran::pack_delta(int n)
820 {
821 int offset = cpgl_->offset_delta();
822
823 for (int i = 0; i < nchoose; i++) {
824 buf[n] = cpgl_->get_data()[i][offset];
825 n += size_one;
826 }
827 }
828
pack_heat(int n)829 void DumpLocalGran::pack_heat(int n)
830 {
831 int offset = cpgl_->offset_heat();
832
833 for (int i = 0; i < nchoose; i++) {
834 buf[n] = cpgl_->get_data()[i][offset];
835 n += size_one;
836 }
837 }
pack_contact_point(int n)838 void DumpLocalGran::pack_contact_point(int n)
839 {
840 int offset = cpgl_->offset_contact_point();
841
842 for (int i = 0; i < nchoose; i++) {
843 vectorCopy3D(&cpgl_->get_data()[i][offset],&buf[n]);
844 n += size_one;
845 }
846 }
847
pack_ms_id1(int n)848 void DumpLocalGran::pack_ms_id1(int n)
849 {
850 int offset = cpgl_->offset_ms_id1();
851
852 for (int i = 0; i < nchoose; i++) {
853 buf[n] = cpgl_->get_data()[i][offset];
854 n += size_one;
855 }
856 }
857
pack_ms_id2(int n)858 void DumpLocalGran::pack_ms_id2(int n)
859 {
860 int offset = cpgl_->offset_ms_id2();
861
862 for (int i = 0; i < nchoose; i++) {
863 buf[n] = cpgl_->get_data()[i][offset];
864 n += size_one;
865 }
866 }
867
868 #endif
869