]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
update step-18 to write .vtu files
authorheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 25 Dec 2013 18:09:10 +0000 (18:09 +0000)
committerheister <heister@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 25 Dec 2013 18:09:10 +0000 (18:09 +0000)
git-svn-id: https://svn.dealii.org/trunk@32122 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/examples/step-18/doc/intro.dox
deal.II/examples/step-18/doc/results.dox
deal.II/examples/step-18/step-18.cc

index 2d54bda4a979548e2cadc943cc74c7fc4992e4fa..5e494a505a9240bdfd5b59d5679bdb6ee7890e4d 100644 (file)
@@ -338,63 +338,23 @@ deferred to a later version.
 This functionality has been implemented in the meantime, and this is the time
 to explain its use. Basically, what we need to do is let every process
 generate graphical output for that subset of cells that it owns, write them
-into separate files and have a way to merge them later on. At this point, it
-should be noted that none of the graphical output formats known to the author
-of this program allows for a simple way to later re-read it and merge it with
-other files corresponding to the same simulation. What deal.II therefore
-offers is the following: When you call the <code>DataOut::build_patches</code>
-function, an intermediate format is generated that contains all the
-information for the data on each cell. Usually, this intermediate format is
-then further processed and converted into one of the graphical formats that we
-can presently write, such as gmv, eps, ucd, gnuplot, or a number of other
-ones. Once written in these formats, there is no way to reconstruct the
-necessary information to merge multiple blocks of output. However, the base
-classes of DataOut also allow to simply dump the intermediate format
-to a file, from which it can later be recovered without loss of information.
-
-This has two advantages: first, simulations may just dump the intermediate
-format data during run-time, and the user may later decide which particular
-graphics format she wants to have. This way, she does not have to re-run the
-entire simulation if graphical output is requested in a different format. One
-typical case is that one would like to take a quick look at the data with
-gnuplot, and then create high-quality pictures using GMV or OpenDX. Since both
-can be generated out of the intermediate format without problem, there is no
-need to re-run the simulation.
-
-In the present context, of more interest is the fact that in contrast to any
-of the other formats, it is simple to merge multiple files of intermediate
-format, if they belong to the same simulation. This is what we will do here:
-we will generate one output file in intermediate format for each processor
-that belongs to this computation (in the sequential case, this will simply be
-a single file). They may then later be read in and merged so that we can
-output a single file in whatever graphical format is requested.
-
-The way to do this is to first instruct the DataOutBase class to
-write intermediate format rather than in gmv or any other graphical
-format. This is simple: just use
-<code>data_out.write_deal_II_intermediate</code>. We will write to a file
-called <code>solution-TTTT.TTTT.d2</code> if there is only one processor, or
-files <code>solution-TTTT.TTTT.NNN.d2</code> if this is really a parallel
-job. Here, <code>TTTT.TTTT</code> denotes the time for which this output has
-been generated, and <code>NNN</code> the number of the MPI process that did this.
-
-The next step is to convert this file or these files into whatever
-format you like. The program that does this is the step-19 tutorial program:
-for example, for the first time step, call it through
-@code
-    ../\step-19/\step-19 solution-0001.0000.*.d2 solution-0001.0000.gmv
-@endcode
-to merge all the intermediate format files into a single file in GMV
-format. More details on the parameters of this program and what it can do for
-you can be found in the documentation of the step-19 tutorial program.
-
-@note In the years since the paragraphs above were written, it has
-also become possible to not only compute solutions in parallel, but
-also to visualize them with some programs (for example with the
-Paraview viewer). To make this efficient, one can not store the entire
-solution in a single file, but instead needs to have a single data
-file for each processor to visualize individually. We discuss this
-concept and how to use it in step-40.
+into separate files and have a way to display all files for a certain timestep
+at the same time. This way the code produces one .vtu file per processor per
+time step. The two common VTK file viewers ParaView and VisIt both support
+opening more than one .vtu file at once. To simplifiy the process of picking
+the correct files and allow moving around in time, both support record files
+that reference all files for a given timestep. Sadly, the record files have a
+different format between VisIt and Paraview, so we write out both formats.
+
+The code will generate the files solution-TTTT.NNN.vtu, where TTTT is the
+timestep number (starting from 1) and NNN is the processor id (starting from
+0). These files contain the locally owned cells for the timestep and
+processor. The files solution-TTTT.visit is the visit record for timestep
+TTTT, while solution-TTTT.pvtu is the same for ParaView. Finally, the file
+solution.pvd is a special record only supported by ParaView that references
+all time steps. So in ParaView, only solution.pvd needs to be opened, while
+one needs to select the group of all .visit files in VisIt for the same
+effect.
 
 
 <h3>Overall structure of the program</h3>
