]> https://gitweb.dealii.org/ - dealii.git/commitdiff
1. Bug fix for wrong sign of curl integrators in integrators/maxwell.h: nitsche_curl_... 2540/head
authorJihuan Tian <jihuan_tian@hotmail.com>
Thu, 21 Apr 2016 23:05:08 +0000 (07:05 +0800)
committerJihuan Tian <jihuan_tian@hotmail.com>
Wed, 11 May 2016 14:07:19 +0000 (22:07 +0800)
2. Add test code 'tests/integrators/maxwell_curl.cc' and related output file for verifying curl integrators in integrators/maxwell.h.
3. Add basis function support checking on faces in the function nitsche_curl_matrix in integrators/maxwell.h.

doc/news/changes.h
include/deal.II/integrators/maxwell.h
tests/integrators/maxwell_curl.cc [new file with mode: 0644]
tests/integrators/maxwell_curl.output [new file with mode: 0644]

index c62ea5a10de492e82b4912af45f7057d752e3c30..c3780ac9c32cd0b914d5b05455671897bd768769 100644 (file)
@@ -241,6 +241,14 @@ inconvenience this causes.
 <h3>Specific improvements</h3>
 
 <ol>
+ <li> Fixed: Corrected the sign of curl calculated in the functions:
+ LocalIntegrators::curl_curl_matrix, LocalIntegrators::curl_matrix,
+ LocalIntegrators::nitsche_curl_matrix and LocalIntegrators::ip_curl_matrix in
+ integrators/maxwell.h.
+ <br>
+ (Jihuan Tian, 2016/05/09)
+ </li>
+
  <li> Fixed: Bug in the RelaxationBlock class function do_step. Before, the 
  corrections were not added together, which leads to a wrong update whenever the
  Jacobi blocks are overlapping. For SOR, SSOR and non-overlapping Jacobi this was 
@@ -287,7 +295,7 @@ inconvenience this causes.
  <br>
  (Martin Kronbichler, Daniel Jodlbauer, 2016/04/21)
  </li>
-
  <li> New: Added an optional string parameter to the ParameterHandler::read_input ()
  and ParameterHandler::read_input_from_string() functions.
  When a line which equals this string is encountered, the parsing of parameters
index 64432a9ee8a065acc3bafa7be2e15647c8020096..e7e08589c5fcbb416912dc433d1e06040e28fa9f 100644 (file)
@@ -37,9 +37,9 @@ namespace LocalIntegrators
    *
    * @f[
    * \nabla\times \mathbf u = \begin{pmatrix}
-   *   \partial_3 u_2 - \partial_2 u_3 \\
-   *   \partial_1 u_3 - \partial_3 u_1 \\
-   *   \partial_2 u_1 - \partial_1 u_2
+   *   \partial_2 u_3 - \partial_3 u_2 \\
+   *   \partial_3 u_1 - \partial_1 u_3 \\
+   *   \partial_1 u_2 - \partial_2 u_1
    * \end{pmatrix}.
    * @f]
    *
@@ -49,7 +49,7 @@ namespace LocalIntegrators
    * function and the vector curl of a scalar function. The current
    * implementation exchanges the sign and we have:
    * @f[
-   *  \nabla \times \mathbf u = \partial_1 u_2 - \partial 2 u_1,
+   *  \nabla \times \mathbf u = \partial_1 u_2 - \partial_2 u_1,
    *  \qquad
    *  \nabla \times p = \begin{pmatrix}
    *    \partial_2 p \\ -\partial_1 p
@@ -197,8 +197,8 @@ namespace LocalIntegrators
                   const unsigned int d1 = (d+1)%dim;
                   const unsigned int d2 = (d+2)%dim;
 
-                  const double cv = fe.shape_grad_component(i,k,d1)[d2] - fe.shape_grad_component(i,k,d2)[d1];
-                  const double cu = fe.shape_grad_component(j,k,d1)[d2] - fe.shape_grad_component(j,k,d2)[d1];
+                  const double cv = fe.shape_grad_component(i,k,d2)[d1] - fe.shape_grad_component(i,k,d1)[d2];
+                  const double cu = fe.shape_grad_component(j,k,d2)[d1] - fe.shape_grad_component(j,k,d1)[d2];
 
                   M(i,j) += dx * cu * cv;
                 }
