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

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

Parent Directory Parent Directory | Revision Log Revision Log | View Patch Patch

revision 2141 by caltinay, Tue Dec 9 04:32:32 2008 UTC revision 2748 by gross, Tue Nov 17 07:32:59 2009 UTC
# Line 1  Line 1 
1    
2  /*******************************************************  /*******************************************************
3  *  *
4  * Copyright (c) 2003-2008 by University of Queensland  * Copyright (c) 2003-2009 by University of Queensland
5  * Earth Systems Science Computational Center (ESSCC)  * Earth Systems Science Computational Center (ESSCC)
6  * http://www.uq.edu.au/esscc  * http://www.uq.edu.au/esscc
7  *  *
# Line 22  Line 22 
22  #include "paso/PasoUtil.h"  #include "paso/PasoUtil.h"
23    
24  #define INT_FORMAT "%d "  #define INT_FORMAT "%d "
25  #define LEN_INT_FORMAT (9+1)  #define LEN_INT_FORMAT (unsigned int)(9+2)
26  #define INT_NEWLINE_FORMAT "%d\n"  #define INT_NEWLINE_FORMAT "%d\n"
27  #define SCALAR_FORMAT "%12.6e\n"  #define SCALAR_FORMAT "%12.6e\n"
28  #define VECTOR_FORMAT "%12.6e %12.6e %12.6e\n"  #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"  #define TENSOR_FORMAT "%12.6e %12.6e %12.6e %12.6e %12.6e %12.6e %12.6e %12.6e %12.6e\n"
30  #define LEN_TENSOR_FORMAT (9*(12+1)+1)  /* strlen("-1.234567e+789 ") == 15 */
31    #define LEN_TENSOR_FORMAT (unsigned int)(9*15+2)
32  #define NEWLINE "\n"  #define NEWLINE "\n"
33  #define LEN_TMP_BUFFER LEN_TENSOR_FORMAT+(MAX_numNodes*LEN_INT_FORMAT+1)+2  #define LEN_TMP_BUFFER LEN_TENSOR_FORMAT+(MAX_numNodes*LEN_INT_FORMAT+1)+2
34  #define NCOMP_MAX 9  #define NCOMP_MAX (unsigned int)9
35    
36  #define __STRCAT(dest, chunk, dest_in_use) \  #define __STRCAT(dest, chunk, dest_in_use) \
37  do {\  do {\
# Line 39  do {\ Line 40  do {\
40  } while(0)  } while(0)
41    
42  #ifdef PASO_MPI  #ifdef PASO_MPI
   
