]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add another test for inhomogeneous constraints.
authorkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 25 Feb 2009 08:52:46 +0000 (08:52 +0000)
committerkronbichler <kronbichler@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 25 Feb 2009 08:52:46 +0000 (08:52 +0000)
git-svn-id: https://svn.dealii.org/trunk@18424 0785d39b-7218-0410-832d-ea1e28bc413d

tests/deal.II/inhomogeneous_constraints.cc
tests/deal.II/inhomogeneous_constraints/cmp/generic
tests/deal.II/inhomogeneous_constraints_nonsymmetric.cc [new file with mode: 0644]
tests/deal.II/inhomogeneous_constraints_nonsymmetric/cmp/generic [new file with mode: 0644]

index ea796cca55f939995352ec5f5888446441a66c25..d0d58c7474dc4af3705bbc050dcd45e1f53ea179 100644 (file)
@@ -41,7 +41,6 @@
 #include <fe/fe_q.h>
 #include <numerics/vectors.h>
 #include <numerics/matrices.h>
-#include <numerics/data_out.h>
 #include <numerics/error_estimator.h>
 #include <lac/compressed_simple_sparsity_pattern.h>
 #include <hp/dof_handler.h>
@@ -214,20 +213,21 @@ void LaplaceProblem<dim>::test_equality ()
          test->value() = 0;
     }
 
-  deallog << "Matrix difference norm: " 
+  deallog << "  Matrix difference norm: " 
          << test_matrix.frobenius_norm() << std::endl;
   Assert (test_matrix.frobenius_norm() < 1e-13, ExcInternalError());
 
                                   // same here -- Dirichlet lines will have
-                                  // nonzero rhs, whereas we will have
-                                  // nonzero one.
+                                  // nonzero rhs, whereas we will have zero
+                                  // rhs when using inhomogeneous
+                                  // constraints.
   for (unsigned int i=0; i<reference_matrix.m(); ++i)
     if (test_all_constraints.is_constrained(i) == false)
       test_rhs(i) -= reference_rhs(i);
     else
       test_rhs(i) = 0;
 
-  deallog << "rhs difference norm: " 
+  deallog << "  RHS difference norm: " 
          << test_rhs.l2_norm() << std::endl;
 
   Assert (test_rhs.l2_norm() < 1e-14, ExcInternalError());
@@ -299,10 +299,6 @@ void LaplaceProblem<dim>::assemble_reference ()
     }
 
   hanging_nodes_only.condense (reference_matrix, reference_rhs);
-  deallog << "Reference matrix nonzeros: " << reference_matrix.n_nonzero_elements() 
-         << ", actually: " << reference_matrix.n_actually_nonzero_elements () 
-         << std::endl;
-
   std::map<unsigned int,double> boundary_values;
   VectorTools::interpolate_boundary_values (dof_handler,
                                            0,
@@ -312,6 +308,10 @@ void LaplaceProblem<dim>::assemble_reference ()
                                      reference_matrix,
                                      solution,
                                      reference_rhs);
+
+  deallog << "  Reference matrix nonzeros: " << reference_matrix.n_nonzero_elements() 
+         << ", actually: " << reference_matrix.n_actually_nonzero_elements () 
+         << std::endl;
 }
 
 
@@ -380,7 +380,7 @@ void LaplaceProblem<dim>::assemble_test_1 ()
     }
 
   test_all_constraints.condense (test_matrix, test_rhs);
-  deallog << "Test matrix 1 nonzeros: " << test_matrix.n_nonzero_elements() 
+  deallog << "  Test matrix 1 nonzeros: " << test_matrix.n_nonzero_elements() 
          << ", actually: " << test_matrix.n_actually_nonzero_elements () 
          << std::endl;
 
@@ -452,7 +452,7 @@ void LaplaceProblem<dim>::assemble_test_2 ()
                                                       test_matrix,
                                                       test_rhs);
     }