@@ -246,7 +246,7 @@ namespace LocalIntegrators
                   const unsigned int d2 = (d+2)%dim;
 
                   const double vv = fetest.shape_value_component(i,k,d);
-                  const double cu = fe.shape_grad_component(j,k,d1)[d2] - fe.shape_grad_component(j,k,d2)[d1];
+                  const double cu = fe.shape_grad_component(j,k,d2)[d1] - fe.shape_grad_component(j,k,d1)[d2];
                   M(i,j) += dx * cu * vv;
                 }
         }
@@ -272,6 +272,7 @@ namespace LocalIntegrators
     void nitsche_curl_matrix (
       FullMatrix<double> &M,
       const FEValuesBase<dim> &fe,
+      const unsigned int face_no,
       double penalty,
       double factor = 1.)
     {
@@ -299,17 +300,20 @@ namespace LocalIntegrators
           const Tensor<1,dim> n = fe.normal_vector(k);
           for (unsigned int i=0; i<n_dofs; ++i)
             for (unsigned int j=0; j<n_dofs; ++j)
-              for (unsigned int d=0; d<d_max; ++d)
+              if (fe.get_fe().has_support_on_face(i, face_no) && fe.get_fe().has_support_on_face(j, face_no))
                 {
-                  const unsigned int d1 = (d+1)%dim;
-                  const unsigned int d2 = (d+2)%dim;
-
-                  const double cv = fe.shape_grad_component(i,k,d1)[d2] - fe.shape_grad_component(i,k,d2)[d1];
-                  const double cu = fe.shape_grad_component(j,k,d1)[d2] - fe.shape_grad_component(j,k,d2)[d1];
-                  const double v= fe.shape_value_component(i,k,d1)*n(d2) - fe.shape_value_component(i,k,d2)*n(d1);
-                  const double u= fe.shape_value_component(j,k,d1)*n(d2) - fe.shape_value_component(j,k,d2)*n(d1);
-
-                  M(i,j) += dx*(2.*penalty*u*v - cv*u - cu*v);
+                  for (unsigned int d=0; d<d_max; ++d)
+                    {
+                      const unsigned int d1 = (d+1)%dim;
+                      const unsigned int d2 = (d+2)%dim;
+
+                      const double cv = fe.shape_grad_component(i,k,d2)[d1] - fe.shape_grad_component(i,k,d1)[d2];
+                      const double cu = fe.shape_grad_component(j,k,d2)[d1] - fe.shape_grad_component(j,k,d1)[d2];
+                      const double v= fe.shape_value_component(i,k,d1)*n[d2] - fe.shape_value_component(i,k,d2)*n[d1];
+                      const double u= fe.shape_value_component(j,k,d1)*n[d2] - fe.shape_value_component(j,k,d2)*n[d1];
+
+                      M(i,j) += dx*(2.*penalty*u*v - cv*u - cu*v);
+                    }
                 }
         }
     }
@@ -433,10 +437,10 @@ namespace LocalIntegrators
                   const unsigned int d1 = (d+1)%dim;
                   const unsigned int d2 = (d+2)%dim;
                   // curl u, curl v
