# ICON4Py Dycore Domain Scientist Readability

How to work in this document: Only work in the most recent subsection (first `##` after this intro). Discussions previously are considered archived, if you want to further discuss something, pull u a summary of the point and continue discussing there.

Goal: Improve the dycore to be domain science friendly. (Not to be confused with the GT4Py feedback workshop.)

## Discussions from 2022-10-18

### "Implementation from paper"

Summary: HV mentioned "as close as possible to paper" as a possible goal for the refactoring. MB doesn't want to emphasize "implementing from paper" to manage expectations. WS states that the e.g. the ICON paper is far from the ICON implementation.

HV: I would be interested in better high-level goals to set the direction in which we want to go. MB brought one which I summarize below as "Software Architecture". But I would be interested to hear other ideas (from scientist point of view) what we want to optimize for in a domain science readability refactoring. The best measure for me is still as close as possible to paper and that does **not** imply formulas before discretization etc, but as close as possible within the boundaries of current GT4Py. (Note: I am asking for what the high-level goal of this refactoring is, not what the goal of GT4Py or somethign else is.)

### Software Architecture

Summary: MB wants the dycore to be architected, for instance by using some classes, putting the stencils in contexts (and changing their names).

HV: I agree, but for this particular exercise, I would focus on domain science readability and than combine that with a software architecture perspective.

### Scope

#### Merge kernels into fewer files

HV: We wanted to do that anyway before starting this refactoring exercise. There is no reason, we have the different components (nh_diffusion, solve_nonhydro, velocity_advection) in separate files. (Small comment to CM: putting in the same file, is possible without limitations; we can still generate as many separate kernels as we like.)

#### Renaming of variables/functions/stencils

HV: This is the minimal task we want to do here and is without consequences (apart from having to rename in all places).

#### Restructuring code

This topic covers AD comment "... rename or just put them together in the file like in the orignal code".

HV: In this exercise, I would go for the best possible implementation ("best" as we are trying to define above) and not limit ourselves by integration constraints (maybe not even performance constraints, i.e. the example from CM). Therefore, I would limit it to one component for now and then discuss how we bring the result back.

#### Unused Fortran code

#### Changing GT4Py and Syntactic Sugar on top of GT4Py

HV: From my point of view adding syntactic sugar on top of GT4Py is absolutely in scope. Changing GT4Py is not in scope short term (because that might spawn longer discussions), but giving feedback how to improve GT4Py (and kicking of the afforementioned discussions) is in scope. Or slightly rephrased, don't let the refactoring be blocked by waiting for possible changes in GT4Py.

### TODO: feel free to add another discussion point here

## Discussions and comments until 2022-10-18

Hannes' original idea was to ask some people with domain experience (e.g. Anurag and Abishek) to refactor the dycore to be more readable, i.e. work on the code and bring it as close as possible to how you would implement it from a paper\*. In the discussion we decided maybe it makes sense to start with having a meeting first and discuss how to approach it, with people interested, especially in the dycore.

\(\*\) Assuming *close as possible to paper* is a good measure of *domain scientist readability*. (One discussion point would be this: what do we want to optimize for in this refactoring.)

MB: I'm not sure the objective of "implementing from paper" is what we need to emphasize. I would like to see the dycore to be "architected", for instance by using some classes, putting the stencils in contexts (and changing their names). I feel that saying "from paper" the expectations are set too high, people would assume there is automatic differentiation, for instance, something that, I guess, would be for a higher level.

OF: There are several takes on this, e.g. I) renaming of variables / functions / stencils II) restructuring of the code III) changes in GT4Py frontend (or add syntactic sugar on top of frontend) IV) ... . Scope of this activity is not entirely clear to me.

AD: I took a look at mo_nh_diffusion.f90 and its GT4Py version. Remarks:
1. The fact that I had to go open different files to read the kernels made it very difficult to follow the code. I still went through the first 8 kernels.
2. Kernel names (*kernel1, *kernel2) are not intuitive. Need to rename them appropriately or just put them together in the file like in the orignal code.
3. There are unused fortran codes that we dont use. What is needed to get mo_nh_diffusion in all python? Understand the issue with halo exhange. We can relax 3 with regards to halo exchange for now.
4. Kernels on its own were readable. Improvements are possible (e.g., kernel 1 line 67, the RHS is not kh_smag_1; its a wind) but needs careful reading. But they aren't many at least as far as I saw.

When I mentioned about readability, I only meant for 1-3 to begin with, which itself is a big step forward. This was my very first careful reading of GT4Py kernels and I found it very much acceptable. Great job!

WS: when you say "from a paper", are you thinking about something like https://rmets.onlinelibrary.wiley.com/doi/10.1002/qj.2378 for the dynamical core? You would find the paper very far from the actual code version. Mauro's suggestion is closer to realistic.

CM: Biggest hurdle to me personally is the split in ICON of the mathematical formulae between ICON setup time (when everything is pre-computed which is runtime constant, for performance reasons) and runtime. Example, divergence (where $c = Cell, e = Edge$):
($\nabla f)_c \approx \frac{1}{A_c} \sum_{e=1}^{3} \ell_e f_e$
Or in gt4py:

