]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
boundary residuals
authorkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 8 Aug 2012 14:28:23 +0000 (14:28 +0000)
committerkanschat <kanschat@0785d39b-7218-0410-832d-ea1e28bc413d>
Wed, 8 Aug 2012 14:28:23 +0000 (14:28 +0000)
git-svn-id: https://svn.dealii.org/trunk@25772 0785d39b-7218-0410-832d-ea1e28bc413d

tests/integrators/laplacian_01.cc
tests/integrators/laplacian_01/cmp/generic

index 4836a589d41fd4cbf96a76b0050018cb873ea698..6ce7df48931ebd248672751e8b91ecbe5d1f2837 100644 (file)
@@ -40,23 +40,37 @@ template <int dim>
 void test_boundary(const FEValuesBase<dim>& fev)
 {
   const unsigned int n = fev.dofs_per_cell;
+  unsigned int d=fev.get_fe().n_components();
+  
   FullMatrix<double> M(n,n);
+  
   Vector<double> u(n), v(n), w(n);
+  std::vector<std::vector<double> >
+    uval    (d,std::vector<double>(fev.n_quadrature_points)),
+    null_val(d,std::vector<double>(fev.n_quadrature_points, 0.));
+  std::vector<std::vector<Tensor<1,dim> > >
+    ugrad   (d,std::vector<Tensor<1,dim> >(fev.n_quadrature_points));
+  
+  std::vector<unsigned int> indices(n);
+  for (unsigned int i=0;i<n;++i)
+    indices[i] = i;
   
   nitsche_matrix(M, fev, 17);
   M.print(deallog,8);
-
-  // deallog << "Residuals";
-  // for (unsigned int i=0;i<n;++i)
-  //   {
-  //     u = 0.;
-  //     u(i) = 1.;
-  //     M.vmult(v,u);
-  //     nitsche_residual(w, fev, 17);
-  //     w.add(-1., v);
-  //     deallog << ' ' << w.l2_norm();
-  //   }
-  // deallog << std::endl;
+  
+  deallog << "Residuals";
+  for (unsigned int i=0;i<n;++i)
+    {
+      u = 0.;
+      u(i) = 1.;
+      fev.get_function_values(u, indices, uval, true);
+      fev.get_function_gradients(u, indices, ugrad, true);
+      nitsche_residual(w, fev, make_slice(uval), make_slice(ugrad), make_slice(null_val), 17);
+      M.vmult(v,u);
+      w.add(-1., v);
+      deallog << ' ' << w.l2_norm();
+    }
+  deallog << std::endl;
 }
 
 
index 476f291320b8d34c3e7e633a367cd9354524c717..0208202d151bd6efa649473d81aa670fa49448d9 100644 (file)
@@ -11,21 +11,25 @@ DEAL::11.     0.33    5.3     0.17
 DEAL::0.33    0       0.17    0       
 DEAL::5.3     0.17    11.     0.33    
 DEAL::0.17    0       0.33    0       
+DEAL::Residuals 0 0 0 0
 DEAL::boundary_matrix 1
 DEAL::0       0.33    0       0.17    
 DEAL::0.33    11.     0.17    5.3     
 DEAL::0       0.17    0       0.33    
 DEAL::0.17    5.3     0.33    11.     
+DEAL::Residuals 0 0 0 0
 DEAL::boundary_matrix 2
 DEAL::11.     5.3     0.33    0.17    
 DEAL::5.3     11.     0.17    0.33    
 DEAL::0.33    0.17    0       0       
 DEAL::0.17    0.33    0       0       
+DEAL::Residuals 0 0 0 0
 DEAL::boundary_matrix 3
 DEAL::0       0       0.33    0.17    
 DEAL::0       0       0.17    0.33    
 DEAL::0.33    0.17    11.     5.3     
 DEAL::0.17    0.33    5.3     11.     
+DEAL::Residuals 0 0 0 0
 DEAL::face_matrix 0
 DEAL::M11
 DEAL::0       0.17    0       0.083   

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.