Fixed codimension one bug in manifold case.

git-svn-id: https://svn.dealii.org/branches/branch_manifold_id_intermediate@32449 0785d39b-7218-0410-832d-ea1e28bc413d

diff --git a/deal.II/source/grid/tria.cc b/deal.II/source/grid/tria.cc
index ac23fb319f6afaf5581f925a458a46c56ccb7e9a..11188685fabb27c88667e798ae6c8be650bb05a4 100644 (file)
--- a/deal.II/source/grid/tria.cc
+++ b/deal.II/source/grid/tria.cc
@@ -4054,8 +4054,13 @@ namespace internal
                 // to the information
                 // stored in the
                 // boundary class
+               const Boundary<dim,spacedim> &manifold = 
+                 (cell->manifold_id() != numbers::invalid_manifold_id ?
+                  cell->get_manifold() :
+                  triangulation.get_boundary(static_cast<types::boundary_id>(cell->material_id())));
+               
                 triangulation.vertices[next_unused_vertex] =
-                  cell->get_boundary().get_new_point_on_quad(cell);
+                  manifold.get_new_point_on_quad (cell);
               }
           }
 
@@ -4441,14 +4446,20 @@ namespace internal
                   // midpoint; otherwise we
                   // have to ask the manifold
                   // object
-                  if (dim == spacedim)
-                    triangulation.vertices[next_unused_vertex] =
-                      (cell->vertex(0) + cell->vertex(1)) / 2;
-                  else
-                    triangulation.vertices[next_unused_vertex] =
-                     cell->get_boundary().get_new_point_on_line(cell);
-                  triangulation.vertices_used[next_unused_vertex] = true;
-
+                 
+                 // Now we always ask the manifold where to put the
+                 // new point. The get_boundary function will return
+                 // a flat boundary if invalid_manifold_id is set,
+                 // behaving as before in the flat case.  Backward
+                 // compatibility requires us to use the material id
+                 // of the cell in the codimension one case, however
+                 // we only do this if the manifold_id is the invalid
+                 // one, otherwise use the material id. This is done
+                 // internally in the get_boundary() function.
+                 triangulation.vertices[next_unused_vertex] =
+                   cell->get_boundary().get_new_point_on_line(cell);
+                 triangulation.vertices_used[next_unused_vertex] = true;
+                    
                   // search for next two
                   // unused cell (++ takes
                   // care of the end of the
@@ -4865,21 +4876,23 @@ namespace internal
                       // lines we can compute the
                       // midpoint as the mean of
                       // the two vertices:
-                      if (line->at_boundary())
+                      // if (line->at_boundary())
                         triangulation.vertices[next_unused_vertex]
                           = line->get_boundary().get_new_point_on_line (line);
-                      else
-                        triangulation.vertices[next_unused_vertex]
-                          = (line->vertex(0) + line->vertex(1)) / 2;
                     }
                   else
-                    // however, if spacedim>dim, we
-                    // always have to ask the
-                    // boundary object for its
-                    // answer
-                    triangulation.vertices[next_unused_vertex]
-                      = triangulation.get_boundary(line->user_index()).get_new_point_on_line (line);
-
+                    // however, if spacedim>dim, we always have to ask
+                    // the boundary object for its answer. We use the
+                    // same object of the cell (which was stored in
+                    // line->user_index() before) unless a manifold_id
+                    // has been set on this very line.
+                   if(line->manifold_id() == numbers::invalid_manifold_id)
+                     triangulation.vertices[next_unused_vertex]
+                       = triangulation.get_boundary(line->user_index()).get_new_point_on_line (line);
+                   else
+                     triangulation.vertices[next_unused_vertex]
+                       = line->get_boundary().get_new_point_on_line (line);
+                     
                   // now that we created
                   // the right point, make
                   // up the two child
@@ -6014,7 +6027,7 @@ namespace internal
                             // care about boundary quads
                             // here
                             triangulation.vertices[next_unused_vertex]