-  deallog << "Test matrix 2 nonzeros: " << test_matrix.n_nonzero_elements() 
+  deallog << "  Test matrix 2 nonzeros: " << test_matrix.n_nonzero_elements() 
          << ", actually: " << test_matrix.n_actually_nonzero_elements () 
          << std::endl;
   test_equality();
@@ -612,7 +612,7 @@ void LaplaceProblem<2>::create_coarse_grid ()
   triangulation.create_triangulation (vertices,
                                     cells,
                                     SubCellData());
-  triangulation.refine_global (3);
+  triangulation.refine_global (2);
 }
 
 
index 91ffd93473b9d603e0bdf4f3beac6d163ee1f37b..d85cc6921550fc1027819075fc4780940d3d766d 100644 (file)
@@ -1,46 +1,46 @@
 
 DEAL::
 DEAL::
-DEAL::   Number of active cells:       768
-DEAL::   Number of degrees of freedom: 3264
+DEAL::   Number of active cells:       192
+DEAL::   Number of degrees of freedom: 864
 DEAL::   Number of constraints       : 0
-DEAL::Reference matrix nonzeros: 49920, actually: 49766
-DEAL::Test matrix 1 nonzeros: 49920, actually: 42514
-DEAL::Matrix difference norm: 0
-DEAL::rhs difference norm: 0
-DEAL::Test matrix 2 nonzeros: 49920, actually: 42514
-DEAL::Matrix difference norm: 0
-DEAL::rhs difference norm: 0
-DEAL:cg::Starting value 28.
-DEAL:cg::Convergence step 48 value 3.5e-07
+DEAL::  Reference matrix nonzeros: 12672, actually: 8974
+DEAL::  Test matrix 1 nonzeros: 12672, actually: 8974
+DEAL::  Matrix difference norm: 0
+DEAL::  RHS difference norm: 0
+DEAL::  Test matrix 2 nonzeros: 12672, actually: 8974
+DEAL::  Matrix difference norm: 0
+DEAL::  RHS difference norm: 0
+DEAL:cg::Starting value 20.
+DEAL:cg::Convergence step 25 value 2.5e-07
 DEAL::Distribute error: 0
 DEAL::
 DEAL::
-DEAL::   Number of active cells:       1173
-DEAL::   Number of degrees of freedom: 5732
-DEAL::   Number of constraints       : 492
-DEAL::Reference matrix nonzeros: 98166, actually: 87212
-DEAL::Test matrix 1 nonzeros: 98166, actually: 78146
-DEAL::Matrix difference norm: 0
-DEAL::rhs difference norm: 0
-DEAL::Test matrix 2 nonzeros: 98166, actually: 78146
-DEAL::Matrix difference norm: 0
-DEAL::rhs difference norm: 0
-DEAL:cg::Starting value 31.
-DEAL:cg::Convergence step 84 value 4.4e-07
+DEAL::   Number of active cells:       309
+DEAL::   Number of degrees of freedom: 1575
+DEAL::   Number of constraints       : 187
+DEAL::  Reference matrix nonzeros: 25901, actually: 17333
+DEAL::  Test matrix 1 nonzeros: 25901, actually: 17333
+DEAL::  Matrix difference norm: 0
+DEAL::  RHS difference norm: 0
+DEAL::  Test matrix 2 nonzeros: 25901, actually: 17333
+DEAL::  Matrix difference norm: 0
+DEAL::  RHS difference norm: 0
+DEAL:cg::Starting value 22.
+DEAL:cg::Convergence step 42 value 3.2e-07
 DEAL::Distribute error: 0
 DEAL::
 DEAL::
