This group describes the simplex support in deal.II.
Simplex and mixed meshes in deal.II are work in progress. Large parts of the library have been ported to be able to operate on such kind of meshes. However, there are still many functions that need to be generalized.
Important Simplex Functionality
Here is an incomplete list of functionality related to simplex computations:
- Mesh generation: GridGenerator::implicit_function(), GridGenerator::convert_hypercube_to_simplex_mesh(), GridGenerator::subdivided_hyper_rectangle_with_simplices(), GridGenerator::subdivided_hyper_cube_with_simplices()
- Quadratures: QGaussWedge, QGaussSimplex, QWitherdenVincentSimplex
- FiniteElements: FE_SimplexP, FE_SimplexDGP, FE_SimplexP_Bubbles FE_PyramidP, FE_PyramidDGP, FE_WedgeP, FE_WedgeDGP
- Mapping: MappingFE
- Other: GridIn::read_vtk(), GridIn::read_msh(), GridIn::read_comsol_mphtxt()
Examples
You can get a good overview of the ported functionalities by taking a look at the tests in the folder "tests/simplex". In the following, we provide two very basic examples to get you started and to provide some implementation details.
Example: simplex mesh
The following code shows how to work with simplex meshes:
#include <fstream>
#include <iostream>
{
public:
void run();
private:
void solve();
};
, fe(2)
{}
{
std::cout <<
"Number of active cells: " <<
triangulation.n_active_cells()
<< std::endl;
}
{
dof_handler.distribute_dofs(fe);
std::cout << "Number of degrees of freedom: " << dof_handler.n_dofs()
<< std::endl;
sparsity_pattern.copy_from(
dsp);
system_matrix.reinit(sparsity_pattern);
solution.reinit(dof_handler.n_dofs());
}
{
fe,
const unsigned int dofs_per_cell = fe.n_dofs_per_cell();
std::vector<types::global_dof_index> local_dof_indices(dofs_per_cell);
for (const auto &cell : dof_handler.active_cell_iterators())
{
fe_values.reinit(cell);
for (const unsigned int q_index : fe_values.quadrature_point_indices())
{
for (const unsigned int i : fe_values.dof_indices())
(fe_values.shape_grad(i, q_index) *
fe_values.shape_grad(
j, q_index) *
fe_values.JxW(q_index));
for (const unsigned int i : fe_values.dof_indices())
cell_rhs(i) += (fe_values.shape_value(i, q_index) *
1. *
fe_values.JxW(q_index));
}
cell->get_dof_indices(local_dof_indices);
for (const unsigned int i : fe_values.dof_indices())
system_matrix.add(local_dof_indices[i],
for (const unsigned int i : fe_values.dof_indices())
}
system_matrix,
solution,
}
void Step3::solve()
{
std::cout << solver_control.last_step()
<< " CG iterations needed to obtain convergence." << std::endl;
}
{
data_out.set_flags(flags);
data_out.attach_dof_handler(dof_handler);
data_out.add_data_vector(solution, "solution");
data_out.build_patches(mapping, 2);
std::ofstream output("solution.vtk");
data_out.write_vtk(output);
}
void Step3::run()
{
solve();
}
{
return 0;
}
void make_sparsity_pattern(const DoFHandler< dim, spacedim > &dof_handler, SparsityPatternBase &sparsity_pattern, const AffineConstraints< number > &constraints={}, const bool keep_constrained_dofs=true, const types::subdomain_id subdomain_id=numbers::invalid_subdomain_id)
@ update_values
Shape function values.
@ update_JxW_values
Transformed quadrature weights.
@ update_gradients
Shape function gradients.
void cell_matrix(FullMatrix< double > &M, const FEValuesBase< dim > &fe, const FEValuesBase< dim > &fetest, const ArrayView< const std::vector< double > > &velocity, const double factor=1.)
const ::parallel::distributed::Triangulation< dim, spacedim > * triangulation
bool write_higher_order_cells
Example: mixed mesh
The following code shows how to work with mixed meshes:
#include <fstream>
#include <iostream>
{
public:
void run();
private:
void solve();
};
{}
{
std::cout <<
"Number of active cells: " <<
triangulation.n_active_cells()
<< std::endl;
}
{
for (const auto &cell : dof_handler.active_cell_iterators())
{
cell->set_active_fe_index(0);
cell->set_active_fe_index(1);
else
}
dof_handler.distribute_dofs(fe);
std::cout << "Number of degrees of freedom: " << dof_handler.n_dofs()
<< std::endl;
sparsity_pattern.copy_from(
dsp);
system_matrix.reinit(sparsity_pattern);
solution.reinit(dof_handler.n_dofs());
}
{
fe,
std::vector<types::global_dof_index> local_dof_indices;
for (const auto &cell : dof_handler.active_cell_iterators())
{
hp_fe_values.reinit(cell);
const auto &fe_values = hp_fe_values.get_present_fe_values();
const unsigned int dofs_per_cell = cell->get_fe().n_dofs_per_cell();
local_dof_indices.resize(dofs_per_cell);
for (const unsigned int q_index : fe_values.quadrature_point_indices())
{
for (const unsigned int i : fe_values.dof_indices())
(fe_values.shape_grad(i, q_index) *
fe_values.shape_grad(
j, q_index) *
fe_values.JxW(q_index));
for (const unsigned int i : fe_values.dof_indices())
cell_rhs(i) += (fe_values.shape_value(i, q_index) *
1. *
fe_values.JxW(q_index));
}
cell->get_dof_indices(local_dof_indices);
for (const unsigned int i : fe_values.dof_indices())
system_matrix.add(local_dof_indices[i],
for (const unsigned int i : fe_values.dof_indices())
}
system_matrix,
solution,
}
void Step3::solve()
{
std::cout << solver_control.last_step()
<< " CG iterations needed to obtain convergence." << std::endl;
}
{
data_out.set_flags(flags);
data_out.attach_dof_handler(dof_handler);
data_out.add_data_vector(solution, "solution");
data_out.build_patches(mapping, 2);
std::ofstream output("solution.vtk");
data_out.write_vtk(output);
}
void Step3::run()
{
solve();
}
{
return 0;
}
#define DEAL_II_NOT_IMPLEMENTED()
constexpr const ReferenceCell Quadrilateral
constexpr const ReferenceCell Triangle
Reference cells
In 2D, we provide triangles and quadrilaterals with the following possible orientations in 3D:
2D: triangle and quadrilateral
Possible orientations of triangles and quadrilaterals in 3D
In 3D, tetrahedra, pyramids, wedges, and hexahedra are available:
Each surface of a 3D reference cell consists of 2D reference cells. The documentation of the enumeration of the numbering of their vertices and lines are given in the right columns.