1 #include <petsc/private/dmpleximpl.h> 2 #include <../src/sys/classes/viewer/impls/vtk/vtkvimpl.h> 3 4 typedef struct { 5 PetscInt nvertices; 6 PetscInt ncells; 7 PetscInt nconn; /* number of entries in cell->vertex connectivity array */ 8 } PieceInfo; 9 10 #if defined(PETSC_USE_REAL_SINGLE) || defined(PETSC_USE_REAL___FP16) 11 /* output in float if single or half precision in memory */ 12 static const char precision[] = "Float32"; 13 typedef float PetscVTUReal; 14 #define MPIU_VTUREAL MPI_FLOAT 15 #elif defined(PETSC_USE_REAL_DOUBLE) || defined(PETSC_USE_REAL___FLOAT128) 16 /* output in double if double or quad precision in memory */ 17 static const char precision[] = "Float64"; 18 typedef double PetscVTUReal; 19 #define MPIU_VTUREAL MPI_DOUBLE 20 #else 21 static const char precision[] = "UnknownPrecision"; 22 typedef PetscReal PetscVTUReal; 23 #define MPIU_VTUREAL MPIU_REAL 24 #endif 25 26 static PetscErrorCode TransferWrite(PetscViewer viewer,FILE *fp,PetscMPIInt srank,PetscMPIInt root,const void *send,void *recv,PetscMPIInt count,MPI_Datatype mpidatatype,PetscMPIInt tag) 27 { 28 PetscMPIInt rank; 29 PetscErrorCode ierr; 30 MPI_Comm comm; 31 32 PetscFunctionBegin; 33 ierr = PetscObjectGetComm((PetscObject)viewer,&comm);CHKERRQ(ierr); 34 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 35 36 if (rank == srank && rank != root) { 37 ierr = MPI_Send((void*)send,count,mpidatatype,root,tag,comm);CHKERRQ(ierr); 38 } else if (rank == root) { 39 const void *buffer; 40 if (root == srank) { /* self */ 41 buffer = send; 42 } else { 43 MPI_Status status; 44 PetscMPIInt nrecv; 45 ierr = MPI_Recv(recv,count,mpidatatype,srank,tag,comm,&status);CHKERRQ(ierr); 46 ierr = MPI_Get_count(&status,mpidatatype,&nrecv);CHKERRQ(ierr); 47 if (count != nrecv) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Array size mismatch"); 48 buffer = recv; 49 } 50 ierr = PetscViewerVTKFWrite(viewer,fp,buffer,count,mpidatatype);CHKERRQ(ierr); 51 } 52 PetscFunctionReturn(0); 53 } 54 55 static PetscErrorCode DMPlexGetVTKConnectivity(DM dm, PetscBool localized, PieceInfo *piece,PetscVTKInt **oconn,PetscVTKInt **ooffsets,PetscVTKType **otypes) 56 { 57 PetscErrorCode ierr; 58 PetscSection coordSection; 59 PetscVTKInt *conn,*offsets; 60 PetscVTKType *types; 61 PetscInt dim,vStart,vEnd,cStart,cEnd,pStart,pEnd,cellHeight,numLabelCells,hasLabel,c,v,countcell,countconn; 62 63 PetscFunctionBegin; 64 ierr = PetscMalloc3(piece->nconn,&conn,piece->ncells,&offsets,piece->ncells,&types);CHKERRQ(ierr); 65 ierr = DMGetCoordinateSection(dm, &coordSection);CHKERRQ(ierr); 66 ierr = DMGetDimension(dm,&dim);CHKERRQ(ierr); 67 ierr = DMPlexGetChart(dm,&pStart,&pEnd);CHKERRQ(ierr); 68 ierr = DMPlexGetVTKCellHeight(dm, &cellHeight);CHKERRQ(ierr); 69 ierr = DMPlexGetHeightStratum(dm, cellHeight, &cStart, &cEnd);CHKERRQ(ierr); 70 ierr = DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd);CHKERRQ(ierr); 71 ierr = DMGetStratumSize(dm, "vtk", 1, &numLabelCells);CHKERRQ(ierr); 72 hasLabel = numLabelCells > 0 ? PETSC_TRUE : PETSC_FALSE; 73 74 countcell = 0; 75 countconn = 0; 76 for (c = cStart; c < cEnd; ++c) { 77 PetscInt nverts,dof = 0,celltype,startoffset,nC=0; 78 79 if (hasLabel) { 80 PetscInt value; 81 82 ierr = DMGetLabelValue(dm, "vtk", c, &value);CHKERRQ(ierr); 83 if (value != 1) continue; 84 } 85 startoffset = countconn; 86 if (localized) { 87 ierr = PetscSectionGetDof(coordSection, c, &dof);CHKERRQ(ierr); 88 } 89 if (!dof) { 90 PetscInt *closure = NULL; 91 PetscInt closureSize; 92 93 ierr = DMPlexGetTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr); 94 for (v = 0; v < closureSize*2; v += 2) { 95 if ((closure[v] >= vStart) && (closure[v] < vEnd)) { 96 if (!localized) conn[countconn++] = closure[v] - vStart; 97 else conn[countconn++] = startoffset + nC; 98 ++nC; 99 } 100 } 101 ierr = DMPlexRestoreTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr); 102 } else { 103 for (nC = 0; nC < dof/dim; nC++) conn[countconn++] = startoffset + nC; 104 } 105 ierr = DMPlexInvertCell(dim, nC, &conn[countconn-nC]);CHKERRQ(ierr); 106 107 offsets[countcell] = countconn; 108 109 nverts = countconn - startoffset; 110 ierr = DMPlexVTKGetCellType_Internal(dm,dim,nverts,&celltype);CHKERRQ(ierr); 111 112 types[countcell] = celltype; 113 countcell++; 114 } 115 if (countcell != piece->ncells) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Inconsistent cell count"); 116 if (countconn != piece->nconn) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Inconsistent connectivity count"); 117 *oconn = conn; 118 *ooffsets = offsets; 119 *otypes = types; 120 PetscFunctionReturn(0); 121 } 122 123 /* 124 Write all fields that have been provided to the viewer 125 Multi-block XML format with binary appended data. 126 */ 127 PetscErrorCode DMPlexVTKWriteAll_VTU(DM dm,PetscViewer viewer) 128 { 129 MPI_Comm comm; 130 PetscSection coordSection; 131 PetscViewer_VTK *vtk = (PetscViewer_VTK*)viewer->data; 132 PetscViewerVTKObjectLink link; 133 FILE *fp; 134 PetscMPIInt rank,size,tag; 135 PetscErrorCode ierr; 136 PetscInt dimEmbed,cellHeight,cStart,cEnd,vStart,vEnd,numLabelCells,hasLabel,c,v,r,i; 137 PetscBool localized; 138 PieceInfo piece,*gpiece = NULL; 139 void *buffer = NULL; 140 const char *byte_order = PetscBinaryBigEndian() ? "BigEndian" : "LittleEndian"; 141 PetscInt loops_per_scalar; 142 143 PetscFunctionBegin; 144 ierr = PetscObjectGetComm((PetscObject)dm,&comm);CHKERRQ(ierr); 145 #if defined(PETSC_USE_COMPLEX) 146 loops_per_scalar = 2; 147 #else 148 loops_per_scalar = 1; 149 #endif 150 ierr = MPI_Comm_size(comm,&size);CHKERRQ(ierr); 151 ierr = MPI_Comm_rank(comm,&rank);CHKERRQ(ierr); 152 tag = ((PetscObject)viewer)->tag; 153 154 ierr = PetscFOpen(comm,vtk->filename,"wb",&fp);CHKERRQ(ierr); 155 ierr = PetscFPrintf(comm,fp,"<?xml version=\"1.0\"?>\n");CHKERRQ(ierr); 156 ierr = PetscFPrintf(comm,fp,"<VTKFile type=\"UnstructuredGrid\" version=\"0.1\" byte_order=\"%s\">\n", byte_order);CHKERRQ(ierr); 157 ierr = PetscFPrintf(comm,fp," <UnstructuredGrid>\n");CHKERRQ(ierr); 158 159 ierr = DMGetCoordinateDim(dm, &dimEmbed);CHKERRQ(ierr); 160 ierr = DMPlexGetVTKCellHeight(dm, &cellHeight);CHKERRQ(ierr); 161 ierr = DMPlexGetHeightStratum(dm, cellHeight, &cStart, &cEnd);CHKERRQ(ierr); 162 ierr = DMPlexGetDepthStratum(dm, 0, &vStart, &vEnd);CHKERRQ(ierr); 163 ierr = DMGetStratumSize(dm, "vtk", 1, &numLabelCells);CHKERRQ(ierr); 164 ierr = DMGetCoordinatesLocalized(dm, &localized);CHKERRQ(ierr); 165 ierr = DMGetCoordinateSection(dm, &coordSection);CHKERRQ(ierr); 166 167 hasLabel = numLabelCells > 0 ? PETSC_TRUE : PETSC_FALSE; 168 piece.nvertices = 0; 169 piece.ncells = 0; 170 piece.nconn = 0; 171 if (!localized) piece.nvertices = vEnd - vStart; 172 for (c = cStart; c < cEnd; ++c) { 173 PetscInt dof = 0; 174 if (hasLabel) { 175 PetscInt value; 176 177 ierr = DMGetLabelValue(dm, "vtk", c, &value);CHKERRQ(ierr); 178 if (value != 1) continue; 179 } 180 if (localized) { 181 ierr = PetscSectionGetDof(coordSection, c, &dof);CHKERRQ(ierr); 182 } 183 if (!dof) { 184 PetscInt *closure = NULL; 185 PetscInt closureSize; 186 187 ierr = DMPlexGetTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr); 188 for (v = 0; v < closureSize*2; v += 2) { 189 if ((closure[v] >= vStart) && (closure[v] < vEnd)) { 190 piece.nconn++; 191 if (localized) piece.nvertices++; 192 } 193 } 194 ierr = DMPlexRestoreTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr); 195 } else { 196 piece.nvertices += dof/dimEmbed; 197 piece.nconn += dof/dimEmbed; 198 } 199 piece.ncells++; 200 } 201 if (!rank) {ierr = PetscMalloc1(size,&gpiece);CHKERRQ(ierr);} 202 ierr = MPI_Gather((PetscInt*)&piece,sizeof(piece)/sizeof(PetscInt),MPIU_INT,(PetscInt*)gpiece,sizeof(piece)/sizeof(PetscInt),MPIU_INT,0,comm);CHKERRQ(ierr); 203 204 /* 205 * Write file header 206 */ 207 if (!rank) { 208 PetscInt boffset = 0; 209 210 for (r=0; r<size; r++) { 211 ierr = PetscFPrintf(PETSC_COMM_SELF,fp," <Piece NumberOfPoints=\"%D\" NumberOfCells=\"%D\">\n",gpiece[r].nvertices,gpiece[r].ncells);CHKERRQ(ierr); 212 /* Coordinate positions */ 213 ierr = PetscFPrintf(PETSC_COMM_SELF,fp," <Points>\n");CHKERRQ(ierr); 214 ierr = PetscFPrintf(PETSC_COMM_SELF,fp," <DataArray type=\"%s\" Name=\"Position\" NumberOfComponents=\"3\" format=\"appended\" offset=\"%D\" />\n",precision,boffset);CHKERRQ(ierr); 215 boffset += gpiece[r].nvertices*3*sizeof(PetscVTUReal) + sizeof(int); 216 ierr = PetscFPrintf(PETSC_COMM_SELF,fp," </Points>\n");CHKERRQ(ierr); 217 /* Cell connectivity */ 218 ierr = PetscFPrintf(PETSC_COMM_SELF,fp," <Cells>\n");CHKERRQ(ierr); 219 ierr = PetscFPrintf(PETSC_COMM_SELF,fp," <DataArray type=\"Int32\" Name=\"connectivity\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",boffset);CHKERRQ(ierr); 220 boffset += gpiece[r].nconn*sizeof(PetscInt) + sizeof(int); 221 ierr = PetscFPrintf(PETSC_COMM_SELF,fp," <DataArray type=\"Int32\" Name=\"offsets\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",boffset);CHKERRQ(ierr); 222 boffset += gpiece[r].ncells*sizeof(PetscInt) + sizeof(int); 223 ierr = PetscFPrintf(PETSC_COMM_SELF,fp," <DataArray type=\"UInt8\" Name=\"types\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",boffset);CHKERRQ(ierr); 224 boffset += gpiece[r].ncells*sizeof(unsigned char) + sizeof(int); 225 ierr = PetscFPrintf(PETSC_COMM_SELF,fp," </Cells>\n");CHKERRQ(ierr); 226 227 /* 228 * Cell Data headers 229 */ 230 ierr = PetscFPrintf(PETSC_COMM_SELF,fp," <CellData>\n");CHKERRQ(ierr); 231 ierr = PetscFPrintf(PETSC_COMM_SELF,fp," <DataArray type=\"Int32\" Name=\"Rank\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",boffset);CHKERRQ(ierr); 232 boffset += gpiece[r].ncells*sizeof(int) + sizeof(int); 233 /* all the vectors */ 234 for (link=vtk->link; link; link=link->next) { 235 Vec X = (Vec)link->vec; 236 DM dmX = NULL; 237 PetscInt bs,nfields,field; 238 const char *vecname = ""; 239 PetscSection section; 240 if ((link->ft != PETSC_VTK_CELL_FIELD) && (link->ft != PETSC_VTK_CELL_VECTOR_FIELD)) continue; 241 if (((PetscObject)X)->name || link != vtk->link) { /* If the object is already named, use it. If it is past the first link, name it to disambiguate. */ 242 ierr = PetscObjectGetName((PetscObject)X,&vecname);CHKERRQ(ierr); 243 } 244 ierr = VecGetDM(X, &dmX);CHKERRQ(ierr); 245 if (!dmX) dmX = dm; 246 ierr = PetscObjectQuery(link->vec, "section", (PetscObject*) §ion);CHKERRQ(ierr); 247 if (!section) {ierr = DMGetLocalSection(dmX, §ion);CHKERRQ(ierr);} 248 ierr = PetscSectionGetDof(section,cStart,&bs);CHKERRQ(ierr); 249 ierr = PetscSectionGetNumFields(section,&nfields);CHKERRQ(ierr); 250 field = 0; 251 if (link->field >= 0) { 252 field = link->field; 253 nfields = field + 1; 254 } 255 for (i=0; field<(nfields?nfields:1); field++) { 256 PetscInt fbs,j; 257 PetscFV fv = NULL; 258 PetscObject f; 259 PetscClassId fClass; 260 const char *fieldname = NULL; 261 char buf[256]; 262 PetscBool vector; 263 if (nfields) { /* We have user-defined fields/components */ 264 ierr = PetscSectionGetFieldDof(section,cStart,field,&fbs);CHKERRQ(ierr); 265 ierr = PetscSectionGetFieldName(section,field,&fieldname);CHKERRQ(ierr); 266 } else fbs = bs; /* Say we have one field with 'bs' components */ 267 ierr = DMGetField(dmX,field,NULL,&f);CHKERRQ(ierr); 268 ierr = PetscObjectGetClassId(f,&fClass);CHKERRQ(ierr); 269 if (fClass == PETSCFV_CLASSID) { 270 fv = (PetscFV) f; 271 } 272 if (nfields && !fieldname) { 273 ierr = PetscSNPrintf(buf,sizeof(buf),"CellField%D",field);CHKERRQ(ierr); 274 fieldname = buf; 275 } 276 vector = PETSC_FALSE; 277 if (link->ft == PETSC_VTK_CELL_VECTOR_FIELD) { 278 vector = PETSC_TRUE; 279 if (fbs > 3) SETERRQ1(PetscObjectComm((PetscObject)viewer),PETSC_ERR_ARG_SIZ,"Cell vector fields can have at most 3 components, %D given\n", fbs); 280 for (j = 0; j < fbs; j++) { 281 const char *compName = NULL; 282 if (fv) { 283 ierr = PetscFVGetComponentName(fv,j,&compName);CHKERRQ(ierr); 284 if (compName) break; 285 } 286 } 287 if (j < fbs) vector = PETSC_FALSE; 288 } 289 if (vector) { 290 #if defined(PETSC_USE_COMPLEX) 291 ierr = PetscFPrintf(comm,fp," <DataArray type=\"%s\" Name=\"%s%s.Re\" NumberOfComponents=\"3\" format=\"appended\" offset=\"%D\" />\n",precision,vecname,fieldname,boffset);CHKERRQ(ierr); 292 boffset += gpiece[r].ncells*3*sizeof(PetscVTUReal) + sizeof(int); 293 ierr = PetscFPrintf(comm,fp," <DataArray type=\"%s\" Name=\"%s%s.Im\" NumberOfComponents=\"3\" format=\"appended\" offset=\"%D\" />\n",precision,vecname,fieldname,boffset);CHKERRQ(ierr); 294 boffset += gpiece[r].ncells*3*sizeof(PetscVTUReal) + sizeof(int); 295 #else 296 ierr = PetscFPrintf(comm,fp," <DataArray type=\"%s\" Name=\"%s%s\" NumberOfComponents=\"3\" format=\"appended\" offset=\"%D\" />\n",precision,vecname,fieldname,boffset);CHKERRQ(ierr); 297 boffset += gpiece[r].ncells*3*sizeof(PetscVTUReal) + sizeof(int); 298 #endif 299 i+=fbs; 300 } else { 301 for (j=0; j<fbs; j++) { 302 const char *compName = NULL; 303 char finalname[256]; 304 if (fv) { 305 ierr = PetscFVGetComponentName(fv,j,&compName);CHKERRQ(ierr); 306 } 307 if (compName) { 308 ierr = PetscSNPrintf(finalname,255,"%s%s.%s",vecname,fieldname,compName);CHKERRQ(ierr); 309 } 310 else if (fbs > 1) { 311 ierr = PetscSNPrintf(finalname,255,"%s%s.%D",vecname,fieldname,j);CHKERRQ(ierr); 312 } else { 313 ierr = PetscSNPrintf(finalname,255,"%s%s",vecname,fieldname);CHKERRQ(ierr); 314 } 315 #if defined(PETSC_USE_COMPLEX) 316 ierr = PetscFPrintf(comm,fp," <DataArray type=\"%s\" Name=\"%s.Re\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",precision,finalname,boffset);CHKERRQ(ierr); 317 boffset += gpiece[r].ncells*sizeof(PetscVTUReal) + sizeof(int); 318 ierr = PetscFPrintf(comm,fp," <DataArray type=\"%s\" Name=\"%s.Im\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",precision,finalname,boffset);CHKERRQ(ierr); 319 boffset += gpiece[r].ncells*sizeof(PetscVTUReal) + sizeof(int); 320 #else 321 ierr = PetscFPrintf(comm,fp," <DataArray type=\"%s\" Name=\"%s\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",precision,finalname,boffset);CHKERRQ(ierr); 322 boffset += gpiece[r].ncells*sizeof(PetscVTUReal) + sizeof(int); 323 #endif 324 i++; 325 } 326 } 327 } 328 if (i != bs) SETERRQ2(comm,PETSC_ERR_PLIB,"Total number of field components %D != block size %D",i,bs); 329 } 330 ierr = PetscFPrintf(PETSC_COMM_SELF,fp," </CellData>\n");CHKERRQ(ierr); 331 332 /* 333 * Point Data headers 334 */ 335 ierr = PetscFPrintf(PETSC_COMM_SELF,fp," <PointData>\n");CHKERRQ(ierr); 336 for (link=vtk->link; link; link=link->next) { 337 Vec X = (Vec)link->vec; 338 DM dmX; 339 PetscInt bs,nfields,field; 340 const char *vecname = ""; 341 PetscSection section; 342 if ((link->ft != PETSC_VTK_POINT_FIELD) && (link->ft != PETSC_VTK_POINT_VECTOR_FIELD)) continue; 343 if (((PetscObject)X)->name || link != vtk->link) { /* If the object is already named, use it. If it is past the first link, name it to disambiguate. */ 344 ierr = PetscObjectGetName((PetscObject)X,&vecname);CHKERRQ(ierr); 345 } 346 ierr = VecGetDM(X, &dmX);CHKERRQ(ierr); 347 if (!dmX) dmX = dm; 348 ierr = PetscObjectQuery(link->vec, "section", (PetscObject*) §ion);CHKERRQ(ierr); 349 if (!section) {ierr = DMGetLocalSection(dmX, §ion);CHKERRQ(ierr);} 350 ierr = PetscSectionGetDof(section,vStart,&bs);CHKERRQ(ierr); 351 ierr = PetscSectionGetNumFields(section,&nfields);CHKERRQ(ierr); 352 field = 0; 353 if (link->field >= 0) { 354 field = link->field; 355 nfields = field + 1; 356 } 357 for (i=0; field<(nfields?nfields:1); field++) { 358 PetscInt fbs,j; 359 const char *fieldname = NULL; 360 char buf[256]; 361 if (nfields) { /* We have user-defined fields/components */ 362 ierr = PetscSectionGetFieldDof(section,vStart,field,&fbs);CHKERRQ(ierr); 363 ierr = PetscSectionGetFieldName(section,field,&fieldname);CHKERRQ(ierr); 364 } else fbs = bs; /* Say we have one field with 'bs' components */ 365 if (nfields && !fieldname) { 366 ierr = PetscSNPrintf(buf,sizeof(buf),"PointField%D",field);CHKERRQ(ierr); 367 fieldname = buf; 368 } 369 if (link->ft == PETSC_VTK_POINT_VECTOR_FIELD) { 370 if (fbs > 3) SETERRQ1(PetscObjectComm((PetscObject)viewer),PETSC_ERR_ARG_SIZ,"Point vector fields can have at most 3 components, %D given\n", fbs); 371 #if defined(PETSC_USE_COMPLEX) 372 ierr = PetscFPrintf(comm,fp," <DataArray type=\"%s\" Name=\"%s%s.Re\" NumberOfComponents=\"3\" format=\"appended\" offset=\"%D\" />\n",precision,vecname,fieldname,boffset);CHKERRQ(ierr); 373 boffset += gpiece[r].nvertices*3*sizeof(PetscVTUReal) + sizeof(int); 374 ierr = PetscFPrintf(comm,fp," <DataArray type=\"%s\" Name=\"%s%s.Im\" NumberOfComponents=\"3\" format=\"appended\" offset=\"%D\" />\n",precision,vecname,fieldname,boffset);CHKERRQ(ierr); 375 boffset += gpiece[r].nvertices*3*sizeof(PetscVTUReal) + sizeof(int); 376 #else 377 ierr = PetscFPrintf(comm,fp," <DataArray type=\"%s\" Name=\"%s%s\" NumberOfComponents=\"3\" format=\"appended\" offset=\"%D\" />\n",precision,vecname,fieldname,boffset);CHKERRQ(ierr); 378 boffset += gpiece[r].nvertices*3*sizeof(PetscVTUReal) + sizeof(int); 379 #endif 380 } else { 381 for (j=0; j<fbs; j++) { 382 char finalname[256]; 383 if (fbs > 1) { 384 ierr = PetscSNPrintf(finalname,255,"%s%s.%D",vecname,fieldname,j);CHKERRQ(ierr); 385 } else { 386 ierr = PetscSNPrintf(finalname,255,"%s%s",vecname,fieldname);CHKERRQ(ierr); 387 } 388 #if defined(PETSC_USE_COMPLEX) 389 ierr = PetscFPrintf(comm,fp," <DataArray type=\"%s\" Name=\"%s.Re\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",precision,finalname,boffset);CHKERRQ(ierr); 390 boffset += gpiece[r].nvertices*sizeof(PetscVTUReal) + sizeof(int); 391 ierr = PetscFPrintf(comm,fp," <DataArray type=\"%s\" Name=\"%s.Im\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",precision,finalname,boffset);CHKERRQ(ierr); 392 boffset += gpiece[r].nvertices*sizeof(PetscVTUReal) + sizeof(int); 393 #else 394 ierr = PetscFPrintf(comm,fp," <DataArray type=\"%s\" Name=\"%s\" NumberOfComponents=\"1\" format=\"appended\" offset=\"%D\" />\n",precision,finalname,boffset);CHKERRQ(ierr); 395 boffset += gpiece[r].nvertices*sizeof(PetscVTUReal) + sizeof(int); 396 #endif 397 } 398 } 399 } 400 } 401 ierr = PetscFPrintf(PETSC_COMM_SELF,fp," </PointData>\n");CHKERRQ(ierr); 402 ierr = PetscFPrintf(PETSC_COMM_SELF,fp," </Piece>\n");CHKERRQ(ierr); 403 } 404 } 405 406 ierr = PetscFPrintf(comm,fp," </UnstructuredGrid>\n");CHKERRQ(ierr); 407 ierr = PetscFPrintf(comm,fp," <AppendedData encoding=\"raw\">\n");CHKERRQ(ierr); 408 ierr = PetscFPrintf(comm,fp,"_");CHKERRQ(ierr); 409 410 if (!rank) { 411 PetscInt maxsize = 0; 412 for (r=0; r<size; r++) { 413 maxsize = PetscMax(maxsize, (PetscInt) (gpiece[r].nvertices*3*sizeof(PetscVTUReal))); 414 maxsize = PetscMax(maxsize, (PetscInt) (gpiece[r].ncells*3*sizeof(PetscVTUReal))); 415 maxsize = PetscMax(maxsize, (PetscInt) (gpiece[r].nconn*sizeof(PetscVTKInt))); 416 } 417 ierr = PetscMalloc(maxsize,&buffer);CHKERRQ(ierr); 418 } 419 for (r=0; r<size; r++) { 420 if (r == rank) { 421 PetscInt nsend; 422 { /* Position */ 423 const PetscScalar *x; 424 PetscVTUReal *y = NULL; 425 Vec coords; 426 PetscBool copy; 427 428 ierr = DMGetCoordinatesLocal(dm,&coords);CHKERRQ(ierr); 429 ierr = VecGetArrayRead(coords,&x);CHKERRQ(ierr); 430 #if defined(PETSC_USE_COMPLEX) 431 copy = PETSC_TRUE; 432 #else 433 copy = (PetscBool) (dimEmbed != 3 || localized || (sizeof(PetscReal) != sizeof(PetscVTUReal))); 434 #endif 435 if (copy) { 436 ierr = PetscMalloc1(piece.nvertices*3,&y);CHKERRQ(ierr); 437 if (localized) { 438 PetscInt cnt; 439 for (c=cStart,cnt=0; c<cEnd; c++) { 440 PetscInt off, dof; 441 442 ierr = PetscSectionGetDof(coordSection, c, &dof);CHKERRQ(ierr); 443 if (!dof) { 444 PetscInt *closure = NULL; 445 PetscInt closureSize; 446 447 ierr = DMPlexGetTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr); 448 for (v = 0; v < closureSize*2; v += 2) { 449 if ((closure[v] >= vStart) && (closure[v] < vEnd)) { 450 ierr = PetscSectionGetOffset(coordSection, closure[v], &off);CHKERRQ(ierr); 451 if (dimEmbed != 3) { 452 y[cnt*3+0] = (PetscVTUReal) PetscRealPart(x[off+0]); 453 y[cnt*3+1] = (PetscVTUReal) ((dimEmbed > 1) ? PetscRealPart(x[off+1]) : 0.0); 454 y[cnt*3+2] = (PetscVTUReal) 0.0; 455 } else { 456 y[cnt*3+0] = (PetscVTUReal) PetscRealPart(x[off+0]); 457 y[cnt*3+1] = (PetscVTUReal) PetscRealPart(x[off+1]); 458 y[cnt*3+2] = (PetscVTUReal) PetscRealPart(x[off+2]); 459 } 460 cnt++; 461 } 462 } 463 ierr = DMPlexRestoreTransitiveClosure(dm, c, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr); 464 } else { 465 ierr = PetscSectionGetOffset(coordSection, c, &off);CHKERRQ(ierr); 466 if (dimEmbed != 3) { 467 for (i=0; i<dof/dimEmbed; i++) { 468 y[cnt*3+0] = (PetscVTUReal) PetscRealPart(x[off + i*dimEmbed + 0]); 469 y[cnt*3+1] = (PetscVTUReal) ((dimEmbed > 1) ? PetscRealPart(x[off + i*dimEmbed + 1]) : 0.0); 470 y[cnt*3+2] = (PetscVTUReal) 0.0; 471 cnt++; 472 } 473 } else { 474 for (i=0; i<dof; i ++) { 475 y[cnt*3+i] = (PetscVTUReal) PetscRealPart(x[off + i]); 476 } 477 cnt += dof/dimEmbed; 478 } 479 } 480 } 481 if (cnt != piece.nvertices) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Count does not match"); 482 } else { 483 for (i=0; i<piece.nvertices; i++) { 484 y[i*3+0] = (PetscVTUReal) PetscRealPart(x[i*dimEmbed+0]); 485 y[i*3+1] = (PetscVTUReal) ((dimEmbed > 1) ? PetscRealPart(x[i*dimEmbed+1]) : 0.); 486 y[i*3+2] = (PetscVTUReal) ((dimEmbed > 2) ? PetscRealPart(x[i*dimEmbed+2]) : 0.); 487 } 488 } 489 } 490 nsend = piece.nvertices*3; 491 ierr = TransferWrite(viewer,fp,r,0,copy ? (const void *) y : (const void *) x,buffer,nsend,MPIU_VTUREAL,tag);CHKERRQ(ierr); 492 ierr = PetscFree(y);CHKERRQ(ierr); 493 ierr = VecRestoreArrayRead(coords,&x);CHKERRQ(ierr); 494 } 495 { /* Connectivity, offsets, types */ 496 PetscVTKInt *connectivity = NULL, *offsets = NULL; 497 PetscVTKType *types = NULL; 498 ierr = DMPlexGetVTKConnectivity(dm,localized,&piece,&connectivity,&offsets,&types);CHKERRQ(ierr); 499 ierr = TransferWrite(viewer,fp,r,0,connectivity,buffer,piece.nconn,MPI_INT,tag);CHKERRQ(ierr); 500 ierr = TransferWrite(viewer,fp,r,0,offsets,buffer,piece.ncells,MPI_INT,tag);CHKERRQ(ierr); 501 ierr = TransferWrite(viewer,fp,r,0,types,buffer,piece.ncells,MPI_CHAR,tag);CHKERRQ(ierr); 502 ierr = PetscFree3(connectivity,offsets,types);CHKERRQ(ierr); 503 } 504 { /* Owners (cell data) */ 505 PetscVTKInt *owners; 506 ierr = PetscMalloc1(piece.ncells,&owners);CHKERRQ(ierr); 507 for (i=0; i<piece.ncells; i++) owners[i] = rank; 508 ierr = TransferWrite(viewer,fp,r,0,owners,buffer,piece.ncells,MPI_INT,tag);CHKERRQ(ierr); 509 ierr = PetscFree(owners);CHKERRQ(ierr); 510 } 511 /* Cell data */ 512 for (link=vtk->link; link; link=link->next) { 513 Vec X = (Vec)link->vec; 514 DM dmX; 515 const PetscScalar *x; 516 PetscVTUReal *y; 517 PetscInt bs, nfields, field; 518 PetscSection section = NULL; 519 520 if ((link->ft != PETSC_VTK_CELL_FIELD) && (link->ft != PETSC_VTK_CELL_VECTOR_FIELD)) continue; 521 ierr = VecGetDM(X, &dmX);CHKERRQ(ierr); 522 if (!dmX) dmX = dm; 523 ierr = PetscObjectQuery(link->vec, "section", (PetscObject*) §ion);CHKERRQ(ierr); 524 if (!section) {ierr = DMGetLocalSection(dmX, §ion);CHKERRQ(ierr);} 525 ierr = PetscSectionGetDof(section,cStart,&bs);CHKERRQ(ierr); 526 ierr = PetscSectionGetNumFields(section,&nfields);CHKERRQ(ierr); 527 field = 0; 528 if (link->field >= 0) { 529 field = link->field; 530 nfields = field + 1; 531 } 532 ierr = VecGetArrayRead(X,&x);CHKERRQ(ierr); 533 ierr = PetscMalloc1(piece.ncells*3,&y);CHKERRQ(ierr); 534 for (i=0; field<(nfields?nfields:1); field++) { 535 PetscInt fbs,j; 536 PetscFV fv = NULL; 537 PetscObject f; 538 PetscClassId fClass; 539 PetscBool vector; 540 if (nfields) { /* We have user-defined fields/components */ 541 ierr = PetscSectionGetFieldDof(section,cStart,field,&fbs);CHKERRQ(ierr); 542 } else fbs = bs; /* Say we have one field with 'bs' components */ 543 ierr = DMGetField(dmX,field,NULL,&f);CHKERRQ(ierr); 544 ierr = PetscObjectGetClassId(f,&fClass);CHKERRQ(ierr); 545 if (fClass == PETSCFV_CLASSID) { 546 fv = (PetscFV) f; 547 } 548 vector = PETSC_FALSE; 549 if (link->ft == PETSC_VTK_CELL_VECTOR_FIELD) { 550 vector = PETSC_TRUE; 551 for (j = 0; j < fbs; j++) { 552 const char *compName = NULL; 553 if (fv) { 554 ierr = PetscFVGetComponentName(fv,j,&compName);CHKERRQ(ierr); 555 if (compName) break; 556 } 557 } 558 if (j < fbs) vector = PETSC_FALSE; 559 } 560 if (vector) { 561 PetscInt cnt, l; 562 for (l = 0; l < loops_per_scalar; l++) { 563 for (c=cStart,cnt=0; c<cEnd; c++) { 564 const PetscScalar *xpoint; 565 PetscInt off, j; 566 567 if (hasLabel) { /* Ignore some cells */ 568 PetscInt value; 569 ierr = DMGetLabelValue(dmX, "vtk", c, &value);CHKERRQ(ierr); 570 if (value != 1) continue; 571 } 572 if (nfields) { 573 ierr = PetscSectionGetFieldOffset(section, c, field, &off);CHKERRQ(ierr); 574 } else { 575 ierr = PetscSectionGetOffset(section, c, &off);CHKERRQ(ierr); 576 } 577 xpoint = &x[off]; 578 for (j = 0; j < fbs; j++) { 579 y[cnt++] = (PetscVTUReal) (l ? PetscImaginaryPart(xpoint[j]) : PetscRealPart(xpoint[j])); 580 } 581 for (; j < 3; j++) y[cnt++] = 0.; 582 } 583 if (cnt != piece.ncells*3) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Count does not match"); 584 ierr = TransferWrite(viewer,fp,r,0,y,buffer,piece.ncells*3,MPIU_VTUREAL,tag);CHKERRQ(ierr); 585 } 586 } else { 587 for (i=0; i<fbs; i++) { 588 PetscInt cnt, l; 589 for (l = 0; l < loops_per_scalar; l++) { 590 for (c=cStart,cnt=0; c<cEnd; c++) { 591 const PetscScalar *xpoint; 592 PetscInt off; 593 594 if (hasLabel) { /* Ignore some cells */ 595 PetscInt value; 596 ierr = DMGetLabelValue(dmX, "vtk", c, &value);CHKERRQ(ierr); 597 if (value != 1) continue; 598 } 599 if (nfields) { 600 ierr = PetscSectionGetFieldOffset(section, c, field, &off);CHKERRQ(ierr); 601 } else { 602 ierr = PetscSectionGetOffset(section, c, &off);CHKERRQ(ierr); 603 } 604 xpoint = &x[off]; 605 y[cnt++] = (PetscVTUReal) (l ? PetscImaginaryPart(xpoint[i]) : PetscRealPart(xpoint[i])); 606 } 607 if (cnt != piece.ncells) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Count does not match"); 608 ierr = TransferWrite(viewer,fp,r,0,y,buffer,piece.ncells,MPIU_VTUREAL,tag);CHKERRQ(ierr); 609 } 610 } 611 } 612 } 613 ierr = PetscFree(y);CHKERRQ(ierr); 614 ierr = VecRestoreArrayRead(X,&x);CHKERRQ(ierr); 615 } 616 /* point data */ 617 for (link=vtk->link; link; link=link->next) { 618 Vec X = (Vec)link->vec; 619 DM dmX; 620 const PetscScalar *x; 621 PetscVTUReal *y; 622 PetscInt bs, nfields, field; 623 PetscSection section = NULL; 624 625 if ((link->ft != PETSC_VTK_POINT_FIELD) && (link->ft != PETSC_VTK_POINT_VECTOR_FIELD)) continue; 626 ierr = VecGetDM(X, &dmX);CHKERRQ(ierr); 627 if (!dmX) dmX = dm; 628 ierr = PetscObjectQuery(link->vec, "section", (PetscObject*) §ion);CHKERRQ(ierr); 629 if (!section) {ierr = DMGetLocalSection(dmX, §ion);CHKERRQ(ierr);} 630 ierr = PetscSectionGetDof(section,vStart,&bs);CHKERRQ(ierr); 631 ierr = PetscSectionGetNumFields(section,&nfields);CHKERRQ(ierr); 632 field = 0; 633 if (link->field >= 0) { 634 field = link->field; 635 nfields = field + 1; 636 } 637 ierr = VecGetArrayRead(X,&x);CHKERRQ(ierr); 638 ierr = PetscMalloc1(piece.nvertices*3,&y);CHKERRQ(ierr); 639 for (i=0; field<(nfields?nfields:1); field++) { 640 PetscInt fbs,j; 641 if (nfields) { /* We have user-defined fields/components */ 642 ierr = PetscSectionGetFieldDof(section,vStart,field,&fbs);CHKERRQ(ierr); 643 } else fbs = bs; /* Say we have one field with 'bs' components */ 644 if (link->ft == PETSC_VTK_POINT_VECTOR_FIELD) { 645 PetscInt cnt, l; 646 for (l = 0; l < loops_per_scalar; l++) { 647 if (!localized) { 648 for (v=vStart,cnt=0; v<vEnd; v++) { 649 PetscInt off; 650 const PetscScalar *xpoint; 651 652 if (nfields) { 653 ierr = PetscSectionGetFieldOffset(section,v,field,&off);CHKERRQ(ierr); 654 } else { 655 ierr = PetscSectionGetOffset(section,v,&off);CHKERRQ(ierr); 656 } 657 xpoint = &x[off]; 658 for (j = 0; j < fbs; j++) { 659 y[cnt++] = (PetscVTUReal) (l ? PetscImaginaryPart(xpoint[j]) : PetscRealPart(xpoint[j])); 660 } 661 for (; j < 3; j++) y[cnt++] = 0.; 662 } 663 } else { 664 for (c=cStart,cnt=0; c<cEnd; c++) { 665 PetscInt *closure = NULL; 666 PetscInt closureSize, off; 667 668 ierr = DMPlexGetTransitiveClosure(dmX, c, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr); 669 for (v = 0, off = 0; v < closureSize*2; v += 2) { 670 if ((closure[v] >= vStart) && (closure[v] < vEnd)) { 671 PetscInt voff; 672 const PetscScalar *xpoint; 673 674 if (nfields) { 675 ierr = PetscSectionGetFieldOffset(section,closure[v],field,&voff);CHKERRQ(ierr); 676 } else { 677 ierr = PetscSectionGetOffset(section,closure[v],&voff);CHKERRQ(ierr); 678 } 679 xpoint = &x[voff]; 680 for (j = 0; j < fbs; j++) { 681 y[cnt + off++] = (PetscVTUReal) (l ? PetscImaginaryPart(xpoint[j]) : PetscRealPart(xpoint[j])); 682 } 683 for (; j < 3; j++) y[cnt + off++] = 0.; 684 } 685 } 686 cnt += off; 687 ierr = DMPlexRestoreTransitiveClosure(dmX, c, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr); 688 } 689 } 690 if (cnt != piece.nvertices*3) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Count does not match"); 691 ierr = TransferWrite(viewer,fp,r,0,y,buffer,piece.nvertices*3,MPIU_VTUREAL,tag);CHKERRQ(ierr); 692 } 693 } else { 694 for (i=0; i<fbs; i++) { 695 PetscInt cnt, l; 696 for (l = 0; l < loops_per_scalar; l++) { 697 if (!localized) { 698 for (v=vStart,cnt=0; v<vEnd; v++) { 699 PetscInt off; 700 const PetscScalar *xpoint; 701 702 if (nfields) { 703 ierr = PetscSectionGetFieldOffset(section,v,field,&off);CHKERRQ(ierr); 704 } else { 705 ierr = PetscSectionGetOffset(section,v,&off);CHKERRQ(ierr); 706 } 707 xpoint = &x[off]; 708 y[cnt++] = (PetscVTUReal) (l ? PetscImaginaryPart(xpoint[i]) : PetscRealPart(xpoint[i])); 709 } 710 } else { 711 for (c=cStart,cnt=0; c<cEnd; c++) { 712 PetscInt *closure = NULL; 713 PetscInt closureSize, off; 714 715 ierr = DMPlexGetTransitiveClosure(dmX, c, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr); 716 for (v = 0, off = 0; v < closureSize*2; v += 2) { 717 if ((closure[v] >= vStart) && (closure[v] < vEnd)) { 718 PetscInt voff; 719 const PetscScalar *xpoint; 720 721 if (nfields) { 722 ierr = PetscSectionGetFieldOffset(section,closure[v],field,&voff);CHKERRQ(ierr); 723 } else { 724 ierr = PetscSectionGetOffset(section,closure[v],&voff);CHKERRQ(ierr); 725 } 726 xpoint = &x[voff]; 727 y[cnt + off++] = (l ? PetscImaginaryPart(xpoint[i]) : PetscRealPart(xpoint[i])); 728 } 729 } 730 cnt += off; 731 ierr = DMPlexRestoreTransitiveClosure(dmX, c, PETSC_TRUE, &closureSize, &closure);CHKERRQ(ierr); 732 } 733 } 734 if (cnt != piece.nvertices) SETERRQ(PETSC_COMM_SELF,PETSC_ERR_PLIB,"Count does not match"); 735 ierr = TransferWrite(viewer,fp,r,0,y,buffer,piece.nvertices,MPIU_VTUREAL,tag);CHKERRQ(ierr); 736 } 737 } 738 } 739 } 740 ierr = PetscFree(y);CHKERRQ(ierr); 741 ierr = VecRestoreArrayRead(X,&x);CHKERRQ(ierr); 742 } 743 } else if (!rank) { 744 PetscInt l; 745 746 ierr = TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].nvertices*3,MPIU_VTUREAL,tag);CHKERRQ(ierr); /* positions */ 747 ierr = TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].nconn,MPI_INT,tag);CHKERRQ(ierr); /* connectivity */ 748 ierr = TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].ncells,MPI_INT,tag);CHKERRQ(ierr); /* offsets */ 749 ierr = TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].ncells,MPI_CHAR,tag);CHKERRQ(ierr); /* types */ 750 ierr = TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].ncells,MPI_INT,tag);CHKERRQ(ierr); /* owner rank (cells) */ 751 /* all cell data */ 752 for (link=vtk->link; link; link=link->next) { 753 Vec X = (Vec)link->vec; 754 PetscInt bs, nfields, field; 755 DM dmX; 756 PetscSection section = NULL; 757 758 if ((link->ft != PETSC_VTK_CELL_FIELD) && (link->ft != PETSC_VTK_CELL_VECTOR_FIELD)) continue; 759 ierr = VecGetDM(X, &dmX);CHKERRQ(ierr); 760 if (!dmX) dmX = dm; 761 ierr = PetscObjectQuery(link->vec, "section", (PetscObject*) §ion);CHKERRQ(ierr); 762 if (!section) {ierr = DMGetLocalSection(dmX, §ion);CHKERRQ(ierr);} 763 ierr = PetscSectionGetDof(section,cStart,&bs);CHKERRQ(ierr); 764 ierr = PetscSectionGetNumFields(section,&nfields);CHKERRQ(ierr); 765 field = 0; 766 if (link->field >= 0) { 767 field = link->field; 768 nfields = field + 1; 769 } 770 for (i=0; field<(nfields?nfields:1); field++) { 771 PetscInt fbs,j; 772 PetscFV fv = NULL; 773 PetscObject f; 774 PetscClassId fClass; 775 PetscBool vector; 776 if (nfields) { /* We have user-defined fields/components */ 777 ierr = PetscSectionGetFieldDof(section,cStart,field,&fbs);CHKERRQ(ierr); 778 } else fbs = bs; /* Say we have one field with 'bs' components */ 779 ierr = DMGetField(dmX,field,NULL,&f);CHKERRQ(ierr); 780 ierr = PetscObjectGetClassId(f,&fClass);CHKERRQ(ierr); 781 if (fClass == PETSCFV_CLASSID) { 782 fv = (PetscFV) f; 783 } 784 vector = PETSC_FALSE; 785 if (link->ft == PETSC_VTK_CELL_VECTOR_FIELD) { 786 vector = PETSC_TRUE; 787 for (j = 0; j < fbs; j++) { 788 const char *compName = NULL; 789 if (fv) { 790 ierr = PetscFVGetComponentName(fv,j,&compName);CHKERRQ(ierr); 791 if (compName) break; 792 } 793 } 794 if (j < fbs) vector = PETSC_FALSE; 795 } 796 if (vector) { 797 for (l = 0; l < loops_per_scalar; l++) { 798 ierr = TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].ncells*3,MPIU_VTUREAL,tag);CHKERRQ(ierr); 799 } 800 } else { 801 for (i=0; i<fbs; i++) { 802 for (l = 0; l < loops_per_scalar; l++) { 803 ierr = TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].ncells,MPIU_VTUREAL,tag);CHKERRQ(ierr); 804 } 805 } 806 } 807 } 808 } 809 /* all point data */ 810 for (link=vtk->link; link; link=link->next) { 811 Vec X = (Vec)link->vec; 812 DM dmX; 813 PetscInt bs, nfields, field; 814 PetscSection section = NULL; 815 816 if ((link->ft != PETSC_VTK_POINT_FIELD) && (link->ft != PETSC_VTK_POINT_VECTOR_FIELD)) continue; 817 ierr = VecGetDM(X, &dmX);CHKERRQ(ierr); 818 if (!dmX) dmX = dm; 819 ierr = PetscObjectQuery(link->vec, "section", (PetscObject*) §ion);CHKERRQ(ierr); 820 if (!section) {ierr = DMGetLocalSection(dmX, §ion);CHKERRQ(ierr);} 821 ierr = PetscSectionGetDof(section,vStart,&bs);CHKERRQ(ierr); 822 ierr = PetscSectionGetNumFields(section,&nfields);CHKERRQ(ierr); 823 field = 0; 824 if (link->field >= 0) { 825 field = link->field; 826 nfields = field + 1; 827 } 828 for (i=0; field<(nfields?nfields:1); field++) { 829 PetscInt fbs; 830 if (nfields) { /* We have user-defined fields/components */ 831 ierr = PetscSectionGetFieldDof(section,vStart,field,&fbs);CHKERRQ(ierr); 832 } else fbs = bs; /* Say we have one field with 'bs' components */ 833 if (link->ft == PETSC_VTK_POINT_VECTOR_FIELD) { 834 for (l = 0; l < loops_per_scalar; l++) { 835 ierr = TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].nvertices*3,MPIU_VTUREAL,tag);CHKERRQ(ierr); 836 } 837 } else { 838 for (i=0; i<fbs; i++) { 839 for (l = 0; l < loops_per_scalar; l++) { 840 ierr = TransferWrite(viewer,fp,r,0,NULL,buffer,gpiece[r].nvertices,MPIU_VTUREAL,tag);CHKERRQ(ierr); 841 } 842 } 843 } 844 } 845 } 846 } 847 } 848 ierr = PetscFree(gpiece);CHKERRQ(ierr); 849 ierr = PetscFree(buffer);CHKERRQ(ierr); 850 ierr = PetscFPrintf(comm,fp,"\n </AppendedData>\n");CHKERRQ(ierr); 851 ierr = PetscFPrintf(comm,fp,"</VTKFile>\n");CHKERRQ(ierr); 852 ierr = PetscFClose(comm,fp);CHKERRQ(ierr); 853 PetscFunctionReturn(0); 854 } 855