-                  const double cv1 = nu1*fe1.shape_grad_component(i,k,d1)[d2] - fe1.shape_grad_component(i,k,d2)[d1];
-                  const double cv2 = nu2*fe2.shape_grad_component(i,k,d1)[d2] - fe2.shape_grad_component(i,k,d2)[d1];
-                  const double cu1 = nu1*fe1.shape_grad_component(j,k,d1)[d2] - fe1.shape_grad_component(j,k,d2)[d1];
-                  const double cu2 = nu2*fe2.shape_grad_component(j,k,d1)[d2] - fe2.shape_grad_component(j,k,d2)[d1];
+                  const double cv1 = nu1*fe1.shape_grad_component(i,k,d2)[d1] - fe1.shape_grad_component(i,k,d1)[d2];
+                  const double cv2 = nu2*fe2.shape_grad_component(i,k,d2)[d1] - fe2.shape_grad_component(i,k,d1)[d2];
+                  const double cu1 = nu1*fe1.shape_grad_component(j,k,d2)[d1] - fe1.shape_grad_component(j,k,d1)[d2];
+                  const double cu2 = nu2*fe2.shape_grad_component(j,k,d2)[d1] - fe2.shape_grad_component(j,k,d1)[d2];
 
                   // u x n, v x n
                   const double u1= fe1.shape_value_component(j,k,d1)*n(d2) - fe1.shape_value_component(j,k,d2)*n(d1);
