Skip to content
Snippets Groups Projects

Compare revisions

Changes are shown as if the source revision was being merged into the target revision. Learn more about comparing revisions.

Source

Select target project
No results found

Target

Select target project
No results found
Show changes
Commits on Source (8)
Hide whitespace changes
Inline Side-by-side
.gitlab-ci.yml
View file @ 5600b6b6
...@@ -29,7 +29,7 @@ stages: ...@@ -29,7 +29,7 @@ stages:
tests-and-coverage: tests-and-coverage:
stage: pretest stage: pretest
extends: .every-commit extends: .every-commit
image: i10git.cs.fau.de:5005/pycodegen/pycodegen/full image: i10git.cs.fau.de:5005/pycodegen/pycodegen/full:cupy12.3
before_script: before_script:
- pip install -e . - pip install -e .
script: script:
...@@ -62,7 +62,7 @@ tests-and-coverage-with-longrun: ...@@ -62,7 +62,7 @@ tests-and-coverage-with-longrun:
stage: test stage: test
when: manual when: manual
allow_failure: true allow_failure: true
image: i10git.cs.fau.de:5005/pycodegen/pycodegen/full image: i10git.cs.fau.de:5005/pycodegen/pycodegen/full:cupy12.3
before_script: before_script:
- pip install sympy --upgrade - pip install sympy --upgrade
- pip install -e . - pip install -e .
...@@ -145,6 +145,7 @@ ubuntu: ...@@ -145,6 +145,7 @@ ubuntu:
.multiarch_template: .multiarch_template:
stage: test stage: test
extends: .every-commit extends: .every-commit
allow_failure: true
before_script: &multiarch_before_script before_script: &multiarch_before_script
# - pip3 install -v . # - pip3 install -v .
- export PYTHONPATH=src - export PYTHONPATH=src
......
pyproject.toml
View file @ 5600b6b6
...@@ -12,7 +12,7 @@ authors = [ ...@@ -12,7 +12,7 @@ authors = [
] ]
license = { file = "COPYING.txt" } license = { file = "COPYING.txt" }
requires-python = ">=3.10" requires-python = ">=3.10"
dependencies = ["sympy>=1.9,<=1.12.1", "numpy>=1.8.0", "appdirs", "joblib", "pyyaml"] dependencies = ["sympy>=1.9,<=1.12.1", "numpy>=1.8.0", "appdirs", "joblib", "pyyaml", "fasteners"]
classifiers = [ classifiers = [
"Development Status :: 4 - Beta", "Development Status :: 4 - Beta",
"Framework :: Jupyter", "Framework :: Jupyter",
......
src/pystencils/backends/simd_instruction_sets.py
View file @ 5600b6b6
import os import os
import platform import platform
from ctypes import CDLL from ctypes import CDLL, c_int, c_size_t, sizeof, byref
from warnings import warn from warnings import warn
import numpy as np import numpy as np
...@@ -38,7 +38,14 @@ def get_supported_instruction_sets(): ...@@ -38,7 +38,14 @@ def get_supported_instruction_sets():
if 'PYSTENCILS_SIMD' in os.environ: if 'PYSTENCILS_SIMD' in os.environ:
return os.environ['PYSTENCILS_SIMD'].split(',') return os.environ['PYSTENCILS_SIMD'].split(',')
if platform.system() == 'Darwin' and platform.machine() == 'arm64': if platform.system() == 'Darwin' and platform.machine() == 'arm64':
return ['neon'] result = ['neon']
libc = CDLL('/usr/lib/libc.dylib')
value = c_int(0)
size = c_size_t(sizeof(value))
status = libc.sysctlbyname(b"hw.optional.arm.FEAT_SME", byref(value), byref(size), None, 0)
if status == 0 and value.value == 1:
result.insert(0, "sme")
return result
elif platform.system() == 'Windows' and platform.machine() == 'ARM64': elif platform.system() == 'Windows' and platform.machine() == 'ARM64':
return ['neon'] return ['neon']
elif platform.system() == 'Linux' and platform.machine() == 'aarch64': elif platform.system() == 'Linux' and platform.machine() == 'aarch64':
...@@ -59,7 +66,7 @@ def get_supported_instruction_sets(): ...@@ -59,7 +66,7 @@ def get_supported_instruction_sets():
length //= 2 length //= 2
result.append(name) result.append(name)
if hwcap2 & (1 << 23): # HWCAP2_SME if hwcap2 & (1 << 23): # HWCAP2_SME
result.append("sme") result.insert(0, "sme") # prepend to list so it is not automatically chosen as best instruction set
return result return result
elif platform.system() == 'Linux' and platform.machine().startswith('riscv'): elif platform.system() == 'Linux' and platform.machine().startswith('riscv'):
libc = CDLL('libc.so.6') libc = CDLL('libc.so.6')
......
src/pystencils/cpu/cpujit.py
View file @ 5600b6b6
...@@ -57,12 +57,14 @@ import tempfile ...@@ -57,12 +57,14 @@ import tempfile
import textwrap import textwrap
import time import time
import warnings import warnings
import pathlib
import numpy as np import numpy as np
from pystencils import FieldType from pystencils import FieldType
from pystencils.astnodes import LoopOverCoordinate from pystencils.astnodes import LoopOverCoordinate
from pystencils.backends.cbackend import generate_c, get_headers from pystencils.backends.cbackend import generate_c, get_headers
from pystencils.backends.simd_instruction_sets import get_supported_instruction_sets
from pystencils.cpu.msvc_detection import get_environment from pystencils.cpu.msvc_detection import get_environment
from pystencils.include import get_pystencils_include_path from pystencils.include import get_pystencils_include_path
from pystencils.kernel_wrapper import KernelWrapper from pystencils.kernel_wrapper import KernelWrapper
...@@ -122,15 +124,15 @@ def get_configuration_file_path(): ...@@ -122,15 +124,15 @@ def get_configuration_file_path():
# 1) Read path from environment variable if found # 1) Read path from environment variable if found
if 'PYSTENCILS_CONFIG' in os.environ: if 'PYSTENCILS_CONFIG' in os.environ:
return os.environ['PYSTENCILS_CONFIG'], True return os.environ['PYSTENCILS_CONFIG']
# 2) Look in current directory for pystencils.json # 2) Look in current directory for pystencils.json
elif os.path.exists("pystencils.json"): elif os.path.exists("pystencils.json"):
return "pystencils.json", True return "pystencils.json"
# 3) Try ~/.pystencils.json # 3) Try ~/.pystencils.json
elif os.path.exists(config_path_in_home): elif os.path.exists(config_path_in_home):
return config_path_in_home, True return config_path_in_home
else: else:
return config_path_in_home, False return config_path_in_home
def create_folder(path, is_file): def create_folder(path, is_file):
...@@ -172,7 +174,11 @@ def read_config(): ...@@ -172,7 +174,11 @@ def read_config():
('restrict_qualifier', '__restrict__') ('restrict_qualifier', '__restrict__')
]) ])
if platform.machine() == 'arm64': if platform.machine() == 'arm64':
default_compiler_config['flags'] = default_compiler_config['flags'].replace('-march=native ', '') if 'sme' in get_supported_instruction_sets():
flag = '-march=armv8.7-a+sme '
else:
flag = ''
default_compiler_config['flags'] = default_compiler_config['flags'].replace('-march=native ', flag)
for libomp in ['/opt/local/lib/libomp/libomp.dylib', '/usr/local/lib/libomp.dylib', for libomp in ['/opt/local/lib/libomp/libomp.dylib', '/usr/local/lib/libomp.dylib',
'/opt/homebrew/lib/libomp.dylib']: '/opt/homebrew/lib/libomp.dylib']:
if os.path.exists(libomp): if os.path.exists(libomp):
...@@ -190,16 +196,22 @@ def read_config(): ...@@ -190,16 +196,22 @@ def read_config():
default_config = OrderedDict([('compiler', default_compiler_config), default_config = OrderedDict([('compiler', default_compiler_config),
('cache', default_cache_config)]) ('cache', default_cache_config)])
config_path, config_exists = get_configuration_file_path() from fasteners import InterProcessLock
config_path = pathlib.Path(get_configuration_file_path())
config_path.parent.mkdir(parents=True, exist_ok=True)
config = default_config.copy() config = default_config.copy()
if config_exists:
with open(config_path, 'r') as json_config_file: lockfile = config_path.with_suffix(config_path.suffix + ".lock")
loaded_config = json.load(json_config_file) with InterProcessLock(lockfile):
config = recursive_dict_update(config, loaded_config) if config_path.exists():
else: with open(config_path, 'r') as json_config_file:
create_folder(config_path, True) loaded_config = json.load(json_config_file)
with open(config_path, 'w') as f: config = recursive_dict_update(config, loaded_config)
json.dump(config, f, indent=4) else:
with open(config_path, 'w') as f:
json.dump(config, f, indent=4)
if config['cache']['object_cache'] is not False: if config['cache']['object_cache'] is not False:
config['cache']['object_cache'] = os.path.expanduser(config['cache']['object_cache']).format(pid=os.getpid()) config['cache']['object_cache'] = os.path.expanduser(config['cache']['object_cache']).format(pid=os.getpid())
......
src/pystencils/field.py
View file @ 5600b6b6
...@@ -948,24 +948,35 @@ def create_numpy_array_with_layout(shape, layout, alignment=False, byte_offset=0 ...@@ -948,24 +948,35 @@ def create_numpy_array_with_layout(shape, layout, alignment=False, byte_offset=0
def spatial_layout_string_to_tuple(layout_str: str, dim: int) -> Tuple[int, ...]: def spatial_layout_string_to_tuple(layout_str: str, dim: int) -> Tuple[int, ...]:
if layout_str in ('fzyx', 'zyxf'): if dim <= 0:
assert dim <= 3 raise ValueError("Dimensionality must be positive")
return tuple(reversed(range(dim)))
layout_str = layout_str.lower()
if layout_str in ('fzyx', 'f', 'reverse_numpy', 'SoA'): if layout_str in ('fzyx', 'zyxf', 'soa', 'aos'):
if dim > 3:
raise ValueError(f"Invalid spatial dimensionality for layout descriptor {layout_str}: May be at most 3.")
return tuple(reversed(range(dim)))
if layout_str in ('f', 'reverse_numpy'):
return tuple(reversed(range(dim))) return tuple(reversed(range(dim)))
elif layout_str in ('c', 'numpy', 'AoS'): elif layout_str in ('c', 'numpy'):
return tuple(range(dim)) return tuple(range(dim))
raise ValueError("Unknown layout descriptor " + layout_str) raise ValueError("Unknown layout descriptor " + layout_str)
def layout_string_to_tuple(layout_str, dim): def layout_string_to_tuple(layout_str, dim):
if dim <= 0:
raise ValueError("Dimensionality must be positive")
layout_str = layout_str.lower() layout_str = layout_str.lower()
if layout_str == 'fzyx' or layout_str == 'soa': if layout_str == 'fzyx' or layout_str == 'soa':
assert dim <= 4 if dim > 4:
raise ValueError(f"Invalid total dimensionality for layout descriptor {layout_str}: May be at most 4.")
return tuple(reversed(range(dim))) return tuple(reversed(range(dim)))
elif layout_str == 'zyxf' or layout_str == 'aos': elif layout_str == 'zyxf' or layout_str == 'aos':
assert dim <= 4 if dim > 4:
raise ValueError(f"Invalid total dimensionality for layout descriptor {layout_str}: May be at most 4.")
return tuple(reversed(range(dim - 1))) + (dim - 1,) return tuple(reversed(range(dim - 1))) + (dim - 1,)
elif layout_str == 'f' or layout_str == 'reverse_numpy': elif layout_str == 'f' or layout_str == 'reverse_numpy':
return tuple(reversed(range(dim))) return tuple(reversed(range(dim)))
......
src/pystencils/include/philox_rand.h
View file @ 5600b6b6
...@@ -74,13 +74,17 @@ OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. ...@@ -74,13 +74,17 @@ OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
#define QUALIFIERS static __forceinline__ __device__ #define QUALIFIERS static __forceinline__ __device__
#elif defined(__OPENCL_VERSION__) #elif defined(__OPENCL_VERSION__)
#define QUALIFIERS static inline #define QUALIFIERS static inline
#elif defined(__ARM_FEATURE_SME)
#define QUALIFIERS __attribute__((arm_streaming_compatible))
#else #else
#define QUALIFIERS inline #define QUALIFIERS inline
#include "myintrin.h" #include "myintrin.h"
#endif #endif
#if defined(__ARM_FEATURE_SME)
#define SVE_QUALIFIERS __attribute__((arm_streaming_compatible)) QUALIFIERS
#else
#define SVE_QUALIFIERS QUALIFIERS
#endif
#define PHILOX_W32_0 (0x9E3779B9) #define PHILOX_W32_0 (0x9E3779B9)
#define PHILOX_W32_1 (0xBB67AE85) #define PHILOX_W32_1 (0xBB67AE85)
#define PHILOX_M4x32_0 (0xD2511F53) #define PHILOX_M4x32_0 (0xD2511F53)
...@@ -749,7 +753,7 @@ QUALIFIERS void philox_double2(uint32 ctr0, int32x4_t ctr1, uint32 ctr2, uint32 ...@@ -749,7 +753,7 @@ QUALIFIERS void philox_double2(uint32 ctr0, int32x4_t ctr1, uint32 ctr2, uint32
#if defined(__ARM_FEATURE_SVE) || defined(__ARM_FEATURE_SME) #if defined(__ARM_FEATURE_SVE) || defined(__ARM_FEATURE_SME)
QUALIFIERS void _philox4x32round(svuint32x4_t & ctr, svuint32x2_t & key) SVE_QUALIFIERS void _philox4x32round(svuint32x4_t & ctr, svuint32x2_t & key)
{ {
svuint32_t lo0 = svmul_u32_x(svptrue_b32(), svget4_u32(ctr, 0), svdup_u32(PHILOX_M4x32_0)); svuint32_t lo0 = svmul_u32_x(svptrue_b32(), svget4_u32(ctr, 0), svdup_u32(PHILOX_M4x32_0));
svuint32_t hi0 = svmulh_u32_x(svptrue_b32(), svget4_u32(ctr, 0), svdup_u32(PHILOX_M4x32_0)); svuint32_t hi0 = svmulh_u32_x(svptrue_b32(), svget4_u32(ctr, 0), svdup_u32(PHILOX_M4x32_0));
...@@ -762,14 +766,14 @@ QUALIFIERS void _philox4x32round(svuint32x4_t & ctr, svuint32x2_t & key) ...@@ -762,14 +766,14 @@ QUALIFIERS void _philox4x32round(svuint32x4_t & ctr, svuint32x2_t & key)
ctr = svset4_u32(ctr, 3, lo0); ctr = svset4_u32(ctr, 3, lo0);
} }
QUALIFIERS void _philox4x32bumpkey(svuint32x2_t & key) SVE_QUALIFIERS void _philox4x32bumpkey(svuint32x2_t & key)
{ {
key = svset2_u32(key, 0, svadd_u32_x(svptrue_b32(), svget2_u32(key, 0), svdup_u32(PHILOX_W32_0))); key = svset2_u32(key, 0, svadd_u32_x(svptrue_b32(), svget2_u32(key, 0), svdup_u32(PHILOX_W32_0)));
key = svset2_u32(key, 1, svadd_u32_x(svptrue_b32(), svget2_u32(key, 1), svdup_u32(PHILOX_W32_1))); key = svset2_u32(key, 1, svadd_u32_x(svptrue_b32(), svget2_u32(key, 1), svdup_u32(PHILOX_W32_1)));
} }
template<bool high> template<bool high>
QUALIFIERS svfloat64_t _uniform_double_hq(svuint32_t x, svuint32_t y) SVE_QUALIFIERS svfloat64_t _uniform_double_hq(svuint32_t x, svuint32_t y)
{ {
// convert 32 to 64 bit // convert 32 to 64 bit
if (high) if (high)
...@@ -796,9 +800,9 @@ QUALIFIERS svfloat64_t _uniform_double_hq(svuint32_t x, svuint32_t y) ...@@ -796,9 +800,9 @@ QUALIFIERS svfloat64_t _uniform_double_hq(svuint32_t x, svuint32_t y)
} }
QUALIFIERS void philox_float4(svuint32_t ctr0, svuint32_t ctr1, svuint32_t ctr2, svuint32_t ctr3, SVE_QUALIFIERS void philox_float4(svuint32_t ctr0, svuint32_t ctr1, svuint32_t ctr2, svuint32_t ctr3,
uint32 key0, uint32 key1, uint32 key0, uint32 key1,
svfloat32_st & rnd1, svfloat32_st & rnd2, svfloat32_st & rnd3, svfloat32_st & rnd4) svfloat32_st & rnd1, svfloat32_st & rnd2, svfloat32_st & rnd3, svfloat32_st & rnd4)
{ {
svuint32x2_t key = svcreate2_u32(svdup_u32(key0), svdup_u32(key1)); svuint32x2_t key = svcreate2_u32(svdup_u32(key0), svdup_u32(key1));
svuint32x4_t ctr = svcreate4_u32(ctr0, ctr1, ctr2, ctr3); svuint32x4_t ctr = svcreate4_u32(ctr0, ctr1, ctr2, ctr3);
...@@ -826,9 +830,9 @@ QUALIFIERS void philox_float4(svuint32_t ctr0, svuint32_t ctr1, svuint32_t ctr2, ...@@ -826,9 +830,9 @@ QUALIFIERS void philox_float4(svuint32_t ctr0, svuint32_t ctr1, svuint32_t ctr2,
} }
QUALIFIERS void philox_double2(svuint32_t ctr0, svuint32_t ctr1, svuint32_t ctr2, svuint32_t ctr3, SVE_QUALIFIERS void philox_double2(svuint32_t ctr0, svuint32_t ctr1, svuint32_t ctr2, svuint32_t ctr3,
uint32 key0, uint32 key1, uint32 key0, uint32 key1,
svfloat64_st & rnd1lo, svfloat64_st & rnd1hi, svfloat64_st & rnd2lo, svfloat64_st & rnd2hi) svfloat64_st & rnd1lo, svfloat64_st & rnd1hi, svfloat64_st & rnd2lo, svfloat64_st & rnd2hi)
{ {
svuint32x2_t key = svcreate2_u32(svdup_u32(key0), svdup_u32(key1)); svuint32x2_t key = svcreate2_u32(svdup_u32(key0), svdup_u32(key1));
svuint32x4_t ctr = svcreate4_u32(ctr0, ctr1, ctr2, ctr3); svuint32x4_t ctr = svcreate4_u32(ctr0, ctr1, ctr2, ctr3);
...@@ -849,9 +853,9 @@ QUALIFIERS void philox_double2(svuint32_t ctr0, svuint32_t ctr1, svuint32_t ctr2 ...@@ -849,9 +853,9 @@ QUALIFIERS void philox_double2(svuint32_t ctr0, svuint32_t ctr1, svuint32_t ctr2
rnd2hi = _uniform_double_hq<true>(svget4_u32(ctr, 2), svget4_u32(ctr, 3)); rnd2hi = _uniform_double_hq<true>(svget4_u32(ctr, 2), svget4_u32(ctr, 3));
} }
QUALIFIERS void philox_float4(uint32 ctr0, svuint32_t ctr1, uint32 ctr2, uint32 ctr3, SVE_QUALIFIERS void philox_float4(uint32 ctr0, svuint32_t ctr1, uint32 ctr2, uint32 ctr3,
uint32 key0, uint32 key1, uint32 key0, uint32 key1,
svfloat32_st & rnd1, svfloat32_st & rnd2, svfloat32_st & rnd3, svfloat32_st & rnd4) svfloat32_st & rnd1, svfloat32_st & rnd2, svfloat32_st & rnd3, svfloat32_st & rnd4)
{ {
svuint32_t ctr0v = svdup_u32(ctr0); svuint32_t ctr0v = svdup_u32(ctr0);
svuint32_t ctr2v = svdup_u32(ctr2); svuint32_t ctr2v = svdup_u32(ctr2);
...@@ -860,16 +864,16 @@ QUALIFIERS void philox_float4(uint32 ctr0, svuint32_t ctr1, uint32 ctr2, uint32 ...@@ -860,16 +864,16 @@ QUALIFIERS void philox_float4(uint32 ctr0, svuint32_t ctr1, uint32 ctr2, uint32
philox_float4(ctr0v, ctr1, ctr2v, ctr3v, key0, key1, rnd1, rnd2, rnd3, rnd4); philox_float4(ctr0v, ctr1, ctr2v, ctr3v, key0, key1, rnd1, rnd2, rnd3, rnd4);
} }
QUALIFIERS void philox_float4(uint32 ctr0, svint32_t ctr1, uint32 ctr2, uint32 ctr3, SVE_QUALIFIERS void philox_float4(uint32 ctr0, svint32_t ctr1, uint32 ctr2, uint32 ctr3,
uint32 key0, uint32 key1, uint32 key0, uint32 key1,
svfloat32_st & rnd1, svfloat32_st & rnd2, svfloat32_st & rnd3, svfloat32_st & rnd4) svfloat32_st & rnd1, svfloat32_st & rnd2, svfloat32_st & rnd3, svfloat32_st & rnd4)
{ {
philox_float4(ctr0, svreinterpret_u32_s32(ctr1), ctr2, ctr3, key0, key1, rnd1, rnd2, rnd3, rnd4); philox_float4(ctr0, svreinterpret_u32_s32(ctr1), ctr2, ctr3, key0, key1, rnd1, rnd2, rnd3, rnd4);
} }
QUALIFIERS void philox_double2(uint32 ctr0, svuint32_t ctr1, uint32 ctr2, uint32 ctr3, SVE_QUALIFIERS void philox_double2(uint32 ctr0, svuint32_t ctr1, uint32 ctr2, uint32 ctr3,
uint32 key0, uint32 key1, uint32 key0, uint32 key1,
svfloat64_st & rnd1lo, svfloat64_st & rnd1hi, svfloat64_st & rnd2lo, svfloat64_st & rnd2hi) svfloat64_st & rnd1lo, svfloat64_st & rnd1hi, svfloat64_st & rnd2lo, svfloat64_st & rnd2hi)
{ {
svuint32_t ctr0v = svdup_u32(ctr0); svuint32_t ctr0v = svdup_u32(ctr0);
svuint32_t ctr2v = svdup_u32(ctr2); svuint32_t ctr2v = svdup_u32(ctr2);
...@@ -878,9 +882,9 @@ QUALIFIERS void philox_double2(uint32 ctr0, svuint32_t ctr1, uint32 ctr2, uint32 ...@@ -878,9 +882,9 @@ QUALIFIERS void philox_double2(uint32 ctr0, svuint32_t ctr1, uint32 ctr2, uint32
philox_double2(ctr0v, ctr1, ctr2v, ctr3v, key0, key1, rnd1lo, rnd1hi, rnd2lo, rnd2hi); philox_double2(ctr0v, ctr1, ctr2v, ctr3v, key0, key1, rnd1lo, rnd1hi, rnd2lo, rnd2hi);
} }
QUALIFIERS void philox_double2(uint32 ctr0, svuint32_t ctr1, uint32 ctr2, uint32 ctr3, SVE_QUALIFIERS void philox_double2(uint32 ctr0, svuint32_t ctr1, uint32 ctr2, uint32 ctr3,
uint32 key0, uint32 key1, uint32 key0, uint32 key1,
svfloat64_st & rnd1, svfloat64_st & rnd2) svfloat64_st & rnd1, svfloat64_st & rnd2)
{ {
svuint32_t ctr0v = svdup_u32(ctr0); svuint32_t ctr0v = svdup_u32(ctr0);
svuint32_t ctr2v = svdup_u32(ctr2); svuint32_t ctr2v = svdup_u32(ctr2);
...@@ -890,9 +894,9 @@ QUALIFIERS void philox_double2(uint32 ctr0, svuint32_t ctr1, uint32 ctr2, uint32 ...@@ -890,9 +894,9 @@ QUALIFIERS void philox_double2(uint32 ctr0, svuint32_t ctr1, uint32 ctr2, uint32
philox_double2(ctr0v, ctr1, ctr2v, ctr3v, key0, key1, rnd1, ignore, rnd2, ignore); philox_double2(ctr0v, ctr1, ctr2v, ctr3v, key0, key1, rnd1, ignore, rnd2, ignore);
} }
QUALIFIERS void philox_double2(uint32 ctr0, svint32_t ctr1, uint32 ctr2, uint32 ctr3, SVE_QUALIFIERS void philox_double2(uint32 ctr0, svint32_t ctr1, uint32 ctr2, uint32 ctr3,
uint32 key0, uint32 key1, uint32 key0, uint32 key1,
svfloat64_st & rnd1, svfloat64_st & rnd2) svfloat64_st & rnd1, svfloat64_st & rnd2)
{ {
philox_double2(ctr0, svreinterpret_u32_s32(ctr1), ctr2, ctr3, key0, key1, rnd1, rnd2); philox_double2(ctr0, svreinterpret_u32_s32(ctr1), ctr2, ctr3, key0, key1, rnd1, rnd2);
} }
......
tests/test_field.py
View file @ 5600b6b6
...@@ -5,21 +5,33 @@ import sympy as sp ...@@ -5,21 +5,33 @@ import sympy as sp
import pystencils as ps import pystencils as ps
from pystencils import TypedSymbol from pystencils import TypedSymbol
from pystencils.typing import create_type from pystencils.typing import create_type
from pystencils.field import Field, FieldType, layout_string_to_tuple from pystencils.field import Field, FieldType, layout_string_to_tuple, spatial_layout_string_to_tuple
def test_field_basic(): def test_field_basic():
f = Field.create_generic('f', spatial_dimensions=2) f = Field.create_generic("f", spatial_dimensions=2)
assert FieldType.is_generic(f) assert FieldType.is_generic(f)
assert f['E'] == f[1, 0] assert f["E"] == f[1, 0]
assert f['N'] == f[0, 1] assert f["N"] == f[0, 1]
assert '_' in f.center._latex('dummy') assert "_" in f.center._latex("dummy")
assert f.index_to_physical(index_coordinates=sp.Matrix([0, 0]), staggered=False)[0] == 0 assert (
assert f.index_to_physical(index_coordinates=sp.Matrix([0, 0]), staggered=False)[1] == 0 f.index_to_physical(index_coordinates=sp.Matrix([0, 0]), staggered=False)[0]
== 0
assert f.physical_to_index(physical_coordinates=sp.Matrix([0, 0]), staggered=False)[0] == 0 )
assert f.physical_to_index(physical_coordinates=sp.Matrix([0, 0]), staggered=False)[1] == 0 assert (
f.index_to_physical(index_coordinates=sp.Matrix([0, 0]), staggered=False)[1]
== 0
)
assert (
f.physical_to_index(physical_coordinates=sp.Matrix([0, 0]), staggered=False)[0]
== 0
)
assert (
f.physical_to_index(physical_coordinates=sp.Matrix([0, 0]), staggered=False)[1]
== 0
)
f1 = f.new_field_with_different_name("f1") f1 = f.new_field_with_different_name("f1")
assert f1.ndim == f.ndim assert f1.ndim == f.ndim
...@@ -28,7 +40,7 @@ def test_field_basic(): ...@@ -28,7 +40,7 @@ def test_field_basic():
fixed = ps.fields("f(5, 5) : double[20, 20]") fixed = ps.fields("f(5, 5) : double[20, 20]")
assert fixed.neighbor_vector((1, 1)).shape == (5, 5) assert fixed.neighbor_vector((1, 1)).shape == (5, 5)
f = Field.create_fixed_size('f', (10, 10), strides=(80, 8), dtype=np.float64) f = Field.create_fixed_size("f", (10, 10), strides=(80, 8), dtype=np.float64)
assert f.spatial_strides == (10, 1) assert f.spatial_strides == (10, 1)
assert f.index_strides == () assert f.index_strides == ()
assert f.center_vector == sp.Matrix([f.center]) assert f.center_vector == sp.Matrix([f.center])
...@@ -37,20 +49,21 @@ def test_field_basic(): ...@@ -37,20 +49,21 @@ def test_field_basic():
assert f1.ndim == f.ndim assert f1.ndim == f.ndim
assert f1.values_per_cell() == f.values_per_cell() assert f1.values_per_cell() == f.values_per_cell()
f = Field.create_fixed_size('f', (8, 8, 2, 2), index_dimensions=2) f = Field.create_fixed_size("f", (8, 8, 2, 2), index_dimensions=2)
assert f.center_vector == sp.Matrix([[f(0, 0), f(0, 1)], assert f.center_vector == sp.Matrix([[f(0, 0), f(0, 1)], [f(1, 0), f(1, 1)]])
[f(1, 0), f(1, 1)]])
field_access = f[1, 1] field_access = f[1, 1]
assert field_access.nr_of_coordinates == 2 assert field_access.nr_of_coordinates == 2
assert field_access.offset_name == 'NE' assert field_access.offset_name == "NE"
neighbor = field_access.neighbor(coord_id=0, offset=-2) neighbor = field_access.neighbor(coord_id=0, offset=-2)
assert neighbor.offsets == (-1, 1) assert neighbor.offsets == (-1, 1)
assert '_' in neighbor._latex('dummy') assert "_" in neighbor._latex("dummy")
f = Field.create_fixed_size('f', (8, 8, 2, 2, 2), index_dimensions=3) f = Field.create_fixed_size("f", (8, 8, 2, 2, 2), index_dimensions=3)
assert f.center_vector == sp.Array([[[f(i, j, k) for k in range(2)] for j in range(2)] for i in range(2)]) assert f.center_vector == sp.Array(
[[[f(i, j, k) for k in range(2)] for j in range(2)] for i in range(2)]
)
f = Field.create_generic('f', spatial_dimensions=5, index_dimensions=2) f = Field.create_generic("f", spatial_dimensions=5, index_dimensions=2)
field_access = f[1, -1, 2, -3, 0](1, 0) field_access = f[1, -1, 2, -3, 0](1, 0)
assert field_access.offsets == (1, -1, 2, -3, 0) assert field_access.offsets == (1, -1, 2, -3, 0)
assert field_access.index == (1, 0) assert field_access.index == (1, 0)
...@@ -60,61 +73,71 @@ def test_error_handling(): ...@@ -60,61 +73,71 @@ def test_error_handling():
struct_dtype = np.dtype([('a', np.int32), ('b', np.float64), ('c', np.uint32)]) struct_dtype = np.dtype([('a', np.int32), ('b', np.float64), ('c', np.uint32)])
Field.create_generic('f', spatial_dimensions=2, index_dimensions=0, dtype=struct_dtype) Field.create_generic('f', spatial_dimensions=2, index_dimensions=0, dtype=struct_dtype)
with pytest.raises(ValueError) as e: with pytest.raises(ValueError) as e:
Field.create_generic('f', spatial_dimensions=2, index_dimensions=1, dtype=struct_dtype) Field.create_generic(
assert 'index dimension' in str(e.value) "f", spatial_dimensions=2, index_dimensions=1, dtype=struct_dtype
)
assert "index dimension" in str(e.value)
arr = np.array([[[(1,)*3, (2,)*3, (3,)*3]]*2], dtype=struct_dtype) arr = np.array([[[(1,) * 3, (2,) * 3, (3,) * 3]] * 2], dtype=struct_dtype)
Field.create_from_numpy_array('f', arr, index_dimensions=0) Field.create_from_numpy_array("f", arr, index_dimensions=0)
with pytest.raises(ValueError) as e: with pytest.raises(ValueError) as e:
Field.create_from_numpy_array('f', arr, index_dimensions=1) Field.create_from_numpy_array("f", arr, index_dimensions=1)
assert 'Structured arrays' in str(e.value) assert "Structured arrays" in str(e.value)
arr = np.zeros([3, 3, 3]) arr = np.zeros([3, 3, 3])
Field.create_from_numpy_array('f', arr, index_dimensions=2) Field.create_from_numpy_array("f", arr, index_dimensions=2)
with pytest.raises(ValueError) as e: with pytest.raises(ValueError) as e:
Field.create_from_numpy_array('f', arr, index_dimensions=3) Field.create_from_numpy_array("f", arr, index_dimensions=3)
assert 'Too many' in str(e.value) assert "Too many" in str(e.value)
Field.create_fixed_size('f', (3, 2, 4), index_dimensions=0, dtype=struct_dtype, layout='reverse_numpy') Field.create_fixed_size(
"f", (3, 2, 4), index_dimensions=0, dtype=struct_dtype, layout="reverse_numpy"
)
with pytest.raises(ValueError) as e: with pytest.raises(ValueError) as e:
Field.create_fixed_size('f', (3, 2, 4), index_dimensions=1, dtype=struct_dtype, layout='reverse_numpy') Field.create_fixed_size(
assert 'Structured arrays' in str(e.value) "f",
(3, 2, 4),
f = Field.create_fixed_size('f', (10, 10)) index_dimensions=1,
dtype=struct_dtype,
layout="reverse_numpy",
)
assert "Structured arrays" in str(e.value)
f = Field.create_fixed_size("f", (10, 10))
with pytest.raises(ValueError) as e: with pytest.raises(ValueError) as e:
f[1] f[1]
assert 'Wrong number of spatial indices' in str(e.value) assert "Wrong number of spatial indices" in str(e.value)
f = Field.create_generic('f', spatial_dimensions=2, index_shape=(3,)) f = Field.create_generic("f", spatial_dimensions=2, index_shape=(3,))
with pytest.raises(ValueError) as e: with pytest.raises(ValueError) as e:
f(3) f(3)
assert 'out of bounds' in str(e.value) assert "out of bounds" in str(e.value)
f = Field.create_fixed_size('f', (10, 10, 3, 4), index_dimensions=2) f = Field.create_fixed_size("f", (10, 10, 3, 4), index_dimensions=2)
with pytest.raises(ValueError) as e: with pytest.raises(ValueError) as e:
f(3, 0) f(3, 0)
assert 'out of bounds' in str(e.value) assert "out of bounds" in str(e.value)
with pytest.raises(ValueError) as e: with pytest.raises(ValueError) as e:
f(1, 0)(1, 0) f(1, 0)(1, 0)
assert 'Indexing an already indexed' in str(e.value) assert "Indexing an already indexed" in str(e.value)
with pytest.raises(ValueError) as e: with pytest.raises(ValueError) as e:
f(1) f(1)
assert 'Wrong number of indices' in str(e.value) assert "Wrong number of indices" in str(e.value)
with pytest.raises(ValueError) as e: with pytest.raises(ValueError) as e:
Field.create_generic('f', spatial_dimensions=2, layout='wrong') Field.create_generic("f", spatial_dimensions=2, layout="wrong")
assert 'Unknown layout descriptor' in str(e.value) assert "Unknown layout descriptor" in str(e.value)
assert layout_string_to_tuple('fzyx', dim=4) == (3, 2, 1, 0) assert layout_string_to_tuple("fzyx", dim=4) == (3, 2, 1, 0)
with pytest.raises(ValueError) as e: with pytest.raises(ValueError) as e:
layout_string_to_tuple('wrong', dim=4) layout_string_to_tuple("wrong", dim=4)
assert 'Unknown layout descriptor' in str(e.value) assert "Unknown layout descriptor" in str(e.value)
def test_decorator_scoping(): def test_decorator_scoping():
dst = ps.fields('dst : double[2D]') dst = ps.fields("dst : double[2D]")
def f1(): def f1():
a = sp.Symbol("a") a = sp.Symbol("a")
...@@ -134,7 +157,7 @@ def test_decorator_scoping(): ...@@ -134,7 +157,7 @@ def test_decorator_scoping():
def test_string_creation(): def test_string_creation():
x, y, z = ps.fields(' x(4), y(3,5) z : double[ 3, 47]') x, y, z = ps.fields(" x(4), y(3,5) z : double[ 3, 47]")
assert x.index_shape == (4,) assert x.index_shape == (4,)
assert y.index_shape == (3, 5) assert y.index_shape == (3, 5)
assert z.spatial_shape == (3, 47) assert z.spatial_shape == (3, 47)
...@@ -142,19 +165,85 @@ def test_string_creation(): ...@@ -142,19 +165,85 @@ def test_string_creation():
def test_itemsize(): def test_itemsize():
x = ps.fields('x: float32[1d]') x = ps.fields("x: float32[1d]")
y = ps.fields('y: float64[2d]') y = ps.fields("y: float64[2d]")
i = ps.fields('i: int16[1d]') i = ps.fields("i: int16[1d]")
assert x.itemsize == 4 assert x.itemsize == 4
assert y.itemsize == 8 assert y.itemsize == 8
assert i.itemsize == 2 assert i.itemsize == 2
def test_spatial_memory_layout_descriptors():
assert (
spatial_layout_string_to_tuple("AoS", 3)
== spatial_layout_string_to_tuple("aos", 3)
== spatial_layout_string_to_tuple("ZYXF", 3)
== spatial_layout_string_to_tuple("zyxf", 3)
== (2, 1, 0)
)
assert (
spatial_layout_string_to_tuple("SoA", 3)
== spatial_layout_string_to_tuple("soa", 3)
== spatial_layout_string_to_tuple("FZYX", 3)
== spatial_layout_string_to_tuple("fzyx", 3)
== spatial_layout_string_to_tuple("f", 3)
== spatial_layout_string_to_tuple("F", 3)
== (2, 1, 0)
)
assert (
spatial_layout_string_to_tuple("c", 3)
== spatial_layout_string_to_tuple("C", 3)
== (0, 1, 2)
)
assert spatial_layout_string_to_tuple("C", 5) == (0, 1, 2, 3, 4)
with pytest.raises(ValueError):
spatial_layout_string_to_tuple("aos", -1)
with pytest.raises(ValueError):
spatial_layout_string_to_tuple("aos", 4)
def test_memory_layout_descriptors():
assert (
layout_string_to_tuple("AoS", 4)
== layout_string_to_tuple("aos", 4)
== layout_string_to_tuple("ZYXF", 4)
== layout_string_to_tuple("zyxf", 4)
== (2, 1, 0, 3)
)
assert (
layout_string_to_tuple("SoA", 4)
== layout_string_to_tuple("soa", 4)
== layout_string_to_tuple("FZYX", 4)
== layout_string_to_tuple("fzyx", 4)
== layout_string_to_tuple("f", 4)
== layout_string_to_tuple("F", 4)
== (3, 2, 1, 0)
)
assert (
layout_string_to_tuple("c", 4)
== layout_string_to_tuple("C", 4)
== (0, 1, 2, 3)
)
assert layout_string_to_tuple("C", 5) == (0, 1, 2, 3, 4)
with pytest.raises(ValueError):
layout_string_to_tuple("aos", -1)
with pytest.raises(ValueError):
layout_string_to_tuple("aos", 5)
def test_staggered(): def test_staggered():
# D2Q5 # D2Q5
j1, j2, j3 = ps.fields('j1(2), j2(2,2), j3(2,2,2) : double[2D]', field_type=FieldType.STAGGERED) j1, j2, j3 = ps.fields(
"j1(2), j2(2,2), j3(2,2,2) : double[2D]", field_type=FieldType.STAGGERED
)
assert j1[0, 1](1) == j1.staggered_access((0, sp.Rational(1, 2))) assert j1[0, 1](1) == j1.staggered_access((0, sp.Rational(1, 2)))
assert j1[0, 1](1) == j1.staggered_access(np.array((0, sp.Rational(1, 2)))) assert j1[0, 1](1) == j1.staggered_access(np.array((0, sp.Rational(1, 2))))
...@@ -163,7 +252,7 @@ def test_staggered(): ...@@ -163,7 +252,7 @@ def test_staggered():
assert j1[0, 1](1) == j1.staggered_access("N") assert j1[0, 1](1) == j1.staggered_access("N")
assert j1[0, 0](1) == j1.staggered_access("S") assert j1[0, 0](1) == j1.staggered_access("S")
assert j1.staggered_vector_access("N") == sp.Matrix([j1.staggered_access("N")]) assert j1.staggered_vector_access("N") == sp.Matrix([j1.staggered_access("N")])
assert j1.staggered_stencil_name == 'D2Q5' assert j1.staggered_stencil_name == "D2Q5"
assert j1.physical_coordinates[0] == TypedSymbol("ctr_0", create_type("int"), nonnegative=True) assert j1.physical_coordinates[0] == TypedSymbol("ctr_0", create_type("int"), nonnegative=True)
assert j1.physical_coordinates[1] == TypedSymbol("ctr_1", create_type("int"), nonnegative=True) assert j1.physical_coordinates[1] == TypedSymbol("ctr_1", create_type("int"), nonnegative=True)
...@@ -176,28 +265,40 @@ def test_staggered(): ...@@ -176,28 +265,40 @@ def test_staggered():
assert j2[0, 1](1, 1) == j2.staggered_access((0, sp.Rational(1, 2)), 1) assert j2[0, 1](1, 1) == j2.staggered_access((0, sp.Rational(1, 2)), 1)
assert j2[0, 1](1, 1) == j2.staggered_access("N", 1) assert j2[0, 1](1, 1) == j2.staggered_access("N", 1)
assert j2.staggered_vector_access("N") == sp.Matrix([j2.staggered_access("N", 0), j2.staggered_access("N", 1)]) assert j2.staggered_vector_access("N") == sp.Matrix(
[j2.staggered_access("N", 0), j2.staggered_access("N", 1)]
)
assert j3[0, 1](1, 1, 1) == j3.staggered_access((0, sp.Rational(1, 2)), (1, 1)) assert j3[0, 1](1, 1, 1) == j3.staggered_access((0, sp.Rational(1, 2)), (1, 1))
assert j3[0, 1](1, 1, 1) == j3.staggered_access("N", (1, 1)) assert j3[0, 1](1, 1, 1) == j3.staggered_access("N", (1, 1))
assert j3.staggered_vector_access("N") == sp.Matrix([[j3.staggered_access("N", (i, j)) assert j3.staggered_vector_access("N") == sp.Matrix(
for j in range(2)] for i in range(2)]) [[j3.staggered_access("N", (i, j)) for j in range(2)] for i in range(2)]
)
# D2Q9 # D2Q9
k1, k2 = ps.fields('k1(4), k2(2) : double[2D]', field_type=FieldType.STAGGERED) k1, k2 = ps.fields("k1(4), k2(2) : double[2D]", field_type=FieldType.STAGGERED)
assert k1[1, 1](2) == k1.staggered_access("NE") assert k1[1, 1](2) == k1.staggered_access("NE")
assert k1[0, 0](2) == k1.staggered_access("SW") assert k1[0, 0](2) == k1.staggered_access("SW")
assert k1[0, 0](3) == k1.staggered_access("NW") assert k1[0, 0](3) == k1.staggered_access("NW")
a = k1.staggered_access("NE") a = k1.staggered_access("NE")
assert a._staggered_offset(a.offsets, a.index[0]) == [sp.Rational(1, 2), sp.Rational(1, 2)] assert a._staggered_offset(a.offsets, a.index[0]) == [
sp.Rational(1, 2),
sp.Rational(1, 2),
]
a = k1.staggered_access("SW") a = k1.staggered_access("SW")
assert a._staggered_offset(a.offsets, a.index[0]) == [sp.Rational(-1, 2), sp.Rational(-1, 2)] assert a._staggered_offset(a.offsets, a.index[0]) == [
sp.Rational(-1, 2),
sp.Rational(-1, 2),
]
a = k1.staggered_access("NW") a = k1.staggered_access("NW")
assert a._staggered_offset(a.offsets, a.index[0]) == [sp.Rational(-1, 2), sp.Rational(1, 2)] assert a._staggered_offset(a.offsets, a.index[0]) == [
sp.Rational(-1, 2),
sp.Rational(1, 2),
]
# sign reversed when using as flux field # sign reversed when using as flux field
r = ps.fields('r(2) : double[2D]', field_type=FieldType.STAGGERED_FLUX) r = ps.fields("r(2) : double[2D]", field_type=FieldType.STAGGERED_FLUX)
assert r[0, 0](0) == r.staggered_access("W") assert r[0, 0](0) == r.staggered_access("W")
assert -r[1, 0](0) == r.staggered_access("E") assert -r[1, 0](0) == r.staggered_access("E")