]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Test ScratchData get_values
authorLuca Heltai <luca.heltai@sissa.it>
Tue, 9 Apr 2019 21:05:17 +0000 (23:05 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Tue, 9 Apr 2019 21:05:17 +0000 (23:05 +0200)
tests/meshworker/scratch_data_04.cc [new file with mode: 0644]
tests/meshworker/scratch_data_04.with_muparser=on.output [new file with mode: 0644]

diff --git a/tests/meshworker/scratch_data_04.cc b/tests/meshworker/scratch_data_04.cc
new file mode 100644 (file)
index 0000000..4887755
--- /dev/null
@@ -0,0 +1,131 @@
+/* ---------------------------------------------------------------------
+ *
+ * Copyright (C) 2018 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.
+ *
+ * ---------------------------------------------------------------------
+ */
+
+// Compute integral of known function using mesh_loop and ScratchData
+
+#include <deal.II/base/function_parser.h>
+#include <deal.II/base/patterns.h>
+#include <deal.II/base/thread_management.h>
+
+#include <deal.II/dofs/dof_tools.h>
+
+#include <deal.II/fe/fe_dgq.h>
+#include <deal.II/fe/fe_q.h>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/grid_tools.h>
+#include <deal.II/grid/grid_tools_cache.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_accessor.h>
+#include <deal.II/grid/tria_iterator.h>
+
+#include <deal.II/lac/affine_constraints.h>
+#include <deal.II/lac/sparse_direct.h>
+
+#include <deal.II/meshworker/copy_data.h>
+#include <deal.II/meshworker/mesh_loop.h>
+#include <deal.II/meshworker/scratch_data.h>
+
+#include <deal.II/numerics/data_out.h>
+#include <deal.II/numerics/vector_tools.h>
+
+#include <cmath>
+#include <fstream>
+#include <iostream>
+#include <unordered_map>
+
+#include "../tests.h"
+
+using namespace MeshWorker;
+
+template <int dim, int spacedim>
+void
+test()
+{
+  Triangulation<dim, spacedim> tria;
+  FE_Q<dim, spacedim>          fe(1);
+  DoFHandler<dim, spacedim>    dh(tria);
+
+  FunctionParser<spacedim> integral_function("x");
+  GridGenerator::hyper_cube(tria);
+  tria.refine_global(5);
+  tria.execute_coarsening_and_refinement();
+  dh.distribute_dofs(fe);
+
+  Vector<double> solution(dh.n_dofs());
+
+  VectorTools::interpolate(dh, integral_function, solution);
+
+  QGauss<dim> quad(3);
+
+  UpdateFlags cell_flags = update_values | update_gradients |
+                           update_quadrature_points | update_JxW_values;
+
+  using ScratchData = ScratchData<dim, spacedim>;
+  struct CopyData
+  {};
+
+
+  double                     H1_norm = 0;
+  FEValuesExtractors::Scalar scalar(0);
+
+  ScratchData scratch(fe, quad, cell_flags);
+  CopyData    copy;
+
+  auto cell = dh.begin_active();
+  auto endc = dh.end();
+
+  typedef decltype(cell) Iterator;
+
+
+  auto cell_integrator = [&H1_norm, &solution, &scalar](const Iterator &cell,
+                                                        ScratchData &   s,
+                                                        CopyData &      c) {
+    const auto &fev = s.reinit(cell);
+    const auto &JxW = s.get_JxW_values();
+
+    s.extract_local_dof_values("solution", solution);
+
+    const auto &u      = s.get_values("solution", scalar);
+    const auto &grad_u = s.get_gradients("solution", scalar);
+
+    for (unsigned int q = 0; q < u.size(); ++q)
+      H1_norm += (u[q] * u[q] + grad_u[q] * grad_u[q]) * JxW[q];
+  };
+
+  const auto empty_copyer = [](const CopyData &) {};
+
+  mesh_loop(cell,
+            endc,
+            cell_integrator,
+            empty_copyer,
+            scratch,
+            copy,
+            assemble_own_cells);
+
+  deallog << "H1 norm: " << std::sqrt(H1_norm) << std::endl;
+}
+
+
+int
+main()
+{
+  initlog(1);
+  MultithreadInfo::set_thread_limit(1); // to make output deterministic
+
+  test<2, 2>();
+}
diff --git a/tests/meshworker/scratch_data_04.with_muparser=on.output b/tests/meshworker/scratch_data_04.with_muparser=on.output
new file mode 100644 (file)
index 0000000..b906801
--- /dev/null
@@ -0,0 +1,2 @@
+
+DEAL::H1 norm: 1.15470

In the beginning the Universe was created. This has made a lot of people very angry and has been widely regarded as a bad move.

Douglas Adams


Typeset in Trocchi and Trocchi Bold Sans Serif.