43  /* writes buffer to file catching the empty buffer case which causes problems  /* writes buffer to file catching the empty buffer case which causes problems
44   * with some MPI versions */   * with some MPI versions */
45  #define MPI_WRITE_ORDERED(BUF, LEN) \  #define MPI_WRITE_ORDERED(BUF) \
46  do {\  do {\
47      if (LEN==0) { strcpy(BUF, " "); LEN=1; }\      int LLEN=0; \
48      MPI_File_write_ordered(mpi_fileHandle_p, BUF, LEN, MPI_CHAR, &mpi_status);\      LLEN=(int) strlen(BUF); \
49        if (LLEN==0) { strcpy(BUF, ""); LLEN=0; }\
50        MPI_File_write_ordered(mpi_fileHandle_p, BUF, LLEN, MPI_CHAR, &mpi_status);\
51  } while(0)  } while(0)
52    
53  /* writes buffer to file on master only */  /* writes buffer to file on master only */
54  #define MPI_RANK0_WRITE_SHARED(BUF) \  #define MPI_RANK0_WRITE_SHARED(BUF) \
55  do {\  do {\
56        int LLEN=0; \
57      if (my_mpi_rank == 0) {\      if (my_mpi_rank == 0) {\
58          MPI_File_iwrite_shared(mpi_fileHandle_p, BUF, strlen(BUF), MPI_CHAR, &mpi_req);\          LLEN=(int) strlen(BUF); \
59            if (LLEN==0) { strcpy(BUF,""); LLEN=0; }\
60            MPI_File_iwrite_shared(mpi_fileHandle_p, BUF, LLEN, MPI_CHAR, &mpi_req);\
61          MPI_Wait(&mpi_req, &mpi_status);\          MPI_Wait(&mpi_req, &mpi_status);\
62      }\      }\
63  } while(0)  } while(0)
# Line 75  void create_MPIInfo(MPI_Info& info) Line 80  void create_MPIInfo(MPI_Info& info)
80    
81  #else  #else
82    
83  #define MPI_WRITE_ORDERED(A, B)  #define MPI_WRITE_ORDERED(A)
84  #define MPI_RANK0_WRITE_SHARED(A)  #define MPI_RANK0_WRITE_SHARED(A)
85    
86  #endif /* PASO_MPI */  #endif /* PASO_MPI */
# Line 90  int nodeInQuadrant(const double *coords, Line 95  int nodeInQuadrant(const double *coords,
95  #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_) )  #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_) )
96    
97      int ret;      int ret;
98      if (type == Rec9) {      if ( (type == Rec9) || (type == Rec9Macro) ) {
99          if (q==0)          if (q==0)
100              ret = INSIDE_2D(coords[2*idx], coords[2*idx+1], 0.25, 0.25, 0.25);              ret = INSIDE_2D(coords[2*idx], coords[2*idx+1], 0.25, 0.25, 0.25);
101          else if (q==1)          else if (q==1)
# Line 101  int nodeInQuadrant(const double *coords, Line 106  int nodeInQuadrant(const double *coords,
106              ret = INSIDE_2D(coords[2*idx], coords[2*idx+1], 0.75, 0.75, 0.25);              ret = INSIDE_2D(coords[2*idx], coords[2*idx+1], 0.75, 0.75, 0.25);
107          else          else
108              ret = 0;              ret = 0;
109      } else if (type == Hex27) {      } else if ((type == Hex27) || (type == Hex27Macro) ){
110          if (q==0)          if (q==0)
111              ret = INSIDE_3D(coords[3*idx], coords[3*idx+1], coords[3*idx+2],              ret = INSIDE_3D(coords[3*idx], coords[3*idx+1], coords[3*idx+2],
112                      0.25, 0.25, 0.25, 0.25);                      0.25, 0.25, 0.25, 0.25);
# Line 138  void Finley_Mesh_saveVTK(const char *fil Line 143  void Finley_Mesh_saveVTK(const char *fil
143                           Finley_Mesh *mesh_p,                           Finley_Mesh *mesh_p,
144                           const dim_t num_data,                           const dim_t num_data,
145                           char **names_p,                           char **names_p,
146                           escriptDataC **data_pp)                           escriptDataC **data_pp,
147                             const char* metadata,
148                             const char*metadata_schema)
149  {  {
150  #ifdef PASO_MPI  #ifdef PASO_MPI
151      MPI_File mpi_fileHandle_p;      MPI_File mpi_fileHandle_p;
# Line 163  void Finley_Mesh_saveVTK(const char *fil Line 170  void Finley_Mesh_saveVTK(const char *fil
170      int mpi_size, i, j, l;      int mpi_size, i, j, l;
171      int cellType=0, nodeType=FINLEY_NODES, elementType=FINLEY_UNKNOWN;      int cellType=0, nodeType=FINLEY_NODES, elementType=FINLEY_UNKNOWN;
172      Finley_ElementFile *elements = NULL;      Finley_ElementFile *elements = NULL;
173      ElementTypeId typeId = NoType;      ElementTypeId typeId = NoRef;
174    
175      const char *vtkHeader = \      const char *vtkHeader = \
176        "<?xml version=\"1.0\"?>\n" \        "<?xml version=\"1.0\"?>\n" \
177        "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\">\n" \        "<VTKFile type=\"UnstructuredGrid\" version=\"0.1\"%s%s>\n%s%s" \
178        "<UnstructuredGrid>\n" \        "<UnstructuredGrid>\n" \
179        "<Piece NumberOfPoints=\"%d\" NumberOfCells=\"%d\">\n" \        "<Piece NumberOfPoints=\"%d\" NumberOfCells=\"%d\">\n" \
180        "<Points>\n" \        "<Points>\n" \
# Line 205  void Finley_Mesh_saveVTK(const char *fil Line 212  void Finley_Mesh_saveVTK(const char *fil
212      my_mpi_rank = mesh_p->Nodes->MPIInfo->rank;      my_mpi_rank = mesh_p->Nodes->MPIInfo->rank;
213      mpi_size = mesh_p->Nodes->MPIInfo->size;      mpi_size = mesh_p->Nodes->MPIInfo->size;
214    
215      /************************************************************************/      /************************************************************************
216      /* open the file and check handle */       * open the file and check handle *
217         */
218      if (mpi_size > 1) {      if (mpi_size > 1) {
219  #ifdef PASO_MPI  #ifdef PASO_MPI
220          const int amode = MPI_MODE_CREATE|MPI_MODE_WRONLY|MPI_MODE_UNIQUE_OPEN;          const int amode = MPI_MODE_CREATE|MPI_MODE_WRONLY|MPI_MODE_UNIQUE_OPEN;
# Line 221  void Finley_Mesh_saveVTK(const char *fil Line 228  void Finley_Mesh_saveVTK(const char *fil
228              sprintf(errorMsg, "saveVTK: File %s could not be opened for writing in parallel.", filename_p);              sprintf(errorMsg, "saveVTK: File %s could not be opened for writing in parallel.", filename_p);
229              Finley_setError(IO_ERROR, errorMsg);              Finley_setError(IO_ERROR, errorMsg);
230          } else {          } else {
231              MPI_File_set_view(mpi_fileHandle_p, MPI_DISPLACEMENT_CURRENT,              ierr=MPI_File_set_view(mpi_fileHandle_p,MPI_DISPLACEMENT_CURRENT,
232                      MPI_CHAR, MPI_CHAR, "native", mpi_info);                      MPI_CHAR, MPI_CHAR, "native", mpi_info);
233          }          }
234  #endif /* PASO_MPI */  #endif /* PASO_MPI */
# Line 368  void Finley_Mesh_saveVTK(const char *fil Line 375  void Finley_Mesh_saveVTK(const char *fil
375              myFirstCell = Finley_ElementFile_getFirstElement(elements);              myFirstCell = Finley_ElementFile_getFirstElement(elements);
376              NN = elements->numNodes;              NN = elements->numNodes;
377              if (nodeType==FINLEY_REDUCED_NODES) {              if (nodeType==FINLEY_REDUCED_NODES) {
378                  typeId = elements->LinearReferenceElement->Type->TypeId;          typeId = elements->referenceElementSet->referenceElement->LinearType->TypeId;
379                  if (hasReducedElements) {                  if (hasReducedElements) {
380                      quadNodes_p=elements->LinearReferenceElementReducedOrder->QuadNodes;                      quadNodes_p=elements->referenceElementSet->referenceElementReducedQuadrature->Parametrization->QuadNodes;
381                  } else {                  } else {
382                      quadNodes_p=elements->LinearReferenceElement->QuadNodes;                      quadNodes_p=elements->referenceElementSet->referenceElement->Parametrization->QuadNodes;
383                  }                  }
384              } else {              } else {
385                  typeId = elements->ReferenceElement->Type->TypeId;                  typeId = elements->referenceElementSet->referenceElement->Type->TypeId;
386                  if (hasReducedElements)                  if (hasReducedElements) {
387                      quadNodes_p=elements->ReferenceElementReducedOrder->QuadNodes;                      quadNodes_p=elements->referenceElementSet->referenceElementReducedQuadrature->Parametrization->QuadNodes;
388                  else                  } else {
389                      quadNodes_p=elements->ReferenceElement->QuadNodes;                      quadNodes_p=elements->referenceElementSet->referenceElement->Parametrization->QuadNodes;
390                    }
391              }              }
392              switch (typeId) {              switch (typeId) {
393                  case Point1:                  case Point1:
# Line 418  void Finley_Mesh_saveVTK(const char *fil Line 426  void Finley_Mesh_saveVTK(const char *fil
426                      numVTKNodesPerElement = 4;                      numVTKNodesPerElement = 4;
427                  break;                  break;
428    
429                  case Rec9:                  case Rec9Macro:
430            case Rec9:
431                      numCellFactor = 4;                      numCellFactor = 4;
432                      cellType = VTK_QUAD;                      cellType = VTK_QUAD;
433                      numVTKNodesPerElement = 4;                      numVTKNodesPerElement = 4;
# Line 435  void Finley_Mesh_saveVTK(const char *fil Line 444  void Finley_Mesh_saveVTK(const char *fil
444                  break;                  break;
445    
446                  case Line3:                  case Line3:
447                    case Line3Macro:
448                  case Tri6Face:                  case Tri6Face:
449                  case Rec8Face:                  case Rec8Face:
450                  case Line3_Contact:                  case Line3_Contact:
# Line 445  void Finley_Mesh_saveVTK(const char *fil Line 455  void Finley_Mesh_saveVTK(const char *fil
455                  break;                  break;
456    
457                  case Tri6:                  case Tri6:
458            case Tri6Macro:    
459                  case Tet10Face:                  case Tet10Face:
460                  case Tri6_Contact:                  case Tri6_Contact:
461                  case Tet10Face_Contact:                  case Tet10Face_Contact:
# Line 460  void Finley_Mesh_saveVTK(const char *fil Line 471  void Finley_Mesh_saveVTK(const char *fil
471                      numVTKNodesPerElement = 8;                      numVTKNodesPerElement = 8;
472                  break;                  break;
473    
474            case Tet10Macro:
475                  case Tet10:                  case Tet10:
476                      cellType = VTK_QUADRATIC_TETRA;                      cellType = VTK_QUADRATIC_TETRA;
477                      numVTKNodesPerElement = 10;                      numVTKNodesPerElement = 10;
# Line 470  void Finley_Mesh_saveVTK(const char *fil Line 482  void Finley_Mesh_saveVTK(const char *fil
482                      numVTKNodesPerElement = 20;                      numVTKNodesPerElement = 20;
483                  break;                  break;
484    
485            case Hex27Macro:
486                  case Hex27:                  case Hex27:
487                      numCellFactor = 8;                      numCellFactor = 8;
488                      cellType = VTK_HEXAHEDRON;                      cellType = VTK_HEXAHEDRON;
# Line 477  void Finley_Mesh_saveVTK(const char *fil Line 490  void Finley_Mesh_saveVTK(const char *fil
490                  break;                  break;
491    
492                  default:                  default:
493                      sprintf(errorMsg, "saveVTK: Element type %s is not supported by VTK.", elements->ReferenceElement->Type->Name);                      sprintf(errorMsg, "saveVTK: Element type %s is not supported by VTK.", elements->referenceElementSet->referenceElement->Type->Name);
494                      Finley_setError(VALUE_ERROR, errorMsg);                      Finley_setError(VALUE_ERROR, errorMsg);
495              }              }
496          }          }
# Line 485  void Finley_Mesh_saveVTK(const char *fil Line 498  void Finley_Mesh_saveVTK(const char *fil
498    
499      /* allocate enough memory for text buffer */      /* allocate enough memory for text buffer */
500    
501      txtBufferSize = strlen(vtkHeader) + 3*LEN_INT_FORMAT + (30+3*maxNameLen);      txtBufferSize = strlen(vtkHeader) + 3*LEN_INT_FORMAT + (30+3*maxNameLen)+strlen(metadata)+strlen(metadata_schema);
   
