Reference documentation for deal.II version GIT relicensing-136-gb80d0be4af 2024-03-18 08:20:02+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\}}\)
Loading...
Searching...
No Matches
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 - 2023 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
22#include <deal.II/base/point.h>
24
25#include <iostream>
26#include <string>
27#include <vector>
28
30
31// Forward declarations
32#ifndef DOXYGEN
33template <int dim, int spacedim>
35class Triangulation;
36template <int dim>
37struct CellData;
38#endif
39
310template <int dim, int spacedim = dim>
312{
313public:
345
349 GridIn();
350
355
359 void
361
366 void
367 read(std::istream &in, Format format = Default);
368
374 void
375 read(const std::string &in, Format format = Default);
376
417 void
418 read_vtk(std::istream &in);
419
434 void
435 read_vtu(std::istream &in);
436
437
445 void
446 read_unv(std::istream &in);
447
457 void
458 read_ucd(std::istream &in,
459 const bool apply_all_indicators_to_manifolds = false);
460
500 void
501 read_abaqus(std::istream &in,
502 const bool apply_all_indicators_to_manifolds = false);
503
507 void
508 read_dbmesh(std::istream &in);
509
513 void
514 read_xda(std::istream &in);
515
523 void
524 read_msh(std::istream &in);
525
526#ifdef DEAL_II_GMSH_WITH_API
584 void
585 read_msh(const std::string &filename);
586#endif
587
624 void
625 read_comsol_mphtxt(std::istream &in);
626
631 void
632 read_tecplot(std::istream &in);
633
668 void
669 read_assimp(const std::string &filename,
670 const unsigned int mesh_index = numbers::invalid_unsigned_int,
671 const bool remove_duplicates = true,
672 const double tol = 1e-12,
673 const bool ignore_unsupported_element_types = true);
674
684 {
689 std::vector<std::vector<int>> id_to_sideset_ids;
690 };
691
749 read_exodusii(const std::string &filename,
750 const bool apply_all_indicators_to_manifolds = false);
751
755 static std::string
756 default_suffix(const Format format);
757
761 static Format
762 parse_format(const std::string &format_name);
763
769 static std::string
771
776 int,
777 << "The section type <" << arg1 << "> in an UNV "
778 << "input file is not implemented.");
779
784 int,
785 << "The element type <" << arg1 << "> in an UNV "
786 << "input file is not implemented.");
787
792 std::string,
793 << "The identifier <" << arg1 << "> as name of a "
794 << "part in an UCD input file is unknown or the "
795 << "respective input routine is not implemented."
796 << "(Maybe the space dimension of triangulation and "
797 << "input file do not match?");
802 "No Triangulation has been attached to this GridIn object "
803 "so that nothing can be filled during any read function "
804 "calls. Please pass a reference to the Triangulation tria "
805 "to be filled in the constructor GridIn(tria) or attach "
806 "it with the function call GridIn::attach_triangulation().");
812 int,
813 int,
814 << "While creating cell " << arg1
815 << ", you are referencing a vertex with index " << arg2
816 << " but no vertex with this index has been described in the input file.");
822 int,
823 int,
824 int,
825 << "While creating cell " << arg1 << " (which is numbered as " << arg2
826 << " in the input file), you are referencing a vertex with index " << arg3
827 << " but no vertex with this index has been described in the input file.");
836 std::string,
837 << "The string <" << arg1
838 << "> is not recognized at the present"
839 << " position of a DB Mesh file.");
840
846 int,
847 << "The specified dimension " << arg1
848 << " is not the same as that of the triangulation to be created.");
849
851 std::string,
852 << "The string <" << arg1
853 << "> is not recognized at the present"
854 << " position of a Gmsh Mesh file.");
855
857 int,
858 << "The Element Identifier <" << arg1 << "> is not "
859 << "supported in the deal.II library when "
860 << "reading meshes in " << dim << " dimensions.\n"
861 << "Supported elements are: \n"
862 << "ELM-TYPE\n"
863 << "1 Line (2 nodes, 1 edge).\n"
864 << "3 Quadrilateral (4 nodes, 4 edges).\n"
865 << "5 Hexahedron (8 nodes, 12 edges, 6 faces) when in 3d.\n"
866 << "15 Point (1 node, ignored when read)");
867
868
871 unsigned int,
872 unsigned int,
873 "While reading a gmsh file, the reader function did not find "
874 "any cells. This sometimes happens if the file only contains a "
875 "surface mesh, but not a volume mesh."
876 "\n\n"
877 "The reader function did find " +
878 std::to_string(arg1) + " lines and " + std::to_string(arg2) +
879 " facets (surface triangles or quadrilaterals).");
880
881protected:
886
907 static void
908 debug_output_grid(const std::vector<CellData<dim>> &cells,
909 const std::vector<Point<spacedim>> &vertices,
910 std::ostream &out);
911
912private:
917 static void
918 skip_empty_lines(std::istream &in);
919
927 static void
928 skip_comment_lines(std::istream &in, const char comment_start);
929
937 static void
938 parse_tecplot_header(std::string &header,
939 std::vector<unsigned int> &tecplot2deal,
940 unsigned int &n_vars,
941 unsigned int &n_vertices,
942 unsigned int &n_cells,
943 std::vector<unsigned int> &IJK,
944 bool &structured,
945 bool &blocked);
946
951};
952
953/* -------------- declaration of explicit specializations ------------- */
954
955#ifndef DOXYGEN
956
957template <>
958void
959GridIn<2>::debug_output_grid(const std::vector<CellData<2>> &cells,
960 const std::vector<Point<2>> &vertices,
961 std::ostream &out);
962
963
964template <>
965void
966GridIn<2, 3>::debug_output_grid(const std::vector<CellData<2>> &cells,
967 const std::vector<Point<3>> &vertices,
968 std::ostream &out);
969template <>
970void
971GridIn<3>::debug_output_grid(const std::vector<CellData<3>> &cells,
972 const std::vector<Point<3>> &vertices,
973 std::ostream &out);
974#endif // DOXYGEN
975
977
978#endif
@ Default
Use GridIn::default_format stored in this object.
Definition grid_in.h:321
@ vtu
Use read_vtu()
Definition grid_in.h:339
@ xda
Use read_xda()
Definition grid_in.h:331
@ tecplot
Use read_tecplot()
Definition grid_in.h:335
@ msh
Use read_msh()
Definition grid_in.h:333
@ ucd
Use read_ucd()
Definition grid_in.h:325
@ dbmesh
Use read_dbmesh()
Definition grid_in.h:329
@ abaqus
Use read_abaqus()
Definition grid_in.h:327
@ assimp
Use read_assimp()
Definition grid_in.h:341
@ unv
Use read_unv()
Definition grid_in.h:323
@ exodusii
Use read_exodusii()
Definition grid_in.h:343
@ vtk
Use read_vtk()
Definition grid_in.h:337
SmartPointer< Triangulation< dim, spacedim >, GridIn< dim, spacedim > > tria
Definition grid_in.h:885
void read_vtk(std::istream &in)
Definition grid_in.cc:163
GridIn()
Definition grid_in.cc:137
static void skip_empty_lines(std::istream &in)
Definition grid_in.cc:3938
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:3403
void read_abaqus(std::istream &in, const bool apply_all_indicators_to_manifolds=false)
Definition grid_in.cc:1169
static std::string default_suffix(const Format format)
Definition grid_in.cc:4249
void read_xda(std::istream &in)
Definition grid_in.cc:1386
void read_comsol_mphtxt(std::istream &in)
Definition grid_in.cc:1450
static void skip_comment_lines(std::istream &in, const char comment_start)
Definition grid_in.cc:3968
void attach_triangulation(Triangulation< dim, spacedim > &tria)
Definition grid_in.cc:154
void read_msh(std::istream &in)
Definition grid_in.cc:2006
static Format parse_format(const std::string &format_name)
Definition grid_in.cc:4284
void read_vtu(std::istream &in)
Definition grid_in.cc:588
void read_tecplot(std::istream &in)
Definition grid_in.cc:3394
Format default_format
Definition grid_in.h:950
ExodusIIData read_exodusii(const std::string &filename, const bool apply_all_indicators_to_manifolds=false)
Definition grid_in.cc:3773
void read_dbmesh(std::istream &in)
Definition grid_in.cc:1225
void read_ucd(std::istream &in, const bool apply_all_indicators_to_manifolds=false)
Definition grid_in.cc:915
void read(std::istream &in, Format format=Default)
Definition grid_in.cc:4182
static void debug_output_grid(const std::vector< CellData< dim > > &cells, const std::vector< Point< spacedim > > &vertices, std::ostream &out)
Definition grid_in.cc:3994
static std::string get_format_names()
Definition grid_in.cc:4341
void read_unv(std::istream &in)
Definition grid_in.cc:615
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:2929
Definition point.h:111
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:502
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:177
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:503
Point< 3 > vertices[4]
#define DeclException0(Exception0)
Definition exceptions.h:471
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:539
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:494
static ::ExceptionBase & ExcUnknownSectionType(int arg1)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition exceptions.h:562
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:516
static ::ExceptionBase & ExcGmshUnsupportedGeometry(int arg1)
static const unsigned int invalid_unsigned_int
Definition types.h:220
std::vector< std::vector< int > > id_to_sideset_ids
Definition grid_in.h:689