Skip to content
GitLab
Projects
Groups
Snippets
/
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Menu
Open sidebar
itischler
waLBerla
Commits
7751b2df
Commit
7751b2df
authored
Apr 12, 2021
by
Helen Schottenhamml
Browse files
Merge branch 'Fix_SP_conversion' into 'master'
FIX: SP/DP conversion warnings Closes #144 See merge request
walberla/walberla!442
parents
8f1f6b17
05b33755
Changes
12
Hide whitespace changes
Inline
Side-by-side
apps/benchmarks/FluidParticleCoupling/GeneratedLBM.py
View file @
7751b2df
...
...
@@ -14,19 +14,21 @@ from collections import OrderedDict
with
CodeGeneration
()
as
ctx
:
generatedMethod
=
"
TRTlike
"
#generatedMethod =
"
D3Q27TRTlike
"
#generatedMethod =
"
cumulant
"
generatedMethod
=
'
TRTlike
'
#generatedMethod =
'
D3Q27TRTlike
'
#generatedMethod =
'
cumulant
'
clear_cache
()
dtype_string
=
'float64'
if
ctx
.
double_accuracy
else
'float32'
cpu_vectorize_info
=
{
'instruction_set'
:
get_vectorize_instruction_set
(
ctx
)}
if
generatedMethod
==
"
TRTlike
"
:
omegaVisc
=
sp
.
Symbol
(
"
omega_visc
"
)
omegaBulk
=
ps
.
fields
(
"
omega_bulk: [3D]
"
,
layout
=
'fzyx'
)
omegaMagic
=
sp
.
Symbol
(
"
omega_magic
"
)
stencil
=
get_stencil
(
"
D3Q19
"
,
'walberla'
)
if
generatedMethod
==
'
TRTlike
'
:
omegaVisc
=
sp
.
Symbol
(
'
omega_visc
'
)
omegaBulk
=
ps
.
fields
(
f
'
omega_bulk:
{
dtype_string
}
[3D]
'
,
layout
=
'fzyx'
)
omegaMagic
=
sp
.
Symbol
(
'
omega_magic
'
)
stencil
=
get_stencil
(
'
D3Q19
'
,
'walberla'
)
x
,
y
,
z
=
MOMENT_SYMBOLS
one
=
sp
.
Rational
(
1
,
1
)
...
...
@@ -53,61 +55,62 @@ with CodeGeneration() as ctx:
generate_lattice_model
(
ctx
,
'GeneratedLBM'
,
collision_rule
,
field_layout
=
'fzyx'
,
cpu_vectorize_info
=
cpu_vectorize_info
)
if
generatedMethod
==
"
D3Q27TRTlike
"
:
if
generatedMethod
==
'
D3Q27TRTlike
'
:
omegaVisc
=
sp
.
Symbol
(
"
omega_visc
"
)
omegaBulk
=
ps
.
fields
(
"
omega_bulk: [3D]
"
,
layout
=
'fzyx'
)
omegaMagic
=
sp
.
Symbol
(
"
omega_magic
"
)
stencil
=
get_stencil
(
"
D3Q27
"
,
'walberla'
)
omegaVisc
=
sp
.
Symbol
(
'
omega_visc
'
)
omegaBulk
=
ps
.
fields
(
f
'
omega_bulk:
{
dtype_string
}
[3D]
'
,
layout
=
'fzyx'
)
omegaMagic
=
sp
.
Symbol
(
'
omega_magic
'
)
stencil
=
get_stencil
(
'
D3Q27
'
,
'walberla'
)
relaxation_rates
=
[
omegaVisc
,
omegaBulk
.
center_vector
,
omegaMagic
,
omegaVisc
,
omegaMagic
,
omegaVisc
]
relaxation_rates
=
[
omegaVisc
,
omegaBulk
.
center_vector
,
omegaMagic
,
omegaVisc
,
omegaMagic
,
omegaVisc
]
methodWithForce
=
create_lb_method
(
stencil
=
stencil
,
method
=
'mrt'
,
maxwellian_moments
=
False
,
weighted
=
True
,
compressible
=
False
,
relaxation_rates
=
relaxation_rates
)
methodWithForce
=
create_lb_method
(
stencil
=
stencil
,
method
=
'mrt'
,
maxwellian_moments
=
False
,
weighted
=
True
,
compressible
=
False
,
relaxation_rates
=
relaxation_rates
)
collision_rule
=
create_lb_collision_rule
(
lb_method
=
methodWithForce
,
optimization
=
{
'cse_global'
:
True
})
generate_lattice_model
(
ctx
,
'GeneratedLBM'
,
collision_rule
,
field_layout
=
'fzyx'
,
cpu_vectorize_info
=
cpu_vectorize_info
)
generate_lattice_model
(
ctx
,
'GeneratedLBM'
,
collision_rule
,
field_layout
=
'fzyx'
,
cpu_vectorize_info
=
cpu_vectorize_info
)
if
generatedMethod
==
"
cumulant
"
:
if
generatedMethod
==
'
cumulant
'
:
x
,
y
,
z
=
MOMENT_SYMBOLS
x
,
y
,
z
=
MOMENT_SYMBOLS
cumulants
=
[
0
]
*
27
cumulants
[
0
]
=
sp
.
sympify
(
1
)
#
000
cumulants
[
0
]
=
sp
.
sympify
(
1
)
#
000
cumulants
[
1
]
=
x
#
100
cumulants
[
2
]
=
y
#
010
cumulants
[
3
]
=
z
#
001
cumulants
[
1
]
=
x
#
100
cumulants
[
2
]
=
y
#
010
cumulants
[
3
]
=
z
#
001
cumulants
[
4
]
=
x
*
y
#
110
cumulants
[
5
]
=
x
*
z
#
101
cumulants
[
6
]
=
y
*
z
#
011
cumulants
[
4
]
=
x
*
y
#
110
cumulants
[
5
]
=
x
*
z
#
101
cumulants
[
6
]
=
y
*
z
#
011
cumulants
[
7
]
=
x
**
2
-
y
**
2
#
200 - 020
cumulants
[
8
]
=
x
**
2
-
z
**
2
#
200 - 002
cumulants
[
9
]
=
x
**
2
+
y
**
2
+
z
**
2
#
200 + 020 + 002
cumulants
[
7
]
=
x
**
2
-
y
**
2
#
200 - 020
cumulants
[
8
]
=
x
**
2
-
z
**
2
#
200 - 002
cumulants
[
9
]
=
x
**
2
+
y
**
2
+
z
**
2
#
200 + 020 + 002
cumulants
[
10
]
=
x
*
y
**
2
+
x
*
z
**
2
#
120 + 102
cumulants
[
11
]
=
x
**
2
*
y
+
y
*
z
**
2
#
210 + 012
cumulants
[
12
]
=
x
**
2
*
z
+
y
**
2
*
z
#
201 + 021
cumulants
[
13
]
=
x
*
y
**
2
-
x
*
z
**
2
#
120 - 102
cumulants
[
14
]
=
x
**
2
*
y
-
y
*
z
**
2
#
210 - 012
cumulants
[
15
]
=
x
**
2
*
z
-
y
**
2
*
z
#
201 - 021
cumulants
[
10
]
=
x
*
y
**
2
+
x
*
z
**
2
#
120 + 102
cumulants
[
11
]
=
x
**
2
*
y
+
y
*
z
**
2
#
210 + 012
cumulants
[
12
]
=
x
**
2
*
z
+
y
**
2
*
z
#
201 + 021
cumulants
[
13
]
=
x
*
y
**
2
-
x
*
z
**
2
#
120 - 102
cumulants
[
14
]
=
x
**
2
*
y
-
y
*
z
**
2
#
210 - 012
cumulants
[
15
]
=
x
**
2
*
z
-
y
**
2
*
z
#
201 - 021
cumulants
[
16
]
=
x
*
y
*
z
#
111
cumulants
[
16
]
=
x
*
y
*
z
#
111
cumulants
[
17
]
=
x
**
2
*
y
**
2
-
2
*
x
**
2
*
z
**
2
+
y
**
2
*
z
**
2
# 220- 2*202 +022
cumulants
[
18
]
=
x
**
2
*
y
**
2
+
x
**
2
*
z
**
2
-
2
*
y
**
2
*
z
**
2
# 220 + 202 - 2*002
cumulants
[
19
]
=
x
**
2
*
y
**
2
+
x
**
2
*
z
**
2
+
y
**
2
*
z
**
2
# 220 + 202 + 022
cumulants
[
17
]
=
x
**
2
*
y
**
2
-
2
*
x
**
2
*
z
**
2
+
y
**
2
*
z
**
2
# 220- 2*202 +022
cumulants
[
18
]
=
x
**
2
*
y
**
2
+
x
**
2
*
z
**
2
-
2
*
y
**
2
*
z
**
2
# 220 + 202 - 2*002
cumulants
[
19
]
=
x
**
2
*
y
**
2
+
x
**
2
*
z
**
2
+
y
**
2
*
z
**
2
# 220 + 202 + 022
cumulants
[
20
]
=
x
**
2
*
y
*
z
# 211
cumulants
[
21
]
=
x
*
y
**
2
*
z
# 121
cumulants
[
22
]
=
x
*
y
*
z
**
2
# 112
cumulants
[
20
]
=
x
**
2
*
y
*
z
# 211
cumulants
[
21
]
=
x
*
y
**
2
*
z
# 121
cumulants
[
22
]
=
x
*
y
*
z
**
2
# 112
cumulants
[
23
]
=
x
**
2
*
y
**
2
*
z
# 221
cumulants
[
24
]
=
x
**
2
*
y
*
z
**
2
# 212
cumulants
[
25
]
=
x
*
y
**
2
*
z
**
2
# 122
cumulants
[
23
]
=
x
**
2
*
y
**
2
*
z
# 221
cumulants
[
24
]
=
x
**
2
*
y
*
z
**
2
# 212
cumulants
[
25
]
=
x
*
y
**
2
*
z
**
2
# 122
cumulants
[
26
]
=
x
**
2
*
y
**
2
*
z
**
2
# 222
...
...
@@ -119,9 +122,9 @@ with CodeGeneration() as ctx:
else
:
return
1
stencil
=
get_stencil
(
"
D3Q27
"
,
'walberla'
)
stencil
=
get_stencil
(
'
D3Q27
'
,
'walberla'
)
omega
=
sp
.
Symbol
(
"
omega
"
)
omega
=
sp
.
Symbol
(
'
omega
'
)
rr_dict
=
OrderedDict
((
c
,
get_relaxation_rate
(
c
,
omega
))
for
c
in
cumulants
)
...
...
@@ -129,12 +132,13 @@ with CodeGeneration() as ctx:
my_method
=
create_with_continuous_maxwellian_eq_moments
(
stencil
,
rr_dict
,
cumulant
=
True
,
compressible
=
True
,)
collision_rule
=
create_lb_collision_rule
(
lb_method
=
my_method
,
optimization
=
{
"
cse_global
"
:
True
,
"
cse_pdfs
"
:
False
})
optimization
=
{
'
cse_global
'
:
True
,
'
cse_pdfs
'
:
False
})
print
(
my_method
.
relaxation_rates
)
generate_lattice_model
(
ctx
,
'GeneratedLBM'
,
collision_rule
,
field_layout
=
'fzyx'
,
cpu_vectorize_info
=
cpu_vectorize_info
)
generate_lattice_model
(
ctx
,
'GeneratedLBM'
,
collision_rule
,
field_layout
=
'fzyx'
,
cpu_vectorize_info
=
cpu_vectorize_info
)
apps/benchmarks/FluidParticleCoupling/GeneratedLBMWithForce.py
View file @
7751b2df
...
...
@@ -15,25 +15,27 @@ from collections import OrderedDict
with
CodeGeneration
()
as
ctx
:
forcing
=
(
sp
.
symbols
(
"
fx
"
),
0
,
0
)
forcemodel
=
Luo
(
forcing
)
forcing
=
(
sp
.
symbols
(
'
fx
'
),
0
,
0
)
forcemodel
=
Luo
(
forcing
)
generatedMethod
=
"
TRTlike
"
#generatedMethod =
"
D3Q27TRTlike
"
#generatedMethod =
"
cumulant
"
#generatedMethod =
"
cumulantTRT
"
generatedMethod
=
'
TRTlike
'
#
generatedMethod =
'
D3Q27TRTlike
'
#
generatedMethod =
'
cumulant
'
#
generatedMethod =
'
cumulantTRT
'
print
(
"
Generating
"
+
generatedMethod
+
"
LBM method
"
)
print
(
'
Generating
'
+
generatedMethod
+
'
LBM method
'
)
clear_cache
()
dtype_string
=
'float64'
if
ctx
.
double_accuracy
else
'float32'
cpu_vectorize_info
=
{
'instruction_set'
:
get_vectorize_instruction_set
(
ctx
)}
if
generatedMethod
==
"
TRTlike
"
:
omegaVisc
=
sp
.
Symbol
(
"
omega_visc
"
)
omegaBulk
=
ps
.
fields
(
"
omega_bulk: [3D]
"
,
layout
=
'fzyx'
)
omegaMagic
=
sp
.
Symbol
(
"
omega_magic
"
)
stencil
=
get_stencil
(
"
D3Q19
"
,
'walberla'
)
if
generatedMethod
==
'
TRTlike
'
:
omegaVisc
=
sp
.
Symbol
(
'
omega_visc
'
)
omegaBulk
=
ps
.
fields
(
f
'
omega_bulk:
{
dtype_string
}
[3D]
'
,
layout
=
'fzyx'
)
omegaMagic
=
sp
.
Symbol
(
'
omega_magic
'
)
stencil
=
get_stencil
(
'
D3Q19
'
,
'walberla'
)
x
,
y
,
z
=
MOMENT_SYMBOLS
one
=
sp
.
Rational
(
1
,
1
)
...
...
@@ -59,63 +61,64 @@ with CodeGeneration() as ctx:
collision_rule
=
create_lb_collision_rule
(
lb_method
=
methodWithForce
,
optimization
=
{
'cse_global'
:
True
})
generate_lattice_model
(
ctx
,
'GeneratedLBMWithForce'
,
collision_rule
,
field_layout
=
'fzyx'
,
cpu_vectorize_info
=
cpu_vectorize_info
)
if
generatedMethod
==
"
D3Q27TRTlike
"
:
if
generatedMethod
==
'
D3Q27TRTlike
'
:
omegaVisc
=
sp
.
Symbol
(
"
omega_visc
"
)
omegaBulk
=
ps
.
fields
(
"
omega_bulk: [3D]
"
,
layout
=
'fzyx'
)
omegaMagic
=
sp
.
Symbol
(
"
omega_magic
"
)
stencil
=
get_stencil
(
"
D3Q27
"
,
'walberla'
)
omegaVisc
=
sp
.
Symbol
(
'
omega_visc
'
)
omegaBulk
=
ps
.
fields
(
f
'
omega_bulk:
{
dtype_string
}
[3D]
'
,
layout
=
'fzyx'
)
omegaMagic
=
sp
.
Symbol
(
'
omega_magic
'
)
stencil
=
get_stencil
(
'
D3Q27
'
,
'walberla'
)
relaxation_rates
=
[
omegaVisc
,
omegaBulk
.
center_vector
,
omegaMagic
,
omegaVisc
,
omegaMagic
,
omegaVisc
]
relaxation_rates
=
[
omegaVisc
,
omegaBulk
.
center_vector
,
omegaMagic
,
omegaVisc
,
omegaMagic
,
omegaVisc
]
methodWithForce
=
create_lb_method
(
stencil
=
stencil
,
method
=
'mrt'
,
maxwellian_moments
=
False
,
weighted
=
True
,
compressible
=
False
,
force_model
=
forcemodel
,
relaxation_rates
=
relaxation_rates
)
methodWithForce
=
create_lb_method
(
stencil
=
stencil
,
method
=
'mrt'
,
maxwellian_moments
=
False
,
weighted
=
True
,
compressible
=
False
,
force_model
=
forcemodel
,
relaxation_rates
=
relaxation_rates
)
collision_rule
=
create_lb_collision_rule
(
lb_method
=
methodWithForce
,
optimization
=
{
'cse_global'
:
True
})
generate_lattice_model
(
ctx
,
'GeneratedLBMWithForce'
,
collision_rule
,
field_layout
=
'fzyx'
,
cpu_vectorize_info
=
cpu_vectorize_info
)
generate_lattice_model
(
ctx
,
'GeneratedLBMWithForce'
,
collision_rule
,
field_layout
=
'fzyx'
,
cpu_vectorize_info
=
cpu_vectorize_info
)
if
generatedMethod
==
"
cumulant
"
:
if
generatedMethod
==
'
cumulant
'
:
x
,
y
,
z
=
MOMENT_SYMBOLS
x
,
y
,
z
=
MOMENT_SYMBOLS
cumulants
=
[
0
]
*
27
cumulants
[
0
]
=
sp
.
sympify
(
1
)
#
000
cumulants
[
0
]
=
sp
.
sympify
(
1
)
#
000
cumulants
[
1
]
=
x
#
100
cumulants
[
2
]
=
y
#
010
cumulants
[
3
]
=
z
#
001
cumulants
[
1
]
=
x
#
100
cumulants
[
2
]
=
y
#
010
cumulants
[
3
]
=
z
#
001
cumulants
[
4
]
=
x
*
y
#
110
cumulants
[
5
]
=
x
*
z
#
101
cumulants
[
6
]
=
y
*
z
#
011
cumulants
[
4
]
=
x
*
y
#
110
cumulants
[
5
]
=
x
*
z
#
101
cumulants
[
6
]
=
y
*
z
#
011
cumulants
[
7
]
=
x
**
2
-
y
**
2
#
200 - 020
cumulants
[
8
]
=
x
**
2
-
z
**
2
#
200 - 002
cumulants
[
9
]
=
x
**
2
+
y
**
2
+
z
**
2
#
200 + 020 + 002
cumulants
[
7
]
=
x
**
2
-
y
**
2
#
200 - 020
cumulants
[
8
]
=
x
**
2
-
z
**
2
#
200 - 002
cumulants
[
9
]
=
x
**
2
+
y
**
2
+
z
**
2
#
200 + 020 + 002
cumulants
[
10
]
=
x
*
y
**
2
+
x
*
z
**
2
#
120 + 102
cumulants
[
11
]
=
x
**
2
*
y
+
y
*
z
**
2
#
210 + 012
cumulants
[
12
]
=
x
**
2
*
z
+
y
**
2
*
z
#
201 + 021
cumulants
[
13
]
=
x
*
y
**
2
-
x
*
z
**
2
#
120 - 102
cumulants
[
14
]
=
x
**
2
*
y
-
y
*
z
**
2
#
210 - 012
cumulants
[
15
]
=
x
**
2
*
z
-
y
**
2
*
z
#
201 - 021
cumulants
[
10
]
=
x
*
y
**
2
+
x
*
z
**
2
#
120 + 102
cumulants
[
11
]
=
x
**
2
*
y
+
y
*
z
**
2
#
210 + 012
cumulants
[
12
]
=
x
**
2
*
z
+
y
**
2
*
z
#
201 + 021
cumulants
[
13
]
=
x
*
y
**
2
-
x
*
z
**
2
#
120 - 102
cumulants
[
14
]
=
x
**
2
*
y
-
y
*
z
**
2
#
210 - 012
cumulants
[
15
]
=
x
**
2
*
z
-
y
**
2
*
z
#
201 - 021
cumulants
[
16
]
=
x
*
y
*
z
#
111
cumulants
[
16
]
=
x
*
y
*
z
#
111
cumulants
[
17
]
=
x
**
2
*
y
**
2
-
2
*
x
**
2
*
z
**
2
+
y
**
2
*
z
**
2
# 220- 2*202 +022
cumulants
[
18
]
=
x
**
2
*
y
**
2
+
x
**
2
*
z
**
2
-
2
*
y
**
2
*
z
**
2
# 220 + 202 - 2*002
cumulants
[
19
]
=
x
**
2
*
y
**
2
+
x
**
2
*
z
**
2
+
y
**
2
*
z
**
2
# 220 + 202 + 022
cumulants
[
17
]
=
x
**
2
*
y
**
2
-
2
*
x
**
2
*
z
**
2
+
y
**
2
*
z
**
2
# 220- 2*202 +022
cumulants
[
18
]
=
x
**
2
*
y
**
2
+
x
**
2
*
z
**
2
-
2
*
y
**
2
*
z
**
2
# 220 + 202 - 2*002
cumulants
[
19
]
=
x
**
2
*
y
**
2
+
x
**
2
*
z
**
2
+
y
**
2
*
z
**
2
# 220 + 202 + 022
cumulants
[
20
]
=
x
**
2
*
y
*
z
# 211
cumulants
[
21
]
=
x
*
y
**
2
*
z
# 121
cumulants
[
22
]
=
x
*
y
*
z
**
2
# 112
cumulants
[
20
]
=
x
**
2
*
y
*
z
# 211
cumulants
[
21
]
=
x
*
y
**
2
*
z
# 121
cumulants
[
22
]
=
x
*
y
*
z
**
2
# 112
cumulants
[
23
]
=
x
**
2
*
y
**
2
*
z
# 221
cumulants
[
24
]
=
x
**
2
*
y
*
z
**
2
# 212
cumulants
[
25
]
=
x
*
y
**
2
*
z
**
2
# 122
cumulants
[
23
]
=
x
**
2
*
y
**
2
*
z
# 221
cumulants
[
24
]
=
x
**
2
*
y
*
z
**
2
# 212
cumulants
[
25
]
=
x
*
y
**
2
*
z
**
2
# 122
cumulants
[
26
]
=
x
**
2
*
y
**
2
*
z
**
2
# 222
cumulants
[
26
]
=
x
**
2
*
y
**
2
*
z
**
2
# 222
def
get_relaxation_rate
(
cumulant
,
omega
):
if
get_order
(
cumulant
)
<=
1
:
...
...
@@ -125,66 +128,68 @@ with CodeGeneration() as ctx:
else
:
return
1
stencil
=
get_stencil
(
"
D3Q27
"
,
'walberla'
)
stencil
=
get_stencil
(
'
D3Q27
'
,
'walberla'
)
omega
=
sp
.
Symbol
(
"
omega
"
)
omega
=
sp
.
Symbol
(
'
omega
'
)
rr_dict
=
OrderedDict
((
c
,
get_relaxation_rate
(
c
,
omega
))
for
c
in
cumulants
)
from
lbmpy.methods
import
create_with_continuous_maxwellian_eq_moments
my_method
=
create_with_continuous_maxwellian_eq_moments
(
stencil
,
rr_dict
,
cumulant
=
True
,
compressible
=
True
,
force_model
=
forcemodel
)
my_method
=
create_with_continuous_maxwellian_eq_moments
(
stencil
,
rr_dict
,
cumulant
=
True
,
compressible
=
True
,
force_model
=
forcemodel
)
collision_rule
=
create_lb_collision_rule
(
lb_method
=
my_method
,
optimization
=
{
"
cse_global
"
:
True
,
"
cse_pdfs
"
:
False
})
optimization
=
{
'
cse_global
'
:
True
,
'
cse_pdfs
'
:
False
})
print
(
my_method
.
relaxation_rates
)
generate_lattice_model
(
ctx
,
'GeneratedLBMWithForce'
,
collision_rule
,
field_layout
=
'fzyx'
,
cpu_vectorize_info
=
cpu_vectorize_info
)
generate_lattice_model
(
ctx
,
'GeneratedLBMWithForce'
,
collision_rule
,
field_layout
=
'fzyx'
,
cpu_vectorize_info
=
cpu_vectorize_info
)
if
generatedMethod
==
"
cumulantTRT
"
:
if
generatedMethod
==
'
cumulantTRT
'
:
x
,
y
,
z
=
MOMENT_SYMBOLS
x
,
y
,
z
=
MOMENT_SYMBOLS
cumulants
=
[
0
]
*
27
cumulants
[
0
]
=
sp
.
sympify
(
1
)
#
000
cumulants
[
0
]
=
sp
.
sympify
(
1
)
#
000
cumulants
[
1
]
=
x
#
100
cumulants
[
2
]
=
y
#
010
cumulants
[
3
]
=
z
#
001
cumulants
[
1
]
=
x
#
100
cumulants
[
2
]
=
y
#
010
cumulants
[
3
]
=
z
#
001
cumulants
[
4
]
=
x
*
y
#
110
cumulants
[
5
]
=
x
*
z
#
101
cumulants
[
6
]
=
y
*
z
#
011
cumulants
[
4
]
=
x
*
y
#
110
cumulants
[
5
]
=
x
*
z
#
101
cumulants
[
6
]
=
y
*
z
#
011
cumulants
[
7
]
=
x
**
2
-
y
**
2
#
200 - 020
cumulants
[
8
]
=
x
**
2
-
z
**
2
#
200 - 002
cumulants
[
9
]
=
x
**
2
+
y
**
2
+
z
**
2
#
200 + 020 + 002
cumulants
[
7
]
=
x
**
2
-
y
**
2
#
200 - 020
cumulants
[
8
]
=
x
**
2
-
z
**
2
#
200 - 002
cumulants
[
9
]
=
x
**
2
+
y
**
2
+
z
**
2
#
200 + 020 + 002
cumulants
[
10
]
=
x
*
y
**
2
+
x
*
z
**
2
#
120 + 102
cumulants
[
11
]
=
x
**
2
*
y
+
y
*
z
**
2
#
210 + 012
cumulants
[
12
]
=
x
**
2
*
z
+
y
**
2
*
z
#
201 + 021
cumulants
[
13
]
=
x
*
y
**
2
-
x
*
z
**
2
#
120 - 102
cumulants
[
14
]
=
x
**
2
*
y
-
y
*
z
**
2
#
210 - 012
cumulants
[
15
]
=
x
**
2
*
z
-
y
**
2
*
z
#
201 - 021
cumulants
[
10
]
=
x
*
y
**
2
+
x
*
z
**
2
#
120 + 102
cumulants
[
11
]
=
x
**
2
*
y
+
y
*
z
**
2
#
210 + 012
cumulants
[
12
]
=
x
**
2
*
z
+
y
**
2
*
z
#
201 + 021
cumulants
[
13
]
=
x
*
y
**
2
-
x
*
z
**
2
#
120 - 102
cumulants
[
14
]
=
x
**
2
*
y
-
y
*
z
**
2
#
210 - 012
cumulants
[
15
]
=
x
**
2
*
z
-
y
**
2
*
z
#
201 - 021
cumulants
[
16
]
=
x
*
y
*
z
#
111
cumulants
[
16
]
=
x
*
y
*
z
#
111
cumulants
[
17
]
=
x
**
2
*
y
**
2
-
2
*
x
**
2
*
z
**
2
+
y
**
2
*
z
**
2
# 220- 2*202 +022
cumulants
[
18
]
=
x
**
2
*
y
**
2
+
x
**
2
*
z
**
2
-
2
*
y
**
2
*
z
**
2
# 220 + 202 - 2*002
cumulants
[
19
]
=
x
**
2
*
y
**
2
+
x
**
2
*
z
**
2
+
y
**
2
*
z
**
2
# 220 + 202 + 022
cumulants
[
17
]
=
x
**
2
*
y
**
2
-
2
*
x
**
2
*
z
**
2
+
y
**
2
*
z
**
2
# 220- 2*202 +022
cumulants
[
18
]
=
x
**
2
*
y
**
2
+
x
**
2
*
z
**
2
-
2
*
y
**
2
*
z
**
2
# 220 + 202 - 2*002
cumulants
[
19
]
=
x
**
2
*
y
**
2
+
x
**
2
*
z
**
2
+
y
**
2
*
z
**
2
# 220 + 202 + 022
cumulants
[
20
]
=
x
**
2
*
y
*
z
# 211
cumulants
[
21
]
=
x
*
y
**
2
*
z
# 121
cumulants
[
22
]
=
x
*
y
*
z
**
2
# 112
cumulants
[
20
]
=
x
**
2
*
y
*
z
# 211
cumulants
[
21
]
=
x
*
y
**
2
*
z
# 121
cumulants
[
22
]
=
x
*
y
*
z
**
2
# 112
cumulants
[
23
]
=
x
**
2
*
y
**
2
*
z
# 221
cumulants
[
24
]
=
x
**
2
*
y
*
z
**
2
# 212
cumulants
[
25
]
=
x
*
y
**
2
*
z
**
2
# 122
cumulants
[
23
]
=
x
**
2
*
y
**
2
*
z
# 221
cumulants
[
24
]
=
x
**
2
*
y
*
z
**
2
# 212
cumulants
[
25
]
=
x
*
y
**
2
*
z
**
2
# 122
cumulants
[
26
]
=
x
**
2
*
y
**
2
*
z
**
2
# 222
cumulants
[
26
]
=
x
**
2
*
y
**
2
*
z
**
2
# 222
def
get_relaxation_rate
(
cumulant
,
omegaVisc
,
omegaMagic
):
if
get_order
(
cumulant
)
<=
1
:
...
...
@@ -194,22 +199,24 @@ with CodeGeneration() as ctx:
else
:
return
omegaMagic
stencil
=
get_stencil
(
"
D3Q27
"
,
'walberla'
)
stencil
=
get_stencil
(
'
D3Q27
'
,
'walberla'
)
omegaVisc
=
sp
.
Symbol
(
"
omega_visc
"
)
omegaMagic
=
sp
.
Symbol
(
"
omega_magic
"
)
omegaVisc
=
sp
.
Symbol
(
'
omega_visc
'
)
omegaMagic
=
sp
.
Symbol
(
'
omega_magic
'
)
rr_dict
=
OrderedDict
((
c
,
get_relaxation_rate
(
c
,
omegaVisc
,
omegaMagic
))
for
c
in
cumulants
)
from
lbmpy.methods
import
create_with_continuous_maxwellian_eq_moments
my_method
=
create_with_continuous_maxwellian_eq_moments
(
stencil
,
rr_dict
,
cumulant
=
True
,
compressible
=
True
,
force_model
=
forcemodel
)
my_method
=
create_with_continuous_maxwellian_eq_moments
(
stencil
,
rr_dict
,
cumulant
=
True
,
compressible
=
True
,
force_model
=
forcemodel
)
collision_rule
=
create_lb_collision_rule
(
lb_method
=
my_method
,
optimization
=
{
"
cse_global
"
:
True
,
"
cse_pdfs
"
:
False
})
optimization
=
{
'
cse_global
'
:
True
,
'
cse_pdfs
'
:
False
})
print
(
my_method
.
relaxation_rates
)
generate_lattice_model
(
ctx
,
'GeneratedLBMWithForce'
,
collision_rule
,
field_layout
=
'fzyx'
,
cpu_vectorize_info
=
cpu_vectorize_info
)
generate_lattice_model
(
ctx
,
'GeneratedLBMWithForce'
,
collision_rule
,
field_layout
=
'fzyx'
,
cpu_vectorize_info
=
cpu_vectorize_info
)
apps/benchmarks/PhaseFieldAllenCahn/InitializerFunctions.cpp
View file @
7751b2df
...
...
@@ -40,10 +40,10 @@ void initPhaseField_bubble(const shared_ptr< StructuredBlockStorage >& blocks, B
// clang-format off
WALBERLA_FOR_ALL_CELLS_INCLUDING_GHOST_LAYER_XYZ
(
phaseField
,
Cell
globalCell
;
blocks
->
transformBlockLocalToGlobalCell
(
globalCell
,
block
,
Cell
(
x
,
y
,
z
));
real_t
Ri
=
s
qrt
(
(
globalCell
[
0
]
-
bubbleMidPoint
[
0
])
*
(
globalCell
[
0
]
-
bubbleMidPoint
[
0
])
+
(
globalCell
[
1
]
-
bubbleMidPoint
[
1
])
*
(
globalCell
[
1
]
-
bubbleMidPoint
[
1
])
+
(
globalCell
[
2
]
-
bubbleMidPoint
[
2
])
*
(
globalCell
[
2
]
-
bubbleMidPoint
[
2
]));
phaseField
->
get
(
x
,
y
,
z
)
=
0.5
+
0.5
*
tanh
(
2.0
*
(
Ri
-
R
)
/
W
);
real_t
Ri
=
s
td
::
sqrt
((
real_c
(
globalCell
[
0
]
)
-
bubbleMidPoint
[
0
])
*
(
real_c
(
globalCell
[
0
]
)
-
bubbleMidPoint
[
0
])
+
(
real_c
(
globalCell
[
1
]
)
-
bubbleMidPoint
[
1
])
*
(
real_c
(
globalCell
[
1
]
)
-
bubbleMidPoint
[
1
])
+
(
real_c
(
globalCell
[
2
]
)
-
bubbleMidPoint
[
2
])
*
(
real_c
(
globalCell
[
2
]
)
-
bubbleMidPoint
[
2
]));
phaseField
->
get
(
x
,
y
,
z
)
=
real_t
(
0.5
)
+
real_t
(
0.5
)
*
std
::
tanh
(
real_t
(
2.0
)
*
(
Ri
-
R
)
/
W
);
)
// clang-format on
}
...
...
@@ -52,9 +52,9 @@ void initPhaseField_bubble(const shared_ptr< StructuredBlockStorage >& blocks, B
void
initPhaseField_RTI
(
const
shared_ptr
<
StructuredBlockStorage
>&
blocks
,
BlockDataID
phaseFieldID
,
const
real_t
W
=
5
)
{
auto
X
=
blocks
->
getDomainCellBB
().
xMax
();
auto
halfY
=
(
blocks
->
getDomainCellBB
().
yMax
())
/
2.0
;
double
perturbation
=
0.05
;
const
real_t
X
=
real_c
(
blocks
->
getDomainCellBB
().
xMax
()
)
;
const
real_t
halfY
=
real_c
(
blocks
->
getDomainCellBB
().
yMax
())
/
real_t
(
2.0
)
;
const
real_t
perturbation
=
real_t
(
0.05
)
;
for
(
auto
&
block
:
*
blocks
)
{
...
...
@@ -62,8 +62,9 @@ void initPhaseField_RTI(const shared_ptr< StructuredBlockStorage >& blocks, Bloc
// clang-format off
WALBERLA_FOR_ALL_CELLS_INCLUDING_GHOST_LAYER_XYZ
(
phaseField
,
Cell
globalCell
;
blocks
->
transformBlockLocalToGlobalCell
(
globalCell
,
block
,
Cell
(
x
,
y
,
z
));
real_t
tmp
=
perturbation
*
X
*
(
cos
((
2.0
*
math
::
pi
*
globalCell
[
0
])
/
X
)
+
cos
((
2.0
*
math
::
pi
*
globalCell
[
2
])
/
X
));
phaseField
->
get
(
x
,
y
,
z
)
=
0.5
+
0.5
*
tanh
(((
globalCell
[
1
]
-
halfY
)
-
tmp
)
/
(
W
/
2.0
));
real_t
tmp
=
perturbation
*
X
*
(
std
::
cos
((
real_t
(
2.0
)
*
math
::
pi
*
real_c
(
globalCell
[
0
]))
/
X
)
+
std
::
cos
((
real_t
(
2.0
)
*
math
::
pi
*
real_c
(
globalCell
[
2
]))
/
X
));
phaseField
->
get
(
x
,
y
,
z
)
=
real_t
(
0.5
)
+
real_t
(
0.5
)
*
std
::
tanh
(((
real_t
(
globalCell
[
1
])
-
halfY
)
-
tmp
)
/
(
W
/
real_t
(
2.0
)));
)
// clang-format on
}
...
...
apps/benchmarks/PhaseFieldAllenCahn/benchmark_multiphase.cpp
View file @
7751b2df
...
...
@@ -84,7 +84,7 @@ using flag_t = walberla::uint8_t;
using
FlagField_T
=
FlagField
<
flag_t
>
;
#if defined(WALBERLA_BUILD_WITH_CUDA)
typedef
cuda
::
GPUField
<
double
>
GPUField
;
typedef
cuda
::
GPUField
<
real_t
>
GPUField
;
#endif
// using CommScheme_T = cuda::communication::UniformGPUScheme<stencil::D2Q9>;
...
...
apps/benchmarks/UniformGridGPU/InitShearVelocity.h
View file @
7751b2df
...
...
@@ -6,7 +6,7 @@ namespace walberla {
inline
void
initShearVelocity
(
const
shared_ptr
<
StructuredBlockStorage
>
&
blocks
,
BlockDataID
velFieldID
,
const
real_t
xMagnitude
=
0.005
,
const
real_t
fluctuationMagnitude
=
0.05
)
const
real_t
xMagnitude
=
real_t
(
0.005
)
,
const
real_t
fluctuationMagnitude
=
real_t
(
0.05
)
)
{
math
::
seedRandomGenerator
(
0
);
auto
halfZ
=
blocks
->
getDomainCellBB
().
zMax
()
/
2
;
...
...
apps/benchmarks/UniformGridGPU/UniformGridGPU.cpp
View file @
7751b2df
...
...
@@ -55,7 +55,7 @@ using FlagField_T = FlagField<flag_t>;
void
initShearVelocity
(
const
shared_ptr
<
StructuredBlockStorage
>
&
blocks
,
BlockDataID
velFieldID
,
const
real_t
xMagnitude
=
0.1
,
const
real_t
fluctuationMagnitude
=
0.05
)
const
real_t
xMagnitude
=
real_t
(
0.1
)
,
const
real_t
fluctuationMagnitude
=
real_t
(
0.05
)
)
{
math
::
seedRandomGenerator
(
0
);
auto
halfZ
=
blocks
->
getDomainCellBB
().
zMax
()
/
2
;
...
...
apps/benchmarks/UniformGridGPU/UniformGridGPU.py
View file @
7751b2df
...
...
@@ -135,7 +135,7 @@ with CodeGeneration() as ctx:
]
lb_method
=
create_lb_method
(
**
options
)
update_rule
=
create_lb_update_rule
(
lb_method
=
lb_method
,
**
options
)
if
not
noopt
:
update_rule
=
insert_fast_divisions
(
update_rule
)
update_rule
=
insert_fast_sqrts
(
update_rule
)
...
...
apps/benchmarks/UniformGridGenerated/InitShearVelocity.h
View file @
7751b2df
...
...
@@ -6,7 +6,7 @@ namespace walberla {
inline
void
initShearVelocity
(
const
shared_ptr
<
StructuredBlockStorage
>
&
blocks
,
BlockDataID
velFieldID
,
const
real_t
xMagnitude
=
0.005
,
const
real_t
fluctuationMagnitude
=
0.05
)
const
real_t
xMagnitude
=
real_t
(
0.005
)
,
const
real_t
fluctuationMagnitude
=
real_t
(
0.05
)
)
{
math
::
seedRandomGenerator
(
0
);
auto
halfZ
=
blocks
->
getDomainCellBB
().
zMax
()
/
2
;
...
...
apps/benchmarks/UniformGridGenerated/ManualKernels.h
View file @
7751b2df
...
...
@@ -30,20 +30,20 @@ public:
for
(
auto
d
=
LM
::
Stencil
::
begin
();
d
!=
LM
::
Stencil
::
end
();
++
d
)
{
const
real_t
pdf
=
src
->
get
(
x
-
d
.
cx
(),
y
-
d
.
cy
(),
z
-
d
.
cz
(),
d
.
toIdx
());
u
[
0
]
+=
pdf
*
d
.
cx
();
u
[
1
]
+=
pdf
*
d
.
cy
();
u
[
2
]
+=
pdf
*
d
.
cz
();
u
[
0
]
+=
pdf
*
real_c
(
d
.
cx
()
)
;
u
[
1
]
+=
pdf
*
real_c
(
d
.
cy
()
)
;
u
[
2
]
+=
pdf
*
real_c
(
d
.
cz
()
)
;
rho
+=
pdf
;