]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added GridOut::write_msh()
authorLuca Heltai <luca.heltai@sissa.it>
Thu, 28 Jan 2021 11:32:49 +0000 (12:32 +0100)
committerLuca Heltai <luca.heltai@sissa.it>
Thu, 28 Jan 2021 11:33:12 +0000 (12:33 +0100)
include/deal.II/grid/grid_out.h
source/grid/grid_out.cc
source/grid/grid_out.inst.in
tests/gmsh/gmsh_api_01.cc [new file with mode: 0644]
tests/gmsh/gmsh_api_01.with_gmsh_with_api=on.output [new file with mode: 0644]
tests/gmsh/gmsh_api_02.cc [new file with mode: 0644]
tests/gmsh/gmsh_api_02.with_gmsh_with_api=on.output [new file with mode: 0644]

index 6ae4408c9d923c0a2527294653a6e6d24e125b8b..09c947817820cb6273244033e0ec8dfd552ae661 100644 (file)
@@ -1098,6 +1098,53 @@ public:
   void
   write_msh(const Triangulation<dim, spacedim> &tria, std::ostream &out) const;
 
+#ifdef DEAL_II_GMSH_WITH_API
+  /**
+   * Write the triangulation in any format supported by gmsh API.
+   *
+   * Gmsh API allows writing its output in several formats through their C++
+   * API. This function translates a Triangulation object into a gmsh collection
+   * of entities and calls the gmsh::write() method with the filename passed as
+   * argument. This method generates a different entity for each unique pair of
+   * non default manifold id and boundary id, and writes a gmsh physical group
+   * for each unique combination, allowing you to read back the triangulation
+   * using the GridIn::read_msh() method that takes a string as argument.
+   *
+   * In particular, all cell objects with non default boundary id or non
+   * default manifold id are grouped in a unique physical tag, whose name
+   * contains the boundary and manifold indicators. The names are constructed
+   * using Patterns::Tools::to_value() applied to a `std::map<std::string, int>`
+   * where the keys are either `MaterialID`, `BoundaryID`, or `ManifoldID`,
+   * i.e., a cell with material id 1, and manifold id 3 would be grouped in a
+   * physical tag (whose numbering is not specified), named
+   * `MaterialID:1, ManifoldID:3`.
+   *
+   * For example, calling the method with a hyper ball grid refined once,
+   * would results in the following physical tags defined in the output file:
+   * @code
+   * MeshFormat
+   * 4.1 0 8
+   * $EndMeshFormat
+   * $PhysicalNames
+   * 3
+   * 1 1 "ManifoldID:0"
+   * 1 2 "BoundaryID:-1, ManifoldID:1"
+   * 2 3 "ManifoldID:1"
+   * $EndPhysicalNames
+   * $Entities
+   * ...
+   * @endcode
+   *
+   * The special boundary id `-1` is used to indicate internal boundaries. The
+   * internal boundaries must be specified whenever it is necessary to specify
+   * a non-flat manifold id.
+   */
+  template <int dim, int spacedim>
+  void
+  write_msh(const Triangulation<dim, spacedim> &tria,
+            const std::string &                 filename) const;
+#endif
+
   /**
    * Write the triangulation in the ucd format.
    *
index 3d39cc3952c6fc85b2f7d408d5a83b81de7d7486..74d95b4812c23178fa9dab91ddb88665d7c8e936 100644 (file)
 #include <boost/algorithm/string.hpp>
 #include <boost/archive/binary_oarchive.hpp>
 
+#ifdef DEAL_II_GMSH_WITH_API
+#  include <gmsh.h>
+#endif
+
 #include <algorithm>
 #include <cmath>
 #include <cstring>
@@ -1470,6 +1474,211 @@ GridOut::write_xfig(const Triangulation<2> &tria,
 
 
 
+#ifdef DEAL_II_GMSH_WITH_API
+template <int dim, int spacedim>
+void
+GridOut::write_msh(const Triangulation<dim, spacedim> &tria,
+                   const std::string &                 filename) const
+{
+  // mesh Type renumbering
+  const std::array<int, 8> dealii_to_gmsh_type = {{15, 1, 2, 3, 4, 7, 6, 5}};
+
+  // Vertex renumbering, by dealii type
+  const std::vector<std::vector<unsigned int>> dealii_to_gmsh = {
+    {{{0}},
+     {{0, 1}},
+     {{0, 1, 2}},
+     {{0, 1, 3, 2}},
+     {{0, 1, 2, 3}},
+     {{0, 1, 3, 2, 4}},
+     {{0, 1, 2, 3, 4, 5}},
+     {{0, 1, 3, 2, 4, 5, 7, 6}}}};
+
+  // Extract all vertices (nodes in gmsh terminology), and store their three
+  // dimensional coordinates (regardless of dim).
+  const auto &             vertices = tria.get_vertices();
+  std::vector<double>      coords(3 * vertices.size());
+  std::vector<std::size_t> nodes(vertices.size());
+
+  // Each node has a strictly positive tag. We assign simply its index+1.
+  std::size_t i = 0;
+  for (const auto &p : vertices)
+    {
+      for (unsigned int d = 0; d < spacedim; ++d)
+        coords[i * 3 + d] = p[d];
+      nodes[i] = i + 1;
+      ++i;
+    }
+
+  // Construct one entity tag per boundary and manifold id pair.
+  // We need to be smart here, in order to save some disk space. All cells need
+  // to be written, but only faces and lines that have non default boundary ids
+  // and/or manifold ids. We collect them into pairs, and for each unique pair,
+  // we create a gmsh entity where we store the elements. Pre-count all the
+  // entities, and make sure we know which pair refers to what entity and
+  // vice-versa.
+  using IdPair = std::pair<types::material_id, types::manifold_id>;
+  std::map<IdPair, int> id_pair_to_entity_tag;
+  std::vector<IdPair>   all_pairs;
+  {
+    std::set<IdPair> set_of_pairs;
+    for (const auto &cell : tria.active_cell_iterators())
+      {
+        set_of_pairs.insert({cell->material_id(), cell->manifold_id()});
+        for (const auto &f : cell->face_iterators())
+          if (f->manifold_id() != numbers::flat_manifold_id ||
+              (f->boundary_id() != 0 &&
+               f->boundary_id() != numbers::internal_face_boundary_id))
+            set_of_pairs.insert({f->boundary_id(), f->manifold_id()});
+        if (dim > 2)
+          for (const auto &l : cell->line_indices())
+            {
+              const auto &f = cell->line(l);
+              if (f->manifold_id() != numbers::flat_manifold_id ||
+                  (f->boundary_id() != 0 &&
+                   f->boundary_id() != numbers::internal_face_boundary_id))
+                set_of_pairs.insert({f->boundary_id(), f->manifold_id()});
+            }
+      }
+    all_pairs = {set_of_pairs.begin(), set_of_pairs.end()};
+
+    int entity = 1;
+    for (const auto &p : set_of_pairs)
+      id_pair_to_entity_tag[p] = entity++;
+  }
+
+  const auto n_entity_tags = id_pair_to_entity_tag.size();
+
+  // All elements in the mesh, by entity tag, and by dealii type.
+  std::vector<std::vector<std::vector<std::size_t>>> element_ids(
+    n_entity_tags, std::vector<std::vector<std::size_t>>(8));
+  std::vector<std::vector<std::vector<std::size_t>>> element_nodes(
+    n_entity_tags, std::vector<std::vector<std::size_t>>(8));
+
+  // One elment id counter for all dimensions.
+  std::size_t element_id = 1;
+
+  const auto add_element = [&](const auto &element, const int &entity_tag) {
+    const auto type = element->reference_cell_type();
+
+    Assert(entity_tag > 0, ExcInternalError());
+    // Add all vertex ids. Make sure we renumber to gmsh, and we add 1 to the
+    // global index.
+    for (const auto &v : element->vertex_indices())
+      element_nodes[entity_tag - 1][type].emplace_back(
+        element->vertex_index(dealii_to_gmsh[type][v]) + 1);
+
+    // Save the element id.
+    element_ids[entity_tag - 1][type].emplace_back(element_id);
+    ++element_id;
+  };
+
+  // Will create a separate gmsh entity, only  if it's a cell, or if the
+  // boundary and/or the manifold ids are not the default ones.
+  // In the meanwhile, also store each pair of dimension and entity tag that was
+  // requested.
+  std::set<std::pair<int, int>> dim_entity_tag;
+
+  auto maybe_add_element =
+    [&](const auto &              element,
+        const types::boundary_id &boundary_or_material_id) {
+      const auto struct_dim  = element->structure_dimension;
+      const auto manifold_id = element->manifold_id();
+
+      // Exclude default boundary/manifold id or invalid/flag
+      const bool non_default_boundary_or_material_id =
+        (boundary_or_material_id != 0 &&
+         boundary_or_material_id != numbers::internal_face_boundary_id);
+      const bool non_default_manifold =
+        manifold_id != numbers::flat_manifold_id;
+      if (struct_dim == dim || non_default_boundary_or_material_id ||
+          non_default_manifold)
+        {
+          const auto entity_tag =
+            id_pair_to_entity_tag[{boundary_or_material_id, manifold_id}];
+          add_element(element, entity_tag);
+          dim_entity_tag.insert({struct_dim, entity_tag});
+        }
+    };
+
+  // Loop recursively over all cells, faces, and possibly lines.
+  for (const auto &cell : tria.active_cell_iterators())
+    {
+      maybe_add_element(cell, cell->material_id());
+      for (const auto &face : cell->face_iterators())
+        maybe_add_element(face, face->boundary_id());
+      if (dim > 2)
+        for (const auto &l : cell->line_indices())
+          maybe_add_element(cell->line(l), cell->line(l)->boundary_id());
+    }
+
+  // Now that we collected everything, plug them into gmsh
+  gmsh::initialize();
+  gmsh::option::setNumber("General.Verbosity", 0);
+  gmsh::model::add("Grid generated in deal.II");
+  for (const auto p : dim_entity_tag)
+    {
+      gmsh::model::addDiscreteEntity(p.first, p.second);
+      gmsh::model::mesh::addNodes(p.first, p.second, nodes, coords);
+    }
+
+  for (unsigned int entity_tag = 0; entity_tag < n_entity_tags; ++entity_tag)
+    for (unsigned int t = 1; t < 8; ++t)
+      {
+        const auto all_element_ids   = element_ids[entity_tag][t];
+        const auto all_element_nodes = element_nodes[entity_tag][t];
+        const auto gmsh_t            = dealii_to_gmsh_type[t];
+        if (all_element_ids.size() > 0)
+          gmsh::model::mesh::addElementsByType(entity_tag + 1,
+                                               gmsh_t,
+                                               all_element_ids,
+                                               all_element_nodes);
+      }
+
+  // Make sure nodes belong to the right entities.
+  gmsh::model::mesh::reclassifyNodes();
+  gmsh::model::mesh::removeDuplicateNodes();
+
+  // Now for each individual pair of dim and entry, add a physical group, if
+  // necessary
+  for (const auto &it : dim_entity_tag)
+    {
+      const auto &d           = it.first;
+      const auto &entity_tag  = it.second;
+      const auto &boundary_id = all_pairs[entity_tag - 1].first;
+      const auto &manifold_id = all_pairs[entity_tag - 1].second;
+
+      std::string physical_name;
+      if (d == dim && boundary_id != 0)
+        physical_name += "MaterialID:" + Utilities::int_to_string(
+                                           static_cast<int>(boundary_id));
+      else if (d < dim && boundary_id != 0)
+        physical_name +=
+          "BoundaryID:" +
+          (boundary_id == numbers::internal_face_boundary_id ?
+             "-1" :
+             Utilities::int_to_string(static_cast<int>(boundary_id)));
+
+      std::string sep = physical_name != "" ? ", " : "";
+      if (manifold_id != numbers::flat_manifold_id)
+        physical_name +=
+          sep + "ManifoldID:" +
+          Utilities::int_to_string(static_cast<int>(manifold_id));
+      const auto physical_tag =
+        gmsh::model::addPhysicalGroup(d, {entity_tag}, -1);
+      if (physical_name != "")
+        gmsh::model::setPhysicalName(d, physical_tag, physical_name);
+    }
+
+
+  gmsh::write(filename);
+  gmsh::clear();
+  gmsh::finalize();
+}
+#endif
+
+
+
 namespace
 {
   /**
index dc59139e57e48929a2529b10c82ba5117ec034f3..68baf3875a352dafebaa64704b31e70e8a08d117 100644 (file)
@@ -26,6 +26,10 @@ for (deal_II_dimension : DIMENSIONS)
 
     template void GridOut::write_msh(const Triangulation<deal_II_dimension> &,
                                      std::ostream &) const;
+#ifdef DEAL_II_GMSH_WITH_API
+    template void GridOut::write_msh(const Triangulation<deal_II_dimension> &,
+                                     const std::string &) const;
+#endif
 
 #if deal_II_dimension != 2
     template void GridOut::write_xfig(
@@ -83,6 +87,11 @@ for (deal_II_dimension : DIMENSIONS; deal_II_space_dimension : SPACE_DIMENSIONS)
     template void GridOut::write_msh(
       const Triangulation<deal_II_dimension, deal_II_space_dimension> &,
       std::ostream &) const;
+#  ifdef DEAL_II_GMSH_WITH_API
+    template void GridOut::write_msh(
+      const Triangulation<deal_II_dimension, deal_II_space_dimension> &,
+      const std::string &) const;
+#  endif
     template void GridOut::write_ucd(
       const Triangulation<deal_II_dimension, deal_II_space_dimension> &,
       std::ostream &) const;
diff --git a/tests/gmsh/gmsh_api_01.cc b/tests/gmsh/gmsh_api_01.cc
new file mode 100644 (file)
index 0000000..e092c79
--- /dev/null
@@ -0,0 +1,40 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2021 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+// Create a hyper ball, and use gmsh to output it.
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/tria.h>
+
+#include "../tests.h"
+
+int
+main()
+{
+  initlog();
+
+  const unsigned int dim      = 2;
+  const unsigned int spacedim = 2;
+
+  Triangulation<dim, spacedim> tria;
+  GridGenerator::hyper_ball(tria);
+  tria.refine_global(1);
+
+  GridOut go;
+  go.write_msh(tria, "output.msh");
+
+  cat_file("output.msh");
+}
diff --git a/tests/gmsh/gmsh_api_01.with_gmsh_with_api=on.output b/tests/gmsh/gmsh_api_01.with_gmsh_with_api=on.output
new file mode 100644 (file)
index 0000000..3d0a3e1
--- /dev/null
@@ -0,0 +1,142 @@
+
+$MeshFormat
+4.1 0 8
+$EndMeshFormat
+$PhysicalNames
+3
+1 1 "ManifoldID:0"
+1 2 "BoundaryID:-1, ManifoldID:1"
+2 3 "ManifoldID:1"
+$EndPhysicalNames
+$Entities
+0 2 2 0
+1 -0.9999999999999999 -0.9999999999999999 0 0.9999999999999999 0.9999999999999999 0 1 1 0 
+4 -0.9999999999999999 -0.9999999999999999 0 0.9999999999999999 0.9999999999999999 0 1 2 0 
+2 -0.9999999999999999 -0.9999999999999999 0 0.9999999999999999 0.9999999999999999 0 1 3 0 
+3 -0.2928932188134525 -0.2928932188134525 0 0.2928932188134525 0.2928932188134525 0 1 4 0 
+$EndEntities
+$Nodes
+4 25 1 25
+1 1 0 8
+1
+2
+7
+8
+9
+11
+13
+19
+-0.7071067811865475 -0.7071067811865475 0
+0.7071067811865475 -0.7071067811865475 0
+-0.7071067811865475 0.7071067811865475 0
+0.7071067811865475 0.7071067811865475 0
+-1.110223024625156e-16 -0.9999999999999999 0
+-0.9999999999999999 -1.110223024625156e-16 0
+0.9999999999999999 -1.110223024625156e-16 0
+-1.110223024625156e-16 0.9999999999999999 0
+1 4 0 12
+10
+12
+14
+15
+16
+17
+18
+20
+21
+22
+24
+25
+-0.5 -0.5 0
+0.5 -0.5 0
+0 -0.2928932188134525 0
+-0.2928932188134525 0 0
+0.2928932188134525 0 0
+0 0.2928932188134525 0
+-0.5 0.5 0
+0.5 0.5 0
+-5.551115123125782e-17 -0.6464466094067262 0
+-0.6464466094067262 -5.551115123125783e-17 0
+0.6464466094067262 -5.551115123125782e-17 0
+-5.551115123125783e-17 0.6464466094067262 0
+2 2 0 4
+3
+4
+5
+6
+-0.2928932188134525 -0.2928932188134525 0
+0.2928932188134525 -0.2928932188134525 0
+-0.2928932188134525 0.2928932188134525 0
+0.2928932188134525 0.2928932188134525 0
+2 3 0 1
+23
+0 0 0
+$EndNodes
+$Elements
+4 60 1 60
+1 1 1 8
+3 1 9 
+7 9 2 
+16 1 11 
+23 11 7 
+35 2 13 
+39 13 8 
+48 7 19 
+55 19 8 
+1 4 1 32
+2 9 21 
+4 10 21 
+6 9 21 
+8 21 12 
+10 21 14 
+11 10 21 
+13 21 14 
+14 21 12 
+17 10 22 
+18 11 22 
+20 10 22 
+21 22 15 
+24 22 18 
+25 11 22 
+27 22 18 
+28 22 15 
+34 13 24 
+36 12 24 
+38 13 24 
+40 24 20 
+42 24 16 
+43 12 24 
+45 24 16 
+46 24 20 
+49 18 25 
+50 19 25 
+52 18 25 
+53 25 17 
+56 25 20 
+57 19 25 
+59 25 20 
+60 25 17 
+2 2 3 16
+1 1 9 21 10 
+5 9 2 12 21 
+9 10 21 14 3 
+12 21 12 4 14 
+15 1 10 22 11 
+19 10 3 15 22 
+22 11 22 18 7 
+26 22 15 5 18 
+33 2 13 24 12 
+37 13 8 20 24 
+41 12 24 16 4 
+44 24 20 6 16 
+47 7 18 25 19 
+51 18 5 17 25 
+54 19 25 20 8 
+58 25 17 6 20 
+2 3 3 4
+29 3 14 23 15 
+30 14 4 16 23 
+31 15 23 17 5 
+32 23 16 6 17 
+$EndElements
+
diff --git a/tests/gmsh/gmsh_api_02.cc b/tests/gmsh/gmsh_api_02.cc
new file mode 100644 (file)
index 0000000..7ab7324
--- /dev/null
@@ -0,0 +1,72 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2021 by the deal.II authors
+//
+// This file is part of the deal.II library.
+//
+// The deal.II library is free software; you can use it, redistribute
+// it, and/or modify it under the terms of the GNU Lesser General
+// Public License as published by the Free Software Foundation; either
+// version 2.1 of the License, or (at your option) any later version.
+// The full text of the license can be found in the file LICENSE.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+// Create all reference cells, and output them with gmsh.
+
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/tria.h>
+
+#include "../tests.h"
+
+template <int dim, int spacedim>
+void
+test(const std::uint8_t kind, const std::string out = "")
+{
+  deallog << "Testing kind(" << (int)kind << ") in dimensions "
+          << "<" << dim << "," << spacedim << ">" << std::endl;
+
+  Triangulation<dim, spacedim> tria;
+  ReferenceCell::make_triangulation(
+    ReferenceCell::internal::make_reference_cell_from_int(kind), tria);
+  GridOut go;
+  go.write_msh(tria, "output.msh");
+  cat_file("output.msh");
+  // Uncomment to inspect each file separately
+  // if (out != "")
+  //   go.write_msh(tria, out);
+  (void)out;
+}
+
+int
+main()
+{
+  initlog();
+
+  // Generate and print all reference cells
+  // Lines
+  test<1, 1>(1, "line_11.msh");
+  test<1, 2>(1, "line_12.msh");
+  test<1, 3>(1, "line_13.msh");
+
+  // Triangles
+  test<2, 2>(2, "tri_22.msh");
+  test<2, 3>(2, "tri_23.msh");
+
+  // Quads
+  test<2, 2>(3, "quad_22.msh");
+  test<2, 3>(3, "quad_23.msh");
+
+  // Tet
+  test<3, 3>(4, "tet_3.msh");
+
+  // Pyramid
+  test<3, 3>(5, "pyramid_3.msh");
+
+  // Wedge
+  test<3, 3>(6, "wedge_3.msh");
+
+  // Hex
+  test<3, 3>(7, "hex_3.msh");
+}
diff --git a/tests/gmsh/gmsh_api_02.with_gmsh_with_api=on.output b/tests/gmsh/gmsh_api_02.with_gmsh_with_api=on.output
new file mode 100644 (file)
index 0000000..a2d7795
--- /dev/null
@@ -0,0 +1,303 @@
+
+DEAL::Testing kind(1) in dimensions <1,1>
+$MeshFormat
+4.1 0 8
+$EndMeshFormat
+$PhysicalNames
+1
+0 1 "BoundaryID:1"
+$EndPhysicalNames
+$Entities
+1 1 0 0
+2 0 0 0 1 1 
+1 0 0 0 1 0 0 1 2 0 
+$EndEntities
+$Nodes
+2 2 1 2
+0 2 0 0
+1 1 0 2
+1
+2
+0 0 0
+1 0 0
+$EndNodes
+$Elements
+1 1 1 1
+1 1 1 1
+1 1 2 
+$EndElements
+
+DEAL::Testing kind(1) in dimensions <1,2>
+$MeshFormat
+4.1 0 8
+$EndMeshFormat
+$PhysicalNames
+1
+0 1 "BoundaryID:1"
+$EndPhysicalNames
+$Entities
+1 1 0 0
+2 0 0 0 1 1 
+1 0 0 0 1 0 0 1 2 0 
+$EndEntities
+$Nodes
+2 2 1 2
+0 2 0 0
+1 1 0 2
+1
+2
+0 0 0
+1 0 0
+$EndNodes
+$Elements
+1 1 1 1
+1 1 1 1
+1 1 2 
+$EndElements
+
+DEAL::Testing kind(1) in dimensions <1,3>
+$MeshFormat
+4.1 0 8
+$EndMeshFormat
+$PhysicalNames
+1
+0 1 "BoundaryID:1"
+$EndPhysicalNames
+$Entities
+1 1 0 0
+2 0 0 0 1 1 
+1 0 0 0 1 0 0 1 2 0 
+$EndEntities
+$Nodes
+2 2 1 2
+0 2 0 0
+1 1 0 2
+1
+2
+0 0 0
+1 0 0
+$EndNodes
+$Elements
+1 1 1 1
+1 1 1 1
+1 1 2 
+$EndElements
+
+DEAL::Testing kind(2) in dimensions <2,2>
+$MeshFormat
+4.1 0 8
+$EndMeshFormat
+$Entities
+0 0 1 0
+1 0 0 0 1 1 0 1 1 0 
+$EndEntities
+$Nodes
+1 3 1 3
+2 1 0 3
+1
+2
+3
+0 0 0
+1 0 0
+0 1 0
+$EndNodes
+$Elements
+1 1 1 1
+2 1 2 1
+1 1 2 3 
+$EndElements
+
+DEAL::Testing kind(2) in dimensions <2,3>
+$MeshFormat
+4.1 0 8
+$EndMeshFormat
+$Entities
+0 0 1 0
+1 0 0 0 1 1 0 1 1 0 
+$EndEntities
+$Nodes
+1 3 1 3
+2 1 0 3
+1
+2
+3
+0 0 0
+1 0 0
+0 1 0
+$EndNodes
+$Elements
+1 1 1 1
+2 1 2 1
+1 1 2 3 
+$EndElements
+
+DEAL::Testing kind(3) in dimensions <2,2>
+$MeshFormat
+4.1 0 8
+$EndMeshFormat
+$Entities
+0 0 1 0
+1 0 0 0 1 1 0 1 1 0 
+$EndEntities
+$Nodes
+1 4 1 4
+2 1 0 4
+1
+2
+3
+4
+0 0 0
+1 0 0
+0 1 0
+1 1 0
+$EndNodes
+$Elements
+1 1 1 1
+2 1 3 1
+1 1 2 4 3 
+$EndElements
+
+DEAL::Testing kind(3) in dimensions <2,3>
+$MeshFormat
+4.1 0 8
+$EndMeshFormat
+$Entities
+0 0 1 0
+1 0 0 0 1 1 0 1 1 0 
+$EndEntities
+$Nodes
+1 4 1 4
+2 1 0 4
+1
+2
+3
+4
+0 0 0
+1 0 0
+0 1 0
+1 1 0
+$EndNodes
+$Elements
+1 1 1 1
+2 1 3 1
+1 1 2 4 3 
+$EndElements
+
+DEAL::Testing kind(4) in dimensions <3,3>
+$MeshFormat
+4.1 0 8
+$EndMeshFormat
+$Entities
+0 0 0 1
+1 0 0 0 1 1 1 1 1 0 
+$EndEntities
+$Nodes
+1 4 1 4
+3 1 0 4
+1
+2
+3
+4
+0 0 0
+1 0 0
+0 1 0
+0 0 1
+$EndNodes
+$Elements
+1 1 1 1
+3 1 4 1
+1 1 2 3 4 
+$EndElements
+
+DEAL::Testing kind(5) in dimensions <3,3>
+$MeshFormat
+4.1 0 8
+$EndMeshFormat
+$Entities
+0 0 0 1
+1 -1 -1 0 1 1 1 1 1 0 
+$EndEntities
+$Nodes
+1 5 1 5
+3 1 0 5
+1
+2
+3
+4
+5
+-1 -1 0
+1 -1 0
+-1 1 0
+1 1 0
+0 0 1
+$EndNodes
+$Elements
+1 1 1 1
+3 1 7 1
+1 1 2 4 3 5 
+$EndElements
+
+DEAL::Testing kind(6) in dimensions <3,3>
+$MeshFormat
+4.1 0 8
+$EndMeshFormat
+$Entities
+0 0 0 1
+1 0 0 0 1 1 1 1 1 0 
+$EndEntities
+$Nodes
+1 6 1 6
+3 1 0 6
+1
+2
+3
+4
+5
+6
+1 0 0
+0 1 0
+0 0 0
+1 0 1
+0 1 1
+0 0 1
+$EndNodes
+$Elements
+1 1 1 1
+3 1 6 1
+1 1 2 3 4 5 6 
+$EndElements
+
+DEAL::Testing kind(7) in dimensions <3,3>
+$MeshFormat
+4.1 0 8
+$EndMeshFormat
+$Entities
+0 0 0 1
+1 0 0 0 1 1 1 1 1 0 
+$EndEntities
+$Nodes
+1 8 1 8
+3 1 0 8
+1
+2
+3
+4
+5
+6
+7
+8
+0 0 0
+1 0 0
+0 1 0
+1 1 0
+0 0 1
+1 0 1
+0 1 1
+1 1 1
+$EndNodes
+$Elements
+1 1 1 1
+3 1 5 1
+1 1 2 4 3 5 6 8 7 
+$EndElements
+

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.