From df331222e5b477030f291dec2c60a7b932dccb6c Mon Sep 17 00:00:00 2001 From: turcksin Date: Wed, 21 May 2014 15:53:04 +0000 Subject: [PATCH] Improve documentation for step-52. git-svn-id: https://svn.dealii.org/trunk@32958 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/examples/step-52/doc/intro.dox | 13 +++++------- deal.II/examples/step-52/doc/results.dox | 2 +- deal.II/examples/step-52/step-52.cc | 25 ++++++++++++------------ 3 files changed, 19 insertions(+), 21 deletions(-) diff --git a/deal.II/examples/step-52/doc/intro.dox b/deal.II/examples/step-52/doc/intro.dox index ac7c4e0329..4a7ba77555 100644 --- a/deal.II/examples/step-52/doc/intro.dox +++ b/deal.II/examples/step-52/doc/intro.dox @@ -21,13 +21,12 @@ fissible and therefore, the neutron flux satisfies the following equation: augmented by appropriate boundary conditions. Here, $v$ is the velocity of neutrons, $D$ is the diffusion coefficient, $\Sigma_a$ is the absorption cross section, and $S$ is a source. Because we are only interested in the -time dependence, we assume that $D$ and $\Sigma_a$ are constant. In this -example, we are only interested in the error in time. The domain is square +time dependence, we assume that $D$ and $\Sigma_a$ are constant. The domain is square $[0,b]\times[0,b]$ and we are looking for a solution of the form: @f{eqnarray*} \phi(x,t) = A\sin(\omega t)(bx-x^2). @f} -By using quadratic finite elements, we will not have any spatial error and all +By using quadratic finite elements, there will not have any spatial error and all the error will come from the time discretization. We impose the following boundary conditions: homogeneous Dirichlet fo $x=0$ and $x=b$ and homogeneous Neumann conditions for $y=0$ and $y=b$. The source is @@ -36,10 +35,8 @@ given by: S=A\left(\frac{1}{v}\omega \cos(\omega t)(bx -x^2) + \sin(\omega t) \left(\Sigma_a (bx-x^2)+2D\right) \right). @f} -Because the solution is a sine, we know that -$\phi\left(x,\frac{\pi}{\omega}\right) = 0$. Therefore, we can easily -compute the error at this time since it is simply the norm of the solution -found. +Because the solution is a sine, we know that $\phi\left(x,\pi\right) = 0$. +Therefore, the error at this time is simply the norm of the numerical solution.

Runge-Kutta

@@ -96,7 +93,7 @@ To simplify the problem, the domain is two dimensional and the mesh is uniform (there is no need to adapt the mesh since we use quadratic finite elements and the exact solution is quadratic). Going from a two dimensional domain to a three dimensional domain is not very challenging. However if the -mesh must be adapted, we cannot forget to: +mesh must be adapted, it is important to note to remember to:
  1. project the solution to the new mesh when the mesh is changed. The mesh used should be the same at the beginning and at the end of the time step. diff --git a/deal.II/examples/step-52/doc/results.dox b/deal.II/examples/step-52/doc/results.dox index 4a70bbe4bc..29fe0f6022 100644 --- a/deal.II/examples/step-52/doc/results.dox +++ b/deal.II/examples/step-52/doc/results.dox @@ -1,6 +1,6 @@

    Results

    -The output of this program consist of the console output and solutions given in +The output of this program consists of the console output and the solutions given in vtu format. The console output is: diff --git a/deal.II/examples/step-52/step-52.cc b/deal.II/examples/step-52/step-52.cc index f49dfb75d1..372727ee82 100644 --- a/deal.II/examples/step-52/step-52.cc +++ b/deal.II/examples/step-52/step-52.cc @@ -91,8 +91,9 @@ namespace Step52 // Evaluate $\left(I-\tau M^{-1} \frac{\partial f(t,y)}{\partial y}\right)^{-1}$ or // equivalently $\left(M-\tau \frac{\partial f}{\partial y}\right)^{-1} M$ at a given // time, for a given $\tau$ and y. - Vector id_minus_tau_J_inverse(const double time, const double tau, - const Vector &y); + Vector id_minus_tau_J_inverse(const double time, + const double tau, + const Vector &y); // Output the results as vtu files. void output_results(unsigned int time_step,TimeStepping::runge_kutta_method method) const; @@ -260,7 +261,7 @@ namespace Step52 // @sect5{Diffusion:evaluate_diffusion} // - // Evaluate the diffusion weak form give a time t and a vector y. + // Evaluate the weak form of the diffusion equation at a given time t and for a given vector y. Vector Diffusion::evaluate_diffusion(const double time, const Vector &y) const { Vector tmp(dof_handler.n_dofs()); @@ -327,13 +328,13 @@ namespace Step52 mass_minus_tau_Jacobian.copy_from(mass_matrix); mass_minus_tau_Jacobian.add(-tau,system_matrix); - // Inverse the matrix to get $\left(M-\tau \frac{\partial f}{\partial y}\right)^{-1}$ + // Inverse the matrix to get $\left(M-\tau \frac{\partial f}{\partial y}\right)^{-1}$. inverse_mass_minus_tau_Jacobian.initialize(mass_minus_tau_Jacobian); // Compute $tmp=My$. mass_matrix.vmult(tmp,y); - // Compute $\left(M-\tau \frac{\partial f}{\partial y}\right)^{-1} tmp$ + // Compute $\left(M-\tau \frac{\partial f}{\partial y}\right)^{-1} tmp$. inverse_mass_minus_tau_Jacobian.vmult(result,tmp); return result; @@ -424,7 +425,7 @@ namespace Step52 } - // sect5{Diffusion::explicit_method} + // @sect5{Diffusion::explicit_method} void Diffusion::explicit_method(TimeStepping::runge_kutta_method method, const unsigned int n_time_steps, const double initial_time, @@ -452,7 +453,7 @@ namespace Step52 - // sect5{Diffusion::implicit_method} + // @sect5{Diffusion::implicit_method} void Diffusion::implicit_method(TimeStepping::runge_kutta_method method, const unsigned int n_time_steps, const double initial_time, @@ -482,7 +483,7 @@ namespace Step52 - // sect5{Diffusion::embedded_explicit_method} + // @sect5{Diffusion::embedded_explicit_method} unsigned int Diffusion::embedded_explicit_method(TimeStepping::runge_kutta_method method, const unsigned int n_time_steps, const double initial_time, @@ -512,7 +513,7 @@ namespace Step52 unsigned int n_steps=0; while (timefinal_time) time_step = final_time-time; @@ -526,7 +527,7 @@ namespace Step52 if ((n_steps+1)%10==0) output_results(n_steps+1,method); - // Update the time step + // Update the time step. time_step = embedded_explicit_runge_kutta.get_status().delta_t_guess; ++n_steps; } @@ -536,7 +537,7 @@ namespace Step52 - // sect5{Diffusion::run} + // @sect5{Diffusion::run} void Diffusion::run() { // Create the grid (a square [0,5]x[0,5]) and refine the mesh four times. @@ -585,7 +586,7 @@ namespace Step52 // Use implicit midpoint. implicit_method(TimeStepping::IMPLICIT_MIDPOINT,n_time_steps,initial_time,final_time); std::cout<<"Implicit Midpoint error: "<