diff --git a/tests/integrators/maxwell_curl.cc b/tests/integrators/maxwell_curl.cc
new file mode 100644 (file)
index 0000000..d435a0e
--- /dev/null
@@ -0,0 +1,109 @@
+// ---------------------------------------------------------------------
+//
+// Copyright (C) 2012 - 2015 by the deal.II authors
+//
+// Author: Jihuan Tian <jihuan_tian@hotmail.com>
+//
+// 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.
+//
+// ---------------------------------------------------------------------
+
+
+// Test curl related functions in integrators/maxwell.h
+
+#include <deal.II/base/logstream.h>
+#include <cstring>
+#include <fstream>
+#include <iostream>
+
+#include <deal.II/grid/grid_generator.h>
+#include <deal.II/grid/tria.h>
+#include <deal.II/grid/tria_iterator.h>
+#include <deal.II/grid/tria_accessor.h>
+
+#include <deal.II/dofs/dof_handler.h>
+#include <deal.II/dofs/dof_accessor.h>
+#include <deal.II/dofs/dof_tools.h>
+
+#include <deal.II/fe/fe_q.h>
+#include <deal.II/fe/fe_system.h>
+#include <deal.II/fe/fe_values.h>
+
+#include <deal.II/lac/full_matrix.h>
+
+#include <deal.II/integrators/maxwell.h>
+
+using namespace dealii;
+using namespace LocalIntegrators::Maxwell;
+
+template <int dim>
+void make_grid(Triangulation<dim> &tr)
+{
+  GridGenerator::subdivided_hyper_cube(tr, 1, -1., 1.);
+}
+
+template <int dim>
+void TestMaxwellCurl(Triangulation<dim> &tr)
+{
+  int element_order = 1;
+  int quadrature_order = 3;
+
+  DoFHandler<dim> dof_handler(tr);
+  FESystem<dim> fe(FE_Q<dim>(element_order), dim);
+  FEValues<dim> fe_values(fe, QGauss<dim>(quadrature_order), update_values | update_JxW_values | update_gradients | update_quadrature_points);
+  FEFaceValues<dim> fe_face_values(fe, QGauss<dim-1>(quadrature_order), update_values | update_JxW_values | update_gradients | update_normal_vectors | update_quadrature_points);
+
+  dof_handler.distribute_dofs(fe);
+
+  const unsigned int dofs_per_cell = fe.dofs_per_cell;
+
+  FullMatrix<double> curl_curl_check(dofs_per_cell, dofs_per_cell);
+  FullMatrix<double> curl_check(dofs_per_cell, dofs_per_cell);
+  FullMatrix<double> nitsche_curl_check(dofs_per_cell, dofs_per_cell);
+
+  curl_curl_check = 0;
+  curl_check = 0;
+  nitsche_curl_check = 0;
+
+  typename::DoFHandler<dim>::cell_iterator cell = dof_handler.begin(0);
+
+  fe_values.reinit(cell);
+
+  curl_curl_matrix<dim>(curl_curl_check, fe_values, 1.);
+  deallog << "curl_curl_matrix" << std::endl;
+  curl_curl_check.print(deallog, 10);
+
+  curl_matrix(curl_check, fe_values, fe_values, 1.);
+  deallog << "curl_matrix" << std::endl;
+  curl_check.print(deallog, 10);
+
+  for (unsigned int face = 0; face < GeometryInfo<dim>::faces_per_cell; ++face)
+    {
+      fe_face_values.reinit(cell, face);
+      nitsche_curl_matrix<dim>(nitsche_curl_check, fe_face_values, face, 1., 1.);
+    }
+
+  deallog << "nitsche_curl_matrix" << std::endl;
+  nitsche_curl_check.print(deallog, 10);
+
+  dof_handler.clear();
+}
+
+int main()
+{
+  const std::string logname = "output";
+  std::ofstream logfile(logname.c_str());
+  deallog.attach(logfile);
+  deallog.threshold_double(1.e-10);
+
+  Triangulation<3> tr;
+  make_grid(tr);
+  TestMaxwellCurl(tr);
+}
diff --git a/tests/integrators/maxwell_curl.output b/tests/integrators/maxwell_curl.output
new file mode 100644 (file)
index 0000000..30f6ffd
--- /dev/null
@@ -0,0 +1,76 @@
+
+DEAL::curl_curl_matrix
+DEAL::0.44      -0.17     -0.17     0.22      0.17      0.17      -0.11     -0.17     -0.083    -0.056    0.17      0.083     -0.11     -0.083    -0.17     -0.056    0.083     0.17      -0.22     -0.083    -0.083    -0.11     0.083     0.083     
+DEAL::-0.17     0.44      -0.17     -0.17     -0.11     -0.083    0.17      0.22      0.17      0.17      -0.056    0.083     -0.083    -0.11     -0.17     -0.083    -0.22     -0.083    0.083     -0.056    0.17      0.083     -0.11     0.083     
+DEAL::-0.17     -0.17     0.44      -0.17     -0.083    -0.11     -0.083    -0.17     -0.11     -0.083    -0.083    -0.22     0.17      0.17      0.22      0.17      0.083     -0.056    0.083     0.17      -0.056    0.083     0.083     -0.11     
+DEAL::0.22      -0.17     -0.17     0.44      0.17      0.17      -0.056    -0.17     -0.083    -0.11     0.17      0.083     -0.056    -0.083    -0.17     -0.11     0.083     0.17      -0.11     -0.083    -0.083    -0.22     0.083     0.083     
+DEAL::0.17      -0.11     -0.083    0.17      0.44      -0.17     -0.17     -0.056    0.083     -0.17     0.22      0.17      0.083     -0.22     -0.083    0.083     -0.11     -0.17     -0.083    -0.11     0.083     -0.083    -0.056    0.17      
+DEAL::0.17      -0.083    -0.11     0.17      -0.17     0.44      0.083     -0.083    -0.22     0.083     -0.17     -0.11     -0.17     0.083     -0.056    -0.17     0.17      0.22      -0.083    0.083     -0.11     -0.083    0.17      -0.056    
+DEAL::-0.11     0.17      -0.083    -0.056    -0.17     0.083     0.44      0.17      -0.17     0.22      -0.17     0.17      -0.22     0.083     -0.083    -0.11     -0.083    0.083     -0.11     0.083     -0.17     -0.056    -0.083    0.17      
+DEAL::-0.17     0.22      -0.17     -0.17     -0.056    -0.083    0.17      0.44      0.17      0.17      -0.11     0.083     -0.083    -0.056    -0.17     -0.083    -0.11     -0.083    0.083     -0.11     0.17      0.083     -0.22     0.083     
+DEAL::-0.083    0.17      -0.11     -0.083    0.083     -0.22     -0.17     0.17      0.44      -0.17     0.083     -0.11     0.083     -0.17     -0.056    0.083     -0.083    -0.11     0.17      -0.17     0.22      0.17      -0.083    -0.056    
+DEAL::-0.056    0.17      -0.083    -0.11     -0.17     0.083     0.22      0.17      -0.17     0.44      -0.17     0.17      -0.11     0.083     -0.083    -0.22     -0.083    0.083     -0.056    0.083     -0.17     -0.11     -0.083    0.17      
+DEAL::0.17      -0.056    -0.083    0.17      0.22      -0.17     -0.17     -0.11     0.083     -0.17     0.44      0.17      0.083     -0.11     -0.083    0.083     -0.056    -0.17     -0.083    -0.22     0.083     -0.083    -0.11     0.17      
+DEAL::0.083     0.083     -0.22     0.083     0.17      -0.11     0.17      0.083     -0.11     0.17      0.17      0.44      -0.083    -0.083    -0.11     -0.083    -0.17     -0.056    -0.17     -0.083    -0.056    -0.17     -0.17     0.22      
+DEAL::-0.11     -0.083    0.17      -0.056    0.083     -0.17     -0.22     -0.083    0.083     -0.11     0.083     -0.083    0.44      -0.17     0.17      0.22      0.17      -0.17     -0.11     -0.17     0.083     -0.056    0.17      -0.083    
+DEAL::-0.083    -0.11     0.17      -0.083    -0.22     0.083     0.083     -0.056    -0.17     0.083     -0.11     -0.083    -0.17     0.44      0.17      -0.17     -0.11     0.083     0.17      0.22      -0.17     0.17      -0.056    -0.083    
+DEAL::-0.17     -0.17     0.22      -0.17     -0.083    -0.056    -0.083    -0.17     -0.056    -0.083    -0.083    -0.11     0.17      0.17      0.44      0.17      0.083     -0.11     0.083     0.17      -0.11     0.083     0.083     -0.22     
+DEAL::-0.056    -0.083    0.17      -0.11     0.083     -0.17     -0.11     -0.083    0.083     -0.22     0.083     -0.083    0.22      -0.17     0.17      0.44      0.17      -0.17     -0.056    -0.17     0.083     -0.11     0.17      -0.083    
+DEAL::0.083     -0.22     0.083     0.083     -0.11     0.17      -0.083    -0.11     -0.083    -0.083    -0.056    -0.17     0.17      -0.11     0.083     0.17      0.44      0.17      -0.17     -0.056    -0.083    -0.17     0.22      -0.17     
+DEAL::0.17      -0.083    -0.056    0.17      -0.17     0.22      0.083     -0.083    -0.11     0.083     -0.17     -0.056    -0.17     0.083     -0.11     -0.17     0.17      0.44      -0.083    0.083     -0.22     -0.083    0.17      -0.11     
+DEAL::-0.22     0.083     0.083     -0.11     -0.083    -0.083    -0.11     0.083     0.17      -0.056    -0.083    -0.17     -0.11     0.17      0.083     -0.056    -0.17     -0.083    0.44      0.17      0.17      0.22      -0.17     -0.17     
+DEAL::-0.083    -0.056    0.17      -0.083    -0.11     0.083     0.083     -0.11     -0.17     0.083     -0.22     -0.083    -0.17     0.22      0.17      -0.17     -0.056    0.083     0.17      0.44      -0.17     0.17      -0.11     -0.083    
+DEAL::-0.083    0.17      -0.056    -0.083    0.083     -0.11     -0.17     0.17      0.22      -0.17     0.083     -0.056    0.083     -0.17     -0.11     0.083     -0.083    -0.22     0.17      -0.17     0.44      0.17      -0.083    -0.11     
+DEAL::-0.11     0.083     0.083     -0.22     -0.083    -0.083    -0.056    0.083     0.17      -0.11     -0.083    -0.17     -0.056    0.17      0.083     -0.11     -0.17     -0.083    0.22      0.17      0.17      0.44      -0.17     -0.17     
+DEAL::0.083     -0.11     0.083     0.083     -0.056    0.17      -0.083    -0.22     -0.083    -0.083    -0.11     -0.17     0.17      -0.056    0.083     0.17      0.22      0.17      -0.17     -0.11     -0.083    -0.17     0.44      -0.17     
+DEAL::0.083     0.083     -0.11     0.083     0.17      -0.056    0.17      0.083     -0.056    0.17      0.17      0.22      -0.083    -0.083    -0.22     -0.083    -0.17     -0.11     -0.17     -0.083    -0.11     -0.17     -0.17     0.44      
+DEAL::curl_matrix
+DEAL::0         0.22      -0.22     0         0.11      -0.11     0         0.11      0.22      0         0.056     0.11      0         -0.22     -0.11     0         -0.11     -0.056    0         -0.11     0.11      0         -0.056    0.056     
+DEAL::-0.22     0         0.22      -0.11     0         -0.22     -0.11     0         0.11      -0.056    0         -0.11     0.22      0         0.11      0.11      0         -0.11     0.11      0         0.056     0.056     0         -0.056    
+DEAL::0.22      -0.22     0         0.11      0.22      0         -0.22     -0.11     0         -0.11     0.11      0         0.11      -0.11     0         0.056     0.11      0         -0.11     -0.056    0         -0.056    0.056     0         
+DEAL::0         0.11      -0.11     0         0.22      -0.22     0         0.056     0.11      0         0.11      0.22      0         -0.11     -0.056    0         -0.22     -0.11     0         -0.056    0.056     0         -0.11     0.11      
+DEAL::-0.11     0         0.22      -0.22     0         -0.22     -0.056    0         0.11      -0.11     0         -0.11     0.11      0         0.11      0.22      0         -0.11     0.056     0         0.056     0.11      0         -0.056    
+DEAL::0.11      -0.22     0         0.22      0.22      0         -0.11     -0.11     0         -0.22     0.11      0         0.056     -0.11     0         0.11      0.11      0         -0.056    -0.056    0         -0.11     0.056     0         
+DEAL::0         0.11      -0.22     0         0.056     -0.11     0         0.22      0.22      0         0.11      0.11      0         -0.11     -0.11     0         -0.056    -0.056    0         -0.22     0.11      0         -0.11     0.056     
+DEAL::-0.11     0         0.11      -0.056    0         -0.11     -0.22     0         0.22      -0.11     0         -0.22     0.11      0         0.056     0.056     0         -0.056    0.22      0         0.11      0.11      0         -0.11     
+DEAL::0.22      -0.11     0         0.11      0.11      0         -0.22     -0.22     0         -0.11     0.22      0         0.11      -0.056    0         0.056     0.056     0         -0.11     -0.11     0         -0.056    0.11      0         
+DEAL::0         0.056     -0.11     0         0.11      -0.22     0         0.11      0.11      0         0.22      0.22      0         -0.056    -0.056    0         -0.11     -0.11     0         -0.11     0.056     0         -0.22     0.11      
+DEAL::-0.056    0         0.11      -0.11     0         -0.11     -0.11     0         0.22      -0.22     0         -0.22     0.056     0         0.056     0.11      0         -0.056    0.11      0         0.11      0.22      0         -0.11     
+DEAL::0.11      -0.11     0         0.22      0.11      0         -0.11     -0.22     0         -0.22     0.22      0         0.056     -0.056    0         0.11      0.056     0         -0.056    -0.11     0         -0.11     0.11      0         
+DEAL::0         0.22      -0.11     0         0.11      -0.056    0         0.11      0.11      0         0.056     0.056     0         -0.22     -0.22     0         -0.11     -0.11     0         -0.11     0.22      0         -0.056    0.11      
+DEAL::-0.22     0         0.11      -0.11     0         -0.11     -0.11     0         0.056     -0.056    0         -0.056    0.22      0         0.22      0.11      0         -0.22     0.11      0         0.11      0.056     0         -0.11     
+DEAL::0.11      -0.11     0         0.056     0.11      0         -0.11     -0.056    0         -0.056    0.056     0         0.22      -0.22     0         0.11      0.22      0         -0.22     -0.11     0         -0.11     0.11      0         
+DEAL::0         0.11      -0.056    0         0.22      -0.11     0         0.056     0.056     0         0.11      0.11      0         -0.11     -0.11     0         -0.22     -0.22     0         -0.056    0.11      0         -0.11     0.22      
+DEAL::-0.11     0         0.11      -0.22     0         -0.11     -0.056    0         0.056     -0.11     0         -0.056    0.11      0         0.22      0.22      0         -0.22     0.056     0         0.11      0.11      0         -0.11     
+DEAL::0.056     -0.11     0         0.11      0.11      0         -0.056    -0.056    0         -0.11     0.056     0         0.11      -0.22     0         0.22      0.22      0         -0.11     -0.11     0         -0.22     0.11      0         
+DEAL::0         0.11      -0.11     0         0.056     -0.056    0         0.22      0.11      0         0.11      0.056     0         -0.11     -0.22     0         -0.056    -0.11     0         -0.22     0.22      0         -0.11     0.11      
+DEAL::-0.11     0         0.056     -0.056    0         -0.056    -0.22     0         0.11      -0.11     0         -0.11     0.11      0         0.11      0.056     0         -0.11     0.22      0         0.22      0.11      0         -0.22     
+DEAL::0.11      -0.056    0         0.056     0.056     0         -0.11     -0.11     0         -0.056    0.11      0         0.22      -0.11     0         0.11      0.11      0         -0.22     -0.22     0         -0.11     0.22      0         
+DEAL::0         0.056     -0.056    0         0.11      -0.11     0         0.11      0.056     0         0.22      0.11      0         -0.056    -0.11     0         -0.11     -0.22     0         -0.11     0.11      0         -0.22     0.22      
+DEAL::-0.056    0         0.056     -0.11     0         -0.056    -0.11     0         0.11      -0.22     0         -0.11     0.056     0         0.11      0.11      0         -0.11     0.11      0         0.22      0.22      0         -0.22     
+DEAL::0.056     -0.056    0         0.11      0.056     0         -0.056    -0.11     0         -0.11     0.11      0         0.11      -0.11     0         0.22      0.11      0         -0.11     -0.22     0         -0.22     0.22      0         
+DEAL::nitsche_curl_matrix
+DEAL::2.7       -0.67     -0.67     1.3       0.33      0.33      0.67      -0.33     -0.33     0.33      0         0.17      0.67      -0.33     -0.33     0.33      0.17      0         0         -0.17     -0.17     0         0         0         
+DEAL::-0.67     2.7       -0.67     -0.33     0.67      -0.33     0.33      1.3       0.33      0         0.33      0.17      -0.33     0.67      -0.33     -0.17     0         -0.17     0.17      0.33      0         0         0         0         
+DEAL::-0.67     -0.67     2.7       -0.33     -0.33     0.67      -0.33     -0.33     0.67      -0.17     -0.17     0         0.33      0.33      1.3       0         0.17      0.33      0.17      0         0.33      0         0         0         
+DEAL::1.3       -0.33     -0.33     2.7       0.67      0.67      0.33      0         -0.17     0.67      0.33      0.33      0.33      -0.17     0         0.67      0.33      0.33      0         0         0         0         0.17      0.17      
+DEAL::0.33      0.67      -0.33     0.67      2.7       -0.67     0         0.33      0.17      -0.33     1.3       0.33      0.17      0         -0.17     0.33      0.67      -0.33     0         0         0         -0.17     0.33      0         
+DEAL::0.33      -0.33     0.67      0.67      -0.67     2.7       0.17      -0.17     0         0.33      -0.33     0.67      0         0.17      0.33      -0.33     0.33      1.3       0         0         0         -0.17     0         0.33      
+DEAL::0.67      0.33      -0.33     0.33      0         0.17      2.7       0.67      -0.67     1.3       -0.33     0.33      0         0.17      -0.17     0         0         0         0.67      0.33      -0.33     0.33      -0.17     0         
+DEAL::-0.33     1.3       -0.33     0         0.33      -0.17     0.67      2.7       0.67      0.33      0.67      0.33      -0.17     0.33      0         0         0         0         0.33      0.67      0.33      0.17      0         0.17      
+DEAL::-0.33     0.33      0.67      -0.17     0.17      0         -0.67     0.67      2.7       -0.33     0.33      0.67      0.17      0         0.33      0         0         0         0.33      -0.33     1.3       0         -0.17     0.33      
+DEAL::0.33      0         -0.17     0.67      -0.33     0.33      1.3       0.33      -0.33     2.7       -0.67     0.67      0         0         0         0         -0.17     0.17      0.33      0.17      0         0.67      -0.33     0.33      
+DEAL::0         0.33      -0.17     0.33      1.3       -0.33     -0.33     0.67      0.33      -0.67     2.7       0.67      0         0         0         0.17      0.33      0         -0.17     0         0.17      -0.33     0.67      0.33      
+DEAL::0.17      0.17      0         0.33      0.33      0.67      0.33      0.33      0.67      0.67      0.67      2.7       0         0         0         -0.17     0         0.33      0         -0.17     0.33      -0.33     -0.33     1.3       
+DEAL::0.67      -0.33     0.33      0.33      0.17      0         0         -0.17     0.17      0         0         0         2.7       -0.67     0.67      1.3       0.33      -0.33     0.67      -0.33     0.33      0.33      0         -0.17     
+DEAL::-0.33     0.67      0.33      -0.17     0         0.17      0.17      0.33      0         0         0         0         -0.67     2.7       0.67      -0.33     0.67      0.33      0.33      1.3       -0.33     0         0.33      -0.17     
+DEAL::-0.33     -0.33     1.3       0         -0.17     0.33      -0.17     0         0.33      0         0         0         0.67      0.67      2.7       0.33      0.33      0.67      0.33      0.33      0.67      0.17      0.17      0         
+DEAL::0.33      -0.17     0         0.67      0.33      -0.33     0         0         0         0         0.17      -0.17     1.3       -0.33     0.33      2.7       0.67      -0.67     0.33      0         0.17      0.67      0.33      -0.33     
+DEAL::0.17      0         0.17      0.33      0.67      0.33      0         0         0         -0.17     0.33      0         0.33      0.67      0.33      0.67      2.7       0.67      0         0.33      -0.17     -0.33     1.3       -0.33     
+DEAL::0         -0.17     0.33      0.33      -0.33     1.3       0         0         0         0.17      0         0.33      -0.33     0.33      0.67      -0.67     0.67      2.7       -0.17     0.17      0         -0.33     0.33      0.67      
+DEAL::0         0.17      0.17      0         0         0         0.67      0.33      0.33      0.33      -0.17     0         0.67      0.33      0.33      0.33      0         -0.17     2.7       0.67      0.67      1.3       -0.33     -0.33     
+DEAL::-0.17     0.33      0         0         0         0         0.33      0.67      -0.33     0.17      0         -0.17     -0.33     1.3       0.33      0         0.33      0.17      0.67      2.7       -0.67     0.33      0.67      -0.33     
+DEAL::-0.17     0         0.33      0         0         0         -0.33     0.33      1.3       0         0.17      0.33      0.33      -0.33     0.67      0.17      -0.17     0         0.67      -0.67     2.7       0.33      -0.33     0.67      
+DEAL::0         0         0         0         -0.17     -0.17     0.33      0.17      0         0.67      -0.33     -0.33     0.33      0         0.17      0.67      -0.33     -0.33     1.3       0.33      0.33      2.7       -0.67     -0.67     
+DEAL::0         0         0         0.17      0.33      0         -0.17     0         -0.17     -0.33     0.67      -0.33     0         0.33      0.17      0.33      1.3       0.33      -0.33     0.67      -0.33     -0.67     2.7       -0.67     
+DEAL::0         0         0         0.17      0         0.33      0         0.17      0.33      0.33      0.33      1.3       -0.17     -0.17     0         -0.33     -0.33     0.67      -0.33     -0.33     0.67      -0.67     -0.67     2.7       

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.