[cig-commits] r13124 - in long/3D/Gale/trunk: . src/Underworld/plugins/Output/VTKOutput
walter at geodynamics.org
walter at geodynamics.org
Thu Oct 23 17:36:30 PDT 2008
Author: walter
Date: 2008-10-23 17:36:30 -0700 (Thu, 23 Oct 2008)
New Revision: 13124
Modified:
long/3D/Gale/trunk/
long/3D/Gale/trunk/src/Underworld/plugins/Output/VTKOutput/VTKOutput.c
Log:
r2358 at earth: boo | 2008-10-23 16:37:30 -0700
Make VTKOutput no longer output pressure separately, since it was meaningless unless averaged over nodes.
Property changes on: long/3D/Gale/trunk
___________________________________________________________________
Name: svk:merge
- 3a629746-de10-0410-b17b-fd6ecaaa963e:/cig:2353
+ 3a629746-de10-0410-b17b-fd6ecaaa963e:/cig:2358
Modified: long/3D/Gale/trunk/src/Underworld/plugins/Output/VTKOutput/VTKOutput.c
===================================================================
--- long/3D/Gale/trunk/src/Underworld/plugins/Output/VTKOutput/VTKOutput.c 2008-10-23 23:35:27 UTC (rev 13123)
+++ long/3D/Gale/trunk/src/Underworld/plugins/Output/VTKOutput/VTKOutput.c 2008-10-24 00:36:30 UTC (rev 13124)
@@ -402,9 +402,7 @@
}
-/* Pressure is stored on cell centers, while everything else is stored
- on cell vertices. So we have to make two different files, one for
- pressure, and one for everything else. */
+/* Everything is stored on cell vertices. */
void VTKOutput_fields(void *context, int myRank, int nprocs) {
FiniteElementContext* self = (FiniteElementContext*) context;
@@ -413,39 +411,28 @@
int nDims, pn, i;
int lower[3], upper[3], p_lower[3], p_upper[3];
- Name field_filename, pressure_filename;
- FILE *fp, *field_fp, *pressure_fp, *pfp, *pfield_fp, *ppressure_fp;
+ Name field_filename;
+ FILE *field_fp, *pfield_fp;
- /* We need to save the grids to be used later when printing out the
- pressure coordinates and to map between 1D and 3D indices for the
- values. */
+ /* We need to save the grids to map between 1D and 3D indices for
+ the values. */
Grid *elGrid, *vertGrid;
/* Open the file */
- Stg_asprintf( &pressure_filename, "%s/pressure.%d.%05d.vts",
- self->outputPath,myRank,
- self->timeStep);
Stg_asprintf( &field_filename, "%s/fields.%d.%05d.vts", self->outputPath,
myRank, self->timeStep);
field_fp=fopen(field_filename,"w");
- pressure_fp=fopen(pressure_filename,"w");
- Memory_Free( pressure_filename );
Memory_Free( field_filename );
/* Print out the parallel control files if rank==0 */
if(myRank==0)
{
- Stg_asprintf( &pressure_filename, "%s/pressure.%05d.pvts",
- self->outputPath,
- self->timeStep);
Stg_asprintf( &field_filename, "%s/fields.%05d.pvts", self->outputPath,
self->timeStep);
pfield_fp=fopen(field_filename,"w");
- ppressure_fp=fopen(pressure_filename,"w");
- Memory_Free( pressure_filename );
Memory_Free( field_filename );
}
@@ -457,7 +444,8 @@
FieldVariable* fieldVar;
fieldVar = FieldVariable_Register_GetByIndex( self->fieldVariable_Register,
var_I );
- if (Stg_Class_IsInstance( fieldVar, FeVariable_Type )) {
+ if (Stg_Class_IsInstance( fieldVar, FeVariable_Type )
+ && strcmp(((FeVariable*)fieldVar)->name,"PressureField")) {
FeVariable* feVar;
Node_LocalIndex lNode_I;
Dof_Index dofAtEachNodeCount;
@@ -492,7 +480,7 @@
p_lower[i]=gen->origin[i];
p_upper[i]=p_lower[i]+gen->range[i];
- /* The grid is split up by elements, so for the pressure
+ /* The grid is split up by elements, so for the element
grid, we simply add the ghost zones. Ghost zones are
only added at the top, not the bottom. */
if(p_upper[i]!=gen->elGrid->sizes[i])
@@ -515,17 +503,8 @@
<CellData></CellData>\n",
lower[0],upper[0]-1,lower[1],upper[1]-1,lower[2],upper[2]-1,
lower[0],upper[0]-1,lower[1],upper[1]-1,lower[2],upper[2]-1);
- fprintf(pressure_fp,"<?xml version=\"1.0\"?>\n\
-<VTKFile type=\"StructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">\n\
- <StructuredGrid WholeExtent=\"%d %d %d %d %d %d\">\n\
- <Piece Extent=\"%d %d %d %d %d %d\">\n\
- <CellData></CellData>\n",
- p_lower[0],p_upper[0]-1,p_lower[1],p_upper[1]-1,
- p_lower[2],p_upper[2]-1,
- p_lower[0],p_upper[0]-1,p_lower[1],p_upper[1]-1,
- p_lower[2],p_upper[2]-1);
- /* Write the coordinates for the fields, but not the pressure */
+ /* Write the coordinates for the fields */
VTKOutput_print_coords(field_fp, feVar->feMesh, gen->vertGrid, nDims,
lower, upper);
fprintf(field_fp," <PointData Scalars=\"StrainRateInvariantField\" Vectors=\"VelocityField\" Tensors=\"Stress\">\n");
@@ -553,37 +532,9 @@
<PDataArray type=\"Float64\" NumberOfComponents=\"3\"/>\n\
</PPoints>\n\
<PPointData Scalars=\"StrainRateInvariantField\" Vectors=\"VelocityField\" Tensors=\"Stress\">\n",global[0]-1,global[1]-1,global[2]-1);
- fprintf(ppressure_fp,"<?xml version=\"1.0\"?>\n\
-<VTKFile type=\"PStructuredGrid\" version=\"0.1\" byte_order=\"LittleEndian\">\n \
- <PStructuredGrid GhostLevel=\"0\" WholeExtent=\"0 %d 0 %d 0 %d\">\n\
- <PCellData></PCellData>\n\
- <PPoints>\n\
- <PDataArray type=\"Float64\" NumberOfComponents=\"3\"/>\n\
- </PPoints>\n\
- <PPointData Scalars=\"PressureField\">\n",
- global[0]-2,global[1]-2,p_global[2]-1);
}
header_printed=1;
}
-
- /* Write the coordinates for the pressure, and set the
- appropriate file pointer to output for this variable. */
- if(!strcmp(feVar->name,"PressureField")){
- VTKOutput_print_coords(pressure_fp, feVar->feMesh, elGrid, nDims,
- p_lower, p_upper);
- fprintf(pressure_fp," <PointData Scalars=\"PressureField\">\n");
- fp=pressure_fp;
- pfp=ppressure_fp;
- low=p_lower;
- up=p_upper;
- grid=elGrid;
- } else {
- fp=field_fp;
- pfp=pfield_fp;
- low=lower;
- up=upper;
- grid=vertGrid;
- }
/* Finally, output the fields. For now, just output every field */
@@ -594,25 +545,25 @@
/* Scalars */
case 2:
case 3:
- fprintf(fp," <DataArray type=\"Float64\" Name=\"%s\" format=\"ascii\">\n",feVar->name);
+ fprintf(field_fp," <DataArray type=\"Float64\" Name=\"%s\" format=\"ascii\">\n",feVar->name);
if(myRank==0)
- fprintf(pfp," <PDataArray type=\"Float64\" Name=\"%s\" format=\"ascii\"/>\n",feVar->name);
+ fprintf(pfield_fp," <PDataArray type=\"Float64\" Name=\"%s\" format=\"ascii\"/>\n",feVar->name);
break;
/* Vectors */
case 4:
case 9:
- fprintf(fp," <DataArray type=\"Float64\" Name=\"%s\" format=\"ascii\" NumberOfComponents=\"3\">\n",feVar->name);
+ fprintf(field_fp," <DataArray type=\"Float64\" Name=\"%s\" format=\"ascii\" NumberOfComponents=\"3\">\n",feVar->name);
if(myRank==0)
- fprintf(pfp," <PDataArray type=\"Float64\" Name=\"%s\" format=\"ascii\" NumberOfComponents=\"3\"/>\n",feVar->name);
+ fprintf(pfield_fp," <PDataArray type=\"Float64\" Name=\"%s\" format=\"ascii\" NumberOfComponents=\"3\"/>\n",feVar->name);
break;
/* Rank 2 Tensors */
case 6:
case 8:
case 18:
case 27:
- fprintf(fp," <DataArray type=\"Float64\" Name=\"%s\" format=\"ascii\" NumberOfComponents=\"9\">\n",feVar->name);
+ fprintf(field_fp," <DataArray type=\"Float64\" Name=\"%s\" format=\"ascii\" NumberOfComponents=\"9\">\n",feVar->name);
if(myRank==0)
- fprintf(pfp," <PDataArray type=\"Float64\" Name=\"%s\" format=\"ascii\" NumberOfComponents=\"9\"/>\n",feVar->name);
+ fprintf(pfield_fp," <PDataArray type=\"Float64\" Name=\"%s\" format=\"ascii\" NumberOfComponents=\"9\"/>\n",feVar->name);
break;
/* Unknown */
default:
@@ -622,15 +573,15 @@
abort();
}
- for(ijk[2]=low[2];ijk[2]<up[2];++ijk[2])
- for(ijk[1]=low[1];ijk[1]<up[1];++ijk[1])
- for(ijk[0]=low[0];ijk[0]<up[0];++ijk[0])
+ for(ijk[2]=lower[2];ijk[2]<upper[2];++ijk[2])
+ for(ijk[1]=lower[1];ijk[1]<upper[1];++ijk[1])
+ for(ijk[0]=lower[0];ijk[0]<upper[0];++ijk[0])
{
double variableValues[MAX_FIELD_COMPONENTS];
Dof_Index dof_I;
unsigned local;
Mesh_GlobalToDomain(feVar->feMesh,MT_VERTEX,
- Grid_Project(grid,ijk),&local);
+ Grid_Project(vertGrid,ijk),&local);
FeVariable_GetValueAtNode(feVar,local,variableValues);
switch(dofAtEachNodeCount*nDims)
@@ -643,27 +594,27 @@
case 9:
case 27:
for ( dof_I = 0; dof_I < dofAtEachNodeCount; dof_I++ ) {
- fprintf(fp, "%.15g ", variableValues[dof_I] );
+ fprintf(field_fp, "%.15g ", variableValues[dof_I] );
}
break;
case 4:
- fprintf(fp, "%.15g %.15g 0", variableValues[0],
+ fprintf(field_fp, "%.15g %.15g 0", variableValues[0],
variableValues[1] );
break;
case 6:
/* Ordering is xx, yy, xy */
- fprintf(fp, "%.15g %.15g 0 %.15g %.15g 0 0 0 0 ",
+ fprintf(field_fp, "%.15g %.15g 0 %.15g %.15g 0 0 0 0 ",
variableValues[0], variableValues[2],
variableValues[2], variableValues[1]);
break;
case 8:
- fprintf(fp, "%.15g %.15g 0 %.15g %.15g 0 0 0 0 ",
+ fprintf(field_fp, "%.15g %.15g 0 %.15g %.15g 0 0 0 0 ",
variableValues[0], variableValues[1],
variableValues[2], variableValues[3]);
break;
case 18:
/* Ordering is xx, yy, zz, xy, xz, yz */
- fprintf(fp, "%.15g %.15g %.15g %.15g %.15g %.15g %.15g %.15g %.15g ",
+ fprintf(field_fp, "%.15g %.15g %.15g %.15g %.15g %.15g %.15g %.15g %.15g ",
variableValues[0], variableValues[3],
variableValues[4], variableValues[3],
variableValues[1], variableValues[5],
@@ -671,23 +622,18 @@
variableValues[2]);
break;
}
- fprintf(fp, "\n" );
+ fprintf(field_fp, "\n" );
}
- fprintf(fp," </DataArray>\n");
+ fprintf(field_fp," </DataArray>\n");
}
}
- fprintf(pressure_fp," </PointData>\n </Piece>\n\
- </StructuredGrid>\n\
-</VTKFile>\n");
fprintf(field_fp," </PointData>\n </Piece>\n\
</StructuredGrid>\n\
</VTKFile>\n");
- fclose(pressure_fp);
fclose(field_fp);
if(myRank==0)
{
- fprintf(ppressure_fp," </PPointData>\n");
fprintf(pfield_fp," </PPointData>\n");
for(i=0;i<nprocs;++i)
{
@@ -695,16 +641,9 @@
Source=\"fields.%d.%05d.vts\"/>\n",lower[0],upper[0]-1,
lower[1],upper[1]-1,lower[2],upper[2]-1,
i,self->timeStep);
- fprintf(ppressure_fp," <Piece Extent=\"%d %d %d %d %d %d\"\n\
- Source=\"pressure.%d.%05d.vts\"/>\n",p_lower[0],p_upper[0]-1,
- p_lower[1],p_upper[1]-1,p_lower[2],p_upper[2]-1,
- i,self->timeStep);
}
- fprintf(ppressure_fp," </PStructuredGrid>\n\
-</VTKFile>\n");
fprintf(pfield_fp," </PStructuredGrid>\n\
</VTKFile>\n");
- fclose(ppressure_fp);
fclose(pfield_fp);
}
More information about the CIG-COMMITS
mailing list