1 /*=========================================================================
2
3 Program: Visualization Toolkit
4 Module: vtkEnSightWriter.cxx
5
6 Copyright (c) Ken Martin, Will Schroeder, Bill Lorensen
7 All rights reserved.
8 See Copyright.txt or http://www.kitware.com/Copyright.htm for details.
9
10 This software is distributed WITHOUT ANY WARRANTY; without even
11 the implied warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
12 PURPOSE. See the above copyright notice for more information.
13
14 =========================================================================*/
15
16 /*----------------------------------------------------------------------------
17 Copyright (c) Sandia Corporation
18 See Copyright.txt or http://www.paraview.org/HTML/Copyright.html for details.
19 ----------------------------------------------------------------------------*/
20
21 /* TODO
22 *
23 *
24 * check that data was written
25 *
26 */
27
28 #include "vtkEnSightWriter.h"
29
30 #include "vtkToolkits.h" // for VTK_USE_PARALLEL
31 #ifdef VTK_USE_PARALLEL
32 # include "vtkMultiProcessController.h"
33 #endif
34
35 #include "vtkBitArray.h"
36 #include "vtkByteSwap.h"
37 #include "vtkCellArray.h"
38 #include "vtkCellData.h"
39 #include "vtkCharArray.h"
40 #include "vtkDataSet.h"
41 #include "vtkDoubleArray.h"
42 #include "vtkErrorCode.h"
43 #include "vtkFieldData.h"
44 #include "vtkFloatArray.h"
45 #include "vtkInformation.h"
46 #include "vtkIntArray.h"
47 #include "vtkLongArray.h"
48 #include "vtkLookupTable.h"
49 #include "vtkMath.h"
50 #include "vtkObjectFactory.h"
51 #include "vtkPointData.h"
52 #include "vtkPoints.h"
53 #include "vtkShortArray.h"
54 #include "vtkStreamingDemandDrivenPipeline.h"
55 #include "vtkUnsignedCharArray.h"
56 #include "vtkUnsignedIntArray.h"
57 #include "vtkUnsignedLongArray.h"
58 #include "vtkUnsignedShortArray.h"
59
60 #include "vtkPriorityQueue.h"
61
62 #include "vtkUnstructuredGrid.h"
63
64 #include <cerrno>
65 #include <cmath>
66 #include <cctype>
67
68 #include <vector>
69 #include <list>
70 #include <map>
71 #include <algorithm>
72
73 // If we are building against a slightly older VTK version,
74 // these cell types are not defined, and won't occur in the input
75
76 #ifndef VTK_QUADRATIC_WEDGE
77 # define VTK_QUADRATIC_WEDGE 26
78 # define VTK_QUADRATIC_PYRAMID 27
79 #endif
80
81 // this undef is required on the hp. vtkMutexLock ends up including
82 // /usr/include/dce/cma_ux.h which has the gall to #define write as cma_write
83
84 #ifdef write
85 # undef write
86 #endif
87
88 //----------------------------------------------------------------------------
89 vtkStandardNewMacro(vtkEnSightWriter);
90
91 //----------------------------------------------------------------------------
92 // Created object with no filename and timestep 0
vtkEnSightWriter()93 vtkEnSightWriter::vtkEnSightWriter()
94 {
95
96 this->BaseName = nullptr;
97 this->FileName = nullptr;
98 this->TimeStep = 0;
99 this->Path=nullptr;
100 this->GhostLevelMultiplier=10000;
101 this->GhostLevel = 0;
102 this->TransientGeometry=false;
103 this->ProcessNumber=0;
104 this->NumberOfProcesses=1;
105 this->NumberOfBlocks=0;
106 this->BlockIDs=nullptr;
107 this->TmpInput = nullptr;
108 }
109
110 //----------------------------------------------------------------------------
~vtkEnSightWriter()111 vtkEnSightWriter::~vtkEnSightWriter()
112 {
113 this->SetBaseName(nullptr);
114 this->SetFileName(nullptr);
115 this->SetPath(nullptr);
116 }
117
118 //----------------------------------------------------------------------------
PrintSelf(ostream & os,vtkIndent indent)119 void vtkEnSightWriter::PrintSelf(ostream& os, vtkIndent indent)
120 {
121 this->Superclass::PrintSelf(os,indent);
122
123 os << indent << "FileName: "
124 << (this->FileName ? this->FileName : "(none)") << "\n";
125 os << indent << "Path: "
126 << (this->Path ? this->Path : "(none)") << "\n";
127 os << indent << "BaseName: "
128 << (this->BaseName ? this->BaseName : "(none)") << "\n";
129
130 os << indent << "TimeStep: " << this->TimeStep << "\n";
131 os << indent << "TransientGeometry: " << this->TransientGeometry << "\n";
132 os << indent << "ProcessNumber: " << this->ProcessNumber << endl;
133 os << indent << "NumberOfProcesses: " << this->NumberOfProcesses << endl;
134 os << indent << "NumberOfBlocks: " << this->NumberOfBlocks << endl;
135 os << indent << "BlockIDs: " << this->BlockIDs << endl;
136 os << indent << "GhostLevel: " << this->GhostLevel << endl;
137 }
138
FillInputPortInformation(int vtkNotUsed (port),vtkInformation * info)139 int vtkEnSightWriter::FillInputPortInformation( int vtkNotUsed(port), vtkInformation* info )
140 {
141 info->Set(vtkAlgorithm::INPUT_REQUIRED_DATA_TYPE(), "vtkUnstructuredGrid");
142 return 1;
143 }
144
145 //----------------------------------------------------------------------------
146 // Specify the input data or filter.
SetInputData(vtkUnstructuredGrid * input)147 void vtkEnSightWriter::SetInputData(vtkUnstructuredGrid *input)
148 {
149 this->SetInputDataInternal(0, input);
150 }
151
152
153 //----------------------------------------------------------------------------
154 // Specify the input data or filter.
GetInput()155 vtkUnstructuredGrid *vtkEnSightWriter::GetInput()
156 {
157 if (this->GetNumberOfInputConnections(0) < 1)
158 {
159 return nullptr;
160 }
161 else if (this->TmpInput)
162 {
163 return this->TmpInput;
164 }
165 else
166 {
167 return (vtkUnstructuredGrid *)(this->Superclass::GetInput());
168 }
169 }
170
171 //----------------------------------------------------------------------------
WriteData()172 void vtkEnSightWriter::WriteData()
173 {
174 int i;
175 unsigned int ui;
176 int blockCount=0;
177 std::list<int>::iterator iter;
178
179 if (this->TmpInput)
180 {
181 this->TmpInput->Delete();
182 this->TmpInput = nullptr;
183 }
184
185 //figure out process ID
186
187 this->ProcessNumber = 0;
188 this->NumberOfProcesses = 1;
189
190 #ifdef VTK_USE_PARALLEL
191 vtkMultiProcessController *c = vtkMultiProcessController::GetGlobalController();
192
193 if (c != nullptr)
194 {
195 this->ProcessNumber=c->GetLocalProcessId();
196 this->NumberOfProcesses = c->GetNumberOfProcesses();
197 }
198 #endif
199
200 vtkUnstructuredGrid *input=this->GetInput();
201 vtkInformation* inInfo = this->GetInputInformation();
202
203 if (this->GhostLevel >
204 inInfo->Get(vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS()))
205 {
206 // re-execute pipeline if necessary to obtain ghost cells
207
208 this->GetInputAlgorithm()->UpdateInformation();
209 inInfo->Set(
210 vtkStreamingDemandDrivenPipeline::UPDATE_NUMBER_OF_GHOST_LEVELS(),
211 this->GhostLevel);
212 this->GetInputAlgorithm()->Update();
213 }
214
215 //get the BlockID Cell Array
216 vtkDataArray *BlockData=input->GetCellData()->GetScalars("BlockId");
217
218 if (BlockData==nullptr || strcmp(BlockData->GetName(),"BlockId"))
219 {
220 BlockData=nullptr;
221 }
222
223 this->ComputeNames();
224
225 if (!this->BaseName)
226 {
227 vtkErrorMacro("A FileName or Path/BaseName must be specified.");
228 return;
229 }
230
231 this->SanitizeFileName(this->BaseName);
232
233 char** blockNames=nullptr;
234 int * elementIDs=nullptr;
235 char charBuffer[1024];
236 char fileBuffer[512];
237 snprintf(charBuffer,sizeof(charBuffer),"%s/%s.%d.%05d.geo",
238 this->Path,this->BaseName,this->ProcessNumber,
239 this->TimeStep);
240
241 //open the geometry file
242 //only if timestep 0 and not transient geometry or transient geometry
243 FILE *fd=nullptr;
244 if (this->ShouldWriteGeometry())
245 {
246 if (!(fd=OpenFile(charBuffer)))
247 return;
248 }
249
250 //Get the FILE's for Point Data Fields
251 std::vector<FILE*> pointArrayFiles;
252 int NumPointArrays=input->GetPointData()->GetNumberOfArrays();
253 for (i=0;i<NumPointArrays;i++)
254 {
255 strcpy(fileBuffer,input->GetPointData()->GetArray(i)->GetName());
256 this->SanitizeFileName(fileBuffer);
257 snprintf(charBuffer,sizeof(charBuffer),"%s/%s.%d.%05d_n.%s",this->Path,this->BaseName,this->ProcessNumber,
258 this->TimeStep,fileBuffer);
259 FILE* ftemp=OpenFile(charBuffer);
260 if (!ftemp)
261 {
262 fclose(fd);
263 return;
264 }
265 pointArrayFiles.push_back(ftemp);
266
267 //write the description line to the file
268 WriteStringToFile(fileBuffer,ftemp);
269 }
270
271 //Get the FILE's for Cell Data Fields
272 std::vector<FILE*> cellArrayFiles;
273 int NumCellArrays=input->GetCellData()->GetNumberOfArrays();
274 for (i=0;i<NumCellArrays;i++)
275 {
276
277 strcpy(fileBuffer,input->GetCellData()->GetArray(i)->GetName());
278 this->SanitizeFileName(fileBuffer);
279 snprintf(charBuffer,sizeof(charBuffer),"%s/%s.%d.%05d_c.%s",this->Path,this->BaseName,this->ProcessNumber,
280 this->TimeStep,fileBuffer);
281 FILE* ftemp=OpenFile(charBuffer);
282 if (!ftemp)
283 {
284 fclose(fd);
285 return;
286 }
287 cellArrayFiles.push_back(ftemp);
288
289 //write the description line to the file
290 WriteStringToFile(fileBuffer,ftemp);
291 }
292
293 //write the header information
294 if (this->ShouldWriteGeometry())
295 {
296 this->WriteStringToFile("C Binary",fd);
297 this->WriteStringToFile("Written by VTK EnSight Writer",fd);
298 //if (this->Title)
299 // this->WriteStringToFile(this->Title,fd);
300 //else
301 this->WriteStringToFile("No Title was Specified",fd);
302 //we will specify node and element ID's
303 this->WriteStringToFile("node id given\n",fd);
304 this->WriteStringToFile("element id given\n",fd);
305 }
306
307 //get the Ghost Cell Array if it exists
308 vtkDataArray *GhostData=input->GetCellData()->GetScalars(vtkDataSetAttributes::GhostArrayName());
309 //if the strings are not the same then we did not get the ghostData array
310 if (GhostData==nullptr || strcmp(GhostData->GetName(), vtkDataSetAttributes::GhostArrayName()))
311 {
312 GhostData=nullptr;
313 }
314
315
316 //data structure to get all the cells for a certain part
317 //basically sort by part# and cell type
318 std::map <int, std::vector <int> > CellsByPart;
319
320 //just a list of part numbers
321 std::list<int> partNumbers;
322
323 //get all the part numbers in the unstructured grid and sort the cells
324 //by part number
325 for (i=0;i<input->GetNumberOfCells();i++)
326 {
327 int key=1;
328 if (BlockData)
329 key=(int)(BlockData->GetTuple(i)[0]);
330 else
331 cout << "No BlockID was found\n";
332 if (CellsByPart.count(key)==0)
333 {
334 CellsByPart[key]=std::vector < int >() ;
335 }
336 CellsByPart[key].push_back(i);
337 partNumbers.push_back(key);
338
339 }
340
341 //remove the duplicates from the partNumbers
342 partNumbers.sort();
343 partNumbers.unique();
344
345 //write out each part
346 for (iter=partNumbers.begin();iter!=partNumbers.end();++iter)
347 {
348 unsigned int j;
349 std::list<int>::iterator iter2;
350 int part=*iter;
351
352 //write the part Header
353 if (this->ShouldWriteGeometry())
354 {
355 blockCount+=1;
356 this->WriteStringToFile("part",fd);
357 this->WriteIntToFile(part,fd);
358 //cout << "part is " << part << endl;
359 this->WriteStringToFile("VTK Part",fd);
360 this->WriteStringToFile("coordinates",fd);
361 }
362
363 //write the part header for data files
364 for (j=0;j<pointArrayFiles.size();j++)
365 {
366 this->WriteStringToFile("part",pointArrayFiles[j]);
367 this->WriteIntToFile(part,pointArrayFiles[j]);
368 this->WriteStringToFile("coordinates",pointArrayFiles[j]);
369 }
370 for (j=0;j<cellArrayFiles.size();j++)
371 {
372 this->WriteStringToFile("part",cellArrayFiles[j]);
373 this->WriteIntToFile(part,cellArrayFiles[j]);
374 }
375
376 //list of VTK Node Indices per part
377 std::list<int> NodesPerPart;
378
379 //map that goes from NodeID to the order, used for element connectivity
380 std::map<int, int> NodeIdToOrder;
381
382 //get a list of all the nodes used for a particular part
383 for (j=0;j<CellsByPart[part].size();j++)
384 {
385 vtkCell* cell=input->GetCell(CellsByPart[part][j]);
386 vtkIdList* points=cell->GetPointIds();
387 for (int k=0;k<points->GetNumberOfIds();k++)
388 {
389 NodesPerPart.push_back(points->GetId(k));
390 }
391 }
392
393 //remove the duplicate Node ID's
394 NodesPerPart.sort();
395 NodesPerPart.unique();
396
397 //write the number of nodes
398 if (this->ShouldWriteGeometry())
399 {
400 this->WriteIntToFile(static_cast<int>(NodesPerPart.size()),fd);
401
402
403 //write the Node ID's to the file
404 //also set up the NodeID->order map
405 int NodeCount=0;
406 for (iter2=NodesPerPart.begin();iter2!=NodesPerPart.end();++iter2)
407 {
408 this->WriteIntToFile(*iter2,fd);
409 NodeIdToOrder[*iter2]=NodeCount+1;
410 NodeCount++;
411 }
412
413
414 //EnSight format requires all the X's, then all the Y's, then all the Z's
415 //write the X Coordinates
416
417 vtkPoints* inputPoints=input->GetPoints();
418 for (iter2=NodesPerPart.begin();iter2!=NodesPerPart.end();++iter2)
419 {
420 this->WriteFloatToFile((float)(inputPoints->GetPoint(*iter2)[0]),fd);
421 }
422 for (iter2=NodesPerPart.begin();iter2!=NodesPerPart.end();++iter2)
423 {
424 this->WriteFloatToFile((float)(inputPoints->GetPoint(*iter2)[1]),fd);
425 }
426 for (iter2=NodesPerPart.begin();iter2!=NodesPerPart.end();++iter2)
427 {
428 this->WriteFloatToFile((float)(inputPoints->GetPoint(*iter2)[2]),fd);
429 }
430
431 }
432
433 //write the Node Data for this part
434 for (j=0;j<pointArrayFiles.size();j++)
435 {
436 vtkDataArray* DataArray=input->GetPointData()->GetArray(j);
437 //figure out what type of data it is
438 int DataSize=DataArray->GetNumberOfComponents();
439
440 for (int CurrentDimension=0;
441 CurrentDimension<DataSize;
442 CurrentDimension++)
443 {
444 for (std::list<int>::iterator k=NodesPerPart.begin();
445 k!=NodesPerPart.end();++k)
446 {
447 this->WriteFloatToFile((float)
448 (DataArray->
449 GetTuple(*k)[CurrentDimension]),
450 pointArrayFiles[j]);
451 }
452 }
453 }
454
455
456
457 //now we need to sort the cell list by element type
458 //map is indexed by cell type has a vector of cell ID's
459 std::map<int, std::vector<int> > CellsByElement;
460 for (j=0;j<CellsByPart[part].size();j++)
461 {
462 int CellType=input->GetCell(CellsByPart[part][j])->GetCellType();
463 int ghostLevel=0;
464 if (GhostData)
465 {
466 ghostLevel=(int)(GhostData->GetTuple(CellsByPart[part][j])[0]);
467 if (ghostLevel & vtkDataSetAttributes::DUPLICATECELL)
468 {
469 ghostLevel=1;
470 }
471 }
472 //we want to sort out the ghost cells from the normal cells
473 //so the element type will be ghostMultiplier*ghostLevel+elementType
474 CellType+=ghostLevel*GhostLevelMultiplier;
475 if (CellsByElement.count(CellType)==0)
476 {
477 CellsByElement[CellType]=std::vector < int >() ;
478 }
479 CellsByElement[CellType].push_back(CellsByPart[part][j]);
480 }
481
482 //now we need to go through each element type that EnSight understands
483 std::vector<int> elementTypes;
484
485 //list the types that EnSight understands
486 //the noticeable absences are the ones without a fixed number of Nodes
487 //for the ghost cell types
488 elementTypes.push_back(VTK_VERTEX);
489 elementTypes.push_back(VTK_LINE);
490 elementTypes.push_back(VTK_TRIANGLE);
491 elementTypes.push_back(VTK_QUAD);
492 elementTypes.push_back(VTK_POLYGON);
493 elementTypes.push_back(VTK_TETRA);
494 elementTypes.push_back(VTK_HEXAHEDRON);
495 elementTypes.push_back(VTK_WEDGE);
496 elementTypes.push_back(VTK_PYRAMID);
497 elementTypes.push_back(VTK_CONVEX_POINT_SET);
498 elementTypes.push_back(VTK_QUADRATIC_EDGE);
499 elementTypes.push_back(VTK_QUADRATIC_TRIANGLE);
500 elementTypes.push_back(VTK_QUADRATIC_QUAD);
501 elementTypes.push_back(VTK_QUADRATIC_TETRA);
502 elementTypes.push_back(VTK_QUADRATIC_HEXAHEDRON);
503 elementTypes.push_back(VTK_QUADRATIC_WEDGE);
504 elementTypes.push_back(VTK_QUADRATIC_PYRAMID);
505 elementTypes.push_back(this->GhostLevelMultiplier+VTK_VERTEX);
506 elementTypes.push_back(this->GhostLevelMultiplier+VTK_LINE);
507 elementTypes.push_back(this->GhostLevelMultiplier+VTK_TRIANGLE);
508 elementTypes.push_back(this->GhostLevelMultiplier+VTK_QUAD);
509 elementTypes.push_back(this->GhostLevelMultiplier+VTK_POLYGON);
510 elementTypes.push_back(this->GhostLevelMultiplier+VTK_TETRA);
511 elementTypes.push_back(this->GhostLevelMultiplier+VTK_HEXAHEDRON);
512 elementTypes.push_back(this->GhostLevelMultiplier+VTK_WEDGE);
513 elementTypes.push_back(this->GhostLevelMultiplier+VTK_PYRAMID);
514 elementTypes.push_back(this->GhostLevelMultiplier+VTK_CONVEX_POINT_SET);
515 elementTypes.push_back(this->GhostLevelMultiplier+VTK_QUADRATIC_EDGE);
516 elementTypes.push_back(this->GhostLevelMultiplier+VTK_QUADRATIC_TRIANGLE);
517 elementTypes.push_back(this->GhostLevelMultiplier+VTK_QUADRATIC_QUAD);
518 elementTypes.push_back(this->GhostLevelMultiplier+VTK_QUADRATIC_TETRA);
519 elementTypes.push_back(this->GhostLevelMultiplier+VTK_QUADRATIC_HEXAHEDRON);
520 elementTypes.push_back(this->GhostLevelMultiplier+VTK_QUADRATIC_WEDGE);
521 elementTypes.push_back(this->GhostLevelMultiplier+VTK_QUADRATIC_PYRAMID);
522
523 //write out each type of element
524 if (this->ShouldWriteGeometry())
525 {
526 for (j=0;j<elementTypes.size();j++)
527 {
528 unsigned int k;
529 int elementType=elementTypes[j];
530 if (CellsByElement.count(elementType)>0)
531 {
532 //switch on element type to write correct type to file
533 this->WriteElementTypeToFile(elementType,fd);
534
535 //number of elements
536 this->WriteIntToFile(
537 static_cast<int>(CellsByElement[elementType].size()),fd);
538
539 //element ID's
540 for (k=0;k<CellsByElement[elementType].size();k++)
541 {
542 int CellId=CellsByElement[elementType][k];
543 this->WriteIntToFile(CellId,fd);
544 }
545
546 //element conenctivity information
547 for (k=0;k<CellsByElement[elementType].size();k++)
548 {
549 int CellId=CellsByElement[elementType][k];
550 vtkIdList *PointIds=input->GetCell(CellId)->GetPointIds();
551 for (int m=0;m<PointIds->GetNumberOfIds();m++)
552 {
553 int PointId=PointIds->GetId(m);
554 this->WriteIntToFile(NodeIdToOrder[PointId],fd);
555 }
556 }
557 }
558 }
559 }
560
561 //write the Cell Data for this part
562 for (j=0;j<cellArrayFiles.size();j++)
563 {
564 vtkDataArray* DataArray=input->GetCellData()->GetArray(j);
565 //figure out what type of data it is
566 int DataSize=DataArray->GetNumberOfComponents();
567
568 for (unsigned int k=0;k<elementTypes.size();k++)
569 {
570 if (!CellsByElement[elementTypes[k]].empty())
571 {
572 this->WriteElementTypeToFile(elementTypes[k],
573 cellArrayFiles[j]);
574 for (unsigned int m=0;m<CellsByElement[elementTypes[k]].size();m++)
575 {
576 for ( int CurrentDimension = 0; CurrentDimension < DataSize;
577 CurrentDimension ++ )
578 {
579 this->WriteFloatToFile
580 ( (float)
581 ( DataArray
582 ->GetTuple( CellsByElement[ elementTypes[k] ][m] )
583 [CurrentDimension]
584 ),
585 cellArrayFiles[j]
586 );
587 }
588 }
589 }
590 }
591 }
592 }
593
594 //now write the empty blocks
595 //use the block list in the exodus model if it exists, otherwise
596 //use the BlockID list if that exists.
597
598 if (this->BlockIDs)
599 {
600 elementIDs=this->BlockIDs;
601 }
602
603 if (elementIDs)
604 {
605 //cout << "have " << this->NumberOfBlocks << " blocks " << endl;
606 for (i=0;i<this->NumberOfBlocks;i++)
607 {
608 unsigned int j;
609 //figure out if the part was already written
610 int part=elementIDs[i];
611 if ( std::find(partNumbers.begin(), partNumbers.end(), part)
612 == partNumbers.end() )
613 {
614 //no information about the part was written to the output files
615 //so write some empty information
616 if (this->ShouldWriteGeometry())
617 {
618 blockCount+=1;
619 this->WriteStringToFile("part",fd);
620 this->WriteIntToFile(part,fd);
621
622 int exodusIndex =
623 this->GetExodusModelIndex(elementIDs,this->NumberOfBlocks,part);
624
625 if (exodusIndex!=-1 && blockNames)
626 {
627 snprintf(charBuffer,sizeof(charBuffer),"Exodus-%s-%d",blockNames[exodusIndex],part);
628 this->WriteStringToFile(charBuffer,fd);
629 }
630 else
631 {
632 this->WriteStringToFile("VTK Part",fd);
633 }
634 }
635
636 //write the part header for data files
637 for (j=0;j<pointArrayFiles.size();j++)
638 {
639 this->WriteStringToFile("part",pointArrayFiles[j]);
640 this->WriteIntToFile(part,pointArrayFiles[j]);
641 }
642 for (j=0;j<cellArrayFiles.size();j++)
643 {
644 this->WriteStringToFile("part",cellArrayFiles[j]);
645 this->WriteIntToFile(part,cellArrayFiles[j]);
646 }
647 }
648 }
649 }
650 //cout << "wrote " << blockCount << "parts\n";
651 if (this->TmpInput)
652 {
653 this->TmpInput->Delete();
654 this->TmpInput = nullptr;
655 }
656
657 //close all the files
658 if (fd)
659 {
660 fclose(fd);
661 }
662
663 for (ui=0;ui<cellArrayFiles.size();ui++)
664 {
665 fclose(cellArrayFiles[ui]);
666 }
667 for (ui=0;ui<pointArrayFiles.size();ui++)
668 {
669 fclose(pointArrayFiles[ui]);
670 }
671
672 }
673
674 //----------------------------------------------------------------------------
WriteCaseFile(int TotalTimeSteps)675 void vtkEnSightWriter::WriteCaseFile(int TotalTimeSteps)
676 {
677
678 vtkUnstructuredGrid *input=this->GetInput();
679 int i;
680
681 this->ComputeNames();
682
683 if (!this->BaseName)
684 {
685 vtkErrorMacro("A FileName or Path/BaseName must be specified.");
686 return;
687 }
688
689 char charBuffer[1024];
690 snprintf(charBuffer,sizeof(charBuffer),"%s/%s.%d.case",this->Path,this->BaseName,this->ProcessNumber);
691
692 //open the geometry file
693 FILE *fd=nullptr;
694 if (!(fd=OpenFile(charBuffer)))
695 {
696 return;
697 }
698
699 this->WriteTerminatedStringToFile("FORMAT\n",fd);
700 this->WriteTerminatedStringToFile("type: ensight gold\n\n",fd);
701 this->WriteTerminatedStringToFile("\nGEOMETRY\n",fd);
702
703 //write the geometry file
704 if (!this->TransientGeometry)
705 {
706 snprintf(charBuffer,sizeof(charBuffer),"model: %s.%d.00000.geo\n",this->BaseName,this->ProcessNumber);
707 this->WriteTerminatedStringToFile(charBuffer,fd);
708 }
709 else
710 {
711 snprintf(charBuffer,sizeof(charBuffer),"model: 1 %s.%d.*****.geo\n",this->BaseName,this->ProcessNumber);
712 this->WriteTerminatedStringToFile(charBuffer,fd);
713 }
714
715 this->WriteTerminatedStringToFile("\nVARIABLE\n",fd);
716
717 char fileBuffer[256];
718
719 //write the Node variable files
720 for (i=0;i<input->GetPointData()->GetNumberOfArrays();i++)
721 {
722
723 strcpy(fileBuffer,input->GetPointData()->GetArray(i)->GetName());
724 // skip arrays that were not written
725 if (strcmp(fileBuffer,"GlobalElementId")==0)
726 {
727 continue;
728 }
729 if (strcmp(fileBuffer,"GlobalNodeId")==0)
730 {
731 continue;
732 }
733 if (strcmp(fileBuffer,"BlockId")==0)
734 {
735 continue;
736 }
737 this->SanitizeFileName(fileBuffer);
738 //figure out what kind of data it is
739 char SmallBuffer[16];
740 switch(input->GetPointData()->GetArray(i)->GetNumberOfComponents())
741 {
742 case(1):
743 strcpy(SmallBuffer,"scalar");
744 break;
745 case(3):
746 strcpy(SmallBuffer,"vector");
747 break;
748 case(6):
749 strcpy(SmallBuffer,"tensor");
750 break;
751 case(9):
752 strcpy(SmallBuffer,"tensor9");
753 break;
754 }
755 if (TotalTimeSteps<=1)
756 {
757 snprintf(charBuffer,sizeof(charBuffer),"%s per node: %s_n %s.%d.00000_n.%s\n",
758 SmallBuffer,
759 fileBuffer,
760 this->BaseName,
761 this->ProcessNumber,
762 fileBuffer);
763 }
764 else
765 {
766 snprintf(charBuffer,sizeof(charBuffer),"%s per node: 1 %s_n %s.%d.*****_n.%s\n",
767 SmallBuffer,
768 fileBuffer,
769 this->BaseName,
770 this->ProcessNumber,
771 fileBuffer);
772
773
774 }
775 this->WriteTerminatedStringToFile(charBuffer,fd);
776 }
777
778 //write the cell variable files
779 for (i=0;i<input->GetCellData()->GetNumberOfArrays();i++)
780 {
781 //figure out what kind of data it is
782 char SmallBuffer[16];
783
784 strcpy(fileBuffer,input->GetCellData()->GetArray(i)->GetName());
785 // skip arrays that were not written
786 if (strcmp(fileBuffer,"GlobalElementId")==0)
787 {
788 continue;
789 }
790 if (strcmp(fileBuffer,"GlobalNodeId")==0)
791 {
792 continue;
793 }
794 if (strcmp(fileBuffer,"BlockId")==0)
795 {
796 continue;
797 }
798 this->SanitizeFileName(fileBuffer);
799 switch(input->GetCellData()->GetArray(i)->GetNumberOfComponents())
800 {
801 case(1):
802 strcpy(SmallBuffer,"scalar");
803 break;
804 case(3):
805 strcpy(SmallBuffer,"vector");
806 break;
807 case(6):
808 strcpy(SmallBuffer,"tensor");
809 break;
810 case(9):
811 strcpy(SmallBuffer,"tensor9");
812 break;
813 }
814 if (TotalTimeSteps<=1)
815 {
816 snprintf(charBuffer,sizeof(charBuffer),"%s per element: %s_c %s.%d.00000_c.%s\n",
817 SmallBuffer,
818 fileBuffer,
819 this->BaseName,
820 this->ProcessNumber,
821 fileBuffer);
822 }
823 else
824 {
825 snprintf(charBuffer,sizeof(charBuffer),"%s per element: 1 %s_c %s.%d.*****_c.%s\n",
826 SmallBuffer,
827 fileBuffer,
828 this->BaseName,
829 this->ProcessNumber,
830 fileBuffer);
831 }
832 this->WriteTerminatedStringToFile(charBuffer,fd);
833 }
834
835 //write time information if we have multiple timesteps
836 if (TotalTimeSteps>1)
837 {
838 this->WriteTerminatedStringToFile("\nTIME\n",fd);
839 this->WriteTerminatedStringToFile("time set: 1\n",fd);
840 snprintf(charBuffer,sizeof(charBuffer),"number of steps: %d\n",TotalTimeSteps);
841 this->WriteTerminatedStringToFile(charBuffer,fd);
842 this->WriteTerminatedStringToFile("filename start number: 00000\n",fd);
843 this->WriteTerminatedStringToFile("filename increment: 00001\n",fd);
844 this->WriteTerminatedStringToFile("time values: \n",fd);
845 for (i=0;i<TotalTimeSteps;i++)
846 {
847 double timestep=i;
848
849 snprintf(charBuffer,sizeof(charBuffer),"%f ",timestep);
850 this->WriteTerminatedStringToFile(charBuffer,fd);
851 if (i%6==0 && i>0)
852 {
853 this->WriteTerminatedStringToFile("\n",fd);
854 }
855 }
856 }
857
858
859 }
860
861 //----------------------------------------------------------------------------
WriteSOSCaseFile(int numProcs)862 void vtkEnSightWriter::WriteSOSCaseFile(int numProcs)
863 {
864 this->ComputeNames();
865
866 if (!this->BaseName)
867 {
868 vtkErrorMacro("A FileName or Path/BaseName must be specified.");
869 return;
870 }
871
872 this->SanitizeFileName(this->BaseName);
873
874 char charBuffer[512];
875 snprintf(charBuffer,sizeof(charBuffer),"%s/%s.case.sos",this->Path,this->BaseName);
876
877 FILE *fd=nullptr;
878 if (!(fd=OpenFile(charBuffer)))
879 return;
880
881 this->WriteTerminatedStringToFile("FORMAT\n",fd);
882 this->WriteTerminatedStringToFile("type: master_server gold\n\n",fd);
883
884 this->WriteTerminatedStringToFile("SERVERS\n",fd);
885 snprintf(charBuffer,sizeof(charBuffer),"number of servers: %d\n\n",numProcs);
886 this->WriteTerminatedStringToFile(charBuffer,fd);
887
888 //write the servers section with placeholders for the ensight server
889 //location and server name
890 int i=0;
891 for (i=0;i<numProcs;i++)
892 {
893 snprintf(charBuffer,sizeof(charBuffer), "#Server %d\n",i);
894 this->WriteTerminatedStringToFile(charBuffer,fd);
895 this->WriteTerminatedStringToFile("#-------\n",fd);
896 snprintf(charBuffer,sizeof(charBuffer), "machine id: MID%05d\n",i);
897 this->WriteTerminatedStringToFile(charBuffer,fd);
898
899 this->WriteTerminatedStringToFile("executable: MEX\n",fd);
900 snprintf(charBuffer,sizeof(charBuffer), "data_path: %s\n",this->Path);
901 this->WriteTerminatedStringToFile(charBuffer,fd);
902 snprintf(charBuffer,sizeof(charBuffer),"casefile: %s.%d.case\n\n",this->BaseName,i);
903 this->WriteTerminatedStringToFile(charBuffer,fd);
904 }
905
906
907 }
908
909 //----------------------------------------------------------------------------
WriteStringToFile(const char * cstring,FILE * file)910 void vtkEnSightWriter::WriteStringToFile(const char* cstring, FILE* file)
911 {
912 char cbuffer[81];
913 unsigned long cstringLength = static_cast<unsigned long>(strlen(cstring));
914 memcpy(cbuffer,cstring,vtkMath::Min(cstringLength,80ul));
915 for (int i = cstringLength; i < 81; ++i)
916 {
917 cbuffer[i] = '\0';
918 }
919
920 // Write a constant 80 bytes to the file.
921 fwrite(cbuffer, sizeof(char),80,file);
922 }
923
924 //----------------------------------------------------------------------------
WriteTerminatedStringToFile(const char * cstring,FILE * file)925 void vtkEnSightWriter::WriteTerminatedStringToFile(const char* cstring, FILE* file)
926 {
927 fwrite(cstring, sizeof(char),std::min(strlen(cstring), static_cast<size_t>(512)),file);
928 }
929
930 //----------------------------------------------------------------------------
WriteIntToFile(const int i,FILE * file)931 void vtkEnSightWriter::WriteIntToFile(const int i,FILE* file)
932 {
933 fwrite(&i, sizeof(int),1,file);
934 }
935
936 //----------------------------------------------------------------------------
WriteFloatToFile(const float f,FILE * file)937 void vtkEnSightWriter::WriteFloatToFile(const float f,FILE* file)
938 {
939 fwrite(&f, sizeof(float),1,file);
940 }
941
942 //----------------------------------------------------------------------------
WriteElementTypeToFile(int elementType,FILE * fd)943 void vtkEnSightWriter::WriteElementTypeToFile(int elementType,FILE* fd)
944 {
945 int ghostLevel=elementType/GhostLevelMultiplier;
946 elementType=elementType%GhostLevelMultiplier;
947 if (ghostLevel==0)
948 {
949 switch(elementType)
950 {
951 case(VTK_VERTEX):
952 this->WriteStringToFile("point",fd);
953 break;
954 case(VTK_LINE):
955 this->WriteStringToFile("bar2",fd);
956 break;
957 case(VTK_TRIANGLE):
958 this->WriteStringToFile("tria3",fd);
959 break;
960 case(VTK_QUAD):
961 this->WriteStringToFile("quad4",fd);
962 break;
963 case(VTK_POLYGON):
964 this->WriteStringToFile("nsided",fd);
965 break;
966 case(VTK_TETRA):
967 this->WriteStringToFile("tetra4",fd);
968 break;
969 case(VTK_HEXAHEDRON):
970 this->WriteStringToFile("hexa8",fd);
971 break;
972 case(VTK_WEDGE):
973 this->WriteStringToFile("penta6",fd);
974 break;
975 case(VTK_PYRAMID):
976 this->WriteStringToFile("pyramid5",fd);
977 break;
978 case(VTK_CONVEX_POINT_SET):
979 this->WriteStringToFile("nfaced",fd);
980 break;
981 case(VTK_QUADRATIC_EDGE):
982 this->WriteStringToFile("bar3",fd);
983 break;
984 case(VTK_QUADRATIC_TRIANGLE):
985 this->WriteStringToFile("tria6",fd);
986 break;
987 case(VTK_QUADRATIC_QUAD):
988 this->WriteStringToFile("quad8",fd);
989 break;
990 case(VTK_QUADRATIC_TETRA):
991 this->WriteStringToFile("tetra10",fd);
992 break;
993 case(VTK_QUADRATIC_HEXAHEDRON):
994 this->WriteStringToFile("hexa20",fd);
995 break;
996 case(VTK_QUADRATIC_WEDGE):
997 this->WriteStringToFile("penta15",fd);
998 break;
999 case(VTK_QUADRATIC_PYRAMID):
1000 this->WriteStringToFile("pyramid13",fd);
1001 break;
1002 }
1003 }
1004 else
1005 {
1006 switch(elementType)
1007 {
1008 case(VTK_VERTEX):
1009 this->WriteStringToFile("g_point",fd);
1010 break;
1011 case(VTK_LINE):
1012 this->WriteStringToFile("g_bar2",fd);
1013 break;
1014 case(VTK_TRIANGLE):
1015 this->WriteStringToFile("g_tria3",fd);
1016 break;
1017 case(VTK_QUAD):
1018 this->WriteStringToFile("g_quad4",fd);
1019 break;
1020 case(VTK_POLYGON):
1021 this->WriteStringToFile("g_nsided",fd);
1022 break;
1023 case(VTK_TETRA):
1024 this->WriteStringToFile("g_tetra4",fd);
1025 break;
1026 case(VTK_HEXAHEDRON):
1027 this->WriteStringToFile("g_hexa8",fd);
1028 break;
1029 case(VTK_WEDGE):
1030 this->WriteStringToFile("g_penta6",fd);
1031 break;
1032 case(VTK_PYRAMID):
1033 this->WriteStringToFile("g_pyramid5",fd);
1034 break;
1035 case(VTK_CONVEX_POINT_SET):
1036 this->WriteStringToFile("g_nfaced",fd);
1037 break;
1038 case(VTK_QUADRATIC_EDGE):
1039 this->WriteStringToFile("g_bar3",fd);
1040 break;
1041 case(VTK_QUADRATIC_TRIANGLE):
1042 this->WriteStringToFile("g_tria6",fd);
1043 break;
1044 case(VTK_QUADRATIC_QUAD):
1045 this->WriteStringToFile("g_quad8",fd);
1046 break;
1047 case(VTK_QUADRATIC_TETRA):
1048 this->WriteStringToFile("g_tetra10",fd);
1049 break;
1050 case(VTK_QUADRATIC_HEXAHEDRON):
1051 this->WriteStringToFile("g_hexa20",fd);
1052 break;
1053 case(VTK_QUADRATIC_WEDGE):
1054 this->WriteStringToFile("g_penta15",fd);
1055 break;
1056 case(VTK_QUADRATIC_PYRAMID):
1057 this->WriteStringToFile("g_pyramid13",fd);
1058 break;
1059 }
1060 }
1061 }
1062
1063 //----------------------------------------------------------------------------
ShouldWriteGeometry()1064 bool vtkEnSightWriter::ShouldWriteGeometry()
1065 {
1066 return (this->TransientGeometry || (this->TimeStep==0));
1067 }
1068
1069 //----------------------------------------------------------------------------
SanitizeFileName(char * name)1070 void vtkEnSightWriter::SanitizeFileName(char* name)
1071 {
1072
1073 char buffer[512];
1074 unsigned int i;
1075 int BufferPosition=0;
1076 for (i=0;i<strlen(name);i++)
1077 {
1078 if (name[i]!='/')
1079 {
1080 buffer[BufferPosition]=name[i];
1081 BufferPosition++;
1082 }
1083 }
1084 buffer[BufferPosition]=0;
1085 for (i=0;i<strlen(buffer);i++)
1086 {
1087 name[i]=buffer[i];
1088 }
1089 name[strlen(buffer)]=0;
1090
1091 }
1092
1093 //----------------------------------------------------------------------------
OpenFile(char * name)1094 FILE* vtkEnSightWriter::OpenFile(char* name)
1095 {
1096 FILE * fd=fopen(name,"wb");
1097
1098 if (fd == nullptr)
1099 {
1100 vtkErrorMacro("Error opening " << name
1101 << ": " << strerror(errno));
1102 return nullptr;
1103 }
1104 return fd;
1105 }
1106
1107 //----------------------------------------------------------------------------
GetExodusModelIndex(int * elementArray,int numberElements,int partID)1108 int vtkEnSightWriter::GetExodusModelIndex(int *elementArray,int numberElements,int partID)
1109 {
1110 int i;
1111 for (i=0;i<numberElements;i++)
1112 {
1113 if (elementArray[i]==partID)
1114 return i;
1115 }
1116 return -1;
1117 }
1118
1119 //----------------------------------------------------------------------------
DefaultNames()1120 void vtkEnSightWriter::DefaultNames()
1121 {
1122 char *path = new char [4];
1123 char *base = new char [20];
1124 strcpy(path, "./");
1125 strcpy(base, "EnSightWriter.out");
1126
1127 this->SetPath(path);
1128 this->SetBaseName(base);
1129 }
1130
1131 //----------------------------------------------------------------------------
ComputeNames()1132 void vtkEnSightWriter::ComputeNames()
1133 {
1134 if (this->Path && this->BaseName)
1135 {
1136 return;
1137 }
1138
1139 if (!this->FileName)
1140 {
1141 this->DefaultNames();
1142 return;
1143 }
1144
1145 // FileName = Path/BaseName.digits.digits
1146
1147 char *path = nullptr;
1148 char *base = nullptr;
1149
1150 char *f = this->FileName;
1151
1152 while (!isgraph(*f)) f++; // find first printable character
1153
1154 if (!*f)
1155 {
1156 // FileName is garbage
1157 DefaultNames();
1158 return;
1159 }
1160
1161 char *buf = new char [strlen(f) + 1];
1162 strcpy(buf, f);
1163
1164 char *slash = strrchr(buf, '/'); // final slash
1165
1166 if (slash)
1167 {
1168 *slash = 0;
1169 path = new char [strlen(buf) + 1];
1170 strcpy(path, buf);
1171 f = slash + 1;
1172 }
1173 else
1174 {
1175 path = new char [4];
1176 strcpy(path, "./");
1177
1178 f = buf;
1179 }
1180
1181 char *firstChar = f;
1182 while (*f && (*f != '.')) f++;
1183 *f = 0;
1184
1185 base = new char [strlen(firstChar) + 1];
1186 strcpy(base, firstChar);
1187
1188 this->SetPath(path);
1189 this->SetBaseName(base);
1190
1191 delete [] buf;
1192 }
1193