-                              = (middle_line->vertex(0) + middle_line->vertex(1)) / 2;
+                              = middle_line->get_boundary().get_new_point_on_line(middle_line);
                             triangulation.vertices_used[next_unused_vertex] = true;
 
                             // now search a slot for the two
@@ -6080,7 +6093,8 @@ namespace internal
 
                     // set the middle vertex
                     // appropriately
-                    if (quad->at_boundary())
+                    if (quad->at_boundary() || 
+                       (quad->manifold_id() != numbers::invalid_manifold_id) )
                       triangulation.vertices[next_unused_vertex]
                        = quad->get_boundary().get_new_point_on_quad (quad);
                     else

diff --git a/tests/manifold/flat_manifold_01.cc b/tests/manifold/flat_manifold_01.cc
deleted file mode 100644 (file)
index 9d463df..0000000
--- a/tests/manifold/flat_manifold_01.cc
+++ /dev/null
@@ -1,55 +0,0 @@
-//----------------------------  manifold_id_01.cc  ---------------------------
-//    Copyright (C) 2011, 2013 by the mathLab team.
-//
-//    This file is subject to LGPL 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.
-//
-//----------------------------  flat_manifold_01.cc  ---------------------------
-
-
-// Test periodicity conditions using FlatManifold
-
-#include "../tests.h"
-#include <fstream>
-#include <base/logstream.h>
-
-
-// all include files you need here
-#include <deal.II/grid/manifold.h>
-#include <deal.II/grid/manifold_lib.h>
-
-template <int spacedim> 
-void test() {
-  Point<spacedim> period;
-  for(unsigned int i=0; i<spacedim; ++i)
-    period[i] = 1.0;
-  FlatManifold<spacedim> manifold(period);
-  deallog << "Testing spacedim = " << spacedim << std::endl;
-  deallog << "Period: " << period << std::endl;
-  std::vector<Point<spacedim> > p(2);
-  std::vector<double> w(2, 0.5);
-  double eps = .1;
-  p[0] = period*eps*2;
-  p[1] = period*(1.0-eps);
-  Point<spacedim> middle = manifold.get_new_point(p, w);
-  deallog << "P(0): " << p[0] << std::endl;
-  deallog << "P(1): " << p[1] << std::endl;
-  deallog << "Middle: " << middle << std::endl;
-}
-
-int main () 
-{
-  std::ofstream logfile("output");
-  deallog.attach(logfile);
-  deallog.depth_console(0);
-  deallog.threshold_double(1.e-5);
-  
-  test<1>();
-  test<2>();
-  test<3>();
-
-  return 0;
-}
-                  

diff --git a/tests/manifold/flat_manifold_01.output b/tests/manifold/flat_manifold_01.output
deleted file mode 100644 (file)
index f6f885f..0000000
--- a/tests/manifold/flat_manifold_01.output
+++ /dev/null
@@ -1,16 +0,0 @@
-
-DEAL::Testing spacedim = 1
-DEAL::Period: 1.00000
-DEAL::P(0): 0.200000
-DEAL::P(1): 0.900000
-DEAL::Middle: 0.0500000
-DEAL::Testing spacedim = 2
-DEAL::Period: 1.00000 1.00000
-DEAL::P(0): 0.200000 0.200000
-DEAL::P(1): 0.900000 0.900000
-DEAL::Middle: 0.0500000 0.0500000
-DEAL::Testing spacedim = 3
-DEAL::Period: 1.00000 1.00000 1.00000
-DEAL::P(0): 0.200000 0.200000 0.200000
-DEAL::P(1): 0.900000 0.900000 0.900000
-DEAL::Middle: 0.0500000 0.0500000 0.0500000

diff --git a/tests/manifold/manifold_id_05.cc b/tests/manifold/manifold_id_05.cc
index 3d94f71390355a91cdde09508216f6c91e01abe8..e6a5a52b558d2182ee19ab55555dee206e474dfb 100644 (file)
--- a/tests/manifold/manifold_id_05.cc
+++ b/tests/manifold/manifold_id_05.cc
@@ -22,7 +22,7 @@
 #include <deal.II/grid/tria_accessor.h>
 #include <deal.II/grid/tria_iterator.h>
 #include <deal.II/grid/grid_generator.h>
