In this tutorial we will show how to access and navigate the Iteration/Expression Tree (IET) rooted in an `Operator`.

# Part I - Top Down

Let’s start with a fairly trivial example. First of all, we disable all performance-related optimizations, to maximize the simplicity of the created IET as well as the readability of the generated code.

``````from devito import configuration
configuration['opt'] = 'noop'
configuration['language'] = 'C'``````

Then, we create a `TimeFunction` with 3 points in each of the space `Dimension`s x and y.

``````from devito import Grid, TimeFunction

grid = Grid(shape=(3, 3))
u = TimeFunction(name='u', grid=grid)``````

We now create an `Operator` that increments by 1 all points in the computational domain.

``````from devito import Eq, Operator

eq = Eq(u.forward, u+1)
op = Operator(eq)``````

An `Operator` is an IET node that can generate, JIT-compile, and run low-level code (e.g., C). Just like all other types of IET nodes, it’s got a number of metadata attached. For example, we can query an `Operator` to retrieve the input/output `Function`s.

``op.input``
``(u(t, x, y),)``
``op.writes``
``(u(t, x, y),)``

If we print `op`, we can see how the generated code looks like.

``print(op)``
``````#define _POSIX_C_SOURCE 200809L
#define START(S) struct timeval start_ ## S , end_ ## S ; gettimeofday(&start_ ## S , NULL);
#define STOP(S,T) gettimeofday(&end_ ## S, NULL); T->S += (double)(end_ ## S .tv_sec-start_ ## S.tv_sec)+(double)(end_ ## S .tv_usec-start_ ## S .tv_usec)/1000000;

#include "stdlib.h"
#include "math.h"
#include "sys/time.h"

struct dataobj
{
void *restrict data;
unsigned long * size;
unsigned long * npsize;
unsigned long * dsize;
int * hsize;
int * hofs;
int * oofs;
void * dmap;
} ;

struct profiler
{
double section0;
} ;

int Kernel(struct dataobj *restrict u_vec, const int time_M, const int time_m, const int x_M, const int x_m, const int y_M, const int y_m, struct profiler * timers)
{
float (*restrict u)[u_vec->size[1]][u_vec->size[2]] __attribute__ ((aligned (64))) = (float (*)[u_vec->size[1]][u_vec->size[2]]) u_vec->data;

for (int time = time_m, t0 = (time)%(2), t1 = (time + 1)%(2); time <= time_M; time += 1, t0 = (time)%(2), t1 = (time + 1)%(2))
{
START(section0)
for (int x = x_m; x <= x_M; x += 1)
{
for (int y = y_m; y <= y_M; y += 1)
{
u[t1][x + 1][y + 1] = u[t0][x + 1][y + 1] + 1;
}
}
STOP(section0,timers)
}

return 0;
}
``````

An `Operator` is the root of an IET that typically consists of several nested `Iteration`s and `Expression`s – two other fundamental IET node types. The user-provided SymPy equations are wrapped within `Expressions`. Loop nest embedding such expressions are constructed by suitably nesting `Iterations`.

The Devito compiler constructs the IET from a collection of `Cluster`s, which represent a higher-level intermediate representation (not covered in this tutorial).

The Devito compiler also attaches to the IET key computational properties, such as sequential, parallel, and affine, which are derived through data dependence analysis.

We can print the IET structure of an `Operator`, as well as the attached computational properties, using the utility function `pprint`.

``````from devito.tools import pprint
pprint(op)``````
``````<Callable Kernel>
<CallableBody <unpacks=0, allocs=0, casts=1, maps=0, objs=0> <unmaps=0, frees=0>>
<List (0, 1, 0)>

<[affine,sequential] Iteration time::time::(time_m, time_M, 1)>
<Section (section0)>

<TimedList (1, 1, 1)>
<[affine,parallel,parallel=] Iteration x::x::(x_m, x_M, 1)>
<[affine,parallel,parallel=] Iteration y::y::(y_m, y_M, 1)>
<ExpressionBundle (1)>

<Expression u[t1, x + 1, y + 1] = u[t0, x + 1, y + 1] + 1>

``````

In this example, `op` is represented as a `<Callable Kernel>`. Attached to it are metadata, such as `_headers` and `_includes`, as well as the `body`, which includes the children IET nodes. Here, the body is the concatenation of an `PointerCast` and a `List` object.

``op._headers``
``OrderedSet([('_POSIX_C_SOURCE', '200809L'), ('START(S)', 'struct timeval start_ ## S , end_ ## S ; gettimeofday(&start_ ## S , NULL);'), ('STOP(S,T)', 'gettimeofday(&end_ ## S, NULL); T->S += (double)(end_ ## S .tv_sec-start_ ## S.tv_sec)+(double)(end_ ## S .tv_usec-start_ ## S .tv_usec)/1000000;')])``
``op._includes``
``OrderedSet(['stdlib.h', 'math.h', 'sys/time.h'])``
``op.body``
``<CallableBody <unpacks=0, allocs=0, casts=1, maps=0, objs=0> <unmaps=0, frees=0>>``

We can explicitly traverse the `body` until we locate the user-provided `SymPy` equations.

``print(op.body.casts[0])  # Printing the PointerCast``
``float (*restrict u)[u_vec->size[1]][u_vec->size[2]] __attribute__ ((aligned (64))) = (float (*)[u_vec->size[1]][u_vec->size[2]]) u_vec->data;``
``print(op.body.body[0])  # Printing the actual body``
``````for (int time = time_m, t0 = (time)%(2), t1 = (time + 1)%(2); time <= time_M; time += 1, t0 = (time)%(2), t1 = (time + 1)%(2))
{
START(section0)
for (int x = x_m; x <= x_M; x += 1)
{
for (int y = y_m; y <= y_M; y += 1)
{
u[t1][x + 1][y + 1] = u[t0][x + 1][y + 1] + 1;
}
}
STOP(section0,timers)
}``````

Below we access the `Iteration` representing the time loop.

``````t_iter = op.body.body[0].body[0]
t_iter``````
``<WithProperties[affine,sequential]::Iteration time[t0,t1]; (time_m, time_M, 1)>``

We can for example inspect the `Iteration` to discover what its iteration bounds are.

``t_iter.limits``
``(time_m, time_M, 1)``

And as we keep going down through the IET, we can eventually reach the `Expression` wrapping the user-provided SymPy equation.

``````expr = t_iter.nodes[0].body[0].body[0].nodes[0].nodes[0].body[0]
expr.view``````
``'<Expression u[t1, x + 1, y + 1] = u[t0, x + 1, y + 1] + 1>'``

Of course, there are mechanisms in place to, for example, find all `Expression`s in a given IET. The Devito compiler has a number of IET visitors, among which `FindNodes`, usable to retrieve all nodes of a particular type. So we easily can get all `Expression`s within `op` as follows

``````from devito.ir.iet import Expression, FindNodes
exprs = FindNodes(Expression).visit(op)
exprs[0].view``````
``'<Expression u[t1, x + 1, y + 1] = u[t0, x + 1, y + 1] + 1>'``