Skip to content

Commit 0752092

Browse files
authored
Update TGV example in documentation (#1363)
Description This PR updates the parameters of both the matrix-based and matrix-free TGV simulations according to the findings in the matrix-free article. To summarize: BDF3 scheme was changed to BDF2. Linear solver parameters were already tuned for the article, so we updated tolerances, number of iterations and krylov vectors accordingly. The example was written with constant time step, which works well so it won't be changed. However, trying constant CFL was added as a suggestion and the script to calculate the kinetic energy rate for that case was already available in the example folder. I do not think it is necessary to re run all the cases, so the part that talks about approximate times was left as it is.
1 parent a7b5f8d commit 0752092

File tree

3 files changed

+70
-92
lines changed

3 files changed

+70
-92
lines changed

doc/source/examples/incompressible-flow/3d-taylor-green-vortex/3d-taylor-green-vortex.rst

+32-44
Original file line numberDiff line numberDiff line change
@@ -10,7 +10,7 @@ Features
1010
---------
1111

1212
- Solvers: ``lethe-fluid`` (with Q2-Q2) or ``lethe-fluid-matrix-free`` (with Q2-Q2 or Q3-Q3)
13-
- Transient problem using ``bdf3`` time integrator
13+
- Transient problem using ``bdf2`` time integrator
1414
- Displays the calculation of enstrophy and total kinetic energy
1515

1616

@@ -34,13 +34,12 @@ The three velocity components :math:`[u_x,u_y,u_z]^T` and the pressure :math:`p`
3434

3535
.. math::
3636
37-
u_{x} &= \sin(x)*\cos(y)*\cos(z) \\
38-
u_{y} &= -\cos(x)*\sin(y)*\cos(z)\\
37+
u_{x} &= \sin(x)\cos(y)\cos(z) \\
38+
u_{y} &= -\cos(x)\sin(y)\cos(z)\\
3939
u_{z} &= 0 \\
40-
p &= \frac{1}{16}*\left[\cos(2x)+\cos(2y)\right]\left[\cos(2z)+2\right]
41-
42-
In this case, the vortex, which is initially 2D, will decay by generating smaller 3D turbulent structures (vortex tubes, rings and sheets). This decay can be monitored through the total kinetic energy of the system. Since the simulation domain is periodic, it can be demontrated that the time derative of the total kinetic energy :math:`E_\mathrm{k}` is directly related to the enstrophy :math:`\mathcal{E}` such that:
40+
p &= \frac{1}{16}\left[\cos(2x)+\cos(2y)\right]\left[\cos(2z)+2\right]
4341
42+
In this case, the vortex, which is initially 2D, will decay by generating smaller 3D turbulent structures (vortex tubes, rings and sheets). This decay can be monitored through the total kinetic energy of the system. Since the simulation domain is periodic, it can be demonstrated that the time derative of the total kinetic energy :math:`E_\mathrm{k}` is directly related to the enstrophy :math:`\mathcal{E}` such that:
4443

4544

4645
.. math::
@@ -70,7 +69,7 @@ The ``mesh`` subsection specifies the computational grid:
7069
set initial refinement = 5
7170
end
7271
73-
The ``type`` specifies the mesh format used. We use the ``hyper_cube`` mesh generated from the deal.II `GridGenerator <https://www.dealii.org/current/doxygen/deal.II/namespaceGridGenerator.html>`_ . We set ``colorize = true`` to be able to adequately set-up the periodic boundary conditions.
72+
The ``type`` specifies the mesh format used. We use the ``hyper_cube`` mesh generated from the deal.II `GridGenerator <https://www.dealii.org/current/doxygen/deal.II/namespaceGridGenerator.html>`_ . We set ``colorize = true`` (the last parameter of the grid arguments) to assign boundary condition ids to each of the walls of the cube.
7473

7574

7675
The last parameter specifies the ``initial refinement`` of the grid. Indicating an ``initial refinement = 5`` implies that the initial mesh is refined 5 times. In 3D, each cell is divided by 8 per refinement. Consequently, the final grid is made of 32768 cells.
@@ -104,7 +103,7 @@ The ``boundary conditions`` subsection establishes the constraints on different
104103
end
105104
end
106105
107-
First, the ``number`` of boundary conditions to be applied must be specified. For each boundary condition, the ``id`` of the boundary as well as its ``type`` must be specified. All boundaries are ``periodic``. The ``x-`` (id=0) is periodic with the ``x+`` boundary (id=1), the ``y-`` (id=2) is periodic with the ``y+`` boundary (id=3) and so on and so forth. For each periodic boundary condition, the periodic direction must be specified. A periodic direction of ``0`` implies that the normal direction of the wall is the :math:`\mathbf{e}_x` vector, ``1`` implies that it's the :math:`\mathbf{e}_y`.
106+
First, the ``number`` of boundary conditions to be applied must be specified. For each boundary condition, the ``id`` of the boundary as well as its ``type`` must be specified. All boundaries are ``periodic``. The ``x-`` boundary (id=0) is periodic with the ``x+`` boundary (id=1), the ``y-`` boundary (id=2) is periodic with the ``y+`` boundary (id=3) and so on and so forth. For each periodic boundary condition, the periodic direction must be specified. A periodic direction of ``0`` implies that the normal direction of the wall is the :math:`\mathbf{e}_x` vector, ``1`` implies that it's the :math:`\mathbf{e}_y`.
108107

109108
Physical Properties
110109
~~~~~~~~~~~~~~~~~~~
@@ -124,7 +123,7 @@ The Reynolds number of 1600 is set solely using the kinematic viscosity since th
124123
FEM Interpolation
125124
~~~~~~~~~~~~~~~~~
126125

127-
The results obtained for the Taylor-Green vortex are highly dependent on the numerical dissipation that occurs within the CFD scheme. Generally, high-order methods outperform traditional second-order accurate methods for this type of flow. In the present case, we will investigate the usage of both second and third degree polynomial.
126+
The results obtained for the Taylor-Green vortex are highly dependent on the numerical dissipation that occurs within the CFD scheme. Generally, high-order methods outperform traditional second-order accurate methods for this type of flow. In the present case, we will investigate the usage of both second and third degree polynomials.
128127

129128
.. code-block:: text
130129
@@ -149,21 +148,19 @@ To monitor the kinetic energy and the enstrophy, we set both calculation to ``tr
149148
Simulation Control
150149
~~~~~~~~~~~~~~~~~~
151150

152-
The ``simulation control`` subsection controls the flow of the simulation. To maximize the temporal accuracy of the simulation, we use a third order ``bdf3`` scheme. Results are written every 2 time-steps. To ensure a more adequate visualization of the high-order elements, we set ``subdivision = 3``. This will allow Paraview to render the high-order solutions with more fidelity.
151+
The ``simulation control`` subsection controls the flow of the simulation. To maximize the temporal accuracy of the simulation, we use a second order ``bdf2`` scheme. Results are written every 2 time-steps. To ensure a more adequate visualization of the high-order elements, we set ``subdivision = 3``. This will allow Paraview to render the high-order solutions with more fidelity.
153152

154153
.. code-block:: text
155154
156155
subsection simulation control
157-
set method = bdf3
156+
set method = bdf2
158157
set time step = 0.05
159-
set number mesh adapt = 0
160158
set time end = 20
161159
set output frequency = 2
162160
set subdivision = 3
163161
end
164162
165163
166-
167164
Matrix-based - Non-linear Solver
168165
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
169166

@@ -192,12 +189,13 @@ Since this is a transient problem, the linear solver can be relatively simple. W
192189
subsection fluid dynamics
193190
set verbosity = verbose
194191
set method = gmres
195-
set max iters = 200
192+
set max iters = 100
196193
set max krylov vectors = 200
197194
set relative residual = 1e-4
198-
set minimum residual = 1e-12
195+
set minimum residual = 1e-7
196+
set preconditioner = ilu
199197
set ilu preconditioner fill = 0
200-
set ilu preconditioner absolute tolerance = 1e-12
198+
set ilu preconditioner absolute tolerance = 1e-10
201199
set ilu preconditioner relative tolerance = 1.00
202200
end
203201
end
@@ -225,45 +223,33 @@ The ``lethe-fluid-matrix-free`` has significantly more parameters for its linear
225223
226224
subsection linear solver
227225
subsection fluid dynamics
228-
set method = gmres
229-
set max iters = 100
230-
set relative residual = 1e-4
231-
set minimum residual = 1e-7
232-
set preconditioner = gcmg
233-
set verbosity = verbos
226+
set verbosity = verbose
227+
set method = gmres
228+
set max iters = 100
229+
set max krylov vectors = 200
230+
set relative residual = 1e-4
231+
set minimum residual = 1e-7
232+
set preconditioner = gcmg
234233
235234
# MG parameters
236-
set mg verbosity = quiet
237-
set mg min level = -1
238-
set mg level min cells = 16
235+
set mg verbosity = quiet
236+
set mg enable hessians in jacobian = false
239237
240238
# Smoother
241-
set mg smoother iterations = 10
239+
set mg smoother iterations = 5
242240
set mg smoother eig estimation = true
243-
241+
244242
# Eigenvalue estimation parameters
245243
set eig estimation smoothing range = 5
246244
set eig estimation cg n iterations = 20
247-
set eig estimation verbosity = verbose
245+
set eig estimation verbosity = quiet
248246
249247
# Coarse-grid solver
250-
set mg coarse grid solver = gmres
251-
set mg gmres max iterations = 2000
252-
set mg gmres tolerance = 1e-7
253-
set mg gmres reduce = 1e-4
254-
set mg gmres max krylov vectors = 30
255-
set mg gmres preconditioner = ilu
256-
set ilu preconditioner fill = 1
257-
set ilu preconditioner absolute tolerance = 1e-10
258-
set ilu preconditioner relative tolerance = 1.00
248+
set mg coarse grid solver = direct
259249
end
260250
end
261251
262-
We set ``mg verbosity = quiet`` to prevent logging of the multigrid parameters during the simulation. Setting ``mg min level = -1`` ensures that the ``mg level min cells = 16`` parameter is used to determine the coarsest level. It is important to ensure that the Taylor-Green vortex has sufficient cells on the coarsest level since periodic boundary conditions are used. Indeed, using a coarsest level with a single cell can lead to a problematic situation where too few degrees of freedom are available on the coarsest level.
263-
264-
The ``smoother``, ``Eigenvalue estimation parameters`` and ``coarse-grid solver`` subsections are explained in the **Theory Guide** (under construction).
265-
266-
252+
We set ``mg verbosity = quiet`` to prevent logging of the multigrid parameters during the simulation. The ``smoother``, ``Eigenvalue estimation parameters`` and ``coarse-grid solver`` subsections are explained in the :doc:`../../../parameters/cfd/linear_solver_control` section.
267253

268254
----------------------
269255
Running the Simulation
@@ -304,14 +290,14 @@ Using the ``enstrophy.dat`` and ``kinetic_energy.dat`` files generated by Lethe,
304290
.. code-block:: text
305291
:class: copy-button
306292
307-
python3 calculate_dissipation_rate.py -i kinetic_energy.dat -o output.dat
293+
python3 calculate_dissipation_rate.py -i output/kinetic_energy.dat
308294
309295
Then, by invoking the second script present in the example, a plot comparing the kinetic energy decay with the enstrophy is generated:
310296

311297
.. code-block:: text
312298
:class: copy-button
313299
314-
python3 plot_dissipation_rate.py -ke kinetic_energy_decay.dat -ens enstrophy.dat -v 0.000625
300+
python3 plot_dissipation_rate.py -ke ke_rate.dat -ens output/enstrophy.dat -v 0.000625
315301
316302
.. tip::
317303

@@ -356,6 +342,8 @@ Possibilities for Extension
356342

357343
- This case is very interesting to postprocess. Try to postprocess this case using other quantities (vorticity, q-criterion) and use the results to generate interesting animations. Feel free to share them with us!
358344

345+
- This case can also be used to experiment with adaptive time step. In the simulation control section add ``adapt = true`` and ``set max cfl = 1``, similar results should be obtained but with significantly less iterations as larger time steps are taken. To postprocess the results use the additional script ``calculate_dissipation_rate_constant_cfl.py`` given in the same folder to calculate the kinetic energy rate.
346+
359347

360348
------------
361349
References

examples/incompressible-flow/3d-taylor-green-vortex/tgv-matrix-based.prm

+13-13
Original file line numberDiff line numberDiff line change
@@ -11,14 +11,13 @@ set dimension = 3
1111
#---------------------------------------------------
1212

1313
subsection simulation control
14-
set method = bdf3
15-
set time step = 0.05
16-
set number mesh adapt = 0
17-
set time end = 20
18-
set output name = tgv
19-
set output frequency = 2
20-
set output path = ./output/
21-
set subdivision = 3
14+
set method = bdf2
15+
set time step = 0.05
16+
set time end = 20
17+
set output name = tgv
18+
set output frequency = 2
19+
set output path = ./output/
20+
set subdivision = 3
2221
end
2322

2423
#---------------------------------------------------
@@ -35,7 +34,7 @@ end
3534
#---------------------------------------------------
3635

3736
subsection timer
38-
set type = iteration
37+
set type = iteration
3938
end
4039

4140
#---------------------------------------------------
@@ -86,7 +85,7 @@ subsection mesh
8685
set type = dealii
8786
set grid type = hyper_cube
8887
set grid arguments = -3.14159265359 : 3.14159265359 : true
89-
set initial refinement = 5
88+
set initial refinement = 5
9089
end
9190

9291
# --------------------------------------------------
@@ -137,12 +136,13 @@ subsection linear solver
137136
subsection fluid dynamics
138137
set verbosity = verbose
139138
set method = gmres
140-
set max iters = 200
139+
set max iters = 100
141140
set max krylov vectors = 200
142141
set relative residual = 1e-4
143-
set minimum residual = 1e-12
142+
set minimum residual = 1e-7
143+
set preconditioner = ilu
144144
set ilu preconditioner fill = 0
145-
set ilu preconditioner absolute tolerance = 1e-12
145+
set ilu preconditioner absolute tolerance = 1e-10
146146
set ilu preconditioner relative tolerance = 1.00
147147
end
148148
end

examples/incompressible-flow/3d-taylor-green-vortex/tgv-matrix-free.prm

+25-35
Original file line numberDiff line numberDiff line change
@@ -11,14 +11,13 @@ set dimension = 3
1111
#---------------------------------------------------
1212

1313
subsection simulation control
14-
set method = bdf3
15-
set time step = 0.05
16-
set number mesh adapt = 0
17-
set time end = 20
18-
set output name = tgv
19-
set output frequency = 2
20-
set output path = ./output/
21-
set subdivision = 3
14+
set method = bdf2
15+
set time step = 0.05
16+
set time end = 20
17+
set output name = tgv
18+
set output frequency = 2
19+
set output path = ./output/
20+
set subdivision = 3
2221
end
2322

2423
#---------------------------------------------------
@@ -35,7 +34,7 @@ end
3534
#---------------------------------------------------
3635

3736
subsection timer
38-
set type = iteration
37+
set type = iteration
3938
end
4039

4140
#---------------------------------------------------
@@ -66,7 +65,7 @@ end
6665

6766
subsection source term
6867
subsection fluid dynamics
69-
set enable = false
68+
set enable = false
7069
end
7170
end
7271

@@ -88,7 +87,7 @@ subsection mesh
8887
set type = dealii
8988
set grid type = hyper_cube
9089
set grid arguments = -3.14159265359 : 3.14159265359 : true
91-
set initial refinement = 5
90+
set initial refinement = 5
9291
end
9392

9493
#---------------------------------------------------
@@ -123,8 +122,8 @@ end
123122

124123
subsection non-linear solver
125124
subsection fluid dynamics
126-
set tolerance = 1e-3
127-
set verbosity = verbose
125+
set tolerance = 1e-3
126+
set verbosity = verbose
128127
end
129128
end
130129

@@ -134,37 +133,28 @@ end
134133

135134
subsection linear solver
136135
subsection fluid dynamics
137-
set method = gmres
138-
set max iters = 100
139-
set relative residual = 1e-4
140-
set minimum residual = 1e-7
141-
set preconditioner = gcmg
142-
set verbosity = verbose
136+
set verbosity = verbose
137+
set method = gmres
138+
set max iters = 100
139+
set max krylov vectors = 200
140+
set relative residual = 1e-4
141+
set minimum residual = 1e-7
142+
set preconditioner = gcmg
143143

144144
# MG parameters
145-
set mg verbosity = verbose
146-
set mg min level = -1
147-
set mg level min cells = 16
145+
set mg verbosity = quiet
146+
set mg enable hessians in jacobian = false
148147

149148
# Smoother
150-
set mg smoother iterations = 10
149+
set mg smoother iterations = 5
151150
set mg smoother eig estimation = true
152151

153152
# Eigenvalue estimation parameters
154153
set eig estimation smoothing range = 5
155-
set eig estimation cg n iterations = 10
156-
set eig estimation verbosity = verbose
154+
set eig estimation cg n iterations = 20
155+
set eig estimation verbosity = quiet
157156

158157
# Coarse-grid solver
159-
set mg coarse grid solver = gmres
160-
set mg gmres max iterations = 2000
161-
set mg gmres tolerance = 1e-7
162-
set mg gmres reduce = 1e-4
163-
set mg gmres max krylov vectors = 30
164-
set mg gmres preconditioner = ilu
165-
166-
set ilu preconditioner fill = 1
167-
set ilu preconditioner absolute tolerance = 1e-10
168-
set ilu preconditioner relative tolerance = 1.00
158+
set mg coarse grid solver = direct
169159
end
170160
end

0 commit comments

Comments
 (0)