/**
* Return the value of a finite element function at quadrature point
- * number @p q_point after a call to @p evalue(true,...).
+ * number @p q_point after a call to @p evaluate(true,...).
*/
__device__ value_type
get_value(const unsigned int q_point) const;
+ /**
+ * Return the value of a finite element function at degree of freedom
+ * @p dof after a call to integrate() or before a call to evaluate().
+ */
+ __device__ value_type
+ get_dof_value(const unsigned int dof) const;
+
/**
* Write a value to the field containing the values on quadrature points
* with component @p q_point. Access to the same fields as through @p
- * get_value(), This specifies the value which is tested by all basis
+ * get_value(). This specifies the value which is tested by all basis
* function on the current cell and integrated over.
*/
__device__ void
submit_value(const value_type &val_in, const unsigned int q_point);
+ /**
+ * Write a value to the field containing the values for the degree of
+ * freedom with index @p dof after a call to integrate() or before
+ * calling evaluate(). Access through the same fields as through
+ * get_dof_value().
+ */
+ __device__ void
+ submit_dof_value(const value_type &val_in, const unsigned int dof);
+
/**
* Return the gradient of a finite element function at quadrature point
* number @p q_point after a call to @p evaluate(...,true).
+ template <int dim,
+ int fe_degree,
+ int n_q_points_1d,
+ int n_components_,
+ typename Number>
+ __device__ typename FEEvaluation<dim,
+ fe_degree,
+ n_q_points_1d,
+ n_components_,
+ Number>::value_type
+ FEEvaluation<dim, fe_degree, n_q_points_1d, n_components_, Number>::
+ get_dof_value(const unsigned int dof) const
+ {
+ return values[dof];
+ }
+
+
+
template <int dim,
int fe_degree,
int n_q_points_1d,
+ template <int dim,
+ int fe_degree,
+ int n_q_points_1d,
+ int n_components_,
+ typename Number>
+ __device__ void
+ FEEvaluation<dim, fe_degree, n_q_points_1d, n_components_, Number>::
+ submit_dof_value(const value_type &val_in, const unsigned int dof)
+ {
+ values[dof] = val_in;
+ }
+
+
+
template <int dim,
int fe_degree,
int n_q_points_1d,
--- /dev/null
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2019 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.md at
+// the top level directory of deal.II.
+//
+// ---------------------------------------------------------------------
+
+
+// check that CUDAWrappers::FEEvaluation::submit_dof_value/get_dof_value
+// works correctly.
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/manifold_lib.h>
+
+#include <deal.II/matrix_free/cuda_fe_evaluation.h>
+
+#include "../tests.h"
+
+// forward declare this function
+template <int dim, int fe_degree>
+void
+test();
+
+template <int dim,
+ int fe_degree,
+ int n_q_points_1d = fe_degree + 1,
+ typename Number = double>
+class MatrixFreeTest
+{
+public:
+ static const unsigned int n_dofs_1d = fe_degree + 1;
+ static const unsigned int n_local_dofs = Utilities::pow(n_dofs_1d, dim);
+ static const unsigned int n_q_points = Utilities::pow(n_q_points_1d, dim);
+
+ MatrixFreeTest(const CUDAWrappers::MatrixFree<dim, Number> &data_in)
+ : data(data_in){};
+
+ __device__ void
+ operator()(
+ const unsigned int cell,
+ const typename CUDAWrappers::MatrixFree<dim, Number>::Data *gpu_data,
+ CUDAWrappers::SharedData<dim, Number> * shared_data,
+ const Number * src,
+ Number * dst) const
+ {
+ CUDAWrappers::FEEvaluation<dim, fe_degree, n_q_points_1d, 1, Number>
+ fe_eval(cell, gpu_data, shared_data);
+
+ const unsigned int tid =
+ (threadIdx.x % n_q_points_1d) +
+ (dim > 1 ? threadIdx.y : 0) * n_q_points_1d +
+ (dim > 2 ? threadIdx.z : 0) * n_q_points_1d * n_q_points_1d;
+
+ bool const evaluate_values = true;
+ bool const evaluate_gradients = true;
+
+ // set to unit vector
+ fe_eval.submit_dof_value(1., tid);
+ __syncthreads();
+ fe_eval.evaluate(evaluate_values, evaluate_gradients);
+
+ // values should evaluate to one, derivatives to zero
+ assert(fe_eval.get_dof_value(tid) == 1.);
+ assert(fe_eval.get_value(tid) == 1.);
+ for (unsigned int e = 0; e < dim; ++e)
+ assert(fe_eval.get_gradient(tid)[e] == 0.);
+ }
+
+
+
+ void
+ test() const
+ {
+ LinearAlgebra::distributed::Vector<Number, MemorySpace::CUDA> dst_dummy;
+ LinearAlgebra::distributed::Vector<Number, MemorySpace::CUDA> src_dummy;
+
+ data.cell_loop(*this, src_dummy, dst_dummy);
+
+ // Check that the kernel was launched correctly
+ AssertCuda(cudaGetLastError());
+ // Check that there was no problem during the execution of the kernel
+ AssertCuda(cudaDeviceSynchronize());
+ deallog << "OK" << std::endl;
+ };
+
+protected:
+ const CUDAWrappers::MatrixFree<dim, Number> &data;
+};
+
+
+
+template <int dim, int fe_degree, typename number>
+void
+do_test(const DoFHandler<dim> & dof,
+ const AffineConstraints<double> &constraints)
+{
+ CUDAWrappers::MatrixFree<dim, number> mf_data;
+ {
+ const QGauss<1> quad(fe_degree + 1);
+ typename CUDAWrappers::MatrixFree<dim, number>::AdditionalData data;
+ data.mapping_update_flags = update_values | update_gradients |
+ update_JxW_values | update_quadrature_points;
+ mf_data.reinit(dof, constraints, quad, data);
+ }
+
+ deallog << "Testing " << dof.get_fe().get_name() << std::endl;
+ MatrixFreeTest<dim, fe_degree, fe_degree + 1, number> mf(mf_data);
+ mf.test();
+}
+
+
+int
+main()
+{
+ initlog();
+ deal_II_exceptions::disable_abort_on_exception();
+
+ test<2, 1>();
+}
+
+
+template <int dim, int fe_degree>
+void
+test()
+{
+ const SphericalManifold<dim> manifold;
+ Triangulation<dim> tria;
+ GridGenerator::hyper_ball(tria);
+ for (const auto &cell : tria.active_cell_iterators())
+ for (unsigned int f = 0; f < GeometryInfo<dim>::faces_per_cell; ++f)
+ if (cell->at_boundary(f))
+ cell->face(f)->set_all_manifold_ids(0);
+ tria.set_manifold(0, manifold);
+
+ // refine first and last cell
+ tria.begin(tria.n_levels() - 1)->set_refine_flag();
+ tria.last()->set_refine_flag();
+ tria.execute_coarsening_and_refinement();
+ tria.refine_global(4 - dim);
+
+ FE_Q<dim> fe(fe_degree);
+ DoFHandler<dim> dof(tria);
+ dof.distribute_dofs(fe);
+
+ AffineConstraints<double> constraints;
+ DoFTools::make_hanging_node_constraints(dof, constraints);
+ constraints.close();
+
+ do_test<dim, fe_degree, double>(dof, constraints);
+}
projects / dealii.git / commitdiff