]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Fixed a few last spelling mistakes and syntax errors
authorBruno <blais.bruno@gmail.com>
Thu, 16 Jul 2020 04:48:46 +0000 (00:48 -0400)
committerBruno Blais <blais.bruno@gmail.com>
Tue, 15 Sep 2020 19:47:35 +0000 (15:47 -0400)
examples/step-68/doc/intro.dox
examples/step-68/step-68.cc

index f3e5dd61a83698ebd339dc790ac9fe7cda15c902..871e8483fb79bb5991e38851c4b5eafb882fe856 100644 (file)
@@ -25,7 +25,7 @@ Particles play an important part in numerical models for a large
  these models can be built using deal.II's particle handling capabilities.
 
 In the present step, we use particles as massless tracer to illustrate
-the dynamic of a vortical flows. Since the particles are massless tracers,
+the dynamic of a vortical flow. Since the particles are massless tracers,
 the position of each particle $i$ is described by the
 following ordinary differential equation (ODE):
 @f[
@@ -46,8 +46,8 @@ is obtained in two different fashions:
 - By evaluating the velocity function on a background triangulation and, using
 a  finite element support, interpolating at the position of the particle.
 
-The first approach is generally not practical, since the velocity profile
-is not known analytically. The second approach, based on interpolating a solution
+The first approach is not practical, since the velocity profile
+is generally not known analytically. The second approach, based on interpolating a solution
 at the position of the particles, mimics exactly what would be done in a
 realistic computational fluid dynamic simulation. In this step, we illustrate both strategies.
 
@@ -55,14 +55,14 @@ realistic computational fluid dynamic simulation. In this step, we illustrate bo
 We note that much greater accuracy could be obtained by using a fourth
 order Runge-Kutta method or another appropriate scheme for the time integration
 of the motion of the particles.  Implementing a more advanced advection scheme
-would be a straightforward extension of this example.
+would be a straightforward extension of this step.
 
 <h3>Particles in deal.II</h3>
 
 In deal.II, Particles::Particle are very simple and flexible entities that can be used
 to build PIC, DEM or any type of particle-based models. Particles have a location
 in real space, a location in the reference space of the element in which they
-lie and a unique ID. In the majority of cases, simulations that include
+are located and a unique ID. In the majority of cases, simulations that include
 particles require a significant number of them. Thus, it becomes interesting
 to handle all particles through an entity which agglomerates all particles.
 In deal.II, this is achieved through the use of the Particles::ParticleHandler class.
@@ -70,9 +70,9 @@ In deal.II, this is achieved through the use of the Particles::ParticleHandler c
 By default, particles do not have a diameter,
 a mass or any other physical properties which we would generally expect of physical particles. However, through
 a ParticleHandler, particles have access to a Particles::PropertyPool. This PropertyPool is
-an array which can be used to store any arbitrary number of properties
+an array which can be used to store an arbitrary number of properties
 associated with the particles. Consequently, users can build their own
-particle solver and attribute the desired properties to the particles (e.g. mass, charge
+particle solver and attribute the desired properties to the particles (e.g. mass, charge,
 diameter, temperature, etc.). In the present tutorial, this is used to
 store the value of the fluid velocity and the process id to which the particles
 belong.
@@ -104,7 +104,7 @@ density function to determine how many and where to generate particles locally.
 
 In other situations, such as the present step, particles must be generated at
 specific locations on cells that may be owned only by a subset of the processors.
-In  most of these situations, the insertion of the particle is done for a very
+In  most of these situations, the insertion of the particles is done for a very
 limited number of time-steps and, consequently, does not constitute a large
 portion of the computational cost. For these occasions, deal.II provides
 convenient Particles::Generators that can globally insert the particles even if
@@ -177,12 +177,12 @@ finite-element problems already discussed in other examples.
 <h3>The testcase</h3>
 
 In the present step, we use particles as massless tracers to illustrate
-the dynamics of a particular vortical flow : the Rayleigh-Kotte Vortex. This flow pattern
+the dynamics of a particular vortical flow : the Rayleigh-Kothe Vortex. This flow pattern
 is generally used as a complex test case for interface tracking methods
 (e.g. volume-of-fluid and level set approches) since
-it leads to strong rotation and elongation of the fluid (Blais, 2013).
+it leads to strong rotation and elongation of the fluid (Blais et al. 2013).
 
-The stream function $\Psi$ of this Rayleigh-Kotte vortex is defined as:
+The stream function $\Psi$ of this Rayleigh-Kothe vortex is defined as:
 
 @f[
 \Psi = \frac{1}{\pi} sin^2 (\pi x) \sin^2 (\pi y) \cos \left( \pi \frac{t}{T} \right)
@@ -193,14 +193,12 @@ where $T$ is the period of the flow. The velocity profile in 2D ($\textbf{u}=[u,
    v &=&  \frac{\partial\Psi}{\partial x} = 2 \cos(\pi x) \sin(\pi x) \sin^2 (\pi y) \cos \left( \pi \frac{t}{T} \right)
 @f}
 
-It can be seen that this velocity reverses periodically due to the term 
+It can be seen that this velocity reverses periodically due to the term
 $\cos \left( \pi \frac{t}{T} \right)$ and that material will end up at its
 starting position after every period of length $t=2T$. We will run this tutorial
 program for exactly one period and compare the final particle location to the
-initial location to illustrate this flow property.
-
-
-This example uses the testcase to produce two models that handle the particles
+initial location to illustrate this flow property. This example uses the testcase
+to produce two models that handle the particles
 slightly differently. The first model prescribes the exact analytical velocity
 solution as the velocity for each particle. Therefore in this model there is no
 error in the assigned velocity to the particles, and any deviation of particle
index 1c186ac509bcf3c1beea2f9f69fe33691a331b45..38da7a344750caf848fb6211548e7d0a434715d0 100644 (file)
@@ -103,7 +103,7 @@ namespace Step68
     // describe the details of the particle tracking simulation and its
     // discretization. The following parameters are about where output should
     // land, the spatial discretization of the velocity (the default is $Q_1$),
-    // the time step, finally,  the output frequency (how many time steps should
+    // the time step and the output frequency (how many time steps should
     // elapse before we generate graphical output again):
     std::string output_directory = "./";
 
@@ -157,7 +157,7 @@ namespace Step68
 
   // The velocity profile is provided as a Function object. We provide the
   // velocity profile. In the present step, this function is hard-coded within
-  // the example. However, it could have been easily made using a ParsedFunction
+  // the example.
   template <int dim>
   class Vortex : public Function<dim>
   {
@@ -214,12 +214,12 @@ namespace Step68
 
     void interpolate_function_to_field();
 
-    // The next two functions are responsible for carrying out explicit Euler
-    // time integration for the cases where the velocity field is interpolated
-    // at the positions of the particles or calculated analytically,
-    // respectively
-    void euler_interpolated(double dt);
-    void euler_analytical(double dt);
+    // The next two functions are responsible for carrying out step of explicit
+    // Euler time integration for the cases where the velocity field is
+    // interpolated at the positions of the particles or calculated
+    // analytically, respectively
+    void euler_step_interpolated(double dt);
+    void euler_step_analytical(double dt);
 
     // The cell_weight() function indicates to the triangulation how much
     // computational work is expected to happen on this cell, and consequently
@@ -492,7 +492,7 @@ namespace Step68
   // using an analytically defined velocity field. This is a relatively trivial
   // usage of the particles.
   template <int dim>
-  void ParticleTracking<dim>::euler_analytical(double dt)
+  void ParticleTracking<dim>::euler_step_analytical(double dt)
   {
     Vector<double> particle_velocity(dim);
 
@@ -524,7 +524,7 @@ namespace Step68
   // We integrate the particle trajectories by interpolating the value of the
   // velocity field at the degrees of freedom to the position of the particles.
   template <int dim>
-  void ParticleTracking<dim>::euler_interpolated(double dt)
+  void ParticleTracking<dim>::euler_step_interpolated(double dt)
   {
     std::vector<types::global_dof_index> dof_indices(fluid_fe.dofs_per_cell);
     Vector<double> dof_data_per_cell(fluid_fe.dofs_per_cell);
@@ -691,10 +691,10 @@ namespace Step68
     if (interpolated_velocity)
       {
         interpolate_function_to_field();
-        euler_interpolated(0.);
+        euler_step_interpolated(0.);
       }
     else
-      euler_analytical(0.);
+      euler_step_analytical(0.);
 
     output_particles(discrete_time.get_step_number());
     output_background(discrete_time.get_step_number());
@@ -714,10 +714,10 @@ namespace Step68
         if (interpolated_velocity)
           {
             interpolate_function_to_field();
-            euler_interpolated(discrete_time.get_previous_step_size());
+            euler_step_interpolated(discrete_time.get_previous_step_size());
           }
         else
-          euler_analytical(discrete_time.get_previous_step_size());
+          euler_step_analytical(discrete_time.get_previous_step_size());
 
         particle_handler.sort_particles_into_subdomains_and_cells();
 

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.