]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add convert_generalized_support_point_values_to_dof_values
authorPraveen C <cpraveen@gmail.com>
Sun, 19 Nov 2017 04:44:46 +0000 (10:14 +0530)
committerPraveen C <cpraveen@gmail.com>
Sun, 19 Nov 2017 04:44:46 +0000 (10:14 +0530)
This function is missing for FE_FaceQ and FE_TraceQ which is now
implemented by copying from FE_DGQ.

See https://github.com/dealii/dealii/issues/5449

include/deal.II/fe/fe_face.h
include/deal.II/fe/fe_trace.h
source/fe/fe_face.cc
source/fe/fe_trace.cc
tests/fe/fe_faceq_interpolate.cc [new file with mode: 0644]
tests/fe/fe_faceq_interpolate.debug.output [new file with mode: 0644]
tests/fe/fe_traceq_interpolate.cc [new file with mode: 0644]
tests/fe/fe_traceq_interpolate.debug.output [new file with mode: 0644]

index 98d45b16491a051d1d188d0480491f7b5da71573..b795fc1fcb1c9472bec6af7cc80ae48b69aefb1c 100644 (file)
@@ -70,6 +70,18 @@ public:
    */
   virtual std::string get_name () const;
 
+   /**
+    * Implementation of the corresponding function in the FiniteElement
+    * class.  Since the current element is interpolatory, the nodal
+    * values are exactly the support point values. Furthermore, since
+    * the current element is scalar, the support point values need to
+    * be vectors of length 1.
+    */
+   virtual
+   void
+   convert_generalized_support_point_values_to_dof_values (const std::vector<Vector<double> > &support_point_values,
+                                                           std::vector<double>                &nodal_values) const;
+
   /**
    * Return the matrix interpolating from a face of of one element to the face
    * of the neighboring element.  The size of the matrix is then
index f4cc2237605ef3418580ef71afee30ea5710aa26..fbadc0ffa4c710d1db9ec1acfecece2461b194eb 100644 (file)
@@ -62,6 +62,18 @@ public:
   std::unique_ptr<FiniteElement<dim,spacedim> >
   clone() const;
 
+   /**
+    * Implementation of the corresponding function in the FiniteElement
+    * class.  Since the current element is interpolatory, the nodal
+    * values are exactly the support point values. Furthermore, since
+    * the current element is scalar, the support point values need to
+    * be vectors of length 1.
+    */
+  virtual
+  void
+  convert_generalized_support_point_values_to_dof_values (const std::vector<Vector<double> > &support_point_values,
+                                                          std::vector<double>                &nodal_values) const;
+
   /**
    * This function returns @p true, if the shape function @p shape_index has
    * non-zero function values somewhere on the face @p face_index.
index 448d2b3b6c1951f9fcd853055affbb6f5b49b91f..3c1cb96a296fee6d410e47bfa316bad70c9fdf32 100644 (file)
@@ -315,7 +315,26 @@ FE_FaceQ<dim,spacedim>::get_constant_modes () const
          (constant_modes, std::vector<unsigned int>(1, 0));
 }
 
+template <int dim, int spacedim>
+void
+FE_FaceQ<dim,spacedim>::
+convert_generalized_support_point_values_to_dof_values (const std::vector<Vector<double> > &support_point_values,
+                                                        std::vector<double>                &nodal_values) const
+{
+   AssertDimension (support_point_values.size(),
+                    this->get_unit_support_points().size());
+   AssertDimension (support_point_values.size(),
+                    nodal_values.size());
+   AssertDimension (this->dofs_per_cell,
+                    nodal_values.size());
+
+   for (unsigned int i=0; i<this->dofs_per_cell; ++i)
+   {
+      AssertDimension (support_point_values[i].size(), 1);
 
+      nodal_values[i] = support_point_values[i](0);
+   }
+}
 
 // ----------------------------- FE_FaceQ<1,spacedim> ------------------------
 
index dee1e7f8699906f5c4af6842016293eb1ad28ef4..4bc1bd659404eddde42678b1f4354389e22c15ac 100644 (file)
@@ -122,6 +122,26 @@ FE_TraceQ<dim,spacedim>::get_constant_modes () const
          (constant_modes, std::vector<unsigned int>(1, 0));
 }
 
+template <int dim, int spacedim>
+void
+FE_TraceQ<dim,spacedim>::
+convert_generalized_support_point_values_to_dof_values (const std::vector<Vector<double> > &support_point_values,
+                                                        std::vector<double>                &nodal_values) const
+{
+   AssertDimension (support_point_values.size(),
+                    this->get_unit_support_points().size());
+   AssertDimension (support_point_values.size(),
+                    nodal_values.size());
+   AssertDimension (this->dofs_per_cell,
+                    nodal_values.size());
+
+   for (unsigned int i=0; i<this->dofs_per_cell; ++i)
+   {
+      AssertDimension (support_point_values[i].size(), 1);
+
+      nodal_values[i] = support_point_values[i](0);
+   }
+}
 
 
 template <int dim, int spacedim>
diff --git a/tests/fe/fe_faceq_interpolate.cc b/tests/fe/fe_faceq_interpolate.cc
new file mode 100644 (file)
index 0000000..2f795a6
--- /dev/null
@@ -0,0 +1,55 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2015 - 2017 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 at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// Call VectorTools::interpolate for a function in FE_TraceQ. The purpose is
+// to check that all functions needed for interpolate to work exist.
+
+#include "../tests.h"
+
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/fe/fe_face.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/numerics/vector_tools.h>
+
+using namespace dealii;
+
+template <int dim>
+void test()
+{
+   Triangulation<dim>    triangulation;
+   FE_FaceQ<dim>         fe(2);
+   DoFHandler<dim>       dof_handler(triangulation);
+
+   GridGenerator::hyper_cube (triangulation, 0, 1);
+   triangulation.refine_global(6);
+
+   dof_handler.distribute_dofs (fe);
+   Vector<double> solution(dof_handler.n_dofs());
+
+   VectorTools::interpolate(dof_handler,
+                            ZeroFunction<dim>(),
+                            solution);
+   deallog << "Success, dim = " << dim << std::endl;
+}
+
+int main()
+{
+   initlog();
+
+   test<2>();
+   test<3>();
+   return 0;
+}
diff --git a/tests/fe/fe_faceq_interpolate.debug.output b/tests/fe/fe_faceq_interpolate.debug.output
new file mode 100644 (file)
index 0000000..472a129
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL::Success, dim = 2
+DEAL::Success, dim = 3
diff --git a/tests/fe/fe_traceq_interpolate.cc b/tests/fe/fe_traceq_interpolate.cc
new file mode 100644 (file)
index 0000000..bfca5bf
--- /dev/null
@@ -0,0 +1,55 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2015 - 2017 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 at
+// the top level of the deal.II distribution.
+//
+// ---------------------------------------------------------------------
+
+// Call VectorTools::interpolate for a function in FE_TraceQ. The purpose is
+// to check that all functions needed for interpolate to work exist.
+
+#include "../tests.h"
+
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/fe/fe_trace.h>
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/numerics/vector_tools.h>
+
+using namespace dealii;
+
+template <int dim>
+void test()
+{
+   Triangulation<dim>    triangulation;
+   FE_TraceQ<dim>         fe(2);
+   DoFHandler<dim>       dof_handler(triangulation);
+
+   GridGenerator::hyper_cube (triangulation, 0, 1);
+   triangulation.refine_global(6);
+
+   dof_handler.distribute_dofs (fe);
+   Vector<double> solution(dof_handler.n_dofs());
+
+   VectorTools::interpolate(dof_handler,
+                            ZeroFunction<dim>(),
+                            solution);
+   deallog << "Success, dim = " << dim << std::endl;
+}
+
+int main()
+{
+   initlog();
+
+   test<2>();
+   test<3>();
+   return 0;
+}
diff --git a/tests/fe/fe_traceq_interpolate.debug.output b/tests/fe/fe_traceq_interpolate.debug.output
new file mode 100644 (file)
index 0000000..472a129
--- /dev/null
@@ -0,0 +1,3 @@
+
+DEAL::Success, dim = 2
+DEAL::Success, dim = 3

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.