Reference documentation for deal.II version 9.5.0
\(\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\}}\)
Loading...
Searching...
No Matches
grid_in.h
Go to the documentation of this file.
1// ---------------------------------------------------------------------
2//
3// Copyright (C) 1999 - 2023 by the deal.II authors
4//
5// This file is part of the deal.II library.
6//
7// The deal.II library is free software; you can use it, redistribute
8// it, and/or modify it under the terms of the GNU Lesser General
9// Public License as published by the Free Software Foundation; either
10// version 2.1 of the License, or (at your option) any later version.
11// The full text of the license can be found in the file LICENSE.md at
12// the top level directory of deal.II.
13//
14// ---------------------------------------------------------------------
15
16#ifndef dealii_grid_in_h
17#define dealii_grid_in_h
18
19
20#include <deal.II/base/config.h>
21
23#include <deal.II/base/point.h>
25
26#include <iostream>
27#include <string>
28#include <vector>
29
31
32// Forward declarations
33#ifndef DOXYGEN
34template <int dim, int spacedim>
36class Triangulation;
37template <int dim>
38struct CellData;
39#endif
40
311template <int dim, int spacedim = dim>
313{
314public:
320 {
345 };
346
350 GridIn();
351
356
360 void
362
367 void
368 read(std::istream &in, Format format = Default);
369
375 void
376 read(const std::string &in, Format format = Default);
377
418 void
419 read_vtk(std::istream &in);
420
435 void
436 read_vtu(std::istream &in);
437
438
446 void
447 read_unv(std::istream &in);
448
458 void
459 read_ucd(std::istream &in,
460 const bool apply_all_indicators_to_manifolds = false);
461
501 void
502 read_abaqus(std::istream &in,
503 const bool apply_all_indicators_to_manifolds = false);
504
508 void
509 read_dbmesh(std::istream &in);
510
514 void
515 read_xda(std::istream &in);
516
524 void
525 read_msh(std::istream &in);
526
527#ifdef DEAL_II_GMSH_WITH_API
585 void
586 read_msh(const std::string &filename);
587#endif
588
625 void
626 read_comsol_mphtxt(std::istream &in);
627
632 void
633 read_tecplot(std::istream &in);
634
669 void
670 read_assimp(const std::string &filename,
671 const unsigned int mesh_index = numbers::invalid_unsigned_int,
672 const bool remove_duplicates = true,
673 const double tol = 1e-12,
674 const bool ignore_unsupported_element_types = true);
675
685 {
690 std::vector<std::vector<int>> id_to_sideset_ids;
691 };
692
750 read_exodusii(const std::string &filename,
751 const bool apply_all_indicators_to_manifolds = false);
752
756 static std::string
757 default_suffix(const Format format);
758
762 static Format
763 parse_format(const std::string &format_name);
764
770 static std::string
772
777 int,
778 << "The section type <" << arg1 << "> in an UNV "
779 << "input file is not implemented.");
780
785 int,
786 << "The element type <" << arg1 << "> in an UNV "
787 << "input file is not implemented.");
788
793 std::string,
794 << "The identifier <" << arg1 << "> as name of a "
795 << "part in an UCD input file is unknown or the "
796 << "respective input routine is not implemented."
797 << "(Maybe the space dimension of triangulation and "
798 << "input file do not match?");
803 "No Triangulation has been attached to this GridIn object "
804 "so that nothing can be filled during any read function "
805 "calls. Please pass a reference to the Triangulation tria "
806 "to be filled in the constructor GridIn(tria) or attach "
807 "it with the function call GridIn::attach_triangulation().");
813 int,
814 int,
815 << "While creating cell " << arg1
816 << ", you are referencing a vertex with index " << arg2
817 << " but no vertex with this index has been described in the input file.");
823 int,
824 int,
825 int,
826 << "While creating cell " << arg1 << " (which is numbered as " << arg2
827 << " in the input file), you are referencing a vertex with index " << arg3
828 << " but no vertex with this index has been described in the input file.");
837 std::string,
838 << "The string <" << arg1
839 << "> is not recognized at the present"
840 << " position of a DB Mesh file.");
841
847 int,
848 << "The specified dimension " << arg1
849 << " is not the same as that of the triangulation to be created.");
850
852 std::string,
853 << "The string <" << arg1
854 << "> is not recognized at the present"
855 << " position of a Gmsh Mesh file.");
856
858 int,
859 << "The Element Identifier <" << arg1 << "> is not "
860 << "supported in the deal.II library when "
861 << "reading meshes in " << dim << " dimensions.\n"
862 << "Supported elements are: \n"
863 << "ELM-TYPE\n"
864 << "1 Line (2 nodes, 1 edge).\n"
865 << "3 Quadrilateral (4 nodes, 4 edges).\n"
866 << "5 Hexahedron (8 nodes, 12 edges, 6 faces) when in 3d.\n"
867 << "15 Point (1 node, ignored when read)");
868
869
872 unsigned int,
873 unsigned int,
874 "While reading a gmsh file, the reader function did not find "
875 "any cells. This sometimes happens if the file only contains a "
876 "surface mesh, but not a volume mesh."
877 "\n\n"
878 "The reader function did find " +
879 std::to_string(arg1) + " lines and " + std::to_string(arg2) +
880 " facets (surface triangles or quadrilaterals).");
881
882protected:
887
908 static void
909 debug_output_grid(const std::vector<CellData<dim>> & cells,
910 const std::vector<Point<spacedim>> &vertices,
911 std::ostream & out);
912
913private:
918 static void
919 skip_empty_lines(std::istream &in);
920
928 static void
929 skip_comment_lines(std::istream &in, const char comment_start);
930
938 static void
939 parse_tecplot_header(std::string & header,
940 std::vector<unsigned int> &tecplot2deal,
941 unsigned int & n_vars,
942 unsigned int & n_vertices,
943 unsigned int & n_cells,
944 std::vector<unsigned int> &IJK,
945 bool & structured,
946 bool & blocked);
947
952};
953
954/* -------------- declaration of explicit specializations ------------- */
955
956#ifndef DOXYGEN
957
958template <>
959void
960GridIn<2>::debug_output_grid(const std::vector<CellData<2>> &cells,
961 const std::vector<Point<2>> & vertices,
962 std::ostream & out);
963
964
965template <>
966void
967GridIn<2, 3>::debug_output_grid(const std::vector<CellData<2>> &cells,
968 const std::vector<Point<3>> & vertices,
969 std::ostream & out);
970template <>
971void
972GridIn<3>::debug_output_grid(const std::vector<CellData<3>> &cells,
973 const std::vector<Point<3>> & vertices,
974 std::ostream & out);
975#endif // DOXYGEN
976
978
979#endif
@ Default
Use GridIn::default_format stored in this object.
Definition grid_in.h:322
@ vtu
Use read_vtu()
Definition grid_in.h:340
@ xda
Use read_xda()
Definition grid_in.h:332
@ tecplot
Use read_tecplot()
Definition grid_in.h:336
@ msh
Use read_msh()
Definition grid_in.h:334
@ ucd
Use read_ucd()
Definition grid_in.h:326
@ dbmesh
Use read_dbmesh()
Definition grid_in.h:330
@ abaqus
Use read_abaqus()
Definition grid_in.h:328
@ assimp
Use read_assimp()
Definition grid_in.h:342
@ unv
Use read_unv()
Definition grid_in.h:324
@ exodusii
Use read_exodusii()
Definition grid_in.h:344
@ vtk
Use read_vtk()
Definition grid_in.h:338
SmartPointer< Triangulation< dim, spacedim >, GridIn< dim, spacedim > > tria
Definition grid_in.h:886
void read_vtk(std::istream &in)
Definition grid_in.cc:165
GridIn()
Definition grid_in.cc:139
static void skip_empty_lines(std::istream &in)
Definition grid_in.cc:3937
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:3402
void read_abaqus(std::istream &in, const bool apply_all_indicators_to_manifolds=false)
Definition grid_in.cc:1171
static std::string default_suffix(const Format format)
Definition grid_in.cc:4248
void read_xda(std::istream &in)
Definition grid_in.cc:1388
void read_comsol_mphtxt(std::istream &in)
Definition grid_in.cc:1452
static void skip_comment_lines(std::istream &in, const char comment_start)
Definition grid_in.cc:3967
void attach_triangulation(Triangulation< dim, spacedim > &tria)
Definition grid_in.cc:156
void read_msh(std::istream &in)
Definition grid_in.cc:2008
static Format parse_format(const std::string &format_name)
Definition grid_in.cc:4283
void read_vtu(std::istream &in)
Definition grid_in.cc:590
void read_tecplot(std::istream &in)
Definition grid_in.cc:3393
Format default_format
Definition grid_in.h:951
ExodusIIData read_exodusii(const std::string &filename, const bool apply_all_indicators_to_manifolds=false)
Definition grid_in.cc:3772
void read_dbmesh(std::istream &in)
Definition grid_in.cc:1227
void read_ucd(std::istream &in, const bool apply_all_indicators_to_manifolds=false)
Definition grid_in.cc:917
void read(std::istream &in, Format format=Default)
Definition grid_in.cc:4181
static void debug_output_grid(const std::vector< CellData< dim > > &cells, const std::vector< Point< spacedim > > &vertices, std::ostream &out)
Definition grid_in.cc:3993
static std::string get_format_names()
Definition grid_in.cc:4340
void read_unv(std::istream &in)
Definition grid_in.cc:617
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:2928
Definition point.h:112
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:472
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:160
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:473
Point< 3 > vertices[4]
#define DeclException0(Exception0)
Definition exceptions.h:465
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)
Definition exceptions.h:533
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)
Definition exceptions.h:488
static ::ExceptionBase & ExcUnknownSectionType(int arg1)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition exceptions.h:556
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:510
static ::ExceptionBase & ExcGmshUnsupportedGeometry(int arg1)
static const unsigned int invalid_unsigned_int
Definition types.h:213
std::vector< std::vector< int > > id_to_sideset_ids
Definition grid_in.h:690