Loading [MathJax]/extensions/TeX/newcommand.js
 deal.II version GIT relicensing-2848-g5241f990fb 2025-03-16 19:30:00+00:00
\newcommand{\dealvcentcolon}{\mathrel{\mathop{:}}} \newcommand{\dealcoloneq}{\dealvcentcolon\mathrel{\mkern-1.2mu}=} \newcommand{\jump}[1]{\left[\!\left[ #1 \right]\!\right]} \newcommand{\average}[1]{\left\{\!\left\{ #1 \right\}\!\right\}}
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages Concepts
grid_in.h
Go to the documentation of this file.
1// ------------------------------------------------------------------------
2//
3// SPDX-License-Identifier: LGPL-2.1-or-later
4// Copyright (C) 1999 - 2024 by the deal.II authors
5//
6// This file is part of the deal.II library.
7//
8// Part of the source code is dual licensed under Apache-2.0 WITH
9// LLVM-exception OR LGPL-2.1-or-later. Detailed license information
10// governing the source code and code contributions can be found in
11// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II.
12//
13// ------------------------------------------------------------------------
14
15#ifndef dealii_grid_in_h
16#define dealii_grid_in_h
17
18
19#include <deal.II/base/config.h>
20
23#include <deal.II/base/point.h>
24
25#include <deal.II/lac/vector.h>
26
27#include <iostream>
28#include <map>
29#include <string>
30#include <vector>
31
33
34// Forward declarations
35#ifndef DOXYGEN
36template <int dim, int spacedim>
38class Triangulation;
39template <int dim>
40struct CellData;
41#endif
42
314template <int dim, int spacedim = dim>
316{
317public:
349
353 GridIn();
354
359
363 void
365
370 void
371 read(std::istream &in, Format format = Default);
372
377 void
378 read(const std::string &in, Format format = Default);
379
420 void
421 read_vtk(std::istream &in);
422
437 void
438 read_vtu(std::istream &in);
439
440
448 void
449 read_unv(std::istream &in);
450
460 void
461 read_ucd(std::istream &in,
462 const bool apply_all_indicators_to_manifolds = false);
463
503 void
504 read_abaqus(std::istream &in,
505 const bool apply_all_indicators_to_manifolds = false);
506
510 void
511 read_dbmesh(std::istream &in);
512
516 void
517 read_xda(std::istream &in);
518
526 void
527 read_msh(std::istream &in);
528
529#ifdef DEAL_II_GMSH_WITH_API
587 void
588 read_msh(const std::string &filename);
589#endif
590
627 void
628 read_comsol_mphtxt(std::istream &in);
629
634 void
635 read_tecplot(std::istream &in);
636
671 void
672 read_assimp(const std::string &filename,
673 const unsigned int mesh_index = numbers::invalid_unsigned_int,
674 const bool remove_duplicates = true,
675 const double tol = 1e-12,
676 const bool ignore_unsupported_element_types = true);
677
687 {
692 std::vector<std::vector<int>> id_to_sideset_ids;
693 };
694
752 read_exodusii(const std::string &filename,
753 const bool apply_all_indicators_to_manifolds = false);
754
758 static std::string
759 default_suffix(const Format format);
760
764 static Format
765 parse_format(const std::string &format_name);
766
772 static std::string
774
791 const std::map<std::string, Vector<double>> &
792 get_cell_data() const;
793
798 int,
799 << "The section type <" << arg1 << "> in an UNV "
800 << "input file is not implemented.");
801
806 int,
807 << "The element type <" << arg1 << "> in an UNV "
808 << "input file is not implemented.");
809
814 std::string,
815 << "The identifier <" << arg1 << "> as name of a "
816 << "part in an UCD input file is unknown or the "
817 << "respective input routine is not implemented."
818 << "(Maybe the space dimension of triangulation and "
819 << "input file do not match?");
824 "No Triangulation has been attached to this GridIn object "
825 "so that nothing can be filled during any read function "
826 "calls. Please pass a reference to the Triangulation tria "
827 "to be filled in the constructor GridIn(tria) or attach "
828 "it with the function call GridIn::attach_triangulation().");
834 int,
835 int,
836 << "While creating cell " << arg1
837 << ", you are referencing a vertex with index " << arg2
838 << " but no vertex with this index has been described in the input file.");
844 int,
845 int,
846 int,
847 << "While creating cell " << arg1 << " (which is numbered as " << arg2
848 << " in the input file), you are referencing a vertex with index " << arg3
849 << " but no vertex with this index has been described in the input file.");
858 std::string,
859 << "The string <" << arg1
860 << "> is not recognized at the present"
861 << " position of a DB Mesh file.");
862
868 int,
869 << "The specified dimension " << arg1
870 << " is not the same as that of the triangulation to be created.");
871
873 std::string,
874 << "The string <" << arg1
875 << "> is not recognized at the present"
876 << " position of a Gmsh Mesh file.");
877
879 int,
880 << "The Element Identifier <" << arg1 << "> is not "
881 << "supported in the deal.II library when "
882 << "reading meshes in " << dim << " dimensions.\n"
883 << "Supported elements are: \n"
884 << "ELM-TYPE\n"
885 << "1 Line (2 nodes, 1 edge).\n"
886 << "2 Triangle (3 nodes, 3 edges).\n"
887 << "3 Quadrilateral (4 nodes, 4 edges).\n"
888 << "4 Tetrahedron (4 nodes, 6 edges, 4 faces) when in 3d.\n"
889 << "5 Hexahedron (8 nodes, 12 edges, 6 faces) when in 3d.\n"
890 << "15 Point (1 node, ignored when read).");
891
892
895 unsigned int,
896 unsigned int,
897 "While reading a gmsh file, the reader function did not find "
898 "any cells. This sometimes happens if the file only contains a "
899 "surface mesh, but not a volume mesh."
900 "\n\n"
901 "The reader function did find " +
902 std::to_string(arg1) + " lines and " + std::to_string(arg2) +
903 " facets (surface triangles or quadrilaterals).");
904
905protected:
910
931 static void
932 debug_output_grid(const std::vector<CellData<dim>> &cells,
933 const std::vector<Point<spacedim>> &vertices,
934 std::ostream &out);
935
936private:
941 static void
942 skip_empty_lines(std::istream &in);
943
951 static void
952 skip_comment_lines(std::istream &in, const char comment_start);
953
961 static void
962 parse_tecplot_header(std::string &header,
963 std::vector<unsigned int> &tecplot2deal,
964 unsigned int &n_vars,
965 unsigned int &n_vertices,
966 unsigned int &n_cells,
967 std::vector<unsigned int> &IJK,
968 bool &structured,
969 bool &blocked);
970
975
985 std::map<std::string, Vector<double>> cell_data;
986};
987
988/* -------------- declaration of explicit specializations ------------- */
989
990#ifndef DOXYGEN
991
992template <>
993void
994GridIn<2>::debug_output_grid(const std::vector<CellData<2>> &cells,
995 const std::vector<Point<2>> &vertices,
996 std::ostream &out);
997
998
999template <>
1000void
1001GridIn<2, 3>::debug_output_grid(const std::vector<CellData<2>> &cells,
1002 const std::vector<Point<3>> &vertices,
1003 std::ostream &out);
1004template <>
1005void
1006GridIn<3>::debug_output_grid(const std::vector<CellData<3>> &cells,
1007 const std::vector<Point<3>> &vertices,
1008 std::ostream &out);
1009#endif // DOXYGEN
1010
1012
1013#endif
@ Default
Use GridIn::default_format stored in this object.
Definition grid_in.h:325
@ vtu
Use read_vtu()
Definition grid_in.h:343
@ xda
Use read_xda()
Definition grid_in.h:335
@ tecplot
Use read_tecplot()
Definition grid_in.h:339
@ msh
Use read_msh()
Definition grid_in.h:337
@ ucd
Use read_ucd()
Definition grid_in.h:329
@ dbmesh
Use read_dbmesh()
Definition grid_in.h:333
@ abaqus
Use read_abaqus()
Definition grid_in.h:331
@ assimp
Use read_assimp()
Definition grid_in.h:345
@ unv
Use read_unv()
Definition grid_in.h:327
@ exodusii
Use read_exodusii()
Definition grid_in.h:347
@ vtk
Use read_vtk()
Definition grid_in.h:341
void read_vtk(std::istream &in)
Definition grid_in.cc:162
GridIn()
Definition grid_in.cc:136
static void skip_empty_lines(std::istream &in)
Definition grid_in.cc:4045
void read_assimp(const std::string &filename, const unsigned int mesh_index=numbers::invalid_unsigned_int, const bool remove_duplicates=true, const double tol=1e-12, const bool ignore_unsupported_element_types=true)
Definition grid_in.cc:3509
std::map< std::string, Vector< double > > cell_data
Definition grid_in.h:985
void read_abaqus(std::istream &in, const bool apply_all_indicators_to_manifolds=false)
Definition grid_in.cc:1246
static std::string default_suffix(const Format format)
Definition grid_in.cc:4349
void read_xda(std::istream &in)
Definition grid_in.cc:1463
void read_comsol_mphtxt(std::istream &in)
Definition grid_in.cc:1527
static void skip_comment_lines(std::istream &in, const char comment_start)
Definition grid_in.cc:4075
void attach_triangulation(Triangulation< dim, spacedim > &tria)
Definition grid_in.cc:153
const std::map< std::string, Vector< double > > & get_cell_data() const
Definition grid_in.cc:658
void read_msh(std::istream &in)
Definition grid_in.cc:2083
static Format parse_format(const std::string &format_name)
Definition grid_in.cc:4384
void read_vtu(std::istream &in)
Definition grid_in.cc:665
void read_tecplot(std::istream &in)
Definition grid_in.cc:3500
Format default_format
Definition grid_in.h:974
ExodusIIData read_exodusii(const std::string &filename, const bool apply_all_indicators_to_manifolds=false)
Definition grid_in.cc:3880
void read_dbmesh(std::istream &in)
Definition grid_in.cc:1302
void read_ucd(std::istream &in, const bool apply_all_indicators_to_manifolds=false)
Definition grid_in.cc:992
ObserverPointer< Triangulation< dim, spacedim >, GridIn< dim, spacedim > > tria
Definition grid_in.h:909
void read(std::istream &in, Format format=Default)
Definition grid_in.cc:4282
static void debug_output_grid(const std::vector< CellData< dim > > &cells, const std::vector< Point< spacedim > > &vertices, std::ostream &out)
Definition grid_in.cc:4101
static std::string get_format_names()
Definition grid_in.cc:4441
void read_unv(std::istream &in)
Definition grid_in.cc:692
static void parse_tecplot_header(std::string &header, std::vector< unsigned int > &tecplot2deal, unsigned int &n_vars, unsigned int &n_vertices, unsigned int &n_cells, std::vector< unsigned int > &IJK, bool &structured, bool &blocked)
Definition grid_in.cc:3032
Definition point.h:113
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:40
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:245
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:41
#define DeclException0(Exception0)
static ::ExceptionBase & ExcInvalidDBMeshFormat()
static ::ExceptionBase & ExcUnknownIdentifier(std::string arg1)
static ::ExceptionBase & ExcInvalidGMSHInput(std::string arg1)
static ::ExceptionBase & ExcInvalidDBMESHInput(std::string arg1)
static ::ExceptionBase & ExcDBMESHWrongDimension(int arg1)
static ::ExceptionBase & ExcInvalidVertexIndex(int arg1, int arg2)
static ::ExceptionBase & ExcNoTriangulationSelected()
#define DeclException2(Exception2, type1, type2, outsequence)
static ::ExceptionBase & ExcUnknownElementType(int arg1)
static ::ExceptionBase & ExcGmshNoCellInformation(unsigned int arg1, unsigned int arg2)
static ::ExceptionBase & ExcInvalidVertexIndexGmsh(int arg1, int arg2, int arg3)
#define DeclExceptionMsg(Exception, defaulttext)
static ::ExceptionBase & ExcUnknownSectionType(int arg1)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
#define DeclException1(Exception1, type1, outsequence)
static ::ExceptionBase & ExcGmshUnsupportedGeometry(int arg1)
constexpr unsigned int invalid_unsigned_int
Definition types.h:232
std::vector< std::vector< int > > id_to_sideset_ids
Definition grid_in.h:692