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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1743 - (show annotations)
Tue Sep 2 00:54:05 2008 UTC (11 years ago) by ksteube
File MIME type: text/plain
File size: 66966 byte(s)
Avoid crash due to zero-length writes in MPI_File_write_ordered.
Write a space instead of empty string.
The space will be at the end of a line and cause no problem.

1
2 /* $Id$ */
3
4 /*******************************************************
5 *
6 * Copyright 2003-2007 by ACceSS MNRF
7 * Copyright 2007 by University of Queensland
8 *
9 * http://esscc.uq.edu.au
10 * Primary Business: Queensland, Australia
11 * Licensed under the Open Software License version 3.0
12 * http://www.opensource.org/licenses/osl-3.0.php
13 *
14 *******************************************************/
15
16 /**************************************************************/
17
18 /* writes data and mesh in a vtk file */
19 /* nodal data needs to be given on FINLEY_NODES or FINLEY_REDUCED_NODES */
20
21 /**************************************************************/
22
23
24 #include "Mesh.h"
25 #include "Assemble.h"
26 #include "vtkCellType.h" /* copied from vtk source directory !!! */
27 #include "paso/PasoUtil.h"
28
29 #define LEN_PRINTED_INT_FORMAT (9+1)
30 #define INT_FORMAT "%d "
31 #define INT_NEWLINE_FORMAT "%d\n"
32 #define FLOAT_SCALAR_FORMAT "%12.6e\n"
33 #define FLOAT_VECTOR_FORMAT "%12.6e %12.6e %12.6e\n"
34 #define FLOAT_TENSOR_FORMAT "%12.6e %12.6e %12.6e %12.6e %12.6e %12.6e %12.6e %12.6e %12.6e\n"
35 #define LEN_PRINTED_FLOAT_SCALAR_FORMAT (12+1)
36 #define LEN_PRINTED_FLOAT_VECTOR_FORMAT (3*(12+1)+1)
37 #define LEN_PRINTED_FLOAT_TENSOR_FORMAT (9*(12+1)+1)
38 #define NEWLINE "\n"
39 #define LEN_TMP_BUFFER LEN_PRINTED_FLOAT_TENSOR_FORMAT+(MAX_numNodes*LEN_PRINTED_INT_FORMAT+1)+2
40 #define NCOMP_MAX 9
41 #define __STRCAT(dest,chunk,dest_in_use) \
42 { \
43 strcpy(&dest[dest_in_use], chunk); \
44 dest_in_use+=strlen(chunk); \
45 }
46 #define INSIDE_1D(_X_,_C_,_R_) ( ABS((_X_)-(_C_)) <= (_R_) )
47 #define INSIDE_2D(_X_,_Y_,_CX_,_CY_,_R_) ( INSIDE_1D(_X_,_CX_,_R_) && INSIDE_1D(_Y_,_CY_,_R_))
48 #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_) )
49
50 void Finley_Mesh_saveVTK(const char * filename_p,
51 Finley_Mesh *mesh_p,
52 const dim_t num_data,
53 char* *names_p,
54 escriptDataC* *data_pp)
55 {
56 #ifdef USE_VTK
57 char error_msg[LenErrorMsg_MAX], *txt_buffer=NULL, tmp_buffer[LEN_TMP_BUFFER];
58 double sampleAvg[NCOMP_MAX], *values, rtmp, *QuadNodes;
59 size_t txt_buffer_in_use;
60 dim_t len_txt_buffer, max_len_names;
61 FILE * fileHandle_p = NULL;
62 int mpi_size, i, j, l, cellType;
63 dim_t i_data, hits, hits_old;
64 dim_t nDim, globalNumPoints, numCells, globalNumCells, numVTKNodesPerElement;
65 dim_t myNumPoints, numPointsPerSample, rank, nComp, nCompReqd;
66 dim_t shape, NN, numCellFactor, myNumCells, max_name_len;
67 bool_t *isCellCentered=NULL,write_celldata=FALSE,write_pointdata=FALSE, reduced_elements=FALSE;
68 bool_t set_scalar=FALSE,set_vector=FALSE, set_tensor=FALSE;
69 index_t myFirstNode, myLastNode, *globalNodeIndex, k, *node_index, myFirstCell;
70 #ifdef PASO_MPI
71 int ierr;
72 /* int amode = MPI_MODE_CREATE | MPI_MODE_WRONLY | MPI_MODE_SEQUENTIAL; */
73 const int amode = MPI_MODE_CREATE | MPI_MODE_WRONLY | MPI_MODE_UNIQUE_OPEN;
74 MPI_File mpi_fileHandle_p;
75 MPI_Status mpi_status;
76 MPI_Request mpi_req;
77 MPI_Info mpi_info=MPI_INFO_NULL;
78 #endif
79 Paso_MPI_rank my_mpi_rank;
80 int nodetype=FINLEY_NODES;
81 int elementtype=FINLEY_UNKNOWN;
82 char elemTypeStr[32];
83 Finley_NodeMapping *nodeMapping=NULL;
84 Finley_ElementFile* elements=NULL;
85 ElementTypeId TypeId;
86
87
88 /****************************************/
89 /* */
90 /* tags in the vtk file */
91
92 char* tags_header="<?xml version=\"1.0\"?>\n" \
93 "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\">\n" \
94 "<UnstructuredGrid>\n" \
95 "<Piece NumberOfPoints=\"%d\" NumberOfCells=\"%d\">\n" \
96 "<Points>\n" \
97 "<DataArray NumberOfComponents=\"%d\" type=\"Float64\" format=\"ascii\">\n";
98 char *tag_End_DataArray = "</DataArray>\n";
99 char* tag_End_PointData = "</PointData>\n";
100 char* tag_End_CellData = "</CellData>\n";
101 char *footer = "</Piece>\n</UnstructuredGrid>\n</VTKFile>\n";
102 char* tags_End_Points_and_Start_Conn = "</DataArray>\n</Points>\n<Cells>\n<DataArray Name=\"connectivity\" type=\"Int32\" format=\"ascii\">\n" ;
103 char* tags_End_Conn_and_Start_Offset = "</DataArray>\n<DataArray Name=\"offsets\" type=\"Int32\" format=\"ascii\">\n";
104 char* tags_End_Offset_and_Start_Type = "</DataArray>\n<DataArray Name=\"types\" type=\"UInt8\" format=\"ascii\">\n";
105 char* tag_Float_DataArray="<DataArray Name=\"%s\" type=\"Float64\" NumberOfComponents=\"%d\" format=\"ascii\">\n";
106 char *tags_End_Type_And_Cells = "</DataArray>\n</Cells>\n";
107
108 int VTK_QUADRATIC_HEXAHEDRON_INDEX[] = { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 16, 17, 18, 19, 12, 13, 14, 15 };
109 /* if there is no mesh we just return */
110 if (mesh_p==NULL) return;
111
112 my_mpi_rank = mesh_p->Nodes->MPIInfo->rank;
113 mpi_size = mesh_p->Nodes->MPIInfo->size;
114 nDim = mesh_p->Nodes->numDim;
115
116 if (! ( (nDim ==2) || (nDim == 3) ) ) {
117 Finley_setError(IO_ERROR, "saveVTK: spatial dimension 2 or 3 is supported only.");
118 return;
119 }
120 /*************************************************************************************/
121
122 /* open the file and check handle */
123
124 if (mpi_size > 1) {
125 #ifdef PASO_MPI
126 /* Collective Call */
127 #ifdef MPIO_HINTS
128 MPI_Info_create(&mpi_info);
129 /* MPI_Info_set(mpi_info, "striping_unit", "424288"); */
130 /* MPI_Info_set(mpi_info, "striping_factor", "16"); */
131 /* MPI_Info_set(mpi_info, "collective_buffering", "true"); */
132 /* MPI_Info_set(mpi_info, "cb_block_size", "131072"); */
133 /* MPI_Info_set(mpi_info, "cb_buffer_size", "1048567"); */
134 /* MPI_Info_set(mpi_info, "cb_nodes", "8"); */
135 /* MPI_Info_set(mpi_info, "access_style", "write_once, sequential"); */
136
137 /*XFS only */
138 /* MPI_Info_set(mpi_info, "direct_write", "true"); */
139 #endif
140 if ( my_mpi_rank == 0) {
141 if (Paso_fileExists(filename_p)) remove(filename_p);
142 }
143 ierr=MPI_File_open(mesh_p->Nodes->MPIInfo->comm, (char*)filename_p, amode,mpi_info, &mpi_fileHandle_p);
144 if (ierr != MPI_SUCCESS) {
145 perror(filename_p);
146 sprintf(error_msg, "saveVTK: File %s could not be opened for writing in parallel.", filename_p);
147 Finley_setError(IO_ERROR,error_msg);
148 } else {
149 MPI_File_set_view(mpi_fileHandle_p,MPI_DISPLACEMENT_CURRENT,MPI_CHAR, MPI_CHAR, "native" , mpi_info);
150 }
151 #endif
152 } else {
153 fileHandle_p = fopen(filename_p, "w");
154 if (fileHandle_p==NULL) {
155 sprintf(error_msg, "saveVTK: File %s could not be opened for writing.", filename_p);
156 Finley_setError(IO_ERROR,error_msg);
157 }
158 }
159 if (! Paso_MPIInfo_noError(mesh_p->Nodes->MPIInfo) ) return;
160 /*************************************************************************************/
161
162 /* find the mesh type to be written */
163
164 isCellCentered=TMPMEMALLOC(num_data,bool_t);
165 max_len_names=0;
166 if (!Finley_checkPtr(isCellCentered)) {
167 reduced_elements=FALSE;
168 nodetype=FINLEY_UNKNOWN;
169 elementtype=FINLEY_UNKNOWN;
170 for (i_data=0;i_data<num_data;++i_data) {
171 if (! isEmpty(data_pp[i_data])) {
172 switch(getFunctionSpaceType(data_pp[i_data]) ) {
173 case FINLEY_NODES:
174 nodetype = (nodetype == FINLEY_REDUCED_NODES) ? FINLEY_REDUCED_NODES : FINLEY_NODES;
175 isCellCentered[i_data]=FALSE;
176 if (elementtype==FINLEY_UNKNOWN || elementtype==FINLEY_ELEMENTS) {
177 elementtype=FINLEY_ELEMENTS;
178 } else {
179 Finley_setError(TYPE_ERROR,"saveVTK: cannot write given data in single file.");
180 }
181 break;
182 case FINLEY_REDUCED_NODES:
183 nodetype = FINLEY_REDUCED_NODES;
184 isCellCentered[i_data]=FALSE;
185 if (elementtype==FINLEY_UNKNOWN || elementtype==FINLEY_ELEMENTS) {
186 elementtype=FINLEY_ELEMENTS;
187 } else {
188 Finley_setError(TYPE_ERROR,"saveVTK: cannot write given data in single file.");
189 }
190 break;
191 case FINLEY_REDUCED_ELEMENTS:
192 reduced_elements=TRUE;
193 case FINLEY_ELEMENTS:
194 isCellCentered[i_data]=TRUE;
195 if (elementtype==FINLEY_UNKNOWN || elementtype==FINLEY_ELEMENTS) {
196 elementtype=FINLEY_ELEMENTS;
197 } else {
198 Finley_setError(TYPE_ERROR,"saveVTK: cannot write given data in single file.");
199 }
200 break;
201 case FINLEY_REDUCED_FACE_ELEMENTS:
202 reduced_elements=TRUE;
203 case FINLEY_FACE_ELEMENTS:
204 isCellCentered[i_data]=TRUE;
205 if (elementtype==FINLEY_UNKNOWN || elementtype==FINLEY_FACE_ELEMENTS) {
206 elementtype=FINLEY_FACE_ELEMENTS;
207 } else {
208 Finley_setError(TYPE_ERROR,"saveVTK: cannot write given data in single file.");
209 }
210 break;
211 case FINLEY_POINTS:
212 isCellCentered[i_data]=TRUE;
213 if (elementtype==FINLEY_UNKNOWN || elementtype==FINLEY_POINTS) {
214 elementtype=FINLEY_POINTS;
215 } else {
216 Finley_setError(TYPE_ERROR,"saveVTK: cannot write given data in single file.");
217 }
218 break;
219 case FINLEY_REDUCED_CONTACT_ELEMENTS_1:
220 reduced_elements=TRUE;
221 case FINLEY_CONTACT_ELEMENTS_1:
222 isCellCentered[i_data]=TRUE;
223 if (elementtype==FINLEY_UNKNOWN || elementtype==FINLEY_CONTACT_ELEMENTS_1) {
224 elementtype=FINLEY_CONTACT_ELEMENTS_1;
225 } else {
226 Finley_setError(TYPE_ERROR,"saveVTK: cannot write given data in single file.");
227 }
228 break;
229 case FINLEY_REDUCED_CONTACT_ELEMENTS_2:
230 reduced_elements=TRUE;
231 case FINLEY_CONTACT_ELEMENTS_2:
232 isCellCentered[i_data]=TRUE;
233 if (elementtype==FINLEY_UNKNOWN || elementtype==FINLEY_CONTACT_ELEMENTS_1) {
234 elementtype=FINLEY_CONTACT_ELEMENTS_1;
235 } else {
236 Finley_setError(TYPE_ERROR,"saveVTK: cannot write given data in single file.");
237 }
238 break;
239 default:
240 sprintf(error_msg,"saveVTK: unknown function space type %d",getFunctionSpaceType(data_pp[i_data]));
241 Finley_setError(TYPE_ERROR,error_msg);
242 }
243 if (isCellCentered[i_data]) {
244 write_celldata=TRUE;
245 } else {
246 write_pointdata=TRUE;
247 }
248 max_len_names =MAX(max_len_names,(dim_t)strlen(names_p[i_data]));
249 }
250 }
251 nodetype = (nodetype == FINLEY_UNKNOWN) ? FINLEY_NODES : nodetype;
252 }
253 if (Finley_noError()) {
254
255 /***************************************/
256
257 /* select number of points and the mesh component */
258
259 if (nodetype == FINLEY_REDUCED_NODES) {
260 myFirstNode = Finley_NodeFile_getFirstReducedNode(mesh_p->Nodes);
261 myLastNode = Finley_NodeFile_getLastReducedNode(mesh_p->Nodes);
262 globalNumPoints = Finley_NodeFile_getGlobalNumReducedNodes(mesh_p->Nodes);
263 globalNodeIndex= Finley_NodeFile_borrowGlobalReducedNodesIndex(mesh_p->Nodes);
264 } else {
265 myFirstNode = Finley_NodeFile_getFirstNode(mesh_p->Nodes);
266 myLastNode = Finley_NodeFile_getLastNode(mesh_p->Nodes);
267 globalNumPoints = Finley_NodeFile_getGlobalNumNodes(mesh_p->Nodes);
268 globalNodeIndex= Finley_NodeFile_borrowGlobalNodesIndex(mesh_p->Nodes);
269 }
270 myNumPoints = myLastNode - myFirstNode;
271 if (elementtype==FINLEY_UNKNOWN) elementtype=FINLEY_ELEMENTS;
272 switch(elementtype) {
273 case FINLEY_ELEMENTS:
274 elements=mesh_p->Elements;
275 break;
276 case FINLEY_FACE_ELEMENTS:
277 elements=mesh_p->FaceElements;
278 break;
279 case FINLEY_POINTS:
280 elements=mesh_p->Points;
281 break;
282 case FINLEY_CONTACT_ELEMENTS_1:
283 elements=mesh_p->ContactElements;
284 break;
285 }
286 if (elements==NULL) {
287 Finley_setError(SYSTEM_ERROR,"saveVTK: undefined element file");
288 } else {
289 /* map finley element type to VTK element type */
290 numCells = elements->numElements;
291 globalNumCells = Finley_ElementFile_getGlobalNumElements(elements);
292 myNumCells= Finley_ElementFile_getMyNumElements(elements);
293 myFirstCell= Finley_ElementFile_getFirstElement(elements);
294 NN = elements->numNodes;
295 if (nodetype==FINLEY_REDUCED_NODES) {
296 TypeId = elements->LinearReferenceElement->Type->TypeId;
297 if (reduced_elements) {
298 QuadNodes=elements->LinearReferenceElementReducedOrder->QuadNodes;
299 } else {
300 QuadNodes=elements->LinearReferenceElement->QuadNodes;
301 }
302 } else {
303 TypeId = elements->ReferenceElement->Type->TypeId;
304 if (reduced_elements) {
305 QuadNodes=elements->ReferenceElementReducedOrder->QuadNodes;
306 } else {
307 QuadNodes=elements->ReferenceElement->QuadNodes;
308 }
309 }
310 switch(TypeId) {
311 case Point1:
312 case Line2Face:
313 case Line3Face:
314 case Point1_Contact:
315 case Line2Face_Contact:
316 case Line3Face_Contact:
317 numCellFactor=1;
318 cellType = VTK_VERTEX;
319 numVTKNodesPerElement = 1;
320 strcpy(elemTypeStr, "VTK_VERTEX");
321 break;
322
323 case Line2:
324 case Tri3Face:
325 case Rec4Face:
326 case Line2_Contact:
327 case Tri3_Contact:
328 case Tri3Face_Contact:
329 case Rec4Face_Contact:
330 numCellFactor=1;
331 cellType = VTK_LINE;
332 numVTKNodesPerElement = 2;
333 strcpy(elemTypeStr, "VTK_LINE");
334 break;
335
336 case Tri3:
337 case Tet4Face:
338 case Tet4Face_Contact:
339 numCellFactor=1;
340 cellType = VTK_TRIANGLE;
341 numVTKNodesPerElement = 3;
342 strcpy(elemTypeStr, "VTK_TRIANGLE");
343 break;
344
345 case Rec4:
346 case Hex8Face:
347 case Rec4_Contact:
348 case Hex8Face_Contact:
349 numCellFactor=1;
350 cellType = VTK_QUAD;
351 numVTKNodesPerElement = 4;
352 strcpy(elemTypeStr, "VTK_QUAD");
353 break;
354
355 case Rec9:
356 numCellFactor=4;
357 cellType = VTK_QUAD;
358 numVTKNodesPerElement = 4;
359 strcpy(elemTypeStr, "VTK_QUAD");
360 break;
361
362 case Tet4:
363 numCellFactor=1;
364 cellType = VTK_TETRA;
365 numVTKNodesPerElement = 4;
366 strcpy(elemTypeStr, "VTK_TETRA");
367 break;
368
369 case Hex8:
370 numCellFactor=1;
371 cellType = VTK_HEXAHEDRON;
372 numVTKNodesPerElement = 8;
373 strcpy(elemTypeStr, "VTK_HEXAHEDRON");
374 break;
375
376 case Line3:
377 case Tri6Face:
378 case Rec8Face:
379 case Line3_Contact:
380 case Tri6Face_Contact:
381 case Rec8Face_Contact:
382 numCellFactor=1;
383 cellType = VTK_QUADRATIC_EDGE;
384 numVTKNodesPerElement = 3;
385 strcpy(elemTypeStr, "VTK_QUADRATIC_EDGE");
386 break;
387
388 case Tri6:
389 case Tet10Face:
390 case Tri6_Contact:
391 case Tet10Face_Contact:
392 numCellFactor=1;
393 cellType = VTK_QUADRATIC_TRIANGLE;
394 numVTKNodesPerElement = 6;
395 strcpy(elemTypeStr, "VTK_QUADRATIC_TRIANGLE");
396 break;
397
398 case Rec8:
399 case Hex20Face:
400 case Rec8_Contact:
401 case Hex20Face_Contact:
402 numCellFactor=1;
403 cellType = VTK_QUADRATIC_QUAD;
404 numVTKNodesPerElement = 8;
405 strcpy(elemTypeStr, "VTK_QUADRATIC_QUAD");
406 break;
407
408 case Tet10:
409 numCellFactor=1;
410 cellType = VTK_QUADRATIC_TETRA;
411 numVTKNodesPerElement = 10;
412 strcpy(elemTypeStr, "VTK_QUADRATIC_TETRA");
413 break;
414
415 case Hex20:
416 numCellFactor=1;
417 cellType = VTK_QUADRATIC_HEXAHEDRON;
418 numVTKNodesPerElement = 20;
419 strcpy(elemTypeStr, "VTK_QUADRATIC_HEXAHEDRON");
420 break;
421
422 case Hex27:
423 numCellFactor=8;
424 cellType = VTK_HEXAHEDRON;
425 numVTKNodesPerElement = 8;
426 strcpy(elemTypeStr, "VTK_HEXAHEDRON");
427 break;
428
429 default:
430 sprintf(error_msg, "saveVTK: Element type %s is not supported by VTK",elements->ReferenceElement->Type->Name);
431 Finley_setError(VALUE_ERROR,error_msg);
432 }
433 }
434 }
435 /***************************************/
436
437 /***************************************/
438 /* */
439 /* allocate text buffer */
440 /* */
441 max_name_len=0;
442 for (i_data =0 ;i_data<num_data;++i_data) max_name_len=MAX(max_name_len,(dim_t)strlen(names_p[i_data]));
443 len_txt_buffer= strlen(tags_header) + 3 * LEN_PRINTED_INT_FORMAT + (30+3*max_name_len); /* header */
444 if (mpi_size > 1) len_txt_buffer=MAX(len_txt_buffer, myNumPoints * LEN_TMP_BUFFER);
445 if (mpi_size > 1) len_txt_buffer=MAX(len_txt_buffer, numCellFactor*myNumCells*(LEN_PRINTED_INT_FORMAT*numVTKNodesPerElement+1));
446 len_txt_buffer=MAX(len_txt_buffer,200+3*max_len_names);
447 len_txt_buffer=MAX(len_txt_buffer, (dim_t)strlen(tag_Float_DataArray) + LEN_PRINTED_INT_FORMAT + max_len_names);
448 if (mpi_size > 1) len_txt_buffer=MAX(len_txt_buffer, numCellFactor*myNumCells*LEN_PRINTED_FLOAT_TENSOR_FORMAT);
449 if (mpi_size > 1) len_txt_buffer=MAX(len_txt_buffer, myNumPoints*LEN_PRINTED_FLOAT_TENSOR_FORMAT);
450 txt_buffer=TMPMEMALLOC(len_txt_buffer+1,char);
451 Finley_checkPtr(txt_buffer);
452
453 if (Finley_noError()) {
454
455 /* select number of points and the mesh component */
456
457 sprintf(txt_buffer,tags_header,globalNumPoints,numCellFactor*globalNumCells,3);
458
459 if (mpi_size > 1) {
460 if ( my_mpi_rank == 0) {
461 #ifdef PASO_MPI
462 MPI_File_iwrite_shared(mpi_fileHandle_p,txt_buffer,strlen(txt_buffer),MPI_CHAR,&mpi_req);
463 MPI_Wait(&mpi_req,&mpi_status);
464 #endif
465 }
466 } else {
467 fprintf(fileHandle_p,txt_buffer);
468 }
469
470 /* write the nodes */
471
472 if (mpi_size > 1) {
473
474 txt_buffer[0] = '\0';
475 txt_buffer_in_use=0;
476 if (nDim==2) {
477 for (i = 0; i < mesh_p->Nodes->numNodes; i++) {
478 if ( (myFirstNode <= globalNodeIndex[i]) && (globalNodeIndex[i] < myLastNode) ) {
479 sprintf(tmp_buffer,FLOAT_VECTOR_FORMAT,
480 mesh_p->Nodes->Coordinates[INDEX2(0, i, nDim)],
481 mesh_p->Nodes->Coordinates[INDEX2(1, i, nDim)],
482 0.);
483 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use);
484 }
485 }
486 } else {
487 for (i = 0; i < mesh_p->Nodes->numNodes; i++) {
488 if ( (myFirstNode <= globalNodeIndex[i]) && (globalNodeIndex[i] < myLastNode) ) {
489 sprintf(tmp_buffer,FLOAT_VECTOR_FORMAT,
490 mesh_p->Nodes->Coordinates[INDEX2(0, i, nDim)],
491 mesh_p->Nodes->Coordinates[INDEX2(1, i, nDim)],
492 mesh_p->Nodes->Coordinates[INDEX2(2, i, nDim)]);
493 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use);
494 }
495 }
496
497 }
498 #ifdef PASO_MPI
499 if (txt_buffer_in_use==0) { strcpy(txt_buffer, " "); txt_buffer_in_use = 1; } /* avoid zero-length writes */
500 MPI_File_write_ordered(mpi_fileHandle_p, txt_buffer,txt_buffer_in_use, MPI_CHAR, &mpi_status);
501 #endif
502 } else {
503 if (nDim==2) {
504 for (i = 0; i < mesh_p->Nodes->numNodes; i++) {
505 if ( (myFirstNode <= globalNodeIndex[i]) && (globalNodeIndex[i] < myLastNode) ) {
506 fprintf(fileHandle_p,FLOAT_VECTOR_FORMAT,
507 mesh_p->Nodes->Coordinates[INDEX2(0, i, nDim)],
508 mesh_p->Nodes->Coordinates[INDEX2(1, i, nDim)],
509 0.);
510 }
511 }
512 } else {
513 for (i = 0; i < mesh_p->Nodes->numNodes; i++) {
514 if ( (myFirstNode <= globalNodeIndex[i]) && (globalNodeIndex[i] < myLastNode) ) {
515 fprintf(fileHandle_p,FLOAT_VECTOR_FORMAT,
516 mesh_p->Nodes->Coordinates[INDEX2(0, i, nDim)],
517 mesh_p->Nodes->Coordinates[INDEX2(1, i, nDim)],
518 mesh_p->Nodes->Coordinates[INDEX2(2, i, nDim)]);
519 }
520 }
521
522 }
523 }
524
525 /* close the Points and open connectivity */
526
527 if (mpi_size > 1) {
528 if ( my_mpi_rank == 0) {
529 #ifdef PASO_MPI
530 MPI_File_iwrite_shared(mpi_fileHandle_p, tags_End_Points_and_Start_Conn, strlen(tags_End_Points_and_Start_Conn), MPI_CHAR, &mpi_req);
531 MPI_Wait(&mpi_req,&mpi_status);
532 #endif
533 }
534 } else {
535 fprintf(fileHandle_p,tags_End_Points_and_Start_Conn);
536 }
537
538 /* write the cells */
539 if (nodetype == FINLEY_REDUCED_NODES) {
540 node_index=elements->ReferenceElement->Type->linearNodes;
541 } else if (VTK_QUADRATIC_HEXAHEDRON==cellType) {
542 node_index=VTK_QUADRATIC_HEXAHEDRON_INDEX;
543 } else if ( (numVTKNodesPerElement!=NN) && (TypeId!=Rec9) && (TypeId!=Hex27) ) {
544 node_index=elements->ReferenceElement->Type->geoNodes;
545 } else {
546 node_index=NULL;
547 }
548
549 if ( mpi_size > 1) {
550 txt_buffer[0] = '\0';
551 txt_buffer_in_use=0;
552 if (node_index == NULL) {
553 if (TypeId==Rec9) {
554 for (i = 0; i < numCells; i++) {
555 if (elements->Owner[i] == my_mpi_rank) {
556 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(0, i, NN)]]);
557 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(4, i, NN)]]);
558 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(8, i, NN)]]);
559 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(7, i, NN)]]);
560 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use)
561 __STRCAT(txt_buffer,NEWLINE,txt_buffer_in_use)
562
563 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(4, i, NN)]]);
564 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(1, i, NN)]]);
565 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(5, i, NN)]]);
566 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(8, i, NN)]]);
567 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use)
568 __STRCAT(txt_buffer,NEWLINE,txt_buffer_in_use)
569
570 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(7, i, NN)]]);
571 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(8, i, NN)]]);
572 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(6, i, NN)]]);
573 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(3, i, NN)]]);
574 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use)
575 __STRCAT(txt_buffer,NEWLINE,txt_buffer_in_use)
576
577 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(8, i, NN)]]);
578 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(5, i, NN)]]);
579 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(2, i, NN)]]);
580 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(6, i, NN)]]);
581 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use)
582 __STRCAT(txt_buffer,NEWLINE,txt_buffer_in_use)
583 }
584 }
585 } else if (TypeId==Hex27) {
586 for (i = 0; i < numCells; i++) {
587 if (elements->Owner[i] == my_mpi_rank) {
588 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 0, i, NN)]]);
589 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 8, i, NN)]]);
590 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(20, i, NN)]]);
591 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(11, i, NN)]]);
592 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(12, i, NN)]]);
593 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(21, i, NN)]]);
594 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
595 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(24, i, NN)]]);
596 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use)
597 __STRCAT(txt_buffer,NEWLINE,txt_buffer_in_use)
598
599 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 8, i, NN)]]);
600 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 1, i, NN)]]);
601 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 9, i, NN)]]);
602 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(20, i, NN)]]);
603 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(21, i, NN)]]);
604 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(13, i, NN)]]);
605 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(22, i, NN)]]);
606 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
607 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use)
608 __STRCAT(txt_buffer,NEWLINE,txt_buffer_in_use)
609
610 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(11, i, NN)]]);
611 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(20, i, NN)]]);
612 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(10, i, NN)]]);
613 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 3, i, NN)]]);
614 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(24, i, NN)]]);
615 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
616 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(23, i, NN)]]);
617 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(15, i, NN)]]);
618 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use)
619 __STRCAT(txt_buffer,NEWLINE,txt_buffer_in_use)
620
621 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(20, i, NN)]]);
622 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 9, i, NN)]]);
623 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 2, i, NN)]]);
624 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(10, i, NN)]]);
625 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
626 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(22, i, NN)]]);
627 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(14, i, NN)]]);
628 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(23, i, NN)]]);
629 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use)
630 __STRCAT(txt_buffer,NEWLINE,txt_buffer_in_use)
631
632 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(12, i, NN)]]);
633 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(21, i, NN)]]);
634 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
635 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(24, i, NN)]]);
636 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 4, i, NN)]]);
637 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(16, i, NN)]]);
638 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(25, i, NN)]]);
639 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(19, i, NN)]]);
640 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use)
641 __STRCAT(txt_buffer,NEWLINE,txt_buffer_in_use)
642
643 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(21, i, NN)]]);
644 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(13, i, NN)]]);
645 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(22, i, NN)]]);
646 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
647 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(16, i, NN)]]);
648 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 5, i, NN)]]);
649 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(17, i, NN)]]);
650 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(25, i, NN)]]);
651 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use)
652 __STRCAT(txt_buffer,NEWLINE,txt_buffer_in_use)
653
654 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(24, i, NN)]]);
655 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
656 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(23, i, NN)]]);
657 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(15, i, NN)]]);
658 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(19, i, NN)]]);
659 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(25, i, NN)]]);
660 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(18, i, NN)]]);
661 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 7, i, NN)]]);
662 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use)
663 __STRCAT(txt_buffer,NEWLINE,txt_buffer_in_use)
664
665 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
666 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(22, i, NN)]]);
667 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(14, i, NN)]]);
668 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(23, i, NN)]]);
669 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(25, i, NN)]]);
670 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(17, i, NN)]]);
671 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 6, i, NN)]]);
672 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(18, i, NN)]]);
673 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use)
674 __STRCAT(txt_buffer,NEWLINE,txt_buffer_in_use)
675 }
676 }
677 } else {
678 for (i = 0; i < numCells; i++) {
679 if (elements->Owner[i] == my_mpi_rank) {
680 for (j = 0; j < numVTKNodesPerElement; j++) {
681 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(j, i, NN)]]);
682 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use)
683 }
684 __STRCAT(txt_buffer,NEWLINE,txt_buffer_in_use)
685 }
686 }
687 }
688 } else {
689 for (i = 0; i < numCells; i++) {
690 if (elements->Owner[i] == my_mpi_rank) {
691 for (j = 0; j < numVTKNodesPerElement; j++) {
692 sprintf(tmp_buffer,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(node_index[j], i, NN)]]);
693 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use)
694 }
695 __STRCAT(txt_buffer,NEWLINE,txt_buffer_in_use)
696 }
697 }
698 }
699 #ifdef PASO_MPI
700 if (txt_buffer_in_use==0) { strcpy(txt_buffer, " "); txt_buffer_in_use = 1; } /* avoid zero-length writes */
701 MPI_File_write_ordered(mpi_fileHandle_p,txt_buffer,txt_buffer_in_use, MPI_CHAR, &mpi_status);
702 #endif
703 } else {
704 if (node_index == NULL) {
705 if (TypeId==Rec9) {
706 for (i = 0; i < numCells; i++) {
707 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(0, i, NN)]]);
708 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(4, i, NN)]]);
709 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(8, i, NN)]]);
710 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(7, i, NN)]]);
711 fprintf(fileHandle_p,NEWLINE);
712
713 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(4, i, NN)]]);
714 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(1, i, NN)]]);
715 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(5, i, NN)]]);
716 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(8, i, NN)]]);
717 fprintf(fileHandle_p,NEWLINE);
718
719 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(7, i, NN)]]);
720 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(8, i, NN)]]);
721 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(6, i, NN)]]);
722 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(3, i, NN)]]);
723 fprintf(fileHandle_p,NEWLINE);
724
725 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(8, i, NN)]]);
726 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(5, i, NN)]]);
727 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(2, i, NN)]]);
728 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(6, i, NN)]]);
729 fprintf(fileHandle_p,NEWLINE);
730 }
731
732 } else if (TypeId==Hex27) {
733 for (i = 0; i < numCells; i++) {
734 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 0, i, NN)]]);
735 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 8, i, NN)]]);
736 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(20, i, NN)]]);
737 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(11, i, NN)]]);
738 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(12, i, NN)]]);
739 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(21, i, NN)]]);
740 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
741 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(24, i, NN)]]);
742 fprintf(fileHandle_p,NEWLINE);
743
744 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 8, i, NN)]]);
745 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 1, i, NN)]]);
746 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 9, i, NN)]]);
747 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(20, i, NN)]]);
748 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(21, i, NN)]]);
749 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(13, i, NN)]]);
750 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(22, i, NN)]]);
751 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
752 fprintf(fileHandle_p,NEWLINE);
753
754 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(11, i, NN)]]);
755 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(20, i, NN)]]);
756 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(10, i, NN)]]);
757 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 3, i, NN)]]);
758 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(24, i, NN)]]);
759 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
760 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(23, i, NN)]]);
761 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(15, i, NN)]]);
762 fprintf(fileHandle_p,NEWLINE);
763
764 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(20, i, NN)]]);
765 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 9, i, NN)]]);
766 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 2, i, NN)]]);
767 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(10, i, NN)]]);
768 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
769 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(22, i, NN)]]);
770 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(14, i, NN)]]);
771 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(23, i, NN)]]);
772 fprintf(fileHandle_p,NEWLINE);
773
774 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(12, i, NN)]]);
775 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(21, i, NN)]]);
776 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
777 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(24, i, NN)]]);
778 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 4, i, NN)]]);
779 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(16, i, NN)]]);
780 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(25, i, NN)]]);
781 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(19, i, NN)]]);
782 fprintf(fileHandle_p,NEWLINE);
783
784 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(21, i, NN)]]);
785 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(13, i, NN)]]);
786 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(22, i, NN)]]);
787 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
788 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(16, i, NN)]]);
789 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 5, i, NN)]]);
790 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(17, i, NN)]]);
791 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(25, i, NN)]]);
792 fprintf(fileHandle_p,NEWLINE);
793
794 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(24, i, NN)]]);
795 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
796 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(23, i, NN)]]);
797 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(15, i, NN)]]);
798 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(19, i, NN)]]);
799 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(25, i, NN)]]);
800 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(18, i, NN)]]);
801 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 7, i, NN)]]);
802 fprintf(fileHandle_p,NEWLINE);
803
804 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(26, i, NN)]]);
805 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(22, i, NN)]]);
806 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(14, i, NN)]]);
807 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(23, i, NN)]]);
808 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(25, i, NN)]]);
809 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(17, i, NN)]]);
810 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2( 6, i, NN)]]);
811 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(18, i, NN)]]);
812 fprintf(fileHandle_p,NEWLINE);
813 }
814 } else {
815 for (i = 0; i < numCells; i++) {
816 for (j = 0; j < numVTKNodesPerElement; j++) {
817 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(j, i, NN)]]);
818 }
819 fprintf(fileHandle_p,NEWLINE);
820 }
821 }
822 } else {
823 for (i = 0; i < numCells; i++) {
824 for (j = 0; j < numVTKNodesPerElement; j++) {
825 fprintf(fileHandle_p,INT_FORMAT,globalNodeIndex[elements->Nodes[INDEX2(node_index[j], i, NN)]]);
826 }
827 fprintf(fileHandle_p,NEWLINE);
828 }
829 }
830
831 }
832 /* finalize the connection and start the offset section */
833 if (mpi_size > 1) {
834 if( my_mpi_rank == 0) {
835 #ifdef PASO_MPI
836 MPI_File_iwrite_shared(mpi_fileHandle_p,tags_End_Conn_and_Start_Offset,strlen(tags_End_Conn_and_Start_Offset),MPI_CHAR,&mpi_req);
837 MPI_Wait(&mpi_req,&mpi_status);
838 #endif
839 }
840 } else {
841 fprintf(fileHandle_p,tags_End_Conn_and_Start_Offset);
842 }
843
844 /* write the offsets */
845
846 if ( mpi_size > 1) {
847 txt_buffer[0] = '\0';
848 txt_buffer_in_use=0;
849 for (i=numVTKNodesPerElement*(myFirstCell*numCellFactor+1); i<=(myFirstCell+myNumCells)*numVTKNodesPerElement*numCellFactor; i+=numVTKNodesPerElement) {
850 sprintf(tmp_buffer, INT_NEWLINE_FORMAT, i);
851 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use);
852 }
853 #ifdef PASO_MPI
854 if (txt_buffer_in_use==0) { strcpy(txt_buffer, " "); txt_buffer_in_use = 1; } /* avoid zero-length writes */
855 MPI_File_write_ordered(mpi_fileHandle_p,txt_buffer,txt_buffer_in_use, MPI_CHAR, &mpi_status);
856 #endif
857 } else {
858 for (i=numVTKNodesPerElement; i<=numCells*numVTKNodesPerElement*numCellFactor; i+=numVTKNodesPerElement) {
859 fprintf(fileHandle_p, INT_NEWLINE_FORMAT, i);
860 }
861
862 }
863 /* finalize the offset section and start the type section */
864 if ( mpi_size > 1) {
865 if ( my_mpi_rank == 0) {
866 #ifdef PASO_MPI
867 MPI_File_iwrite_shared(mpi_fileHandle_p,tags_End_Offset_and_Start_Type,strlen(tags_End_Offset_and_Start_Type),MPI_CHAR,&mpi_req);
868 MPI_Wait(&mpi_req,&mpi_status);
869 #endif
870 }
871 } else {
872 fprintf(fileHandle_p,tags_End_Offset_and_Start_Type);
873 }
874 /* write element type */
875 sprintf(tmp_buffer, INT_NEWLINE_FORMAT, cellType);
876 if ( mpi_size > 1) {
877 txt_buffer[0] = '\0';
878 txt_buffer_in_use=0;
879 for (i=0; i<numCells*numCellFactor; i++) __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use);
880 #ifdef PASO_MPI
881 if (txt_buffer_in_use==0) { strcpy(txt_buffer, " "); txt_buffer_in_use = 1; } /* avoid zero-length writes */
882 MPI_File_write_ordered(mpi_fileHandle_p,txt_buffer,txt_buffer_in_use, MPI_CHAR, &mpi_status);
883 #endif
884 } else {
885 for (i=0; i<numCells*numCellFactor; i++) fprintf(fileHandle_p, tmp_buffer);
886 }
887 /* finalize cell information */
888 if ( mpi_size > 1) {
889 if ( my_mpi_rank == 0) {
890 #ifdef PASO_MPI
891 MPI_File_iwrite_shared(mpi_fileHandle_p,tags_End_Type_And_Cells,strlen(tags_End_Type_And_Cells),MPI_CHAR,&mpi_req);
892 MPI_Wait(&mpi_req,&mpi_status);
893 #endif
894 }
895 } else {
896 fprintf(fileHandle_p,tags_End_Type_And_Cells);
897 }
898 }
899
900 /* Write cell data */
901 if (write_celldata && Finley_noError()) {
902 /* mark the active data arrays */
903 txt_buffer[0] = '\0';
904 set_scalar=FALSE,set_vector=FALSE, set_tensor=FALSE;
905 strcat(txt_buffer, "<CellData");
906 for (i_data =0 ;i_data<num_data;++i_data) {
907 if (! isEmpty(data_pp[i_data]) && isCellCentered[i_data]) {
908 /* if the rank == 0: --> scalar data */
909 /* if the rank == 1: --> vector data */
910 /* if the rank == 2: --> tensor data */
911
912 switch(getDataPointRank(data_pp[i_data])) {
913 case 0:
914 if (! set_scalar) {
915 strcat(txt_buffer," Scalars=\"");
916 strcat(txt_buffer,names_p[i_data]);
917 strcat(txt_buffer,"\"");
918 set_scalar=TRUE;
919 }
920 break;
921 case 1:
922 if (! set_vector) {
923 strcat(txt_buffer," Vectors=\"");
924 strcat(txt_buffer,names_p[i_data]);
925 strcat(txt_buffer,"\"");
926 set_vector=TRUE;
927 }
928 break;
929 case 2:
930 if (! set_tensor) {
931 strcat(txt_buffer," Tensors=\"");
932 strcat(txt_buffer,names_p[i_data]);
933 strcat(txt_buffer,"\"");
934 set_tensor=TRUE;
935 }
936 break;
937 default:
938 sprintf(error_msg, "saveVTK: data %s: Vtk can't handle objects with rank greater than 2.",names_p[i_data]);
939 Finley_setError(VALUE_ERROR,error_msg);
940 return;
941 }
942 }
943 }
944 strcat(txt_buffer, ">\n");
945 if ( mpi_size > 1) {
946 if ( my_mpi_rank == 0) {
947 #ifdef PASO_MPI
948 MPI_File_iwrite_shared(mpi_fileHandle_p,txt_buffer,strlen(txt_buffer),MPI_CHAR,&mpi_req);
949 MPI_Wait(&mpi_req,&mpi_status);
950 #endif
951 }
952 } else {
953 fprintf(fileHandle_p,txt_buffer);
954 }
955 /* write the arrays */
956 for (i_data =0 ;i_data<num_data;++i_data) {
957 if (! isEmpty(data_pp[i_data]) && isCellCentered[i_data]) {
958 txt_buffer[0] = '\0';
959 txt_buffer_in_use=0;
960 numPointsPerSample=getNumDataPointsPerSample(data_pp[i_data]);
961 rank = getDataPointRank(data_pp[i_data]);
962 nComp = getDataPointSize(data_pp[i_data]);
963 nCompReqd=1; /* the number of components mpi_required by vtk */
964 shape=0;
965 if (rank == 0) {
966 nCompReqd = 1;
967 } else if (rank == 1) {
968 shape=getDataPointShape(data_pp[i_data], 0);
969 if (shape>3) {
970 Finley_setError(VALUE_ERROR, "saveVTK: rank 1 object must have less then 4 components");
971 }
972 nCompReqd = 3;
973 } else {
974 shape=getDataPointShape(data_pp[i_data], 0);
975 if (shape>3 || shape != getDataPointShape(data_pp[i_data], 1)) {
976 Finley_setError(VALUE_ERROR, "saveVTK: rank 2 object must have less then 4x4 components and must have a square shape");
977 }
978 nCompReqd = 9;
979 }
980 if (Finley_noError()) {
981 sprintf(txt_buffer,tag_Float_DataArray,names_p[i_data], nCompReqd);
982 if ( mpi_size > 1) {
983 if ( my_mpi_rank == 0) {
984 #ifdef PASO_MPI
985 MPI_File_iwrite_shared(mpi_fileHandle_p,txt_buffer,strlen(txt_buffer),MPI_CHAR,&mpi_req);
986 MPI_Wait(&mpi_req,&mpi_status);
987 #endif
988 }
989 } else {
990 fprintf(fileHandle_p,txt_buffer);
991 }
992
993 for (i=0; i<numCells; i++) {
994 if (elements->Owner[i] == my_mpi_rank) {
995 values = getSampleData(data_pp[i_data], i);
996 for (l=0; l< numCellFactor;++l) {
997 /* averaging over the number of points in the sample */
998 if (isExpanded(data_pp[i_data])) {
999 for (k=0; k<MIN(nComp,NCOMP_MAX); k++) sampleAvg[k]=0;
1000 hits=0;
1001 for (j=0; j<numPointsPerSample; j++) {
1002 hits_old=hits;
1003 if (TypeId==Rec9) {
1004 switch(l) {
1005 case 0:
1006 if (INSIDE_2D(QuadNodes[2*j],QuadNodes[2*j+1],0.25,0.25,0.25)) hits++;
1007 break;
1008 case 1:
1009 if (INSIDE_2D(QuadNodes[2*j],QuadNodes[2*j+1],0.75,0.25,0.25)) hits++;
1010 break;
1011 case 2:
1012 if (INSIDE_2D(QuadNodes[2*j],QuadNodes[2*j+1],0.25,0.75,0.25)) hits++;
1013 break;
1014 case 3:
1015 if (INSIDE_2D(QuadNodes[2*j],QuadNodes[2*j+1],0.75,0.75,0.25)) hits++;
1016 break;
1017 }
1018 } else if (TypeId==Hex27) {
1019 switch(l) {
1020 case 0:
1021 if (INSIDE_3D(QuadNodes[3*j],QuadNodes[3*j+1],QuadNodes[3*j+2],0.25,0.25,0.25,0.25)) hits++;
1022 break;
1023 case 1:
1024 if (INSIDE_3D(QuadNodes[3*j],QuadNodes[3*j+1],QuadNodes[3*j+2],0.75,0.25,0.25,0.25)) hits++;
1025 break;
1026 case 2:
1027 if (INSIDE_3D(QuadNodes[3*j],QuadNodes[3*j+1],QuadNodes[3*j+2],0.25,0.75,0.25,0.25)) hits++;
1028 break;
1029 case 3:
1030 if (INSIDE_3D(QuadNodes[3*j],QuadNodes[3*j+1],QuadNodes[3*j+2],0.75,0.75,0.25,0.25)) hits++;
1031 break;
1032 case 4:
1033 if (INSIDE_3D(QuadNodes[3*j],QuadNodes[3*j+1],QuadNodes[3*j+2],0.25,0.25,0.75,0.25)) hits++;
1034 break;
1035 case 5:
1036 if (INSIDE_3D(QuadNodes[3*j],QuadNodes[3*j+1],QuadNodes[3*j+2],0.75,0.25,0.75,0.25)) hits++;
1037 break;
1038 case 6:
1039 if (INSIDE_3D(QuadNodes[3*j],QuadNodes[3*j+1],QuadNodes[3*j+2],0.25,0.75,0.75,0.25)) hits++;
1040 break;
1041 case 7:
1042 if (INSIDE_3D(QuadNodes[3*j],QuadNodes[3*j+1],QuadNodes[3*j+2],0.75,0.75,0.75,0.25)) hits++;
1043 break;
1044 }
1045 } else {
1046 hits++;
1047 }
1048 if (hits_old<hits) for (k=0; k<MIN(nComp,NCOMP_MAX); k++) {
1049 sampleAvg[k] += values[INDEX2(k,j,nComp)];
1050 }
1051 }
1052 for (k=0; k<MIN(nComp,NCOMP_MAX); k++) sampleAvg[k] /=MAX(hits,1);
1053 } else {
1054 for (k=0; k<MIN(nComp,NCOMP_MAX); k++) sampleAvg[k] = values[k];
1055 }
1056 /* if the number of required components is more than the number
1057 * of actual components, pad with zeros
1058 */
1059 /* probably only need to get shape of first element */
1060 /* write the data different ways for scalar, vector and tensor */
1061 if (nCompReqd == 1) {
1062 sprintf(tmp_buffer,FLOAT_SCALAR_FORMAT,sampleAvg[0]);
1063 } else if (nCompReqd == 3) {
1064 if (shape==1) {
1065 sprintf(tmp_buffer,FLOAT_VECTOR_FORMAT,sampleAvg[0],0.,0.);
1066 } else if (shape==2) {
1067 sprintf(tmp_buffer,FLOAT_VECTOR_FORMAT,sampleAvg[0],sampleAvg[1],0.);
1068 } else if (shape==3) {
1069 sprintf(tmp_buffer,FLOAT_VECTOR_FORMAT,sampleAvg[0],sampleAvg[1],sampleAvg[2]);
1070 }
1071 } else if (nCompReqd == 9) {
1072 if (shape==1) {
1073 sprintf(tmp_buffer,FLOAT_TENSOR_FORMAT,sampleAvg[0],0.,0.,
1074 0.,0.,0.,
1075 0.,0.,0.);
1076 } else if (shape==2) {
1077 sprintf(tmp_buffer,FLOAT_TENSOR_FORMAT,sampleAvg[0],sampleAvg[1],0.,
1078 sampleAvg[2],sampleAvg[3],0.,
1079 0.,0.,0.);
1080 } else if (shape==3) {
1081 sprintf(tmp_buffer,FLOAT_TENSOR_FORMAT,sampleAvg[0],sampleAvg[1],sampleAvg[2],
1082 sampleAvg[3],sampleAvg[4],sampleAvg[5],
1083 sampleAvg[6],sampleAvg[7],sampleAvg[8]);
1084 }
1085 }
1086 /* this needs a bit mor work!!! */
1087 if ( mpi_size > 1) {
1088 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use);
1089 } else {
1090 fprintf(fileHandle_p,tmp_buffer);
1091 }
1092 }
1093 }
1094 }
1095 if ( mpi_size > 1) {
1096 #ifdef PASO_MPI
1097 if (txt_buffer_in_use==0) { strcpy(txt_buffer, " "); txt_buffer_in_use = 1; } /* avoid zero-length writes */
1098 MPI_File_write_ordered(mpi_fileHandle_p,txt_buffer,txt_buffer_in_use, MPI_CHAR, &mpi_status);
1099 #endif
1100 if ( my_mpi_rank == 0) {
1101 #ifdef PASO_MPI
1102 MPI_File_iwrite_shared(mpi_fileHandle_p,tag_End_DataArray,strlen(tag_End_DataArray),MPI_CHAR,&mpi_req);
1103 MPI_Wait(&mpi_req,&mpi_status);
1104 #endif
1105 }
1106 } else {
1107 fprintf(fileHandle_p,tag_End_DataArray);
1108 }
1109 }
1110 }
1111 }
1112 if ( mpi_size > 1) {
1113 if ( my_mpi_rank == 0) {
1114 #ifdef PASO_MPI
1115 MPI_File_iwrite_shared(mpi_fileHandle_p,tag_End_CellData,strlen(tag_End_CellData),MPI_CHAR,&mpi_req);
1116 MPI_Wait(&mpi_req,&mpi_status);
1117 #endif
1118 }
1119 } else {
1120 fprintf(fileHandle_p,tag_End_CellData);
1121 }
1122 }
1123 /* point data */
1124 if (write_pointdata && Finley_noError()) {
1125 /* mark the active data arrays */
1126 set_scalar=FALSE,set_vector=FALSE, set_tensor=FALSE;
1127 txt_buffer[0] = '\0';
1128 strcat(txt_buffer, "<PointData");
1129 for (i_data =0 ;i_data<num_data;++i_data) {
1130 if (! isEmpty(data_pp[i_data]) && !isCellCentered[i_data]) {
1131 /* if the rank == 0: --> scalar data */
1132 /* if the rank == 1: --> vector data */
1133 /* if the rank == 2: --> tensor data */
1134
1135 switch(getDataPointRank(data_pp[i_data])) {
1136 case 0:
1137 if (! set_scalar) {
1138 strcat(txt_buffer," Scalars=\"");
1139 strcat(txt_buffer,names_p[i_data]);
1140 strcat(txt_buffer,"\"");
1141 set_scalar=TRUE;
1142 }
1143 break;
1144 case 1:
1145 if (! set_vector) {
1146 strcat(txt_buffer," Vectors=\"");
1147 strcat(txt_buffer,names_p[i_data]);
1148 strcat(txt_buffer,"\"");
1149 set_vector=TRUE;
1150 }
1151 break;
1152 case 2:
1153 if (! set_tensor) {
1154 strcat(txt_buffer," Tensors=\"");
1155 strcat(txt_buffer,names_p[i_data]);
1156 strcat(txt_buffer,"\"");
1157 set_tensor=TRUE;
1158 }
1159 break;
1160 default:
1161 sprintf(error_msg, "saveVTK: data %s: Vtk can't handle objects with rank greater than 2.",names_p[i_data]);
1162 Finley_setError(VALUE_ERROR,error_msg);
1163 return;
1164 }
1165 }
1166 }
1167 strcat(txt_buffer, ">\n");
1168 if ( mpi_size > 1) {
1169 if ( my_mpi_rank == 0) {
1170 #ifdef PASO_MPI
1171 MPI_File_iwrite_shared(mpi_fileHandle_p,txt_buffer,strlen(txt_buffer),MPI_CHAR,&mpi_req);
1172 MPI_Wait(&mpi_req,&mpi_status);
1173 #endif
1174 }
1175 } else {
1176 fprintf(fileHandle_p,txt_buffer);
1177 }
1178 /* write the arrays */
1179 for (i_data =0 ;i_data<num_data;++i_data) {
1180 if (! isEmpty(data_pp[i_data]) && !isCellCentered[i_data]) {
1181 txt_buffer[0] = '\0';
1182 txt_buffer_in_use=0;
1183 numPointsPerSample=getNumDataPointsPerSample(data_pp[i_data]);
1184 rank = getDataPointRank(data_pp[i_data]);
1185 nComp = getDataPointSize(data_pp[i_data]);
1186 if (getFunctionSpaceType(data_pp[i_data]) == FINLEY_REDUCED_NODES) {
1187 nodeMapping=mesh_p->Nodes->reducedNodesMapping;
1188 } else {
1189 nodeMapping=mesh_p->Nodes->nodesMapping;
1190 }
1191 nCompReqd=1; /* the number of components mpi_required by vtk */
1192 shape=0;
1193 if (rank == 0) {
1194 nCompReqd = 1;
1195 } else if (rank == 1) {
1196 shape=getDataPointShape(data_pp[i_data], 0);
1197 if (shape>3) {
1198 Finley_setError(VALUE_ERROR, "saveVTK: rank 1 object must have less then 4 components");
1199 }
1200 nCompReqd = 3;
1201 } else {
1202 shape=getDataPointShape(data_pp[i_data], 0);
1203 if (shape>3 || shape != getDataPointShape(data_pp[i_data], 1)) {
1204 Finley_setError(VALUE_ERROR, "saveVTK: rank 2 object must have less then 4x4 components and must have a square shape");
1205 }
1206 nCompReqd = 9;
1207 }
1208 if (Finley_noError()) {
1209 sprintf(txt_buffer,tag_Float_DataArray,names_p[i_data], nCompReqd);
1210 if ( mpi_size > 1) {
1211 if ( my_mpi_rank == 0) {
1212 #ifdef PASO_MPI
1213 MPI_File_iwrite_shared(mpi_fileHandle_p,txt_buffer,strlen(txt_buffer),MPI_CHAR,&mpi_req);
1214 MPI_Wait(&mpi_req,&mpi_status);
1215 #endif
1216 }
1217 } else {
1218 fprintf(fileHandle_p,txt_buffer);
1219 }
1220 for (i=0; i<mesh_p->Nodes->numNodes; i++) {
1221 k=globalNodeIndex[i];
1222 if ( (myFirstNode <= k) && (k < myLastNode) ) {
1223 values = getSampleData(data_pp[i_data], nodeMapping->target[i]);
1224 /* if the number of mpi_required components is more than the number
1225 * of actual components, pad with zeros
1226 */
1227 /* probably only need to get shape of first element */
1228 /* write the data different ways for scalar, vector and tensor */
1229 if (nCompReqd == 1) {
1230 sprintf(tmp_buffer,FLOAT_SCALAR_FORMAT,values[0]);
1231 } else if (nCompReqd == 3) {
1232 if (shape==1) {
1233 sprintf(tmp_buffer,FLOAT_VECTOR_FORMAT,values[0],0.,0.);
1234 } else if (shape==2) {
1235 sprintf(tmp_buffer,FLOAT_VECTOR_FORMAT,values[0],values[1],0.);
1236 } else if (shape==3) {
1237 sprintf(tmp_buffer,FLOAT_VECTOR_FORMAT,values[0],values[1],values[2]);
1238 }
1239 } else if (nCompReqd == 9) {
1240 if (shape==1) {
1241 sprintf(tmp_buffer,FLOAT_TENSOR_FORMAT,values[0],0.,0.,
1242 0.,0.,0.,
1243 0.,0.,0.);
1244 } else if (shape==2) {
1245 sprintf(tmp_buffer,FLOAT_TENSOR_FORMAT,values[0],values[1],0.,
1246 values[2],values[3],0.,
1247 0.,0.,0.);
1248 } else if (shape==3) {
1249 sprintf(tmp_buffer,FLOAT_TENSOR_FORMAT,values[0],values[1],values[2],
1250 values[3],values[4],values[5],
1251 values[6],values[7],values[8]);
1252 }
1253 }
1254 if ( mpi_size > 1) {
1255 __STRCAT(txt_buffer,tmp_buffer,txt_buffer_in_use);
1256 } else {
1257 fprintf(fileHandle_p,tmp_buffer);
1258 }
1259 }
1260 }
1261 if ( mpi_size > 1) {
1262 #ifdef PASO_MPI
1263 if (txt_buffer_in_use==0) { strcpy(txt_buffer, " "); txt_buffer_in_use = 1; } /* avoid zero-length writes */
1264 MPI_File_write_ordered(mpi_fileHandle_p,txt_buffer,txt_buffer_in_use, MPI_CHAR, &mpi_status);
1265 #endif
1266 if ( my_mpi_rank == 0) {
1267 #ifdef PASO_MPI
1268 MPI_File_iwrite_shared(mpi_fileHandle_p,tag_End_DataArray,strlen(tag_End_DataArray),MPI_CHAR,&mpi_req);
1269 MPI_Wait(&mpi_req,&mpi_status);
1270 #endif
1271 }
1272 } else {
1273 fprintf(fileHandle_p,tag_End_DataArray);
1274 }
1275 }
1276 }
1277 }
1278 if ( mpi_size > 1) {
1279 if ( my_mpi_rank == 0) {
1280 #ifdef PASO_MPI
1281 MPI_File_iwrite_shared(mpi_fileHandle_p,tag_End_PointData,strlen(tag_End_PointData),MPI_CHAR,&mpi_req);
1282 MPI_Wait(&mpi_req,&mpi_status);
1283 #endif
1284 }
1285 } else {
1286 fprintf(fileHandle_p,tag_End_PointData);
1287 }
1288 }
1289 if (Finley_noError()) {
1290 if ( mpi_size > 1) {
1291 if ( my_mpi_rank == 0) {
1292 #ifdef PASO_MPI
1293 MPI_File_iwrite_shared(mpi_fileHandle_p,footer,strlen(footer),MPI_CHAR,&mpi_req);
1294 MPI_Wait(&mpi_req,&mpi_status);
1295 #ifdef MPIO_HINTS
1296 MPI_Info_free(&mpi_info);
1297 #undef MPIO_HINTS
1298 #endif
1299 #endif
1300 }
1301 #ifdef PASO_MPI
1302 MPI_File_close(&mpi_fileHandle_p);
1303 #endif
1304 } else {
1305 fprintf(fileHandle_p,footer);
1306 fclose(fileHandle_p);
1307 }
1308 }
1309 TMPMEMFREE(isCellCentered);
1310 TMPMEMFREE(txt_buffer);
1311 return;
1312 #else
1313 /* Don't kill the job if saveVTK() doesn't work */
1314 fprintf(stderr, "\n\nsaveVTK warning: VTK is not available\n\n\n");
1315 #endif
1316 }

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26