--- /dev/null
+//---------------------------- non_primitive_1.cc ---------------------------
+// $Id$
+// Version: $Name$
+//
+// Copyright (C) 2002 by the deal.II authors and Anna Schneebeli
+//
+// 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.
+//
+//---------------------------- non_primitive_1.cc ---------------------------
+
+
+// check some things about Nedelec elements, basically also that the
+// DoFRenumbering::component_wise function also works for
+// non_primitive elements, for which it did not work previously since
+// there is no component to associate a non-primitive shape function
+// with
+//
+// this program is a modified version of one by Anna Schneebeli,
+// University of Basel
+
+#include <base/logstream.h>
+#include <grid/tria.h>
+#include <dofs/dof_handler.h>
+#include <grid/grid_generator.h>
+#include <grid/tria_accessor.h>
+#include <grid/tria_iterator.h>
+#include <dofs/dof_accessor.h>
+#include <dofs/dof_tools.h>
+#include <fe/fe_system.h>
+#include <fe/fe_q.h>
+#include <fe/fe_nedelec.h>
+#include <fe/fe_base.h>
+#include <numerics/dof_renumbering.h>
+#include <iostream>
+#include <fstream>
+
+
+class SystemTest
+{
+ public:
+ SystemTest ();
+ void run ();
+
+ private:
+ void make_grid_and_dofs ();
+ void shape_to_components ();
+ void check_numbering ();
+
+
+ Triangulation<2> triangulation;
+ FESystem<2> fe;
+ DoFHandler<2> dof_handler;
+
+
+};
+
+SystemTest::SystemTest () :
+ fe (FE_Nedelec<2>(1), 2,
+ FE_Q<2>(1), 1),
+ dof_handler (triangulation)
+{};
+
+
+void SystemTest::make_grid_and_dofs ()
+{
+
+ GridGenerator::hyper_cube (triangulation, -1, 1);
+ triangulation.refine_global (0);
+ deallog << "Number of active cells: " << triangulation.n_active_cells()
+ << std::endl;
+ deallog << "Total number of cells: " << triangulation.n_cells()
+ << std::endl;
+
+ dof_handler.distribute_dofs (fe);
+ deallog << "Number of degrees of freedom: " << dof_handler.n_dofs()
+ << std::endl;
+
+};
+
+void SystemTest::shape_to_components ()
+{
+ // testing, if the shape function
+ // with index i is of type Nedelec:
+ // (i.e. the first component of the FESystem)
+ // 1 for yes, 0 for no.
+
+ for(unsigned int i = 0; i<fe.dofs_per_cell; i++)
+ deallog <<" shapefunction "<< i << " is Nedelec: "
+ << (fe.is_primitive(i) ? "false" : "true") << std::endl;
+};
+
+
+
+void SystemTest::check_numbering ()
+{
+ DoFHandler<2>::active_cell_iterator cell = dof_handler.begin_active(),
+ endc = dof_handler.end();
+ std::vector<unsigned int> local_dof_indices(fe.dofs_per_cell);
+
+ for (; cell!=endc; ++cell)
+ {
+ cell->get_dof_indices (local_dof_indices);
+ for (unsigned int i=0; i<fe.dofs_per_cell; i++)
+ deallog <<" DoF "<< local_dof_indices[i] << " belongs to base element "
+ << fe.system_to_base_index(i).first.first
+ << ", instance " << fe.system_to_base_index(i).first.second
+ << std::endl;
+ deallog<< std::endl;
+ };
+
+
+ //Now: Componentwise reodering of the dofs
+
+ deallog << " Now we renumber the DoFs component-wise:" << std::endl;
+ deallog << " ****************************************" << std::endl;
+ DoFRenumbering::component_wise (dof_handler);
+
+ cell = dof_handler.begin_active();
+ endc = dof_handler.end();
+
+ for (; cell!=endc; ++cell)
+ {
+ cell->get_dof_indices (local_dof_indices);
+ for(unsigned int i=0; i<fe.dofs_per_cell; i++)
+ deallog <<" DoF "<< local_dof_indices[i] << " belongs to base "
+ << fe.system_to_base_index(i).first.first
+ << ", instance " << fe.system_to_base_index(i).first.second
+ << std::endl;
+ deallog << std::endl;
+ };
+
+};
+
+
+void SystemTest::run ()
+{
+ make_grid_and_dofs ();
+ shape_to_components ();
+ check_numbering();
+};
+
+
+
+int main ()
+{
+ std::ofstream logfile("anna_1.output");
+ deallog.attach(logfile);
+ deallog.depth_console(0);
+
+ SystemTest system_test;
+ system_test.run ();
+ return 0;
+};
--- /dev/null
+
+DEAL::Number of active cells: 1
+DEAL::Total number of cells: 1
+DEAL::Number of degrees of freedom: 12
+DEAL:: shapefunction 0 is Nedelec: false
+DEAL:: shapefunction 1 is Nedelec: false
+DEAL:: shapefunction 2 is Nedelec: false
+DEAL:: shapefunction 3 is Nedelec: false
+DEAL:: shapefunction 4 is Nedelec: true
+DEAL:: shapefunction 5 is Nedelec: true
+DEAL:: shapefunction 6 is Nedelec: true
+DEAL:: shapefunction 7 is Nedelec: true
+DEAL:: shapefunction 8 is Nedelec: true
+DEAL:: shapefunction 9 is Nedelec: true
+DEAL:: shapefunction 10 is Nedelec: true
+DEAL:: shapefunction 11 is Nedelec: true
+DEAL:: DoF 0 belongs to base element 1, instance 0
+DEAL:: DoF 1 belongs to base element 1, instance 0
+DEAL:: DoF 2 belongs to base element 1, instance 0
+DEAL:: DoF 3 belongs to base element 1, instance 0
+DEAL:: DoF 4 belongs to base element 0, instance 0
+DEAL:: DoF 5 belongs to base element 0, instance 1
+DEAL:: DoF 6 belongs to base element 0, instance 0
+DEAL:: DoF 7 belongs to base element 0, instance 1
+DEAL:: DoF 8 belongs to base element 0, instance 0
+DEAL:: DoF 9 belongs to base element 0, instance 1
+DEAL:: DoF 10 belongs to base element 0, instance 0
+DEAL:: DoF 11 belongs to base element 0, instance 1
+DEAL::
+DEAL:: Now we renumber the DoFs component-wise:
+DEAL:: ****************************************
+DEAL:: DoF 8 belongs to base 1, instance 0
+DEAL:: DoF 9 belongs to base 1, instance 0
+DEAL:: DoF 10 belongs to base 1, instance 0
+DEAL:: DoF 11 belongs to base 1, instance 0
+DEAL:: DoF 0 belongs to base 0, instance 0
+DEAL:: DoF 4 belongs to base 0, instance 1
+DEAL:: DoF 1 belongs to base 0, instance 0
+DEAL:: DoF 5 belongs to base 0, instance 1
+DEAL:: DoF 2 belongs to base 0, instance 0
+DEAL:: DoF 6 belongs to base 0, instance 1
+DEAL:: DoF 3 belongs to base 0, instance 0
+DEAL:: DoF 7 belongs to base 0, instance 1
+DEAL::