@@ -553,8 +513,7 @@ for (unsigned int i=0; i<dofs_per_cell; ++i)
 
 <li> <code>output_results ()</code>: This function simply outputs the solution
   based on what we have said above, i.e. every processor computes output only
-  for its own portion of the domain, and this can then be later merged by an
-  external program. In addition to the solution, we also compute the norm of
+  for its own portion of the domain. In addition to the solution, we also compute the norm of
   the stress averaged over all the quadrature points on each cell.
 </ul>
 
index b59460a793083c1defcc64028f2e7f19de2252cf..95e369d4c588bf54d30390dac1095785a581472d 100644 (file)
@@ -90,40 +90,20 @@ problem to keep a computer busy for a while. At the end of the day,
 this is what we have for output:
 @code
 examples/\step-18> ls -l *.d2
--rw-r--r--  1 bangerth wheeler 8797414 May 25 09:10 solution-0001.0000.d2
--rw-r--r--  1 bangerth wheeler 8788500 May 25 09:32 solution-0002.0000.d2
--rw-r--r--  1 bangerth wheeler 8763718 May 25 09:55 solution-0003.0000.d2
--rw-r--r--  1 bangerth wheeler 8738940 May 25 10:17 solution-0004.0000.d2
--rw-r--r--  1 bangerth wheeler 8710104 May 25 10:39 solution-0005.0000.d2
--rw-r--r--  1 bangerth wheeler 8685388 May 25 11:01 solution-0006.0000.d2
--rw-r--r--  1 bangerth wheeler 8649088 May 25 11:23 solution-0007.0000.d2
--rw-r--r--  1 bangerth wheeler 8585146 May 25 11:45 solution-0008.0000.d2
--rw-r--r--  1 bangerth wheeler 8489764 May 25 12:07 solution-0009.0000.d2
--rw-r--r--  1 bangerth wheeler 8405388 May 25 12:29 solution-0010.0000.d2
+-rw-r--r--  1 bangerth wheeler 8797414 May 25 09:10 solution-0001.0000.vtu
+-rw-r--r--  1 bangerth wheeler 8788500 May 25 09:32 solution-0002.0000.vtu
+-rw-r--r--  1 bangerth wheeler 8763718 May 25 09:55 solution-0003.0000.vtu
+-rw-r--r--  1 bangerth wheeler 8738940 May 25 10:17 solution-0004.0000.vtu
+-rw-r--r--  1 bangerth wheeler 8710104 May 25 10:39 solution-0005.0000.vtu
+-rw-r--r--  1 bangerth wheeler 8685388 May 25 11:01 solution-0006.0000.vtu
+-rw-r--r--  1 bangerth wheeler 8649088 May 25 11:23 solution-0007.0000.vtu
+-rw-r--r--  1 bangerth wheeler 8585146 May 25 11:45 solution-0008.0000.vtu
+-rw-r--r--  1 bangerth wheeler 8489764 May 25 12:07 solution-0009.0000.vtu
+-rw-r--r--  1 bangerth wheeler 8405388 May 25 12:29 solution-0010.0000.vtu
 @endcode
 
 
