Express finite difference operators as convolutions

When #2539 is done we should start the work of figuring out and expressing existing finite difference operators as a convolution operation.

Since the assumption is that the 'convolution operation' will work on raw numbers, for each operator we need to determine:
- A correct basis in which the RHS operand needs to be and appropriate, space-invariant, convolution kernel for the operation
- We need to define function that applied to the operator will return the expected basis of the operand
- Some operators (e.g. divergence), may require additional operations to be done on input and output. We need another function that will communicate extra operations necessary.

The ultimate goal is to have all the information in the operator that would allow to translate a operation:
```
Op(field)
```
into
```
field |> field -> transform_to_basis(Op, field) |> field -> preprocess(Op, field) |> field -> convolve(Op, field |> field -> postprocess(Op, field)
```
I am abusing the notation/julia syntax here, but basically we need to be able to split application of the operator into multiple operations, some of them potentially noops. It should be done on during the evaluation of the `Broadcasted` (`StencilBroadcasted`) expression in a follow-up.