Skip to content

Commit 457c221

Browse files
oguevremontblaisbmivaiahepap
authored
Added homogeneous and heterogeneous first order reactions for tracer (chaos-polymtl#1411)
Description The tracer physics needed one new assembler and two new reaction constant models to take into account simple reactions where only one reactive species is considered. This PR adds this feature, with documentation and application tests of homogeneous and heterogeneous first order reactions. Testing 3 new tests are added, that test the features both with and without solid, and with the discontinuous Galerkin formulation as well. [lethe-fluid/tracer_first_order_reaction_homogeneous] [lethe-fluid/tracer_first_order_reaction_homogeneous_dg] [lethe-fluid-sharp/tracer_first_order_reaction_heterogeneous] Documentation The new parameters are documented. Co-authored-by: Bruno Blais <[email protected]> Co-authored-by: mivaia <[email protected]> Co-authored-by: hepap <[email protected]>
1 parent 8edc05b commit 457c221

29 files changed

+1402
-79
lines changed

CHANGELOG.md

Lines changed: 6 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -3,6 +3,12 @@
33
All notable changes to the Lethe project will be documented in this file.
44
The format is based on [Keep a Changelog](http://keepachangelog.com/).
55

6+
## [Master] - 2025-02-22
7+
8+
### Added
9+
10+
- MINOR In the tracer physics, a new feature is added that allows to take into account single species reaction of arbitrary reaction constant and order. This feature works in CG or DG, with or without immersed solids. [#1411](https://github.com/chaos-polymtl/lethe/pull/1411)
11+
612
## [Master] - 2025-02-20
713

814
### Added
Lines changed: 45 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,45 @@
1+
Running on 1 MPI rank(s)...
2+
Number of active cells: 1024
3+
Number of degrees of freedom: 3267
4+
Volume of triangulation: 1
5+
Number of tracer degrees of freedom: 1089
6+
7+
*****************************
8+
Steady iteration: 1/3
9+
*****************************
10+
div u : 3.22831e-09
11+
error at the weak BC : 0.000732218
12+
L2 error velocity : 0.00107872 L2 error pressure: 0.0119908
13+
L2 error tracer: 2.46686e-06
14+
15+
*****************************
16+
Steady iteration: 2/3
17+
*****************************
18+
Number of active cells: 4096
19+
Number of degrees of freedom: 12675
20+
Volume of triangulation: 1
21+
Number of tracer degrees of freedom: 4225
22+
div u : -3.03625e-08
23+
error at the weak BC : 0.000187088
24+
L2 error velocity : 0.000275274 L2 error pressure: 0.00423168
25+
L2 error tracer: 4.47348e-07
26+
27+
*****************************
28+
Steady iteration: 3/3
29+
*****************************
30+
Number of active cells: 16384
31+
Number of degrees of freedom: 49923
32+
Volume of triangulation: 1
33+
Number of tracer degrees of freedom: 16641
34+
div u : 2.93362e-10
35+
error at the weak BC : 4.74152e-05
36+
L2 error velocity : 7.0998e-05 L2 error pressure: 0.00140056
37+
L2 error tracer: 9.70565e-08
38+
cells error_velocity error_pressure
39+
1024 1.078724e-03 - 1.199081e-02 -
40+
4096 2.752737e-04 1.97 4.231680e-03 1.50
41+
16384 7.099801e-05 1.96 1.400565e-03 1.60
42+
cells error_tracer
43+
1024 2.466861e-06 -
44+
4096 4.473484e-07 2.46
45+
16384 9.705649e-08 2.20
Lines changed: 173 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,173 @@
1+
# SPDX-FileCopyrightText: Copyright (c) 2025 The Lethe Authors
2+
# SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception OR LGPL-2.1-or-later
3+
# This MMS tests a 1.5-order heterogeneous reaction across a disk.
4+
5+
set dimension = 2
6+
7+
subsection simulation control
8+
set method = steady
9+
set output frequency = 0
10+
set number mesh adapt = 2
11+
end
12+
13+
subsection multiphysics
14+
set fluid dynamics = true
15+
set tracer = true
16+
end
17+
18+
subsection FEM
19+
set velocity order = 1
20+
set pressure order = 1
21+
set tracer order = 1
22+
end
23+
24+
subsection mesh
25+
set type = dealii
26+
set grid type = subdivided_hyper_rectangle
27+
set grid arguments = 1, 1: 0, 0: 1, 1 : true
28+
set initial refinement = 5
29+
end
30+
31+
subsection mesh adaptation
32+
set type = uniform
33+
end
34+
35+
subsection boundary conditions
36+
set number = 4
37+
subsection bc 0
38+
set id = 0
39+
set type = function weak
40+
end
41+
subsection bc 1
42+
set id = 1
43+
set type = function weak
44+
end
45+
subsection bc 2
46+
set id = 2
47+
set type = function weak
48+
end
49+
subsection bc 3
50+
set id = 3
51+
set type = function weak
52+
end
53+
end
54+
55+
subsection boundary conditions tracer
56+
set number = 4
57+
set time dependent = false
58+
subsection bc 0
59+
set id = 0
60+
set type = dirichlet
61+
subsection dirichlet
62+
set Function expression = (sin(pi*x)*cos(pi*y) + 1)*(0.5*tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0) + 0.5)/2
63+
end
64+
end
65+
subsection bc 1
66+
set id = 1
67+
set type = dirichlet
68+
subsection dirichlet
69+
set Function expression = (sin(pi*x)*cos(pi*y) + 1)*(0.5*tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0) + 0.5)/2
70+
end
71+
end
72+
subsection bc 2
73+
set id = 2
74+
set type = dirichlet
75+
subsection dirichlet
76+
set Function expression = (sin(pi*x)*cos(pi*y) + 1)*(0.5*tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0) + 0.5)/2
77+
end
78+
end
79+
subsection bc 3
80+
set id = 3
81+
set type = dirichlet
82+
subsection dirichlet
83+
set Function expression = (sin(pi*x)*cos(pi*y) + 1)*(0.5*tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0) + 0.5)/2
84+
end
85+
end
86+
end
87+
88+
subsection initial conditions
89+
set type = nodal
90+
subsection uvwp
91+
set Function expression = (((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*x)^2*sin(pi*y)*cos(pi*y); -(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*x)*sin(pi*y)^2*cos(pi*x); sin(pi*x) + sin(pi*y)
92+
end
93+
subsection tracer
94+
set Function expression = (sin(pi*x)*cos(pi*y) + 1)*(0.5*tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0) + 0.5)/2
95+
end
96+
end
97+
98+
subsection particles
99+
set assemble Navier-Stokes inside particles = true
100+
set number of particles = 1
101+
subsection particle info 0
102+
set type = sphere
103+
set shape arguments = 0.1
104+
subsection position
105+
set Function expression = 0.5; 0.5
106+
end
107+
end
108+
end
109+
110+
subsection physical properties
111+
set number of fluids = 1
112+
subsection fluid 0
113+
set tracer diffusivity model = immersed solid tanh
114+
set tracer reaction constant model = immersed solid tanh
115+
set tracer reaction order = 1.5
116+
subsection immersed solid tanh
117+
set tracer diffusivity inside = 1
118+
set tracer diffusivity outside = 1
119+
set tracer reaction constant inside = 9.86960440109
120+
set tracer reaction constant outside = 0
121+
set thickness = 0.1
122+
end
123+
end
124+
end
125+
126+
subsection analytical solution
127+
set enable = true
128+
set verbosity = verbose
129+
subsection uvwp
130+
set Function expression = (((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*x)^2*sin(pi*y)*cos(pi*y); -(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*x)*sin(pi*y)^2*cos(pi*x); sin(pi*x) + sin(pi*y)
131+
end
132+
subsection tracer
133+
set Function expression = (sin(pi*x)*cos(pi*y) + 1)*(0.5*tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0) + 0.5)/2
134+
end
135+
end
136+
137+
subsection source term
138+
subsection tracer
139+
set Function expression = 0.353553390593274*((sin(pi*x)*cos(pi*y) + 1)*(0.5*tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0) + 0.5))^1.5*(-pi^2*(0.5*tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0) + 0.5) + pi^2) + (2.5*(1 - tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0)^2)*(1.0*x - 0.5)*(sin(pi*x)*cos(pi*y) + 1)/((x - 0.5)^2 + (y - 0.5)^2)^0.5 + pi*(0.5*tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0) + 0.5)*cos(pi*x)*cos(pi*y)/2)*(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*x)^2*sin(pi*y)*cos(pi*y) - (2.5*(1 - tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0)^2)*(1.0*y - 0.5)*(sin(pi*x)*cos(pi*y) + 1)/((x - 0.5)^2 + (y - 0.5)^2)^0.5 - pi*(0.5*tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0) + 0.5)*sin(pi*x)*sin(pi*y)/2)*(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*x)*sin(pi*y)^2*cos(pi*x) - (-10.0*pi*(1.0*x - 0.5)*(tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0)^2 - 1)*cos(pi*x)*cos(pi*y)/((x - 0.5)^2 + (y - 0.5)^2)^0.5 + (sin(pi*x)*cos(pi*y) + 1)*(tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0)^2 - 1)*(5.0*(1.0*x - 0.5)^2/((x - 0.5)^2 + (y - 0.5)^2)^1.5 + 100.0*(1.0*x - 0.5)^2*tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0)/((x - 0.5)^2 + (y - 0.5)^2)^1.0 - 5.0/((x - 0.5)^2 + (y - 0.5)^2)^0.5) - 0.5*pi^2*(tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0) + 1)*sin(pi*x)*cos(pi*y))/2 - (10.0*pi*(1.0*y - 0.5)*(tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0)^2 - 1)*sin(pi*x)*sin(pi*y)/((x - 0.5)^2 + (y - 0.5)^2)^0.5 + (sin(pi*x)*cos(pi*y) + 1)*(tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0)^2 - 1)*(5.0*(1.0*y - 0.5)^2/((x - 0.5)^2 + (y - 0.5)^2)^1.5 + 100.0*(1.0*y - 0.5)^2*tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0)/((x - 0.5)^2 + (y - 0.5)^2)^1.0 - 5.0/((x - 0.5)^2 + (y - 0.5)^2)^0.5) - 0.5*pi^2*(tanh(10.0*((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 1.0) + 1)*sin(pi*x)*cos(pi*y))/2
140+
end
141+
subsection fluid dynamics
142+
set Function expression = ((1.0*x - 0.5)*sin(pi*x)^2*sin(pi*y)*cos(pi*y)/((x - 0.5)^2 + (y - 0.5)^2)^0.5 + 2*pi*(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*x)*sin(pi*y)*cos(pi*x)*cos(pi*y))*(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*x)^2*sin(pi*y)*cos(pi*y) - (((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*((1.0*y - 0.5)*sin(pi*x)^2*sin(pi*y)*cos(pi*y)/((x - 0.5)^2 + (y - 0.5)^2)^0.5 - pi*(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*x)^2*sin(pi*y)^2 + pi*(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*x)^2*cos(pi*y)^2)*sin(pi*x)*sin(pi*y)^2*cos(pi*x) - (4*pi*(1.0*x - 0.5)*sin(pi*x)*cos(pi*x)/((x - 0.5)^2 + (y - 0.5)^2)^0.5 - ((1.0*x - 0.5)^2/((x - 0.5)^2 + (y - 0.5)^2)^1.5 - 1.0/((x - 0.5)^2 + (y - 0.5)^2)^0.5)*sin(pi*x)^2 - 2*pi^2*(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*(sin(pi*x)^2 - cos(pi*x)^2))*sin(pi*y)*cos(pi*y) - (-2*pi*(1.0*y - 0.5)*sin(pi*y)^2/((x - 0.5)^2 + (y - 0.5)^2)^0.5 + 2*pi*(1.0*y - 0.5)*cos(pi*y)^2/((x - 0.5)^2 + (y - 0.5)^2)^0.5 - ((1.0*y - 0.5)^2/((x - 0.5)^2 + (y - 0.5)^2)^1.5 - 1.0/((x - 0.5)^2 + (y - 0.5)^2)^0.5)*sin(pi*y)*cos(pi*y) - 4*pi^2*(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*y)*cos(pi*y))*sin(pi*x)^2 + pi*cos(pi*x); -(-(1.0*y - 0.5)*sin(pi*x)*sin(pi*y)^2*cos(pi*x)/((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 2*pi*(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*x)*sin(pi*y)*cos(pi*x)*cos(pi*y))*(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*x)*sin(pi*y)^2*cos(pi*x) + (((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*(-(1.0*x - 0.5)*sin(pi*x)*sin(pi*y)^2*cos(pi*x)/((x - 0.5)^2 + (y - 0.5)^2)^0.5 + pi*(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*x)^2*sin(pi*y)^2 - pi*(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*y)^2*cos(pi*x)^2)*sin(pi*x)^2*sin(pi*y)*cos(pi*y) - (-4*pi*(1.0*y - 0.5)*sin(pi*y)*cos(pi*y)/((x - 0.5)^2 + (y - 0.5)^2)^0.5 + ((1.0*y - 0.5)^2/((x - 0.5)^2 + (y - 0.5)^2)^1.5 - 1.0/((x - 0.5)^2 + (y - 0.5)^2)^0.5)*sin(pi*y)^2 + 2*pi^2*(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*(sin(pi*y)^2 - cos(pi*y)^2))*sin(pi*x)*cos(pi*x) - (2*pi*(1.0*x - 0.5)*sin(pi*x)^2/((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 2*pi*(1.0*x - 0.5)*cos(pi*x)^2/((x - 0.5)^2 + (y - 0.5)^2)^0.5 + ((1.0*x - 0.5)^2/((x - 0.5)^2 + (y - 0.5)^2)^1.5 - 1.0/((x - 0.5)^2 + (y - 0.5)^2)^0.5)*sin(pi*x)*cos(pi*x) + 4*pi^2*(((x - 0.5)^2 + (y - 0.5)^2)^0.5 - 0.1)*sin(pi*x)*cos(pi*x))*sin(pi*y)^2 + pi*cos(pi*y); (1.0*x - 0.5)*sin(pi*x)^2*sin(pi*y)*cos(pi*y)/((x - 0.5)^2 + (y - 0.5)^2)^0.5 - (1.0*y - 0.5)*sin(pi*x)*sin(pi*y)^2*cos(pi*x)/((x - 0.5)^2 + (y - 0.5)^2)^0.5
143+
end
144+
end
145+
146+
subsection non-linear solver
147+
subsection fluid dynamics
148+
set verbosity = quiet
149+
set tolerance = 1e-6
150+
set max iterations = 10
151+
end
152+
subsection tracer
153+
set verbosity = quiet
154+
set tolerance = 1e-6
155+
set max iterations = 10
156+
end
157+
end
158+
159+
subsection linear solver
160+
subsection fluid dynamics
161+
set verbosity = quiet
162+
set method = gmres
163+
set relative residual = 1e-3
164+
set minimum residual = 1e-10
165+
set max krylov vectors = 200
166+
end
167+
subsection tracer
168+
set verbosity = quiet
169+
set method = gmres
170+
set relative residual = 1e-6
171+
set minimum residual = 1e-10
172+
end
173+
end
Lines changed: 39 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,39 @@
1+
Running on 1 MPI rank(s)...
2+
Number of active cells: 256
3+
Number of degrees of freedom: 867
4+
Volume of triangulation: 1
5+
Number of tracer degrees of freedom: 289
6+
7+
*****************************
8+
Steady iteration: 1/3
9+
*****************************
10+
L2 error velocity: 0.00380898
11+
L2 error tracer: 1.57375e-06
12+
13+
*****************************
14+
Steady iteration: 2/3
15+
*****************************
16+
Number of active cells: 1024
17+
Number of degrees of freedom: 3267
18+
Volume of triangulation: 1
19+
Number of tracer degrees of freedom: 1089
20+
L2 error velocity: 0.000962316
21+
L2 error tracer: 9.81605e-08
22+
23+
*****************************
24+
Steady iteration: 3/3
25+
*****************************
26+
Number of active cells: 4096
27+
Number of degrees of freedom: 12675
28+
Volume of triangulation: 1
29+
Number of tracer degrees of freedom: 4225
30+
L2 error velocity: 0.000241505
31+
L2 error tracer: 6.13192e-09
32+
cells error_velocity error_pressure
33+
256 3.808975e-03 - 1.772614e-02 -
34+
1024 9.623156e-04 1.98 5.122861e-03 1.79
35+
4096 2.415048e-04 1.99 1.490588e-03 1.78
36+
cells error_tracer
37+
256 1.573750e-06 -
38+
1024 9.816048e-08 4.00
39+
4096 6.131918e-09 4.00

0 commit comments

Comments
 (0)