[GT4Py] Embedded Field-View - skip_values and boundaries

# Embedded Field-View - skip_values and boundaries - Shaped by: @egparedes, @havogt - Appetite (FTEs, weeks): 2-3 devs, full cycle - Developers: ## Problem This is a continuation of the embedded field view implementation. The following features are missing: - support for connectivities with _skip values_ - support for boundary conditions (semantic change in `where`) The following features are not supported with non-numpy fields: - `scan`: we need to pre-allocate a field for the result which is currently done with the numpy allocator; it should be deduced from the arguments. - broadcasting a scalar: we cannot infer the buffer, instead we would need the FunctionField which is not yet merged. The following cleanups where discussed in the previous cycle: - use the Standard Array API for slicing to a 0d array - implement the synatx for indexing as discussed in [GT4Py team workshop](https://hackmd.io/tE9eaKoOTL-apBRLtEUlJg) ## Solution Sub-projects in order of priority ### Connectivities with skip values The feature can be implemented in the `Field` implementation. As the connectivity encodes the mask of its local dimension, we don't need to store an extra mask. For this to work we need a 1to1 mapping of LocalDimension to connectivity (not only offset to connectivity). This is already a requirement in the gtfn backend for a similar reason. 3 steps are needed: - Take skip values into account when computing the inverse image: they need to be ignored. - Remap (`take`) needs to handle skip values (no action is needed as long as the field contains at least 1 value, then `-1` is implicitly wrapping around). This requires good documentation. Think of possible corner cases. - Reductions: need to take into account which values to ignore, e.g. `np.sum(..., initial=..., where=connectivity>=0)`. Prototype: https://github.com/havogt/gt4py/tree/embedded_skip_value_connectivity ### Shaping for Boundary conditions Look at boundary condition cases and decide which implementation is more useful in the short-term (eventually we want to implement both): - non-overlapping `concat` - with precedence `.at.set` Use icon4py examples as a starting point, e.g. `fused_velocity_advection_stencil_1_to_7` from https://github.com/C2SM/icon4py/pull/261. If we conclude quickly that the non-overlapping `concat` is the way to go, we can implement it with the strategy outlined in the previous cycles (TODO add reference). ``` 1_k_level: Field[[K], float] # Domain(K:(0,1)) interior: Field[[K], float] # Domain(K:(1,nlev)) # old semantics k_index: Field[[K], int] # Domain(K:(0,nlev)) where(k_index > 1 == 0, 1_k_level, interior) # not doable with new semantics # new way to express that concat(1_k_level, interior) # other example f1 = Field[[K], float] # Domain(K:(0,nlev)) f2 = Field[[K], float] # Domain(K:(0,nlev)) where(k_index < n, f1, f2) #or concat(f1[K(0):K(n)], f2[K(n):K(nlev)]) ``` ### 0d field slicing: Nikki Always keep fields objects, e.g. if slicing on a 1d field with a scalar value, return obejct should not be a scalar, rather a 0d field. Changing how slicing works here. Instead of slicing to a scalar, we should return a 0d array. First idea: This is the standard behavior of the Array API. Additionally, we add a `.item()` method to `Field` to extract the scalar value. Agreed idea: Instead a fbuiltin `as_scalar` since it's used also in the DSL. It is an error to call `as_scalar` on a non-0d field. Finalized idea: `as_scalar` is pointless ### Indexing: Nikki For domain slices, use the syntax sketched in the [GT4Py team workshop](https://hackmd.io/tE9eaKoOTL-apBRLtEUlJg): `f[I(-1):I(5)]`. The more compact alternative `f[I[-1:5]]` was discarded (for now) because it can also be interpreted as remapping a field. E.g. going from: ```python! f[(named_range(I, unit_range(5,10)), ...)] ``` To: ```python! f[I(5):I(10), ...] ``` `(slice(I(5), I(10)),)` -> `((I, UnitRange(5,10)))` `(slice(I(5), None),)` -> `((I, UnitRange(5,Infinity)))` `(slice(I(5), J(10))` -> raise ... `(slice(None), ...)` -> `(slice(None), ...)` f[:] Dimensions' indexes can also be interchangeable. If field has `[IDim, JDim]`, indexing can be expressed as: `f[J(5):J(10), I(0):I(6)]` #### Rabbit hole We need to investigate if we need to implement the same indexing when we parse (instead of embedded execution). If this is needed, we should implement it first for embedded, then shape implementation for parsing. ### Optional: Allocator in fields Store the allocator in the field, and use it for `scan` to decide how to allocate the result buffer. ### Optional: Shaping of how to proceed with FunctionField The result should be a clear idea if we should proceed with the FunctionField implementation or if we should implement it based with an NDArrayLike namespace. The _NDArrayLike namespace_ approach would be to have a single `Field` class (for which the domain computations are implemented) and dispatch the buffer operations to the respective buffer or function implementation. This requires understanding how the origin of infinity Fields fit into this picture ## Rabbit holes see hints in the sub-projects ## No-gos A full implementation of `concat` is not reachable. ## Progress  - [ ] Connectivities with skip values (Hannes) - [x] implemented - [ ] merged - [ ] Boundary conditions - [x] explore limited concat_where (1D mask, with automatic concat) - [ ] explore lazy field - [x] 0-d Fields (Nikki) - [ ] Index syntax (Nikki) - [ ] ... - [x] cupy support (Hannes): allocator in fields not needed for scan, see https://github.com/GridTools/gt4py/pull/1372 - [x] implementation - [x] merged - [ ] jax support (Hannes) - [ ] draft PR: blocked by proper implementation of creating jax arrays (not via allocators) - [ ] from previous cycle: as_offset - [x] extended remap (+ tests) (Enrique)