Write data to disk using hdf5 format in single or double precision.
56 hid_t dataspace, memspace, dataset;
57 hid_t strspace, stratt, string_type;
59 hid_t file_identifier, group, timestep;
60 hid_t file_access = 0;
62 hid_t plist_id_mpiio = 0;
73 static float ****node_coords, ****cell_coords;
74 char filename[512], filenamexmf[512], tstepname[32];
75 char *coords =
"/cell_coords/X /cell_coords/Y /cell_coords/Z ";
76 char *cname[] = {
"X",
"Y",
"Z"};
78 int rank, nd, nv, ns, nc, ngh, ii, jj, kk;
79 int n1p, n2p, n3p, nprec;
87 if (node_coords == NULL) {
94 node_coords =
ARRAY_4D(3, n3p, n2p, n1p,
float);
97 for (kk = 0; kk < n3p; kk++) { x3 = grid[
KDIR].
xl[
KBEG+kk];
98 for (jj = 0; jj < n2p; jj++) { x2 = grid[
JDIR].
xl[
JBEG+jj];
99 for (ii = 0; ii < n1p; ii++) { x1 = grid[
IDIR].
xl[
IBEG+ii];
101 node_coords[
JDIR][kk][jj][ii] = 0.0;
102 node_coords[
KDIR][kk][jj][ii] = 0.0;
104 #if GEOMETRY == CARTESIAN || GEOMETRY == CYLINDRICAL
106 node_coords[
JDIR][kk][jj][ii] = (float)x2; ,
107 node_coords[
KDIR][kk][jj][ii] = (float)x3;)
109 D_EXPAND(node_coords[
IDIR][kk][jj][ii] = (
float)(x1*cos(x2)); ,
110 node_coords[
JDIR][kk][jj][ii] = (float)(x1*sin(x2)); ,
111 node_coords[
KDIR][kk][jj][ii] = (float)(x3);)
114 D_EXPAND(node_coords[
IDIR][kk][jj][ii] = (
float)(x1*sin(x2)); ,
115 node_coords[
JDIR][kk][jj][ii] = (float)(x1*cos(x2)); ,
116 node_coords[
KDIR][kk][jj][ii] = 0.0;)
118 D_EXPAND(node_coords[
IDIR][kk][jj][ii] = (
float)(x1*sin(x2)*cos(x3)); ,
119 node_coords[
JDIR][kk][jj][ii] = (float)(x1*sin(x2)*sin(x3)); ,
120 node_coords[
KDIR][kk][jj][ii] = (float)(x1*cos(x2));)
123 print1 (
"! HDF5_IO: Unknown geometry\n");
128 for (kk = 0; kk <
NX3; kk++) { x3 = grid[
KDIR].
x[
KBEG+kk];
129 for (jj = 0; jj <
NX2; jj++) { x2 = grid[
JDIR].
x[
JBEG+jj];
130 for (ii = 0; ii <
NX1; ii++) { x1 = grid[
IDIR].
x[
IBEG+ii];
132 cell_coords[
JDIR][kk][jj][ii] = 0.0;
133 cell_coords[
KDIR][kk][jj][ii] = 0.0;
135 #if GEOMETRY == CARTESIAN || GEOMETRY == CYLINDRICAL
137 cell_coords[
JDIR][kk][jj][ii] = (float)x2; ,
138 cell_coords[
KDIR][kk][jj][ii] = (float)x3;)
140 D_EXPAND(cell_coords[
IDIR][kk][jj][ii] = (
float)(x1*cos(x2)); ,
141 cell_coords[
JDIR][kk][jj][ii] = (float)(x1*sin(x2)); ,
142 cell_coords[
KDIR][kk][jj][ii] = (float)(x3);)
145 D_EXPAND(cell_coords[
IDIR][kk][jj][ii] = (
float)(x1*sin(x2)); ,
146 cell_coords[
JDIR][kk][jj][ii] = (float)(x1*cos(x2)); ,
147 cell_coords[
KDIR][kk][jj][ii] = 0.0;)
149 D_EXPAND(cell_coords[
IDIR][kk][jj][ii] = (
float)(x1*sin(x2)*cos(x3)); ,
150 cell_coords[
JDIR][kk][jj][ii] = (float)(x1*sin(x2)*sin(x3)); ,
151 cell_coords[
KDIR][kk][jj][ii] = (float)(x1*cos(x2));)
154 print1 (
"! HDF5_IO: Unknown geometry\n");
166 for (nd = 0; nd <
DIMENSIONS; nd++) wgrid[nd] = grid + DIMENSIONS - nd - 1;
169 MPI_Barrier (MPI_COMM_WORLD);
170 if (
prank == 0)time(&tbeg);
173 sprintf (filename,
"%s/data.%04d.%s", output->
dir,output->
nfile, output->
ext);
179 file_access = H5Pcreate(H5P_FILE_ACCESS);
181 H5Pset_fapl_mpiposix(file_access, MPI_COMM_WORLD, 1);
183 H5Pset_fapl_mpio(file_access, MPI_COMM_WORLD, MPI_INFO_NULL);
185 file_identifier = H5Fcreate(filename, H5F_ACC_TRUNC, H5P_DEFAULT, file_access);
186 H5Pclose(file_access);
188 file_identifier = H5Fcreate(filename, H5F_ACC_TRUNC, H5P_DEFAULT, H5P_DEFAULT);
191 sprintf (tstepname,
"Timestep_%d", output->
nfile);
192 timestep = H5Gcreate(file_identifier, tstepname, 0);
194 tspace = H5Screate(H5S_SCALAR);
195 tattr = H5Acreate(timestep,
"Time", H5T_NATIVE_DOUBLE, tspace, H5P_DEFAULT);
196 err = H5Awrite(tattr, H5T_NATIVE_DOUBLE, &
g_time);
200 group = H5Gcreate(timestep,
"vars", 0);
204 strspace = H5Screate_simple(1, &dimstr, NULL);
205 string_type = H5Tcopy(H5T_C_S1);
206 H5Tset_size(string_type, strlen(
"/cell_coords/X "));
207 H5Tset_strpad(string_type,H5T_STR_SPACEPAD);
208 stratt = H5Acreate(group,
"coords",string_type,strspace, H5P_DEFAULT);
209 err = H5Awrite(stratt, string_type, coords);
213 for (nd = 0; nd <
DIMENSIONS; nd++) dimens[nd] = wgrid[nd]->np_int_glob;
215 dataspace = H5Screate_simple(rank, dimens, NULL);
219 start[nd] = wgrid[nd]->
beg - wgrid[nd]->
nghost;
221 count[nd] = wgrid[nd]->
np_int;
224 err = H5Sselect_hyperslab(dataspace, H5S_SELECT_SET, start, stride, count, NULL);
227 for (nd = 0; nd <
DIMENSIONS; nd++) dimens[nd] = wgrid[nd]->np_tot;
229 memspace = H5Screate_simple(rank,dimens,NULL);
232 start[nd] = wgrid[nd]->
nghost;
234 count[nd] = wgrid[nd]->
np_int;
236 err = H5Sselect_hyperslab(memspace,H5S_SELECT_SET, start, stride, count, NULL);
243 plist_id_mpiio = H5Pcreate(H5P_DATASET_XFER);
244 H5Pset_dxpl_mpio(plist_id_mpiio,H5FD_MPIO_COLLECTIVE);
247 for (nv = 0; nv < output->
nvar; nv++) {
254 dataset = H5Dcreate(group, output->
var_name[nv], H5T_NATIVE_DOUBLE,
255 dataspace, H5P_DEFAULT);
257 err = H5Dwrite(dataset, H5T_NATIVE_DOUBLE, memspace, dataspace,
258 plist_id_mpiio, output->
V[nv][0][0]);
260 err = H5Dwrite(dataset, H5T_NATIVE_DOUBLE, memspace, dataspace,
261 H5P_DEFAULT, output->
V[nv][0][0]);
267 dataset = H5Dcreate(group, output->
var_name[nv], H5T_NATIVE_FLOAT,
268 dataspace, H5P_DEFAULT);
271 err = H5Dwrite(dataset, H5T_NATIVE_FLOAT, memspace,
272 dataspace, plist_id_mpiio, Vpt);
274 err = H5Dwrite(dataset, H5T_NATIVE_FLOAT, memspace,
275 dataspace, H5P_DEFAULT, Vpt);
282 H5Pclose(plist_id_mpiio);
302 group = H5Gcreate(timestep,
"stag_vars", 0);
305 plist_id_mpiio = H5Pcreate(H5P_DATASET_XFER);
306 H5Pset_dxpl_mpio(plist_id_mpiio,H5FD_MPIO_COLLECTIVE);
315 dimens[nd] = wgrid[nd]->
np_int_glob + (ns == (DIMENSIONS-1-nd));
317 dataspace = H5Screate_simple(rank, dimens, NULL);
321 start[nd] = wgrid[nd]->
beg - wgrid[nd]->
nghost;
323 count[nd] = wgrid[nd]->
np_int;
324 if (ns == DIMENSIONS-1-nd){
325 if (grid[ns].lbound != 0) count[nd] += 1;
329 err = H5Sselect_hyperslab(dataspace, H5S_SELECT_SET,
330 start, stride, count, NULL);
334 dimens[nd] = wgrid[nd]->
np_tot + (ns==(DIMENSIONS-1-nd));
336 memspace = H5Screate_simple(rank,dimens,NULL);
339 start[nd] = wgrid[nd]->
nghost;
341 count[nd] = wgrid[nd]->
np_int;
342 if (ns == (DIMENSIONS-1-nd) && grid[ns].lbound != 0) {
348 err = H5Sselect_hyperslab(memspace,H5S_SELECT_SET, start, stride, count, NULL);
351 dataset = H5Dcreate(group, output->
var_name[
NVAR+ns], H5T_NATIVE_DOUBLE,
352 dataspace, H5P_DEFAULT);
354 err = H5Dwrite(dataset, H5T_NATIVE_DOUBLE, memspace, dataspace,
355 plist_id_mpiio, output->
V[
NVAR+ns][0][0]);
357 err = H5Dwrite(dataset, H5T_NATIVE_DOUBLE, memspace, dataspace,
358 H5P_DEFAULT, output->
V[
NVAR+ns][0][0]);
363 dataset = H5Dcreate(group, output->
var_name[
NVAR+ns], H5T_NATIVE_FLOAT,
364 dataspace, H5P_DEFAULT);
367 err = H5Dwrite(dataset, H5T_NATIVE_FLOAT, memspace,
368 dataspace, plist_id_mpiio, Vpt);
370 err = H5Dwrite(dataset, H5T_NATIVE_FLOAT, memspace,
371 dataspace, H5P_DEFAULT, Vpt);
381 H5Pclose(plist_id_mpiio);
388 group = H5Gcreate(file_identifier,
"cell_coords", 0);
390 for (nd = 0; nd <
DIMENSIONS; nd++) dimens[nd] = wgrid[nd]->np_int_glob;
391 dataspace = H5Screate_simple(rank, dimens, NULL);
395 start[nd] = wgrid[nd]->
beg - wgrid[nd]->
nghost;
397 count[nd] = wgrid[nd]->
np_int;
399 err = H5Sselect_hyperslab(dataspace, H5S_SELECT_SET,
400 start, stride, count, NULL);
403 for (nd = 0; nd <
DIMENSIONS; nd++) dimens[nd] = wgrid[nd]->np_int;
404 memspace = H5Screate_simple(rank,dimens,NULL);
411 plist_id_mpiio = H5Pcreate(H5P_DATASET_XFER);
412 H5Pset_dxpl_mpio(plist_id_mpiio,H5FD_MPIO_COLLECTIVE);
415 for (nc = 0; nc < 3; nc++) {
416 dataset = H5Dcreate(group, cname[nc], H5T_NATIVE_FLOAT, dataspace, H5P_DEFAULT);
419 err = H5Dwrite(dataset, H5T_NATIVE_FLOAT, memspace,
420 dataspace, plist_id_mpiio, cell_coords[nc][0][0]);
422 err = H5Dwrite(dataset, H5T_NATIVE_FLOAT, memspace,
423 dataspace, H5P_DEFAULT, cell_coords[nc][0][0]);
430 H5Pclose(plist_id_mpiio);
436 group = H5Gcreate(file_identifier,
"node_coords", 0);
438 for (nd = 0; nd <
DIMENSIONS; nd++) dimens[nd] = wgrid[nd]->np_int_glob+1;
440 dataspace = H5Screate_simple(rank, dimens, NULL);
444 start[nd] = wgrid[nd]->
beg - wgrid[nd]->
nghost;
446 count[nd] = wgrid[nd]->
np_int;
447 if (wgrid[nd]->
rbound != 0) count[nd] += 1;
450 err = H5Sselect_hyperslab(dataspace, H5S_SELECT_SET,
451 start, stride, count, NULL);
455 dimens[nd] = wgrid[nd]->
np_int;
456 if (wgrid[nd]->
rbound != 0) dimens[nd] += 1;
458 memspace = H5Screate_simple(rank,dimens,NULL);
465 plist_id_mpiio = H5Pcreate(H5P_DATASET_XFER);
466 H5Pset_dxpl_mpio(plist_id_mpiio,H5FD_MPIO_COLLECTIVE);
469 for (nc = 0; nc < 3; nc++) {
470 dataset = H5Dcreate(group, cname[nc], H5T_NATIVE_FLOAT, dataspace, H5P_DEFAULT);
473 err = H5Dwrite(dataset, H5T_NATIVE_FLOAT, memspace,
474 dataspace, plist_id_mpiio, node_coords[nc][0][0]);
476 err = H5Dwrite(dataset, H5T_NATIVE_FLOAT, memspace,
477 dataspace, H5P_DEFAULT, node_coords[nc][0][0]);
484 H5Pclose(plist_id_mpiio);
490 H5Fclose(file_identifier);
497 sprintf(xmfext,
"dbl.xmf");
500 sprintf(xmfext,
"flt.xmf");
504 sprintf (filenamexmf,
"%s/data.%04d.%s", output->
dir, output->
nfile, xmfext);
506 fxmf = fopen(filenamexmf,
"w");
507 fprintf(fxmf,
"<?xml version=\"1.0\" ?>\n");
508 fprintf(fxmf,
"<!DOCTYPE Xdmf SYSTEM \"Xdmf.dtd\" []>\n");
509 fprintf(fxmf,
"<Xdmf Version=\"2.0\">\n");
510 fprintf(fxmf,
" <Domain>\n");
511 fprintf(fxmf,
" <Grid Name=\"node_mesh\" GridType=\"Uniform\">\n");
512 fprintf(fxmf,
" <Time Value=\"%12.6e\"/>\n",
g_time);
514 fprintf(fxmf,
" <Topology TopologyType=\"2DSMesh\" NumberOfElements=\"%d %d\"/>\n",
515 wgrid[0]->np_int_glob+1, wgrid[1]->np_int_glob+1);
516 fprintf(fxmf,
" <Geometry GeometryType=\"X_Y\">\n");
518 fprintf(fxmf,
" <DataItem Dimensions=\"%d %d\" NumberType=\"Float\" Precision=\"4\" Format=\"HDF\">\n",
519 wgrid[0]->np_int_glob+1, wgrid[1]->np_int_glob+1 );
520 fprintf(fxmf,
" %s:/node_coords/%s\n",filename,cname[nd]);
521 fprintf(fxmf,
" </DataItem>\n");
523 #elif DIMENSIONS == 3
524 fprintf(fxmf,
" <Topology TopologyType=\"3DSMesh\" NumberOfElements=\"%d %d %d\"/>\n",
525 wgrid[0]->np_int_glob+1, wgrid[1]->np_int_glob+1, wgrid[2]->np_int_glob+1);
526 fprintf(fxmf,
" <Geometry GeometryType=\"X_Y_Z\">\n");
528 fprintf(fxmf,
" <DataItem Dimensions=\"%d %d %d\" NumberType=\"Float\" Precision=\"4\" Format=\"HDF\">\n",
529 wgrid[0]->np_int_glob+1, wgrid[1]->np_int_glob+1, wgrid[2]->np_int_glob+1 );
530 fprintf(fxmf,
" %s:/node_coords/%s\n",filename,cname[nd]);
531 fprintf(fxmf,
" </DataItem>\n");
534 fprintf(fxmf,
" </Geometry>\n");
535 for (nv = 0; nv < output->
nvar; nv++) {
537 fprintf(fxmf,
" <Attribute Name=\"%s\" AttributeType=\"Scalar\" Center=\"Cell\">\n",
540 fprintf(fxmf,
" <DataItem Dimensions=\"%d %d\" NumberType=\"Float\" Precision=\"%d\" Format=\"HDF\">\n",
541 wgrid[0]->np_int_glob, wgrid[1]->np_int_glob, nprec);
542 #elif DIMENSIONS == 3
543 fprintf(fxmf,
" <DataItem Dimensions=\"%d %d %d\" NumberType=\"Float\" Precision=\"%d\" Format=\"HDF\">\n",
544 wgrid[0]->np_int_glob, wgrid[1]->np_int_glob, wgrid[2]->np_int_glob, nprec);
546 fprintf(fxmf,
" %s:%s/vars/%s\n",filename,tstepname,output->
var_name[nv]);
547 fprintf(fxmf,
" </DataItem>\n");
548 fprintf(fxmf,
" </Attribute>\n");
551 fprintf(fxmf,
" <Attribute Name=\"%s\" AttributeType=\"Scalar\" Center=\"Cell\">\n",
554 fprintf(fxmf,
" <DataItem Dimensions=\"%d %d\" NumberType=\"Float\" Precision=\"4\" Format=\"HDF\">\n",
555 wgrid[0]->np_int_glob, wgrid[1]->np_int_glob);
556 #elif DIMENSIONS == 3
557 fprintf(fxmf,
" <DataItem Dimensions=\"%d %d %d\" NumberType=\"Float\" Precision=\"4\" Format=\"HDF\">\n",
558 wgrid[0]->np_int_glob, wgrid[1]->np_int_glob, wgrid[2]->np_int_glob);
560 fprintf(fxmf,
" %s:/cell_coords/%s\n",filename,cname[nd]);
561 fprintf(fxmf,
" </DataItem>\n");
562 fprintf(fxmf,
" </Attribute>\n");
564 fprintf(fxmf,
" </Grid>\n");
565 fprintf(fxmf,
" </Domain>\n");
566 fprintf(fxmf,
"</Xdmf>\n");
573 MPI_Barrier (MPI_COMM_WORLD);
576 print1 (
" [%5.2f sec]",difftime(tend,tbeg));
double *** V[64]
pointer to arrays being written - same for all
long int NX1
Number of interior zones in the X1 directions (boundaries excluded) for the local processor...
void print1(const char *fmt,...)
int nfile
current number being saved - one per output
int np_int_glob
Total number of points in the global domain (boundaries excluded).
int rbound
Same as lbound, but for the right edge of the grid.
char dir[256]
output directory name
float *** Convert_dbl2flt(double ***Vdbl, double unit, int swap_endian)
char ext[8]
output extension
long int NX2
Number of interior zones in the X2 directions (boundaries excluded) for the local processor...
int beg
Global start index for the local array.
#define ARRAY_4D(nx, ny, nz, nv, type)
int nghost
Number of ghost zones.
int * dump_var
select vars being written - one per output
D_EXPAND(tot/[n]=(double) grid[IDIR].np_int_glob;, tot/[n]=(double) grid[JDIR].np_int_glob;, tot/[n]=(double) grid[KDIR].np_int_glob;)
double g_time
The current integration time.
int type
output format (DBL, FLT, ...) - one per output
int * stag_var
centered or staggered variable - same for all
long int KBEG
Lower grid index of the computational domain in the the X3 direction for the local processor...
int np_tot
Total number of points in the local domain (boundaries included).
long int JBEG
Lower grid index of the computational domain in the the X2 direction for the local processor...
#define QUIT_PLUTO(e_code)
long int NX3
Number of interior zones in the X3 directions (boundaries excluded) for the local processor...
char ** var_name
variable names - same for all
long int IBEG
Lower grid index of the computational domain in the the X1 direction for the local processor...
int np_int
Total number of points in the local domain (boundaries excluded).