MACSio  0.9
Multi-purpose, Application-Centric, Scalable I/O Proxy App
 All Data Structures Files Functions Variables Typedefs Enumerations Enumerator Macros Modules Pages
macsio_silo.c
Go to the documentation of this file.
1 /*
2 Copyright (c) 2015, Lawrence Livermore National Security, LLC.
3 Produced at the Lawrence Livermore National Laboratory.
4 Written by Mark C. Miller
5 
6 LLNL-CODE-676051. All rights reserved.
7 
8 This file is part of MACSio
9 
10 Please also read the LICENSE file at the top of the source code directory or
11 folder hierarchy.
12 
13 This program is free software; you can redistribute it and/or modify it under
14 the terms of the GNU General Public License (as published by the Free Software
15 Foundation) version 2, dated June 1991.
16 
17 This program is distributed in the hope that it will be useful, but WITHOUT
18 ANY WARRANTY; without even the IMPLIED WARRANTY OF MERCHANTABILITY or FITNESS
19 FOR A PARTICULAR PURPOSE. See the terms and conditions of the GNU General
20 Public License for more details.
21 
22 You should have received a copy of the GNU General Public License along with
23 this program; if not, write to the Free Software Foundation, Inc., 59 Temple
24 Place, Suite 330, Boston, MA 02111-1307 USA
25 */
26 
27 #include <errno.h>
28 #include <stdio.h>
29 #include <stdlib.h>
30 #include <string.h>
31 #ifdef _WIN32
32  #ifndef WINDOWS_LEAN_AND_MEAN
33  #define WINDOWS_LEAN_AND_MEAN
34  #endif
35  #include <windows.h>
36 #endif
37 
38 #include <macsio_clargs.h>
39 #include <macsio_data.h>
40 #include <macsio_iface.h>
41 #include <macsio_log.h>
42 #include <macsio_main.h>
43 #include <macsio_mif.h>
44 #include <macsio_utils.h>
45 
46 #include <silo.h>
47 
48 #ifdef HAVE_MPI
49 #include <mpi.h>
50 #endif
51 
62 #define NARRVALS(Arr) (sizeof(Arr)/sizeof(Arr[0]))
63 
64 /*
65  * BEGIN StringToDriver utility code from Silo's tests
66  */
67 
68 #define CHECK_SYMBOL(A) if (!strncmp(str, #A, strlen(str))) return A
69 
70 #define CHECK_SYMBOLN_INT(A) \
71 if (!strncmp(tok, #A, strlen(#A))) \
72 { \
73  int n = sscanf(tok, #A"=%d", &driver_ints[driver_nints]);\
74  if (n == 1) \
75  { \
76  DBAddOption(opts, A, &driver_ints[driver_nints]);\
77  driver_nints++; \
78  got_it = 1; \
79  } \
80 }
81 
82 #define CHECK_SYMBOLN_STR(A) \
83 if (!strncmp(tok, #A, strlen(#A))) \
84 { \
85  driver_strs[driver_nstrs] = strdup(&tok[strlen(#A)]+1);\
86  DBAddOption(opts, A, driver_strs[driver_nstrs]); \
87  driver_nstrs++; \
88  got_it = 1; \
89 }
90 
91 #define CHECK_SYMBOLN_SYM(A) \
92 if (!strncmp(tok, #A, strlen(#A))) \
93 { \
94  driver_ints[driver_nints] = StringToDriver(&tok[strlen(#A)]+1);\
95  DBAddOption(opts, A, &driver_ints[driver_nints]); \
96  driver_nints++; \
97  got_it = 1; \
98 }
99 
100 static DBoptlist *driver_opts[] = {0,0,0,0,0,0,0,0,0,0};
101 static int driver_opts_ids[] = {-1,-1,-1,-1,-1,-1,-1,-1,-1,-1};
102 static int driver_ints[100];
103 static int driver_nints = 0;
104 static char *driver_strs[] = {0,0,0,0,0,0,0,0,0,0};
105 static int driver_nstrs = 0;
106 static const int driver_nopts = sizeof(driver_opts)/sizeof(driver_opts[0]);
107 
108 static void CleanupDriverStuff()
109 {
110  int i;
111  for (i = 0; i < driver_nopts; i++)
112  {
113  if (driver_opts_ids[i] != -1) DBUnregisterFileOptionsSet(driver_opts_ids[i]);
114  if (driver_opts[i]) DBFreeOptlist(driver_opts[i]);
115  }
116  for (i = 0; i < sizeof(driver_strs)/sizeof(driver_strs[0]); i++)
117  if (driver_strs[i]) free(driver_strs[i]);
118 }
119 
120 static void MakeDriverOpts(DBoptlist **_opts, int *opts_id)
121 {
122  DBoptlist *opts = DBMakeOptlist(30);
123  int i;
124 
125  for (i = 0; i < driver_nopts; i++)
126  {
127  if (driver_opts[i] == 0)
128  {
129  driver_opts[i] = opts;
130  break;
131  }
132  }
133 
134  *_opts = opts;
135  *opts_id = DBRegisterFileOptionsSet(opts);
136  driver_opts_ids[i] = *opts_id;
137 }
138 
139 static int StringToDriver(const char *str)
140 {
141  DBoptlist *opts = 0;
142  int opts_id = -1;
143 
144  CHECK_SYMBOL(DB_PDB);
145  CHECK_SYMBOL(DB_PDBP);
146  CHECK_SYMBOL(DB_HDF5);
147  CHECK_SYMBOL(DB_HDF5_SEC2);
148  CHECK_SYMBOL(DB_HDF5_STDIO);
149  CHECK_SYMBOL(DB_HDF5_CORE);
150  CHECK_SYMBOL(DB_HDF5_SPLIT);
151  CHECK_SYMBOL(DB_HDF5_MPIO);
152  CHECK_SYMBOL(DB_HDF5_MPIP);
153  CHECK_SYMBOL(DB_HDF5_LOG);
154  CHECK_SYMBOL(DB_HDF5_DIRECT);
155  CHECK_SYMBOL(DB_HDF5_FAMILY);
156  CHECK_SYMBOL(DB_HDF5_SILO);
157 
158  CHECK_SYMBOL(DB_FILE_OPTS_H5_DEFAULT_DEFAULT);
159  CHECK_SYMBOL(DB_FILE_OPTS_H5_DEFAULT_SEC2);
160  CHECK_SYMBOL(DB_FILE_OPTS_H5_DEFAULT_STDIO);
161  CHECK_SYMBOL(DB_FILE_OPTS_H5_DEFAULT_CORE);
162  CHECK_SYMBOL(DB_FILE_OPTS_H5_DEFAULT_LOG);
163  CHECK_SYMBOL(DB_FILE_OPTS_H5_DEFAULT_SPLIT);
164  CHECK_SYMBOL(DB_FILE_OPTS_H5_DEFAULT_DIRECT);
165  CHECK_SYMBOL(DB_FILE_OPTS_H5_DEFAULT_FAMILY);
166  CHECK_SYMBOL(DB_FILE_OPTS_H5_DEFAULT_MPIP);
167  CHECK_SYMBOL(DB_FILE_OPTS_H5_DEFAULT_MPIO);
168  CHECK_SYMBOL(DB_FILE_OPTS_H5_DEFAULT_SILO);
169 
170  CHECK_SYMBOL(DB_H5VFD_DEFAULT);
171  CHECK_SYMBOL(DB_H5VFD_SEC2);
172  CHECK_SYMBOL(DB_H5VFD_STDIO);
173  CHECK_SYMBOL(DB_H5VFD_CORE);
174  CHECK_SYMBOL(DB_H5VFD_LOG);
175  CHECK_SYMBOL(DB_H5VFD_SPLIT);
176  CHECK_SYMBOL(DB_H5VFD_DIRECT);
177  CHECK_SYMBOL(DB_H5VFD_FAMILY);
178  CHECK_SYMBOL(DB_H5VFD_MPIO);
179  CHECK_SYMBOL(DB_H5VFD_MPIP);
180  CHECK_SYMBOL(DB_H5VFD_SILO);
181 
182  if (!strncmp(str, "DB_HDF5_OPTS(", 13))
183  {
184  char *tok, *tmpstr;;
185 
186  MakeDriverOpts(&opts, &opts_id);
187 
188  tmpstr = strdup(&str[13]);
189  errno = 0;
190  tok = strtok(tmpstr, ",)");
191 
192  while (tok)
193  {
194  int got_it = 0;
195 
196  CHECK_SYMBOLN_SYM(DBOPT_H5_VFD)
197  CHECK_SYMBOLN_SYM(DBOPT_H5_RAW_FILE_OPTS)
198  CHECK_SYMBOLN_STR(DBOPT_H5_RAW_EXTENSION)
199  CHECK_SYMBOLN_SYM(DBOPT_H5_META_FILE_OPTS)
200  CHECK_SYMBOLN_STR(DBOPT_H5_META_EXTENSION)
201  CHECK_SYMBOLN_INT(DBOPT_H5_CORE_ALLOC_INC)
202  CHECK_SYMBOLN_INT(DBOPT_H5_CORE_NO_BACK_STORE)
203  CHECK_SYMBOLN_INT(DBOPT_H5_META_BLOCK_SIZE)
204  CHECK_SYMBOLN_INT(DBOPT_H5_SMALL_RAW_SIZE)
205  CHECK_SYMBOLN_INT(DBOPT_H5_ALIGN_MIN)
206  CHECK_SYMBOLN_INT(DBOPT_H5_ALIGN_VAL)
207  CHECK_SYMBOLN_INT(DBOPT_H5_DIRECT_MEM_ALIGN)
208  CHECK_SYMBOLN_INT(DBOPT_H5_DIRECT_BLOCK_SIZE)
209  CHECK_SYMBOLN_INT(DBOPT_H5_DIRECT_BUF_SIZE)
210  CHECK_SYMBOLN_STR(DBOPT_H5_LOG_NAME)
211  CHECK_SYMBOLN_INT(DBOPT_H5_LOG_BUF_SIZE)
212  CHECK_SYMBOLN_INT(DBOPT_H5_SIEVE_BUF_SIZE)
213  CHECK_SYMBOLN_INT(DBOPT_H5_CACHE_NELMTS)
214  CHECK_SYMBOLN_INT(DBOPT_H5_CACHE_NBYTES)
215  CHECK_SYMBOLN_INT(DBOPT_H5_FAM_SIZE)
216  CHECK_SYMBOLN_SYM(DBOPT_H5_FAM_FILE_OPTS)
217  CHECK_SYMBOLN_INT(DBOPT_H5_SILO_BLOCK_SIZE)
218  CHECK_SYMBOLN_INT(DBOPT_H5_SILO_BLOCK_COUNT)
219  CHECK_SYMBOLN_INT(DBOPT_H5_SILO_LOG_STATS)
220  CHECK_SYMBOLN_INT(DBOPT_H5_SILO_USE_DIRECT)
221  CHECK_SYMBOLN_STR(DB_FILE_OPTS_H5_DEFAULT_DEFAULT)
222  CHECK_SYMBOLN_STR(DB_FILE_OPTS_H5_DEFAULT_SEC2)
223  CHECK_SYMBOLN_STR(DB_FILE_OPTS_H5_DEFAULT_STDIO)
224  CHECK_SYMBOLN_STR(DB_FILE_OPTS_H5_DEFAULT_CORE)
225  CHECK_SYMBOLN_STR(DB_FILE_OPTS_H5_DEFAULT_LOG)
226  CHECK_SYMBOLN_STR(DB_FILE_OPTS_H5_DEFAULT_SPLIT)
227  CHECK_SYMBOLN_STR(DB_FILE_OPTS_H5_DEFAULT_DIRECT)
228  CHECK_SYMBOLN_STR(DB_FILE_OPTS_H5_DEFAULT_FAMILY)
229  CHECK_SYMBOLN_STR(DB_FILE_OPTS_H5_DEFAULT_MPIP)
230  CHECK_SYMBOLN_STR(DB_FILE_OPTS_H5_DEFAULT_MPIO)
231  CHECK_SYMBOLN_STR(DB_FILE_OPTS_H5_DEFAULT_SILO)
232 
233  if (!got_it)
234  {
235  fprintf(stderr, "Unable to determine driver from string \"%s\"\n", tok);
236  exit(-1);
237  }
238 
239  tok = strtok(0, ",)");
240  if (errno != 0)
241  {
242  fprintf(stderr, "Unable to determine driver from string \"%s\"\n", tok);
243  exit(-1);
244  }
245  }
246 
247  free(tmpstr);
248 
249  return DB_HDF5_OPTS(opts_id);
250  }
251 
252  fprintf(stderr, "Unable to determine driver from string \"%s\"\n", str);
253  exit(-1);
254 }
255 
256 /*
257  * END StringToDriver utility code from Silo's tests
258  */
259 
260 /* the name you want to assign to the interface */
261 static char const *iface_name = "silo";
262 static char const *iface_ext = "silo";
263 
264 static const char *filename;
265 static int has_mesh = 0;
266 static int driver = DB_HDF5;
267 static int show_all_errors = FALSE;
268 #warning MOVE LOG HANDLE TO IO CONTEXT
269 
270 static int process_args(int argi, int argc, char *argv[])
271 {
272  const MACSIO_CLARGS_ArgvFlags_t argFlags = {
276  char *driver_str = 0;
277  char *compression_str = 0;
278  int cksums = 0;
279  int hdf5friendly = 0;
280  int show_all_errors = 0;
281 
282  MACSIO_CLARGS_ProcessCmdline(0, argFlags, argi, argc, argv,
283  "--driver %s", "DB_HDF5",
284  "Specify Silo's I/O driver (DB_PDB|DB_HDF5 or variants)",
285  &driver_str,
286  "--checksums", "",
287  "Enable checksum checks",
288  &cksums,
289  "--hdf5friendly", "",
290  "Generate HDF5 friendly files",
291  &hdf5friendly,
292  "--compression %s", MACSIO_CLARGS_NODEFAULT,
293  "Specify compression method to use",
294  &compression_str,
295  "--show-all-errors", "",
296  "Show all errors Silo encounters",
297  &show_all_errors,
299 
300  if (driver_str)
301  {
302  driver = StringToDriver(driver_str);
303  free(driver_str);
304  }
305  if (compression_str)
306  {
307  if (*compression_str)
308  DBSetCompression(compression_str);
309  free(compression_str);
310  }
311  DBSetEnableChecksums(cksums);
312  DBSetFriendlyHDF5Names(hdf5friendly);
313  DBShowErrors(show_all_errors?DB_ALL_AND_DRVR:DB_ALL, NULL);
314 
315  return 0;
316 }
317 
318 static void *CreateSiloFile(const char *fname, const char *nsname, void *userData)
319 {
320  int driverId = *((int*) userData);
321  DBfile *siloFile = DBCreate(fname, DB_CLOBBER, DB_LOCAL, "macsio output file", driverId);
322  if (siloFile && nsname)
323  {
324  DBMkDir(siloFile, nsname);
325  DBSetDir(siloFile, nsname);
326  }
327  return (void *) siloFile;
328 }
329 
330 static void *OpenSiloFile(const char *fname, const char *nsname, MACSIO_MIF_ioFlags_t ioFlags,
331  void *userData)
332 {
333  int driverId = *((int*) userData);
334  DBfile *siloFile = DBOpen(fname, driverId, ioFlags.do_wr ? DB_APPEND : DB_READ);
335  if (siloFile && nsname)
336  {
337  if (ioFlags.do_wr)
338  DBMkDir(siloFile, nsname);
339  DBSetDir(siloFile, nsname);
340  }
341  return (void *) siloFile;
342 }
343 
344 static void CloseSiloFile(void *file, void *userData)
345 {
346  DBfile *siloFile = (DBfile *) file;
347  if (siloFile)
348  DBClose(siloFile);
349 }
350 
351 static void write_rect_mesh_part(DBfile *dbfile, json_object *part)
352 {
353  json_object *coordobj;
354  char const *coordnames[] = {"X","Y","Z"};
355  void const *coords[3];
356  int ndims = JsonGetInt(part, "Mesh/GeomDim");
357  int dims[3] = {1,1,1};
358  int dimsz[3] = {1,1,1};
359 
360  coordobj = JsonGetObj(part, "Mesh/Coords/XAxisCoords");
361  coords[0] = json_object_extarr_data(coordobj);
362  dims[0] = JsonGetInt(part, "Mesh/LogDims", 0);
363  dimsz[0] = dims[0]-1;
364  if (ndims > 1)
365  {
366  coordobj = JsonGetObj(part, "Mesh/Coords/YAxisCoords");
367  coords[1] = json_object_extarr_data(coordobj);
368  dims[1] = JsonGetInt(part, "Mesh/LogDims", 1);
369  dimsz[1] = dims[1]-1;
370 
371  }
372  if (ndims > 2)
373  {
374  coordobj = JsonGetObj(part, "Mesh/Coords/ZAxisCoords");
375  coords[2] = json_object_extarr_data(coordobj);
376  dims[2] = JsonGetInt(part, "Mesh/LogDims", 2);
377  dimsz[2] = dims[2]-1;
378  }
379 
380  DBPutQuadmesh(dbfile, "mesh", (char**) coordnames, coords,
381  dims, ndims, DB_DOUBLE, DB_COLLINEAR, 0);
382 
383  json_object *vars_array = JsonGetObj(part, "Vars");
384  for (int i = 0; i < json_object_array_length(vars_array); i++)
385  {
386  json_object *varobj = json_object_array_get_idx(vars_array, i);
387  int cent = strcmp(JsonGetStr(varobj, "centering"),"zone")?DB_NODECENT:DB_ZONECENT;
388  int *d = cent==DB_NODECENT?dims:dimsz;
389  json_object *dataobj = JsonGetObj(varobj, "data");
390  int dtype = json_object_extarr_type(dataobj)==json_extarr_type_flt64?DB_DOUBLE:DB_INT;
391 
392  DBPutQuadvar1(dbfile, JsonGetStr(varobj, "name"), "mesh",
393  (void *)json_object_extarr_data(dataobj), d, ndims, 0, 0, dtype, cent, 0);
394  }
395 }
396 
397 static void write_mesh_part(DBfile *dbfile, json_object *part)
398 {
399  if (!strcmp(JsonGetStr(part, "Mesh/MeshType"), "rectilinear"))
400  write_rect_mesh_part(dbfile, part);
401 }
402 
403 static void WriteMultiXXXObjects(json_object *main_obj, DBfile *siloFile, int dumpn, MACSIO_MIF_baton_t *bat)
404 {
405  int i, j;
406  char const *file_ext = JsonGetStr(main_obj, "clargs/fileext");
407  char const *file_base = JsonGetStr(main_obj, "clargs/filebase");
408  int numChunks = JsonGetInt(main_obj, "problem/global/TotalParts");
409  char **blockNames = (char **) malloc(numChunks * sizeof(char*));
410  int *blockTypes = (int *) malloc(numChunks * sizeof(int));
411 
412  /* Go to root directory in the silo file */
413  DBSetDir(siloFile, "/");
414 
415  /* Construct the lists of individual object names */
416  for (i = 0; i < numChunks; i++)
417  {
418  int rank_owning_chunk = MACSIO_DATA_GetRankOwningPart(main_obj, i);
419  int groupRank = MACSIO_MIF_RankOfGroup(bat, rank_owning_chunk);
420  blockNames[i] = (char *) malloc(1024);
421  if (groupRank == 0)
422  {
423  /* this mesh block is in the file 'root' owns */
424  sprintf(blockNames[i], "/domain_%07d/mesh", i);
425  }
426  else
427  {
428 #warning USE SILO NAMESCHEMES INSTEAD
429  sprintf(blockNames[i], "%s_silo_%05d_%03d.%s:/domain_%07d/mesh",
430  JsonGetStr(main_obj, "clargs/filebase"),
431  groupRank, dumpn,
432  JsonGetStr(main_obj, "clargs/fileext"),
433  i);
434  }
435  blockTypes[i] = DB_QUADMESH;
436  }
437 
438  /* Write the multi-block objects */
439  DBPutMultimesh(siloFile, "multi_mesh", numChunks, blockNames, blockTypes, 0);
440 
441  json_object *first_part = JsonGetObj(main_obj, "problem/parts", 0);
442  json_object *vars_array = JsonGetObj(first_part, "Vars");
443  int numVars = json_object_array_length(vars_array);
444  for (j = 0; j < numVars; j++)
445  {
446  for (i = 0; i < numChunks; i++)
447  {
448  int rank_owning_chunk = MACSIO_DATA_GetRankOwningPart(main_obj, i);
449  int groupRank = MACSIO_MIF_RankOfGroup(bat, rank_owning_chunk);
450  if (groupRank == 0)
451  {
452  /* this mesh block is in the file 'root' owns */
453  sprintf(blockNames[i], "/domain_%07d/%s", i,
454  JsonGetStr(vars_array, "", j, "name"));
455  }
456  else
457  {
458 #warning USE SILO NAMESCHEMES INSTEAD
459  sprintf(blockNames[i], "%s_silo_%05d_%03d.%s:/domain_%07d/%s",
460  JsonGetStr(main_obj, "clargs/filebase"),
461  groupRank,
462  dumpn,
463  JsonGetStr(main_obj, "clargs/fileext"),
464  i,
465  JsonGetStr(vars_array, "", j, "name"));
466  }
467  blockTypes[i] = DB_QUADVAR;
468  }
469 
470  /* Write the multi-block objects */
471  DBPutMultivar(siloFile, JsonGetStr(vars_array, "", j, "name"), numChunks, blockNames, blockTypes, 0);
472 
473 #warning WRITE MULTIBLOCK DOMAIN ASSIGNMENT AS A TINY QUADMESH OF SAME PHYSICAL SIZE OF MESH BUT FEWER ZONES
474 
475  }
476 
477  /* Clean up */
478  for (i = 0; i < numChunks; i++)
479  free(blockNames[i]);
480  free(blockNames);
481  free(blockTypes);
482 }
483 
484 static void WriteDecompMesh(json_object *main_obj, DBfile *siloFile, int dumpn, MACSIO_MIF_baton_t *bat)
485 {
486  int i, ndims = JsonGetInt(main_obj, "clargs/part_dim");
487  int total_parts = JsonGetInt(main_obj, "problem/global/TotalParts");
488  int zdims[3] = {0, 0, 0}, dims[3] = {0, 0, 0};
489  double bounds[6] = {0, 0, 0, 0, 0, 0};
490  char const *coordnames[] = {"X","Y","Z"};
491  double *coords[3];
492  double *color = (double *) malloc(total_parts * sizeof(double));
493 
494  if (ndims >= 1)
495  {
496  double *x, xdelta;
497  dims[0] = JsonGetInt(main_obj, "problem/global/PartsLogDims/0")+1;
498  bounds[0] = JsonGetInt(main_obj, "problem/global/Bounds/0");
499  bounds[3] = JsonGetInt(main_obj, "problem/global/Bounds/3");
500  x = (double*) malloc(dims[0] * sizeof(double));
501  xdelta = MACSIO_UTILS_XDelta(dims, bounds);
502  for (i = 0; i < dims[0]; i++)
503  x[i] = bounds[0] + i * xdelta;
504  coords[0] = x;
505  }
506 
507  if (ndims >= 2)
508  {
509  double *y, ydelta;
510  dims[1] = JsonGetInt(main_obj, "problem/global/PartsLogDims/1")+1;
511  bounds[1] = JsonGetInt(main_obj, "problem/global/Bounds/1");
512  bounds[4] = JsonGetInt(main_obj, "problem/global/Bounds/4");
513  y = (double*) malloc(dims[1] * sizeof(double));
514  ydelta = MACSIO_UTILS_YDelta(dims, bounds);
515  for (i = 0; i < dims[1]; i++)
516  y[i] = bounds[1] + i * ydelta;
517  coords[1] = y;
518  }
519 
520  if (ndims >= 3)
521  {
522  double *z, zdelta;
523  dims[2] = JsonGetInt(main_obj, "problem/global/PartsLogDims/2")+1;
524  bounds[2] = JsonGetInt(main_obj, "problem/global/Bounds/2");
525  bounds[5] = JsonGetInt(main_obj, "problem/global/Bounds/5");
526  z = (double*) malloc(dims[2] * sizeof(double));
527  zdelta = MACSIO_UTILS_ZDelta(dims, bounds);
528  for (i = 0; i < dims[2]; i++)
529  z[i] = bounds[2] + i * zdelta;
530  coords[2] = z;
531  }
532 
533  DBPutQuadmesh(siloFile, "part_mesh", (char**) coordnames, coords,
534  dims, ndims, DB_DOUBLE, DB_COLLINEAR, 0);
535 
536  for (i = 0; i < total_parts; i++)
537  color[i] = MACSIO_DATA_GetRankOwningPart(main_obj, i);
538  for (i = 0; i < ndims; i++)
539  zdims[i] = dims[i]-1;
540 
541  DBPutQuadvar1(siloFile, "decomp", "part_mesh", color,
542  zdims, ndims, NULL, 0, DB_DOUBLE, DB_ZONECENT, 0);
543 }
544 
545 #warning HOW IS A NEW DUMP CLASS HANDLED
546 #warning ADD TIMING LOGIC
547 static void main_dump(int argi, int argc, char **argv, json_object *main_obj, int dumpn, double dumpt)
548 {
549  DBfile *siloFile;
550  int numGroups = -1;
551  int rank, size;
552  char fileName[256];
553  MACSIO_MIF_baton_t *bat;
555  JsonGetInt(main_obj, "clargs/exercise_scr")&0x1};
556 
557  /* Without this barrier, I get strange behavior with Silo's MACSIO_MIF interface */
558 #warning CONFIRM THIS IS STILL NEEDED
559  mpi_errno = MPI_Barrier(MACSIO_MAIN_Comm);
560 
561  /* process cl args */
562  process_args(argi, argc, argv);
563 
564  rank = JsonGetInt(main_obj, "parallel/mpi_rank");
565  size = JsonGetInt(main_obj, "parallel/mpi_size");
566 
567 #warning MOVE TO A FUNCTION
568  /* ensure we're in MIF mode and determine the file count */
569  json_object *parfmode_obj = JsonGetObj(main_obj, "clargs/parallel_file_mode");
570  if (parfmode_obj)
571  {
572  json_object *modestr = json_object_array_get_idx(parfmode_obj, 0);
573  json_object *filecnt = json_object_array_get_idx(parfmode_obj, 1);
574  if (modestr && !strcmp(json_object_get_string(modestr), "SIF"))
575  {
576  MACSIO_LOG_MSG(Die, ("Silo plugin doesn't support SIF mode"));
577  }
578  else if (modestr && strcmp(json_object_get_string(modestr), "MIF"))
579  {
580  MACSIO_LOG_MSG(Warn, ("Ignoring non-standard MIF mode"));
581  }
582  if (filecnt)
583  numGroups = json_object_get_int(filecnt);
584  else
585  numGroups = size;
586  }
587  else
588  {
589  char const * modestr = JsonGetStr(main_obj, "clargs/parallel_file_mode");
590  if (!strcmp(modestr, "SIF"))
591  {
592  MACSIO_LOG_MSG(Die, ("Silo plugin doesn't support SIF mode"));
593  }
594  else if (!strcmp(modestr, "MIFMAX"))
595  numGroups = JsonGetInt(main_obj, "parallel/mpi_size");
596  else if (!strcmp(modestr, "MIFAUTO"))
597  {
598  /* Call utility to determine optimal file count */
599 #warning ADD UTILIT TO DETERMINE OPTIMAL FILE COUNT
600  }
601  }
602 
603  /* Initialize MACSIO_MIF, pass a pointer to the driver type as the user data. */
604  bat = MACSIO_MIF_Init(numGroups, ioFlags, MACSIO_MAIN_Comm, 1,
606 
607  /* Construct name for the silo file */
608 #warning CHANGE NAMING SCHEME SO LS WORKS BETTER
609  sprintf(fileName, "%s_silo_%05d_%03d.%s",
610  JsonGetStr(main_obj, "clargs/filebase"),
611  MACSIO_MIF_RankOfGroup(bat, rank),
612  dumpn,
613  JsonGetStr(main_obj, "clargs/fileext"));
614 
615  /* Wait for write access to the file. All processors call this.
616  * Some processors (the first in each group) return immediately
617  * with write access to the file. Other processors wind up waiting
618  * until they are given control by the preceeding processor in
619  * the group when that processor calls "HandOffBaton" */
620  siloFile = (DBfile *) MACSIO_MIF_WaitForBaton(bat, fileName, 0);
621 
622  json_object *parts = JsonGetObj(main_obj, "problem/parts");
623  int numParts = json_object_array_length(parts);
624 
625  for (int i = 0; i < numParts; i++)
626  {
627  char domain_dir[256];
628  json_object *this_part = JsonGetObj(main_obj, "problem/parts", i);
629 
630  snprintf(domain_dir, sizeof(domain_dir), "domain_%07d", JsonGetInt(this_part, "Mesh/ChunkID"));
631 
632  DBMkDir(siloFile, domain_dir);
633  DBSetDir(siloFile, domain_dir);
634 
635  write_mesh_part(siloFile, this_part);
636 
637  DBSetDir(siloFile, "..");
638  }
639 
640  /* If this is the 'root' processor, also write Silo's multi-XXX objects */
641  if (rank == 0)
642  {
643  WriteMultiXXXObjects(main_obj, siloFile, dumpn, bat);
644 
645 #warning DECOMP MESH SHOULDN'T BE INCLUDED IN PERFORMANCE NUMBERS
646  /* output a top-level quadmesh and vars to indicate processor decomp */
647  if (MACSIO_LOG_DebugLevel >= 2)
648  {
649  WriteDecompMesh(main_obj, siloFile, dumpn, bat);
650  }
651  }
652 
653  /* Hand off the baton to the next processor. This winds up closing
654  * the file so that the next processor that opens it can be assured
655  * of getting a consistent and up to date view of the file's contents. */
656  MACSIO_MIF_HandOffBaton(bat, siloFile);
657 
658  /* We're done using MACSIO_MIF, so finish it off */
659  MACSIO_MIF_Finish(bat);
660 }
661 
662 #warning TO BE MOVED TO SILO LIBRARY
663 static void DBSplitMultiName(char *mname, char **file, char **dir, char **obj)
664 {
665  char *_file, *_dir, *_obj;
666  char *colon = strchr(mname, ':');
667  char *lastslash = strrchr(mname, '/');
668 
669  /* Split the incoming string by inserting null chars */
670  if (colon) *colon = '\0';
671  if (lastslash) *lastslash = '\0';
672 
673  if (colon)
674  {
675  if (file) *file = mname;
676  if (lastslash)
677  {
678  if (dir) *dir = colon+1;
679  if (obj) *obj = lastslash+1;
680  }
681  else
682  {
683  if (dir) *dir = 0;
684  if (obj) *obj = colon+1;
685  }
686  }
687  else
688  {
689  if (file) *file = 0;
690  if (lastslash)
691  {
692  if (dir) *dir = mname;
693  if (obj) *obj = lastslash+1;
694  }
695  else
696  {
697  if (dir) *dir = 0;
698  if (obj) *obj = mname;
699  }
700  }
701 }
702 
703 #warning TO BE MOVED TO SILO LIBRARY
704 static char const *
705 DBGetFilename(DBfile const *f)
706 {
707  return f->pub.name;
708 }
709 
710 #warning SHOULD USE MACSIO_MIF FOR READ TOO BUT INTERFACE IS LACKING
711 static
712 void main_load(int argi, int argc, char **argv, char const *path, json_object *main_obj, json_object **data_read_obj)
713 {
714  int my_rank = JsonGetInt(main_obj, "parallel/mpi_rank");
715  int mpi_size = JsonGetInt(main_obj, "parallel/mpi_rank");
716  int i, num_parts, my_part_cnt, use_ns = 0, maxlen = 0, bcast_data[3];
717  char *all_meshnames, *my_meshnames;
718  char *var_names_list = strdup(JsonGetStr(main_obj, "clargs/read_vars"));
719  char *var_names_list_orig = var_names_list;
720  char *vname;
721  int *my_part_ids, *all_part_cnts;
722  DBfile *partFile = 0;
723  int silo_driver = DB_UNKNOWN;
724 
725  /* Open the root file */
726  if (my_rank == 0)
727  {
728  DBfile *rootFile = DBOpen(path, DB_UNKNOWN, DB_READ);
729  DBmultimesh *mm = DBGetMultimesh(rootFile, JsonGetStr(main_obj, "clargs/read_mesh"));
730 
731  /* Examine multimesh for count of mesh pieces and count of files */
732  num_parts = mm->nblocks;
733  use_ns = mm->block_ns ? 1 : 0;
734 
735  /* Reformat all the meshname strings to a single, long buffer */
736  if (!use_ns)
737  {
738  int i;
739  for (i = 0; i < num_parts; i++)
740  {
741  int len = strlen(mm->meshnames[i]);
742  if (len > maxlen) maxlen = len;
743  }
744  maxlen++; /* for nul char */
745  all_meshnames = (char *) calloc(num_parts * maxlen, sizeof(char));
746  for (i = 0; i < num_parts; i++)
747  strcpy(&all_meshnames[i*maxlen], mm->meshnames[i]);
748  }
749 
750  bcast_data[0] = num_parts;
751  bcast_data[1] = use_ns;
752  bcast_data[2] = maxlen;
753 
754  DBClose(rootFile);
755  }
756 #ifdef HAVE_MPI
757  MPI_Bcast(bcast_data, NARRVALS(bcast_data), MPI_INT, 0, MACSIO_MAIN_Comm);
758 #endif
759  num_parts = bcast_data[0];
760  use_ns = bcast_data[1];
761  maxlen = bcast_data[2];
762 
763 #if 0
764  MACSIO_DATA_SimpleAssignKPartsToNProcs(num_parts, mpi_size, my_rank, &all_part_cnts,
765  &my_part_cnt, &my_part_ids);
766 #endif
767 
768  my_meshnames = (char *) calloc(my_part_cnt * maxlen, sizeof(char));
769 
770  if (use_ns)
771  {
772  }
773 #ifdef HAVE_MPI
774  else
775  {
776  int *displs = 0;
777 
778  if (my_rank == 0)
779  {
780  displs = (int *) malloc(num_parts * sizeof(int));
781  displs[0] = 0;
782  for (i = 1; i < num_parts; i++)
783  displs[i] = displs[i-1] + all_part_cnts[i-1] * maxlen;
784  }
785 
786  /* MPI_scatter the block names or the external arrays for any namescheme */
787  MPI_Scatterv(all_meshnames, all_part_cnts, displs, MPI_CHAR,
788  my_meshnames, my_part_cnt, MPI_CHAR, 0, MACSIO_MAIN_Comm);
789 
790  if (my_rank == 0)
791  free(displs);
792  }
793 #endif
794 
795  /* Iterate finding correct file/dir combo and reading mesh pieces and variables */
796  for (i = 0; i < my_part_cnt; i++)
797  {
798  char *partFileName, *partDirName, *partObjName;
799  DBObjectType silo_objtype;
800 
801  DBSplitMultiName(&my_meshnames[i*maxlen], &partFileName, &partDirName, &partObjName);
802 
803  /* if filename is different from current, close current */
804  if (partFile && strcmp(DBGetFilename(partFile), partFileName))
805  DBClose(partFile);
806 
807  /* Open the file containing this part */
808  partFile = DBOpen(partFileName, silo_driver, DB_READ);
809  silo_driver = DBGetDriverType(partFile);
810 
811  DBSetDir(partFile, partDirName);
812 
813  /* Get the mesh part */
814  silo_objtype = (DBObjectType) DBInqMeshtype(partFile, partObjName);
815 
816  switch (silo_objtype)
817  {
818  case DB_QUADRECT:
819  case DB_QUADCURV:
820  case DB_QUADMESH:
821  {
822  DBquadmesh *qm = DBGetQuadmesh(partFile, partObjName);
823  break;
824  }
825  case DB_UCDMESH:
826  {
827  break;
828  }
829  case DB_POINTMESH:
830  {
831  break;
832  }
833  default:
834  {
835  }
836  }
837 
838  while (vname = strsep(&var_names_list, ", "))
839  {
840  DBObjectType silo_vartype = DBInqVarType(partFile, vname);
841  switch (silo_vartype)
842  {
843  case DB_QUADVAR:
844  {
845  DBquadvar *qv = DBGetQuadvar(partFile, vname);
846  break;
847  }
848  case DB_POINTVAR:
849  {
850  break;
851  }
852  case DB_UCDVAR:
853  {
854  break;
855  }
856  default: continue;
857  }
858  }
859  free(var_names_list_orig);
860 
861  /* Add the mesh part to the returned json object */
862 
863  DBSetDir(partFile, "/");
864  }
865 
866  free(my_meshnames);
867  if (my_rank == 0)
868  free(all_meshnames);
869 }
870 
872 {
873  MACSIO_IFACE_Handle_t iface;
874 
875  if (strlen(iface_name) >= MACSIO_IFACE_MAX_NAME)
876  MACSIO_LOG_MSG(Die, ("Interface name \"%s\" too long",iface_name));
877 
878  /* Take this slot in the map */
879  strcpy(iface.name, iface_name);
880  strcpy(iface.ext, iface_ext);
881 
882  /* Must define at least these two methods */
883  iface.dumpFunc = main_dump;
884  iface.loadFunc = main_load;
886 
887  if (!MACSIO_IFACE_Register(&iface))
888  MACSIO_LOG_MSG(Die, ("Failed to register interface \"%s\"", iface_name));
889 
890  return 0;
891 }
892 
893 /* this one statement is the only statement requiring compilation by
894  a C++ compiler. That is because it involves initialization and non
895  constant expressions (a function call in this case). This function
896  call is guaranteed to occur during *initialization* (that is before
897  even 'main' is called) and so will have the effect of populating the
898  iface_map array merely by virtue of the fact that this code is linked
899  with a main. */
901 
#define JsonGetObj(OBJ,...)
Get the object at specified path.
Definition: json_object.h:771
static int StringToDriver(const char *str)
Definition: macsio_silo.c:139
static void * OpenSiloFile(const char *fname, const char *nsname, MACSIO_MIF_ioFlags_t ioFlags, void *userData)
Definition: macsio_silo.c:330
int MACSIO_MIF_RankOfGroup(MACSIO_MIF_baton_t const *Bat, int rankInComm)
Rank of the group in which a given (global) rank exists.
Definition: macsio_mif.c:311
#define JsonGetStr(OBJ,...)
Get the string value at specified path.
Definition: json_object.h:765
unsigned int do_wr
Definition: macsio_mif.h:155
void * MACSIO_MIF_WaitForBaton(MACSIO_MIF_baton_t *Bat, char const *fname, char const *nsname)
Wait for exclusive access to the group's file.
Definition: macsio_mif.c:210
#define MACSIO_CLARGS_WARN
Definition: macsio_clargs.h:35
static int driver
Definition: macsio_silo.c:266
static int show_all_errors
Definition: macsio_silo.c:267
static int has_mesh
Definition: macsio_silo.c:265
int mpi_errno
Error code returned by most recent MPI call.
Definition: macsio_log.c:43
#define MACSIO_MIF_WRITE
Definition: macsio_mif.h:149
static char const * iface_ext
Definition: macsio_silo.c:262
void const * json_object_extarr_data(struct json_object *jso)
Definition: json_object.c:1206
static DBoptlist * driver_opts[]
Definition: macsio_silo.c:100
static char const * iface_name
Definition: macsio_silo.c:261
static char const * DBGetFilename(DBfile const *f)
Definition: macsio_silo.c:705
int MACSIO_CLARGS_ProcessCmdline(void **retobj, MACSIO_CLARGS_ArgvFlags_t flags, int argi, int argc, char **argv,...)
int MACSIO_IFACE_Register(MACSIO_IFACE_Handle_t const *iface)
Definition: macsio_iface.c:35
static void main_load(int argi, int argc, char **argv, char const *path, json_object *main_obj, json_object **data_read_obj)
Definition: macsio_silo.c:712
#define CHECK_SYMBOLN_STR(A)
Definition: macsio_silo.c:82
static int register_this_interface()
Definition: macsio_silo.c:871
void MACSIO_MIF_Finish(MACSIO_MIF_baton_t *bat)
End a MACSIO_MIF I/O operation and free resources.
Definition: macsio_mif.c:191
int MACSIO_LOG_DebugLevel
Definition: macsio_log.c:44
ProcessArgsFunc processArgsFunc
Definition: macsio_iface.h:59
#define MACSIO_IFACE_MAX_NAME
Definition: macsio_iface.h:33
#define CHECK_SYMBOLN_INT(A)
Definition: macsio_silo.c:70
static void CleanupDriverStuff()
Definition: macsio_silo.c:108
int json_object_array_length(struct json_object *obj)
Definition: json_object.c:854
static void CloseSiloFile(void *file, void *userData)
Definition: macsio_silo.c:344
#define JsonGetInt(OBJ,...)
Get the int value at specified path.
Definition: json_object.h:747
#define NARRVALS(Arr)
Definition: macsio_silo.c:62
double MACSIO_UTILS_YDelta(int const *dims, double const *bounds)
Definition: macsio_utils.c:236
static int process_args(int argi, int argc, char *argv[])
Definition: macsio_silo.c:270
static const int driver_nopts
Definition: macsio_silo.c:106
double MACSIO_UTILS_ZDelta(int const *dims, double const *bounds)
Definition: macsio_utils.c:241
static int driver_opts_ids[]
Definition: macsio_silo.c:101
static void DBSplitMultiName(char *mname, char **file, char **dir, char **obj)
Definition: macsio_silo.c:663
#define MACSIO_CLARGS_END_OF_ARGS
Definition: macsio_clargs.h:53
#define MACSIO_CLARGS_NODEFAULT
Definition: macsio_clargs.h:54
const char * json_object_get_string(struct json_object *obj)
Definition: json_object.c:751
#define MACSIO_CLARGS_ASSIGN_OFF
Definition: macsio_clargs.h:43
char ext[MACSIO_IFACE_MAX_NAME]
Definition: macsio_iface.h:55
#define MACSIO_CLARGS_TOMEM
Definition: macsio_clargs.h:39
struct json_object * json_object_array_get_idx(struct json_object *obj, int idx)
Definition: json_object.c:870
int MACSIO_DATA_GetRankOwningPart(json_object *main_obj, int chunkId)
Definition: macsio_data.c:946
static int dummy
Definition: macsio_silo.c:900
static void * CreateSiloFile(const char *fname, const char *nsname, void *userData)
Definition: macsio_silo.c:318
static int driver_ints[100]
Definition: macsio_silo.c:102
static int driver_nstrs
Definition: macsio_silo.c:105
MACSIO_MIF_baton_t * MACSIO_MIF_Init(int numFiles, MACSIO_MIF_ioFlags_t ioFlags, MPI_Comm mpiComm, int mpiTag, MACSIO_MIF_CreateCB createCb, MACSIO_MIF_OpenCB openCb, MACSIO_MIF_CloseCB closeCb, void *clientData)
Initialize MACSIO_MIF for a MIF I/O operation.
Definition: macsio_mif.c:102
#define CHECK_SYMBOLN_SYM(A)
Definition: macsio_silo.c:91
static const char * filename
Definition: macsio_silo.c:264
char name[MACSIO_IFACE_MAX_NAME]
Definition: macsio_iface.h:54
static int driver_nints
Definition: macsio_silo.c:103
int32_t json_object_get_int(struct json_object *obj)
Definition: json_object.c:504
#define CHECK_SYMBOL(A)
Definition: macsio_silo.c:68
enum json_extarr_type json_object_extarr_type(struct json_object *jso)
Definition: json_object.c:1143
static void WriteMultiXXXObjects(json_object *main_obj, DBfile *siloFile, int dumpn, MACSIO_MIF_baton_t *bat)
Definition: macsio_silo.c:403
static void main_dump(int argi, int argc, char **argv, json_object *main_obj, int dumpn, double dumpt)
Definition: macsio_silo.c:547
double MACSIO_UTILS_XDelta(int const *dims, double const *bounds)
Definition: macsio_utils.c:231
static void write_mesh_part(DBfile *dbfile, json_object *part)
Definition: macsio_silo.c:397
MPI_Comm MACSIO_MAIN_Comm
Definition: macsio_main.c:279
void MACSIO_MIF_HandOffBaton(MACSIO_MIF_baton_t const *Bat, void *file)
Release exclusive access to the group's file.
Definition: macsio_mif.c:284
#define MACSIO_LOG_MSG(SEV, MSG)
Convenience macro for logging a message to the main log.
Definition: macsio_log.h:133
static void WriteDecompMesh(json_object *main_obj, DBfile *siloFile, int dumpn, MACSIO_MIF_baton_t *bat)
Definition: macsio_silo.c:484
static void MakeDriverOpts(DBoptlist **_opts, int *opts_id)
Definition: macsio_silo.c:120
int MACSIO_DATA_SimpleAssignKPartsToNProcs(int k, int n, int my_rank, int *my_part_cnt, int **my_part_ids)
Definition: macsio_data.c:961
static char * driver_strs[]
Definition: macsio_silo.c:104
static void write_rect_mesh_part(DBfile *dbfile, json_object *part)
Definition: macsio_silo.c:351