-DEAL::   Number of active cells:       1644
-DEAL::   Number of degrees of freedom: 9764
-DEAL::   Number of constraints       : 1434
-DEAL::Reference matrix nonzeros: 200326, actually: 161794
-DEAL::Test matrix 1 nonzeros: 200326, actually: 150182
-DEAL::Matrix difference norm: 0
-DEAL::rhs difference norm: 0
-DEAL::Test matrix 2 nonzeros: 200326, actually: 150182
-DEAL::Matrix difference norm: 0
-DEAL::rhs difference norm: 0
-DEAL:cg::Starting value 35.
-DEAL:cg::Convergence step 107 value 5.2e-07
+DEAL::   Number of active cells:       450
+DEAL::   Number of degrees of freedom: 2692
+DEAL::   Number of constraints       : 475
+DEAL::  Reference matrix nonzeros: 51910, actually: 33908
+DEAL::  Test matrix 1 nonzeros: 51910, actually: 33908
+DEAL::  Matrix difference norm: 0
+DEAL::  RHS difference norm: 0
+DEAL::  Test matrix 2 nonzeros: 51910, actually: 33908
+DEAL::  Matrix difference norm: 0
+DEAL::  RHS difference norm: 0
+DEAL:cg::Starting value 25.
+DEAL:cg::Convergence step 61 value 3.8e-07
 DEAL::Distribute error: 0
