pystencils issueshttps://i10git.cs.fau.de/pycodegen/pystencils/-/issues2024-03-27T18:22:08+01:00https://i10git.cs.fau.de/pycodegen/pystencils/-/issues/48Python 3.102024-03-27T18:22:08+01:00Jan HönigPython 3.10We switch to `python3.10` when `conda` and the latest Ubuntu LTS version support that.
Apply pattern matching at important functions.We switch to `python3.10` when `conda` and the latest Ubuntu LTS version support that.
Apply pattern matching at important functions.Markus HolzerMarkus Holzerhttps://i10git.cs.fau.de/pycodegen/pystencils/-/issues/82C++ Code Printer2024-03-12T11:14:20+01:00Frederik HennigC++ Code Printerhttps://i10git.cs.fau.de/pycodegen/pystencils/-/issues/81Improved Expression Manipulation
2024-03-10T17:59:43+01:00Frederik HennigImproved Expression Manipulation
## Freeze
`FreezeExpressions` is currently still quite naive; it should be more intelligent and take SymPy's peculiarities into account:
- Expand `Pow` to multiplications where sensible
- Detect subtractions (e.g. a + (-1 * b)) and m...## Freeze
`FreezeExpressions` is currently still quite naive; it should be more intelligent and take SymPy's peculiarities into account:
- Expand `Pow` to multiplications where sensible
- Detect subtractions (e.g. a + (-1 * b)) and map them to `PsSub`
- [x] Extend `FreezeExpressions`
## Constant Folding
The new AST is completely static, but constant folding will be necessary as a general optimization and during many transformations. Integers can always be folded; we should be more conservative about floating-point constants.
- [x] Implement a constant folding passhttps://i10git.cs.fau.de/pycodegen/pystencils/-/issues/80kernelcreation: Translation from SymPy-frontend2024-02-28T10:40:00+01:00Frederik Hennigkernelcreation: Translation from SymPy-frontendThe machinery behind `create_kernel` mapping a kernel expressed in the SymPy-based frontend to the backend needs to be rebuilt. The following primary operations need to be reimplemented or adapted:
- **Freeze AST:** The input Assignmen...The machinery behind `create_kernel` mapping a kernel expressed in the SymPy-based frontend to the backend needs to be rebuilt. The following primary operations need to be reimplemented or adapted:
- **Freeze AST:** The input AssignmentCollection (or whatever container is used) is transformed to the pymbolic-based representation. Output of the freeze should be a pure but still untyped backend-AST containing pymbolic expressions.
- **Type Inference:** The type inference algorithm should compute all missing expression types on the frozen AST. All `Variables` must be replaced by `PsTypedVariable`, all constants replaced by `PsTypedConstant`. After this, the AST is fully typed and must remain so.
- **Loop Nest/ GPU indexing:** Depending on the type of kernel, the AST is wrapped into a loop nest or prepended with a GPU index calculation
- **Optimization**: Various optimizations can be applied. These include:
- Constant folding
- Loop-invariant code motion
- OpenMP parallelization
- Loop transformations (tiling, blocking, splitting)
- Vectorization
After the types are figured out, any operation on the AST must maintain the invariant that the AST is fully and correctly typed.
## Data Structure Correspondence
The following correspondences and translations between front-end and back-end structures exist and might have to be tracked by the translator:
- `sp.Symbol` --(freeze)--> `pb.Variable` --(typing)--> `PsTypedVariable`
- `ps.TypedSymbol` --(freeze)--> `PsTypedVariable`
- `ps.Field` --(freeze)--> `PsLinearizedArray` + `PsTypedVariable`s for pointers and indexing info
## Concept, Components, and Progress
### Iteration Spaces
The kernel's iteration space is modelled using the `IterationSpace` class.
Two variants exist: full and sparse iteration spaces.
The iteration space controls the translation of field accesses to array accesses
and will later be materialized to an index source (e.g. a loop nest) by the platform-specific code generator.
`FullIterationSpace`:
- [x] Creation from ghost layers
- [x] Creation from iteration slice
`SparseIterationSpace`:
- [x] Creation from an index field (requires: struct data types)
### Supported Functions
- [ ] Set up a collection of all supported functions and define their typing behaviour
### FreezeExpressions
The freeze is performed by `freeze.FreezeExpressions`. It inherits from `pb.SympyToPymbolicMapper` and extends it functionality where necessary.
- [x] `Field.Access`:
- [x] Generic and Absolute accesses
- [x] Index field accesses: Require `SparseIterationSpace`
- [x] Buffer accesses: Require inference of buffer indexing
- [ ] Staggered Accesses
- [ ] Exponentials (`sp.Pow`): Expand to multiplications where sensible
- [ ] Functions and function calls: Check if supported and map onto backend function symbols
### Typification
Untyped expressions are typified by `typification.Typifier`.
- [x] Variables: Assign default type from configuration
- [x] Constants: Assign target type from current expression context, or fail if target type is unknown
- [x] Arithmetic expressions: Typify all operands homogenously
- [x] Struct member accesses (`Lookup`): Use struct type
- [x] Array accesses: Typify index expression with `index_dtype` target type
- [x] Functions: Typify according to function type rules (see above *Supported Functions*)
### Platform-Specific Code Generation
The remaining code generation process is controlled by the platform-specific code generator.
All platform-specific code generators shall be subclasses of `platform.Platform` which defines their common interface. They all shall implement:
- [ ] Materialization of the iteration space (sparse and full)
- [x] Basic CPU
- [x] Materialization of functions (map mathematical, in particular transcendental functions, to target-specific implementations)
- [ ] Target-specific optimizations (including vectorization)
- [x] Provide collection of required header files, preprocessor definitions, etc.
The platform instance shall be created from a `TargetSpec` provided by the user.
The main code generation driver shall then call the `Platform` object to finalize code generation.https://i10git.cs.fau.de/pycodegen/pystencils/-/issues/79C Code Printer2024-02-28T09:49:45+01:00Frederik HennigC Code PrinterThe C code printer should be a simple structure-recursive algorithm and not contain any transformation logic, nor should infer or compute any additional information from the AST.
Code printing is currently implemented in `emission.py`. ...The C code printer should be a simple structure-recursive algorithm and not contain any transformation logic, nor should infer or compute any additional information from the AST.
Code printing is currently implemented in `emission.py`. The printer should be kept minimial.
The `CPrinter` class is implemented as a visitor over the AST nodes and calls an instance of pymbolic's [CCodeMapper](https://documen.tician.de/pymbolic/mappers.html#pymbolic.mapper.c_code.CCodeMapper) to print expressions.
An extension of the printer for specific vector instruction sets, etc., should not be necessary: instead, all intrinsics, function calls, etc. should be modelled using pymbolic's `FunctionSymbol` and `Call` nodes.https://i10git.cs.fau.de/pycodegen/pystencils/-/issues/78Type System2024-02-28T09:49:31+01:00Frederik HennigType SystemWe shall extend and refactor the existing type system of pystencils. While the exact structure of the type system still bears discussion, the following features shall be supported:
- The base class `AbstractType` should model a type's ...We shall extend and refactor the existing type system of pystencils. While the exact structure of the type system still bears discussion, the following features shall be supported:
- The base class `AbstractType` should model a type's name and `const` qualifier, such that every type may be `const`.
- Any numerical types have to model their set of legal arithmetic and logical operations.
- Pointers should be general, they shall be able to point to any type.
- A `CustomType` class should be available to model any types outside the scope of pystencils which are managed only through their name (e.g. any types occuring in integration with C++ frameworks)
## Progress
- [x] Scalar Numeric Types: Signed/Unsigned Integers, IEEE-754 floats
- [x] Vector Numeric Types
- [x] Struct typeshttps://i10git.cs.fau.de/pycodegen/pystencils/-/issues/76Abstract Syntax Tree2024-01-27T14:07:33+01:00Frederik HennigAbstract Syntax TreeThe abstract syntax tree is the core data structure of the backend.
It represents the generated C code on a level of abstraction fit for pystencils.
The range of syntax and structure represented by the AST is from the kernel function (up...The abstract syntax tree is the core data structure of the backend.
It represents the generated C code on a level of abstraction fit for pystencils.
The range of syntax and structure represented by the AST is from the kernel function (upper limit) to just above the arithmetic-logical expressions (lower limit). Above that, code generation is the business of whatever framework pystencils is embedded in. Below, we use `pymbolic`.
## Class Hierarchy and Syntax Modelling
The base class of the AST nodes shall be `PsAstNode`. Its primary tasks are:
- Manage the branching structure
- Manage custom metadata
There shall be a dedicated class for each required syntactic structure, including:
- *Loops:* We model increment-type loops by `PsLoop`
- *Expressions:* Expressions are modelled using pymbolic. Each pymbolic expression in the AST is wrapped in the `PsExpression` leaf node or its subclasses.
## Typing in Expressions
In a well-formed AST, all pymbolic expressions must be fully typed; they may contain no bare `Variable`s and no bare constants, but only their typed variants `PsTypedVariable` and `PsTypedConstants` (see pycodegen/pystencils#77).
## Correspondence to frontend data structures and metadata
The AST is made up of a fixed set of backend data structures, and shall contain no frontend (e.g. SymPy-related) data structures in its primary structure. However, mappings to front-end structures might be managed within the AST metadata.
A metadata dict shall exist on each AST node for any frontend to store arbitrary data in. The backend shall not access that metadata.
## Visitors
As most algorithms on the AST will follow the visitor pattern, a type-based visitor dispatcher is provided to simplify the implementation of recursive visitor-type algorithms.https://i10git.cs.fau.de/pycodegen/pystencils/-/issues/77pymbolic extension: Typed Variables, Constants, and LinearizedArray2024-01-27T13:29:42+01:00Frederik Hennigpymbolic extension: Typed Variables, Constants, and LinearizedArrayThe `pymbolic` class hierarchy needs to be extended with a number of typed expressions, in particular for typed variables and constants.
## PsTypedVariable
The `PsTypedVariable` class shall be a direct subclass of the pymbolic `Variabl...The `pymbolic` class hierarchy needs to be extended with a number of typed expressions, in particular for typed variables and constants.
## PsTypedVariable
The `PsTypedVariable` class shall be a direct subclass of the pymbolic `Variable`, annotating a name with a type. In the backend, it takes the position of the old `TypedSymbol`.
## PsTypedConstant
pymbolic allows the injection of arbitrary types for constants, as long as they implement the arithmetic operators. We shall utilize this by introducing `PsTypedConstant` to model constants annotated with types in the backend, to make sure that all operations (e.g. constant folding) respect the type's behaviour.
## Linearized Arrays
To model linearized multi-dimensional arrays, we introduce the class `PsLinearizedArray` to the backend. We may model accesses to such arrays using a sublcass of `pymbolic.Subscript`, linked to a particular array instance.https://i10git.cs.fau.de/pycodegen/pystencils/-/issues/72Divisions are printed using `pow`2023-09-28T14:49:11+02:00Daniel BauerDivisions are printed using `pow`# Description
pystencils prints divisions like `1/x` as `pow(x, -1.0)`.
The floating point `pow` function is of course undesirable due to performance concerns.
It looks like pystencils intends to [print powers with small integer expon...# Description
pystencils prints divisions like `1/x` as `pow(x, -1.0)`.
The floating point `pow` function is of course undesirable due to performance concerns.
It looks like pystencils intends to [print powers with small integer exponent using multiplications and divisions](https://i10git.cs.fau.de/pycodegen/pystencils/-/blob/master/pystencils/backends/cbackend.py#L447-450).
However, the [types of powers are collated between the base and exponent](https://i10git.cs.fau.de/pycodegen/pystencils/-/blob/master/pystencils/typing/leaf_typing.py#L215-224).
Note that sympy converts divisions to powers automatically.
# MWE
```python
import pystencils as ps
from pystencils import CreateKernelConfig
from pystencils.astnodes import Block, KernelFunction, SympyAssignment
from sympy.abc import x, y
equ = SympyAssignment(y, 1/x)
typed_equ = ps.typing.transformations.add_types(equ, CreateKernelConfig())
print(typed_equ)
kernel = KernelFunction(
Block([typed_equ]),
ps.Target.CPU,
ps.Backend.C,
ps.cpu.cpujit.make_python_function,
None,
)
code = ps.get_code_str(kernel)
print(code)
```
## Actual output
```c
y ← x**CastFunc(-1, double)
FUNC_PREFIX void kernel(double x)
{
const double y = pow(x, -1.0);
}
```
## Expected output
```c
y ← x**CastFunc(-1, int)
FUNC_PREFIX void kernel(double x)
{
const double y = 1.0 / x;
}
```https://i10git.cs.fau.de/pycodegen/pystencils/-/issues/61Remove one of the Assignements2023-08-29T10:29:01+02:00Markus HolzerRemove one of the AssignementsThere exists the `Assignment` provided by SymPy and the `SympyAssignment`. Both do almost the same but are treated differently in the back. It is confusing to use both in pystencils. A better idea would be to remove the `Assignment` by S...There exists the `Assignment` provided by SymPy and the `SympyAssignment`. Both do almost the same but are treated differently in the back. It is confusing to use both in pystencils. A better idea would be to remove the `Assignment` by SymPy and only keep the pystencils version. This version is better adapted to our purposeshttps://i10git.cs.fau.de/pycodegen/pystencils/-/issues/71Installation fails for Python 3.112023-07-21T11:59:59+02:00Daniel BauerInstallation fails for Python 3.11Trying to install on Python 3.11 gives the error
```
gcc -Wsign-compare -DNDEBUG -O2 -Wall -U_FORTIFY_SOURCE -D_FORTIFY_SOURCE=3 -fstack-protector-strong -funwind-tables -fasynchronous-unwind-tables -fstack-clash-protection -Werror...Trying to install on Python 3.11 gives the error
```
gcc -Wsign-compare -DNDEBUG -O2 -Wall -U_FORTIFY_SOURCE -D_FORTIFY_SOURCE=3 -fstack-protector-strong -funwind-tables -fasynchronous-unwind-tables -fstack-clash-protection -Werror=return-type -g -DOPENSSL_LOAD_CONF -fwrapv -fno-semantic-interposition -O2 -Wall -U_FORTIFY_SOURCE -D_FORTIFY_SOURCE=3 -fstack-protector-strong -funwind-tables -fasynchronous-unwind-tables -fstack-clash-protection -Werror=return-type -g -IVendor/ -O2 -Wall -U_FORTIFY_SOURCE -D_FORTIFY_SOURCE=3 -fstack-protector-strong -funwind-tables -fasynchronous-unwind-tables -fstack-clash-protection -Werror=return-type -g -IVendor/ -fPIC -I/usr/include/python3.11 -c pystencils/boundaries/createindexlistcython.c -o build/temp.linux-x86_64-cpython-311/pystencils/boundaries/createindexlistcython.o
pystencils/boundaries/createindexlistcython.c:211:12: fatal error: longintrepr.h: No such file or directory
211 | #include "longintrepr.h"
| ^~~~~~~~~~~~~~~
compilation terminated.
error: command '/usr/bin/gcc' failed with exit code 1
```
Seems like every Python project using Cython eventually runs into this.
I *think* using a recent Cython version fixes the problem.
See https://github.com/cython/cython/issues/4461 for more information.Markus HolzerMarkus Holzerhttps://i10git.cs.fau.de/pycodegen/pystencils/-/issues/70pyCUDA vs CuPy2023-06-22T17:26:18+02:00Markus HolzerpyCUDA vs CuPyShould we replace [PyCuda](https://documen.tician.de/pycuda/) with [CuPy](https://cupy.dev/)?
Advantages of [CuPy](https://cupy.dev/):
- AMD support
- probably higher maintained due to NVIDIA support
- SciPy compatible.Should we replace [PyCuda](https://documen.tician.de/pycuda/) with [CuPy](https://cupy.dev/)?
Advantages of [CuPy](https://cupy.dev/):
- AMD support
- probably higher maintained due to NVIDIA support
- SciPy compatible.Markus HolzerMarkus Holzerhttps://i10git.cs.fau.de/pycodegen/pystencils/-/issues/69sync data to GPU when loade2023-06-22T17:26:18+02:00Markus Holzersync data to GPU when loadeCalling dh.load_all does not sync the GPU arrays. This is not user friendlyCalling dh.load_all does not sync the GPU arrays. This is not user friendlyMarkus HolzerMarkus Holzerhttps://i10git.cs.fau.de/pycodegen/pystencils/-/issues/40SymPy 1.10 is broken2023-06-03T15:05:25+02:00Markus HolzerSymPy 1.10 is brokenThe Latest Sympy master is at the moment not supported by pystencils. Problems occur with `deepcopy`. An example of failure is shown here:
https://i10git.cs.fau.de/holzer/pystencils/-/jobs/658250.
The problem was introduced with this c...The Latest Sympy master is at the moment not supported by pystencils. Problems occur with `deepcopy`. An example of failure is shown here:
https://i10git.cs.fau.de/holzer/pystencils/-/jobs/658250.
The problem was introduced with this commit to SymPy:
https://github.com/sympy/sympy/commit/ef0de0c80ab13501076d9bf611a111250bbc88b0Markus HolzerMarkus Holzerhttps://i10git.cs.fau.de/pycodegen/pystencils/-/issues/68Loosen compiler requirement at startup2023-06-01T12:08:41+02:00Max ElfnerLoosen compiler requirement at startupHi,
when starting pystencils on a Win10 machine with no VS installed, I get the output below which is quite verbose but pretty much the issue is that Visual Studio compilers are not installed. This is a pretty strict requirement since M...Hi,
when starting pystencils on a Win10 machine with no VS installed, I get the output below which is quite verbose but pretty much the issue is that Visual Studio compilers are not installed. This is a pretty strict requirement since MSVC can be hard to license / install in a corporate framework. (Side remark: Installation can be performed hacking mingw in or just using conda)
I would like to ask the following:
- Would it be possible to use (many) of the functions, e.g. assignments and kernels, especially display, without a compiler?
- If (see docs) numba / llvm is supported, cna even the compilation be performed without a C compiler?
- **And so:** Wouldn't it be an option to perform the check for a compiler later on, making the package usable if no MSVC is installed?
and lastly:
- I can't seem to find documentation on where a config file is located or what the syntax might be to explicitly explained to specify the compiler to use (`mingw gcc`?)
```python
import pystencils as ps
Traceback (most recent call last):
File "C:\...\pystencils\cpu\msvc_detection.py", line 21, in get_environment
version_nr = find_latest_msvc_version_using_environment_variables()
File "C:\...\pystencils\cpu\msvc_detection.py", line 45, in find_latest_msvc_version_using_environment_variables
raise ValueError("Visual Studio not found.")
ValueError: Visual Studio not found.
During handling of the above exception, another exception occurred:
Traceback (most recent call last):
... longer traceback ...
ValueError: Visual Studio not found. Write path to VS folder in pystencils config
```https://i10git.cs.fau.de/pycodegen/pystencils/-/issues/64Re-introduce support for complex numbers2023-05-31T12:42:53+02:00Alexander ReinauerRe-introduce support for complex numberssince !292 complex numbers are not supported anymore, which is a feature I'm dependent on.since !292 complex numbers are not supported anymore, which is a feature I'm dependent on.https://i10git.cs.fau.de/pycodegen/pystencils/-/issues/67Matplotlib arrow rendering2023-03-28T09:23:10+02:00Helen SchottenhammlMatplotlib arrow renderingStarting from `matplotlib = 3.5`, pystencils' definition for 3D arrows is deprecated. Therefore, 3D stencils cannot be rendered anymore.
The exact issue is stated, e.g., here : https://github.com/matplotlib/matplotlib/issues/21688
W...Starting from `matplotlib = 3.5`, pystencils' definition for 3D arrows is deprecated. Therefore, 3D stencils cannot be rendered anymore.
The exact issue is stated, e.g., here : https://github.com/matplotlib/matplotlib/issues/21688
We should either put requirements on the matplotlib version or apply the suggested fix, i.e., changing https://i10git.cs.fau.de/pycodegen/pystencils/-/blob/master/pystencils/stencil.py#L419 to
```
class Arrow3D(FancyArrowPatch):
def __init__(self, xs, ys, zs, *args, **kwargs):
super().__init__((0, 0), (0, 0), *args, **kwargs)
self._verts3d = xs, ys, zs
def do_3d_projection(self, renderer=None):
xs3d, ys3d, zs3d = self._verts3d
xs, ys, zs = proj3d.proj_transform(xs3d, ys3d, zs3d, self.axes.M)
self.set_positions((xs[0], ys[0]), (xs[1], ys[1]))
return np.min(zs)
```Markus HolzerMarkus Holzerhttps://i10git.cs.fau.de/pycodegen/pystencils/-/issues/66Absolute access is probably not copied correctly after _eval_subs()2023-03-24T19:55:26+01:00Nils KohlAbsolute access is probably not copied correctly after _eval_subs()I think there is a bug in `Field.Access::_eval_subs(self, old, new)` as the `is_absolute_access` property is not copied.
It should be
```
def _eval_subs(self, old, new):
return Field.Access(self.field,
...I think there is a bug in `Field.Access::_eval_subs(self, old, new)` as the `is_absolute_access` property is not copied.
It should be
```
def _eval_subs(self, old, new):
return Field.Access(self.field,
tuple(sp.sympify(a).subs(old, new) for a in self.offsets),
tuple(sp.sympify(a).subs(old, new) for a in self.index),
is_absolute_access=self.is_absolute_access, # << this is missing, maybe there is more ...
dtype=self.dtype)
```
instead of
```
def _eval_subs(self, old, new):
return Field.Access(self.field,
tuple(sp.sympify(a).subs(old, new) for a in self.offsets),
tuple(sp.sympify(a).subs(old, new) for a in self.index),
dtype=self.dtype)
```
Otherwise when performing sympy substitutions, that property might get lost and subsequent `resolve_field_access()` endeavors may produce garbage.Markus HolzerMarkus Holzerhttps://i10git.cs.fau.de/pycodegen/pystencils/-/issues/58iteration_slice with step-size for GPU-kernels2023-03-17T11:26:56+01:00Alexander Reinaueriteration_slice with step-size for GPU-kernelsStep-size of `iteration_slice` is ignored for GPU-kernels.
```
import pystencils as ps
dh = ps.create_data_handling(domain_size=(5, 5), periodicity=True, default_target="gpu")
iteration_slice = ps.make_slice[1:-1:2, 1:-1:2]
config = p...Step-size of `iteration_slice` is ignored for GPU-kernels.
```
import pystencils as ps
dh = ps.create_data_handling(domain_size=(5, 5), periodicity=True, default_target="gpu")
iteration_slice = ps.make_slice[1:-1:2, 1:-1:2]
config = ps.CreateKernelConfig(target=dh.default_target, iteration_slice=iteration_slice)
field = dh.add_array("a")
assign = ps.Assignment(field.center, 1.0)
kernel = ps.create_kernel(assign, config=config).compile()
dh.fill(field.name, 0, ghost_layers=True)
if config.target == ps.enums.Target.GPU:
dh.to_gpu(field.name)
dh.run_kernel(kernel)
if config.target == ps.enums.Target.GPU:
dh.to_cpu(field.name)
print(dh.gather_array(field.name, ghost_layers=True))
```
results in
```
[[0. 0. 0. 0. 0. 0. 0.]
[0. 1. 1. 1. 1. 1. 0.]
[0. 1. 1. 1. 1. 1. 0.]
[0. 1. 1. 1. 1. 1. 0.]
[0. 1. 1. 1. 1. 1. 0.]
[0. 1. 1. 1. 1. 1. 0.]
[0. 0. 0. 0. 0. 0. 0.]]
```https://i10git.cs.fau.de/pycodegen/pystencils/-/issues/63SimpleExtrapolationOutflow normal to cast to integer2022-12-07T08:52:18+01:00Corentin LothodeSimpleExtrapolationOutflow normal to cast to integerI just had a problem with generating outflow conditions, because I entered a list of floating point numbers for the normal instead of integers. Is it possible to add an error or a cast to integers when generating the code ?I just had a problem with generating outflow conditions, because I entered a list of floating point numbers for the normal instead of integers. Is it possible to add an error or a cast to integers when generating the code ?