"find a valid muparser library on your system and also did "
"not choose the one that comes bundled with deal.II.");
+ /**
+ * This function requires support for the Assimp library.
+ */
+ DeclExceptionMsg (ExcNeedsAssimp,
+ "You are attempting to use functionality that is only available "
+ "if deal.II was configured to use Assimp, but cmake did not "
+ "find a valid Assimp library.");
+
#ifdef DEAL_II_WITH_CUDA
/**
* This exception is raised if an error happened in a CUDA kernel.
/**
* This class implements an input mechanism for grid data. It allows to read a
* grid structure into a triangulation object. At present, UCD (unstructured
- * cell data), DB Mesh, XDA, Gmsh, Tecplot, NetCDF, UNV, VTK, and Cubit are
- * supported as input format for grid data. Any numerical data other than
+ * cell data), DB Mesh, XDA, Gmsh, Tecplot, NetCDF, UNV, VTK, ASSIMP, and Cubit
+ * are supported as input format for grid data. Any numerical data other than
* geometric (vertex locations) and topological (how vertices form cells,
* faces, and edges) information is ignored, but the readers for the various
* formats generally do read information that associates material ids or
- * boundary ids to cells or faces (see
- * @ref GlossMaterialId "this"
- * and
- * @ref GlossBoundaryIndicator "this"
- * glossary entry for more information).
+ * boundary ids to cells or faces (see @ref GlossMaterialId "this" and @ref
+ * GlossBoundaryIndicator "this" glossary entry for more information).
*
* @note Since deal.II only supports line, quadrilateral and hexahedral
* meshes, the functions in this class can only read meshes that consist
/// Use read_tecplot()
tecplot,
/// Use read_vtk()
- vtk
+ vtk,
+ /// Use read_assimp()
+ assimp,
};
/**
*/
void read_tecplot (std::istream &in);
+ /**
+ * Read in a file supported by Assimp, and generate a Triangulation
+ * out of it. If you specify a @p mesh_index, only the mesh with
+ * the given index will be extracted, otherwise all meshes which are
+ * present in the file will be used to generate the Triangulation.
+ *
+ * This function can only be used to read two-dimensional meshes (possibly
+ * embedded in three dimensions), as this is the only format supported
+ * by the supported graphical formats.
+ *
+ * If @p remove_duplicates is set to true (the default), then
+ * duplicated vertices will be removed if their distance is lower
+ * than @p tol.
+ *
+ * Only the elements compatible with the given dimension and spacedimension
+ * will be extracted from the mesh, and only those elements that are
+ * compatible with deal.II are supported. If you set
+ * `ignore_unsupported_element_types`, all the other element types are simply
+ * ignored by this algorithm. If your mesh contains a mixture of triangles
+ * and quadrilaterals, for example, only the quadrilaterals will be
+ * extracted. Careful that your mesh may not make any sense if you are mixing
+ * compatible and incompatible element types. If
+ * `ignore_unsupported_element_types` is set to `false`, then an exception is
+ * thrown when an unsupporte type is encountered.
+ *
+ * @param filename The file to read from
+ * @param mesh_index Index of the mesh within the file
+ * @param remove_duplicates Remove duplicated vertices
+ * @param tol Tolerance to use when removing vertices
+ * @param ignore_unsupported_element_types Don't throw exceptions if we
+ * encounter unsupported types during parsing
+ */
+ 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);
+
/**
* Return the standard suffix for a file in this format.
*/
GridIn<3>::debug_output_grid (const std::vector<CellData<3> > &cells,
const std::vector<Point<3> > &vertices,
std::ostream &out);
-
#endif // DOXYGEN
DEAL_II_NAMESPACE_CLOSE
#include <netcdfcpp.h>
#endif
+#ifdef DEAL_II_WITH_ASSIMP
+DEAL_II_DISABLE_EXTRA_DIAGNOSTICS
+#pragma GCC diagnostic push
+#pragma GCC diagnostic ignored "-Waddress-of-packed-member"
+#include <assimp/Importer.hpp> // C++ importer interface
+#include <assimp/scene.h> // Output data structure
+#include <assimp/postprocess.h> // Post processing flags
+#pragma GCC diagnostic pop
+DEAL_II_ENABLE_EXTRA_DIAGNOSTICS
+
+#undef AI_CONFIG_PP_RVC_FLAGS
+#define AI_CONFIG_PP_RVC_FLAGS \
+ aiComponent_NORMALS | \
+ aiComponent_TANGENTS_AND_BITANGENTS | \
+ aiComponent_COLORS | \
+ aiComponent_TEXCOORDS | \
+ aiComponent_BONEWEIGHTS | \
+ aiComponent_ANIMATIONS | \
+ aiComponent_TEXTURES | \
+ aiComponent_LIGHTS | \
+ aiComponent_CAMERAS | \
+ aiComponent_MATERIALS
+#endif
+
DEAL_II_NAMESPACE_OPEN
}
+
+template <int dim, int spacedim>
+void GridIn<dim, spacedim>::read_assimp(const std::string &filename,
+ const unsigned int mesh_index,
+ const bool remove_duplicates,
+ const double tol,
+ const bool ignore_unsupported_types)
+{
+ // Only good for surface grids.
+ AssertThrow(dim<3, ExcImpossibleInDim(dim));
+
+ // Create an istance of the Importer class
+ Assimp::Importer importer;
+
+ // And have it read the given file with some postprocessing
+ const aiScene *scene = importer.ReadFile( filename.c_str(),
+ aiProcess_RemoveComponent |
+ aiProcess_JoinIdenticalVertices |
+ aiProcess_ImproveCacheLocality |
+ aiProcess_SortByPType |
+ aiProcess_OptimizeGraph |
+ aiProcess_OptimizeMeshes);
+
+ // If the import failed, report it
+ AssertThrow(scene, ExcMessage(importer.GetErrorString()))
+
+ AssertThrow(scene->mNumMeshes, ExcMessage("Input file contains no meshes."));
+
+ AssertThrow((mesh_index == numbers::invalid_unsigned_int) ||
+ (mesh_index < scene->mNumMeshes),
+ ExcMessage("Too few meshes in the file."));
+
+ int start_mesh = (mesh_index == numbers::invalid_unsigned_int ?
+ 0 : mesh_index);
+ int end_mesh = (mesh_index == numbers::invalid_unsigned_int ?
+ scene->mNumMeshes : mesh_index+1);
+
+ // Deal.II objects are created empty, and then filled with imported file.
+ std::vector<Point<spacedim> > vertices;
+ std::vector<CellData<dim> > cells;
+ SubCellData subcelldata;
+
+ // A series of counters to merge cells.
+ unsigned int v_offset=0;
+ unsigned int c_offset=0;
+
+ // The index of the mesh will be used as a material index.
+ for (int m=start_mesh; m<end_mesh; ++m)
+ {
+ const aiMesh *mesh = scene->mMeshes[m];
+
+ // Check that we know what to do with this mesh, otherwise just
+ // ignore it
+ if ( (dim == 2) && mesh->mPrimitiveTypes != aiPrimitiveType_POLYGON)
+ {
+ std::stringstream s;
+ s << "Incompatible mesh " << m
+ << "/" << scene->mNumMeshes;
+ AssertThrow(ignore_unsupported_types, ExcMessage(s.str()));
+ continue;
+ }
+ else if ( (dim == 1) && mesh->mPrimitiveTypes != aiPrimitiveType_LINE)
+ {
+ std::stringstream s;
+ s << "Incompatible mesh " << m
+ << "/" << scene->mNumMeshes;
+ AssertThrow(ignore_unsupported_types, ExcMessage(s.str()));
+ continue;
+ }
+ // Vertices
+ const unsigned int n_vertices = mesh->mNumVertices;
+ const aiVector3D *mVertices = mesh->mVertices;
+
+ // Faces
+ const unsigned int n_faces = mesh->mNumFaces;
+ const aiFace *mFaces = mesh->mFaces;
+
+ vertices.resize(v_offset+n_vertices);
+ cells.resize(c_offset+n_faces);
+
+ for (unsigned int i=0; i<n_vertices; ++i)
+ for (unsigned int d=0; d<spacedim; ++d)
+ vertices[i+v_offset][d] = mVertices[i][d];
+
+ unsigned int valid_cell = c_offset;
+ for (unsigned int i=0; i<n_faces; ++i)
+ {
+ if (mFaces[i].mNumIndices == GeometryInfo<dim>::vertices_per_cell)
+ {
+ for (unsigned int f=0; f<GeometryInfo<dim>::vertices_per_cell; ++f)
+ {
+ cells[valid_cell].vertices[f] = mFaces[i].mIndices[f]+v_offset;
+ }
+ cells[valid_cell].material_id = (types::material_id) m;
+ ++valid_cell;
+ }
+ else
+ {
+ std::stringstream s;
+ s << "Face " << i << " of mesh "
+ << m << " has " << mFaces[i].mNumIndices
+ << " vertices. We expected only "
+ << GeometryInfo<dim>::vertices_per_cell;
+ AssertThrow(ignore_unsupported_types, ExcMessage(s.str()));
+ }
+ }
+ cells.resize(valid_cell);
+
+ // The vertices are added all at once. Cells are check for
+ // validity, so only valid_cells are now present in the deal.II
+ // list of cells.
+ v_offset += n_vertices;
+ c_offset = valid_cell;
+ }
+
+ // No cells were read
+ if (cells.size() == 0)
+ return;
+
+ if (remove_duplicates)
+ {
+ // The function delete_duplicated_vertices() needs to be called more
+ // than once if a vertex is duplicated more than once. So we keep
+ // calling it until the number of vertices does not change any more.
+ unsigned int n_verts = 0;
+ while (n_verts != vertices.size())
+ {
+ n_verts = vertices.size();
+ std::vector<unsigned int> considered_vertices;
+ GridTools::delete_duplicated_vertices(vertices, cells, subcelldata,
+ considered_vertices, tol);
+ }
+ }
+
+ GridTools::delete_unused_vertices(vertices, cells, subcelldata);
+ if (dim == spacedim)
+ GridReordering<dim, spacedim>::invert_all_cells_of_negative_grid
+ (vertices,cells);
+
+ GridReordering<dim, spacedim>::reorder_cells(cells);
+ if (dim == 2)
+ tria->create_triangulation_compatibility(vertices, cells, subcelldata);
+ else
+ tria->create_triangulation(vertices, cells, subcelldata);
+}
+
+
template <int dim, int spacedim>
void GridIn<dim, spacedim>::skip_empty_lines (std::istream &in)
{
case netcdf:
Assert(false, ExcMessage("There is no read_netcdf(istream &) function. "
- "Use the read(_netcdf)(string &filename) "
+ "Use the read_netcdf)(string &filename) "
"functions, instead."));
return;
read_tecplot (in);
return;
+ case assimp:
+ Assert(false, ExcMessage("There is no read_assimp(istream &) function. "
+ "Use the read_netcdf)(string &filename, ...) "
+ "functions, instead."));
+ return;
+
case Default:
break;
}
template <int dim, int spacedim>
std::string GridIn<dim, spacedim>::get_format_names ()
{
- return "dbmesh|msh|unv|vtk|ucd|abaqus|xda|netcdf|tecplot";
+ return "dbmesh|msh|unv|vtk|ucd|abaqus|xda|netcdf|tecplot|assimp";
}
namespace
--- /dev/null
+//-----------------------------------------------------------
+//
+// Copyright (C) 2017 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 at
+// the top level of the deal.II distribution.
+//
+//-----------------------------------------------------------
+
+#include "../tests.h"
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_in.h>
+#include <deal.II/grid/grid_out.h>
+
+
+int main ()
+{
+ initlog();
+ Triangulation<2,3> tria;
+ GridIn<2,3> gridin;
+ gridin.attach_triangulation(tria);
+ gridin.read_assimp(SOURCE_DIR "/grids/torus.obj", -1, true, 1e-3);
+
+ GridOut go;
+ go.write_msh(tria, deallog.get_file_stream());
+ return 0;
+}
--- /dev/null
+
+$NOD
+16
+1 -0.364915 0.735417 -1.86297
+2 0.935085 0.735417 -0.562966
+3 0.635085 0.435417 -0.562966
+4 -0.364915 0.435417 -1.56297
+5 -0.364915 1.03542 -1.56297
+6 0.635085 1.03542 -0.562966
+7 -0.364915 0.735417 -1.26297
+8 0.335085 0.735417 -0.562966
+9 -1.66492 0.735417 -0.562966
+10 -1.36492 0.435417 -0.562966
+11 -1.36492 1.03542 -0.562966
+12 -1.06491 0.735417 -0.562966
+13 -0.364915 0.735417 0.737034
+14 -0.364915 0.435417 0.437034
+15 -0.364915 1.03542 0.437034
+16 -0.364915 0.735417 0.137034
+$ENDNOD
+$ELM
+16
+1 3 0 0 4 3 4 1 2
+2 3 0 0 4 6 2 1 5
+3 3 0 0 4 8 6 5 7
+4 3 0 0 4 3 8 7 4
+5 3 0 0 4 4 10 9 1
+6 3 0 0 4 5 1 9 11
+7 3 0 0 4 7 5 11 12
+8 3 0 0 4 4 7 12 10
+9 3 0 0 4 10 14 13 9
+10 3 0 0 4 11 9 13 15
+11 3 0 0 4 12 11 15 16
+12 3 0 0 4 10 12 16 14
+13 3 0 0 4 14 3 2 13
+14 3 0 0 4 15 13 2 6
+15 3 0 0 4 16 15 6 8
+16 3 0 0 4 14 16 8 3
+$ENDELM
--- /dev/null
+# Blender v2.74 (sub 0) OBJ File: ''
+# www.blender.org
+mtllib torus.mtl
+o SuperToroid
+v 0.935085 0.735417 -0.562966
+v 0.635085 1.035417 -0.562966
+v 0.335085 0.735417 -0.562966
+v 0.635085 0.435417 -0.562966
+v -0.364915 0.735417 -1.862966
+v -0.364915 1.035417 -1.562966
+v -0.364915 0.735417 -1.262966
+v -0.364915 0.435417 -1.562966
+v -1.664915 0.735417 -0.562966
+v -1.364915 1.035417 -0.562966
+v -1.064915 0.735417 -0.562966
+v -1.364915 0.435417 -0.562966
+v -0.364915 0.735417 0.737034
+v -0.364915 1.035417 0.437034
+v -0.364915 0.735417 0.137034
+v -0.364915 0.435417 0.437034
+vn 0.577400 -0.577400 -0.577400
+vn 0.577400 0.577400 -0.577400
+vn -0.577400 0.577400 0.577400
+vn -0.577400 -0.577400 0.577400
+vn -0.577400 -0.577400 -0.577400
+vn -0.577400 0.577400 -0.577400
+vn 0.577400 0.577400 0.577400
+vn 0.577400 -0.577400 0.577400
+usemtl None
+s off
+f 5//1 1//1 4//1 8//1
+f 1//2 5//2 6//2 2//2
+f 2//3 6//3 7//3 3//3
+f 3//4 7//4 8//4 4//4
+f 9//5 5//5 8//5 12//5
+f 5//6 9//6 10//6 6//6
+f 6//7 10//7 11//7 7//7
+f 7//8 11//8 12//8 8//8
+f 13//4 9//4 12//4 16//4
+f 9//3 13//3 14//3 10//3
+f 10//2 14//2 15//2 11//2
+f 11//1 15//1 16//1 12//1
+f 1//8 13//8 16//8 4//8
+f 13//7 1//7 2//7 14//7
+f 14//6 2//6 3//6 15//6
+f 15//5 3//5 4//5 16//5