Skip to content
GitLab
Menu
Projects
Groups
Snippets
Loading...
Help
Help
Support
Community forum
Keyboard shortcuts
?
Submit feedback
Contribute to GitLab
Sign in
Toggle navigation
Menu
Open sidebar
Stephan Seitz
pystencils
Commits
ebfd56e1
Commit
ebfd56e1
authored
Jun 15, 2020
by
Markus Holzer
Browse files
fixed docstring bugs and expm1 optimisation
parent
35c029a0
Changes
5
Show whitespace changes
Inline
Side-by-side
pystencils/assignment.py
View file @
ebfd56e1
...
...
@@ -66,6 +66,7 @@ try:
if
int
(
sympy_version
[
0
])
<=
1
and
int
(
sympy_version
[
1
])
<=
4
:
def
hash_fun
(
self
):
return
hash
((
self
.
lhs
,
self
.
rhs
))
Assignment
.
__hash__
=
hash_fun
except
Exception
:
pass
...
...
@@ -94,16 +95,19 @@ def assignment_from_stencil(stencil_array, input_field, output_field,
... [0, 6, 0]]
By default 'visual ordering is used - i.e. the stencil is applied as the nested lists are written down
>>> assignment_from_stencil(stencil, f, g, order='visual')
Assignment(g_C, 3*f_W + 6*f_S + 4*f_C + 2*f_N + 5*f_E)
>>> expected_output = Assignment(g[0, 0], 3*f[-1, 0] + 6*f[0, -1] + 4*f[0, 0] + 2*f[0, 1] + 5*f[1, 0])
>>> assignment_from_stencil(stencil, f, g, order='visual') == expected_output
True
'numpy' ordering uses the first coordinate of the stencil array for x offset, second for y offset etc.
>>> assignment_from_stencil(stencil, f, g, order='numpy')
Assignment(g_C, 2*f_W + 3*f_S + 4*f_C + 5*f_N + 6*f_E)
>>> expected_output = Assignment(g[0, 0], 2*f[-1, 0] + 3*f[0, -1] + 4*f[0, 0] + 5*f[0, 1] + 6*f[1, 0])
>>> assignment_from_stencil(stencil, f, g, order='numpy') == expected_output
True
You can also pass field accesses to apply the stencil at an already shifted position:
>>> assignment_from_stencil(stencil, f[1, 0], g[2, 0])
Assignment(g_2E, 3*f_C + 6*f_SE + 4*f_E + 2*f_NE + 5*f_2E)
>>> expected_output = Assignment(g[2, 0], 3*f[0, 0] + 6*f[1, -1] + 4*f[1, 0] + 2*f[1, 1] + 5*f[2, 0])
>>> assignment_from_stencil(stencil, f[1, 0], g[2, 0]) == expected_output
True
"""
from
pystencils.field
import
Field
...
...
pystencils/fd/derivation.py
View file @
ebfd56e1
import
itertools
import
warnings
from
collections
import
defaultdict
import
numpy
as
np
...
...
pystencils/fd/finitedifferences.py
View file @
ebfd56e1
...
...
@@ -21,10 +21,13 @@ def diffusion(scalar, diffusion_coeff, idx=None):
Examples:
>>> f = Field.create_generic('f', spatial_dimensions=2)
>>> d = sp.Symbol("d")
>>> dx = sp.Symbol("dx")
>>> diffusion_term = diffusion(scalar=f, diffusion_coeff=sp.Symbol("d"))
>>> discretization = Discretization2ndOrder()
>>> discretization(diffusion_term)
(f_W*d + f_S*d - 4*f_C*d + f_N*d + f_E*d)/dx**2
>>> expected_output = ((f[-1, 0] + f[0, -1] - 4*f[0, 0] + f[0, 1] + f[1, 0]) * d) / dx**2
>>> sp.simplify(discretization(diffusion_term) - expected_output)
0
"""
if
isinstance
(
scalar
,
Field
):
first_arg
=
scalar
.
center
...
...
@@ -313,8 +316,9 @@ def discretize_center(term, symbols_to_field_dict, dx, dim=3):
>>> term
x*x^Delta^0
>>> f = Field.create_generic('f', spatial_dimensions=3)
>>> discretize_center(term, { x: f }, dx=1, dim=3)
f_C*(-f_W/2 + f_E/2)
>>> expected_output = f[0, 0, 0] * (-f[-1, 0, 0]/2 + f[1, 0, 0]/2)
>>> sp.simplify(discretize_center(term, { x: f }, dx=1, dim=3) - expected_output)
0
"""
substitutions
=
{}
for
symbols
,
field
in
symbols_to_field_dict
.
items
():
...
...
@@ -396,8 +400,10 @@ def discretize_divergence(vector_term, symbols_to_field_dict, dx):
>>> x, dx = sp.symbols("x dx")
>>> grad_x = grad(x, dim=3)
>>> f = Field.create_generic('f', spatial_dimensions=3)
>>> sp.simplify(discretize_divergence(grad_x, {x : f}, dx))
(f_W + f_S + f_B - 6*f_C + f_T + f_N + f_E)/dx**2
>>> expected_output = (f[-1, 0, 0] + f[0, -1, 0] + f[0, 0, -1] -
... 6*f[0, 0, 0] + f[0, 0, 1] + f[0, 1, 0] + f[1, 0, 0])/dx**2
>>> sp.simplify(discretize_divergence(grad_x, {x : f}, dx) - expected_output)
0
"""
dim
=
len
(
vector_term
)
result
=
0
...
...
pystencils_tests/test_fd_derivation.ipynb
View file @
ebfd56e1
...
...
@@ -40,11 +40,11 @@
standard_2d_00
.
get_stencil
().
as_array
()
```
%%%% Output: execute_result
$\displaystyle \left[\begin{matrix}0 & 0 & 0\\1 & -2 & 1\\0 & 0 & 0\end{matrix}\right]$
⎡0 0 0⎤
⎢ ⎥
⎢1 -2 1⎥
⎢ ⎥
...
...
@@ -77,11 +77,11 @@
isotropic_2d_00_res
.
as_array
()
```
%%%% Output: execute_result
$\displaystyle \left[\begin{matrix}\frac{1}{12} & - \frac{1}{6} & \frac{1}{12}\\\frac{5}{6} & - \frac{5}{3} & \frac{5}{6}\\\frac{1}{12} & - \frac{1}{6} & \frac{1}{12}\end{matrix}\right]$
⎡1/12 -1/6 1/12⎤
⎢ ⎥
⎢5/6 -5/3 5/6 ⎥
⎢ ⎥
...
...
@@ -99,34 +99,14 @@
%% Cell type:code id: tags:
```
python
expected_result
=
sp
.
Matrix
([[
1
,
-
2
,
1
],
[
10
,
-
20
,
10
],
[
1
,
-
2
,
1
]])
/
12
assert
expected_result
[:]
==
isotropic_2d_00_res
.
as_array
()
[:]
expected_result
=
sp
.
Array
([[
1
,
-
2
,
1
],
[
10
,
-
20
,
10
],
[
1
,
-
2
,
1
]])
/
12
assert
expected_result
==
isotropic_2d_00_res
.
as_array
()
```
%% Cell type:code id: tags:
```
python
type
(
isotropic_2d_00_res
.
as_array
())
```
%%%% Output: execute_result
sympy.tensor.array.dense_ndim_array.MutableDenseNDimArray
%% Cell type:code id: tags:
```
python
type
(
expected_result
)
```
%%%% Output: execute_result
sympy.matrices.dense.MutableDenseMatrix
%% Cell type:markdown id: tags:
## Isotropic laplacian
%% Cell type:code id: tags:
...
...
@@ -138,23 +118,23 @@
iso_laplacian
```
%%%% Output: execute_result
$\displaystyle \left[\begin{matrix}\frac{1}{6} & \frac{2}{3} & \frac{1}{6}\\\frac{2}{3} & - \frac{10}{3} & \frac{2}{3}\\\frac{1}{6} & \frac{2}{3} & \frac{1}{6}\end{matrix}\right]$
⎡1/6 2/3 1/6⎤
⎢ ⎥
⎢2/3 -10/3 2/3⎥
⎢ ⎥
⎣1/6 2/3 1/6⎦
%% Cell type:code id: tags:
```
python
expected_result
=
sp
.
Matrix
([[
1
,
4
,
1
],
[
4
,
-
20
,
4
],
[
1
,
4
,
1
]])
/
6
assert
iso_laplacian
[:]
==
expected_result
[:]
expected_result
=
sp
.
Array
([[
1
,
4
,
1
],
[
4
,
-
20
,
4
],
[
1
,
4
,
1
]])
/
6
assert
iso_laplacian
==
expected_result
```
%% Cell type:markdown id: tags:
# stencils for staggered fields
...
...
pystencils_tests/test_sympy_optimizations.py
View file @
ebfd56e1
...
...
@@ -11,8 +11,11 @@ def test_sympy_optimizations():
x
,
y
,
z
=
pystencils
.
fields
(
'x, y, z: float32[2d]'
)
# Triggers Sympy's expm1 optimization
# Sympy's expm1 optimization is tedious to use and the behaviour is highly depended on the sympy version. In
# some cases the exp expression has to be encapsulated in brackets or multiplied with 1 or 1.0
# for sympy to work properly ...
assignments
=
pystencils
.
AssignmentCollection
({
x
[
0
,
0
]:
sp
.
exp
(
y
[
0
,
0
])
-
1
x
[
0
,
0
]:
1.0
*
(
sp
.
exp
(
y
[
0
,
0
])
-
1
)
})
assignments
=
optimize_assignments
(
assignments
,
optims_pystencils_cpu
)
...
...
@@ -28,7 +31,7 @@ def test_evaluate_constant_terms():
for
target
in
(
'cpu'
,
'gpu'
):
x
,
y
,
z
=
pystencils
.
fields
(
'x, y, z: float32[2d]'
)
# Triggers Sympy's
expm1
optimization
# Triggers Sympy's
cos
optimization
assignments
=
pystencils
.
AssignmentCollection
({
x
[
0
,
0
]:
-
sp
.
cos
(
1
)
+
y
[
0
,
0
]
})
...
...
@@ -53,8 +56,6 @@ def test_do_not_evaluate_constant_terms():
x
[
0
,
0
]:
-
sp
.
cos
(
1
)
+
y
[
0
,
0
]
})
optimize_assignments
(
assignments
,
optimizations
)
ast
=
pystencils
.
create_kernel
(
assignments
,
target
=
target
)
code
=
pystencils
.
get_code_str
(
ast
)
assert
'cos('
in
code
...
...
Write
Preview
Markdown
is supported
0%
Try again
or
attach a new file
.
Attach a file
Cancel
You are about to add
0
people
to the discussion. Proceed with caution.
Finish editing this message first!
Cancel
Please
register
or
sign in
to comment
