diff --git a/.gitignore b/.gitignore index f4dd1f61..cc966aa2 100644 --- a/.gitignore +++ b/.gitignore @@ -56,7 +56,7 @@ main-debug main #Output files -*.vtk +applications/*/*.vtk *.vtu *.vts *.frt diff --git a/include/field_input/IntegrationTools/PField.hh b/include/field_input/IntegrationTools/PField.hh index 8dff31c7..134878c7 100644 --- a/include/field_input/IntegrationTools/PField.hh +++ b/include/field_input/IntegrationTools/PField.hh @@ -1,2 +1,2 @@ -#include "pfield/Body.hh" -#include "pfield/Coordinate.hh" \ No newline at end of file +#include +#include \ No newline at end of file diff --git a/include/field_input/IntegrationTools/pfield/Body.hh b/include/field_input/IntegrationTools/pfield/Body.hh index 9d28d020..649f29f9 100644 --- a/include/field_input/IntegrationTools/pfield/Body.hh +++ b/include/field_input/IntegrationTools/pfield/Body.hh @@ -1,8 +1,8 @@ #ifndef Body_HH #define Body_HH -#include "./Mesh.hh" -#include "./PField.hh" +#include +#include #include #include #include @@ -11,8 +11,9 @@ namespace PRISMS { - /// A class for a Body: a combination of Mesh and Field(s)) - /// + /** + * \brief A class for a Body-a combination of Mesh and Field(s)-used in file read-in. + */ template class Body { @@ -21,109 +22,104 @@ namespace PRISMS std::vector> scalar_field; - // std::vector< PField, DIM > > vector_field; + /** + * \brief Constructor. + */ + Body() = default; - // std::vector< PField, DIM > > tensor_field; - - // ---------------------------------------------------------- - // Constructors - Body() {}; - - /// Read from a 2D vtk file - /// For now: - /// only ASCII files - /// only rectilinear grids (though output as UNSTRUCTURED_GRID) - /// only (2d) Quad elements - /// + /** + * \brief Read from a 2D/3D vtk file with quad/hex elements + */ void read_vtk(const std::string &vtkfile) { - std::cout << "Begin reading unstructured vtk file" << std::endl; - - // read in vtk file here - std::ifstream infile_mesh(vtkfile.c_str()); - - // read mesh info - mesh.read_vtk(infile_mesh); + std::cout << "Begin reading unstructured vtk file\n"; + // Read vtk file std::ifstream infile(vtkfile.c_str()); + if (!infile.is_open()) + { + throw std::runtime_error("Could not open VTK file: " + vtkfile); + } - // read point data - std::istringstream ss; - std::string str, name, type, line; - int numcomp; - unsigned long int Npoints, u, p; + // Read mesh info + mesh.read_vtk(infile); - while (!infile.eof()) + // Read point data + std::istringstream ss; + std::string str; + std::string name; + std::string type; + std::string line; + int numcomp = 0; + unsigned long int Npoints = 0; + + while (std::getline(infile, line)) { - std::getline(infile, line); + // Check for Point data + if (line.rfind("POINT_DATA", 0) == 0) + { + ss.clear(); + ss.str(line); + ss >> str >> Npoints; + } - if (line[0] == 'P') + // Check for Scalar data + if (line.rfind("SCALARS", 0) == 0) { - if (line.size() > 9 && line.substr(0, 10) == "POINT_DATA") + ss.clear(); + ss.str(line); + ss >> str >> name >> type >> numcomp; + + // Skip LOOKUP_TABLE line that follows SCALAR line + std::getline(infile, line); + + // Read data + std::vector gid(Npoints); + for (auto &value : gid) { - // std::cout << line << "\n"; - ss.clear(); - ss.str(line); - ss >> str >> Npoints; + infile >> value; } - } - if (line[0] == 'S') - { - if (line.size() > 6 && line.substr(0, 7) == "SCALARS") + // Construct field descriptions + std::vector var_name(DIM); + std::vector var_description(DIM); + switch (DIM) { - ss.clear(); - ss.str(line); - ss >> str >> name >> type >> numcomp; + case 1: + throw std::runtime_error( + "1D dimensional vtk read-in is not currently supported."); - // read LOOKUP_TABLE line - std::getline(infile, line); + break; - // read data - std::cout << "begin reading data" << std::endl; + case 2: + var_name[0] = "x"; + var_description[0] = "x coordinate"; + var_name[1] = "y"; + var_description[1] = "y coordinate"; - std::vector gid(Npoints); - for (unsigned int i = 0; i < Npoints; i++) - { - infile >> gid[i]; - // std::cout << data[i] << std::endl; - } - std::cout << " done" << std::endl; + break; - // construct field - std::vector var_name(DIM); - std::vector var_description(DIM); + case 3: + var_name[0] = "x"; + var_description[0] = "x coordinate"; + var_name[1] = "y"; + var_description[1] = "y coordinate"; + var_name[2] = "z"; + var_description[2] = "z coordinate"; - if (DIM == 2) - { - var_name[0] = "x"; - var_description[0] = "x coordinate"; - var_name[1] = "y"; - var_description[1] = "y coordinate"; - } - if (DIM > 2) - { - var_name[2] = "z"; - var_description[2] = "z coordinate"; - } + break; - std::cout << "Construct PField '" << name << "'" << std::endl; - scalar_field.push_back(PField(name, - var_name, - var_description, - mesh, - gid, - 0.0)); - std::cout << " done" << std::endl; - - gid.clear(); - // + default: + throw std::runtime_error("Invalid dimension for vtk read-in."); } + + // Construct field + std::cout << "\tConstructing PField '" << name << "'\n"; + scalar_field.emplace_back(name, var_name, var_description, mesh, gid, 0.0); } } - - infile.close(); + std::cout << "Finished reading VTK file." << std::endl; } void @@ -327,11 +323,17 @@ namespace PRISMS for (unsigned int i = 0; i < scalar_field.size(); i++) { if (scalar_field[i].name() == name) - return scalar_field[i]; + { + return scalar_field[i]; + } } throw std::invalid_argument("Could not find scalar_field named '" + name + "'"); } }; + + template class Body; + template class Body; + } // namespace PRISMS #endif diff --git a/include/field_input/IntegrationTools/pfield/Coordinate.hh b/include/field_input/IntegrationTools/pfield/Coordinate.hh index b2fc93c5..1c31188e 100644 --- a/include/field_input/IntegrationTools/pfield/Coordinate.hh +++ b/include/field_input/IntegrationTools/pfield/Coordinate.hh @@ -1,42 +1,56 @@ - #ifndef Coordinate_HH #define Coordinate_HH +#include + namespace PRISMS { - /// A class for a coordinate, templated by dimension - /// This is a possible option anyplace 'Coordinate' class template is used - /// but it is not the only option. Any class that implements - /// 'Coordinate::operator[]()' should work - - template + /** + * \brief A class for a coordinate, templated by dimension. This is a possible option + * anyplace 'Coordinate' class template is used but it is not the only option. Any class + * that implements 'Coordinate::operator[]()' should work. + */ + template class Coordinate { - float _coord[DIM]; - public: + /** + * \brief Get the number of dimensions. + */ [[nodiscard]] int size() const { - return DIM; + return dim; } + /** + * \brief Get the ith coordinate. + */ float & operator[](int i) { - return _coord[i]; + return coordinate[i]; } + /** + * \brief Get a const reference to the ith coordinate. + */ const float & operator[](int i) const { - return _coord[i]; + return coordinate[i]; } + /** + * \brief Friend function for output stream. + */ template friend std::ostream & operator<<(std::ostream &outstream, const Coordinate &coord); + + private: + float coordinate[dim]; }; template diff --git a/include/field_input/IntegrationTools/pfield/Mesh.hh b/include/field_input/IntegrationTools/pfield/Mesh.hh index 413d369e..8d3a879e 100644 --- a/include/field_input/IntegrationTools/pfield/Mesh.hh +++ b/include/field_input/IntegrationTools/pfield/Mesh.hh @@ -1,14 +1,15 @@ - #ifndef Mesh_HH #define Mesh_HH -#include "../datastruc/Bin.hh" -#include "../pfunction/PFuncBase.hh" -#include "./interpolation/Hexahedron.hh" -#include "./interpolation/Interpolator.hh" -#include "./interpolation/Quad.hh" +#include + #include #include +#include +#include +#include +#include +#include #include #include @@ -19,32 +20,22 @@ namespace PRISMS construct_basis_function(PFuncBase>, double> *&bfunc, const std::string &name) { - if (name == "Quad") - { - bfunc = new Quad(); - } - else - { - std::cout << "Error in construct_basis_function (2D): unknown name: " << name - << std::endl; - exit(1); - } + AssertThrow(name == "Quad", + dealii::ExcMessage( + "Error in construct_basis_function (2D): unknown name: " + name)); + + bfunc = new Quad(); } inline void construct_basis_function(PFuncBase>, double> *&bfunc, const std::string &name) { - if (name == "Hexahedron") - { - bfunc = new Hexahedron(); - } - else - { - std::cout << "Error in construct_basis_function (3D): unknown name: " << name - << std::endl; - exit(1); - } + AssertThrow(name == "Hexahedron", + dealii::ExcMessage( + "Error in construct_basis_function (3D): unknown name: " + name)); + + bfunc = new Hexahedron(); } template @@ -57,34 +48,27 @@ namespace PRISMS const std::vector &cell_node, const std::vector> &node) { - if (name == "Quad") - { - Interpolator *interp_ptr; - - // std::cout << "cell nodes: " << cell_node[0] << " " << cell_node[2] << - // std::endl; - - PRISMS::Coordinate<2> dim; - dim[0] = node[cell_node[2]][0] - node[cell_node[0]][0]; - dim[1] = node[cell_node[2]][1] - node[cell_node[0]][1]; - - // QuadValues(const Coordinate &node, const Coordinate &dim, int node_index) - for (int j = 0; j < 4; j++) - { - interp.push_back(interp_ptr); - interp.back() = new PRISMS::QuadValues(cell_node[j], - cell, - bfunc_ptr, - node[cell_node[j]], - dim, - j); - } - } - else + AssertThrow(name == "Quad", + dealii::ExcMessage( + "Error in construct_interpolating_function (2D): unknown name: " + + name)); + + Interpolator *interp_ptr = nullptr; + + PRISMS::Coordinate<2> dim; + dim[0] = node[cell_node[2]][0] - node[cell_node[0]][0]; + dim[1] = node[cell_node[2]][1] - node[cell_node[0]][1]; + + // QuadValues(const Coordinate &node, const Coordinate &dim, int node_index) + for (int j = 0; j < 4; j++) { - std::cout << "Error in construct_interpolating_function (2D): unknown name: " - << name << std::endl; - exit(1); + interp.push_back(interp_ptr); + interp.back() = new PRISMS::QuadValues(cell_node[j], + cell, + bfunc_ptr, + node[cell_node[j]], + dim, + j); } } @@ -98,32 +82,28 @@ namespace PRISMS const std::vector &cell_node, const std::vector> &node) { - if (name == "Hexahedron") - { - Interpolator *interp_ptr; - - PRISMS::Coordinate<3> dim; - dim[0] = node[cell_node[6]][0] - node[cell_node[0]][0]; - dim[1] = node[cell_node[6]][1] - node[cell_node[0]][1]; - dim[2] = node[cell_node[6]][2] - node[cell_node[0]][2]; - - // QuadValues(const Coordinate &node, const Coordinate &dim, int node_index) - for (int j = 0; j < 8; j++) - { - interp.push_back(interp_ptr); - interp.back() = new PRISMS::HexahedronValues(cell_node[j], - cell, - bfunc_ptr, - node[cell_node[j]], - dim, - j); - } - } - else + AssertThrow(name == "Hexahedron", + dealii::ExcMessage( + "Error in construct_interpolating_function (3D): unknown name: " + + name)); + + Interpolator *interp_ptr = nullptr; + + PRISMS::Coordinate<3> dim; + dim[0] = node[cell_node[6]][0] - node[cell_node[0]][0]; + dim[1] = node[cell_node[6]][1] - node[cell_node[0]][1]; + dim[2] = node[cell_node[6]][2] - node[cell_node[0]][2]; + + // QuadValues(const Coordinate &node, const Coordinate &dim, int node_index) + for (int j = 0; j < 8; j++) { - std::cout << "Error in construct_interpolating_function (3D): unknown name: " - << name << std::endl; - exit(1); + interp.push_back(interp_ptr); + interp.back() = new PRISMS::HexahedronValues(cell_node[j], + cell, + bfunc_ptr, + node[cell_node[j]], + dim, + j); } } @@ -161,8 +141,7 @@ namespace PRISMS Bin *, Coordinate> _bin; public: - // still need a constructor - Mesh() {}; + Mesh() = default; ~Mesh() { @@ -180,15 +159,20 @@ namespace PRISMS void read_vtk(std::ifstream &infile) { - std::cout << "Read unstructured mesh" << std::endl; + std::cout << "\tReading unstructured mesh\n"; std::istringstream ss; - std::string line, str, type; + std::string line; + std::string str; + std::string type; unsigned int uli_dummy; double d_dummy; - unsigned long int Npoints, Ncells, Ncell_numbers, u; + unsigned long int Npoints; + unsigned long int Ncells; + unsigned long int Ncell_numbers; + unsigned long int u; std::vector cell_node; PRISMS::Coordinate _coord; @@ -200,7 +184,6 @@ namespace PRISMS while (!infile.eof()) { std::getline(infile, line); - // std::cout << "line: " << line << std::endl; if (line[0] == 'P') { @@ -212,7 +195,6 @@ namespace PRISMS if (line.size() > 5 && line.substr(0, 6) == "POINTS") { // read header line - // std::cout << line << "\n"; ss.clear(); ss.str(line); ss >> str >> Npoints >> type; @@ -221,38 +203,30 @@ namespace PRISMS std::vector> value(DIM); std::vector> hist(DIM); - std::cout << "Read POINTS: " << Npoints << std::endl; + std::cout << "\tReading POINTS: " << Npoints << std::endl; _node.reserve(Npoints); - std::cout << " reserve OK" << std::endl; for (unsigned int i = 0; i < Npoints; i++) { if (DIM == 2) { infile >> _coord[0] >> _coord[1] >> d_dummy; - // std::cout << _coord[0] << " " << _coord[1] << " " << d_dummy - // << std::endl; } else if (DIM == 3) { infile >> _coord[0] >> _coord[1] >> _coord[2]; - // std::cout << _coord[0] << " " << _coord[1] << " " << - // _coord[3] << std::endl; } for (int j = 0; j < DIM; j++) - add_once(value[j], hist[j], _coord[j]); + { + add_once(value[j], hist[j], _coord[j]); + } _node.push_back(_coord); } - std::cout << " done" << std::endl; // create bins - - std::cout << "Determine Body size" << std::endl; for (int j = 0; j < DIM; j++) { std::sort(value[j].begin(), value[j].end()); - // std::cout << "j: " << j << " back(): " << value[j].back() << - // std::endl; min.push_back(value[j][0]); N.push_back(value[j].size()); incr.push_back((value[j].back() - value[j][0]) / (1.0 * N.back())); @@ -265,20 +239,19 @@ namespace PRISMS min[j] -= incr[j]; N[j] += 2; } - std::cout << " Min Coordinate: "; + std::cout << "\tMin Coordinate: "; for (int j = 0; j < DIM; j++) - std::cout << _min[j] << " "; + { + std::cout << _min[j] << " "; + } std::cout << std::endl; - std::cout << " Max Coordinate: "; + std::cout << "\tMax Coordinate: "; for (int j = 0; j < DIM; j++) - std::cout << _max[j] << " "; - std::cout << std::endl; - - std::cout << " done" << std::endl; + { + std::cout << _max[j] << " "; + } - std::cout << "Initialize Bin" << std::endl; _bin = Bin *, Coordinate>(min, incr, N); - std::cout << " done" << std::endl; } } @@ -286,7 +259,6 @@ namespace PRISMS { if (line.size() > 4 && line.substr(0, 5) == "CELLS") { - // std::cout << line << "\n"; ss.clear(); ss.str(line); @@ -309,7 +281,7 @@ namespace PRISMS } bfunc_ptr = _bfunc.back(); - std::cout << "Read CELLS: " << Ncells << std::endl; + std::cout << "\n\tReading CELLS: " << Ncells << std::endl; for (unsigned int i = 0; i < Ncells; i++) { infile >> uli_dummy; @@ -320,9 +292,6 @@ namespace PRISMS infile >> cell_node[j]; } - // std::cout << cell_node[0] << " " << cell_node[1] << " " << - // cell_node[2] << " " << cell_node[3] << std::endl; - // create interpolator if (DIM == 2) { @@ -343,16 +312,12 @@ namespace PRISMS _node); } } - std::cout << " done" << std::endl; // bin interpolators - std::cout << "Bin interpolating functions" << std::endl; - for (unsigned int i = 0; i < _interp.size(); i++) { _bin.add_range(_interp[i], _interp[i]->min(), _interp[i]->max()); } - std::cout << " done max_bin_size: " << _bin.max_size() << std::endl; } else if (line.size() > 9 && line.substr(0, 10) == "CELL_TYPES") { @@ -797,8 +762,6 @@ namespace PRISMS bfunc[i] = 0.0; } node_index[i] = (*bin[i]).node(); - // std::cout << "i: " << i << " bfunc: " << bfunc[i] << " node: " << - // _node[ node_index[i]] << std::endl; } return; } @@ -807,8 +770,6 @@ namespace PRISMS // bfunc[i] = 0.0; // node_index[i] = (*bin[i]).node(); // } - // std::cout << "i: " << i << " bfunc: " << bfunc[i] << " node: " << _node[ - // node_index[i]] << std::endl; } }; @@ -820,7 +781,6 @@ namespace PRISMS std::vector &node_index, int &s) { - // std::cout << "begin Mesh::grad_basis_functions()" << std::endl; std::vector *> &bin = _bin.contents(coord); s = bin.size(); @@ -843,8 +803,6 @@ namespace PRISMS bfunc[i] = 0.0; } node_index[i] = (*bin[i]).node(); - // std::cout << "i: " << i << " bfunc: " << bfunc[i] << " node: " << - // _node[ node_index[i]] << std::endl; } return; } @@ -853,10 +811,7 @@ namespace PRISMS // bfunc[i] = 0.0; // node_index[i] = (*bin[i]).node(); // } - // std::cout << "i: " << i << " bfunc: " << bfunc[i] << " node: " << _node[ - // node_index[i]] << std::endl; } - // std::cout << "finish Mesh::grad_basis_functions()" << std::endl; } // Set 'bfunc' to evaluated hess basis functions at coord, and 's' is the length @@ -890,8 +845,6 @@ namespace PRISMS bfunc[i] = 0.0; } node_index[i] = (*bin[i]).node(); - // std::cout << "i: " << i << " bfunc: " << bfunc[i] << " node: " << - // _node[ node_index[i]] << std::endl; } return; } @@ -900,8 +853,6 @@ namespace PRISMS // bfunc[i] = 0.0; // node_index[i] = (*bin[i]).node(); // } - // std::cout << "i: " << i << " bfunc: " << bfunc[i] << " node: " << _node[ - // node_index[i]] << std::endl; } } @@ -909,8 +860,6 @@ namespace PRISMS void add_once(std::vector &list, std::vector &hist, double val) { - // std::cout << "begin add_once()" << std::endl; - for (unsigned int i = 0; i < list.size(); i++) { if (list[i] == val) @@ -922,8 +871,6 @@ namespace PRISMS list.push_back(val); hist.push_back(1); - - // std::cout << "finish add_once()" << std::endl; } }; } // namespace PRISMS diff --git a/tests/field_input/field_input.cc b/tests/field_input/field_input.cc index e69de29b..6d5ff277 100644 --- a/tests/field_input/field_input.cc +++ b/tests/field_input/field_input.cc @@ -0,0 +1,64 @@ +#include "catch.hpp" + +#include + +TEST_CASE("Unstructured vtk File read-in") +{ + SECTION("2D") + { + using ScalarField = PRISMS::PField; + using Body = PRISMS::Body; + Body body; + + REQUIRE_THROWS_AS(body.read_vtk("invalid"), std::runtime_error); + + body.read_vtk("field_input/test_2D_1degree.vtk"); + + ScalarField n_field = body.find_scalar_field("n"); + ScalarField nx_field = body.find_scalar_field("|nx|"); + + for (auto x : {0.0, 5.0, 10.0}) + { + for (auto y : {0.0, 5.0, 10.0}) + { + double point[2] = {x, y}; + REQUIRE(n_field(point) == x); + REQUIRE(nx_field(point) == 1.0); + } + } + + REQUIRE_THROWS_AS(body.find_scalar_field("invalid"), std::invalid_argument); + } + + SECTION("3D") + { + using ScalarField = PRISMS::PField; + using Body = PRISMS::Body; + Body body; + + REQUIRE_THROWS_AS(body.read_vtk("invalid"), std::runtime_error); + + body.read_vtk("field_input/test_3D_1degree.vtk"); + + ScalarField n_field = body.find_scalar_field("n"); + ScalarField nx_field = body.find_scalar_field("|nx|"); + + for (auto x : {0.0, 5.0, 10.0}) + { + for (auto y : {0.0, 5.0, 10.0}) + { + for (auto z : {0.0, 5.0, 10.0}) + { + double point[3] = {x, y, z}; + REQUIRE(n_field(point) == x); + REQUIRE(nx_field(point) == 1.0); + } + } + } + + REQUIRE_THROWS_AS(body.find_scalar_field("invalid"), std::invalid_argument); + } +} + +TEST_CASE("Rectilinear vtk File read-in") +{} \ No newline at end of file diff --git a/tests/field_input/test_2D_1degree.vtk b/tests/field_input/test_2D_1degree.vtk new file mode 100644 index 00000000..728f097b --- /dev/null +++ b/tests/field_input/test_2D_1degree.vtk @@ -0,0 +1,38 @@ +# vtk DataFile Version 3.0 +#This file was generated by the deal.II library on 2024/12/17 at 7:19:10 +ASCII +DATASET UNSTRUCTURED_GRID + +POINTS 16 double +0 0 0 +5 0 0 +0 5 0 +5 5 0 +5 0 0 +10 0 0 +5 5 0 +10 5 0 +0 5 0 +5 5 0 +0 10 0 +5 10 0 +5 5 0 +10 5 0 +5 10 0 +10 10 0 + +CELLS 4 20 +4 0 1 3 2 +4 4 5 7 6 +4 8 9 11 10 +4 12 13 15 14 + +CELL_TYPES 4 + 9 9 9 9 +POINT_DATA 16 +SCALARS n double 1 +LOOKUP_TABLE default +0 5 0 5 5 10 5 10 0 5 0 5 5 10 5 10 +SCALARS |nx| double 1 +LOOKUP_TABLE default +1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 diff --git a/tests/field_input/test_2D_2degree.vtk b/tests/field_input/test_2D_2degree.vtk new file mode 100644 index 00000000..69499321 --- /dev/null +++ b/tests/field_input/test_2D_2degree.vtk @@ -0,0 +1,70 @@ +# vtk DataFile Version 3.0 +#This file was generated by the deal.II library on 2024/12/17 at 7:20:15 +ASCII +DATASET UNSTRUCTURED_GRID + +POINTS 36 double +0 0 0 +2.5 0 0 +5 0 0 +0 2.5 0 +2.5 2.5 0 +5 2.5 0 +0 5 0 +2.5 5 0 +5 5 0 +5 0 0 +7.5 0 0 +10 0 0 +5 2.5 0 +7.5 2.5 0 +10 2.5 0 +5 5 0 +7.5 5 0 +10 5 0 +0 5 0 +2.5 5 0 +5 5 0 +0 7.5 0 +2.5 7.5 0 +5 7.5 0 +0 10 0 +2.5 10 0 +5 10 0 +5 5 0 +7.5 5 0 +10 5 0 +5 7.5 0 +7.5 7.5 0 +10 7.5 0 +5 10 0 +7.5 10 0 +10 10 0 + +CELLS 16 80 +4 0 1 4 3 +4 1 2 5 4 +4 3 4 7 6 +4 4 5 8 7 +4 9 10 13 12 +4 10 11 14 13 +4 12 13 16 15 +4 13 14 17 16 +4 18 19 22 21 +4 19 20 23 22 +4 21 22 25 24 +4 22 23 26 25 +4 27 28 31 30 +4 28 29 32 31 +4 30 31 34 33 +4 31 32 35 34 + +CELL_TYPES 16 + 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 9 +POINT_DATA 36 +SCALARS n double 1 +LOOKUP_TABLE default +0 2.5 5 0 2.5 5 0 2.5 5 5 7.5 10 5 7.5 10 5 7.5 10 0 2.5 5 0 2.5 5 0 2.5 5 5 7.5 10 5 7.5 10 5 7.5 10 +SCALARS |nx| double 1 +LOOKUP_TABLE default +1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 diff --git a/tests/field_input/test_3D_1degree.vtk b/tests/field_input/test_3D_1degree.vtk new file mode 100644 index 00000000..c6c3ea58 --- /dev/null +++ b/tests/field_input/test_3D_1degree.vtk @@ -0,0 +1,90 @@ +# vtk DataFile Version 3.0 +#This file was generated by the deal.II library on 2024/12/17 at 10:15:23 +ASCII +DATASET UNSTRUCTURED_GRID + +POINTS 64 double +0 0 0 +5 0 0 +0 5 0 +5 5 0 +0 0 5 +5 0 5 +0 5 5 +5 5 5 +5 0 0 +10 0 0 +5 5 0 +10 5 0 +5 0 5 +10 0 5 +5 5 5 +10 5 5 +0 5 0 +5 5 0 +0 10 0 +5 10 0 +0 5 5 +5 5 5 +0 10 5 +5 10 5 +5 5 0 +10 5 0 +5 10 0 +10 10 0 +5 5 5 +10 5 5 +5 10 5 +10 10 5 +0 0 5 +5 0 5 +0 5 5 +5 5 5 +0 0 10 +5 0 10 +0 5 10 +5 5 10 +5 0 5 +10 0 5 +5 5 5 +10 5 5 +5 0 10 +10 0 10 +5 5 10 +10 5 10 +0 5 5 +5 5 5 +0 10 5 +5 10 5 +0 5 10 +5 5 10 +0 10 10 +5 10 10 +5 5 5 +10 5 5 +5 10 5 +10 10 5 +5 5 10 +10 5 10 +5 10 10 +10 10 10 + +CELLS 8 72 +8 0 1 3 2 4 5 7 6 +8 8 9 11 10 12 13 15 14 +8 16 17 19 18 20 21 23 22 +8 24 25 27 26 28 29 31 30 +8 32 33 35 34 36 37 39 38 +8 40 41 43 42 44 45 47 46 +8 48 49 51 50 52 53 55 54 +8 56 57 59 58 60 61 63 62 + +CELL_TYPES 8 + 12 12 12 12 12 12 12 12 +POINT_DATA 64 +SCALARS n double 1 +LOOKUP_TABLE default +0 5 0 5 0 5 0 5 5 10 5 10 5 10 5 10 0 5 0 5 0 5 0 5 5 10 5 10 5 10 5 10 0 5 0 5 0 5 0 5 5 10 5 10 5 10 5 10 0 5 0 5 0 5 0 5 5 10 5 10 5 10 5 10 +SCALARS |nx| double 1 +LOOKUP_TABLE default +1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1 1