deal.II version GIT relicensing-2330-gf6dfc6c370 2025-01-06 13:10: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\}}\)
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 - 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
312template <int dim, int spacedim = dim>
314{
315public:
347
351 GridIn();
352
357
361 void
363
368 void
369 read(std::istream &in, Format format = Default);
370
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
789 const std::map<std::string, Vector<double>> &
790 get_cell_data() const;
791
796 int,
797 << "The section type <" << arg1 << "> in an UNV "
798 << "input file is not implemented.");
799
804 int,
805 << "The element type <" << arg1 << "> in an UNV "
806 << "input file is not implemented.");
807
812 std::string,
813 << "The identifier <" << arg1 << "> as name of a "
814 << "part in an UCD input file is unknown or the "
815 << "respective input routine is not implemented."
816 << "(Maybe the space dimension of triangulation and "
817 << "input file do not match?");
822 "No Triangulation has been attached to this GridIn object "
823 "so that nothing can be filled during any read function "
824 "calls. Please pass a reference to the Triangulation tria "
825 "to be filled in the constructor GridIn(tria) or attach "
826 "it with the function call GridIn::attach_triangulation().");
832 int,
833 int,
834 << "While creating cell " << arg1
835 << ", you are referencing a vertex with index " << arg2
836 << " but no vertex with this index has been described in the input file.");
842 int,
843 int,
844 int,
845 << "While creating cell " << arg1 << " (which is numbered as " << arg2
846 << " in the input file), you are referencing a vertex with index " << arg3
847 << " but no vertex with this index has been described in the input file.");
856 std::string,
857 << "The string <" << arg1
858 << "> is not recognized at the present"
859 << " position of a DB Mesh file.");
860
866 int,
867 << "The specified dimension " << arg1
868 << " is not the same as that of the triangulation to be created.");
869
871 std::string,
872 << "The string <" << arg1
873 << "> is not recognized at the present"
874 << " position of a Gmsh Mesh file.");
875
877 int,
878 << "The Element Identifier <" << arg1 << "> is not "
879 << "supported in the deal.II library when "
880 << "reading meshes in " << dim << " dimensions.\n"
881 << "Supported elements are: \n"
882 << "ELM-TYPE\n"
883 << "1 Line (2 nodes, 1 edge).\n"
884 << "2 Triangle (3 nodes, 3 edges).\n"
885 << "3 Quadrilateral (4 nodes, 4 edges).\n"
886 << "4 Tetrahedron (4 nodes, 6 edges, 4 faces) when in 3d.\n"
887 << "5 Hexahedron (8 nodes, 12 edges, 6 faces) when in 3d.\n"
888 << "15 Point (1 node, ignored when read).");
889
890
893 unsigned int,
894 unsigned int,
895 "While reading a gmsh file, the reader function did not find "
896 "any cells. This sometimes happens if the file only contains a "
897 "surface mesh, but not a volume mesh."
898 "\n\n"
899 "The reader function did find " +
900 std::to_string(arg1) + " lines and " + std::to_string(arg2) +
901 " facets (surface triangles or quadrilaterals).");
902
903protected:
908
929 static void
930 debug_output_grid(const std::vector<CellData<dim>> &cells,
931 const std::vector<Point<spacedim>> &vertices,
932 std::ostream &out);
933
934private:
939 static void
940 skip_empty_lines(std::istream &in);
941
949 static void
950 skip_comment_lines(std::istream &in, const char comment_start);
951
959 static void
960 parse_tecplot_header(std::string &header,
961 std::vector<unsigned int> &tecplot2deal,
962 unsigned int &n_vars,
963 unsigned int &n_vertices,
964 unsigned int &n_cells,
965 std::vector<unsigned int> &IJK,
966 bool &structured,
967 bool &blocked);
968
973
983 std::map<std::string, Vector<double>> cell_data;
984};
985
986/* -------------- declaration of explicit specializations ------------- */
987
988#ifndef DOXYGEN
989
990template <>
991void
992GridIn<2>::debug_output_grid(const std::vector<CellData<2>> &cells,
993 const std::vector<Point<2>> &vertices,
994 std::ostream &out);
995
996
997template <>
998void
999GridIn<2, 3>::debug_output_grid(const std::vector<CellData<2>> &cells,
1000 const std::vector<Point<3>> &vertices,
1001 std::ostream &out);
1002template <>
1003void
1004GridIn<3>::debug_output_grid(const std::vector<CellData<3>> &cells,
1005 const std::vector<Point<3>> &vertices,
1006 std::ostream &out);
1007#endif // DOXYGEN
1008
1010
1011#endif
@ Default
Use GridIn::default_format stored in this object.
Definition grid_in.h:323
@ vtu
Use read_vtu()
Definition grid_in.h:341
@ xda
Use read_xda()
Definition grid_in.h:333
@ tecplot
Use read_tecplot()
Definition grid_in.h:337
@ msh
Use read_msh()
Definition grid_in.h:335
@ ucd
Use read_ucd()
Definition grid_in.h:327
@ dbmesh
Use read_dbmesh()
Definition grid_in.h:331
@ abaqus
Use read_abaqus()
Definition grid_in.h:329
@ assimp
Use read_assimp()
Definition grid_in.h:343
@ unv
Use read_unv()
Definition grid_in.h:325
@ exodusii
Use read_exodusii()
Definition grid_in.h:345
@ vtk
Use read_vtk()
Definition grid_in.h:339
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:4040
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:3504
std::map< std::string, Vector< double > > cell_data
Definition grid_in.h:983
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:4344
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:4070
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:4379
void read_vtu(std::istream &in)
Definition grid_in.cc:665
void read_tecplot(std::istream &in)
Definition grid_in.cc:3495
Format default_format
Definition grid_in.h:972
ExodusIIData read_exodusii(const std::string &filename, const bool apply_all_indicators_to_manifolds=false)
Definition grid_in.cc:3875
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:907
void read(std::istream &in, Format format=Default)
Definition grid_in.cc:4277
static void debug_output_grid(const std::vector< CellData< dim > > &cells, const std::vector< Point< spacedim > > &vertices, std::ostream &out)
Definition grid_in.cc:4096
static std::string get_format_names()
Definition grid_in.cc:4436
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:3030
Definition point.h:111
#define DEAL_II_NAMESPACE_OPEN
Definition config.h:498
#define DEAL_II_CXX20_REQUIRES(condition)
Definition config.h:175
#define DEAL_II_NAMESPACE_CLOSE
Definition config.h:499
#define DeclException0(Exception0)
Definition exceptions.h:466
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:534
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:489
static ::ExceptionBase & ExcUnknownSectionType(int arg1)
#define DeclException3(Exception3, type1, type2, type3, outsequence)
Definition exceptions.h:557
#define DeclException1(Exception1, type1, outsequence)
Definition exceptions.h:511
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:690