]> https://gitweb.dealii.org/ - dealii.git/commitdiff
Improve documentation
authorDaniel Arndt <arndtd@ornl.gov>
Fri, 10 May 2019 21:13:21 +0000 (17:13 -0400)
committerDaniel Arndt <arndtd@ornl.gov>
Fri, 10 May 2019 21:13:21 +0000 (17:13 -0400)
examples/step-64/doc/intro.dox
examples/step-64/step-64.cu

index 5350f84abdd9719e96c45fbea4d3e7b2cae69ef6..a7d69b296f473042f029e062244a7957dfa6cd69 100644 (file)
@@ -29,11 +29,11 @@ amount is fairly small and the use of CUDA is limited to a few keywords.
 <h3>The test case</h3>
 
 In this example, we consider the Helmholtz problem @f{eqnarray*} - \nabla \cdot
-\nabla u + a(\mathbf r) u &=&1,\\ u &=& 0 \quad \text{on} \partial \Omega @f}
+\nabla u + a(\mathbf x) u &=&1,\\ u &=& 0 \quad \text{on } \partial \Omega @f}
 where $a(\mathbf x)$ is a variable coefficient.
 
 We choose as domain $\Omega=[0,1]^3$ and $a(\mathbf x)=\frac{10}{0.05 +
-2\|mathbf r\|^2}$, Since the coefficient is symmetric around the origin but
+2\|\mathbf x\|^2}$, Since the coefficient is symmetric around the origin but
 the domain is not, we will end up with a non-symmetric solution.
 
 <h3>Moving data to and from the device</h3>
@@ -42,28 +42,32 @@ GPUs (we will use device from now on to refer to the GPU) have their own memory
 that is separate from the memory accessible to the CPU (we will use the term
 "host" from now on). A normal calculation on the device can be divided in three
 separate steps:
- 1) the data is moved from the host to the device
- 2) the computation is done on the device
- 3) the result is moved back from the device to the host
+ -# the data is moved fromi the host to the device,
+ -# the computation is done on the device,
+ -# the result is moved back from the device to the host
+
 The data movements can either be done explicitly by the user code or done
 automatically using UVM (Unified Virtual Memory). In deal.II, only the first
 method is supported. While it means an extra burden for the user, it allows a
 better control of data movement and more importantly it avoids to mistakenly run
 important kernels on the host instead of the device.
 
-The data movement in deal.II is done using
-LinearAlgebra::ReadWriteVector<Number>. These vectors can be seen as buffers on
-the host that are used to either store data from the device or to send data to
-the device. There are two types of vectors that can be used on the device:
-LinearAlgebra::CUDAWrappers::Vector<Number>, which is similar to the more common
-Vector<Number>, and LinearAlgebra::distributed::Vector<Number,
+The data movement in deal.II is done using LinearAlgebra::ReadWriteVector. These
+vectors can be seen as buffers on the host that are used to either store data
+from the device or to send data to the device. There are two types of vectors
+that can be used on the device:
+- LinearAlgebra::CUDAWrappers::Vector, which is similar to the more common
+Vector<Number>, and
+- LinearAlgebra::distributed::Vector<Number,
 MemorySpace::CUDA>, which is a regular
-LinearAlgebra::distributed::Vector<Number> where we have specified which memory
-space to use. If no memory space is specified, the default is MemorySpace::Host.
+LinearAlgebra::distributed::Vector where we have specified which memory
+space to use.
+
+If no memory space is specified, the default is MemorySpace::Host.
 
 Next, we show how to move data to/from the device using
-LinearAlgebra::CUDAWrappers::Vector<Number>:
-<code>
+LinearAlgebra::CUDAWrappers::Vector:
+@code
 unsigned int size = 10;
 LinearAlgebra::CUDAWrappers::Vector<double> vector_dev(10);
 LinearAlgebra::ReadWriteVector<double> rw_vector(10);
