From 007a87a2947a49310596335e7fe641f2b999ded4 Mon Sep 17 00:00:00 2001 From: buerg Date: Wed, 3 Jul 2013 22:40:19 +0000 Subject: [PATCH] Add test for grid with nonstandard orientated faces for FE_Nedelec<3>. By Alexander Grayver. git-svn-id: https://svn.dealii.org/trunk@29932 0785d39b-7218-0410-832d-ea1e28bc413d --- tests/fe/fe_face_orientation_nedelec_0.cc | 176 ++++++++++++++++++ .../fe_face_orientation_nedelec_0/cmp/generic | 65 +++++++ 2 files changed, 241 insertions(+) create mode 100644 tests/fe/fe_face_orientation_nedelec_0.cc create mode 100644 tests/fe/fe_face_orientation_nedelec_0/cmp/generic diff --git a/tests/fe/fe_face_orientation_nedelec_0.cc b/tests/fe/fe_face_orientation_nedelec_0.cc new file mode 100644 index 0000000000..287784b0d7 --- /dev/null +++ b/tests/fe/fe_face_orientation_nedelec_0.cc @@ -0,0 +1,176 @@ +//---------------------------------------------------------------------- +// $Id: fe_face_orientation_nedelec.h 29507 2013-05-14 14:11:09Z kronbichler $ +// Version: $Name$ +// +// Copyright (C) 2013 by the deal.II authors +// +// This file is subject to QPL and may not be distributed +// without copyright and license information. Please refer +// to the file deal.II/doc/license.html for the text and +// further information on this license. +// +//---------------------------------------------------------------------- + +// Test FE_Nedelec<3> for meshes with faces with non-standard orientation. + +#include "../tests.h" +#include +#include +#include +#include +#include +#include +#include +#include +#include +#include + +void create_reference_triangulation (Triangulation<3>& triangulation) +{ + GridGenerator::hyper_cube (triangulation, -1.0, 1.0); + triangulation.refine_global(1); +} + +void create_triangulation (Triangulation<3>& triangulation) +{ + static const Point<3> vertices_parallelograms[] = { + Point<3> (-1., -1., -1.), //0 + Point<3> (0., -1., -1.), + Point<3> (1., -1., -1.), + + Point<3> (-1., -1., 0.), //3 + Point<3> (0., -1., 0.), + Point<3> (1., -1., 0.), + + Point<3> (-1., -1., 1.), //6 + Point<3> (0., -1., 1.), + Point<3> (1., -1., 1.), + + Point<3> (-1., 0., -1.), //9 + Point<3> (0., 0., -1.), + Point<3> (1., 0., -1.), + + Point<3> (-1., 0., 0.), //12 + Point<3> (0., 0., 0.), + Point<3> (1., 0., 0.), + + Point<3> (-1., 0., 1.), //15 + Point<3> (0., 0., 1.), + Point<3> (1., 0., 1.), + + Point<3> (-1., 1., -1.), //18 + Point<3> (0., 1., -1.), + Point<3> (1., 1., -1.), + + Point<3> (-1., 1., 0.), //21 + Point<3> (0., 1., 0.), + Point<3> (1., 1., 0.), + + Point<3> (-1., 1., 1.), //24 + Point<3> (0., 1., 1.), + Point<3> (1., 1., 1.) + }; + const unsigned n_vertices = sizeof(vertices_parallelograms) / sizeof(vertices_parallelograms[0]); + + const unsigned n_cells = 8; + + const std::vector > vertices (&vertices_parallelograms[0], + &vertices_parallelograms[n_vertices]); + + // create grid with all possible combintations of face_flip, face_orientation and face_rotation flags + static const int cell_vertices[][GeometryInfo<3>::vertices_per_cell] = { + {0, 1, 9, 10, 3, 4, 12, 13}, // cell 1 standard + {10, 11, 13, 14, 1, 2, 4, 5}, // cell 2 rotated by 270 deg + {9, 10, 18, 19, 12, 13, 21, 22}, // cell 3 standard + {13, 14, 10, 11, 22, 23, 19, 20}, // cell 4 rotated by 90 deg + {3, 4, 12, 13, 6, 7, 15, 16}, // cell 5 standard + {4, 5, 13, 14, 7, 8, 16, 17}, // cell 6 standard + {24, 25, 15, 16, 21, 22, 12, 13}, // cell 7 rotated by 180 deg + {13, 14, 22, 23, 16, 17, 25, 26} // cell 8 standard + }; + + std::vector > cells (n_cells, CellData<3>()); + for (unsigned i = 0; i::vertices_per_cell; ++j) + cells[i].vertices[j] = cell_vertices[i][j]; + cells[i].material_id = 0; + } + + triangulation.create_triangulation (vertices, cells, SubCellData()); +} + +void evaluate (const FE_Nedelec<3>& fe, const DoFHandler<3>& dof_handler_ref, const Vector& u_ref, const DoFHandler<3>& dof_handler, const Vector& u) +{ + const FEValuesExtractors::Vector component (0); + const QGauss<3> quadrature (2); + const unsigned int n_q_points = quadrature.size (); + Functions::FEFieldFunction<3> fe_field_function (dof_handler, u); + FEValues<3> fe_values (fe, quadrature, update_quadrature_points | update_values); + std::vector > values (n_q_points, Vector (3)); + std::vector > values_ref (n_q_points); + + for (DoFHandler<3>::active_cell_iterator cell = dof_handler_ref.begin_active (); cell != dof_handler_ref.end (); ++cell) + { + fe_values.reinit (cell); + fe_values[component].get_function_values (u_ref, values_ref); + fe_field_function.vector_value_list (fe_values.get_quadrature_points (), values); + std::vector dof_indices (fe.dofs_per_cell);cell->get_dof_indices (dof_indices); + for (unsigned int q_point = 0; q_point < n_q_points; ++q_point) + { + for (unsigned int d = 0; d < 3; ++d) + deallog << values_ref[q_point][d] - values[q_point] (d) << " "; + + deallog << std::endl; + } + } +} + +void set_reference_solution (Vector& vector) +{ + for (unsigned int i = 0; i < vector.size (); ++i) + vector (i) = 1.0; +} + +void set_solution (Vector& vector, const DoFHandler<3>& dof_handler, const DoFHandler<3>& dof_handler_ref, const Vector& u_ref) +{ + ConstraintMatrix constraints; + + constraints.close (); + + Functions::FEFieldFunction<3> fe_field_function (dof_handler_ref, u_ref); + + VectorTools::project (dof_handler, constraints, QGauss<3> (2), fe_field_function, vector); +} + +int main () +{ + initlog (__FILE__); + deallog.threshold_double (1.e-10); + + Triangulation<3> tria_ref; + + create_reference_triangulation (tria_ref); + + FE_Nedelec<3> fe (0); + DoFHandler<3> dof_handler_ref (tria_ref); + + dof_handler_ref.distribute_dofs (fe); + + Vector u_ref (dof_handler_ref.n_dofs ()); + + set_reference_solution (u_ref); + + Triangulation<3> tria; + + create_triangulation (tria); + + DoFHandler<3> dof_handler (tria); + + dof_handler.distribute_dofs (fe); + + Vector u (dof_handler.n_dofs ()); + + set_solution (u, dof_handler, dof_handler_ref, u_ref); + evaluate (fe, dof_handler_ref, u_ref, dof_handler, u); +} diff --git a/tests/fe/fe_face_orientation_nedelec_0/cmp/generic b/tests/fe/fe_face_orientation_nedelec_0/cmp/generic new file mode 100644 index 0000000000..6453add255 --- /dev/null +++ b/tests/fe/fe_face_orientation_nedelec_0/cmp/generic @@ -0,0 +1,65 @@ + +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 +DEAL::0 0 0 -- 2.39.5