-#include <deal.II/grid/manifold_lib.h>
+#include <deal.II/grid/tria_boundary_lib.h>
 #include <deal.II/grid/grid_out.h>
 
 // Helper function
@@ -35,13 +35,13 @@ void test(unsigned int ref=1){
        for(unsigned int i=0; i<dim; ++i)
          center[i] = .5;
 
-       PolarManifold<spacedim> boundary(center);
+       HyperBallBoundary<dim,spacedim> boundary(center, .5*std::sqrt((double)dim));
        Triangulation<dim,spacedim> tria;
        GridGenerator::hyper_cube (tria);
        typename Triangulation<dim,spacedim>::active_cell_iterator cell;
        
        tria.begin_active()->face(0)->set_manifold_id(1);
-       tria.set_manifold(1,boundary);
+       tria.set_boundary(1,boundary);
        
        tria.refine_global(2);
 

diff --git a/tests/manifold/manifold_id_06.cc b/tests/manifold/manifold_id_06.cc
index ccf108b4ff663798ea7092e9075e1da18ef3169e..8fa4e892ca284764c533b840c829520b93a1c83f 100644 (file)
--- a/tests/manifold/manifold_id_06.cc
+++ b/tests/manifold/manifold_id_06.cc
@@ -22,8 +22,8 @@
 #include <deal.II/grid/tria_accessor.h>
 #include <deal.II/grid/tria_iterator.h>
 #include <deal.II/grid/grid_generator.h>
-#include <deal.II/grid/manifold_lib.h>
 #include <deal.II/grid/grid_out.h>
+#include <deal.II/grid/tria_boundary_lib.h>
 
 // Helper function
 template <int dim, int spacedim>
@@ -37,7 +37,7 @@ void test(unsigned int ref=1){
 
        double radius=center.norm();
        
-       PolarManifold<spacedim> boundary(center);       
+       HyperBallBoundary<dim,spacedim> boundary(center, .25*std::sqrt((double)spacedim));      
        Triangulation<dim,spacedim> tria;
        GridGenerator::hyper_cube (tria);
        typename Triangulation<dim,spacedim>::active_cell_iterator cell;
@@ -52,7 +52,7 @@ void test(unsigned int ref=1){
              if(cell->face(f)->center().distance(center)< radius)
                cell->face(f)->set_all_manifold_ids(1);
 
-       tria.set_manifold(1,boundary);
+       tria.set_boundary(1,boundary);
        tria.refine_global(2);
 
        GridOut gridout;

diff --git a/tests/manifold/manifold_id_07.cc b/tests/manifold/manifold_id_07.cc
index 943033102516ca51244cf5a4450bdbb9b3df929e..9e0662f45aee616e50daf2a495d6e73abb2af9df 100644 (file)
--- a/tests/manifold/manifold_id_07.cc
+++ b/tests/manifold/manifold_id_07.cc
@@ -24,7 +24,7 @@
 #include <deal.II/grid/tria_accessor.h>
 #include <deal.II/grid/tria_iterator.h>
 #include <deal.II/grid/grid_generator.h>
-#include <deal.II/grid/manifold_lib.h>
+#include <deal.II/grid/tria_boundary_lib.h>
 #include <deal.II/grid/grid_out.h>
 
 // Helper function
@@ -38,12 +38,12 @@ void test(unsigned int ref=1){
          center[i] = .5;
 
        Triangulation<dim,spacedim> tria;
-       PolarManifold<spacedim> boundary(center);       
+       HyperBallBoundary<dim,spacedim> boundary(center,.5*std::sqrt((double)dim));     
        GridGenerator::hyper_cube (tria);
        typename Triangulation<dim,spacedim>::active_cell_iterator cell;
        
        tria.begin_active()->face(0)->set_manifold_id(1);
-       tria.set_manifold(1,boundary);
+       tria.set_boundary(1,boundary);
        
        tria.refine_global(2);