]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fixed step-34.
authorLuca Heltai <luca.heltai@gmail.com>
Tue, 31 Mar 2015 19:07:32 +0000 (21:07 +0200)
committerLuca Heltai <luca.heltai@sissa.it>
Tue, 7 Apr 2015 14:02:16 +0000 (16:02 +0200)
examples/step-34/step-34.cc

index 09641afb9dc6ff34b63e19c4a8f9c75e08c5bb57..84e6c1f8c425576bbef4a0d1249de0213b454b80 100644 (file)
@@ -42,7 +42,7 @@
 #include <deal.II/grid/grid_generator.h>
 #include <deal.II/grid/grid_in.h>
 #include <deal.II/grid/grid_out.h>
-#include <deal.II/grid/tria_boundary_lib.h>
+#include <deal.II/grid/manifold_lib.h>
 
 #include <deal.II/dofs/dof_handler.h>
 #include <deal.II/dofs/dof_accessor.h>
@@ -183,7 +183,7 @@ namespace Step34
     //
     // Experimenting a little with the computation of the angles gives very
     // accurate results for simpler geometries. To verify this you can comment
-    // out, in the read_domain() method, the tria.set_boundary(1, boundary)
+    // out, in the read_domain() method, the tria.set_manifold(1, manifold)
     // line, and check the alpha that is generated by the program. By removing
     // this call, whenever the mesh is refined new nodes will be placed along
     // the straight lines that made up the coarse mesh, rather than be pulled
@@ -481,21 +481,22 @@ namespace Step34
   // normals to the mesh are external to both the circle in 2d or the sphere
   // in 3d.
   //
-  // The other detail that is required for appropriate refinement of the
-  // boundary element mesh, is an accurate description of the manifold that
-  // the mesh is approximating. We already saw this several times for the
-  // boundary of standard finite element meshes (for example in step-5 and
-  // step-6), and here the principle and usage is the same, except that the
-  // HyperBallBoundary class takes an additional template parameter that
-  // specifies the embedding space dimension. The function object still has to
-  // be static to live at least as long as the triangulation object to which
-  // it is attached.
+  // The other detail that is required for appropriate refinement of
+  // the boundary element mesh, is an accurate description of the
+  // manifold that the mesh is approximating. We already saw this
+  // several times for the boundary of standard finite element meshes
+  // (for example in step-5 and step-6), and here the principle and
+  // usage is the same, except that the SphericalManifold class takes
+  // an additional template parameter that specifies the embedding
+  // space dimension. The function object still has to be static to
+  // live at least as long as the triangulation object to which it is
+  // attached.
 
   template <int dim>
   void BEMProblem<dim>::read_domain()
   {
     static const Point<dim> center = Point<dim>();
-    static const HyperBallBoundary<dim-1, dim> boundary(center,1.);
+    static const SphericalManifold<dim-1, dim> manifold(center);
 
     std::ifstream in;
     switch (dim)
@@ -515,8 +516,9 @@ namespace Step34
     GridIn<dim-1, dim> gi;
     gi.attach_triangulation (tria);
     gi.read_ucd (in);
-
-    tria.set_boundary(1, boundary);
+    
+    tria.set_all_manifold_ids(1);
+    tria.set_manifold(1, manifold);
   }
 
 

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.