]> https://gitweb.dealii.org/ - dealii-svn.git/commitdiff
Add a random distortion function to the grid class.
authorwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 1 Oct 1998 16:34:39 +0000 (16:34 +0000)
committerwolf <wolf@0785d39b-7218-0410-832d-ea1e28bc413d>
Thu, 1 Oct 1998 16:34:39 +0000 (16:34 +0000)
git-svn-id: https://svn.dealii.org/trunk@606 0785d39b-7218-0410-832d-ea1e28bc413d

deal.II/deal.II/include/grid/tria.h
deal.II/deal.II/source/grid/tria.cc

index b939231f13e3f3b3fe5a9c4bd70366f04cacb2cc..aaf813716e739c9714b793b240b7763cca0e81e5 100644 (file)
@@ -550,6 +550,13 @@ class TriaDimensionInfo<2> {
  *   coarse level cells so, that the interface between cells is also the
  *   interface between regions of different materials.
  *
+ *   Finally, there is a special function for folks who like bad grids:
+ *   #Triangulation<dim>::distort_random#. It moves all the vertices in the
+ *   grid a bit around by a random value, leaving behind a distorted mesh.
+ *   Note that you should apply this function to the final mesh, since refinement
+ *   smoothes the mesh a bit.
+ *
+ *
  *
  *   \subsection{Refinement and coarsening of a triangulation}
  *
@@ -1281,6 +1288,24 @@ class Triangulation : public TriaDimensionInfo<dim> {
                                      */
     void create_hyper_ball (const Point<dim> &center = Point<dim>(),
                            const double radius = 1.);
+
+                                    /**
+                                     * Distort the grid by randomly moving
+                                     * around all the vertices of the grid.
+                                     * The direction of moving is random,
+                                     * while the length of the shift vector
+                                     * has a value of #factor# times the
+                                     * average length of the active lines
+                                     * adjacent to this vertex. Note that
+                                     * #factor# should obviously be well
+                                     * below #0.5#.
+                                     *
+                                     * If #keep_boundary# is set to #true#
+                                     * (which is the default), then boundary
+                                     * vertices are not moved.
+                                     */
+    void distort_random (const double factor,
+                        const bool   keep_boundary=true);
     
                                      
                                     /**
index 2f927842ebe664b6ebf103277f2e40aa7c5ca64c..2723c83b506baa195f5eb0ad91a9fe5771eaefad 100644 (file)
@@ -668,6 +668,68 @@ void Triangulation<2>::create_hyper_ball (const Point<2> &p, const double radius
 
 
 
+template <int dim>
+void Triangulation<dim>::distort_random (const double factor,
+                                        const bool   keep_boundary) {
+                                  // note number of lines adjacent
+                                  // to a vertex
+  vector<short unsigned int> adjacent_lines (vertices.size(), 0);
+                                  // sum up the length of the lines
+                                  // adjecent to the vertex
+  vector<double>             cumulated_length (vertices.size(), 0);
+                                  // also note if a vertex is at
+                                  // the boundary
+  vector<bool>               at_boundary (vertices.size(), false);
+  
+  for (active_line_iterator line=begin_active_line();
+       line != end_line(); ++line)
+    {
+      if (keep_boundary && line->at_boundary())
+       {
+         at_boundary[line->vertex_index(0)] = true;
+         at_boundary[line->vertex_index(1)] = true;
+       };
+      
+      adjacent_lines[line->vertex_index(0)]++;
+      adjacent_lines[line->vertex_index(1)]++;
+      cumulated_length[line->vertex_index(0)] += line->diameter();
+      cumulated_length[line->vertex_index(1)] += line->diameter();
+    };
+
+
+  const unsigned int n_vertices = vertices.size();
+  Point<dim> shift_vector;
+  
+  for (unsigned int vertex=0; vertex<n_vertices; ++vertex) 
+    {
+                                      // ignore this vertex if we whall keep
+                                      // the boundary and this vertex *is* at
+                                      // the boundary
+      if (keep_boundary && at_boundary[vertex])
+       continue;
+      
+                                      // first compute a random shift vector
+      for (unsigned int d=0; d<dim; ++d)
+       shift_vector(d) = rand()*1.0/RAND_MAX;
+
+                                      // bring the random vector to the
+                                      // specified length by dividing
+                                      // through its original length, then
+                                      // multiplication by the average length
+                                      // of the adjacent lines times the given
+                                      // factor
+      shift_vector *= (factor *
+                      cumulated_length[vertex]) /
+                     (adjacent_lines[vertex] *
+                      sqrt(shift_vector.square()));
+
+                                      // finally move the vertex
+      vertices[vertex] += shift_vector;
+    };
+};
+
+
+
 
 template <int dim>
 void Triangulation<dim>::set_all_refine_flags () {

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.