@@ -73,10 +77,11 @@ vector_dev.import(rw_vector, VectorOperations::insert);
 // Do some computations on the device
 // Move the data to the host
 rw_vector.import(vector_dev, VectorOperations::insert);
-</code>
-Using LinearAlgebra::distributed::Vector<Number, MemorySpace::CUDA> is similar
-but the import() stage may involve MPI communication:
-<code>
+@endcode
+Using LinearAlgebra::distributed::Vector<Number,
+MemorySpace::CUDA> is similar
+but the `import()` stage may involve MPI communication:
+@code
 IndexSet locally_owned_dofs, locally_relevant_dofs;
 // Fill the two IndexSet...
 LinearAlgebra::distributed::Vector<double, MemorySpace::CUDA>
@@ -91,15 +96,15 @@ distributed_vector_dev.import(owned_rw_vector, VectorOperations::insert);
 LinearAlgebra::ReadWriteVector<double> relevant_rw_vector(locally_relevant_dofs);
 // Move the data to the host and do an MPI communication
 relevnt_rw_vector(distributed_vector_dev, VectorOperations::insert);
-</code>
+@endcode
 
-import() supports two kinds of strategies: VectorOperations::insert and
-VectorOperations::add.
+While importing a vector, values can either by inserted (using VectorOperation::insert)
+or added (using VectorOperation::add).
 
 
 <h3>Matrix-vector product implementation</h3>
 
 The code necessary to evaluate the matrix-free operator on the device is very
 similar to the one on the host. However, there are a few differences, the main
-ones being that the local_apply() function in Step-37 and the loop over
+ones being that the `local_apply()` function in Step-37 and the loop over
 quadrature points both need to be encapsulated in their own functors.
index 8c727e52154206d034559571fbb8f8488ca4e9b0..04453d555d52b1c19529caea021d4dad668dc09a 100644 (file)
@@ -53,7 +53,7 @@ namespace Step64
 
 
   // Define a class that implements the varying coefficients we want to use in
-  // the Helmholtzoperator.
+  // the Helmholtz operator.
   // Later, we want to pass an object of this type to a
   // CUDAWrappers::MatrixFree<dim, double> object that expects the class to have
   // an operator that fills the values provided in the constructor for a given
@@ -92,7 +92,7 @@ namespace Step64
   {
     const unsigned int pos = CUDAWrappers::local_q_point_id<dim, double>(
       cell, gpu_data, n_dofs_1d, n_q_points);
-    const auto q_point =
+    const Tensor<1, dim> q_point =
       CUDAWrappers::get_quadrature_point<dim, double>(cell,
                                                       gpu_data,
                                                       n_dofs_1d);
@@ -107,9 +107,9 @@ namespace Step64
   }
 
 
-  // The class HelmholtzOperatorQuad implements the evaluation of the Helmholtz
-  // operator in each quadrature point. It uses a similar mechanism as the
-  // MatrixFree framework introduced in step-37.
+  // The class `HelmholtzOperatorQuad` implements the evaluation of the
+  // Helmholtz operator at each quadrature point. It uses a similar mechanism as
+  // the MatrixFree framework introduced in step-37.
   template <int dim, int fe_degree>
   class HelmholtzOperatorQuad
   {
@@ -127,25 +127,25 @@ namespace Step64
   };
 
 
-  // The Helmholtz operator reads
-  // \begin{align*}
-  // -\Delta u + c\cdot u
-  // \end{align*}
-  // and consists of two parts that are correspond to the two function calls
+  // The Helmholtz problem reads in weak form
+  // @f{eqnarray*}
+  //   (\nabla v, \nabla u)+ (v, a(\mathbf x) u) &=&(v,1) \quad \forall v.
+  // @f}
+  // The two terms on the left-hand side correspond to the two function calls
   // here.
   template <int dim, int fe_degree>
   __device__ void HelmholtzOperatorQuad<dim, fe_degree>::
                   operator()(CUDAWrappers::FEEvaluation<dim, fe_degree> *fe_eval,
              const unsigned int                          q) const
   {
-    fe_eval->submit_value(/*coef **/ fe_eval->get_value(q), q);
+    fe_eval->submit_value(coef * fe_eval->get_value(q), q);
     fe_eval->submit_gradient(fe_eval->get_gradient(q), q);
   }
 
 
   // Finally, we need to define a class that implements the whole operator