edge_length: Field[[EdgeDim], float],
edge_orientation: Field[[EdgeDim], float],
input_field: Field[[EdgeDim, KDim], float],
area: Field[[CellDim, float],

div = neighbor_sum(
    edge_length(C2E)*edge_orientation(C2E)*input_field(C2E), axis=C2EDim
)/area

In ICON this is split into a setup time function that defines `geofac_div`:

          ptr_int%geofac_div(jc,je,jb) = &
            & ptr_patch%edges%primal_edge_length(ile,ibe) * &
            & ptr_patch%cells%edge_orientation(jc,jb,je)  / &
            & ptr_patch%cells%area(jc,jb)

and then the runtime divergence stencil:

geofac_div: Field[[CellDim, C2EDim], float],
input_field: Field[[EdgeDim, KDim], float],

div = neighbor_sum(input_field(C2E) * geofac_div, axis=C2EDim)

This split is done for performance reasons, however, I think it greatly reduces readability. `geofac_div` is not a variable that should exist. The divergence should be defined as a gt4py function once (like above) and then be used where ever it is needed.

CM: Some of the things proposed here will get in the way of our stencil-by-stencil verification, and I am not sure it is a good idea to do them so early:
1. Renaming stencil: not problematic
2. Putting everything between two halo exchanges into one python file: not problematic (thanks to DSL preprocessor)
3. Deleting unused Fortran code paths: not problematic (although makes updating ICON very hard probably)
4. Restructuring: Problematic, if we cannot match between Fortran and gt4py kernel by kernel, we need a new way to verify

## First meeting

or can we discuss this offline?

### Questions:

(no order)

#### Which component do we start with?

HV: Diffusion, if Anurag will be the main contributor from the domain science side?
AD: I can do that and can also contribute in other parts of dycore. But as I said before, the first step in that direction is to do 1-3.
CM: Diffusion or velocity advection are good candidated imo.

#### When do we start?

HV: We could start now or wait until we have a Python driver (beginning of next year) for the module that we want to refactor.
Why is it good to have a Python driver? If we have a Python driver, we can easily test the refactored code against the existing version. Alternatively, developing the Python driver could be part of the refactoring to gain more familarity with the dycore in GT4Py.

AD: Do you think the driver is needed for 1-3? Driver is much more complex and will take some planning. But you are right, it is time to start working on that too. We (I and MB) have discussed about it.

CM: I think everything will be easier if code between two halo exchanges is in one gt4py file. This is part of the current cycle, however, it will bleed into next cycle since I had to deal with probtest issue.

#### How much do we diverge from the Fortran structure?

HV: Opinion: We don't set any limits to make this exercise maximally useful. Even to a degree where we cannot currently integrate the refactored code into Fortran, e.g. because of halo exchanges.

MB: If the interface to the dycore is well redefined in ICON, I do not see why we could not integrate the Python dycore, even with halo-exchanges in Python. I would agree in not setting the limit, I'm wondering if the ability to run ICON with the Python dycore is not a requirement for EXCLAIM.

AD: Yes, this is a good question. Doing 1-3 will most likely lead to significant code divergence such that authomatic merge will not be possible. This is point of discussion in the next EC+BoD meeting. But we do need for quite sometime a dycore that is driven by Fortran.

CM: Things that will make future ICON merges very hard:
1. Renaming stencil: not problematic
2. Putting everything between two halo exchanges into one python file: not problematic (thanks to DSL preprocessor)
3. Deleting unused Fortran code paths: problematic
4. Restructuring: problematic

#### What do we optimize for?

MB: My opinion is that we should optimize for code maintainability and extendability. For instance we could design the code using constructs that can be understood by DaCe, so that we could play with it later without refactoring completely. Not sure it this is realistic (do we know what are the boundaries here? Linus could know this one.)

WS: "optimize for code maintainability and extensibility" is a great goal, but I don't understand what this means in practice. I would rather emphasize "readability" by the domain scientist.

AD: Do 1-3 cause any interference with what you want?

CM: I think it is too early for major restructuring. Renaming and writing bigger python files are doable and do not seem to have negative side effects. Also we could try to eliminate this split of formulae between ICON setup time and runtime, which will also improve readability in my opinion.

#### People?

Who will actively contribute in this refactoring?

HV: Longer term: Could we even make this interesting to external people to contribute?

WS: Getting external people to contribute would be the goal. In the case of the dynamical core this will be a challenge. Possibly we could approach Daniel Reinert w.r.t. the advection. Andreas Jocksch is also now, to an extent, an advection domain scientist.

MB: DO you mean external to EXCLAIM? I'm definitely in for "external to the development team", for "external to EXCLAIM" the issue is to find people willing to invest the time.

AD: 1-3 doesn't need external contribution and is a big step forward.

CM: Short term: Anurag, Andreas Jocksch, ... Long term: DWD, MPI