From: Maximilian Bergbauer Date: Wed, 23 Apr 2025 08:15:38 +0000 (+0200) Subject: Add test X-Git-Url: https://gitweb.dealii.org/cgi-bin/gitweb.cgi?a=commitdiff_plain;h=7816deeb8e8912331abd005730fb32743ba1027b;p=dealii.git Add test --- diff --git a/tests/matrix_free/point_evaluation_37.cc b/tests/matrix_free/point_evaluation_37.cc new file mode 100644 index 0000000000..04a64c4b8a --- /dev/null +++ b/tests/matrix_free/point_evaluation_37.cc @@ -0,0 +1,147 @@ +// ------------------------------------------------------------------------ +// +// SPDX-License-Identifier: LGPL-2.1-or-later +// Copyright (C) 2020 - 2024 by the deal.II authors +// +// This file is part of the deal.II library. +// +// Part of the source code is dual licensed under Apache-2.0 WITH +// LLVM-exception OR LGPL-2.1-or-later. Detailed license information +// governing the source code and code contributions can be found in +// LICENSE.md and CONTRIBUTING.md at the top level directory of deal.II. +// +// ------------------------------------------------------------------------ + + +// check FEPointEvaluation lexicographic numbering + +#include + +#include +#include +#include + +#include + +#include +#include +#include + +#include +#include + +#include + +#include "../tests.h" + +using namespace dealii; + +int +main(int argc, char *argv[]) +{ + initlog(); + + constexpr unsigned int dim = 3; + constexpr unsigned int degree = 2; + constexpr unsigned int components = 2; + + using Number = double; + constexpr unsigned int n_lanes = VectorizedArray::size(); + + const unsigned int cell_batch_index = 0; + + Triangulation triangulation; + MappingQ mapping(1); + + GridGenerator::subdivided_hyper_cube(triangulation, 4, -1, 1); + + auto dof_handler = DoFHandler(triangulation); + + FESystem fe_system(FE_Q(degree), components); + + dof_handler.distribute_dofs(fe_system); + + auto matrix_free = MatrixFree>(); + typename MatrixFree>::AdditionalData + additional_data; + additional_data.mapping_update_flags = + update_gradients | update_values | update_quadrature_points; + additional_data.hold_all_faces_to_owned_cells = true; + + AffineConstraints affine_constraints{}; + QGauss quadrature(degree + 1); + std::vector> normals(quadrature.size()); + for (auto &normal : normals) + { + const Number x = 1. / ((dim == 3) ? std::sqrt(3.) : std::sqrt(2.)); + normal[0] = x; + normal[1] = x; + if (dim == 3) + normal[2] = x; + } + NonMatching::ImmersedSurfaceQuadrature surface_quadrature( + quadrature.get_points(), quadrature.get_weights(), normals); + + // setup matrixfree object + matrix_free.reinit( + mapping, dof_handler, affine_constraints, quadrature, additional_data); + + NonMatching::MappingInfo> mapping_info( + mapping, + additional_data.mapping_update_flags | update_normal_vectors | + update_JxW_values); + + std::vector>> cell_vec( + 1, triangulation.begin_active()); + std::vector> quad_vec( + 1, surface_quadrature); + mapping_info.reinit_surface(cell_vec, quad_vec); + + FEPointEvaluation> + evaluator_surface(mapping_info, fe_system, 0, true); + + FEEvaluation> + evaluator(matrix_free, 0, 0); + + LinearAlgebra::distributed::Vector src; + matrix_free.initialize_dof_vector(src); + for (unsigned int i = 0; i < src.size(); i++) + src[i] = random_value(); + + const auto dofs_per_cell = evaluator.dofs_per_cell; + + evaluator.reinit(cell_batch_index); + evaluator.read_dof_values(src); + + for (unsigned int v = 0; v < 1; ++v) + { + evaluator_surface.reinit(0); + + StridedArrayView strided_view( + &evaluator.begin_dof_values()[0][v], dofs_per_cell); + + evaluator_surface.evaluate(StridedArrayView( + &evaluator.begin_dof_values()[0][v], + dofs_per_cell), + EvaluationFlags::gradients); + + for (const auto q : evaluator_surface.quadrature_point_indices()) + { + evaluator_surface.submit_normal_derivative( + evaluator_surface.get_normal_derivative(q), q); + } + + evaluator_surface.test_and_sum( + StridedArrayView(&evaluator.begin_dof_values()[0][v], + dofs_per_cell), + EvaluationFlags::gradients); + + const unsigned int dofs_per_component = evaluator.dofs_per_component; + for (unsigned int c = 0; c < components; ++c) + { + for (unsigned int i = 0; i < dofs_per_component; ++i) + deallog << strided_view[c * dofs_per_component + i] << ' '; + deallog << std::endl; + } + } +} diff --git a/tests/matrix_free/point_evaluation_37.output b/tests/matrix_free/point_evaluation_37.output new file mode 100644 index 0000000000..ded0e4291a --- /dev/null +++ b/tests/matrix_free/point_evaluation_37.output @@ -0,0 +1,3 @@ + +DEAL::15.3798 -2.49427 -4.00648 10.9579 -23.2549 0.836823 -1.63726 -3.53717 3.13241 6.43517 -8.00859 6.13889 -13.7625 -33.3275 18.8530 -3.74731 14.8151 6.72491 -3.41090 5.34329 -2.13885 -0.168920 8.65233 -10.3875 1.26111 4.07237 7.27893 +DEAL::-6.30453 -17.4819 -2.01155 -9.86062 4.73463 3.69133 -1.32309 4.30278 -1.55771 1.06091 7.60697 3.98305 10.4435 17.9998 8.97276 2.12405 -3.88467 -5.69541 3.30566 2.26495 -0.945841 -0.205502 -21.0919 -6.64979 -2.28229 2.62285 6.18152