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 Anton Gladky(TU Bergakademie Freiberg), gladky.anton@gmail.com
36 ------------------------------------------------------------------------- */
37
38 #ifdef LAMMPS_VTK
39 #include <string.h>
40 #include "dump_atom_vtk.h"
41 #include "atom.h"
42 #include "group.h"
43 #include "error.h"
44 #include "memory.h"
45 #include "comm.h"
46 #include <sstream>
47 #include <vtkVersion.h>
48 #ifndef VTK_MAJOR_VERSION
49 #include <vtkConfigure.h>
50 #endif
51 #include<vtkCellArray.h>
52 #include<vtkDoubleArray.h>
53 #include<vtkIntArray.h>
54 #include<vtkPoints.h>
55 #include<vtkPointData.h>
56 #include<vtkCellData.h>
57 #include<vtkSmartPointer.h>
58 #include<vtkUnstructuredGrid.h>
59 #include<vtkXMLUnstructuredGridWriter.h>
60
61 #ifdef vtkGenericDataArray_h
62 #define InsertNextTupleValue InsertNextTypedTuple
63 #endif
64
65 using namespace LAMMPS_NS;
66
67 /* ---------------------------------------------------------------------- */
68
DumpATOMVTK(LAMMPS * lmp,int narg,char ** arg)69 DumpATOMVTK::DumpATOMVTK(LAMMPS *lmp, int narg, char **arg) :
70 Dump(lmp, narg, arg),
71 tmpEXP(lmp)
72 {
73 if (narg != 5) error->all(FLERR,"Illegal dump command");
74 if (multiproc) error->all(FLERR,"Invalid dump filename");
75
76 sortBuffer = new SortBuffer(lmp, true);
77
78 size_one = 17;
79
80 char *str = (char *) "%d %g %g %g";
81 int n = strlen(str) + 1;
82 format_default = new char[n];
83 strcpy(format_default,str);
84 }
85
86 /* ---------------------------------------------------------------------- */
87
init_style()88 void DumpATOMVTK::init_style()
89 {
90
91 }
92
93 /* ---------------------------------------------------------------------- */
94
write_header(bigint)95 void DumpATOMVTK::write_header(bigint /* n */)
96 {
97 }
98
99 /* ---------------------------------------------------------------------- */
100
count()101 int DumpATOMVTK::count()
102 {
103 n_calls_ = 0;
104
105 return Dump::count();
106 }
107
108 /* ---------------------------------------------------------------------- */
109
pack(int * ids)110 void DumpATOMVTK::pack(int *ids)
111 {
112 int n = 0;
113 int m = 0;
114
115 int *tag = atom->tag;
116 int *type = atom->type;
117 int *mask = atom->mask;
118 double *rmass = atom->rmass;
119 double *mass = atom->mass;
120 double **x = atom->x;
121 double **v = atom->v;
122 double **f = atom->f;
123 double **o = atom->omega;
124 int nlocal = atom->nlocal;
125
126 for (int i = 0; i < nlocal; i++) {
127 if (mask[i] & groupbit) {
128
129 if (ids) ids[n++] = tag[i];
130
131 double massTemp;
132
133 if (rmass) {
134 massTemp=rmass[i];
135 } else {
136 massTemp=mass[type[i]];
137 }
138
139 int me = comm->me;
140
141 buf[m++] = x[i][0];
142 buf[m++] = x[i][1];
143 buf[m++] = x[i][2];
144
145 buf[m++] = atom->radius[i];
146 buf[m++] = massTemp;
147 buf[m++] = static_cast<double>(tag[i]);
148 buf[m++] = static_cast<double>(type[i]);
149
150 buf[m++] = v[i][0];
151 buf[m++] = v[i][1];
152 buf[m++] = v[i][2];
153
154 buf[m++] = o[i][0];
155 buf[m++] = o[i][1];
156 buf[m++] = o[i][2];
157
158 buf[m++] = f[i][0];
159 buf[m++] = f[i][1];
160 buf[m++] = f[i][2];
161
162 buf[m++] = static_cast<double>(me);
163 }
164 }
165
166 setFileCurrent();
167 tmpEXP.setFileName(filecurrent);
168 return;
169 }
170
171 /* ---------------------------------------------------------------------- */
172
write_data(int n,double * mybuf)173 void DumpATOMVTK::write_data(int n, double *mybuf)
174 {
175 if (comm->me != 0) return;
176 n_calls_++;
177 int m = 0;
178 for (int i = 0; i < n; i++) {
179 DumpATOMVTK::DataVTK tmpVTKDat(
180 V3(mybuf[m+0], mybuf[m+1], mybuf[m+2]),
181 mybuf[m+3], mybuf[m+4], static_cast<int> (mybuf[m+5]), static_cast<int> (mybuf[m+6]),
182 V3(mybuf[m+7], mybuf[m+8], mybuf[m+9]),
183 V3(mybuf[m+10], mybuf[m+11], mybuf[m+12]),
184 V3(mybuf[m+13], mybuf[m+14], mybuf[m+15]),
185 static_cast<int> (mybuf[m+16]));
186 tmpEXP.add(tmpVTKDat);
187 m += size_one;
188 }
189
190 if(n_calls_ == comm->nprocs) {
191 tmpEXP.writeSER();
192 tmpEXP.clear();
193 delete [] filecurrent;
194 }
195 }
196
197 /* ---------------------------------------------------------------------- */
198
DataVTK(V3 Pos,double Rad,double Mass,int Id,int Type,V3 VelL,V3 VelA,V3 Force,int proc)199 DumpATOMVTK::DataVTK::DataVTK(V3 Pos, double Rad, double Mass, int Id, int Type,
200 V3 VelL, V3 VelA, V3 Force, int proc) {
201 _Pos = Pos;
202 _Rad = Rad;
203 _Mass = Mass;
204 _Id = Id;
205 _Type = Type;
206 _VelL = VelL;
207 _VelA = VelA;
208 _Force = Force;
209 _proc = proc;
210 }
211 /* ---------------------------------------------------------------------- */
212
serialize()213 std::string DumpATOMVTK::DataVTK::serialize() {
214 std::string tmp;
215 std::ostringstream stringStream;
216
217 stringStream <<_Pos[0]<<' '<<_Pos[1]<<' '<<_Pos[2]<<' '<<_Rad<<' '<<_Mass<<' '<<_Id
218 <<' '<<_Type<<' '
219 <<_VelL[0]<<' '<<_VelL[1]<<' '<<_VelL[2]<<' '
220 <<_VelA[0]<<' '<<_VelA[1]<<' '<<_VelA[2]<<' '
221 <<_Force[0]<<' '<<_Force[1]<<' '<<_Force[2]<<' '<<_proc<<'\n';
222
223 tmp = stringStream.str();
224 return tmp;
225 }
226
227 /* ---------------------------------------------------------------------- */
228
add(DumpATOMVTK::DataVTK & d)229 void DumpATOMVTK::vtkExportData::add(DumpATOMVTK::DataVTK & d) {
230 vtkData.push_back(d);
231 }
232
233 /* ---------------------------------------------------------------------- */
234
setFileName(const char * fileName)235 void DumpATOMVTK::vtkExportData::setFileName(const char * fileName) {
236 _fileName = fileName;
237 _setFileName = true;
238 }
239
240 /* ---------------------------------------------------------------------- */
241
vtkExportData(LAMMPS * lmp)242 DumpATOMVTK::vtkExportData::vtkExportData(LAMMPS *lmp) :
243 DumpVTK(lmp)
244 {
245 _setFileName=false;
246 }
247 /* ---------------------------------------------------------------------- */
248
size()249 int DumpATOMVTK::vtkExportData::size() {
250 return vtkData.size();
251 }
252
253 /* ---------------------------------------------------------------------- */
254
writeSER()255 void DumpATOMVTK::vtkExportData::writeSER() {
256
257 vtkSmartPointer<vtkPoints> spheresPos = vtkSmartPointer<vtkPoints>::New();
258 vtkSmartPointer<vtkCellArray> spheresCells = vtkSmartPointer<vtkCellArray>::New();
259
260 vtkSmartPointer<vtkDoubleArray> radii = vtkSmartPointer<vtkDoubleArray>::New();
261 radii->SetNumberOfComponents(1);
262 radii->SetName("radii");
263
264 vtkSmartPointer<vtkDoubleArray> spheresMass = vtkSmartPointer<vtkDoubleArray>::New();
265 spheresMass->SetNumberOfComponents(1);
266 spheresMass->SetName("mass");
267
268 vtkSmartPointer<vtkIntArray> spheresId = vtkSmartPointer<vtkIntArray>::New();
269 spheresId->SetNumberOfComponents(1);
270 spheresId->SetName("id");
271
272 vtkSmartPointer<vtkIntArray> spheresType = vtkSmartPointer<vtkIntArray>::New();
273 spheresType->SetNumberOfComponents(1);
274 spheresType->SetName("type");
275
276 vtkSmartPointer<vtkIntArray> spheresProc = vtkSmartPointer<vtkIntArray>::New();
277 spheresProc->SetNumberOfComponents(1);
278 spheresProc->SetName("proc");
279
280 vtkSmartPointer<vtkDoubleArray> spheresVelL = vtkSmartPointer<vtkDoubleArray>::New();
281 spheresVelL->SetNumberOfComponents(3);
282 spheresVelL->SetName("velocity_lin");
283
284 vtkSmartPointer<vtkDoubleArray> spheresVelA = vtkSmartPointer<vtkDoubleArray>::New();
285 spheresVelA->SetNumberOfComponents(3);
286 spheresVelA->SetName("velocity_ang");
287
288 vtkSmartPointer<vtkDoubleArray> spheresForce = vtkSmartPointer<vtkDoubleArray>::New();
289 spheresForce->SetNumberOfComponents(3);
290 spheresForce->SetName("force");
291
292 for (unsigned int i=0; i < vtkData.size(); i++) {
293 vtkIdType pid[1];
294 pid[0] = spheresPos->InsertNextPoint(vtkData[i]._Pos[0], vtkData[i]._Pos[1], vtkData[i]._Pos[2]);
295 radii->InsertNextValue(vtkData[i]._Rad);
296
297 double vv[3] = {vtkData[i]._VelL[0], vtkData[i]._VelL[1], vtkData[i]._VelL[2]};
298 spheresVelL->InsertNextTupleValue(vv);
299
300 double oo[3] = {vtkData[i]._VelA[0], vtkData[i]._VelA[1], vtkData[i]._VelA[2]};
301 spheresVelA->InsertNextTupleValue(oo);
302
303 double ff[3] = {vtkData[i]._Force[0], vtkData[i]._Force[1], vtkData[i]._Force[2]};
304 spheresForce->InsertNextTupleValue(ff);
305
306 spheresMass->InsertNextValue(vtkData[i]._Mass);
307
308 spheresId->InsertNextValue(vtkData[i]._Id);
309 spheresType->InsertNextValue(vtkData[i]._Type);
310 spheresProc->InsertNextValue(vtkData[i]._proc);
311
312 spheresCells->InsertNextCell(1,pid);
313 }
314
315 vtkSmartPointer<vtkUnstructuredGrid> spheresUg = vtkSmartPointer<vtkUnstructuredGrid>::New();
316
317 spheresUg->SetPoints(spheresPos);
318 spheresUg->SetCells(VTK_VERTEX, spheresCells);
319 spheresUg->GetPointData()->AddArray(radii);
320 spheresUg->GetPointData()->AddArray(spheresId);
321 spheresUg->GetPointData()->AddArray(spheresType);
322 spheresUg->GetPointData()->AddArray(spheresProc);
323 spheresUg->GetPointData()->AddArray(spheresMass);
324 spheresUg->GetPointData()->AddArray(spheresVelL);
325 spheresUg->GetPointData()->AddArray(spheresVelA);
326 spheresUg->GetPointData()->AddArray(spheresForce);
327
328 vtkSmartPointer<vtkXMLUnstructuredGridWriter> writer = vtkSmartPointer<vtkXMLUnstructuredGridWriter>::New();
329 setVtkWriterOptions(vtkXMLWriter::SafeDownCast(writer));
330 #if VTK_MAJOR_VERSION < 6
331 writer->SetInput(spheresUg);
332 #else
333 writer->SetInputData(spheresUg);
334 #endif
335 writer->SetFileName(_fileName);
336 writer->Write();
337 }
338
339 /* ---------------------------------------------------------------------- */
340
show()341 void DumpATOMVTK::vtkExportData::show() {
342 for (unsigned int i=0; i < vtkData.size(); i++) {
343 std::cerr << vtkData[i].serialize();
344 }
345 }
346
347 /* ---------------------------------------------------------------------- */
348
clear()349 void DumpATOMVTK::vtkExportData::clear() {
350 vtkData.clear();
351 }
352
353 /* ---------------------------------------------------------------------- */
354
setFileCurrent()355 void DumpATOMVTK::setFileCurrent() {
356 if (multifile == 0) filecurrent = filename;
357 else {
358 filecurrent = new char[strlen(filename) + 16];
359 char *ptr = strchr(filename,'*');
360 *ptr = '\0';
361 if (padflag == 0)
362 sprintf(filecurrent,"%s" BIGINT_FORMAT "%s",
363 filename,update->ntimestep,ptr+1);
364 else {
365 char bif[8],pad[16];
366 strcpy(bif,BIGINT_FORMAT);
367 sprintf(pad,"%%s%%0%d_%%d%s%%s",padflag,&bif[1]);
368 sprintf(filecurrent,pad,filename,comm->me,update->ntimestep,ptr+1);
369 }
370 *ptr = '*';
371 }
372 }
373
374 /* ---------------------------------------------------------------------- */
375
modify_param(int narg,char ** arg)376 int DumpATOMVTK::vtkExportData::modify_param(int narg, char **arg)
377 {
378 const int mvtk = DumpVTK::modify_param(narg, arg);
379 if (mvtk > 0)
380 return mvtk;
381
382 return 0;
383 }
384
385 /* ---------------------------------------------------------------------- */
386
modify_param(int narg,char ** arg)387 int DumpATOMVTK::modify_param(int narg, char **arg)
388 {
389 return tmpEXP.modify_param(narg, arg);
390 }
391
392 #endif
393