-  // evaluation that corresponds to matrix-vector product in matrix-based
-  // approaches. It corresponds
+  // evaluation that corresponds to matrix-vector product in matrix-based
+  // approaches.
   template <int dim, int fe_degree>
   class LocalHelmholtzOperator
   {
@@ -174,9 +174,10 @@ namespace Step64
 
 
   // This is the call operator that performs the Helmholtz operator evaluation
-  // on a given cell similar to the MatrixFree framework. In particular, we need
-  // access to both values and gradients of the source vector and we write value
-  // and gradient information to the destincation vector.
+  // on a given cell similar to the MatrixFree framework on the CPU.
+  // In particular, we need access to both values and gradients of the source
+  // vector and we write value and gradient information to the destination
+  // vector.
   template <int dim, int fe_degree>
   __device__ void LocalHelmholtzOperator<dim, fe_degree>::operator()(
     const unsigned int                                          cell,
@@ -272,7 +273,6 @@ namespace Step64
   {
   public:
     HelmholtzProblem();
-    ~HelmholtzProblem();
 
     void run();
 
@@ -289,8 +289,8 @@ namespace Step64
 
     parallel::distributed::Triangulation<dim> triangulation;
 
-    DoFHandler<dim> dof_handler;
     FE_Q<dim>       fe;
+    DoFHandler<dim> dof_handler;
 
     IndexSet locally_owned_dofs;
     IndexSet locally_relevant_dofs;
@@ -326,21 +326,13 @@ namespace Step64
   HelmholtzProblem<dim, fe_degree>::HelmholtzProblem()
     : mpi_communicator(MPI_COMM_WORLD)
     , triangulation(mpi_communicator)
-    , dof_handler(triangulation)
     , fe(fe_degree)
+    , dof_handler(triangulation)
     , pcout(std::cout, Utilities::MPI::this_mpi_process(mpi_communicator) == 0)
   {}
 
 
 
-  template <int dim, int fe_degree>
-  HelmholtzProblem<dim, fe_degree>::~HelmholtzProblem()
-  {
-    dof_handler.clear();
-  }
-
-
-
   template <int dim, int fe_degree>
   void HelmholtzProblem<dim, fe_degree>::setup_system()
   {
@@ -429,11 +421,12 @@ namespace Step64
 
 
   // This solve() function finally contains the calls to the new classes
-  // previously dicussed. Here we don't use any preconditioner, i.e. the
-  // identity, to focus just on the pecuiarities of the CUDA MatrixFree
-  // framework. Of course, in a real application the choice of a suitable
-  // preconditioner is crucial but we have at least the same restructions as in
-  // step-37 since matrix entries are computed on the fly and not stored.
+  // previously dicussed. Here we don't use any preconditioner, i.e.
+  // precondition by the identity, to focus just on the peculiarities of the
+  // CUDAWrappers::MatrixFree framework. Of course, in a real application the
+  // choice of a suitable preconditioner is crucial but we have at least the
+  // same restructions as in step-37 since matrix entries are computed on the
+  // fly and not stored.
   template <int dim, int fe_degree>
   void HelmholtzProblem<dim, fe_degree>::solve()
   {
@@ -497,8 +490,8 @@ namespace Step64
   }
 
 
-  // Nothing surprising in the run function as well. We simply compute the
-  // solution on a series of (globally) refined meshes.
+  // There is nothing surprising in the `run()` function either. We simply
+  // compute the solution on a series of (globally) refined meshes.
   template <int dim, int fe_degree>
   void HelmholtzProblem<dim, fe_degree>::run()
   {

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.