]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Added default constructor to Quadrature 2489/head
authordanshapero <shapero.daniel@gmail.com>
Thu, 7 Apr 2016 19:47:11 +0000 (12:47 -0700)
committerdanshapero <shapero.daniel@gmail.com>
Sat, 9 Apr 2016 00:40:57 +0000 (17:40 -0700)
doc/news/changes.h
include/deal.II/base/quadrature.h
include/deal.II/base/quadrature_lib.h
tests/base/quadrature_move.cc [new file with mode: 0644]
tests/base/quadrature_move.with_cxx11=on.output [new file with mode: 0644]

index 241d096123b6838ea70db1f9a636541bd82102a2..c3b790306ac34434ac8f0168413d01280b8710f5 100644 (file)
@@ -116,7 +116,12 @@ inconvenience this causes.
  <br>
  (Daniel Arndt, 2016/04/08)
  </li>
+
+ <li> New: A move constructor has been added to Quadrature.
+ <br>
+ (Daniel Shapero, 2016/04/08)
+ </li>
+
  <li> Fixed: The multigrid transfer performed invalid data accesses on
  multigrid hierarchies that define the coarse level as a level larger than
  0. This has been fixed.
index 227fb4f75c638ed1f34212ac61b6b85e78c0a146..1a1413a8aac29e6bd5de6cec6a204497327db164 100644 (file)
@@ -125,6 +125,17 @@ public:
    */
   Quadrature (const Quadrature<dim> &q);
 
+#ifdef DEAL_II_WITH_CXX11
+  /**
+   * Move constructor. Construct a new quadrature object by transferring the
+   * internal data of another quadrature object.
+   *
+   * @note this constructor is only available if deal.II is configured with
+   * C++11 support.
+   */
+  Quadrature (Quadrature<dim> &&) = default;
+#endif
+
   /**
    * Construct a quadrature formula from given vectors of quadrature points
    * (which should really be in the unit cell) and the corresponding weights.
index 1c39bbeed11a524e7036a1805f5d396e2adbe97d..96e3950946b93d147c05b544c3a30616397eb03a 100644 (file)
@@ -294,6 +294,14 @@ public:
              const double alpha = 1,
              const bool factor_out_singular_weight=false);
 
+#ifdef DEAL_II_WITH_CXX11
+  /**
+   * Move constructor. We cannot rely on the move constructor for `Quadrature`,
+   * since it does not know about the additional member `fraction` of this class.
+   */
+  QGaussLogR(QGaussLogR<dim> &&) = default;
+#endif
+
 protected:
   /**
    * This is the length of interval $(0,origin)$, or 1 if either of the two
@@ -577,6 +585,14 @@ public:
   QGaussRadauChebyshev(const unsigned int n,
                        EndPoint ep=QGaussRadauChebyshev::left);
 
+#ifdef DEAL_II_WITH_CXX11
+  /**
+   * Move constructor. We cannot rely on the move constructor for `Quadrature`,
+   * since it does not know about the additional member `ep` of this class.
+   */
+  QGaussRadauChebyshev(QGaussRadauChebyshev<dim> &&) = default;
+#endif
+
 private:
   const EndPoint ep;
   /// Computes the points of the quadrature formula.
diff --git a/tests/base/quadrature_move.cc b/tests/base/quadrature_move.cc
new file mode 100644 (file)
index 0000000..eb0fec7
--- /dev/null
@@ -0,0 +1,82 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 1998 - 2015 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.
+//
+// ---------------------------------------------------------------------
+
+
+#include "../tests.h"
+#include <fstream>
+#include <cmath>
+
+#include <deal.II/base/quadrature_lib.h>
+
+
+template <class Quad, typename... Args>
+std::string check_q_move(Args &&...args)
+{
+  Quad quad1(args...);
+  const unsigned int size1 = quad1.size();
+
+  std::vector<double> weights1 = quad1.get_weights();
+
+  Quad quad2(std::move(quad1));
+  const unsigned int size2 = quad2.size();
+
+  std::vector<double> weights2 = quad2.get_weights();
+
+  if (size1 != size2) return "NOPE";
+
+  for (unsigned short i = 0; i < size1; ++i)
+    if (std::fabs(weights1[i] - weights2[i]) > 1.0e-16)
+      return "NOPE";
+
+  return "OK";
+}
+
+
+template <template <int dim> class Quad, typename... Args>
+void check_quadrature_move(Args &&...args)
+{
+  deallog << check_q_move<Quad<1>>(std::forward<Args>(args)...) << 1 << " "
+          << check_q_move<Quad<2>>(std::forward<Args>(args)...) << 2 << " "
+          << check_q_move<Quad<3>>(std::forward<Args>(args)...) << 3
+          << std::endl;
+}
+
+
+int main()
+{
+  initlog();
+
+  check_quadrature_move<QMidpoint>();
+  check_quadrature_move<QTrapez>();
+  check_quadrature_move<QSimpson>();
+  check_quadrature_move<QMilne>();
+  check_quadrature_move<QWeddle>();
+
+  for (unsigned int p = 2; p < 5; ++p)
+    {
+      check_quadrature_move<QGauss>(p);
+      check_quadrature_move<QGaussLobatto>(p);
+    }
+
+  const auto ep = QGaussRadauChebyshev<1>::right;
+  for (unsigned int p = 2; p < 5; ++p)
+    {
+      deallog << "Gauss Log R: " << check_q_move<QGaussLogR<1>>(p) << std::endl;
+      deallog << "Gauss Radau Chebyshev: "
+              << check_q_move<QGaussRadauChebyshev<1>>(p, ep) << std::endl;
+    }
+
+  return 0;
+}
diff --git a/tests/base/quadrature_move.with_cxx11=on.output b/tests/base/quadrature_move.with_cxx11=on.output
new file mode 100644 (file)
index 0000000..5e60b40
--- /dev/null
@@ -0,0 +1,18 @@
+
+DEAL::OK1 OK2 OK3
+DEAL::OK1 OK2 OK3
+DEAL::OK1 OK2 OK3
+DEAL::OK1 OK2 OK3
+DEAL::OK1 OK2 OK3
+DEAL::OK1 OK2 OK3
+DEAL::OK1 OK2 OK3
+DEAL::OK1 OK2 OK3
+DEAL::OK1 OK2 OK3
+DEAL::OK1 OK2 OK3
+DEAL::OK1 OK2 OK3
+DEAL::Gauss Log R: OK
+DEAL::Gauss Radau Chebyshev: OK
+DEAL::Gauss Log R: OK
+DEAL::Gauss Radau Chebyshev: OK
+DEAL::Gauss Log R: OK
+DEAL::Gauss Radau Chebyshev: OK

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.