]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add test denis_1.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 12 Nov 2002 16:08:33 +0000 (16:08 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Tue, 12 Nov 2002 16:08:33 +0000 (16:08 +0000)
git-svn-id: https://svn.dealii.org/trunk@6756 0785d39b-7218-0410-832d-ea1e28bc413d

tests/bits/Makefile
tests/bits/denis_1.cc [new file with mode: 0644]

index ff2ea93248e79fdf7968cf1d0ccb21cdd40c89ce..6d7a678a19cf489ac880d79b4ee118d59b75faa7 100644 (file)
@@ -28,10 +28,12 @@ anna_3.exe              : anna_3.g.$(OBJEXT)              $(libraries)
 anna_5.exe              : anna_5.g.$(OBJEXT)              $(libraries)
 gerold_1.exe            : gerold_1.g.$(OBJEXT)            $(libraries)
 roy_1.exe               : roy_1.g.$(OBJEXT)               $(libraries)
+denis_1.exe             : denis_1.g.$(OBJEXT)             $(libraries)
 
 tests = anna_1 anna_2 anna_3 anna_5 \
        gerold_1 \
-       roy_1
+       roy_1 \
+        denis_1
 
 ############################################################
 
diff --git a/tests/bits/denis_1.cc b/tests/bits/denis_1.cc
new file mode 100644 (file)
index 0000000..457957c
--- /dev/null
@@ -0,0 +1,96 @@
+//----------------------------  denis_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.
+//
+//----------------------------  denis_1.cc  ---------------------------
+
+
+// check for a bug in DerivativeApproximation::approximate_gradient
+//
+// this program is a modified version of one by
+// Denis Danilov <danilovdenis@yandex.ru>, 
+
+
+#include <base/logstream.h>
+#include <lac/vector.h>
+#include <grid/tria.h>
+#include <dofs/dof_handler.h>
+#include <grid/grid_generator.h>
+#include <grid/tria_iterator.h>
+#include <dofs/dof_accessor.h>
+#include <numerics/data_out.h>
+#include <fe/fe_q.h>
+#include <numerics/derivative_approximation.h>
+#include <fstream>
+#include <iostream>
+
+int main() 
+{
+  std::ofstream logfile("denis_1.output");
+  deallog.attach(logfile);
+  deallog.depth_console(0);
+
+  Triangulation<2>   triangulation;
+  FE_Q<2>            fe(2);
+  DoFHandler<2>      dof_handler(triangulation);
+  Vector<double>       phi_solution;
+  Vector<float>        gradient_phi;
+  float                gradient_phi_min, gradient_phi_max;
+  double delta = 0.05;
+
+  GridGenerator::hyper_cube(triangulation, 0, 1);
+  triangulation.refine_global(5);
+
+  dof_handler.distribute_dofs(fe);
+
+  DoFHandler<2>::active_cell_iterator cell = dof_handler.begin_active(), endc = dof_handler.end();
+  
+  phi_solution.reinit(dof_handler.n_dofs());
+
+  for(; cell!=endc; ++cell)
+    {
+      for(unsigned int vertex=0;
+          vertex < GeometryInfo<2>::vertices_per_cell;
+          ++vertex)                               
+       {
+         double x, y, r;
+         x = cell->vertex(vertex)(0);
+         y = cell->vertex(vertex)(1);
+         r = sqrt( x*x + y*y );
+         phi_solution( cell->vertex_dof_index(vertex,0) )
+           = 0.5 * ( 1 - tanh( ( r - 0.5 )
+                               /( 2*M_SQRT2 * delta ) ) ) ;
+       }
+    }
+
+  gradient_phi.reinit(triangulation.n_active_cells());
+  DerivativeApproximation::approximate_gradient(dof_handler, phi_solution, gradient_phi);
+
+  gradient_phi_min = 1e30;
+  gradient_phi_max = -1;
+
+  cell = dof_handler.begin_active();  
+  for(unsigned int cell_no=0; cell!=endc; ++cell, ++cell_no)
+    {
+      if( gradient_phi(cell_no) < gradient_phi_min )
+       gradient_phi_min =  gradient_phi(cell_no);
+      
+      if( gradient_phi(cell_no) > gradient_phi_max )
+       gradient_phi_max =  gradient_phi(cell_no);
+      
+    }
+
+  deallog << "gradient_phi_min: " << gradient_phi_min << std::endl;
+  deallog << "gradient_phi_max: " << gradient_phi_max << std::endl;
+
+  //////////////////////
+  
+  
+};

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.