From ded5e9f49130b45642a660058e7693bc9d4703d5 Mon Sep 17 00:00:00 2001 From: bangerth Date: Thu, 13 May 2010 03:01:41 +0000 Subject: [PATCH] Patch by Scott Miller: Implement writing in VTU format. git-svn-id: https://svn.dealii.org/trunk@21120 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/base/include/base/data_out_base.h | 48 ++- deal.II/base/source/data_out_base.cc | 380 ++++++++++++++++++++++ 2 files changed, 419 insertions(+), 9 deletions(-) diff --git a/deal.II/base/include/base/data_out_base.h b/deal.II/base/include/base/data_out_base.h index 969300459f..7f4d162b0b 100644 --- a/deal.II/base/include/base/data_out_base.h +++ b/deal.II/base/include/base/data_out_base.h @@ -1246,9 +1246,6 @@ class DataOutBase */ unsigned int memory_consumption () const; }; - - - /** * Flags controlling the details * of output in deal.II @@ -1395,6 +1392,12 @@ class DataOutBase */ vtk, + /** + * Output in @ref + * SoftwareVTU format. + */ + vtu, + /** * Output in deal.II * intermediate format. @@ -1674,11 +1677,9 @@ class DataOutBase std::ostream &out); /** - * Write the given list of patches to the output stream in @ref - * SoftwareVTK format. - * - * This is the file format used by @ref SoftwareVTK, as described in - * their manual, section 14.3. + * Write the given list of patches to the output stream in @ref SoftwareVTK + * format. The data is written in the traditional VTK format as opposed to the + * XML-based format that write_vtu() produces. * * The vector_data_ranges argument denotes ranges of components in the * output that are considered a vector, rather than simply a @@ -1693,6 +1694,26 @@ class DataOutBase const std::vector > &vector_data_ranges, const VtkFlags &flags, std::ostream &out); + + +/** + * Write the given list of patches to the output stream in @ref SoftwareVTK + * format. The data is written in the XML-based VTK format as opposed to the + * traditional format that write_vtk() produces. + * + * The vector_data_ranges argument denotes ranges of components in the + * output that are considered a vector, rather than simply a + * collection of scalar fields. The VTK output format has special + * provisions that allow these components to be output by a single + * name rather than having to group several scalar fields into a + * vector later on in the visualization program. + */ + template + static void write_vtu (const std::vector > &patches, + const std::vector &data_names, + const std::vector > &vector_data_ranges, + const VtkFlags &flags, + std::ostream &out); /** * Write the given list of patches to the output stream in deal.II @@ -1816,7 +1837,7 @@ class DataOutBase *
  • povray: .pov *
  • eps: .eps *
  • gmv: .gmv - *
  • vtk: .vtk. + *
  • vtk: .vtk *
  • deal_II_intermediate: .d2. * */ @@ -2105,6 +2126,7 @@ class DataOutInterface : private DataOutBase using DataOutBase::tecplot; using DataOutBase::tecplot_binary; using DataOutBase::vtk; + using DataOutBase::vtu; using DataOutBase::deal_II_intermediate; using DataOutBase::parse_output_format; using DataOutBase::get_output_format_names; @@ -2197,6 +2219,14 @@ class DataOutInterface : private DataOutBase * DataOut::write_vtk. */ void write_vtk (std::ostream &out) const; + + /** + * Obtain data through get_patches() + * and write it to out + * in Vtu (VTK's XML) format. See + * DataOut::write_vtu. + */ + void write_vtu (std::ostream &out) const; /** * Obtain data through get_patches() diff --git a/deal.II/base/source/data_out_base.cc b/deal.II/base/source/data_out_base.cc index 8e62c1eb2b..5ceaa164e9 100644 --- a/deal.II/base/source/data_out_base.cc +++ b/deal.II/base/source/data_out_base.cc @@ -589,6 +589,69 @@ namespace }; + class VtuStream + { + public: + /** + * Constructor, storing + * persistent values for + * later use. + */ + VtuStream (std::ostream& stream, + const DataOutBase::VtkFlags flags); + + /** + * Output operator for points. + */ + template + void write_point (const unsigned int index, + const Point&); + + /** + * Write dim-dimensional cell + * with first vertex at + * number start and further + * vertices offset by the + * specified values. Values + * not needed are ignored. + * + * The order of vertices for + * these cells in different + * dimensions is + *
      + *
    1. [0,1] + *
    2. [] + *
    3. [] + *
    + */ + template + void write_cell(const unsigned int index, + const unsigned int start, + const unsigned int x_offset, + const unsigned int y_offset, + const unsigned int z_offset); + + /** + * Forwarding of output stream + */ + template + std::ostream& operator<< (const T&); + + private: + /** + * The ostream to use. Since + * the life span of these + * objects is small, we use a + * very simple storage + * technique. + */ + std::ostream& stream; + + /** + * The flags controlling the output + */ + const DataOutBase::VtkFlags flags; + }; //----------------------------------------------------------------------// @@ -913,8 +976,52 @@ namespace stream << '\n'; } + VtuStream::VtuStream(std::ostream& out, const DataOutBase::VtkFlags f) + : + stream(out), flags(f) + {} + template + void + VtuStream::write_point (const unsigned int, + const Point& p) + { + // write out coordinates + stream << p; + // fill with zeroes + for (unsigned int i=dim; i<3; ++i) + stream << " 0"; + stream << '\n'; + } + + + template + void + VtuStream::write_cell( + unsigned int, + unsigned int start, + unsigned int d1, + unsigned int d2, + unsigned int d3) + { + stream << start << '\t' + << start+d1; + if (dim>=2) + { + stream << '\t' << start+d2+d1 + << '\t' << start+d2; + if (dim>=3) + { + stream << '\t' << start+d3 + << '\t' << start+d3+d1 + << '\t' << start+d3+d2+d1 + << '\t' << start+d3+d2; + } + } + stream << '\n'; + } + template std::ostream& @@ -3997,6 +4104,268 @@ DataOutBase::write_vtk (const std::vector > &patches, } +template +void +DataOutBase::write_vtu (const std::vector > &patches, + const std::vector &data_names, + const std::vector > &vector_data_ranges, + const VtkFlags &flags, + std::ostream &out) +{ + AssertThrow (out, ExcIO()); + + Assert (patches.size() > 0, ExcNoPatches()); + + VtuStream vtu_out(out, flags); + + const unsigned int n_data_sets = data_names.size(); + // check against # of data sets in + // first patch. checks against all + // other patches are made in + // write_gmv_reorder_data_vectors + Assert ((patches[0].data.n_rows() == n_data_sets && !patches[0].points_are_available) || + (patches[0].data.n_rows() == n_data_sets+spacedim && patches[0].points_are_available), + ExcDimensionMismatch (patches[0].points_are_available + ? + (n_data_sets + spacedim) + : + n_data_sets, + patches[0].data.n_rows())); + + /////////////////////// + // preamble + if (true) + { + std::time_t time1= std::time (0); + std::tm *time = std::localtime(&time1); + out << " \n"; + out << " \n"; + // LittleEndian vs BigEndian should not + // matter; we are writing an ascii + // file. + out << ""; + out << '\n'; + out << ""; + out << '\n'; + } + + + // first count the number of cells + // and cells for later use + unsigned int n_nodes; + unsigned int n_cells; + compute_sizes(patches, n_nodes, n_cells); + // in gmv format the vertex + // coordinates and the data have an + // order that is a bit unpleasant + // (first all x coordinates, then + // all y coordinate, ...; first all + // data of variable 1, then + // variable 2, etc), so we have to + // copy the data vectors a bit around + // + // note that we copy vectors when + // looping over the patches since we + // have to write them one variable + // at a time and don't want to use + // more than one loop + // + // this copying of data vectors can + // be done while we already output + // the vertices, so do this on a + // separate task and when wanting + // to write out the data, we wait + // for that task to finish + Table<2,double> data_vectors (n_data_sets, n_nodes); + + void (*fun_ptr) (const std::vector > &, + Table<2,double> &) + = &DataOutBase::template write_gmv_reorder_data_vectors; + Threads::Task<> reorder_task = Threads::new_task (fun_ptr, patches, data_vectors); + + /////////////////////////////// + // first make up a list of used + // vertices along with their + // coordinates + // + // note that we have to print + // d=1..3 dimensions + out << " \n"; + out << " \n"; + out << " \n"; + write_nodes(patches, vtu_out); + out << " \n"; + out << " \n \n"; + ///////////////////////////////// + // now for the cells + out << " \n"; + out << " \n"; + write_cells(patches, vtu_out); + out << " \n"; + + // XML VTU format uses offsets; this is + // different than the VTK format, which + // puts the number of nodes per cell in + // front of the connectivity list. + out << " \n"; + for(unsigned int i=0; i::vertices_per_cell); + + out << " \n"; + + // next output the types of the + // cells. since all cells are + // the same, this is simple + out << " \n"; + + for (unsigned int i=0; i \n"; + out << " \n"; + + + /////////////////////////////////////// + // data output. + + // now write the data vectors to + // @p{out} first make sure that all + // data is in place + reorder_task.join (); + + // then write data. the + // 'POINT_DATA' means: node data + // (as opposed to cell data, which + // we do not support explicitly + // here). all following data sets + // are point data + out << " \n"; + + // when writing, first write out + // all vector data, then handle the + // scalar data sets that have been + // left over + std::vector data_set_written (n_data_sets, false); + for (unsigned int n_th_vector=0; n_th_vector(vector_data_ranges[n_th_vector]) >= + std_cxx1x::get<0>(vector_data_ranges[n_th_vector]), + ExcLowerRange (std_cxx1x::get<1>(vector_data_ranges[n_th_vector]), + std_cxx1x::get<0>(vector_data_ranges[n_th_vector]))); + AssertThrow (std_cxx1x::get<1>(vector_data_ranges[n_th_vector]) < n_data_sets, + ExcIndexRange (std_cxx1x::get<1>(vector_data_ranges[n_th_vector]), + 0, n_data_sets)); + AssertThrow (std_cxx1x::get<1>(vector_data_ranges[n_th_vector]) + 1 + - std_cxx1x::get<0>(vector_data_ranges[n_th_vector]) <= 3, + ExcMessage ("Can't declare a vector with more than 3 components " + "in VTK")); + + // mark these components as already + // written: + for (unsigned int i=std_cxx1x::get<0>(vector_data_ranges[n_th_vector]); + i<=std_cxx1x::get<1>(vector_data_ranges[n_th_vector]); + ++i) + data_set_written[i] = true; + + // write the + // header. concatenate all the + // component names with double + // underscores unless a vector + // name has been specified + out << " (vector_data_ranges[n_th_vector]) != "") + out << std_cxx1x::get<2>(vector_data_ranges[n_th_vector]); + else + { + for (unsigned int i=std_cxx1x::get<0>(vector_data_ranges[n_th_vector]); + i(vector_data_ranges[n_th_vector]); + ++i) + out << data_names[i] << "__"; + out << data_names[std_cxx1x::get<1>(vector_data_ranges[n_th_vector])]; + } + + out << "\" NumberOfComponents=\"3\" format=\"ascii\"> \n"; + + // now write data. pad all + // vectors to have three + // components + for (unsigned int n=0; n(vector_data_ranges[n_th_vector]) - + std_cxx1x::get<0>(vector_data_ranges[n_th_vector])) + { + case 0: + out << data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector]), n) << " 0 0" + << '\n'; + break; + + case 1: + out << data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector]), n) << ' ' + << data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector])+1, n) << " 0" + << '\n'; + break; + case 2: + out << data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector]), n) << ' ' + << data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector])+1, n) << ' ' + << data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector])+2, n) + << '\n'; + break; + + default: + // VTK doesn't + // support + // anything else + // than vectors + // with 1, 2, or + // 3 components + Assert (false, ExcInternalError()); + } + } + + out << " \n"; + } + + // now do the left over scalar data sets + for (unsigned int data_set=0; data_set \n"; + + std::copy (data_vectors[data_set].begin(), + data_vectors[data_set].end(), + std::ostream_iterator(out, " ")); + out << '\n'; + out << " \n"; + } + + out << " \n"; + + // Finish up writing a valid XML file + out << " \n"; + out << " \n"; + out << " \n"; + // make sure everything now gets to + // disk + out.flush (); + + // assert the stream is still ok + AssertThrow (out, ExcIO()); +} + + template void @@ -4213,6 +4582,13 @@ void DataOutInterface::write_vtk (std::ostream &out) const vtk_flags, out); } +template +void DataOutInterface::write_vtu (std::ostream &out) const +{ + DataOutBase::write_vtu (get_patches(), get_dataset_names(), + get_vector_data_ranges(), + vtk_flags, out); +} template @@ -4275,6 +4651,10 @@ DataOutInterface::write (std::ostream &out, case vtk: write_vtk (out); break; + + case vtu: + write_vtu (out); + break; case deal_II_intermediate: write_deal_II_intermediate (out); -- 2.39.5