From 58f7cede4785b4491bcbdcc517c0c6c89b45dd7a Mon Sep 17 00:00:00 2001 From: Wolfgang Bangerth Date: Thu, 23 Dec 2021 22:20:20 -0700 Subject: [PATCH] Separate 2d povray output from the general case. --- source/base/data_out_base.cc | 608 ++++++++++++++++++----------------- 1 file changed, 321 insertions(+), 287 deletions(-) diff --git a/source/base/data_out_base.cc b/source/base/data_out_base.cc index b1a7a278cf..f924359703 100644 --- a/source/base/data_out_base.cc +++ b/source/base/data_out_base.cc @@ -4068,44 +4068,52 @@ namespace DataOutBase } - - template - void - write_povray( - const std::vector> &patches, - const std::vector & data_names, - const std::vector< - std::tuple> &, - const PovrayFlags &flags, - std::ostream & out) + namespace { - AssertThrow(out, ExcIO()); + template + void + do_write_povray(const std::vector> &, + const std::vector &, + const PovrayFlags &, + std::ostream &) + { + Assert(false, + ExcMessage("Writing files in POVRAY format is only supported " + "for two-dimensional meshes.")); + } + + + + void + do_write_povray(const std::vector> &patches, + const std::vector &data_names, + const PovrayFlags & flags, + std::ostream & out) + { + AssertThrow(out, ExcIO()); #ifndef DEAL_II_WITH_MPI - // verify that there are indeed patches to be written out. most - // of the times, people just forget to call build_patches when there - // are no patches, so a warning is in order. that said, the - // assertion is disabled if we support MPI since then it can - // happen that on the coarsest mesh, a processor simply has no cells it - // actually owns, and in that case it is legit if there are no patches - Assert(patches.size() > 0, ExcNoPatches()); + // verify that there are indeed patches to be written out. most + // of the times, people just forget to call build_patches when there + // are no patches, so a warning is in order. that said, the + // assertion is disabled if we support MPI since then it can + // happen that on the coarsest mesh, a processor simply has no cells it + // actually owns, and in that case it is legit if there are no patches + Assert(patches.size() > 0, ExcNoPatches()); #else - if (patches.size() == 0) - return; + if (patches.size() == 0) + return; #endif - Assert(dim == 2, - ExcNotImplemented()); // only for 2-D surfaces on a 2-D plane - Assert(spacedim == 2, ExcNotImplemented()); + constexpr int dim = 2; + constexpr int spacedim = 2; - const unsigned int n_data_sets = data_names.size(); - (void)n_data_sets; + const unsigned int n_data_sets = data_names.size(); + (void)n_data_sets; - // write preamble - { - out << "/* This file was generated by the deal.II library." << '\n' + // write preamble + { + out + << "/* This file was generated by the deal.II library." << '\n' << " Date = " << Utilities::System::get_date() << '\n' << " Time = " << Utilities::System::get_time() << '\n' << '\n' @@ -4113,285 +4121,311 @@ namespace DataOutBase << '\n' << "*/ " << '\n'; - // include files - out << "#include \"colors.inc\" " << '\n' - << "#include \"textures.inc\" " << '\n'; + // include files + out << "#include \"colors.inc\" " << '\n' + << "#include \"textures.inc\" " << '\n'; + + + // use external include file for textures, camera and light + if (flags.external_data) + out << "#include \"data.inc\" " << '\n'; + else // all definitions in data file + { + // camera + out << '\n' + << '\n' + << "camera {" << '\n' + << " location <1,4,-7>" << '\n' + << " look_at <0,0,0>" << '\n' + << " angle 30" << '\n' + << "}" << '\n'; + + // light + out << '\n' + << "light_source {" << '\n' + << " <1,4,-7>" << '\n' + << " color Grey" << '\n' + << "}" << '\n'; + out << '\n' + << "light_source {" << '\n' + << " <0,20,0>" << '\n' + << " color White" << '\n' + << "}" << '\n'; + } + } + // max. and min. height of solution + Assert(patches.size() > 0, ExcNoPatches()); + double hmin = patches[0].data(0, 0); + double hmax = patches[0].data(0, 0); - // use external include file for textures, camera and light - if (flags.external_data) - out << "#include \"data.inc\" " << '\n'; - else // all definitions in data file + for (const auto &patch : patches) { - // camera - out << '\n' - << '\n' - << "camera {" << '\n' - << " location <1,4,-7>" << '\n' - << " look_at <0,0,0>" << '\n' - << " angle 30" << '\n' - << "}" << '\n'; - - // light - out << '\n' - << "light_source {" << '\n' - << " <1,4,-7>" << '\n' - << " color Grey" << '\n' - << "}" << '\n'; - out << '\n' - << "light_source {" << '\n' - << " <0,20,0>" << '\n' - << " color White" << '\n' - << "}" << '\n'; - } - } + const unsigned int n_subdivisions = patch.n_subdivisions; - // max. and min. height of solution - Assert(patches.size() > 0, ExcNoPatches()); - double hmin = patches[0].data(0, 0); - double hmax = patches[0].data(0, 0); + Assert((patch.data.n_rows() == n_data_sets && + !patch.points_are_available) || + (patch.data.n_rows() == n_data_sets + spacedim && + patch.points_are_available), + ExcDimensionMismatch(patch.points_are_available ? + (n_data_sets + spacedim) : + n_data_sets, + patch.data.n_rows())); + Assert(patch.data.n_cols() == + Utilities::fixed_power(n_subdivisions + 1), + ExcInvalidDatasetSize(patch.data.n_cols(), + n_subdivisions + 1)); - for (const auto &patch : patches) - { - const unsigned int n_subdivisions = patch.n_subdivisions; + for (unsigned int i = 0; i < n_subdivisions + 1; ++i) + for (unsigned int j = 0; j < n_subdivisions + 1; ++j) + { + const int dl = i * (n_subdivisions + 1) + j; + if (patch.data(0, dl) < hmin) + hmin = patch.data(0, dl); + if (patch.data(0, dl) > hmax) + hmax = patch.data(0, dl); + } + } - Assert((patch.data.n_rows() == n_data_sets && - !patch.points_are_available) || - (patch.data.n_rows() == n_data_sets + spacedim && - patch.points_are_available), - ExcDimensionMismatch(patch.points_are_available ? - (n_data_sets + spacedim) : - n_data_sets, - patch.data.n_rows())); - Assert(patch.data.n_cols() == - Utilities::fixed_power(n_subdivisions + 1), - ExcInvalidDatasetSize(patch.data.n_cols(), n_subdivisions + 1)); + out << "#declare HMIN=" << hmin << ";" << '\n' + << "#declare HMAX=" << hmax << ";" << '\n' + << '\n'; - for (unsigned int i = 0; i < n_subdivisions + 1; ++i) - for (unsigned int j = 0; j < n_subdivisions + 1; ++j) - { - const int dl = i * (n_subdivisions + 1) + j; - if (patch.data(0, dl) < hmin) - hmin = patch.data(0, dl); - if (patch.data(0, dl) > hmax) - hmax = patch.data(0, dl); - } - } + if (!flags.external_data) + { + // texture with scaled niveau lines 10 lines in the surface + out << "#declare Tex=texture{" << '\n' + << " pigment {" << '\n' + << " gradient y" << '\n' + << " scale y*(HMAX-HMIN)*" << 0.1 << '\n' + << " color_map {" << '\n' + << " [0.00 color Light_Purple] " << '\n' + << " [0.95 color Light_Purple] " << '\n' + << " [1.00 color White] " << '\n' + << "} } }" << '\n' + << '\n'; + } - out << "#declare HMIN=" << hmin << ";" << '\n' - << "#declare HMAX=" << hmax << ";" << '\n' - << '\n'; + if (!flags.bicubic_patch) + { + // start of mesh header + out << '\n' << "mesh {" << '\n'; + } - if (!flags.external_data) - { - // texture with scaled niveau lines 10 lines in the surface - out << "#declare Tex=texture{" << '\n' - << " pigment {" << '\n' - << " gradient y" << '\n' - << " scale y*(HMAX-HMIN)*" << 0.1 << '\n' - << " color_map {" << '\n' - << " [0.00 color Light_Purple] " << '\n' - << " [0.95 color Light_Purple] " << '\n' - << " [1.00 color White] " << '\n' - << "} } }" << '\n' - << '\n'; - } + // loop over all patches + for (const auto &patch : patches) + { + const unsigned int n_subdivisions = patch.n_subdivisions; + const unsigned int n = n_subdivisions + 1; + const unsigned int d1 = 1; + const unsigned int d2 = n; - if (!flags.bicubic_patch) - { - // start of mesh header - out << '\n' << "mesh {" << '\n'; - } + Assert((patch.data.n_rows() == n_data_sets && + !patch.points_are_available) || + (patch.data.n_rows() == n_data_sets + spacedim && + patch.points_are_available), + ExcDimensionMismatch(patch.points_are_available ? + (n_data_sets + spacedim) : + n_data_sets, + patch.data.n_rows())); + Assert(patch.data.n_cols() == Utilities::fixed_power(n), + ExcInvalidDatasetSize(patch.data.n_cols(), + n_subdivisions + 1)); - // loop over all patches - for (const auto &patch : patches) - { - const unsigned int n_subdivisions = patch.n_subdivisions; - const unsigned int n = n_subdivisions + 1; - const unsigned int d1 = 1; - const unsigned int d2 = n; - Assert((patch.data.n_rows() == n_data_sets && - !patch.points_are_available) || - (patch.data.n_rows() == n_data_sets + spacedim && - patch.points_are_available), - ExcDimensionMismatch(patch.points_are_available ? - (n_data_sets + spacedim) : - n_data_sets, - patch.data.n_rows())); - Assert(patch.data.n_cols() == Utilities::fixed_power(n), - ExcInvalidDatasetSize(patch.data.n_cols(), n_subdivisions + 1)); + std::vector> ver(n * n); + for (unsigned int i2 = 0; i2 < n; ++i2) + for (unsigned int i1 = 0; i1 < n; ++i1) + { + // compute coordinates for this patch point, storing in ver + ver[i1 * d1 + i2 * d2] = + get_equispaced_location(patch, {i1, i2}, n_subdivisions); + } - std::vector> ver(n * n); - for (unsigned int i2 = 0; i2 < n; ++i2) - for (unsigned int i1 = 0; i1 < n; ++i1) + if (!flags.bicubic_patch) { - // compute coordinates for this patch point, storing in ver - ver[i1 * d1 + i2 * d2] = - get_equispaced_location(patch, {i1, i2}, n_subdivisions); + // approximate normal vectors in patch + std::vector> nrml; + // only if smooth triangles are used + if (flags.smooth) + { + nrml.resize(n * n); + // These are difference quotients of the surface + // mapping. We take them symmetric inside the + // patch and one-sided at the edges + Point<3> h1, h2; + // Now compute normals in every point + for (unsigned int i = 0; i < n; ++i) + for (unsigned int j = 0; j < n; ++j) + { + const unsigned int il = (i == 0) ? i : (i - 1); + const unsigned int ir = + (i == n_subdivisions) ? i : (i + 1); + const unsigned int jl = (j == 0) ? j : (j - 1); + const unsigned int jr = + (j == n_subdivisions) ? j : (j + 1); + + h1(0) = + ver[ir * d1 + j * d2](0) - ver[il * d1 + j * d2](0); + h1(1) = patch.data(0, ir * d1 + j * d2) - + patch.data(0, il * d1 + j * d2); + h1(2) = + ver[ir * d1 + j * d2](1) - ver[il * d1 + j * d2](1); + + h2(0) = + ver[i * d1 + jr * d2](0) - ver[i * d1 + jl * d2](0); + h2(1) = patch.data(0, i * d1 + jr * d2) - + patch.data(0, i * d1 + jl * d2); + h2(2) = + ver[i * d1 + jr * d2](1) - ver[i * d1 + jl * d2](1); + + nrml[i * d1 + j * d2](0) = + h1(1) * h2(2) - h1(2) * h2(1); + nrml[i * d1 + j * d2](1) = + h1(2) * h2(0) - h1(0) * h2(2); + nrml[i * d1 + j * d2](2) = + h1(0) * h2(1) - h1(1) * h2(0); + + // normalize Vector + double norm = + std::sqrt(std::pow(nrml[i * d1 + j * d2](0), 2.) + + std::pow(nrml[i * d1 + j * d2](1), 2.) + + std::pow(nrml[i * d1 + j * d2](2), 2.)); + + if (nrml[i * d1 + j * d2](1) < 0) + norm *= -1.; + + for (unsigned int k = 0; k < 3; ++k) + nrml[i * d1 + j * d2](k) /= norm; + } + } + + // setting up triangles + for (unsigned int i = 0; i < n_subdivisions; ++i) + for (unsigned int j = 0; j < n_subdivisions; ++j) + { + // down/left vertex of triangle + const int dl = i * d1 + j * d2; + if (flags.smooth) + { + // writing smooth_triangles + + // down/right triangle + out << "smooth_triangle {" << '\n' + << "\t<" << ver[dl](0) << "," << patch.data(0, dl) + << "," << ver[dl](1) << ">, <" << nrml[dl](0) + << ", " << nrml[dl](1) << ", " << nrml[dl](2) + << ">," << '\n'; + out << " \t<" << ver[dl + d1](0) << "," + << patch.data(0, dl + d1) << "," << ver[dl + d1](1) + << ">, <" << nrml[dl + d1](0) << ", " + << nrml[dl + d1](1) << ", " << nrml[dl + d1](2) + << ">," << '\n'; + out << "\t<" << ver[dl + d1 + d2](0) << "," + << patch.data(0, dl + d1 + d2) << "," + << ver[dl + d1 + d2](1) << ">, <" + << nrml[dl + d1 + d2](0) << ", " + << nrml[dl + d1 + d2](1) << ", " + << nrml[dl + d1 + d2](2) << ">}" << '\n'; + + // upper/left triangle + out << "smooth_triangle {" << '\n' + << "\t<" << ver[dl](0) << "," << patch.data(0, dl) + << "," << ver[dl](1) << ">, <" << nrml[dl](0) + << ", " << nrml[dl](1) << ", " << nrml[dl](2) + << ">," << '\n'; + out << "\t<" << ver[dl + d1 + d2](0) << "," + << patch.data(0, dl + d1 + d2) << "," + << ver[dl + d1 + d2](1) << ">, <" + << nrml[dl + d1 + d2](0) << ", " + << nrml[dl + d1 + d2](1) << ", " + << nrml[dl + d1 + d2](2) << ">," << '\n'; + out << "\t<" << ver[dl + d2](0) << "," + << patch.data(0, dl + d2) << "," << ver[dl + d2](1) + << ">, <" << nrml[dl + d2](0) << ", " + << nrml[dl + d2](1) << ", " << nrml[dl + d2](2) + << ">}" << '\n'; + } + else + { + // writing standard triangles down/right triangle + out << "triangle {" << '\n' + << "\t<" << ver[dl](0) << "," << patch.data(0, dl) + << "," << ver[dl](1) << ">," << '\n'; + out << "\t<" << ver[dl + d1](0) << "," + << patch.data(0, dl + d1) << "," << ver[dl + d1](1) + << ">," << '\n'; + out << "\t<" << ver[dl + d1 + d2](0) << "," + << patch.data(0, dl + d1 + d2) << "," + << ver[dl + d1 + d2](1) << ">}" << '\n'; + + // upper/left triangle + out << "triangle {" << '\n' + << "\t<" << ver[dl](0) << "," << patch.data(0, dl) + << "," << ver[dl](1) << ">," << '\n'; + out << "\t<" << ver[dl + d1 + d2](0) << "," + << patch.data(0, dl + d1 + d2) << "," + << ver[dl + d1 + d2](1) << ">," << '\n'; + out << "\t<" << ver[dl + d2](0) << "," + << patch.data(0, dl + d2) << "," << ver[dl + d2](1) + << ">}" << '\n'; + } + } } + else + { + // writing bicubic_patch + Assert(n_subdivisions == 3, + ExcDimensionMismatch(n_subdivisions, 3)); + out << '\n' + << "bicubic_patch {" << '\n' + << " type 0" << '\n' + << " flatness 0" << '\n' + << " u_steps 0" << '\n' + << " v_steps 0" << '\n'; + for (int i = 0; i < 16; ++i) + { + out << "\t<" << ver[i](0) << "," << patch.data(0, i) << "," + << ver[i](1) << ">"; + if (i != 15) + out << ","; + out << '\n'; + } + out << " texture {Tex}" << '\n' << "}" << '\n'; + } + } + if (!flags.bicubic_patch) + { + // the end of the mesh + out << " texture {Tex}" << '\n' << "}" << '\n' << '\n'; + } - if (!flags.bicubic_patch) - { - // approximate normal vectors in patch - std::vector> nrml; - // only if smooth triangles are used - if (flags.smooth) - { - nrml.resize(n * n); - // These are difference quotients of the surface - // mapping. We take them symmetric inside the - // patch and one-sided at the edges - Point<3> h1, h2; - // Now compute normals in every point - for (unsigned int i = 0; i < n; ++i) - for (unsigned int j = 0; j < n; ++j) - { - const unsigned int il = (i == 0) ? i : (i - 1); - const unsigned int ir = - (i == n_subdivisions) ? i : (i + 1); - const unsigned int jl = (j == 0) ? j : (j - 1); - const unsigned int jr = - (j == n_subdivisions) ? j : (j + 1); - - h1(0) = - ver[ir * d1 + j * d2](0) - ver[il * d1 + j * d2](0); - h1(1) = patch.data(0, ir * d1 + j * d2) - - patch.data(0, il * d1 + j * d2); - h1(2) = - ver[ir * d1 + j * d2](1) - ver[il * d1 + j * d2](1); - - h2(0) = - ver[i * d1 + jr * d2](0) - ver[i * d1 + jl * d2](0); - h2(1) = patch.data(0, i * d1 + jr * d2) - - patch.data(0, i * d1 + jl * d2); - h2(2) = - ver[i * d1 + jr * d2](1) - ver[i * d1 + jl * d2](1); - - nrml[i * d1 + j * d2](0) = h1(1) * h2(2) - h1(2) * h2(1); - nrml[i * d1 + j * d2](1) = h1(2) * h2(0) - h1(0) * h2(2); - nrml[i * d1 + j * d2](2) = h1(0) * h2(1) - h1(1) * h2(0); - - // normalize Vector - double norm = - std::sqrt(std::pow(nrml[i * d1 + j * d2](0), 2.) + - std::pow(nrml[i * d1 + j * d2](1), 2.) + - std::pow(nrml[i * d1 + j * d2](2), 2.)); - - if (nrml[i * d1 + j * d2](1) < 0) - norm *= -1.; - - for (unsigned int k = 0; k < 3; ++k) - nrml[i * d1 + j * d2](k) /= norm; - } - } + // make sure everything now gets to disk + out.flush(); - // setting up triangles - for (unsigned int i = 0; i < n_subdivisions; ++i) - for (unsigned int j = 0; j < n_subdivisions; ++j) - { - // down/left vertex of triangle - const int dl = i * d1 + j * d2; - if (flags.smooth) - { - // writing smooth_triangles - - // down/right triangle - out << "smooth_triangle {" << '\n' - << "\t<" << ver[dl](0) << "," << patch.data(0, dl) - << "," << ver[dl](1) << ">, <" << nrml[dl](0) << ", " - << nrml[dl](1) << ", " << nrml[dl](2) << ">," << '\n'; - out << " \t<" << ver[dl + d1](0) << "," - << patch.data(0, dl + d1) << "," << ver[dl + d1](1) - << ">, <" << nrml[dl + d1](0) << ", " - << nrml[dl + d1](1) << ", " << nrml[dl + d1](2) - << ">," << '\n'; - out << "\t<" << ver[dl + d1 + d2](0) << "," - << patch.data(0, dl + d1 + d2) << "," - << ver[dl + d1 + d2](1) << ">, <" - << nrml[dl + d1 + d2](0) << ", " - << nrml[dl + d1 + d2](1) << ", " - << nrml[dl + d1 + d2](2) << ">}" << '\n'; - - // upper/left triangle - out << "smooth_triangle {" << '\n' - << "\t<" << ver[dl](0) << "," << patch.data(0, dl) - << "," << ver[dl](1) << ">, <" << nrml[dl](0) << ", " - << nrml[dl](1) << ", " << nrml[dl](2) << ">," << '\n'; - out << "\t<" << ver[dl + d1 + d2](0) << "," - << patch.data(0, dl + d1 + d2) << "," - << ver[dl + d1 + d2](1) << ">, <" - << nrml[dl + d1 + d2](0) << ", " - << nrml[dl + d1 + d2](1) << ", " - << nrml[dl + d1 + d2](2) << ">," << '\n'; - out << "\t<" << ver[dl + d2](0) << "," - << patch.data(0, dl + d2) << "," << ver[dl + d2](1) - << ">, <" << nrml[dl + d2](0) << ", " - << nrml[dl + d2](1) << ", " << nrml[dl + d2](2) - << ">}" << '\n'; - } - else - { - // writing standard triangles down/right triangle - out << "triangle {" << '\n' - << "\t<" << ver[dl](0) << "," << patch.data(0, dl) - << "," << ver[dl](1) << ">," << '\n'; - out << "\t<" << ver[dl + d1](0) << "," - << patch.data(0, dl + d1) << "," << ver[dl + d1](1) - << ">," << '\n'; - out << "\t<" << ver[dl + d1 + d2](0) << "," - << patch.data(0, dl + d1 + d2) << "," - << ver[dl + d1 + d2](1) << ">}" << '\n'; - - // upper/left triangle - out << "triangle {" << '\n' - << "\t<" << ver[dl](0) << "," << patch.data(0, dl) - << "," << ver[dl](1) << ">," << '\n'; - out << "\t<" << ver[dl + d1 + d2](0) << "," - << patch.data(0, dl + d1 + d2) << "," - << ver[dl + d1 + d2](1) << ">," << '\n'; - out << "\t<" << ver[dl + d2](0) << "," - << patch.data(0, dl + d2) << "," << ver[dl + d2](1) - << ">}" << '\n'; - } - } - } - else - { - // writing bicubic_patch - Assert(n_subdivisions == 3, - ExcDimensionMismatch(n_subdivisions, 3)); - out << '\n' - << "bicubic_patch {" << '\n' - << " type 0" << '\n' - << " flatness 0" << '\n' - << " u_steps 0" << '\n' - << " v_steps 0" << '\n'; - for (int i = 0; i < 16; ++i) - { - out << "\t<" << ver[i](0) << "," << patch.data(0, i) << "," - << ver[i](1) << ">"; - if (i != 15) - out << ","; - out << '\n'; - } - out << " texture {Tex}" << '\n' << "}" << '\n'; - } - } + AssertThrow(out, ExcIO()); + } + } // namespace - if (!flags.bicubic_patch) - { - // the end of the mesh - out << " texture {Tex}" << '\n' << "}" << '\n' << '\n'; - } - // make sure everything now gets to disk - out.flush(); - AssertThrow(out, ExcIO()); + template + void + write_povray( + const std::vector> &patches, + const std::vector & data_names, + const std::vector< + std::tuple> &, + const PovrayFlags &flags, + std::ostream & out) + { + do_write_povray(patches, data_names, flags, out); } -- 2.39.5