From 1a078acd0f73013944e2d4ccf2317ef65f851888 Mon Sep 17 00:00:00 2001 From: bangerth Date: Fri, 22 Aug 2008 20:06:35 +0000 Subject: [PATCH] Make beta and time step dimension independent. git-svn-id: https://svn.dealii.org/trunk@16657 0785d39b-7218-0410-832d-ea1e28bc413d --- deal.II/examples/step-31/doc/results.dox | 52 +++++++++++++++++++++++- deal.II/examples/step-31/step-31.cc | 4 +- 2 files changed, 52 insertions(+), 4 deletions(-) diff --git a/deal.II/examples/step-31/doc/results.dox b/deal.II/examples/step-31/doc/results.dox index 6083f95148..6a2f055b1a 100644 --- a/deal.II/examples/step-31/doc/results.dox +++ b/deal.II/examples/step-31/doc/results.dox @@ -196,14 +196,62 @@ is no longer than one grid point away (which for $Q_1$ elements is $h_K$, but for $Q_2$ elements is $h_K/2$). +
Results for 3d
+ +One can repeat these experiments in 3d and find the optimal time step +for each value of $\beta$ and find the best value of $\beta$. What one +finds is that for the same $\beta$ already used in 2d, the time steps +needs to be a bit small, by around a factor of 1.2 or so. This is +easily explained: the time step restriction is +$k=\min_K \frac{ch_K}{\|\mathbf{u}\|_{L^\infty(K)}}$ where $h_K$ is +the diameter of the cell. However, what is really needed is the +distance between mesh points, which is $\frac{h_K}{\sqrt{d}}$. So a +more appropriate form would be +$k=\min_K \frac{ch_K}{\|\mathbf{u}\|_{L^\infty(K)}\sqrt{d}}$. + +The second find is that one needs to choose $\beta$ slightly bigger +(about $\beta=0.05$ or so). This then again reduces the time step we +can take. + + + +
Conclusions
Concluding, $\beta=0.03$ appears to be a good choice for the -stabilization parameter, and then +stabilization parameter in 2d, and $\beta=0.05$ in 3d. In a dimension +independent way, we can model this as $\beta=0.015d$. As we have seen +in the sections above, in 2d $k=\frac 14 \frac 1{q_T}\frac{h_K}{\|\mathbf{u}\|_{L^\infty(K)}}$ is an appropriate time step, where $q_T$ is the polynomial degree of the temperature shape functions (in the program, this corresponds to -the variable temperature_degree). +the variable temperature_degree). To reconcile this with +the findings in 3d for the same $\beta$, we could write this as +$k=\frac 1{2\sqrt{2}\sqrt{d}} \frac +1{q_T}\frac{h_K}{\|\mathbf{u}\|_{L^\infty(K)}}$ +but this doesn't take into account that we also have to increase +$\beta$ in 3d. The final form that takes all these factors in reads as +follows: +@f{eqnarray*} + k = + \frac 1{2\sqrt{2}} \frac 1{\sqrt{d}} + \frac 2d + \frac 1{q_T} + \frac{h_K}{\|\mathbf{u}\|_{L^\infty(K)}} + = + \frac 1{d\sqrt{2}\sqrt{d}} + \frac 1{q_T} + \frac{h_K}{\|\mathbf{u}\|_{L^\infty(K)}}. +@f} +In the first form (in the center of the equation), $\frac +1{2\sqrt{2}}$ is a universal constant, $\frac 1{\sqrt{d}}$ +is the factor that accounts for the difference between cell diameter +and grid point separation, +$\frac 2d$ accounts for the increase in $\beta$ with space dimension, +$\frac 1{q_T}$ accounts for the distance between grid points for +higher order elements, and $\frac{h_K}{\|\mathbf{u}\|_{L^\infty(K)}}$ +for the local speed of transport relative to the cell size. This is +the formula that we use in the program. As for the question of whether to use $Q_1$ or $Q_2$ elements for the temperature, the following considerations may be useful: First, diff --git a/deal.II/examples/step-31/step-31.cc b/deal.II/examples/step-31/step-31.cc index d513c327bb..5ced6fe730 100644 --- a/deal.II/examples/step-31/step-31.cc +++ b/deal.II/examples/step-31/step-31.cc @@ -703,7 +703,7 @@ compute_viscosity(const std::vector &old_temperature, const double cell_diameter, const double old_time_step) { - const double beta = 0.03; + const double beta = 0.015 * dim; const double alpha = 1; if (global_u_infty == 0) @@ -1908,7 +1908,7 @@ void BoussinesqFlowProblem::solve () } old_time_step = time_step; - time_step = 1./4. / + time_step = 1./(std::sqrt(2.)*dim*std::sqrt(1.*dim)) / temperature_degree * GridTools::minimal_cell_diameter(triangulation) / std::max (get_maximal_velocity(), .01); -- 2.39.5