-Let us convert these files in deal.II intermediate format to gmv
-format (this assumes that you have already compiled the
-step-19 example program):
-@code
-examples/\step-18> ../\step-19/\step-19
-
-Converter from deal.II intermediate format to other graphics formats.
-
-Usage: ./\step-19 [-p parameter_file] list_of_input_files [-x output_format] -o output_file
-
-examples/\step-18> ../\step-19/\step-19 solution-0001.0000.d2 -x gmv -o solution-0001.0000.gmv
-examples/\step-18> ../\step-19/\step-19 solution-0002.0000.d2 -x gmv -o solution-0002.0000.gmv
-[...]
-@endcode
-Of course, since we have run the program only in sequential mode, we
-do have only one intermediate file for each time step that we have to
-take as input.
-
-
-
-If we visualize these files with GMV, we get to see the full picture
+If we visualize these files with VisIt or Paraview, we get to see the full picture
 of the disaster our forced compression wreaks on the cylinder (colors in the
 images encode the norm of the stress in the material):
 
@@ -233,39 +213,26 @@ Timestep 20 at time 10
 That's quite a good number of unknowns, given that we are in 3d. The output of
 this program are 16 files for each time step:
 @code
-examples/\step-18> ls -l solution-0001.000*
--rw-r--r--    1 bangerth mfw       4325219 Aug 11 09:44 solution-0001.0000-000.d2
--rw-r--r--    1 bangerth mfw       4454460 Aug 11 09:44 solution-0001.0000-001.d2
--rw-r--r--    1 bangerth mfw       4485242 Aug 11 09:43 solution-0001.0000-002.d2
--rw-r--r--    1 bangerth mfw       4517364 Aug 11 09:43 solution-0001.0000-003.d2
--rw-r--r--    1 bangerth mfw       4462829 Aug 11 09:43 solution-0001.0000-004.d2
--rw-r--r--    1 bangerth mfw       4482487 Aug 11 09:43 solution-0001.0000-005.d2
--rw-r--r--    1 bangerth mfw       4548619 Aug 11 09:43 solution-0001.0000-006.d2
--rw-r--r--    1 bangerth mfw       4522421 Aug 11 09:43 solution-0001.0000-007.d2
--rw-r--r--    1 bangerth mfw       4337529 Aug 11 09:43 solution-0001.0000-008.d2
--rw-r--r--    1 bangerth mfw       4163047 Aug 11 09:43 solution-0001.0000-009.d2
--rw-r--r--    1 bangerth mfw       4288247 Aug 11 09:43 solution-0001.0000-010.d2
--rw-r--r--    1 bangerth mfw       4350410 Aug 11 09:43 solution-0001.0000-011.d2
--rw-r--r--    1 bangerth mfw       4458427 Aug 11 09:43 solution-0001.0000-012.d2
--rw-r--r--    1 bangerth mfw       4466037 Aug 11 09:43 solution-0001.0000-013.d2
--rw-r--r--    1 bangerth mfw       4505679 Aug 11 09:44 solution-0001.0000-014.d2
--rw-r--r--    1 bangerth mfw       4340488 Aug 11 09:44 solution-0001.0000-015.d2
-@endcode
-We merge and convert these 16 intermediate files into a single gmv file as
-follows:
-@code
-examples/\step-18> time ../\step-19/\step-19 solution-0001.0000-* -x gmv -o solution-0001.0000.gmv
-
-real    0m45.929s
-user    0m41.290s
-sys     0m0.990s
-examples/\step-18> ls -l solution-0001.0000.gmv
--rw-r--r--    1 bangerth mfw      68925360 Aug 11 17:04 solution-0001.0000.gmv
+examples/\step-18> ls -l solution-0001*
+-rw-r--r--    1 bangerth mfw       4325219 Aug 11 09:44 solution-0001.000.d2
+-rw-r--r--    1 bangerth mfw       4454460 Aug 11 09:44 solution-0001.001.d2
+-rw-r--r--    1 bangerth mfw       4485242 Aug 11 09:43 solution-0001.002.d2
+-rw-r--r--    1 bangerth mfw       4517364 Aug 11 09:43 solution-0001.003.d2
+-rw-r--r--    1 bangerth mfw       4462829 Aug 11 09:43 solution-0001.004.d2
+-rw-r--r--    1 bangerth mfw       4482487 Aug 11 09:43 solution-0001.005.d2
+-rw-r--r--    1 bangerth mfw       4548619 Aug 11 09:43 solution-0001.006.d2
+-rw-r--r--    1 bangerth mfw       4522421 Aug 11 09:43 solution-0001.007.d2
+-rw-r--r--    1 bangerth mfw       4337529 Aug 11 09:43 solution-0001.008.d2
+-rw-r--r--    1 bangerth mfw       4163047 Aug 11 09:43 solution-0001.009.d2
+-rw-r--r--    1 bangerth mfw       4288247 Aug 11 09:43 solution-0001.010.d2
+-rw-r--r--    1 bangerth mfw       4350410 Aug 11 09:43 solution-0001.011.d2
+-rw-r--r--    1 bangerth mfw       4458427 Aug 11 09:43 solution-0001.012.d2
+-rw-r--r--    1 bangerth mfw       4466037 Aug 11 09:43 solution-0001.013.d2
+-rw-r--r--    1 bangerth mfw       4505679 Aug 11 09:44 solution-0001.014.d2
+-rw-r--r--    1 bangerth mfw       4340488 Aug 11 09:44 solution-0001.015.d2
 @endcode
 