502      if (mpi_size > 1) {      if (mpi_size > 1) {
503          txtBufferSize = MAX(txtBufferSize, myNumPoints * LEN_TMP_BUFFER);          txtBufferSize = MAX(txtBufferSize, myNumPoints * LEN_TMP_BUFFER);
504          txtBufferSize = MAX(txtBufferSize, numCellFactor * myNumCells *          txtBufferSize = MAX(txtBufferSize, numCellFactor * myNumCells *
# Line 506  void Finley_Mesh_saveVTK(const char *fil Line 518  void Finley_Mesh_saveVTK(const char *fil
518      if (Finley_noError()) {      if (Finley_noError()) {
519          const index_t *nodeIndex;          const index_t *nodeIndex;
520          if (FINLEY_REDUCED_NODES == nodeType) {          if (FINLEY_REDUCED_NODES == nodeType) {
521              nodeIndex = elements->ReferenceElement->Type->linearNodes;              nodeIndex = elements->referenceElementSet->referenceElement->Type->linearNodes;
522          } else if (Rec9 == typeId) {          } else if ( (Rec9 == typeId) || (Rec9Macro == typeId) ) {
523              nodeIndex = VTK_REC9_INDEX;              nodeIndex = VTK_REC9_INDEX;
524          } else if (Hex20 == typeId) {          } else if (Hex20 == typeId) {
525              nodeIndex = VTK_HEX20_INDEX;              nodeIndex = VTK_HEX20_INDEX;
526          } else if (Hex27 == typeId) {          } else if ( (Hex27 == typeId) || (Hex27Macro == typeId) ){
527              nodeIndex = VTK_HEX27_INDEX;              nodeIndex = VTK_HEX27_INDEX;
528          } else if (numVTKNodesPerElement != NN) {          } else if (numVTKNodesPerElement  !=  elements->referenceElementSet->referenceElement->Type->numNodes) {
529              nodeIndex = elements->ReferenceElement->Type->geoNodes;              nodeIndex = elements->referenceElementSet->referenceElement->Type->relevantGeoNodes;
530          } else {          } else {
531              nodeIndex = NULL;              nodeIndex = NULL;
532          }          }
533    
534          sprintf(txtBuffer, vtkHeader, globalNumPoints,          if (strlen(metadata)>0) {
535                  numCellFactor*globalNumCells, 3);             if (strlen(metadata_schema)>0) {
536                  sprintf(txtBuffer, vtkHeader," ",metadata_schema,metadata,"\n",globalNumPoints, numCellFactor*globalNumCells, 3);
537               } else {
538                  sprintf(txtBuffer, vtkHeader,"","",metadata,"\n",globalNumPoints, numCellFactor*globalNumCells, 3);
539               }
540            } else {
541               if (strlen(metadata_schema)>0) {
542                  sprintf(txtBuffer, vtkHeader," ",metadata_schema,"","",globalNumPoints, numCellFactor*globalNumCells, 3);
543               } else {
544                  sprintf(txtBuffer, vtkHeader,"","","","",globalNumPoints, numCellFactor*globalNumCells, 3);
545               }
546            }
547    
548          if (mpi_size > 1) {          if (mpi_size > 1) {
549              /* write the nodes */              /* write the nodes */
# Line 548  void Finley_Mesh_saveVTK(const char *fil Line 571  void Finley_Mesh_saveVTK(const char *fil
571                      }                      }
572                  }                  }
573              } /* nDim */              } /* nDim */
574              MPI_WRITE_ORDERED(txtBuffer, txtBufferInUse);              MPI_WRITE_ORDERED(txtBuffer);
575    
576              /* write the cells */              /* write the cells */
577              MPI_RANK0_WRITE_SHARED(tags_End_Points_and_Start_Conn);              MPI_RANK0_WRITE_SHARED(tags_End_Points_and_Start_Conn);
# Line 578  void Finley_Mesh_saveVTK(const char *fil Line 601  void Finley_Mesh_saveVTK(const char *fil
601                      }                      }
602                  }                  }
603              } /* nodeIndex */              } /* nodeIndex */
604              MPI_WRITE_ORDERED(txtBuffer, txtBufferInUse);              MPI_WRITE_ORDERED(txtBuffer);
605    
606              /* write the offsets */              /* write the offsets */
607              MPI_RANK0_WRITE_SHARED(tags_End_Conn_and_Start_Offset);              MPI_RANK0_WRITE_SHARED(tags_End_Conn_and_Start_Offset);
# Line 590  void Finley_Mesh_saveVTK(const char *fil Line 613  void Finley_Mesh_saveVTK(const char *fil
613                  sprintf(tmpBuffer, INT_NEWLINE_FORMAT, i);                  sprintf(tmpBuffer, INT_NEWLINE_FORMAT, i);
614                  __STRCAT(txtBuffer, tmpBuffer, txtBufferInUse);                  __STRCAT(txtBuffer, tmpBuffer, txtBufferInUse);
615              }              }
616              MPI_WRITE_ORDERED(txtBuffer, txtBufferInUse);              MPI_WRITE_ORDERED(txtBuffer);
617    
618              /* write element type */              /* write element type */
619              sprintf(tmpBuffer, INT_NEWLINE_FORMAT, cellType);              sprintf(tmpBuffer, INT_NEWLINE_FORMAT, cellType);
# Line 602  void Finley_Mesh_saveVTK(const char *fil Line 625  void Finley_Mesh_saveVTK(const char *fil
625              {              {
626                  __STRCAT(txtBuffer, tmpBuffer, txtBufferInUse);                  __STRCAT(txtBuffer, tmpBuffer, txtBufferInUse);
627              }              }
628              MPI_WRITE_ORDERED(txtBuffer, txtBufferInUse);              MPI_WRITE_ORDERED(txtBuffer);
629              /* finalize cell information */              /* finalize cell information */
630              strcpy(txtBuffer, "</DataArray>\n</Cells>\n");              strcpy(txtBuffer, "</DataArray>\n</Cells>\n");
631              MPI_RANK0_WRITE_SHARED(txtBuffer);              MPI_RANK0_WRITE_SHARED(txtBuffer);
# Line 762  void Finley_Mesh_saveVTK(const char *fil Line 785  void Finley_Mesh_saveVTK(const char *fil
785                  txtBufferInUse = 0;                  txtBufferInUse = 0;
786                  for (i=0; i<numCells; i++) {                  for (i=0; i<numCells; i++) {
787                      if (elements->Owner[i] == my_mpi_rank) {                      if (elements->Owner[i] == my_mpi_rank) {
788                          double *values = getSampleData(data_pp[dataIdx], i);                          void* sampleBuffer=allocSampleBuffer(data_pp[dataIdx]);
789                            __const double *values = getSampleDataRO(data_pp[dataIdx], i,sampleBuffer);
790                          for (l = 0; l < numCellFactor; l++) {                          for (l = 0; l < numCellFactor; l++) {
791                              double sampleAvg[NCOMP_MAX];                              double sampleAvg[NCOMP_MAX];
792                              dim_t nCompUsed = MIN(nComp, NCOMP_MAX);                              dim_t nCompUsed = MIN(nComp, NCOMP_MAX);
793    
794                              /* average over number of points in the sample */                              /* average over number of points in the sample */
795                              if (isExpanded(data_pp[dataIdx])) {                              if (isExpanded(data_pp[dataIdx])) {
796                                  dim_t hits=0, hits_old;                                  dim_t hits=0;
797                                  for (k=0; k<nCompUsed; k++) sampleAvg[k]=0;                                  for (k=0; k<nCompUsed; k++) sampleAvg[k]=0;
798                                  for (j=0; j<numPointsPerSample; j++) {                                  for (j=0; j<numPointsPerSample; j++) {
                                     hits_old=hits;  
799                                      if (nodeInQuadrant(quadNodes_p, typeId, j, l)) {                                      if (nodeInQuadrant(quadNodes_p, typeId, j, l)) {
800                                          hits++;                                          hits++;
801                                          for (k=0; k<nCompUsed; k++) {                                          for (k=0; k<nCompUsed; k++) {
# Line 828  void Finley_Mesh_saveVTK(const char *fil Line 851  void Finley_Mesh_saveVTK(const char *fil
851                                  fputs(tmpBuffer, fileHandle_p);                                  fputs(tmpBuffer, fileHandle_p);
852                              }                              }
853                          } /* for l (numCellFactor) */                          } /* for l (numCellFactor) */
854                            freeSampleBuffer(sampleBuffer);
855                      } /* if I am the owner */                      } /* if I am the owner */
856                  } /* for i (numCells) */                  } /* for i (numCells) */
857    
858                  if ( mpi_size > 1) {                  if ( mpi_size > 1) {
859                      MPI_WRITE_ORDERED(txtBuffer, txtBufferInUse);                      MPI_WRITE_ORDERED(txtBuffer);
860                      MPI_RANK0_WRITE_SHARED(tag_End_DataArray);                      MPI_RANK0_WRITE_SHARED(tag_End_DataArray);
861                  } else {                  } else {
862                      fputs(tag_End_DataArray, fileHandle_p);                      fputs(tag_End_DataArray, fileHandle_p);
# Line 943  void Finley_Mesh_saveVTK(const char *fil Line 967  void Finley_Mesh_saveVTK(const char *fil
967                  for (i=0; i<mesh_p->Nodes->numNodes; i++) {                  for (i=0; i<mesh_p->Nodes->numNodes; i++) {
968                      k = globalNodeIndex[i];                      k = globalNodeIndex[i];
969                      if ( (myFirstNode <= k) && (k < myLastNode) ) {                      if ( (myFirstNode <= k) && (k < myLastNode) ) {
970                          double *values = getSampleData(data_pp[dataIdx], nodeMapping->target[i]);                          void* sampleBuffer=allocSampleBuffer(data_pp[dataIdx]);
971                            __const double *values = getSampleDataRO(data_pp[dataIdx], nodeMapping->target[i], sampleBuffer);
972                          /* if the number of mpi_required components is more than                          /* if the number of mpi_required components is more than
973                           * the number of actual components, pad with zeros.                           * the number of actual components, pad with zeros.
974                           * Probably only need to get shape of first element */                           * Probably only need to get shape of first element */
# Line 983  void Finley_Mesh_saveVTK(const char *fil Line 1008  void Finley_Mesh_saveVTK(const char *fil
1008                          } else {                          } else {
1009                              fputs(tmpBuffer, fileHandle_p);                              fputs(tmpBuffer, fileHandle_p);
1010                          }                          }
1011                            freeSampleBuffer(sampleBuffer);                 /* no-one needs values anymore */
1012                      } /* if this is my node */                      } /* if this is my node */
1013                  } /* for i (numNodes) */                  } /* for i (numNodes) */
1014    
1015                  if ( mpi_size > 1) {                  if ( mpi_size > 1) {
1016                      MPI_WRITE_ORDERED(txtBuffer, txtBufferInUse);                      MPI_WRITE_ORDERED(txtBuffer);
1017                      MPI_RANK0_WRITE_SHARED(tag_End_DataArray);                      MPI_RANK0_WRITE_SHARED(tag_End_DataArray);
1018                  } else {                  } else {
1019                      fputs(tag_End_DataArray, fileHandle_p);                      fputs(tag_End_DataArray, fileHandle_p);
# Line 1015  void Finley_Mesh_saveVTK(const char *fil Line 1041  void Finley_Mesh_saveVTK(const char *fil
1041      if ( mpi_size > 1) {      if ( mpi_size > 1) {
1042  #ifdef PASO_MPI  #ifdef PASO_MPI
1043          MPI_File_close(&mpi_fileHandle_p);          MPI_File_close(&mpi_fileHandle_p);
1044            MPI_Barrier(mesh_p->Nodes->MPIInfo->comm);
1045  #endif  #endif
1046      } else {      } else {
1047          fclose(fileHandle_p);          fclose(fileHandle_p);

Legend:
Removed from v.2141  
changed lines
  Added in v.2748

  ViewVC Help
Powered by ViewVC 1.1.26