]> https://gitweb.dealii.org/ - dealii.git/commitdiff
simplified Grid generation in hp_constraints_neither_dominate_01
authorDenis Davydov <davydden@gmail.com>
Wed, 2 Sep 2015 20:47:09 +0000 (22:47 +0200)
committerDenis Davydov <davydden@gmail.com>
Wed, 2 Sep 2015 21:14:57 +0000 (23:14 +0200)
Also fixed indentation.

tests/hp/hp_constraints_neither_dominate_01.cc

index 36824cecce9b1d772acd48bee9377e22d328c7aa..943696b68a28294408621ee27f945e973aa5a287 100644 (file)
@@ -70,7 +70,8 @@
 #include <fstream>
 #include <iostream>
 
-// #define DEBUG_OUTPUT_VTK
+//#define DEBUG_OUTPUT_VTK
+unsigned int counter=0;
 
 const double eps = 1e-10;
 
@@ -80,16 +81,19 @@ template<int dim>
 void test2cells(const unsigned int p1=2,
                 const unsigned int p2=1)
 {
+  Assert(dim>1,ExcInternalError());
   Triangulation<dim>   triangulation;
   {
-    Triangulation<dim>   triangulationL;
-    Triangulation<dim>   triangulationR;
-    GridGenerator::hyper_cube (triangulationL, -1,0); //create a square [-1,0]^d domain
-    GridGenerator::hyper_cube (triangulationR, -1,0); //create a square [-1,0]^d domain
-    Point<dim> shift_vector;
-    shift_vector[0] = 1.0;
-    GridTools::shift(shift_vector,triangulationR);
-    GridGenerator::merge_triangulations (triangulationL, triangulationR, triangulation);
+    Point<dim> p1,p2;
+    for (unsigned int d=0; d<dim; d++)
+      p1[d] = -1;
+    p2[0] = 1.0;
+    std::vector<unsigned int> repetitoins(dim,1);
+    repetitoins[0]=2;
+    GridGenerator::subdivided_hyper_rectangle(triangulation,
+                                              repetitoins,
+                                              p1,
+                                              p2);
   }
 
   hp::DoFHandler<dim> dof_handler(triangulation);
@@ -123,6 +127,7 @@ void test2cells(const unsigned int p1=2,
 
 #ifdef DEBUG_OUTPUT_VTK
   // output to check if all is good:
+  counter++;
   std::vector<Vector<double>> shape_functions;
   std::vector<std::string> names;
   for (unsigned int s=0; s < dof_handler.n_dofs(); s++)
@@ -164,11 +169,9 @@ void test2cells(const unsigned int p1=2,
   for (unsigned int i = 0; i < shape_functions.size(); i++)
     data_out.add_data_vector (shape_functions[i], names[i]);
 
-  data_out.build_patches(15);
-
-  std::string filename = "2cell-shape_functions_p1="+std::to_string(p1)
-    +"_p2="+std::to_string(p2)
-    +"_"+dealii::Utilities::int_to_string(dim)+"D.vtu";
+  data_out.build_patches(0);
+  std::string filename = "shape_functions_"+dealii::Utilities::int_to_string(counter,1)
+                         +"_"+dealii::Utilities::int_to_string(dim)+"D.vtu";
   std::ofstream output (filename.c_str ());
   data_out.write_vtu (output);
 #endif
@@ -189,30 +192,30 @@ void test2cells(const unsigned int p1=2,
   std::vector<unsigned int> local_face_dof_indices;
   std::vector<Vector<double> > values;
   for (typename hp::DoFHandler<dim>::active_cell_iterator
-         cell = dof_handler.begin_active();
+       cell = dof_handler.begin_active();
        cell != dof_handler.end(); ++cell)
-      {
-        const unsigned int fe_index = cell->active_fe_index();
-        local_face_dof_indices.resize(fe_collection[fe_index].dofs_per_face);
-        for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
-          if (!cell->at_boundary(f)) // that's enough for our simple mesh
-            {
-              deallog << "cell="<<cell<<" face="<<f<<std::endl;
-              fe_face_values_hp.reinit(cell,f);
-              const FEFaceValues<dim> &fe_face_values = fe_face_values_hp.get_present_fe_values();
-
-              const unsigned int n_q_points = fe_face_values.n_quadrature_points;
-              values.resize (n_q_points,
-                             Vector<double>(2));
-
-              fe_face_values.get_function_values (solution,values);
-
-              const std::vector<dealii::Point<dim> > &q_points = fe_face_values.get_quadrature_points();
-
-              for (unsigned int q = 0; q < n_q_points; q++)
-                deallog << "u["<<q_points[q]<<"]={"<<values[q][0]<<","<<values[q][1]<<"}"<<std::endl;
-            }
-      }
+    {
+      const unsigned int fe_index = cell->active_fe_index();
+      local_face_dof_indices.resize(fe_collection[fe_index].dofs_per_face);
+      for (unsigned int f=0; f<GeometryInfo<dim>::faces_per_cell; ++f)
+        if (!cell->at_boundary(f)) // that's enough for our simple mesh
+          {
+            deallog << "cell="<<cell<<" face="<<f<<std::endl;
+            fe_face_values_hp.reinit(cell,f);
+            const FEFaceValues<dim> &fe_face_values = fe_face_values_hp.get_present_fe_values();
+
+            const unsigned int n_q_points = fe_face_values.n_quadrature_points;
+            values.resize (n_q_points,
+                           Vector<double>(2));
+
+            fe_face_values.get_function_values (solution,values);
+
+            const std::vector<dealii::Point<dim> > &q_points = fe_face_values.get_quadrature_points();
+
+            for (unsigned int q = 0; q < n_q_points; q++)
+              deallog << "u["<<q_points[q]<<"]={"<<values[q][0]<<","<<values[q][1]<<"}"<<std::endl;
+          }
+    }
 
   dof_handler.clear();
 }
@@ -227,34 +230,34 @@ int main (int argc,char **argv)
 
 
   try
-  {
-    test2cells<2>(1,2); // Q(1) x FE_Nothing vs Q(2) x Q(1) => common Q(1) x FE_Nothing
-    test2cells<2>(2,1); // Q(2) x FE_Nothing vs Q(1) x Q(1) => common Q(1) x FE_Nothing
-    test2cells<3>(1,2); // Q(1) x FE_Nothing vs Q(2) x Q(1) => common Q(1) x FE_Nothing
-    test2cells<3>(2,1); // Q(2) x FE_Nothing vs Q(1) x Q(1) => common Q(1) x FE_Nothing
-  }
+    {
+      test2cells<2>(1,2); // Q(1) x FE_Nothing vs Q(2) x Q(1) => common Q(1) x FE_Nothing
+      test2cells<2>(2,1); // Q(2) x FE_Nothing vs Q(1) x Q(1) => common Q(1) x FE_Nothing
+      test2cells<3>(1,2); // Q(1) x FE_Nothing vs Q(2) x Q(1) => common Q(1) x FE_Nothing
+      test2cells<3>(2,1); // Q(2) x FE_Nothing vs Q(1) x Q(1) => common Q(1) x FE_Nothing
+    }
   catch (std::exception &exc)
-  {
-    std::cerr << std::endl << std::endl
-              << "----------------------------------------------------"
-              << std::endl;
-    std::cerr << "Exception on processing: " << std::endl
-              << exc.what() << std::endl
-              << "Aborting!" << std::endl
-              << "----------------------------------------------------"
-              << std::endl;
-
-    return 1;
-  }
+    {
+      std::cerr << std::endl << std::endl
+                << "----------------------------------------------------"
+                << std::endl;
+      std::cerr << "Exception on processing: " << std::endl
+                << exc.what() << std::endl
+                << "Aborting!" << std::endl
+                << "----------------------------------------------------"
+                << std::endl;
+
+      return 1;
+    }
   catch (...)
-  {
-    std::cerr << std::endl << std::endl
-              << "----------------------------------------------------"
-              << std::endl;
-    std::cerr << "Unknown exception!" << std::endl
-              << "Aborting!" << std::endl
-              << "----------------------------------------------------"
-              << std::endl;
-    return 1;
-  };
+    {
+      std::cerr << std::endl << std::endl
+                << "----------------------------------------------------"
+                << std::endl;
+      std::cerr << "Unknown exception!" << std::endl
+                << "Aborting!" << std::endl
+                << "----------------------------------------------------"
+                << std::endl;
+      return 1;
+    };
 }

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.