-Doing so for all time steps, we obtain gmv files that we can visualize (albeit
-with some difficulty, due to their size gmv isn't exactly fast when plotting
-them). Here are first the mesh on which we compute as well as the partitioning
+Here are first the mesh on which we compute as well as the partitioning
 for the 16 processors:
 
 
index 1d6347f8ab6bf630694946c8a37c294f91658dfe..ff5f5bb61928b9f4efe1316366305af6bd591d0e 100644 (file)
@@ -15,7 +15,8 @@
  * ---------------------------------------------------------------------
 
  *
- * Author: Wolfgang Bangerth, University of Texas at Austin, 2000, 2004, 2005
+ * Author: Wolfgang Bangerth, University of Texas at Austin, 2000, 2004, 2005,
+ * Timo Heister, 2013
  */
 
 
@@ -1284,11 +1285,10 @@ namespace Step18
 
   // @sect4{TopLevel::output_results}
 
-  // This function generates the graphical output in intermediate format as
-  // explained in the introduction. Each process will only work on the cells
-  // it owns, and then write the result into a file of its own. These files
-  // may later be merged to get a single file in any of the supported output
-  // files, as mentioned in the introduction.
+  // This function generates the graphical output in .vtu format as explained
+  // in the introduction. Each process will only work on the cells it owns,
+  // and then write the result into a file of its own. Additionally, processor
+  // 0 will write the record files the reference all the .vtu files.
   //
   // The crucial part of this function is to give the <code>DataOut</code>
   // class a way to only work on the cells that the present process owns. This
@@ -1468,50 +1468,70 @@ namespace Step18
     data_out.build_patches ();
 
 
