/[escript]/trunk/finley/src/Mesh_saveVTK.c
ViewVC logotype

Annotation of /trunk/finley/src/Mesh_saveVTK.c

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2421 - (hide annotations)
Thu May 14 02:29:06 2009 UTC (10 years, 8 months ago) by gross
File MIME type: text/plain
File size: 45983 byte(s)
vtk writer supports now meta data and meta data schema
1 jgs 110
2 ksteube 1312 /*******************************************************
3 ksteube 1811 *
4     * Copyright (c) 2003-2008 by University of Queensland
5     * Earth Systems Science Computational Center (ESSCC)
6     * http://www.uq.edu.au/esscc
7     *
8     * Primary Business: Queensland, Australia
9     * Licensed under the Open Software License version 3.0
10     * http://www.opensource.org/licenses/osl-3.0.php
11     *
12     *******************************************************/
13 ksteube 1312
14 caltinay 2141 /***************************************************************************/
15     /* Writes data and mesh in VTK XML format to a VTU file. */
16     /* Nodal data needs to be given on FINLEY_NODES or FINLEY_REDUCED_NODES */
17     /***************************************************************************/
18 ksteube 1811
19 jgs 110 #include "Mesh.h"
20 gross 1062 #include "Assemble.h"
21 caltinay 2141 #include "vtkCellType.h" /* copied from vtk source directory */
22 gross 1669 #include "paso/PasoUtil.h"
23 jgs 110
24 ksteube 1312 #define INT_FORMAT "%d "
25 phornby 2172 #define LEN_INT_FORMAT (unsigned int)(9+1)
26 ksteube 1312 #define INT_NEWLINE_FORMAT "%d\n"
27 caltinay 2141 #define SCALAR_FORMAT "%12.6e\n"
28     #define VECTOR_FORMAT "%12.6e %12.6e %12.6e\n"
29     #define TENSOR_FORMAT "%12.6e %12.6e %12.6e %12.6e %12.6e %12.6e %12.6e %12.6e %12.6e\n"
30 phornby 2172 #define LEN_TENSOR_FORMAT (unsigned int)(9*(12+1)+1)
31 ksteube 1312 #define NEWLINE "\n"
32 caltinay 2141 #define LEN_TMP_BUFFER LEN_TENSOR_FORMAT+(MAX_numNodes*LEN_INT_FORMAT+1)+2
33 phornby 2172 #define NCOMP_MAX (unsigned int)9
34 caltinay 2141
35     #define __STRCAT(dest, chunk, dest_in_use) \
36     do {\
37     strcpy(&dest[dest_in_use], chunk);\
38     dest_in_use += strlen(chunk);\
39     } while(0)
40    
41     #ifdef PASO_MPI
42     /* writes buffer to file catching the empty buffer case which causes problems
43     * with some MPI versions */
44 gross 2385 #define MPI_WRITE_ORDERED(BUF) \
45 caltinay 2141 do {\
46 gross 2385 int LLEN=0; \
47     LLEN=(int) strlen(BUF); \
48     if (LLEN==0) { strcpy(BUF, ""); LLEN=0; }\
49     MPI_File_write_ordered(mpi_fileHandle_p, BUF, LLEN, MPI_CHAR, &mpi_status);\
50 caltinay 2141 } while(0)
51    
52     /* writes buffer to file on master only */
53     #define MPI_RANK0_WRITE_SHARED(BUF) \
54     do {\
55 gross 2385 int LLEN=0; \
56 caltinay 2141 if (my_mpi_rank == 0) {\
57 gross 2385 LLEN=(int) strlen(BUF); \
58     if (LLEN==0) { strcpy(BUF,""); LLEN=0; }\
59     MPI_File_iwrite_shared(mpi_fileHandle_p, BUF, LLEN, MPI_CHAR, &mpi_req);\
60 caltinay 2141 MPI_Wait(&mpi_req, &mpi_status);\
61     }\
62     } while(0)
63    
64     /* For reference only. Investigation needed as to which values may improve
65     * performance */
66     #if 0
67     void create_MPIInfo(MPI_Info& info)
68     {
69     MPI_Info_create(&info);
70     MPI_Info_set(info, "access_style", "write_once, sequential");
71     MPI_Info_set(info, "collective_buffering", "true");
72     MPI_Info_set(info, "cb_block_size", "131072");
73     MPI_Info_set(info, "cb_buffer_size", "1048567");
74     MPI_Info_set(info, "cb_nodes", "8");
75     MPI_Info_set(info, "striping_factor", "16");
76     MPI_Info_set(info, "striping_unit", "424288");
77 dhawcroft 793 }
78 caltinay 2141 #endif
79 dhawcroft 793
80 caltinay 2141 #else
81    
82 gross 2385 #define MPI_WRITE_ORDERED(A)
83 caltinay 2141 #define MPI_RANK0_WRITE_SHARED(A)
84    
85     #endif /* PASO_MPI */
86    
87    
88     /* Returns one if the node given by coords and idx is within the quadrant
89     * indexed by q and if the element type is Rec9 or Hex27, zero otherwise */
90     int nodeInQuadrant(const double *coords, ElementTypeId type, int idx, int q)
91     {
92     #define INSIDE_1D(_X_,_C_,_R_) ( ABS((_X_)-(_C_)) <= (_R_) )
93     #define INSIDE_2D(_X_,_Y_,_CX_,_CY_,_R_) ( INSIDE_1D(_X_,_CX_,_R_) && INSIDE_1D(_Y_,_CY_,_R_))
94     #define INSIDE_3D(_X_,_Y_,_Z_,_CX_,_CY_,_CZ_,_R_) ( INSIDE_1D(_X_,_CX_,_R_) && INSIDE_1D(_Y_,_CY_,_R_) && INSIDE_1D(_Z_,_CZ_,_R_) )
95    
96     int ret;
97     if (type == Rec9) {
98     if (q==0)
99     ret = INSIDE_2D(coords[2*idx], coords[2*idx+1], 0.25, 0.25, 0.25);
100     else if (q==1)
101     ret = INSIDE_2D(coords[2*idx], coords[2*idx+1], 0.75, 0.25, 0.25);
102     else if (q==2)
103     ret = INSIDE_2D(coords[2*idx], coords[2*idx+1], 0.25, 0.75, 0.25);
104     else if (q==3)
105     ret = INSIDE_2D(coords[2*idx], coords[2*idx+1], 0.75, 0.75, 0.25);
106     else
107     ret = 0;
108     } else if (type == Hex27) {
109     if (q==0)
110     ret = INSIDE_3D(coords[3*idx], coords[3*idx+1], coords[3*idx+2],
111     0.25, 0.25, 0.25, 0.25);
112     else if (q==1)
113     ret = INSIDE_3D(coords[3*idx], coords[3*idx+1], coords[3*idx+2],
114     0.75, 0.25, 0.25, 0.25);
115     else if (q==2)
116     ret = INSIDE_3D(coords[3*idx], coords[3*idx+1], coords[3*idx+2],
117     0.25, 0.75, 0.25, 0.25);
118     else if (q==3)
119     ret = INSIDE_3D(coords[3*idx], coords[3*idx+1], coords[3*idx+2],
120     0.75, 0.75, 0.25, 0.25);
121     else if (q==4)
122     ret = INSIDE_3D(coords[3*idx], coords[3*idx+1], coords[3*idx+2],
123     0.25, 0.25, 0.75, 0.25);
124     else if (q==5)
125     ret = INSIDE_3D(coords[3*idx], coords[3*idx+1], coords[3*idx+2],
126     0.75, 0.25, 0.75, 0.25);
127     else if (q==6)
128     ret = INSIDE_3D(coords[3*idx], coords[3*idx+1], coords[3*idx+2],
129     0.25, 0.75, 0.75, 0.25);
130     else if (q==7)
131     ret = INSIDE_3D(coords[3*idx], coords[3*idx+1], coords[3*idx+2],
132     0.75, 0.75, 0.75, 0.25);
133     else
134     ret = 0;
135     } else {
136     ret = 1;
137     }
138     return ret;
139     }
140    
141     void Finley_Mesh_saveVTK(const char *filename_p,
142 ksteube 1312 Finley_Mesh *mesh_p,
143     const dim_t num_data,
144 caltinay 2141 char **names_p,
145 gross 2421 escriptDataC **data_pp,
146     const char* metadata,
147     const char*metadata_schema)
148 dhawcroft 793 {
149 caltinay 2141 #ifdef PASO_MPI
150     MPI_File mpi_fileHandle_p;
151     MPI_Status mpi_status;
152     MPI_Request mpi_req;
153     MPI_Info mpi_info = MPI_INFO_NULL;
154     #endif
155     Paso_MPI_rank my_mpi_rank;
156     FILE *fileHandle_p = NULL;
157     char errorMsg[LenErrorMsg_MAX], *txtBuffer;
158     char tmpBuffer[LEN_TMP_BUFFER];
159     size_t txtBufferSize, txtBufferInUse, maxNameLen;
160     double *quadNodes_p = NULL;
161     dim_t dataIdx, nDim;
162     dim_t numCells=0, globalNumCells=0, numVTKNodesPerElement=0;
163     dim_t myNumPoints=0, globalNumPoints=0;
164     dim_t shape, NN=0, numCellFactor=1, myNumCells=0;
165     bool_t *isCellCentered;
166     bool_t writeCellData=FALSE, writePointData=FALSE, hasReducedElements=FALSE;
167     index_t myFirstNode=0, myLastNode=0, *globalNodeIndex=NULL;
168     index_t myFirstCell=0, k;
169     int mpi_size, i, j, l;
170     int cellType=0, nodeType=FINLEY_NODES, elementType=FINLEY_UNKNOWN;
171     Finley_ElementFile *elements = NULL;
172     ElementTypeId typeId = NoType;
173 dhawcroft 793
174 caltinay 2141 const char *vtkHeader = \
175     "<?xml version=\"1.0\"?>\n" \
176 gross 2421 "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\"%s%s>\n%s%s" \
177 caltinay 2141 "<UnstructuredGrid>\n" \
178     "<Piece NumberOfPoints=\"%d\" NumberOfCells=\"%d\">\n" \
179     "<Points>\n" \
180     "<DataArray NumberOfComponents=\"%d\" type=\"Float64\" format=\"ascii\">\n";
181     char *vtkFooter = "</Piece>\n</UnstructuredGrid>\n</VTKFile>\n";
182     const char *tag_Float_DataArray="<DataArray Name=\"%s\" type=\"Float64\" NumberOfComponents=\"%d\" format=\"ascii\">\n";
183     char *tags_End_Points_and_Start_Conn = "</DataArray>\n</Points>\n<Cells>\n<DataArray Name=\"connectivity\" type=\"Int32\" format=\"ascii\">\n" ;
184     char *tags_End_Conn_and_Start_Offset = "</DataArray>\n<DataArray Name=\"offsets\" type=\"Int32\" format=\"ascii\">\n";
185     char *tags_End_Offset_and_Start_Type = "</DataArray>\n<DataArray Name=\"types\" type=\"UInt8\" format=\"ascii\">\n";
186     char *tag_End_DataArray = "</DataArray>\n";
187 dhawcroft 793
188 caltinay 2141 const int VTK_HEX20_INDEX[] =
189     { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 16, 17, 18, 19, 12, 13, 14, 15 };
190     const int VTK_REC9_INDEX[] =
191     { 0, 4, 8, 7, 4, 1, 5, 8, 7, 8, 6, 3, 8, 5, 2, 6 };
192     const int VTK_HEX27_INDEX[] =
193     { 0, 8, 20, 11, 12, 21, 26, 24,
194     8, 1, 9, 20, 21, 13, 22, 26,
195     11, 20, 10, 3, 24, 26, 23, 15,
196     20, 9, 2, 10, 26, 22, 14, 23,
197     12, 21, 26, 24, 4, 16, 25, 19,
198     21, 13, 22, 26, 16, 5, 17, 25,
199     24, 26, 23, 15, 19, 25, 18, 7,
200     26, 22, 14, 23, 25, 17, 6, 18 };
201 caltinay 1985
202 caltinay 2141 /* if there is no mesh we just return */
203     if (mesh_p==NULL) return;
204 dhawcroft 793
205 caltinay 2141 nDim = mesh_p->Nodes->numDim;
206 dhawcroft 796
207 caltinay 2141 if (nDim != 2 && nDim != 3) {
208     Finley_setError(TYPE_ERROR, "saveVTK: spatial dimension 2 or 3 is supported only.");
209     return;
210     }
211     my_mpi_rank = mesh_p->Nodes->MPIInfo->rank;
212     mpi_size = mesh_p->Nodes->MPIInfo->size;
213 dhawcroft 793
214 gross 2308 /************************************************************************
215     * open the file and check handle *
216     */
217 caltinay 2141 if (mpi_size > 1) {
218     #ifdef PASO_MPI
219     const int amode = MPI_MODE_CREATE|MPI_MODE_WRONLY|MPI_MODE_UNIQUE_OPEN;
220     int ierr;
221     if (my_mpi_rank == 0 && Paso_fileExists(filename_p)) {
222     remove(filename_p);
223     }
224     ierr = MPI_File_open(mesh_p->Nodes->MPIInfo->comm, (char*)filename_p,
225     amode, mpi_info, &mpi_fileHandle_p);
226     if (ierr != MPI_SUCCESS) {
227     sprintf(errorMsg, "saveVTK: File %s could not be opened for writing in parallel.", filename_p);
228     Finley_setError(IO_ERROR, errorMsg);
229     } else {
230 gross 2310 ierr=MPI_File_set_view(mpi_fileHandle_p,MPI_DISPLACEMENT_CURRENT,
231 caltinay 2141 MPI_CHAR, MPI_CHAR, "native", mpi_info);
232     }
233     #endif /* PASO_MPI */
234     } else {
235 ksteube 1312 fileHandle_p = fopen(filename_p, "w");
236     if (fileHandle_p==NULL) {
237 caltinay 2141 sprintf(errorMsg, "saveVTK: File %s could not be opened for writing.", filename_p);
238     Finley_setError(IO_ERROR, errorMsg);
239     }
240     }
241     if (!Paso_MPIInfo_noError(mesh_p->Nodes->MPIInfo)) return;
242 dhawcroft 793
243 caltinay 2141 /* General note: From this point if an error occurs Finley_setError is
244     * called and subsequent steps are skipped until the end of this function
245     * where allocated memory is freed and the file is closed. */
246 dhawcroft 793
247 caltinay 2141 /************************************************************************/
248     /* find the mesh type to be written */
249 dhawcroft 793
250 caltinay 2141 isCellCentered = TMPMEMALLOC(num_data, bool_t);
251     maxNameLen = 0;
252     if (!Finley_checkPtr(isCellCentered)) {
253     for (dataIdx=0; dataIdx<num_data; ++dataIdx) {
254     if (! isEmpty(data_pp[dataIdx])) {
255     switch(getFunctionSpaceType(data_pp[dataIdx]) ) {
256     case FINLEY_NODES:
257     nodeType = (nodeType == FINLEY_REDUCED_NODES) ? FINLEY_REDUCED_NODES : FINLEY_NODES;
258     isCellCentered[dataIdx] = FALSE;
259     if (elementType==FINLEY_UNKNOWN || elementType==FINLEY_ELEMENTS) {
260     elementType = FINLEY_ELEMENTS;
261     } else {
262     Finley_setError(TYPE_ERROR, "saveVTK: cannot write given data in single file.");
263     }
264     break;
265     case FINLEY_REDUCED_NODES:
266     nodeType = FINLEY_REDUCED_NODES;
267     isCellCentered[dataIdx] = FALSE;
268     if (elementType==FINLEY_UNKNOWN || elementType==FINLEY_ELEMENTS) {
269     elementType = FINLEY_ELEMENTS;
270     } else {
271     Finley_setError(TYPE_ERROR, "saveVTK: cannot write given data in single file.");
272     }
273     break;
274     case FINLEY_REDUCED_ELEMENTS:
275     hasReducedElements = TRUE;
276     case FINLEY_ELEMENTS:
277     isCellCentered[dataIdx] = TRUE;
278     if (elementType==FINLEY_UNKNOWN || elementType==FINLEY_ELEMENTS) {
279     elementType = FINLEY_ELEMENTS;
280     } else {
281     Finley_setError(TYPE_ERROR, "saveVTK: cannot write given data in single file.");
282     }
283     break;
284     case FINLEY_REDUCED_FACE_ELEMENTS:
285     hasReducedElements = TRUE;
286     case FINLEY_FACE_ELEMENTS:
287     isCellCentered[dataIdx] = TRUE;
288     if (elementType==FINLEY_UNKNOWN || elementType==FINLEY_FACE_ELEMENTS) {
289     elementType = FINLEY_FACE_ELEMENTS;
290     } else {
291     Finley_setError(TYPE_ERROR, "saveVTK: cannot write given data in single file.");
292     }
293     break;
294     case FINLEY_POINTS:
295     isCellCentered[dataIdx]=TRUE;
296     if (elementType==FINLEY_UNKNOWN || elementType==FINLEY_POINTS) {
297     elementType = FINLEY_POINTS;
298     } else {
299     Finley_setError(TYPE_ERROR, "saveVTK: cannot write given data in single file.");
300     }
301     break;
302     case FINLEY_REDUCED_CONTACT_ELEMENTS_1:
303     hasReducedElements = TRUE;
304     case FINLEY_CONTACT_ELEMENTS_1:
305     isCellCentered[dataIdx] = TRUE;
306     if (elementType==FINLEY_UNKNOWN || elementType==FINLEY_CONTACT_ELEMENTS_1) {
307     elementType = FINLEY_CONTACT_ELEMENTS_1;
308     } else {
309     Finley_setError(TYPE_ERROR, "saveVTK: cannot write given data in single file.");
310     }
311     break;
312     case FINLEY_REDUCED_CONTACT_ELEMENTS_2:
313     hasReducedElements = TRUE;
314     case FINLEY_CONTACT_ELEMENTS_2:
315     isCellCentered[dataIdx] = TRUE;
316     if (elementType==FINLEY_UNKNOWN || elementType==FINLEY_CONTACT_ELEMENTS_1) {
317     elementType = FINLEY_CONTACT_ELEMENTS_1;
318     } else {
319     Finley_setError(TYPE_ERROR, "saveVTK: cannot write given data in single file.");
320     }
321     break;
322     default:
323     sprintf(errorMsg, "saveVTK: unknown function space type %d",getFunctionSpaceType(data_pp[dataIdx]));
324     Finley_setError(TYPE_ERROR, errorMsg);
325     }
326     if (isCellCentered[dataIdx]) {
327     writeCellData = TRUE;
328     } else {
329     writePointData = TRUE;
330     }
331     maxNameLen = MAX(maxNameLen, strlen(names_p[dataIdx]));
332     }
333     }
334     }
335 dhawcroft 793
336 caltinay 2141 /************************************************************************/
337     /* select number of points and the mesh component */
338 gross 1741
339 caltinay 2141 if (Finley_noError()) {
340     if (nodeType == FINLEY_REDUCED_NODES) {
341     myFirstNode = Finley_NodeFile_getFirstReducedNode(mesh_p->Nodes);
342     myLastNode = Finley_NodeFile_getLastReducedNode(mesh_p->Nodes);
343     globalNumPoints = Finley_NodeFile_getGlobalNumReducedNodes(mesh_p->Nodes);
344     globalNodeIndex = Finley_NodeFile_borrowGlobalReducedNodesIndex(mesh_p->Nodes);
345     } else {
346     myFirstNode = Finley_NodeFile_getFirstNode(mesh_p->Nodes);
347     myLastNode = Finley_NodeFile_getLastNode(mesh_p->Nodes);
348     globalNumPoints = Finley_NodeFile_getGlobalNumNodes(mesh_p->Nodes);
349     globalNodeIndex = Finley_NodeFile_borrowGlobalNodesIndex(mesh_p->Nodes);
350 ksteube 1312 }
351 caltinay 2141 myNumPoints = myLastNode - myFirstNode;
352     if (elementType==FINLEY_UNKNOWN) elementType=FINLEY_ELEMENTS;
353     switch(elementType) {
354     case FINLEY_ELEMENTS:
355     elements = mesh_p->Elements;
356     break;
357     case FINLEY_FACE_ELEMENTS:
358     elements = mesh_p->FaceElements;
359     break;
360     case FINLEY_POINTS:
361     elements = mesh_p->Points;
362     break;
363     case FINLEY_CONTACT_ELEMENTS_1:
364     elements = mesh_p->ContactElements;
365     break;
366     }
367     if (elements==NULL) {
368     Finley_setError(SYSTEM_ERROR, "saveVTK: undefined element file");
369     } else {
370     /* map finley element type to VTK element type */
371     numCells = elements->numElements;
372     globalNumCells = Finley_ElementFile_getGlobalNumElements(elements);
373     myNumCells = Finley_ElementFile_getMyNumElements(elements);
374     myFirstCell = Finley_ElementFile_getFirstElement(elements);
375     NN = elements->numNodes;
376     if (nodeType==FINLEY_REDUCED_NODES) {
377     typeId = elements->LinearReferenceElement->Type->TypeId;
378     if (hasReducedElements) {
379     quadNodes_p=elements->LinearReferenceElementReducedOrder->QuadNodes;
380     } else {
381     quadNodes_p=elements->LinearReferenceElement->QuadNodes;
382     }
383     } else {
384     typeId = elements->ReferenceElement->Type->TypeId;
385     if (hasReducedElements)
386     quadNodes_p=elements->ReferenceElementReducedOrder->QuadNodes;
387     else
388     quadNodes_p=elements->ReferenceElement->QuadNodes;
389     }
390     switch (typeId) {
391     case Point1:
392     case Line2Face:
393     case Line3Face:
394     case Point1_Contact:
395     case Line2Face_Contact:
396     case Line3Face_Contact:
397     cellType = VTK_VERTEX;
398     numVTKNodesPerElement = 1;
399     break;
400 dhawcroft 793
401 caltinay 2141 case Line2:
402     case Tri3Face:
403     case Rec4Face:
404     case Line2_Contact:
405     case Tri3_Contact:
406     case Tri3Face_Contact:
407     case Rec4Face_Contact:
408     cellType = VTK_LINE;
409     numVTKNodesPerElement = 2;
410     break;
411 dhawcroft 793
412 caltinay 2141 case Tri3:
413     case Tet4Face:
414     case Tet4Face_Contact:
415     cellType = VTK_TRIANGLE;
416     numVTKNodesPerElement = 3;
417     break;
418 dhawcroft 793
419 caltinay 2141 case Rec4:
420     case Hex8Face:
421     case Rec4_Contact:
422     case Hex8Face_Contact:
423     cellType = VTK_QUAD;
424     numVTKNodesPerElement = 4;
425     break;
426 dhawcroft 793
427 caltinay 2141 case Rec9:
428     numCellFactor = 4;
429     cellType = VTK_QUAD;
430     numVTKNodesPerElement = 4;
431     break;
432 dhawcroft 793
433 caltinay 2141 case Tet4:
434     cellType = VTK_TETRA;
435     numVTKNodesPerElement = 4;
436     break;
437 dhawcroft 793
438 caltinay 2141 case Hex8:
439     cellType = VTK_HEXAHEDRON;
440     numVTKNodesPerElement = 8;
441     break;
442 dhawcroft 796
443 caltinay 2141 case Line3:
444     case Tri6Face:
445     case Rec8Face:
446     case Line3_Contact:
447     case Tri6Face_Contact:
448     case Rec8Face_Contact:
449     cellType = VTK_QUADRATIC_EDGE;
450     numVTKNodesPerElement = 3;
451     break;
452 dhawcroft 793
453 caltinay 2141 case Tri6:
454     case Tet10Face:
455     case Tri6_Contact:
456     case Tet10Face_Contact:
457     cellType = VTK_QUADRATIC_TRIANGLE;
458     numVTKNodesPerElement = 6;
459     break;
460 dhawcroft 793
461 caltinay 2141 case Rec8:
462     case Hex20Face:
463     case Rec8_Contact:
464     case Hex20Face_Contact:
465     cellType = VTK_QUADRATIC_QUAD;
466     numVTKNodesPerElement = 8;
467     break;
468 dhawcroft 793
469 caltinay 2141 case Tet10:
470     cellType = VTK_QUADRATIC_TETRA;
471     numVTKNodesPerElement = 10;
472     break;
473 gross 1741
474 caltinay 2141 case Hex20:
475     cellType = VTK_QUADRATIC_HEXAHEDRON;
476     numVTKNodesPerElement = 20;
477     break;
478 gross 1741
479 caltinay 2141 case Hex27:
480     numCellFactor = 8;
481     cellType = VTK_HEXAHEDRON;
482     numVTKNodesPerElement = 8;
483     break;
484 gross 1741
485 caltinay 2141 default:
486     sprintf(errorMsg, "saveVTK: Element type %s is not supported by VTK.", elements->ReferenceElement->Type->Name);
487     Finley_setError(VALUE_ERROR, errorMsg);
488     }
489     }
490     }
491 gross 1741
492 caltinay 2141 /* allocate enough memory for text buffer */
493 gross 1741
494 gross 2421 txtBufferSize = strlen(vtkHeader) + 3*LEN_INT_FORMAT + (30+3*maxNameLen)+strlen(metadata)+strlen(metadata_schema);
495 caltinay 2141 if (mpi_size > 1) {
496 phornby 2172 txtBufferSize = MAX(txtBufferSize, myNumPoints * LEN_TMP_BUFFER);
497 caltinay 2141 txtBufferSize = MAX(txtBufferSize, numCellFactor * myNumCells *
498     (LEN_INT_FORMAT * numVTKNodesPerElement + 1));
499     txtBufferSize = MAX(txtBufferSize,
500     numCellFactor * myNumCells * LEN_TENSOR_FORMAT);
501     txtBufferSize = MAX(txtBufferSize, myNumPoints * LEN_TENSOR_FORMAT);
502     }
503     txtBuffer = TMPMEMALLOC(txtBufferSize+1, char);
504 gross 1741
505 caltinay 2141 /* sets error if memory allocation failed */
506     Finley_checkPtr(txtBuffer);
507 gross 1741
508 caltinay 2141 /************************************************************************/
509     /* write number of points and the mesh component */
510 gross 1741
511 caltinay 2141 if (Finley_noError()) {
512     const index_t *nodeIndex;
513     if (FINLEY_REDUCED_NODES == nodeType) {
514     nodeIndex = elements->ReferenceElement->Type->linearNodes;
515     } else if (Rec9 == typeId) {
516     nodeIndex = VTK_REC9_INDEX;
517     } else if (Hex20 == typeId) {
518     nodeIndex = VTK_HEX20_INDEX;
519     } else if (Hex27 == typeId) {
520     nodeIndex = VTK_HEX27_INDEX;
521     } else if (numVTKNodesPerElement != NN) {
522     nodeIndex = elements->ReferenceElement->Type->geoNodes;
523 ksteube 1312 } else {
524 caltinay 2141 nodeIndex = NULL;
525 ksteube 1312 }
526 gross 1741
527 gross 2421 if (strlen(metadata)>0) {
528     if (strlen(metadata_schema)>0) {
529     sprintf(txtBuffer, vtkHeader," ",metadata_schema,metadata,"\n",globalNumPoints, numCellFactor*globalNumCells, 3);
530     } else {
531     sprintf(txtBuffer, vtkHeader,"","",metadata,"\n",globalNumPoints, numCellFactor*globalNumCells, 3);
532     }
533     } else {
534     if (strlen(metadata_schema)>0) {
535     sprintf(txtBuffer, vtkHeader," ",metadata_schema,"","",globalNumPoints, numCellFactor*globalNumCells, 3);
536     } else {
537     sprintf(txtBuffer, vtkHeader,"","","","",globalNumPoints, numCellFactor*globalNumCells, 3);
538     }
539     }
540 gross 1741
541 caltinay 2141 if (mpi_size > 1) {
542     /* write the nodes */
543     MPI_RANK0_WRITE_SHARED(txtBuffer);
544     txtBuffer[0] = '\0';
545     txtBufferInUse = 0;
546     if (nDim==2) {
547     for (i = 0; i < mesh_p->Nodes->numNodes; i++) {
548     if ( (myFirstNode <= globalNodeIndex[i]) && (globalNodeIndex[i] < myLastNode) ) {
549     sprintf(tmpBuffer, VECTOR_FORMAT,
550     mesh_p->Nodes->Coordinates[INDEX2(0, i, nDim)],
551     mesh_p->Nodes->Coordinates[INDEX2(1, i, nDim)],
552     0.);
553     __STRCAT(txtBuffer, tmpBuffer, txtBufferInUse);
554     }
555     }
556     } else {
557     for (i = 0; i < mesh_p->Nodes->numNodes; i++) {
558     if ( (myFirstNode <= globalNodeIndex[i]) && (globalNodeIndex[i] < myLastNode) ) {
559     sprintf(tmpBuffer, VECTOR_FORMAT,
560     mesh_p->Nodes->Coordinates[INDEX2(0, i, nDim)],
561     mesh_p->Nodes->Coordinates[INDEX2(1, i, nDim)],
562     mesh_p->Nodes->Coordinates[INDEX2(2, i, nDim)]);
563     __STRCAT(txtBuffer, tmpBuffer, txtBufferInUse);
564     }
565     }
566     } /* nDim */
567 gross 2385 MPI_WRITE_ORDERED(txtBuffer);
568 gross 1741
569 caltinay 2141 /* write the cells */
570     MPI_RANK0_WRITE_SHARED(tags_End_Points_and_Start_Conn);
571     txtBuffer[0] = '\0';
572     txtBufferInUse = 0;
573     if (nodeIndex == NULL) {
574     for (i = 0; i < numCells; i++) {
575     if (elements->Owner[i] == my_mpi_rank) {
576     for (j = 0; j < numVTKNodesPerElement; j++) {
577     sprintf(tmpBuffer, INT_FORMAT, globalNodeIndex[elements->Nodes[INDEX2(j, i, NN)]]);
578     __STRCAT(txtBuffer, tmpBuffer, txtBufferInUse);
579     }
580     __STRCAT(txtBuffer, NEWLINE, txtBufferInUse);
581     }
582     }
583     } else {
584     for (i = 0; i < numCells; i++) {
585     if (elements->Owner[i] == my_mpi_rank) {
586     for (l = 0; l < numCellFactor; l++) {
587     const int* idx=&nodeIndex[l*numVTKNodesPerElement];
588     for (j = 0; j < numVTKNodesPerElement; j++) {
589     sprintf(tmpBuffer, INT_FORMAT, globalNodeIndex[elements->Nodes[INDEX2(idx[j], i, NN)]]);
590     __STRCAT(txtBuffer, tmpBuffer, txtBufferInUse);
591     }
592     __STRCAT(txtBuffer, NEWLINE, txtBufferInUse);
593     }
594     }
595     }
596     } /* nodeIndex */
597 gross 2385 MPI_WRITE_ORDERED(txtBuffer);
598 gross 1741
599 caltinay 2141 /* write the offsets */
600     MPI_RANK0_WRITE_SHARED(tags_End_Conn_and_Start_Offset);
601     txtBuffer[0] = '\0';
602     txtBufferInUse = 0;
603     for (i = numVTKNodesPerElement*(myFirstCell*numCellFactor+1);
604     i <= (myFirstCell+myNumCells)*numVTKNodesPerElement*numCellFactor; i += numVTKNodesPerElement)
605     {
606     sprintf(tmpBuffer, INT_NEWLINE_FORMAT, i);
607     __STRCAT(txtBuffer, tmpBuffer, txtBufferInUse);
608     }
609 gross 2385 MPI_WRITE_ORDERED(txtBuffer);
610 gross 1741
611 caltinay 2141 /* write element type */
612     sprintf(tmpBuffer, INT_NEWLINE_FORMAT, cellType);
613     MPI_RANK0_WRITE_SHARED(tags_End_Offset_and_Start_Type);
614     txtBuffer[0] = '\0';
615     txtBufferInUse = 0;
616     for (i = numVTKNodesPerElement*(myFirstCell*numCellFactor+1);
617     i <= (myFirstCell+myNumCells)*numVTKNodesPerElement*numCellFactor; i += numVTKNodesPerElement)
618     {
619     __STRCAT(txtBuffer, tmpBuffer, txtBufferInUse);
620     }
621 gross 2385 MPI_WRITE_ORDERED(txtBuffer);
622 caltinay 2141 /* finalize cell information */
623     strcpy(txtBuffer, "</DataArray>\n</Cells>\n");
624     MPI_RANK0_WRITE_SHARED(txtBuffer);
625 gross 1741
626 caltinay 2141 } else { /***** mpi_size == 1 *****/
627 gross 1741
628 caltinay 2141 /* write the nodes */
629     fputs(txtBuffer, fileHandle_p);
630     if (nDim==2) {
631     for (i = 0; i < mesh_p->Nodes->numNodes; i++) {
632     if ( (myFirstNode <= globalNodeIndex[i]) && (globalNodeIndex[i] < myLastNode) ) {
633     fprintf(fileHandle_p, VECTOR_FORMAT,
634     mesh_p->Nodes->Coordinates[INDEX2(0, i, nDim)],
635     mesh_p->Nodes->Coordinates[INDEX2(1, i, nDim)],
636     0.);
637     }
638     }
639     } else {
640     for (i = 0; i < mesh_p->Nodes->numNodes; i++) {
641     if ( (myFirstNode <= globalNodeIndex[i]) && (globalNodeIndex[i] < myLastNode) ) {
642     fprintf(fileHandle_p, VECTOR_FORMAT,
643     mesh_p->Nodes->Coordinates[INDEX2(0, i, nDim)],
644     mesh_p->Nodes->Coordinates[INDEX2(1, i, nDim)],
645     mesh_p->Nodes->Coordinates[INDEX2(2, i, nDim)]);
646     }
647     }
648     } /* nDim */
649 gross 1741
650 caltinay 2141 /* write the cells */
651     fputs(tags_End_Points_and_Start_Conn, fileHandle_p);
652     if (nodeIndex == NULL) {
653     for (i = 0; i < numCells; i++) {
654     for (j = 0; j < numVTKNodesPerElement; j++) {
655     fprintf(fileHandle_p, INT_FORMAT, globalNodeIndex[elements->Nodes[INDEX2(j, i, NN)]]);
656     }
657     fprintf(fileHandle_p, NEWLINE);
658     }
659     } else {
660     for (i = 0; i < numCells; i++) {
661     for (l = 0; l < numCellFactor; l++) {
662     const int* idx=&nodeIndex[l*numVTKNodesPerElement];
663     for (j = 0; j < numVTKNodesPerElement; j++) {
664     fprintf(fileHandle_p, INT_FORMAT, globalNodeIndex[elements->Nodes[INDEX2(idx[j], i, NN)]]);
665     }
666     fprintf(fileHandle_p, NEWLINE);
667     }
668     }
669     } /* nodeIndex */
670 gross 1741
671 caltinay 2141 /* write the offsets */
672     fputs(tags_End_Conn_and_Start_Offset, fileHandle_p);
673     for (i = numVTKNodesPerElement; i <= numCells*numVTKNodesPerElement*numCellFactor; i += numVTKNodesPerElement) {
674     fprintf(fileHandle_p, INT_NEWLINE_FORMAT, i);
675     }
676 gross 1741
677 caltinay 2141 /* write element type */
678     sprintf(tmpBuffer, INT_NEWLINE_FORMAT, cellType);
679     fputs(tags_End_Offset_and_Start_Type, fileHandle_p);
680     for (i = 0; i < numCells*numCellFactor; i++)
681     fputs(tmpBuffer, fileHandle_p);
682     /* finalize cell information */
683     fputs("</DataArray>\n</Cells>\n", fileHandle_p);
684     } /* mpi_size */
685 gross 1741
686 caltinay 2141 } /* Finley_noError */
687 dhawcroft 796
688 caltinay 2141 /************************************************************************/
689     /* write cell data */
690 dhawcroft 796
691 caltinay 2141 if (writeCellData && Finley_noError()) {
692     bool_t set_scalar=FALSE, set_vector=FALSE, set_tensor=FALSE;
693     /* mark the active data arrays */
694     strcpy(txtBuffer, "<CellData");
695     for (dataIdx = 0; dataIdx < num_data; dataIdx++) {
696     if (!isEmpty(data_pp[dataIdx]) && isCellCentered[dataIdx]) {
697     /* rank == 0 <--> scalar data */
698     /* rank == 1 <--> vector data */
699     /* rank == 2 <--> tensor data */
700     switch (getDataPointRank(data_pp[dataIdx])) {
701     case 0:
702     if (!set_scalar) {
703     strcat(txtBuffer, " Scalars=\"");
704     strcat(txtBuffer, names_p[dataIdx]);
705     strcat(txtBuffer, "\"");
706     set_scalar = TRUE;
707     }
708     break;
709     case 1:
710     if (!set_vector) {
711     strcat(txtBuffer, " Vectors=\"");
712     strcat(txtBuffer, names_p[dataIdx]);
713     strcat(txtBuffer, "\"");
714     set_vector = TRUE;
715     }
716     break;
717     case 2:
718     if (!set_tensor) {
719     strcat(txtBuffer, " Tensors=\"");
720     strcat(txtBuffer, names_p[dataIdx]);
721     strcat(txtBuffer, "\"");
722     set_tensor = TRUE;
723     }
724     break;
725     default:
726     sprintf(errorMsg, "saveVTK: data %s: VTK supports data with rank <= 2 only.", names_p[dataIdx]);
727     Finley_setError(VALUE_ERROR, errorMsg);
728     }
729     }
730     if (!Finley_noError())
731     break;
732 dhawcroft 818 }
733 dhawcroft 793 }
734 caltinay 2141 /* only continue if no error occurred */
735     if (writeCellData && Finley_noError()) {
736     strcat(txtBuffer, ">\n");
737     if ( mpi_size > 1) {
738     MPI_RANK0_WRITE_SHARED(txtBuffer);
739     } else {
740     fputs(txtBuffer, fileHandle_p);
741 dhawcroft 793 }
742 dhawcroft 818
743 caltinay 2141 /* write the arrays */
744     for (dataIdx = 0; dataIdx < num_data; dataIdx++) {
745     if (!isEmpty(data_pp[dataIdx]) && isCellCentered[dataIdx]) {
746     dim_t numPointsPerSample=getNumDataPointsPerSample(data_pp[dataIdx]);
747     dim_t rank = getDataPointRank(data_pp[dataIdx]);
748     dim_t nComp = getDataPointSize(data_pp[dataIdx]);
749     dim_t nCompReqd = 1; /* number of components mpi_required */
750     if (rank == 0) {
751     nCompReqd = 1;
752     shape = 0;
753     } else if (rank == 1) {
754     nCompReqd = 3;
755     shape = getDataPointShape(data_pp[dataIdx], 0);
756     if (shape > 3) {
757     Finley_setError(VALUE_ERROR, "saveVTK: rank 1 objects must have 3 components at most.");
758     }
759     } else {
760     nCompReqd = 9;
761     shape = getDataPointShape(data_pp[dataIdx], 0);
762     if (shape > 3 || shape != getDataPointShape(data_pp[dataIdx], 1)) {
763     Finley_setError(VALUE_ERROR, "saveVTK: rank 2 objects of shape 2x2 or 3x3 supported only.");
764     }
765     }
766     /* bail out if an error occurred */
767     if (!Finley_noError())
768     break;
769 dhawcroft 793
770 caltinay 2141 sprintf(txtBuffer, tag_Float_DataArray, names_p[dataIdx], nCompReqd);
771     if ( mpi_size > 1) {
772     MPI_RANK0_WRITE_SHARED(txtBuffer);
773     } else {
774     fputs(txtBuffer, fileHandle_p);
775     }
776 gross 1741
777 caltinay 2141 txtBuffer[0] = '\0';
778     txtBufferInUse = 0;
779     for (i=0; i<numCells; i++) {
780     if (elements->Owner[i] == my_mpi_rank) {
781 jfenwick 2271 void* sampleBuffer=allocSampleBuffer(data_pp[dataIdx]);
782     __const double *values = getSampleDataRO(data_pp[dataIdx], i,sampleBuffer);
783 caltinay 2141 for (l = 0; l < numCellFactor; l++) {
784     double sampleAvg[NCOMP_MAX];
785     dim_t nCompUsed = MIN(nComp, NCOMP_MAX);
786    
787     /* average over number of points in the sample */
788     if (isExpanded(data_pp[dataIdx])) {
789     dim_t hits=0, hits_old;
790     for (k=0; k<nCompUsed; k++) sampleAvg[k]=0;
791     for (j=0; j<numPointsPerSample; j++) {
792     hits_old=hits;
793     if (nodeInQuadrant(quadNodes_p, typeId, j, l)) {
794     hits++;
795     for (k=0; k<nCompUsed; k++) {
796     sampleAvg[k] += values[INDEX2(k,j,nComp)];
797     }
798 gross 1741 }
799 caltinay 2141 }
800     for (k=0; k<nCompUsed; k++)
801     sampleAvg[k] /= MAX(hits, 1);
802     } else {
803     for (k=0; k<nCompUsed; k++)
804     sampleAvg[k] = values[k];
805     } /* isExpanded */
806    
807     /* if the number of required components is more than
808     * the number of actual components, pad with zeros
809     */
810     /* probably only need to get shape of first element */
811     if (nCompReqd == 1) {
812     sprintf(tmpBuffer, SCALAR_FORMAT, sampleAvg[0]);
813     } else if (nCompReqd == 3) {
814     if (shape==1) {
815     sprintf(tmpBuffer, VECTOR_FORMAT,
816     sampleAvg[0], 0.f, 0.f);
817     } else if (shape==2) {
818     sprintf(tmpBuffer, VECTOR_FORMAT,
819     sampleAvg[0], sampleAvg[1], 0.f);
820     } else if (shape==3) {
821     sprintf(tmpBuffer, VECTOR_FORMAT,
822     sampleAvg[0],sampleAvg[1],sampleAvg[2]);
823     }
824     } else if (nCompReqd == 9) {
825     if (shape==1) {
826     sprintf(tmpBuffer, TENSOR_FORMAT,
827     sampleAvg[0], 0.f, 0.f,
828     0.f, 0.f, 0.f,
829     0.f, 0.f, 0.f);
830     } else if (shape==2) {
831     sprintf(tmpBuffer, TENSOR_FORMAT,
832     sampleAvg[0], sampleAvg[1], 0.f,
833     sampleAvg[2], sampleAvg[3], 0.f,
834     0.f, 0.f, 0.f);
835     } else if (shape==3) {
836     sprintf(tmpBuffer, TENSOR_FORMAT,
837     sampleAvg[0],sampleAvg[1],sampleAvg[2],
838     sampleAvg[3],sampleAvg[4],sampleAvg[5],
839     sampleAvg[6],sampleAvg[7],sampleAvg[8]);
840     }
841     }
842 gross 1741 if ( mpi_size > 1) {
843 caltinay 2141 __STRCAT(txtBuffer, tmpBuffer, txtBufferInUse);
844 gross 1741 } else {
845 caltinay 2141 fputs(tmpBuffer, fileHandle_p);
846 gross 1741 }
847 caltinay 2141 } /* for l (numCellFactor) */
848 jfenwick 2271 freeSampleBuffer(sampleBuffer);
849 caltinay 2141 } /* if I am the owner */
850     } /* for i (numCells) */
851    
852     if ( mpi_size > 1) {
853 gross 2385 MPI_WRITE_ORDERED(txtBuffer);
854 caltinay 2141 MPI_RANK0_WRITE_SHARED(tag_End_DataArray);
855     } else {
856     fputs(tag_End_DataArray, fileHandle_p);
857     }
858     } /* !isEmpty && cellCentered */
859     } /* for dataIdx */
860    
861     strcpy(txtBuffer, "</CellData>\n");
862     if ( mpi_size > 1) {
863     MPI_RANK0_WRITE_SHARED(txtBuffer);
864     } else {
865     fputs(txtBuffer, fileHandle_p);
866 dhawcroft 793 }
867 caltinay 2141 } /* if noError && writeCellData */
868 dhawcroft 793
869 caltinay 2141 /************************************************************************/
870     /* write point data */
871    
872     if (writePointData && Finley_noError()) {
873     /* mark the active data arrays */
874     bool_t set_scalar=FALSE, set_vector=FALSE, set_tensor=FALSE;
875     strcpy(txtBuffer, "<PointData");
876     for (dataIdx = 0; dataIdx<num_data; dataIdx++) {
877     if (!isEmpty(data_pp[dataIdx]) && !isCellCentered[dataIdx]) {
878     switch (getDataPointRank(data_pp[dataIdx])) {
879     case 0:
880     if (!set_scalar) {
881     strcat(txtBuffer, " Scalars=\"");
882     strcat(txtBuffer, names_p[dataIdx]);
883     strcat(txtBuffer, "\"");
884     set_scalar = TRUE;
885     }
886     break;
887     case 1:
888     if (!set_vector) {
889     strcat(txtBuffer, " Vectors=\"");
890     strcat(txtBuffer, names_p[dataIdx]);
891     strcat(txtBuffer, "\"");
892     set_vector = TRUE;
893     }
894     break;
895     case 2:
896     if (!set_tensor) {
897     strcat(txtBuffer, " Tensors=\"");
898     strcat(txtBuffer, names_p[dataIdx]);
899     strcat(txtBuffer, "\"");
900     set_tensor = TRUE;
901     }
902     break;
903     default:
904     sprintf(errorMsg, "saveVTK: data %s: VTK supports data with rank <= 2 only.", names_p[dataIdx]);
905     Finley_setError(VALUE_ERROR, errorMsg);
906     }
907 dhawcroft 793 }
908 caltinay 2141 if (!Finley_noError())
909     break;
910 dhawcroft 793 }
911 caltinay 2141 }
912     /* only continue if no error occurred */
913     if (writePointData && Finley_noError()) {
914     strcat(txtBuffer, ">\n");
915     if ( mpi_size > 1) {
916     MPI_RANK0_WRITE_SHARED(txtBuffer);
917     } else {
918     fputs(txtBuffer, fileHandle_p);
919 gross 1062 }
920 caltinay 2141
921     /* write the arrays */
922     for (dataIdx = 0; dataIdx < num_data; dataIdx++) {
923     if (!isEmpty(data_pp[dataIdx]) && !isCellCentered[dataIdx]) {
924     Finley_NodeMapping* nodeMapping;
925     dim_t rank = getDataPointRank(data_pp[dataIdx]);
926     dim_t nCompReqd = 1; /* number of components mpi_required */
927     if (getFunctionSpaceType(data_pp[dataIdx]) == FINLEY_REDUCED_NODES) {
928     nodeMapping = mesh_p->Nodes->reducedNodesMapping;
929     } else {
930     nodeMapping = mesh_p->Nodes->nodesMapping;
931     }
932     if (rank == 0) {
933     nCompReqd = 1;
934     shape = 0;
935     } else if (rank == 1) {
936     nCompReqd = 3;
937     shape = getDataPointShape(data_pp[dataIdx], 0);
938     if (shape > 3) {
939     Finley_setError(VALUE_ERROR, "saveVTK: rank 1 objects must have 3 components at most.");
940     }
941     } else {
942     nCompReqd = 9;
943     shape=getDataPointShape(data_pp[dataIdx], 0);
944     if (shape > 3 || shape != getDataPointShape(data_pp[dataIdx], 1)) {
945     Finley_setError(VALUE_ERROR, "saveVTK: rank 2 objects of shape 2x2 or 3x3 supported only.");
946     }
947     }
948     /* bail out if an error occurred */
949     if (!Finley_noError())
950     break;
951    
952     sprintf(txtBuffer, tag_Float_DataArray, names_p[dataIdx], nCompReqd);
953     if ( mpi_size > 1) {
954     MPI_RANK0_WRITE_SHARED(txtBuffer);
955     } else {
956     fputs(txtBuffer, fileHandle_p);
957     }
958    
959     txtBuffer[0] = '\0';
960     txtBufferInUse = 0;
961     for (i=0; i<mesh_p->Nodes->numNodes; i++) {
962     k = globalNodeIndex[i];
963     if ( (myFirstNode <= k) && (k < myLastNode) ) {
964 jfenwick 2271 void* sampleBuffer=allocSampleBuffer(data_pp[dataIdx]);
965     __const double *values = getSampleDataRO(data_pp[dataIdx], nodeMapping->target[i], sampleBuffer);
966 caltinay 2141 /* if the number of mpi_required components is more than
967     * the number of actual components, pad with zeros.
968     * Probably only need to get shape of first element */
969     if (nCompReqd == 1) {
970     sprintf(tmpBuffer, SCALAR_FORMAT, values[0]);
971     } else if (nCompReqd == 3) {
972     if (shape==1) {
973     sprintf(tmpBuffer, VECTOR_FORMAT,
974     values[0], 0.f, 0.f);
975     } else if (shape==2) {
976     sprintf(tmpBuffer, VECTOR_FORMAT,
977     values[0], values[1], 0.f);
978     } else if (shape==3) {
979     sprintf(tmpBuffer, VECTOR_FORMAT,
980     values[0], values[1], values[2]);
981     }
982     } else if (nCompReqd == 9) {
983     if (shape==1) {
984     sprintf(tmpBuffer, TENSOR_FORMAT,
985     values[0], 0.f, 0.f,
986     0.f, 0.f, 0.f,
987     0.f, 0.f, 0.f);
988     } else if (shape==2) {
989     sprintf(tmpBuffer, TENSOR_FORMAT,
990     values[0], values[1], 0.f,
991     values[2], values[3], 0.f,
992     0.f, 0.f, 0.f);
993     } else if (shape==3) {
994     sprintf(tmpBuffer, TENSOR_FORMAT,
995     values[0], values[1], values[2],
996     values[3], values[4], values[5],
997     values[6], values[7], values[8]);
998     }
999     }
1000     if ( mpi_size > 1) {
1001     __STRCAT(txtBuffer, tmpBuffer, txtBufferInUse);
1002     } else {
1003     fputs(tmpBuffer, fileHandle_p);
1004     }
1005 jfenwick 2271 freeSampleBuffer(sampleBuffer); /* no-one needs values anymore */
1006 caltinay 2141 } /* if this is my node */
1007     } /* for i (numNodes) */
1008    
1009     if ( mpi_size > 1) {
1010 gross 2385 MPI_WRITE_ORDERED(txtBuffer);
1011 caltinay 2141 MPI_RANK0_WRITE_SHARED(tag_End_DataArray);
1012     } else {
1013     fputs(tag_End_DataArray, fileHandle_p);
1014     }
1015     } /* !isEmpty && !isCellCentered */
1016     } /* for dataIdx */
1017    
1018     strcpy(txtBuffer, "</PointData>\n");
1019 ksteube 1312 if ( mpi_size > 1) {
1020 caltinay 2141 MPI_RANK0_WRITE_SHARED(txtBuffer);
1021 gross 1062 } else {
1022 caltinay 2141 fputs(txtBuffer, fileHandle_p);
1023 gross 1062 }
1024 caltinay 2141 } /* if noError && writePointData */
1025    
1026     /* Final write to VTK file */
1027     if (Finley_noError()) {
1028     if (mpi_size > 1) {
1029     MPI_RANK0_WRITE_SHARED(vtkFooter);
1030     } else {
1031     fputs(vtkFooter, fileHandle_p);
1032 dhawcroft 793 }
1033 caltinay 2141 }
1034    
1035     if ( mpi_size > 1) {
1036 caltinay 2094 #ifdef PASO_MPI
1037 caltinay 2141 MPI_File_close(&mpi_fileHandle_p);
1038 jfenwick 2256 MPI_Barrier(mesh_p->Nodes->MPIInfo->comm);
1039 caltinay 2094 #endif
1040 caltinay 2141 } else {
1041     fclose(fileHandle_p);
1042     }
1043     TMPMEMFREE(isCellCentered);
1044     TMPMEMFREE(txtBuffer);
1045 jgs 110 }
1046 caltinay 2141

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26