diff --git a/tests/deal.II/inhomogeneous_constraints_nonsymmetric.cc b/tests/deal.II/inhomogeneous_constraints_nonsymmetric.cc
new file mode 100644 (file)
index 0000000..1591ca7
--- /dev/null
@@ -0,0 +1,488 @@
+//----------------  inhomogeneous_constraints_nonsymmetric.cc  -------------------
+//    $Id$
+//    Version: $Name$ 
+//
+//    Copyright (C) 2009 by the deal.II authors
+//
+//    This file is subject to QPL and may not be  distributed
+//    without copyright and license information. Please refer
+//    to the file deal.II/doc/license.html for the  text  and
+//    further information on this license.
+//
+//----------------  inhomogeneous_constraints_nonsymmetric.cc  -------------------
+
+
+// this function tests the correctness of the implementation of
+// inhomogeneous constraints on a nonsymmetric matrix that comes from a
+// discretization of the first derivative
+
+#include "../tests.h"
+
+#include <base/quadrature_lib.h>
+#include <base/function.h>
+#include <base/logstream.h>
+#include <base/utilities.h>
+#include <lac/vector.h>
+#include <lac/full_matrix.h>
+#include <lac/sparse_matrix.h>
+#include <grid/tria.h>
+#include <grid/grid_generator.h>
+#include <grid/grid_refinement.h>
+#include <dofs/dof_accessor.h>
+#include <dofs/dof_tools.h>
+#include <dofs/dof_constraints.h>
+#include <fe/fe_q.h>
+#include <fe/fe_values.h>
+#include <numerics/vectors.h>
+#include <numerics/matrices.h>
+#include <numerics/error_estimator.h>
+#include <lac/compressed_simple_sparsity_pattern.h>
+
+#include <fstream>
+#include <iostream>
+#include <complex>
+
+std::ofstream logfile("inhomogeneous_constraints_nonsymmetric/output");
+
+using namespace dealii;
+
+template <int dim>
+class AdvectionProblem 
+{
+  public:
+    AdvectionProblem ();
+    ~AdvectionProblem ();
+
+    void run ();
+    
+  private:
+    void setup_system ();
+    void test_equality ();
+    void assemble_reference ();
+    void assemble_test_1 ();
+    void assemble_test_2 ();
+
+    Triangulation<dim>   triangulation;
+
+    DoFHandler<dim>      dof_handler;
+    FE_Q<dim>            fe;
+
+    ConstraintMatrix     hanging_nodes_only;
+    ConstraintMatrix     test_all_constraints;
+
+    SparsityPattern      sparsity_pattern;
+    SparseMatrix<double> reference_matrix;
+    SparseMatrix<double> test_matrix;
+
+    Vector<double>       reference_rhs;
+    Vector<double>       test_rhs;
+};
+
+
+
+template <int dim>
+class RightHandSide : public Function<dim>
+{
+  public:
+    RightHandSide () : Function<dim> () {}
+    
+    virtual double value (const Point<dim>   &p,
+                         const unsigned int  component) const;
+};
+
+
+template <int dim>
+double
+RightHandSide<dim>::value (const Point<dim>   &p,
+                          const unsigned int  /*component*/) const
+{
+  double product = 1;
+  for (unsigned int d=0; d<dim; ++d)
+    product *= (p[d]+1);
+  return product;
+}
+
+
+template <int dim>
+AdvectionProblem<dim>::AdvectionProblem ()
+               :
+                dof_handler (triangulation),
+               fe (2)
+{}
+
+
+template <int dim>
+AdvectionProblem<dim>::~AdvectionProblem () 
+{
+  dof_handler.clear ();
+}
+
+
+template <int dim>
+void AdvectionProblem<dim>::setup_system ()
+{
+  dof_handler.distribute_dofs (fe);
+
+  reference_rhs.reinit (dof_handler.n_dofs());
+  test_rhs.reinit (dof_handler.n_dofs());
+
+  hanging_nodes_only.clear ();
+  test_all_constraints.clear ();
+
+                                  // add boundary conditions as
+                                  // inhomogeneous constraints here. In
+                                  // contrast to step-27, we choose a
+                                  // constant function with value 1 here.
+  {
+    std::map<unsigned int,double> boundary_values;
+    VectorTools::interpolate_boundary_values (dof_handler,
+                                             0,
+                                             RightHandSide<dim>(),
+                                             boundary_values);
+    std::map<unsigned int,double>::const_iterator boundary_value = 
+      boundary_values.begin();
+    for ( ; boundary_value !=boundary_values.end(); ++boundary_value)
+      {
+       test_all_constraints.add_line(boundary_value->first);
+       test_all_constraints.set_inhomogeneity (boundary_value->first, 
+                                               boundary_value->second);
+      }
+  }
+  DoFTools::make_hanging_node_constraints (dof_handler,
+                                          hanging_nodes_only);
+  DoFTools::make_hanging_node_constraints (dof_handler,
+                                          test_all_constraints);
+  hanging_nodes_only.close ();
+  test_all_constraints.close ();
+
+  CompressedSimpleSparsityPattern csp (dof_handler.n_dofs(),
+                                      dof_handler.n_dofs());
+  DoFTools::make_sparsity_pattern (dof_handler, csp,
+                                  hanging_nodes_only, true);
+  sparsity_pattern.copy_from (csp);
+
+  reference_matrix.reinit (sparsity_pattern);
+  test_matrix.reinit (sparsity_pattern);
+}
+
+
+
+                                  // test whether we are equal with the
+                                  // standard matrix and right hand side
+template <int dim>
+void AdvectionProblem<dim>::test_equality ()
+{
+                                  // need to manually go through the
+                                  // matrix, since we can have different
+                                  // entries in constrained lines.
+  for (unsigned int i=0; i<reference_matrix.m(); ++i)
+    {
+      SparseMatrix<double>::const_iterator reference = reference_matrix.begin(i);
+      SparseMatrix<double>::iterator test = test_matrix.begin(i);
+      if (test_all_constraints.is_constrained(i) == false)
+       {
+         for ( ; test != test_matrix.end(i); ++test, ++reference)
+             test->value() -= reference->value();
+       }
+      else
+       for ( ; test != test_matrix.end(i); ++test)
+         test->value() = 0;
+    }
+
+  deallog << "  Matrix difference norm: " 
+         << test_matrix.frobenius_norm() << std::endl;
+  Assert (test_matrix.frobenius_norm() < 1e-13, ExcInternalError());
+
+                                  // same here -- Dirichlet lines will have
+                                  // nonzero rhs, whereas we will have zero
+                                  // rhs when using inhomogeneous
+                                  // constraints.
+  for (unsigned int i=0; i<reference_matrix.m(); ++i)
+    if (test_all_constraints.is_constrained(i) == false)
+      test_rhs(i) -= reference_rhs(i);
+    else
+      test_rhs(i) = 0;
+
+  deallog << "  RHS difference norm: " 
+         << test_rhs.l2_norm() << std::endl;
+
+  Assert (test_rhs.l2_norm() < 1e-14, ExcInternalError());
+}
+
+
+
+
+template <int dim>
+void AdvectionProblem<dim>::assemble_reference () 
+{
+  reference_matrix = 0;
+  reference_rhs = 0;
+
+  QGauss<dim> quadrature_formula (3);
+  FEValues<dim> fe_values (fe, quadrature_formula, 
+                          update_values    |  update_gradients | 
+                          update_quadrature_points  |  update_JxW_values);
+
+  const RightHandSide<dim> rhs_function;
+  const unsigned int dofs_per_cell = fe.dofs_per_cell;
+  const unsigned int n_q_points = quadrature_formula.size();
+
+  FullMatrix<double>   cell_matrix (dofs_per_cell, dofs_per_cell);
+  Vector<double>       cell_rhs (dofs_per_cell);
+
+  std::vector<unsigned int> local_dof_indices (dofs_per_cell);
+  std::vector<double>  rhs_values (n_q_points);
+  
+  typename DoFHandler<dim>::active_cell_iterator
+    cell = dof_handler.begin_active(),
+    endc = dof_handler.end();
+  for (; cell!=endc; ++cell)
+    {
+      cell_matrix = 0;
+      cell_rhs = 0;
+      fe_values.reinit (cell);
+
+      rhs_function.value_list (fe_values.get_quadrature_points(),
+                              rhs_values);
+
+      Tensor<1,dim> advection_direction;
+      advection_direction[0] = 1;
+      advection_direction[1] = 1;
+      advection_direction[dim-1] = -1;
+      
+      for (unsigned int q_point=0; q_point<n_q_points; ++q_point)
+       for (unsigned int i=0; i<dofs_per_cell; ++i)
+         {
+           for (unsigned int j=0; j<dofs_per_cell; ++j)
+             cell_matrix(i,j) += (fe_values.shape_value(i,q_point) *
+                                  advection_direction *
+                                  fe_values.shape_grad(j,q_point) *
+                                  fe_values.JxW(q_point));
+
+           cell_rhs(i) += (fe_values.shape_value(i,q_point) *
+                           rhs_values[q_point] *
+                           fe_values.JxW(q_point));
+         }
+
+      local_dof_indices.resize (dofs_per_cell);
+      cell->get_dof_indices (local_dof_indices);
+
+      reference_matrix.add(local_dof_indices, cell_matrix);
+      for (unsigned int i=0; i<dofs_per_cell; ++i)
+       reference_rhs(local_dof_indices[i]) += cell_rhs(i);
+    }
+
+  hanging_nodes_only.condense (reference_matrix, reference_rhs);
+
+                                  // use some other rhs vector as dummy for
+                                  // application of Dirichlet conditions
+  std::map<unsigned int,double> boundary_values;
+  VectorTools::interpolate_boundary_values (dof_handler,
+                                           0,
+                                           RightHandSide<dim>(),
+                                           boundary_values);
+  MatrixTools::apply_boundary_values (boundary_values,
+                                     reference_matrix,
+                                     test_rhs,
+                                     reference_rhs);
+
+  deallog << "  Reference matrix nonzeros: " << reference_matrix.n_nonzero_elements() 
+         << ", actually: " << reference_matrix.n_actually_nonzero_elements () 
+         << std::endl;
+}
+
+
+
+template <int dim>
+void AdvectionProblem<dim>::assemble_test_1 () 
+{
+  test_matrix = 0;
+  test_rhs = 0;
+
+
+  QGauss<dim> quadrature_formula (3);
+  FEValues<dim> fe_values (fe, quadrature_formula, 
+                          update_values    |  update_gradients | 
+                          update_quadrature_points  |  update_JxW_values);
+
+  const RightHandSide<dim> rhs_function;
+  const unsigned int dofs_per_cell = fe.dofs_per_cell;
+  const unsigned int n_q_points = quadrature_formula.size();
+
+  FullMatrix<double>   cell_matrix (dofs_per_cell, dofs_per_cell);
+  Vector<double>       cell_rhs (dofs_per_cell);
+
+  std::vector<unsigned int> local_dof_indices (dofs_per_cell);
+  std::vector<double>  rhs_values (n_q_points);
+  
+  typename DoFHandler<dim>::active_cell_iterator
+    cell = dof_handler.begin_active(),
+    endc = dof_handler.end();
+  for (; cell!=endc; ++cell)
+    {
+      cell_matrix = 0;
+      cell_rhs = 0;
+      fe_values.reinit (cell);
+
+      rhs_function.value_list (fe_values.get_quadrature_points(),
+                              rhs_values);
+
+      Tensor<1,dim> advection_direction;
+      advection_direction[0] = 1;
+      advection_direction[1] = 1;
+      advection_direction[dim-1] = -1;
+      
+      for (unsigned int q_point=0; q_point<n_q_points; ++q_point)
+       for (unsigned int i=0; i<dofs_per_cell; ++i)
+         {
+           for (unsigned int j=0; j<dofs_per_cell; ++j)
+             cell_matrix(i,j) += (fe_values.shape_value(i,q_point) *
+                                  advection_direction *
+                                  fe_values.shape_grad(j,q_point) *
+                                  fe_values.JxW(q_point));
+
+           cell_rhs(i) += (fe_values.shape_value(i,q_point) *
+                           rhs_values[q_point] *
+                           fe_values.JxW(q_point));
+         }
+
+      local_dof_indices.resize (dofs_per_cell);
+      cell->get_dof_indices (local_dof_indices);
+
+      test_matrix.add(local_dof_indices, cell_matrix);
+      for (unsigned int i=0; i<dofs_per_cell; ++i)
+       test_rhs(local_dof_indices[i]) += cell_rhs(i);
+       
+    }
+
+  test_all_constraints.condense (test_matrix, test_rhs);
+  deallog << "  Test matrix 1 nonzeros: " << test_matrix.n_nonzero_elements() 
+         << ", actually: " << test_matrix.n_actually_nonzero_elements () 
+         << std::endl;
+
+  test_equality();
+}
+
+
+
+template <int dim>
+void AdvectionProblem<dim>::assemble_test_2 () 
+{
+  test_matrix = 0;
+  test_rhs = 0;
+
+  QGauss<dim> quadrature_formula (3);
+  FEValues<dim> fe_values (fe, quadrature_formula, 
+                          update_values    |  update_gradients | 
+                          update_quadrature_points  |  update_JxW_values);
+
+  const RightHandSide<dim> rhs_function;
+  const unsigned int dofs_per_cell = fe.dofs_per_cell;
+  const unsigned int n_q_points = quadrature_formula.size();
+
+  FullMatrix<double>   cell_matrix (dofs_per_cell, dofs_per_cell);
+  Vector<double>       cell_rhs (dofs_per_cell);
+
+  std::vector<unsigned int> local_dof_indices (dofs_per_cell);
+  std::vector<double>  rhs_values (n_q_points);
+  
+  typename DoFHandler<dim>::active_cell_iterator
+    cell = dof_handler.begin_active(),
+    endc = dof_handler.end();
+  for (; cell!=endc; ++cell)
+    {
+      cell_matrix = 0;
+      cell_rhs = 0;
+      fe_values.reinit (cell);
+
+      rhs_function.value_list (fe_values.get_quadrature_points(),
+                              rhs_values);
+
+      Tensor<1,dim> advection_direction;
+      advection_direction[0] = 1;
+      advection_direction[1] = 1;
+      advection_direction[dim-1] = -1;
+      
+      for (unsigned int q_point=0; q_point<n_q_points; ++q_point)
+       for (unsigned int i=0; i<dofs_per_cell; ++i)
+         {
+           for (unsigned int j=0; j<dofs_per_cell; ++j)
+             cell_matrix(i,j) += (fe_values.shape_value(i,q_point) *
+                                  advection_direction *
+                                  fe_values.shape_grad(j,q_point) *
+                                  fe_values.JxW(q_point));
+
+           cell_rhs(i) += (fe_values.shape_value(i,q_point) *
+                           rhs_values[q_point] *
+                           fe_values.JxW(q_point));
+         }
+
+      local_dof_indices.resize (dofs_per_cell);
+      cell->get_dof_indices (local_dof_indices);
+
+      test_all_constraints.distribute_local_to_global (cell_matrix,
+                                                      cell_rhs,
+                                                      local_dof_indices,
+                                                      test_matrix,
+                                                      test_rhs);
+    }
+  deallog << "  Test matrix 2 nonzeros: " << test_matrix.n_nonzero_elements() 
+         << ", actually: " << test_matrix.n_actually_nonzero_elements () 
+         << std::endl;
+  test_equality();
+}
+
+
+template <int dim>
+void AdvectionProblem<dim>::run () 
+{
+  GridGenerator::hyper_cube (triangulation);
+  triangulation.refine_global (2);
+
+                                  // manually refine the first two cells to
+                                  // create some hanging nodes
+  {
+    typename DoFHandler<dim>::active_cell_iterator cell = dof_handler.begin_active();
+    cell->set_refine_flag();
+    cell++;
+    cell->set_refine_flag();
+  }
+  triangulation.execute_coarsening_and_refinement();  
+
+  setup_system ();
+
+  deallog << std::endl << std::endl
+         << "  Number of active cells:       "
+         << triangulation.n_active_cells()
+         << std::endl
+         << "  Number of degrees of freedom: "
+         << dof_handler.n_dofs()
+         << std::endl
+         << "  Number of constraints       : "
+         << hanging_nodes_only.n_constraints()
+         << std::endl;
+
+  assemble_reference ();
+  assemble_test_1 ();
+  assemble_test_2 ();
+}
+
+
+
+int main () 
+{
+  deallog << std::setprecision (2);
+  logfile << std::setprecision (2);
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+  deallog.threshold_double(1.e-12);
+
+  {
+    AdvectionProblem<2> advection_problem;
+    advection_problem.run ();
+  }
+  {
+    AdvectionProblem<3> advection_problem;
+    advection_problem.run ();
+  }
+}
diff --git a/tests/deal.II/inhomogeneous_constraints_nonsymmetric/cmp/generic b/tests/deal.II/inhomogeneous_constraints_nonsymmetric/cmp/generic
new file mode 100644 (file)
index 0000000..d6a1000
--- /dev/null
@@ -0,0 +1,25 @@
+
+DEAL::
+DEAL::
+DEAL::  Number of active cells:       22
+DEAL::  Number of degrees of freedom: 114
+DEAL::  Number of constraints       : 9
+DEAL::  Reference matrix nonzeros: 1618, actually: 701
+DEAL::  Test matrix 1 nonzeros: 1618, actually: 701
+DEAL::  Matrix difference norm: 0
+DEAL::  RHS difference norm: 0
+DEAL::  Test matrix 2 nonzeros: 1618, actually: 698
+DEAL::  Matrix difference norm: 0
+DEAL::  RHS difference norm: 0
+DEAL::
+DEAL::
+DEAL::  Number of active cells:       78
+DEAL::  Number of degrees of freedom: 928
+DEAL::  Number of constraints       : 87
+DEAL::  Reference matrix nonzeros: 46980, actually: 14420
+DEAL::  Test matrix 1 nonzeros: 46980, actually: 14420
+DEAL::  Matrix difference norm: 0
+DEAL::  RHS difference norm: 0
+DEAL::  Test matrix 2 nonzeros: 46980, actually: 14416
+DEAL::  Matrix difference norm: 0
+DEAL::  RHS difference norm: 0

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.