-    // Now that we have generated the intermediate format, let us determine
-    // the name of the file we will want to write it to. We compose it of the
-    // prefix <code>solution-</code>, followed by a representation of the
-    // present time written as a fixed point number so that file names sort
-    // naturally:
-    std::ostringstream filename;
-    filename << "solution-";
-    filename << std::setfill('0');
-    filename.setf(std::ios::fixed, std::ios::floatfield);
-    filename << std::setw(9) << std::setprecision(4) << present_time;
-
-    // Next, in case there are multiple processes working together, we have to
-    // generate different file names for the output of each process. In our
-    // case, we encode the process number as a three-digit integer, padded
-    // with zeros. The assertion in the first line of the block makes sure
-    // that there are less than 1000 processes (a very conservative check, but
-    // worth having anyway) as our scheme of generating process numbers would
-    // overflow if there were 1000 processes or more. Note that we choose to
-    // use <code>AssertThrow</code> rather than <code>Assert</code> since the
-    // number of processes is a variable that depends on input files or the
-    // way the process is started, rather than static assumptions in the
-    // program code. Therefore, it is inappropriate to use <code>Assert</code>
-    // that is optimized away in optimized mode, whereas here we actually can
-    // assume that users will run the largest computations with the most
-    // processors in optimized mode, and we should check our assumptions in
-    // this particular case, and not only when running in debug mode:
-    if (n_mpi_processes != 1)
-      {
-        AssertThrow (n_mpi_processes < 1000, ExcNotImplemented());
-
-        filename << '-';
-        filename << std::setfill('0');
-        filename << std::setw(3) << this_mpi_process;
-      }
-
-    // To the file name, attach the file name suffix usually used for the
-    // deal.II intermediate format. To determine it, we use the same function
-    // that has already been used in step-13:
-    filename << data_out.default_suffix(DataOut<dim>::deal_II_intermediate);
+    // Let us determine the name of the file we will want to write it to. We
+    // compose it of the prefix <code>solution-</code>, followed by the time
+    // step number, and finally the processor id (encoded as a three digit
+    // number):
+    std::string filename = "solution-" + Utilities::int_to_string(timestep_no,4)
+                          + "." + Utilities::int_to_string(this_mpi_process,3)
+                          + ".vtu";
+    
+    // The following assertion makes sure that there are less than 1000
+    // processes (a very conservative check, but worth having anyway) as our
+    // scheme of generating process numbers would overflow if there were 1000
+    // processes or more. Note that we choose to use <code>AssertThrow</code>
+    // rather than <code>Assert</code> since the number of processes is a
+    // variable that depends on input files or the way the process is started,
+    // rather than static assumptions in the program code. Therefore, it is
+    // inappropriate to use <code>Assert</code> that is optimized away in
+    // optimized mode, whereas here we actually can assume that users will run
+    // the largest computations with the most processors in optimized mode,
+    // and we should check our assumptions in this particular case, and not
+    // only when running in debug mode:
+    AssertThrow (n_mpi_processes < 1000, ExcNotImplemented());    
 
     // With the so-completed filename, let us open a file and write the data
-    // we have generated into it, using the intermediate format:
-    std::ofstream output (filename.str().c_str());
-    data_out.write_deal_II_intermediate (output);
+    // we have generated into it:
+    std::ofstream output (filename.c_str());
+    data_out.write_vtu (output);
+
+    // The record files must be written only once and not by each processor,
+    // so we do this on processor 0:
+    if (this_mpi_process==0)
+      {
+       // Here we collect all filenames of the current timestep (same format as above)
+       std::vector<std::string> filenames;
+        for (unsigned int i=0; i<n_mpi_processes; ++i)
+          filenames.push_back ("solution-" + Utilities::int_to_string(timestep_no,4)
+                              + "." + Utilities::int_to_string(i,3)
+                              + ".vtu");
+
+       // Now we write the .visit file. The naming is similar to the .vtu files, only
+       // that the file obviously doesn't contain a processor id.
+        const std::string
+        visit_master_filename = ("solution-" +
+                                Utilities::int_to_string(timestep_no,4) +
+                                 ".visit");
+        std::ofstream visit_master (visit_master_filename.c_str());
+        data_out.write_visit_record (visit_master, filenames);
+
+       // Similarly, we write the paraview .pvtu:
+        const std::string
+        pvtu_master_filename = ("solution-" +
+                               Utilities::int_to_string(timestep_no,4) +
+                                ".pvtu");
+        std::ofstream pvtu_master (pvtu_master_filename.c_str());
+        data_out.write_pvtu_record (pvtu_master, filenames);
+
+       // Finally, we write the paraview record, that references all .pvtu files and
+       // their respective time. Note that the variable times_and_names is declared
+       // static, so it will retain the entries from the pervious timesteps.
+       static std::vector<std::pair<double,std::string> > times_and_names;
+       times_and_names.push_back (std::pair<double,std::string> (present_time, pvtu_master_filename));
+       std::ofstream pvd_output ("solution.pvd");
+       data_out.write_pvd_record (pvd_output, times_and_names);
+      }
+    
   }
 
 

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.