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