From 14f0c1b8e032414e7f89d01d550b9478f66ef27a Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Fri, 28 Jan 2011 18:48:06 +0000 Subject: [PATCH] Write data in VTU format compressed, if possible. git-svn-id: https://svn.dealii.org/trunk@23270 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/doc/news/changes.h | 9 +- deal.II/source/base/data_out_base.cc | 789 ++++++++++++++++----------- 2 files changed, 493 insertions(+), 305 deletions(-) diff --git a/deal.II/doc/news/changes.h b/deal.II/doc/news/changes.h index 4eaf74addf..f787685f2f 100644 --- a/deal.II/doc/news/changes.h +++ b/deal.II/doc/news/changes.h @@ -40,7 +40,7 @@ on any system.
(Wolfgang Bangerth, 2011/01/22) -
  • Fixed: When using the --enable-mpi to +
  • Fixed: When using the --enable-mpi to ./configure, the script only tried mpiCC as the MPI C++ compiler. However, on some systems, it is called mpicxx. The script now tries that as well. @@ -71,6 +71,13 @@ should be fixed now.

    Specific improvements

      +
    1. Changed: If the libz library was detected during library +configuration, the function DataOutBase::write_vtu now writes data in compressed +format, saving a good fraction of disk space in output files. +
      +(Wolfgang Bangerth, 2011/01/28) +
    2. +
    3. New: Trilinos and PETSc vectors now have a function has_ghost_elements().
      (Timo Heister, 2011/01/26) diff --git a/deal.II/source/base/data_out_base.cc b/deal.II/source/base/data_out_base.cc index b6f36224f3..1a196e4d4a 100644 --- a/deal.II/source/base/data_out_base.cc +++ b/deal.II/source/base/data_out_base.cc @@ -48,6 +48,7 @@ # include #endif + DEAL_II_NAMESPACE_OPEN @@ -64,6 +65,203 @@ namespace } +namespace +{ + // the functions in this namespace are + // taken from the libb64 project, see + // http://sourceforge.net/projects/libb64 + // + // libb64 has been placed in the public + // domain + namespace base64 + { + typedef enum + { + step_A, step_B, step_C + } base64_encodestep; + + typedef struct + { + base64_encodestep step; + char result; + } base64_encodestate; + + void base64_init_encodestate(base64_encodestate* state_in) + { + state_in->step = step_A; + state_in->result = 0; + } + + inline + char base64_encode_value(char value_in) + { + static const char* encoding + = "ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz0123456789+/"; + if (value_in > 63) return '='; + return encoding[(int)value_in]; + } + + int base64_encode_block(const char* plaintext_in, + int length_in, + char* code_out, + base64_encodestate *state_in) + { + const char* plainchar = plaintext_in; + const char* const plaintextend = plaintext_in + length_in; + char* codechar = code_out; + char result; + char fragment; + + result = state_in->result; + + switch (state_in->step) + { + while (1) + { + case step_A: + if (plainchar == plaintextend) + { + state_in->result = result; + state_in->step = step_A; + return codechar - code_out; + } + fragment = *plainchar++; + result = (fragment & 0x0fc) >> 2; + *codechar++ = base64_encode_value(result); + result = (fragment & 0x003) << 4; + case step_B: + if (plainchar == plaintextend) + { + state_in->result = result; + state_in->step = step_B; + return codechar - code_out; + } + fragment = *plainchar++; + result |= (fragment & 0x0f0) >> 4; + *codechar++ = base64_encode_value(result); + result = (fragment & 0x00f) << 2; + case step_C: + if (plainchar == plaintextend) + { + state_in->result = result; + state_in->step = step_C; + return codechar - code_out; + } + fragment = *plainchar++; + result |= (fragment & 0x0c0) >> 6; + *codechar++ = base64_encode_value(result); + result = (fragment & 0x03f) >> 0; + *codechar++ = base64_encode_value(result); + } + } + /* control should not reach here */ + return codechar - code_out; + } + + int base64_encode_blockend(char* code_out, base64_encodestate* state_in) + { + char* codechar = code_out; + + switch (state_in->step) + { + case step_B: + *codechar++ = base64_encode_value(state_in->result); + *codechar++ = '='; + *codechar++ = '='; + break; + case step_C: + *codechar++ = base64_encode_value(state_in->result); + *codechar++ = '='; + break; + case step_A: + break; + } + *codechar++ = '\0'; + + return codechar - code_out; + } + } + + + /** + * Do a base64 encoding of the given data. + * + * The function allocates memory as + * necessary and returns a pointer to + * it. The calling function must release + * this memory again. + */ + char * + encode_block (const char *data, + const int data_size) + { + base64::base64_encodestate state; + base64::base64_init_encodestate(&state); + + char *encoded_data = new char[2*data_size+1]; + + const int encoded_length_data + = base64::base64_encode_block (data, data_size, + encoded_data, &state); + base64::base64_encode_blockend (encoded_data + encoded_length_data, + &state); + + return encoded_data; + } + + + +#ifdef HAVE_LIBZ + /** + * Do a zlib compression followed + * by a base64 encoding of the + * given data. The result is then + * written to the given stream. + */ + template + void write_compressed_block (const std::vector &data, + std::ostream &output_stream) + { + if (data.size() != 0) + { + // allocate a buffer for compressing + // data and do so + uLongf compressed_data_length + = compressBound (data.size() * sizeof(T)); + char *compressed_data = new char[compressed_data_length]; + int err = compress2 ((Bytef *) compressed_data, + &compressed_data_length, + (const Bytef *) &data[0], + data.size() * sizeof(T), + Z_BEST_COMPRESSION); + Assert (err == Z_OK, ExcInternalError()); + + // now encode the compression header + const uint32_t compression_header[4] + = { 1, /* number of blocks */ + 32768, /* size of block */ + (uint32_t)(data.size() * sizeof(T)), /* size of last block */ + (uint32_t)compressed_data_length }; /* list of compressed sizes of blocks */ + + char *encoded_header = encode_block ((char*)&compression_header[0], + 4 * sizeof(compression_header[0])); + output_stream << encoded_header; + delete[] encoded_header; + + // next do the compressed + // data encoding in base64 + char *encoded_data = encode_block (compressed_data, + compressed_data_length); + delete[] compressed_data; + + output_stream << encoded_data; + delete[] encoded_data; + } + } +#endif +} + + //----------------------------------------------------------------------// //Auxiliary data @@ -208,6 +406,12 @@ namespace void write_point (const unsigned int index, const Point&); + /** + * Do whatever is necessary to + * terminate the list of points. + */ + void flush_points (); + /** * Write dim-dimensional cell * with first vertex at @@ -232,6 +436,12 @@ namespace const unsigned int y_offset, const unsigned int z_offset); + /** + * Do whatever is necessary to + * terminate the list of cells. + */ + void flush_cells (); + /** * Write a complete set of * data for a single node. @@ -292,6 +502,12 @@ namespace void write_point (const unsigned int index, const Point&); + /** + * Do whatever is necessary to + * terminate the list of points. + */ + void flush_points (); + /** * Write dim-dimensional cell * with first vertex at @@ -316,6 +532,12 @@ namespace const unsigned int y_offset, const unsigned int z_offset); + /** + * Do whatever is necessary to + * terminate the list of cells. + */ + void flush_cells (); + /** * Forwarding of output stream */ @@ -375,6 +597,12 @@ namespace void write_point (const unsigned int index, const Point&); + /** + * Do whatever is necessary to + * terminate the list of points. + */ + void flush_points (); + /** * Write dim-dimensional cell * with first vertex at @@ -399,6 +627,12 @@ namespace const unsigned int y_offset, const unsigned int z_offset); + /** + * Do whatever is necessary to + * terminate the list of cells. + */ + void flush_cells (); + /** * Forwarding of output stream */ @@ -459,6 +693,12 @@ namespace void write_point (const unsigned int index, const Point&); + /** + * Do whatever is necessary to + * terminate the list of points. + */ + void flush_points (); + /** * Write dim-dimensional cell * with first vertex at @@ -487,6 +727,12 @@ namespace const unsigned int y_offset, const unsigned int z_offset); + /** + * Do whatever is necessary to + * terminate the list of cells. + */ + void flush_cells (); + /** * Write a complete set of * data for a single node. @@ -546,6 +792,12 @@ namespace void write_point (const unsigned int index, const Point&); + /** + * Do whatever is necessary to + * terminate the list of points. + */ + void flush_points (); + /** * Write dim-dimensional cell * with first vertex at @@ -570,6 +822,12 @@ namespace const unsigned int y_offset, const unsigned int z_offset); + /** + * Do whatever is necessary to + * terminate the list of cells. + */ + void flush_cells (); + /** * Forwarding of output stream */ @@ -611,6 +869,12 @@ namespace void write_point (const unsigned int index, const Point&); + /** + * Do whatever is necessary to + * terminate the list of points. + */ + void flush_points (); + /** * Write dim-dimensional cell * with first vertex at @@ -635,12 +899,32 @@ namespace const unsigned int y_offset, const unsigned int z_offset); + /** + * Do whatever is necessary to + * terminate the list of cells. + */ + void flush_cells (); + /** * Forwarding of output stream */ template std::ostream& operator<< (const T&); + /** + * Forwarding of output stream. + * + * If libz was found during + * configuration, this operator + * compresses and encodes the + * entire data + * block. Otherwise, it simply + * writes it element by + * element. + */ + template + std::ostream& operator<< (const std::vector&); + private: /** * The ostream to use. Since @@ -655,6 +939,21 @@ namespace * The flags controlling the output */ const DataOutBase::VtkFlags flags; + + /** + * A list of vertices and + * cells, to be used in case we + * want to compress the data. + * + * The data types of these + * arrays needs to match what + * we print in the XML-preamble + * to the respective parts of + * VTU files (e.g. Float64 and + * Int32) + */ + std::vector vertices; + std::vector cells; }; @@ -689,6 +988,11 @@ namespace } + void + DXStream::flush_points () + {} + + template void DXStream::write_cell( @@ -729,6 +1033,11 @@ namespace } + void + DXStream::flush_cells () + {} + + template inline void @@ -771,6 +1080,11 @@ namespace } + void + GmvStream::flush_points () + {} + + template void GmvStream::write_cell( @@ -804,6 +1118,11 @@ namespace + void + GmvStream::flush_cells () + {} + + //----------------------------------------------------------------------// TecplotStream::TecplotStream(std::ostream& out, const DataOutBase::TecplotFlags f) @@ -824,6 +1143,11 @@ namespace } + void + TecplotStream::flush_points () + {} + + template void TecplotStream::write_cell( @@ -854,6 +1178,12 @@ namespace + void + TecplotStream::flush_cells () + {} + + + //----------------------------------------------------------------------// UcdStream::UcdStream(std::ostream& out, const DataOutBase::UcdFlags f) @@ -879,6 +1209,12 @@ namespace } + + void + UcdStream::flush_points () + {} + + template void UcdStream::write_cell( @@ -917,6 +1253,12 @@ namespace } + + void + UcdStream::flush_cells () + {} + + template inline void @@ -953,6 +1295,12 @@ namespace } + + void + VtkStream::flush_points () + {} + + template void VtkStream::write_cell( @@ -980,6 +1328,13 @@ namespace stream << '\n'; } + + void + VtkStream::flush_cells () + {} + + + VtuStream::VtuStream(std::ostream& out, const DataOutBase::VtkFlags f) : stream(out), flags(f) @@ -991,12 +1346,35 @@ namespace VtuStream::write_point (const unsigned int, const Point& p) { +#if !defined(HAVE_LIBZ) // write out coordinates stream << p; // fill with zeroes for (unsigned int i=dim; i<3; ++i) stream << " 0"; stream << '\n'; +#else + // if we want to compress, then + // first collect all the data in + // an array + for (unsigned int i=0; i=2) @@ -1023,7 +1402,55 @@ namespace << '\t' << start+d3+d2; } } - stream << '\n'; + stream << '\n'; +#else + cells.push_back (start); + cells.push_back (start+d1); + if (dim>=2) + { + cells.push_back (start+d2+d1); + cells.push_back (start+d2); + if (dim>=3) + { + cells.push_back (start+d3); + cells.push_back (start+d3+d1); + cells.push_back (start+d3+d2+d1); + cells.push_back (start+d3+d2); + } + } +#endif + } + + + + void + VtuStream::flush_cells () + { +#ifdef HAVE_LIBZ + // compress the data we have in + // memory and write them to the + // stream. then release the data + *this << cells << '\n'; + cells.clear (); +#endif + } + + + template + std::ostream& + VtuStream::operator<< (const std::vector &data) + { +#ifdef HAVE_LIBZ + // compress the data we have in + // memory and write them to the + // stream. then release the data + write_compressed_block (data, stream); +#else + for (unsigned int i=0; i void -DataOutBase::write_nodes ( - const std::vector >& patches, - STREAM& out) +DataOutBase::write_nodes (const std::vector >& patches, + STREAM& out) { Assert (dim<=3, ExcNotImplemented()); unsigned int count = 0; @@ -1785,7 +2211,8 @@ DataOutBase::write_nodes ( // it here. Point node; - for (typename std::vector >::const_iterator patch=patches.begin(); + for (typename std::vector >::const_iterator + patch=patches.begin(); patch!=patches.end(); ++patch) { const unsigned int n_subdivisions = patch->n_subdivisions; @@ -1810,12 +2237,13 @@ DataOutBase::write_nodes ( out.write_point(count++, node); } } + out.flush_points (); } template -void DataOutBase::write_cells( - const std::vector >& patches, - STREAM& out) +void +DataOutBase::write_cells (const std::vector >& patches, + STREAM& out) { Assert (dim<=3, ExcNotImplemented()); unsigned int count = 0; @@ -1823,7 +2251,8 @@ void DataOutBase::write_cells( // Array to hold all the node // numbers of a cell. 8 is // sufficient for 3D - for (typename std::vector >::const_iterator patch=patches.begin(); + for (typename std::vector >::const_iterator + patch=patches.begin(); patch!=patches.end(); ++patch) { const unsigned int n_subdivisions = patch->n_subdivisions; @@ -1833,9 +2262,6 @@ void DataOutBase::write_cells( const unsigned int n2 = (dim>1) ? n_subdivisions : 1; const unsigned int n3 = (dim>2) ? n_subdivisions : 1; // Offsets of outer loops -// const unsigned int d1 = 1; -// const unsigned int d2 = n_subdivisions+1; -// const unsigned int d3 = d2*d2; unsigned int d1 = 1; unsigned int d2 = n; unsigned int d3 = n*n; @@ -1851,6 +2277,8 @@ void DataOutBase::write_cells( // of the first vertex of this patch first_vertex_of_patch += Utilities::fixed_power(n_subdivisions+1); } + + out.flush_cells (); } @@ -4274,216 +4702,6 @@ DataOutBase::write_vtk (const std::vector > &patches, } -namespace -{ - // the functions in this namespace are - // taken from the libb64 project, see - // http://sourceforge.net/projects/libb64 - // - // libb64 has been placed in the public - // domain - namespace base64 - { - typedef enum - { - step_A, step_B, step_C - } base64_encodestep; - - typedef struct - { - base64_encodestep step; - char result; - } base64_encodestate; - - void base64_init_encodestate(base64_encodestate* state_in) - { - state_in->step = step_A; - state_in->result = 0; - } - - inline - char base64_encode_value(char value_in) - { - static const char* encoding - = "ABCDEFGHIJKLMNOPQRSTUVWXYZabcdefghijklmnopqrstuvwxyz0123456789+/"; - if (value_in > 63) return '='; - return encoding[(int)value_in]; - } - - int base64_encode_block(const char* plaintext_in, - int length_in, - char* code_out, - base64_encodestate *state_in) - { - const char* plainchar = plaintext_in; - const char* const plaintextend = plaintext_in + length_in; - char* codechar = code_out; - char result; - char fragment; - - result = state_in->result; - - switch (state_in->step) - { - while (1) - { - case step_A: - if (plainchar == plaintextend) - { - state_in->result = result; - state_in->step = step_A; - return codechar - code_out; - } - fragment = *plainchar++; - result = (fragment & 0x0fc) >> 2; - *codechar++ = base64_encode_value(result); - result = (fragment & 0x003) << 4; - case step_B: - if (plainchar == plaintextend) - { - state_in->result = result; - state_in->step = step_B; - return codechar - code_out; - } - fragment = *plainchar++; - result |= (fragment & 0x0f0) >> 4; - *codechar++ = base64_encode_value(result); - result = (fragment & 0x00f) << 2; - case step_C: - if (plainchar == plaintextend) - { - state_in->result = result; - state_in->step = step_C; - return codechar - code_out; - } - fragment = *plainchar++; - result |= (fragment & 0x0c0) >> 6; - *codechar++ = base64_encode_value(result); - result = (fragment & 0x03f) >> 0; - *codechar++ = base64_encode_value(result); - } - } - /* control should not reach here */ - return codechar - code_out; - } - - int base64_encode_blockend(char* code_out, base64_encodestate* state_in) - { - char* codechar = code_out; - - switch (state_in->step) - { - case step_B: - *codechar++ = base64_encode_value(state_in->result); - *codechar++ = '='; - *codechar++ = '='; - break; - case step_C: - *codechar++ = base64_encode_value(state_in->result); - *codechar++ = '='; - break; - case step_A: - break; - } - *codechar++ = '\0'; - - return codechar - code_out; - } - } - - - /** - * Do a base64 encoding of the given data. - * - * The function allocates memory as - * necessary and returns a pointer to - * it. The calling function must release - * this memory again. - */ - char * - encode_block (const char *data, - const int data_size) - { - base64::base64_encodestate state; - base64::base64_init_encodestate(&state); - - char *encoded_data = new char[sizeof(uint32_t)+2*data_size+1]; - - uint32_t int_header = data_size; - - const int encoded_length_header - = base64::base64_encode_block ((char*)&int_header, sizeof(int_header), - encoded_data, &state); - - const int encoded_length_data - = base64::base64_encode_block (data, data_size, - encoded_data + encoded_length_header, &state); - base64::base64_encode_blockend (encoded_data + encoded_length_header + encoded_length_data, - &state); - - return encoded_data; - } - - - -#ifdef HAVE_LIBZ - /** - * Do a zlib compression followed by a - * base64 encoding of the given data. - * - * The function allocates memory as - * necessary and returns a pointer to - * it. The calling function must release - * this memory again. - */ - template - std::pair - compress_and_encode_block (const std::vector &data) - { - if (data.size() != 0) - { - // allocate a buffer for compressing - // data and do so - uLongf compressed_data_length - = compressBound (data.size() * sizeof(T)); - char *compressed_data = new char[compressed_data_length]; - int err = compress2 ((Bytef *) compressed_data, - &compressed_data_length, - (const Bytef *) &data[0], - data.size() * sizeof(T), - Z_BEST_COMPRESSION); - Assert (err == Z_OK, ExcInternalError()); - - // now encode the compression header - const uint32_t compression_header[5] - = { 1, /* number of blocks */ - (uint32_t)(data.size() * sizeof(T)), /* size of block */ - 0, /* size of last block */ - (uint32_t)compressed_data_length, - 0}; /* list of compressed sizes of blocks */ - - char *encoded_header = encode_block ((char*)&compression_header[0], - 5 * sizeof(compression_header[0])); - - - // now do the encoding in base64 - char *encoded_data = encode_block (compressed_data, - compressed_data_length); - - // release the buffer for the - // compressed data and return the - // encoded data - delete[] compressed_data; - - return std::make_pair (encoded_header, encoded_data); - } - else - return std::pair(0,0); - } -#endif -} - - template void DataOutBase::write_vtu (const std::vector > &patches, @@ -4550,7 +4768,7 @@ DataOutBase::write_vtu (const std::vector > &patches, << "\n-->\n"; out << " > &patches, } +#ifdef HAVE_LIBZ + const char *ascii_or_binary = "binary"; +#else + const char *ascii_or_binary = "ascii"; +#endif + + // first count the number of cells // and cells for later use unsigned int n_nodes; @@ -4610,72 +4835,16 @@ DataOutBase::write_vtu (const std::vector > &patches, out << "\n"; out << " \n"; -#if !defined(HAVE_LIBZ) || 1 - out << " \n"; + out << " \n"; write_nodes(patches, vtu_out); out << " \n"; -#else - out << " \n"; - { - // collect vertices into one big array, - // then compress and encode it. we don't - // know up front how many vertices we - // need, so only allocate a typical guess - std::vector vertices; - vertices.reserve (patches.size() * 3 * GeometryInfo::vertices_per_cell); - for (typename std::vector >::const_iterator patch=patches.begin(); - patch!=patches.end(); ++patch) - { - const unsigned int n_subdivisions = patch->n_subdivisions; - const unsigned int n = n_subdivisions+1; - // Length of loops in all - // dimensions. If a dimension - // is not used, a loop of - // length one will do the job. - const unsigned int n1 = (dim>0) ? n : 1; - const unsigned int n2 = (dim>1) ? n : 1; - const unsigned int n3 = (dim>2) ? n : 1; - - for (unsigned int i3=0; i3 node; - compute_node(node, &*patch, - i1, - i2, - i3, - n_subdivisions); - for (unsigned int d=0; d - data - = compress_and_encode_block (vertices); - out << data.first << data.second; - - delete[] data.first; - delete[] data.second; - } - out << " \n"; -#endif out << " \n\n"; ///////////////////////////////// // now for the cells out << " \n"; - out << " \n"; + out << " \n"; write_cells(patches, vtu_out); out << " \n"; @@ -4683,20 +4852,27 @@ DataOutBase::write_vtu (const std::vector > &patches, // different than the VTK format, which // puts the number of nodes per cell in // front of the connectivity list. - out << " \n"; + out << " \n"; + + std::vector offsets (n_cells); for(unsigned int i=0; i::vertices_per_cell); + offsets[i] = (i+1)*GeometryInfo::vertices_per_cell; + vtu_out << offsets; out << "\n"; 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"; + { + // this should compress well :-) + std::vector cell_types (n_cells, vtk_cell_type[dim]); + vtu_out << cell_types; + } out << "\n"; out << " \n"; out << " \n"; @@ -4762,31 +4938,35 @@ DataOutBase::write_vtu (const std::vector > &patches, out << data_names[std_cxx1x::get<1>(vector_data_ranges[n_th_vector])]; } - out << "\" NumberOfComponents=\"3\" format=\"ascii\">\n"; + out << "\" NumberOfComponents=\"3\" format=\"" + << ascii_or_binary << "\">\n"; // now write data. pad all // vectors to have three // components + std::vector data; + data.reserve (n_nodes*dim); + 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'; + data.push_back (data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector]), n)); + data.push_back (0); + data.push_back (0); 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'; + data.push_back (data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector]), n)); + data.push_back (data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector])+1, n)); + data.push_back (0); 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'; + data.push_back (data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector]), n)); + data.push_back (data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector])+1, n)); + data.push_back (data_vectors(std_cxx1x::get<0>(vector_data_ranges[n_th_vector])+2, n)); break; default: @@ -4799,7 +4979,7 @@ DataOutBase::write_vtu (const std::vector > &patches, Assert (false, ExcInternalError()); } } - + vtu_out << data; out << " \n"; } @@ -4809,11 +4989,12 @@ DataOutBase::write_vtu (const std::vector > &patches, { out << " \n"; + << "\" format=\"" + << ascii_or_binary << "\">\n"; - std::copy (data_vectors[data_set].begin(), - data_vectors[data_set].end(), - std::ostream_iterator(out, "\n")); + std::vector data (data_vectors[data_set].begin(), + data_vectors[data_set].end()); + vtu_out << data; out << " \n"; } -- 2.39.5