]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add one more test and fix it in fe_values.cc.
authorbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 1 Nov 2010 21:36:23 +0000 (21:36 +0000)
committerbangerth <bangerth@0785d39b-7218-0410-832d-ea1e28bc413d>
Mon, 1 Nov 2010 21:36:23 +0000 (21:36 +0000)
git-svn-id: https://svn.dealii.org/trunk@22580 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/source/fe/fe_values.cc
tests/codim_one/fe_values_extractor_01.cc [new file with mode: 0644]
tests/codim_one/fe_values_extractor_01/cmp/generic [new file with mode: 0644]
tests/codim_one/fe_values_extractor_01/square.msh [new file with mode: 0644]

index 1dd877aeeb4b294d76e63a838b8c360ce5fd6b8c..e89b8a4911bddadf4130feeebbaf8b8b144d0b46 100644 (file)
@@ -131,14 +131,14 @@ namespace FEValuesViews
                  first_vector_component (first_vector_component),
                  shape_function_data (fe_values.fe->dofs_per_cell)
   {
-    Assert (first_vector_component+dim-1 < fe_values.fe->n_components(),
-           ExcIndexRange(first_vector_component+dim-1, 0,
+    Assert (first_vector_component+spacedim-1 < fe_values.fe->n_components(),
+           ExcIndexRange(first_vector_component+spacedim-1, 0,
                          fe_values.fe->n_components()));
 
     const std::vector<unsigned int> shape_function_to_row_table
       = make_shape_function_to_row_table (*fe_values.fe);
 
-    for (unsigned int d=0; d<dim; ++d)
+    for (unsigned int d=0; d<spacedim; ++d)
       {
        const unsigned int component = first_vector_component + d;
 
@@ -180,7 +180,7 @@ namespace FEValuesViews
     for (unsigned int i=0; i<fe_values.fe->dofs_per_cell; ++i)
       {
        unsigned int n_nonzero_components = 0;
-       for (unsigned int d=0; d<dim; ++d)
+       for (unsigned int d=0; d<spacedim; ++d)
          if (shape_function_data[i].is_nonzero_shape_function_component[d]
              == true)
            ++n_nonzero_components;
@@ -191,7 +191,7 @@ namespace FEValuesViews
          shape_function_data[i].single_nonzero_component = -1;
        else
          {
-           for (unsigned int d=0; d<dim; ++d)
+           for (unsigned int d=0; d<spacedim; ++d)
              if (shape_function_data[i].is_nonzero_shape_function_component[d]
                  == true)
                {
@@ -1250,8 +1250,8 @@ namespace internal
                                       // dimensionality 'dim' of the
                                       // manifold, not 'spacedim' of
                                       // the output vector
-      const unsigned int n_vectors = (fe.n_components() >= dim ?
-                                     fe.n_components()-dim+1 :
+      const unsigned int n_vectors = (fe.n_components() >= spacedim ?
+                                     fe.n_components()-spacedim+1 :
                                      0);
       vectors.resize (n_vectors);
       for (unsigned int component=0; component<n_vectors; ++component)
diff --git a/tests/codim_one/fe_values_extractor_01.cc b/tests/codim_one/fe_values_extractor_01.cc
new file mode 100644 (file)
index 0000000..928c8f9
--- /dev/null
@@ -0,0 +1,79 @@
+
+//----------------------------  template.cc  ---------------------------
+//    $Id$
+//    Version: $Name$
+//
+//    Copyright (C) 2010 by the deal.II authors and Sebastian Pauletti
+//
+//    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.
+//
+//----------------------------  template.cc  ---------------------------
+
+
+// DataOut::build_patches appeared to have a problem with outputting
+// lines in 2d where nodes were numbered differently when writing data
+// vectors as opposed to writing node locations. in the end this
+// turned out to be a feature: the mesh was a circle of lines, so
+// there are equally many cells as their were nodes, and consequently
+// DataOut assumed that it had cell_data, rather than
+// dof_data. passing the correct argument fixed the problem, but it
+// won't hurt to have this test anyway.
+
+#include "../tests.h"
+#include <fstream>
+#include <base/logstream.h>
+
+#include <grid/tria.h>
+#include <dofs/dof_handler.h>
+#include <fe/fe_q.h>
+#include <fe/fe_system.h>
+#include <fe/fe_values.h>
+#include <grid/grid_in.h>
+#include <lac/vector.h>
+#include <base/quadrature_lib.h>
+#include <numerics/data_out.h>
+
+std::ofstream logfile("fe_values_extractor_01/output");
+
+
+int main ()
+{
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+   const unsigned int dim = 1;
+
+   Triangulation<dim,dim+1>    triangulation;
+   FESystem <dim,dim+1> fe (FE_Q<dim,dim+1> (1), dim+1);
+   DoFHandler<dim,dim+1>       dof_handler(triangulation);
+
+   Vector<double> soln;
+
+   GridIn<dim,dim+1> grid_in;
+   grid_in.attach_triangulation (triangulation);
+   std::ifstream fname("data_out_02/square.msh");
+   grid_in.read_msh (fname);
+
+   dof_handler.distribute_dofs (fe);
+   soln.reinit (dof_handler.n_dofs());
+   soln = 1;
+
+   std::vector<Tensor<1,dim+1> > local_velocity_values (1);
+   const FEValuesExtractors::Vector velocities (0);
+   QGauss<dim>  quadrature_formula(1);
+   DoFHandler<dim,dim+1>::active_cell_iterator
+     cell     = dof_handler.begin_active(),
+     endc     = dof_handler.end();
+   FEValues<dim,dim+1> fe_v (fe, quadrature_formula, update_values);
+
+   for (; cell!=endc; ++cell) {
+     fe_v.reinit (cell);
+
+     fe_v[velocities].get_function_values (soln, local_velocity_values);
+
+     deallog << local_velocity_values.at(0) << std::endl;
+   }
+}
diff --git a/tests/codim_one/fe_values_extractor_01/cmp/generic b/tests/codim_one/fe_values_extractor_01/cmp/generic
new file mode 100644 (file)
index 0000000..1a437e9
--- /dev/null
@@ -0,0 +1,5 @@
+
+DEAL::1.00000 1.00000
+DEAL::1.00000 1.00000
+DEAL::1.00000 1.00000
+DEAL::1.00000 1.00000
diff --git a/tests/codim_one/fe_values_extractor_01/square.msh b/tests/codim_one/fe_values_extractor_01/square.msh
new file mode 100644 (file)
index 0000000..997b625
--- /dev/null
@@ -0,0 +1,17 @@
+$MeshFormat
+2.1 0 8
+$EndMeshFormat
+$Nodes
+4
+1 0 0 0
+2 1 0 0
+3 1 1 0
+4 0 1 0
+$EndNodes
+$Elements
+4
+1 1 3 0 1 0 1 2
+2 1 3 0 1 0 2 3
+3 1 3 0 1 0 3 4
+4 1 3 0 1 0 4 1
+$EndElements

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.