]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Add another test for complex-valued DerivativeForm objects. 2797/head
authorWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 11 Jul 2016 11:56:23 +0000 (06:56 -0500)
committerWolfgang Bangerth <bangerth@math.tamu.edu>
Mon, 11 Jul 2016 11:56:38 +0000 (06:56 -0500)
tests/base/derivative_forms_03.cc [new file with mode: 0644]
tests/base/derivative_forms_03.output [new file with mode: 0644]

diff --git a/tests/base/derivative_forms_03.cc b/tests/base/derivative_forms_03.cc
new file mode 100644 (file)
index 0000000..8bb3741
--- /dev/null
@@ -0,0 +1,62 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2000 - 2016 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.
+//
+// ---------------------------------------------------------------------
+
+// exercise DerivativeForm::norm() for complex data types
+
+
+#include "../tests.h"
+#include <deal.II/base/derivative_form.h>
+
+
+template<int dim, int spacedim>
+void test()
+{
+  DerivativeForm<1,dim,spacedim,std::complex<double> > dF;
+  double dF_norm_sqr = 0;
+  for (unsigned int i=0; i<spacedim; ++i)
+    for (unsigned int j=0; j<dim; ++j)
+      {
+        dF[i][j] = std::complex<double>(i+2*j+1,i+2*j+1);
+        dF_norm_sqr += (i+2*j+1)*(i+2*j+1)*2;
+      }
+
+  DerivativeForm<2,dim,spacedim,std::complex<double> > ddF;
+  double ddF_norm_sqr = 0;
+  for (unsigned int i=0; i<spacedim; ++i)
+    for (unsigned int j=0; j<dim; ++j)
+      for (unsigned int k=0; k<dim; ++k)
+        {
+          ddF[i][j][k] = std::complex<double>(i+2*j+3*k+1,i+2*j+3*k+1);
+          ddF_norm_sqr += (i+2*j+3*k+1)*(i+2*j+3*k+1)*2;
+        }
+
+  // output the norms of these objects
+  deallog << "||dF||: " << dF.norm() << std::endl;
+  deallog << "||ddF||: " << ddF.norm() << std::endl;
+
+  Assert (std::fabs(dF.norm()-std::sqrt(dF_norm_sqr)) < 1e-12, ExcInternalError());
+  Assert (std::fabs(ddF.norm()-std::sqrt(ddF_norm_sqr)) < 1e-12, ExcInternalError());
+}
+
+int main()
+{
+  initlog();
+  test<1,1>();
+  test<1,2>();
+  test<1,3>();
+  test<2,2>();
+  test<2,3>();
+  test<3,3>();
+}
diff --git a/tests/base/derivative_forms_03.output b/tests/base/derivative_forms_03.output
new file mode 100644 (file)
index 0000000..070de4f
--- /dev/null
@@ -0,0 +1,13 @@
+
+DEAL::||dF||: 1.41421
+DEAL::||ddF||: 1.41421
+DEAL::||dF||: 3.16228
+DEAL::||ddF||: 3.16228
+DEAL::||dF||: 5.29150
+DEAL::||ddF||: 5.29150
+DEAL::||dF||: 7.74597
+DEAL::||ddF||: 17.6635
+DEAL::||dF||: 11.3137
+DEAL::||ddF||: 24.0832
+DEAL::||dF||: 18.6548
+DEAL::||ddF||: 56.1249

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.