From fa053ee3527ee036bbefa98d6ea32fff61abb3d7 Mon Sep 17 00:00:00 2001 From: leicht Date: Thu, 14 Dec 2006 10:58:11 +0000 Subject: [PATCH] add functionality to read in tecplot (grid) files (in 2d) git-svn-id: https://svn.dealii.org/trunk@14238 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/base/source/path_search.cc | 10 +- deal.II/deal.II/include/grid/grid_in.h | 37 +- deal.II/deal.II/source/grid/grid_in.cc | 475 ++++++++++++++++++++++++- deal.II/doc/news/changes.html | 9 + 4 files changed, 526 insertions(+), 5 deletions(-) diff --git a/deal.II/base/source/path_search.cc b/deal.II/base/source/path_search.cc index 5ba5cb88a4..097dbb1533 100644 --- a/deal.II/base/source/path_search.cc +++ b/deal.II/base/source/path_search.cc @@ -13,6 +13,8 @@ #include #include +#include +#include #include #include @@ -40,12 +42,18 @@ PathSearch::initialize_classes() v.push_back(std::string(".prm")); suffix_lists.insert(map_type(std::string("PARAMETER"), v)); -//TODO[GK]: instead of listing by hand here, query the GridIn class which formats it can presently read +//TODO[GK]: instead of listing by hand here, query the GridIn class which formats it can presently read. +//[Tobias Leicht]: the problem with this approach is that we would could not link with a base-lib alone +// (undefined reference problem). we do that in the tests, however. v.clear(); v.push_back(empty); v.push_back(std::string(".inp")); v.push_back(std::string(".xda")); v.push_back(std::string(".dbmesh")); + v.push_back(std::string(".dat")); + v.push_back(std::string(".plt")); + v.push_back(std::string(".nc")); + v.push_back(std::string(".msh")); suffix_lists.insert(map_type(std::string("MESH"), v)); } diff --git a/deal.II/deal.II/include/grid/grid_in.h b/deal.II/deal.II/include/grid/grid_in.h index bf60945689..4c90207c41 100644 --- a/deal.II/deal.II/include/grid/grid_in.h +++ b/deal.II/deal.II/include/grid/grid_in.h @@ -123,7 +123,11 @@ class SubCellData; * mesh generator (see http://www.geuz.org/gmsh/ ). The documentation * in the @p GMSH manual explains how to generate meshes compatible with * the deal.II library (i.e. quads rather than triangles). - * + * + *
  • Tecplot format: this format is used by @p TECPLOT and often + * serves as a basis for data exchange between different applications. Note, + * that currently only the ASCII format is supported, binary data cannot be + * read. * * *

    Structure of input grid data. The GridReordering class

    @@ -208,7 +212,9 @@ class GridIn /// Use read_msh() msh, /// Use read_netcdf() - netcdf + netcdf, + /// Use read_tecplot() + tecplot }; /** @@ -277,6 +283,14 @@ class GridIn */ void read_netcdf (const std::string &filename); + /** + * Read grid data from a file containing + * tecplot ASCII data. This also works in + * the absence of any tecplot + * installation. + */ + void read_tecplot (std::istream &in); + /** * Returns the standard suffix * for a file in this format. @@ -384,6 +398,25 @@ class GridIn static void skip_comment_lines (std::istream &in, const char comment_start); + /** + * This function does the nasty work (due + * to very lax conventions and different + * versions of the tecplot format) of + * extracting the important parameters from + * a tecplot header, contained in the + * string @p header. The other variables + * are output variables, their value has no + * influence on the function execution.. + */ + static void parse_tecplot_header(std::string &header, + unsigned int (&tecplot2deal)[dim], + unsigned int &n_vars, + unsigned int &n_vertices, + unsigned int &n_cells, + unsigned int (&IJK)[dim], + bool &structured, + bool &blocked); + /** * This function can write the * raw cell data objects created diff --git a/deal.II/deal.II/source/grid/grid_in.cc b/deal.II/deal.II/source/grid/grid_in.cc index d7d7edc346..dff3062b06 100644 --- a/deal.II/deal.II/source/grid/grid_in.cc +++ b/deal.II/deal.II/source/grid/grid_in.cc @@ -12,6 +12,7 @@ //--------------------------------------------------------------------------- #include +#include #include #include @@ -21,6 +22,7 @@ #include #include #include +#include #ifdef HAVE_LIBNETCDF #include @@ -28,7 +30,6 @@ DEAL_II_NAMESPACE_OPEN - template GridIn::GridIn () : tria(0), default_format(ucd) @@ -1224,6 +1225,449 @@ void GridIn<3>::read_netcdf (const std::string &filename) #endif +template +void GridIn::parse_tecplot_header(std::string &header, + unsigned int (&tecplot2deal)[dim], + unsigned int &n_vars, + unsigned int &n_vertices, + unsigned int &n_cells, + unsigned int (&IJK)[dim], + bool &structured, + bool &blocked) +{ + // initialize the output variables + n_vars=0; + n_vertices=0; + n_cells=0; + switch(dim) + { + case 3: + IJK[2]=0; + case 2: + IJK[1]=0; + case 1: + IJK[0]=0; + } + structured=true; + blocked=false; + + // convert the string to upper case + transform(header.begin(),header.end(),header.begin(),::toupper); + + // replace all tabs, commas, newlines by + // whitespaces + std::replace(header.begin(),header.end(),'\t',' '); + std::replace(header.begin(),header.end(),',',' '); + std::replace(header.begin(),header.end(),'\n',' '); + + // now remove whitespace in front of and + // after '=' + std::string::size_type pos=header.find("="); + + while(pos!=static_cast(std::string::npos)) + if(header[pos+1]==' ') + header.erase(pos+1,1); + else if (header[pos-1]==' ') + { + header.erase(pos-1,1); + --pos; + } + else + pos=header.find("=",++pos); + + // split the string into individual entries + std::vector entries=Utilities::break_text_into_lines(header,1,' '); + + // now go through the list and try to extract + for (unsigned int i=0; i1, + ExcMessage("Tecplot file contains Y data, which is not possible for 1d plot")); + tecplot2deal[1]=n_vars; + } + else if (entries[i]=="\"Z\"") + switch (dim) + { + case 1: + AssertThrow(false, + ExcMessage("Tecplot file contains Y data, which is not possible for 1d plot")); + break; + case 2: + // we assume, that z + // contains the data + // which is intended + // as y in deal + tecplot2deal[1]=n_vars; + break; + case 3: + tecplot2deal[2]=n_vars; + break; + } + ++n_vars; + ++i; + } + // set i back, so that the next + // string is treated correctly + --i; + + AssertThrow(n_vars>=dim, + ExcMessage("Tecplot file must contain at least one variable for each dimension")); + for (unsigned int i=1; i0, + ExcMessage("Tecplot file must contain at least one variable for each dimension.")); + } + else if (Utilities::match_at_string_start(entries[i],"ZONETYPE=ORDERED")) + structured=true; + else if (Utilities::match_at_string_start(entries[i],"ZONETYPE=FELINESEG") && dim==1) + structured=false; + else if (Utilities::match_at_string_start(entries[i],"ZONETYPE=FEQUADRILATERAL") && dim==2) + structured=false; + else if (Utilities::match_at_string_start(entries[i],"ZONETYPE=FEBRICK") && dim==3) + structured=false; + else if (Utilities::match_at_string_start(entries[i],"ZONETYPE=")) + // unsupported ZONETYPE + { + AssertThrow(false,ExcMessage("The tecplot file contains an unsupported ZONETYPE.")); + } + else if (Utilities::match_at_string_start(entries[i],"DATAPACKING=POINT")) + blocked=false; + else if (Utilities::match_at_string_start(entries[i],"DATAPACKING=BLOCK")) + blocked=true; + else if (Utilities::match_at_string_start(entries[i],"F=POINT")) + { + structured=true; + blocked=false; + } + else if (Utilities::match_at_string_start(entries[i],"F=BLOCK")) + { + structured=true; + blocked=true; + } + else if (Utilities::match_at_string_start(entries[i],"F=FEPOINT")) + { + structured=false; + blocked=false; + } + else if (Utilities::match_at_string_start(entries[i],"F=FEBLOCK")) + { + structured=false; + blocked=true; + } + else if (Utilities::match_at_string_start(entries[i],"ET=QUADRILATERAL") && dim==2) + structured=false; + else if (Utilities::match_at_string_start(entries[i],"ET=BRICK") && dim==3) + structured=false; + else if (Utilities::match_at_string_start(entries[i],"ET=")) + // unsupported ElementType + { + AssertThrow(false,ExcMessage("The tecplot file contains an unsupported ElementType.")); + } + else if (Utilities::match_at_string_start(entries[i],"I=")) + IJK[0]=Utilities::get_integer_at_position(entries[i],2).first; + else if (Utilities::match_at_string_start(entries[i],"J=")) + { + IJK[1]=Utilities::get_integer_at_position(entries[i],2).first; + AssertThrow(dim>1 || IJK[1]==1, + ExcMessage("Parameter 'J=' found in tecplot, although this is only possible for dimensions greater than 1.")); + } + else if (Utilities::match_at_string_start(entries[i],"K=")) + { + IJK[2]=Utilities::get_integer_at_position(entries[i],2).first; + AssertThrow(dim>2 || IJK[2]==1, + ExcMessage("Parameter 'K=' found in tecplot, although this is only possible for dimensions greater than 2.")); + } + else if (Utilities::match_at_string_start(entries[i],"N=")) + n_vertices=Utilities::get_integer_at_position(entries[i],2).first; + else if (Utilities::match_at_string_start(entries[i],"E=")) + n_cells=Utilities::get_integer_at_position(entries[i],2).first; + } + + // now we have read all the fields we are + // interested in. do some checks and + // calculate the variables + if (structured) + { + n_vertices=1; + n_cells=1; + for (unsigned int d=0; d0, + ExcMessage("Tecplot file does not contain a complete and consistent set of parameters")); + n_vertices*=IJK[d]; + n_cells*=(IJK[d]-1); + } + } + else + { + AssertThrow(n_vertices>0, + ExcMessage("Tecplot file does not contain a complete and consistent set of parameters")); + if (n_cells==0) + // this means an error, although + // tecplot itself accepts entries like + // 'J=20' instead of 'E=20'. therefore, + // take the max of IJK + n_cells=*std::max_element(IJK,IJK+dim); + AssertThrow(n_cells>0, + ExcMessage("Tecplot file does not contain a complete and consistent set of parameters")); + } +} + + + +#if deal_II_dimension == 2 + +template <> +void GridIn<2>::read_tecplot (std::istream &in) +{ + const unsigned int dim=2; + Assert (tria != 0, ExcNoTriangulationSelected()); + AssertThrow (in, ExcIO()); + + // skip comments at start of file + skip_comment_lines (in, '#'); + + // some strings for parsing the header + std::string line, header; + + // first, concatenate all header lines + // create a searchstring with almost all + // letters. exclude e and E from the letters + // to search, as they might appear in + // exponential notation + std::string letters ="abcdfghijklmnopqrstuvwxyzABCDFGHIJKLMNOPQRSTUVWXYZ"; + + getline(in,line); + while(line.find_first_of(letters)!=std::string::npos) + { + header+=" "+line; + getline(in,line); + } + + // now create some variables holding + // important information on the mesh, get + // this information from the header string + unsigned int tecplot2deal[dim], + n_vars, + n_vertices, + n_cells, + IJK[dim]; + bool structured, + blocked; + + parse_tecplot_header(header, + tecplot2deal,n_vars,n_vertices,n_cells,IJK, + structured,blocked); + + // reserve space for vertices. note, that in + // tecplot vertices are ordered beginning + // with 1, whereas in deal all indices start + // with 0. in order not to use -1 for all the + // connectivity information, a 0th vertex + // (unused) is inserted at the origin. + std::vector > vertices(n_vertices+1); + vertices[0]=Point(); + // reserve space for cells + std::vector > cells(n_cells); + SubCellData subcelldata; + + if (blocked) + { + // blocked data format. first we get all + // the values of the first variable for + // all points, after that we get all + // values for the second variable and so + // on. + + // dummy variable to read in all the info + // we do not want to use + double dummy; + // which is the first index to read in + // the loop (see below) + unsigned int next_index=0; + + // note, that we have already read the + // first line containing the first variable + if (tecplot2deal[0]==0) + { + // we need the information in this + // line, so extract it + std::vector first_var=Utilities::break_text_into_lines(line,1); + char *endptr; + for (unsigned int i=1; i>vertices[j](next_index); + // now we got all values of the first + // variable, so increase the counter + next_index=1; + } + + // main loop over all variables + for (unsigned int i=1; i>vertices[j](next_index); + ++next_index; + } + else + { + // we do not need this line, read + // it in and discard it + for (unsigned int j=1; j>dummy; + } + } + Assert(next_index==dim, ExcInternalError()); + } + else + { + // the data is not blocked, so we get all + // the variables for one point, then the + // next and so on. create a vector to + // hold these components + std::vector vars(n_vars); + + // now fill the first vertex. note, that we + // have already read the first line + // containing the first vertex + std::vector first_vertex=Utilities::break_text_into_lines(line,1); + char *endptr; + for (unsigned int d=0; d>vars[i]; + // fill the vertex + // coordinates. respect the position + // of coordinates in the list of + // variables + for (unsigned int i=0; i boundary_vertices(2*I+2*J-4); + unsigned int k=0; + for (unsigned int i=1;i::vertices_per_cell; ++j) + in>>cells[i].vertices[j]; + } + // do some clean-up on vertices + GridTools::delete_unused_vertices (vertices, cells, subcelldata); + } + + // check that no forbidden arrays are + // used. as we do not read in any + // subcelldata, nothing should happen here. + Assert (subcelldata.check_consistency(dim), ExcInternalError()); + AssertThrow (in, ExcIO()); + + // do some cleanup on cells + GridReordering::invert_all_cells_of_negative_grid (vertices, cells); + GridReordering::reorder_cells (cells); + tria->create_triangulation_compatibility (vertices, cells, subcelldata); +} + +#endif + +template +void GridIn::read_tecplot(std::istream &) +{ + Assert(false, ExcNotImplemented()); +} template @@ -1269,6 +1713,9 @@ void GridIn::skip_comment_lines (std::istream &in, // put back first character of // first non-comment line in.putback (c); + + // at last: skip additional empty lines, if present + skip_empty_lines(in); } @@ -1494,6 +1941,10 @@ void GridIn::read (std::istream& in, "functions, instead.")); return; + case tecplot: + read_tecplot (in); + return; + case Default: break; } @@ -1518,6 +1969,8 @@ GridIn::default_suffix (const Format format) return ".xda"; case netcdf: return ".nc"; + case tecplot: + return ".dat"; default: Assert (false, ExcNotImplemented()); return ".unknown_format"; @@ -1551,6 +2004,24 @@ GridIn::parse_format (const std::string &format_name) if (format_name == "nc") return netcdf; + if (format_name == "tecplot") + return tecplot; + + if (format_name == "dat") + return tecplot; + + if (format_name == "plt") + // Actually, this is the extension for the + // tecplot binary format, which we do not + // support right now. However, some people + // tend to create tecplot ascii files with + // the extension 'plt' instead of + // 'dat'. Thus, include this extension + // here. If it actually is a binary file, + // the read_tecplot() function will fail + // and throw an exception, anyway. + return tecplot; + AssertThrow (false, ExcInvalidState ()); // return something weird return Format(Default); @@ -1561,7 +2032,7 @@ GridIn::parse_format (const std::string &format_name) template std::string GridIn::get_format_names () { - return "dbmesh|msh|ucd|xda|netcdf"; + return "dbmesh|msh|ucd|xda|netcdf|tecplot"; } diff --git a/deal.II/doc/news/changes.html b/deal.II/doc/news/changes.html index 19cc56e4c1..8752cc5f9c 100644 --- a/deal.II/doc/news/changes.html +++ b/deal.II/doc/news/changes.html @@ -928,6 +928,15 @@ inconvenience this causes.
      +
    1. Extended: The GridIn class can now read in tecplot files + in ASCII format (block and point format, ordered and unstructured grids, + format specifiers acccording to Tecplot 10 and younger versions). At the + moment the implementation is restricted to 2d grids but can easily be + extended to 3d as well. +
      + (Tobias Leicht 2006/12/14) +

      +
    2. Extended: So far, the GridReordering::invert_all_cells_of_negative_grid function did nothing in 2d. Now, it inverts cells from clockwise to -- 2.39.5