---
# System prepended metadata

title: 'Frontend proposal: Field view with local operations'
tags: [frontend]

---

###### tags: `frontend`

# Frontend proposal: Field view with local operations

The goal of this project is to develop a draft specification and prototype implementation of a functional frontend and its connection to the iterator-view model developed in the last cycle ([spec](https://hackmd.io/@gridtools/B1Yh6B8cd), [notes](https://hackmd.io/rT-Uwd2ASq6-h5rdlovZ_g), [impl](https://github.com/fthaler/gt4py_new_model/tree/iterator_v2/src/unstructured), [tests](https://github.com/fthaler/gt4py_new_model/tree/iterator_v2/tests/unstructured_tests)). The central datastructure exposed to the user will be a `Field` mapping indices (e.g. in I, J, K, Vertex etc.) to values (e.g. int, tuple of int etc.) while operations on fields are expressed locally.

## Description

As standard differential operators as well as tridiagonal solvers are inherently functional the primary target of this project is a purely functional frontend, i.e. every call can be replaced by its corresponding value without changing the program behaviour. The Laplacian is the prototypical example and will serve as a bare mimimum here to outline the frontend and its connection to the iterator-view model. 

```python
# laplacian operator corresponding to a fencil with a single closure
def laplacian(field: Field[(I, J), dtype]) -> Field[(I, J), dtype]:
    # stencils are functions with a local view
    @stencil
    def stencil(pos, f):
        """
        pos: tuple of indices
        f: accessor like pointing to current value of field
        """
        # boundary conditions must be enforced
        #  a mask field can be used in the iterator-view
        if pos not in field.domain[1:-1, 1:-1]
            return 0
        # a call to f corresponds to shift of field by its argument
        return -4 * f() + f(I-1) + f(I+1) + f(0, -1) + f(0, 1)

    # application of a stencil to a field returns a new field (field view)
    return apply_stencil(stencil, field)

# double-laplacian corresponding to either a fencil
#  - with two closures
#  - with one closure where stencil is lifted
#  - with one closure where the stencil is inlined
def laplap(field: Field[(I, J), dtype]) -> Field[(I, J), dtype]:
    return laplacian(laplacian(field))
```

__Prior work__
Experiments with a `Field` datastructure have been made in the last cycle with an unpolished but working prototype implementation for Cartesian in pure python. Additionally a functional IR based on the simply typed lambda calculus was developed, successfully lowered into using a tracing approach and an inlining pass was demonstrated. The lowering is partially chaotic and needs cleanup. A [language specification](https://gridtools.github.io/gtc/docs/latest/gt_frontend/language_reference.html#abstract) for the (stateful) unstructured prototype exists of which parts can be reused.

__Partially-functional extension__
Stateful operations will not disappear for the near future or ever. We need to enable users of the current GT4Py toolchain (e.g. Vulcan, ECMWF, Stefano Ubbiali) as well as future users porting existing models (e.g. Exclaim, MeteoSwiss) transitioning to a functional toolchain. For that purpose `apply_stencil` should be extended with an additional `out` argument into which the result is written. To ensure referential transparency with respect to the DSL the passed field must be externally supplied and may only be used once.

```python
apply_stencil(stencil: Callable, *inputs: Field, *, out: Field)
```

### Execution modes

Performance and debuggability are two opposing properties not archievable at the same time. To balance between the two, different execution modes may be chosen by the user depending on his needs in a particular situation. These executions modes rely on increasing fractions of the toolchain and serve as milestones in the implementation.

1. *Direct mode*
   Everything is executed directly in python. No code generation is needed. <br>
   This mode is mainly meant for verification and to allow users to understand the semantics. Users can propose new patterns to be supported by the toolchain by implementing them in this mode. It uses regular python for-loops and is thus prohibitively slow.
2. *Fast debug mode* (default)
   Each `apply_stencil` call is executed using one of the backends, while the remaining code is executed directly in python. <br>
   In this mode the user can set breakpoints inbetween `apply_stencil` calls and introspect or plot temporary fields.
3. *Fast mode*
   Everything is executed using one of the backends. Debugging is only possible in the target language, e.g. C++. <br>
   



## Tasks

- *Direct mode*: Implement the frontend syntax in pure python and test core cartesian as well as unstructured examples collected in the last cycle
- Write a specification of the frontend syntax & semantic
- *Fast debug mode*: Execute single `apply_stencil` calls by
  - Tracing the stencil into Lambda-Let IR (already done)
  - Lower into iterator view IR (`Field.__getitem__` -> `deref`+`shift`, `apply_stencil` -> `fencil`+`closure`
  - Execute using one of the iterator view backends (only python available currently)
- *Fast mode*: Translate an entire function consisting of multiple `apply_stencil` calls (tracing of the frontend partially done)
- [Optional] *Partially-functional extension*: Implement the `apply_stencil` "overload" with `out` argument
- [Optional] Operations on Fields with the same domain, e.g. `a+b`
- [Optional] Python bindings to the prototype implementation of `antonaf` currently being developed
- [Optional] Introduce a lift operation in the frontend and investigate its usefulness

## Goals

- Draft frontend specification
- Frontend prototype connected to the iterator-view prototype

## Appetite

1 FTE Full cycle (Cartesian)

1 FTE Full cycle (Unstructured)

## Risk factors

- The frontend prototype will rely on an advanced form of tracing that is complex (1000 lines of dense code) and only known to a single developer (tehrengruber).
- An insight of the last cycle is that this form of tracing is not sustainble and needs to be replaced (conceptually clear, implementation effort not so much).
- @havogt, developer of python part of the iterator-view prototype on holidays during most of the cycle

## Examples

### Cartesian

```python
def laplacian(field: "Field"):
    def laplacian(pos, f):
        if pos not in field.domain[1:-1, 1:-1]
            return 0
        return -4 * f() + f(I-1) + f(I+1) + f(0, -1) + f(0, 1)

    return apply_stencil(stencil, field)
```

```python
def laplacian(field: "Field"):
    @stencil
    def stencil(pos, f):
        if pos not in field.domain[1:-1, 1:-1]
            return 0
        return -4 * f(0, 0) + f(-1, 0) + f(1, 0) + f(0, -1) + f(0, 1)

    return apply_stencil(stencil, field)
```

__Composition__

```python
def laplap(field):
    lap = laplacian(field)
    laplap = laplacian(lap)
    return laplap
```

### Unstructured

```python
def laplacian(field: Field):
    @stencil
    def stencil(v1):
        weights = (-4, 1, 1, 1, 1)
        return sum(weight[l] * field[v2] for l, v2 in enumerate(v2v[v1]))

    return apply_stencil(stencil, field)
```

```python
def fvm_nabla(mesh, pp: Field[Vertex, dtype], S_MXX: Field[Edge, dtype], S_MYY: Field[Edge, dtype]):
    def stencil(e, s_mxx):
        zavg = 0.5 * sum(pp[v] for v in e2v[e])
        return s_mxx * zavg
        
    zavgS_MXX = apply_stencil(stencil, mesh.edges, S_MXX)
    zavgS_MYY = apply_stencil(stencil, mesh.edges, S_MYY)
        
    def stencil(v):
        pnabla_MXX = sum(S_MXX[e] * sign[v, local_e] for local_e, e in graph(v2e[v, :]))
        pnabla_MYY = sum(S_MYY[e] * sign[v, local_e] for local_e, e in graph(v2e[v, :]))
        
        return (pnabla_MXX / vol, pnabla_MYY / vol)
        
    return apply_stencil(stencil, mesh.vertices)
```

Composition works the same as before.