From 6673f4d22a89bf26bc6b6acb84343b83d12e0dcb Mon Sep 17 00:00:00 2001 From: Martin Kronbichler Date: Tue, 31 May 2022 21:23:47 +0200 Subject: [PATCH] Reduce run time of matrix-free RT tests --- tests/matrix_free/matrix_vector_rt_common.h | 36 ++++++++++--------- .../matrix_vector_rt_face_common.h | 34 +++++++++--------- 2 files changed, 36 insertions(+), 34 deletions(-) diff --git a/tests/matrix_free/matrix_vector_rt_common.h b/tests/matrix_free/matrix_vector_rt_common.h index 9f2ee3d0e2..90ea0b41c8 100644 --- a/tests/matrix_free/matrix_vector_rt_common.h +++ b/tests/matrix_free/matrix_vector_rt_common.h @@ -160,13 +160,14 @@ protected: }; -template +template void do_test(const DoFHandler & dof, const AffineConstraints &constraints, const TestType test_type) { deallog << "Testing " << enum_to_string(test_type) << std::endl; + const unsigned int fe_degree = dof.get_fe().degree - 1; // constraints.distribute(solution); MatrixFree mf_data; @@ -192,7 +193,7 @@ do_test(const DoFHandler & dof, } // MatrixFree solution - MatrixFreeTest mf(mf_data, test_type); + MatrixFreeTest mf(mf_data, test_type); mf.test_functions(solution, initial_condition); @@ -212,6 +213,9 @@ do_test(const DoFHandler & dof, std::vector local_dof_indices(dofs_per_cell); FullMatrix local_matrix(dofs_per_cell, dofs_per_cell); + std::vector> phi_val(dofs_per_cell); + std::vector> phi_grad(dofs_per_cell); + std::vector phi_div(dofs_per_cell); const FEValuesExtractors::Vector velocities(0); @@ -222,29 +226,27 @@ do_test(const DoFHandler & dof, local_matrix = 0; for (const auto q : fe_val.quadrature_point_indices()) - for (unsigned int i = 0; i < dofs_per_cell; ++i) - { - const Tensor<1, dim> phi_i = fe_val[velocities].value(i, q) * 10.; - const Tensor<2, dim> grad_phi_i = fe_val[velocities].gradient(i, q); - const Number div_phi_i = fe_val[velocities].divergence(i, q); - + { + const double JxW = fe_val.JxW(q); + for (unsigned int i = 0; i < dofs_per_cell; ++i) + { + phi_val[i] = fe_val[velocities].value(i, q); + phi_grad[i] = fe_val[velocities].gradient(i, q); + phi_div[i] = fe_val[velocities].divergence(i, q); + } + for (unsigned int i = 0; i < dofs_per_cell; ++i) for (unsigned int j = 0; j < dofs_per_cell; ++j) { - const Tensor<1, dim> phi_j = fe_val[velocities].value(j, q); - const Tensor<2, dim> grad_phi_j = - fe_val[velocities].gradient(j, q); - const Number div_phi_j = fe_val[velocities].divergence(j, q); - if (test_type < TestType::gradients) - local_matrix(i, j) += phi_j * phi_i * fe_val.JxW(q); + local_matrix(i, j) += 10. * (phi_val[j] * phi_val[i]) * JxW; if (test_type == TestType::gradients || test_type == TestType::values_gradients) local_matrix(i, j) += - scalar_product(grad_phi_i, grad_phi_j) * fe_val.JxW(q); + scalar_product(phi_grad[i], phi_grad[j]) * JxW; else if (test_type == TestType::divergence) - local_matrix(i, j) += div_phi_i * div_phi_j * fe_val.JxW(q); + local_matrix(i, j) += phi_div[i] * phi_div[j] * JxW; } - } + } cell->get_dof_indices(local_dof_indices); for (unsigned int i = 0; i < dofs_per_cell; ++i) for (unsigned int j = 0; j < dofs_per_cell; ++j) diff --git a/tests/matrix_free/matrix_vector_rt_face_common.h b/tests/matrix_free/matrix_vector_rt_face_common.h index b9a1cd1cef..dd65b7f391 100644 --- a/tests/matrix_free/matrix_vector_rt_face_common.h +++ b/tests/matrix_free/matrix_vector_rt_face_common.h @@ -252,7 +252,7 @@ do_test(const DoFHandler & dof, initial_condition[i] = random_value(); } - MatrixFreeTest mf(mf_data, test_type); + MatrixFreeTest mf(mf_data, test_type); mf.test_functions(solution, initial_condition); @@ -273,6 +273,9 @@ do_test(const DoFHandler & dof, std::vector local_dof_indices(dofs_per_cell); FullMatrix local_matrix(dofs_per_cell, dofs_per_cell); + std::vector> phi_val(dofs_per_cell); + std::vector> phi_grad(dofs_per_cell); + std::vector phi_div(dofs_per_cell); const FEValuesExtractors::Vector velocities(0); // Assemble matrix @@ -283,30 +286,27 @@ do_test(const DoFHandler & dof, local_matrix = 0; for (const auto q : fe_val.quadrature_point_indices()) - for (unsigned int i = 0; i < dofs_per_cell; ++i) - { - const Tensor<1, dim> phi_i = fe_val[velocities].value(i, q) * 10.; - const Tensor<2, dim> grad_phi_i = - fe_val[velocities].gradient(i, q); - const Number div_phi_i = fe_val[velocities].divergence(i, q); - + { + const double JxW = fe_val.JxW(q); + for (unsigned int i = 0; i < dofs_per_cell; ++i) + { + phi_val[i] = fe_val[velocities].value(i, q); + phi_grad[i] = fe_val[velocities].gradient(i, q); + phi_div[i] = fe_val[velocities].divergence(i, q); + } + for (unsigned int i = 0; i < dofs_per_cell; ++i) for (unsigned int j = 0; j < dofs_per_cell; ++j) { - const Tensor<1, dim> phi_j = fe_val[velocities].value(j, q); - const Tensor<2, dim> grad_phi_j = - fe_val[velocities].gradient(j, q); - const Number div_phi_j = fe_val[velocities].divergence(j, q); - if (test_type < TestType::gradients) - local_matrix(i, j) += phi_j * phi_i * fe_val.JxW(q); + local_matrix(i, j) += 10. * (phi_val[j] * phi_val[i]) * JxW; if (test_type == TestType::gradients || test_type == TestType::values_gradients) local_matrix(i, j) += - scalar_product(grad_phi_i, grad_phi_j) * fe_val.JxW(q); + scalar_product(phi_grad[i], phi_grad[j]) * JxW; else if (test_type == TestType::divergence) - local_matrix(i, j) += div_phi_i * div_phi_j * fe_val.JxW(q); + local_matrix(i, j) += phi_div[i] * phi_div[j] * JxW; } - } + } cell->get_dof_indices(local_dof_indices); for (unsigned int i = 0; i < dofs_per_cell; ++i) for (unsigned int j = 0; j < dofs_per_cell; ++j) -- 2.39.5