Year: 2021

Abstract:
Recent work in equivariant deep learning bears strong similarities to physics. Fields over a base space are fundamental entities in both subjects, as are equivariant maps between these fields. In deep learning, however, these maps are usually defined by convolutions with a kernel, whereas they are partial differential operators (PDOs) in physics. Developing the theory of equivariant PDOs in the context of deep learning could bring these subjects even closer together and lead to a stronger flow of ideas. In this work, we derive a G-steerability constraint that completely characterizes when a PDO between feature vector fields is equivariant, for arbitrary symmetry groups G. We then fully solve this constraint for several important groups. We use our solutions as equivariant drop-in replacements for convolutional layers and benchmark them in that role. Finally, we develop a framework for equivariant maps based on Schwartz distributions that unifies classical convolutions and differential operators and gives insight about the relation between the two.

text: Recent work in equivariant deep learning bears strong similarities to physics. Fields over a base space are fundamental entities in both subjects, as are equivariant maps between these fields. In deep learning, however, these maps are usually defined by convolutions with a kernel, whereas they are partial differential operators (PDOs) in physics. Developing the theory of equivariant PDOs in the context of deep learning could bring these subjects even closer together and lead to a stronger flow of ideas. In this work, we derive a G-steerability constraint that completely characterizes when a PDO between feature vector fields is equivariant, for arbitrary symmetry groups G. We then fully solve this constraint for several important groups. We use our solutions as equivariant drop-in replacements for convolutional layers and benchmark them in that role. Finally, we develop a framework for equivariant maps based on Schwartz distributions that unifies classical convolutions and differential operators and gives insight about the relation between the two.

heading: 1 INTRODUCTION
text: In many machine learning tasks, the data exhibits certain symmetries, such as translation- and sometimes rotation-invariance in image classification. To exploit those symmetries, equivariant neural networks have been widely studied and successfully applied in the past years, beginning with Group convolutional neural networks (Cohen & Welling, 2016; Weiler et al., 2018b). A significant generalization of Group convolutional networks is given by steerable CNNs (Cohen & Welling, 2017; Weiler et al., 2018a; Weiler & Cesa, 2019), which unify many different pre-existing equivariant models (Weiler & Cesa, 2019). They do this by representing features as fields of feature vectors over a base space, such as R2 in the case of two-dimensional images. Layers are then linear equivariant maps between these fields.
This is very reminiscent of physics. There, fields are used to model particles and their interactions, with physical space or spacetime as the base space. The maps between these fields are also equivariant, with the symmetries being part of fundamental physical laws.
It is also noteworthy that these symmetries are largely ones that appear the most often in deep learning, such as translation and rotation equivariance. These similarities have already led to ideas from physics being applied in equivariant deep learning (Lang & Weiler, 2021).
However, one remaining difference is that physics uses equivariant partial differential operators (PDOs) to define maps between fields, such as the gradient or Laplacian. Therefore, using PDOs instead of convolutions in deep learning would complete the analogy to physics and could lead to even more transfer of ideas between subjects.
Equivariant PDO-based networks have already been designed in prior work (Shen et al., 2020; Smets et al., 2020; Sharp et al., 2020). Most relevant for our work are PDO-eConvs (Shen et al., 2020), which can be seen as the PDO-analogon of group convolutions. However, PDO-eConvs are only one instance of equivariant PDOs and do not cover the most common PDOs from physics, such as the gradient, divergence, etc. Very similarly to how steerable CNNs (Cohen & Welling, 2017; Weiler et al., 2018a; Weiler & Cesa, 2019) generalize group convolutions, we generalize PDO-eConvs by characterizing the set of all translation equivariant PDOs between feature fields over Euclidean spaces. Because of this analogy, we dub these equivariant differential operators steerable PDOs.
These steerable PDOs and their similarity to steerable CNNs also raise the question of how equivariant PDOs and kernels relate to each other, and whether they can be unified. We present a framework for equivariant maps that contains both steerable PDOs and convolutions with steerable kernels as special cases. We then prove that this framework defines the most general set of translation equivariant, linear, continuous maps between feature fields, complementing recent work (Aronsson, 2021) that describes when equivariant maps are convolutions. Since formally developing this framework requires the theory of Schwartz distributions, we cover it mainly in Appendix E, and the main paper can be read without any knowledge of distributions. However, we reference the main results from this framework in the paper where appropriate.
In order to make steerable PDOs practically applicable, we describe an approach to find complete bases for vector spaces of equivariant PDOs and then apply this method to the most important cases. We have also implemented steerable PDOs for all subgroups of O(2) (see supplementary material for the code). Our code extends the E2CNN library1 (Weiler & Cesa, 2019), which will allow practitioners to easily use both steerable kernels and steerable PDOs within the same library, and even to combine both inside the same network. Finally, we test our approach empirically by comparing steerable PDOs to steerable CNNs. In particular, we benchmark different discretization methods for the numerical implementation.
In summary, our main contributions are as follows:
• We develop the theory of equivariant PDOs on Euclidean spaces, giving a practical characterization of precisely when a PDO is equivariant under any given symmetry.
• We unify equivariant PDOs and kernels into one framework that provably contains all translation equivariant, linear, continuous maps between feature spaces.
• We describe a method for finding bases of the vector spaces of equivariant PDOs, and provide explicit bases for many important cases.
• We benchmark steerable PDOs using different discretization procedures and provide an implementation of steerable PDOs as an extension of the E2CNN library.

heading: 1.1 RELATED WORK
text: Equivariant convolutional networks Our approach to equivariance follows the one taken by steerable CNNs (Cohen & Welling, 2017; Weiler et al., 2018a; Weiler & Cesa, 2019; Brandstetter et al., 2021). They represent each feature as a map from the base space, such as Rd, to a fiber Rc that is equipped with a representation ρ of the point group G. Compared to vanilla CNNs, which have fiber R, steerable CNNs thus extend the codomain of feature maps.
A different approach is taken by group convolutional networks (Cohen & Welling, 2016; Hoogeboom et al., 2018; Weiler et al., 2018b). They represent each feature as a map from a group H acting on the input space to R. Because the input to the network usually does not lie in H , this requires a lifting map from the input space to H . Compared to vanilla CNNs, group convolutional networks can thus be understood as extending the domain of feature maps.
When H = Rd oG is the semidirect product of the translation group and a pointwise group G, then group convolutions on H are equivalent to G-steerable convolutions with regular representations. For finite G, the group convolution over G simply becomes a finite sum. LieConvs (Finzi et al., 2020) describe a way of implementing group convolutions even for infinite groups by using a Monte Carlo approximation for the convolution integral. Steerable CNNs with regular representations would have to use similar approximations for infinite groups, but they can instead also use (non-regular) finite-dimensional representations. Both the group convolutional and the steerable approach can be applied to non-Euclidean input spaces—LieConvs define group convolutions on arbitrary Lie groups and steerable convolutions can be defined on Riemannian manifolds (Cohen et al., 2019b; Weiler et al., 2021) and homogeneous spaces (Cohen et al., 2019a).
One practical advantage of the group convolutional approach employed by LieConvs is that it doesn’t require solving any equivariance constraints, which tends to make implementation of new groups easier. They also require somewhat less heavy theoretical machinery. On the other hand, steerable
1https://quva-lab.github.io/e2cnn/
CNNs are much more general. This makes them interesting from a theoretical angle and also has more practical advantages; for example, they can naturally represent the symmetries of vector field input or output. Since our focus is developing the theory of equivariant PDOs and the connection to physics, where vector fields are ubiquitous, we are taking the steerable perspective in this paper.
Equivariant PDO-based networks The work most closely related to ours are PDO-eConvs (Shen et al., 2020), which apply the group convolutional perspective to PDOs. Unlike LieConvs, they are not designed to work with infinite groups. The steerable PDOs we introduce generalizes PDO-eConvs, which are obtained as a special case by using regular representations.
A different approach to equivariant PDO-based networks was taken by Smets et al. (2020). Instead of applying a differential operator to input features, they use layers that map an initial condition for a PDE to its solution at a fixed later time. The PDE has a fixed form but several learnable parameters and constraints on these parameters—combined with the form of the PDE—guarantee equivariance.
Sharp et al. (2020) also use a PDE, namely the diffusion equation, as part of their DiffusionNet model, which can learn on 3D surfaces. But instead of exploiting the connection between PDOs and convolutions, they solve the PDE in the spectral domain. Interestingly, the time evolution operator for the diffusion equation is exp(t∆), which can be interpreted as an infinite power series in the Laplacian, very reminiscent of the finite Laplacian polynomials that naturally appear throughout this paper. Studying the equivariance of such infinite series of PDOs might be an interesting direction for future work.
2 STEERABLE PDOS
In this section, we develop the theory of equivariant PDOs. We will represent all features as smooth fields f : Rd → Rc that associate a feature vector f(x) ∈ Rc, called the fiber at x, to each point x ∈ Rd. We write Fi = C∞(Rd,Rci) for the space of these fields f in layer i. Additionally, we have a group of transformations acting on the input space Rd, which describes under which symmetries we want the PDOs to be equivariant. We will always use a group of the form H = (Rd,+)oG, for some G ≤ GL(d,R). Here, (Rd,+) refers to the group of translations of Rd, while G is some group of linear invertible transformations. The full group of symmetries H is the semidi-
rect product of these two, meaning that each element h ∈ H can be uniquely written as h = tg, where t ∈ Rd is a translation and g ∈ G a linear transformation. For example, if G = {e} is the trivial group, we consider only equivariance under translations, as in classical CNNs, while for G = SO(d) we additionally consider rotational equivariance.
Each feature space Fi has an associated group representation ρi : G→ GL(ci,R), which determines how each fiber Rci transforms under transformations of the input space. Briefly, ρi associates an invertible matrix ρi(g) to each group element g, such that ρi(g)ρi(g′) = ρi(gg′); more details on representation theory can be found in Appendix B. To see why these representations are necessary, consider the feature space F = C∞(R2,R2) and the group G = SO(2). The two channels could simply be two independent scalar fields, meaning that rotations of the input move each fiber but do not transform the fibers themselves. Formally, this would mean using trivial representations ρ(g) = 1 for both channels. On the other hand, the two channels could together form a vector field, which means that each fiber would need to itself be rotated in addition to being moved. This would correspond to the representation ρ(g) = g. These two cases are visualized in Fig. 1.
In general, the transformation of a feature f ∈ Fi under an input transformation tg with t ∈ Rd and g ∈ G is given by
((tg) Bi f)(x) := ρi(g)f(g −1(x− t)) . (1)
The g−1(x− t) term moves each fiber spatially, whereas the ρi(g) is responsible for the individual transformation of each fiber.
For a network, we will need maps between adjacent feature spaces Fi and Fi+1. Since during this theory section, we only consider single layers in isolation, we will drop the index i and simply denote the layer map as Φ : Fin → Fout. We are particularly interested in equivariant maps Φ, i.e. maps that commute with the action of H on the feature spaces:
Φ(hBin f) = hBout Φ(f) ∀h ∈ H, f ∈ Fin . (2)
We call Φ translation-equivariant if Eq. (2) holds for h ∈ (Rd,+), i.e. for pure translations. Analogously, Φ is G-equivariant if it holds for linear transformations h ∈ G. Because H is the semidirect product of Rd and G, a map is H-equivariant if and only if it is both translation- and G-equivariant.

heading: 2.1 PDOS AS MAPS BETWEEN FEATURE SPACES
text: We want to use PDOs for the layer map Φ, so we need to introduce some notation for PDOs between multi-dimensional feature fields. As shown in Fig. 2, such a multi-dimensional PDO can be interpreted as a matrix of one-dimensional PDOs. Specifically, a PDO from C∞(Rd,Rcin) to C∞(Rd,Rcout) is described by a cout × cin matrix. For example, the 2D divergence operator, which maps from R2 to R = R1 can be written as the 1 × 2 matrix (∂1 ∂2). This is exactly analogous to convolutional kernels, which can also be interpreted as cout × cin matrices of scalar-valued kernels. To work with the one-dimensional PDOs that make up the entries of this matrix, we use multi-index notation, so for a tuple α = (α1, . . . , αd) ∈ Nd0, we write ∂α := ∂ α1 1 . . . ∂ αd d . A general one-
dimensional PDO is a sum ∑ α cα∂
α, where the coefficients cα are smooth functions cα : Rd → R (so for now no spatial weight sharing is assumed). The sum ranges over all multi-indices α, but we require all but a finite number of coefficients to be zero everywhere, so the sum is effectively finite. As described, a PDO between general feature spaces is then a matrix of these one-dimensional PDOs.

heading: 2.2 EQUIVARIANCE CONSTRAINT FOR PDOS
text: We now derive a complete characterization of the PDOs that are H-equivariant in the sense defined by Eq. (2). Because a map is equivariant under the full symmetries H = (Rd,+) oG if and only if it is both translation equivariant and G-equivariant, we split up our treatment into these two requirements.
First, we note that translation equivariance corresponds to spatial weight sharing, just like in CNNs: Proposition 1. A cout × cin-PDO Φ with matrix entries ∑ α c ij α ∂
α is translation equivariant if and only if all coefficients cijα are constants, i.e. c ij α (x) = c ij α (x ′) for all x, x′ ∈ Rd.
The proof can be found in Appendix F.
So from now on, we restrict our attention to PDOs with constant coefficients and ask under which circumstances they are additionally equivariant under the action of the point group G ≤ GL(d,R). To answer that, we make use of a duality between polynomials and PDOs that will appear throughout this paper: for a PDO ∑ α cα∂ α, with cα ∈ R, there is an associated polynomial simply given by
∑ α cαx α, where xα := xα11 . . . x αd d . Conversely, for any polynomial p = ∑ α cαx
α, we get a PDO by formally plugging in ∂ = (∂1, . . . , ∂d) for x, yielding p(∂) := ∑ α cα∂
α. We will denote the map in this direction by D, so we write D(p) := p(∂).
We can extend this duality to PDOs between multi-dimensional fields: for a matrix P of polynomials, we define D(P ) component-wise, so D(P ) will be a matrix of PDOs given by D(P )ij := D(Pij). To avoid confusion, we will always denote polynomials by lowercase letters, such as p, q, and matrices of polynomials by uppercase letters, like P and Q.
As a simple example of the correspondence between polynomials and PDOs, the Laplacian operator ∆ is given by ∆ = D(|x|2). The gradient (∂1, . . . , ∂d)T is induced by the d×1 matrix (x1, . . . , xd) T and the 3D curl (d = 3) is induced by the 3× 3 matrix
P = ( 0 −x3 x2 x3 0 −x1 −x2 x1 0 ) =⇒ D(P ) = ( 0 −∂3 ∂2 ∂3 0 −∂1 −∂2 ∂1 0 ) . (3)
Finally, note that D is a ring isomorphism, i.e. a bijection with D(p + q) = D(p) + D(q) and D(pq) = D(p) ◦D(q), allowing us to switch viewpoints at will. We are now ready to state our main result on the G-equivariance of PDOs:
Theorem 2. For any matrix of polynomials P , the differential operator D(P ) is G-equivariant if and only if it satisfies the PDO G-steerability constraint,
P ( (g−1)Tx ) = ρout(g)P (x)ρin(g) −1 ∀g ∈ G, x ∈ Rd . (4)
P (x) just means that x is plugged into the polynomials in each entry of P , which results in a real-valued matrix P (x) ∈ Rcout×cin , and similarly for P ( (g−1)Tx ) . Appendix C provides more intuition for the action of group elements on polynomials. Because D is an isomorphism, this result completely characterizes when a PDO with constant coefficients is G-equivariant, and in conjunction with Proposition 1, when any PDO between feature fields is H-equivariant. The proof of Theorem 2 can again be found in Appendix F.
To avoid confusion, we would like to point out that this description of PDOs as polynomials is only a useful trick that lets us express certain operations more easily and will later let us connect steerable PDOs to steerable kernels. Convolving with these polynomials is not a meaningful operation; they only become useful when ∂ is plugged in for x.

heading: 2.3 EXAMPLES OF EQUIVARIANT PDOS
text: To build intuition, we explore the space of equivariant PDOs for two simple cases before covering the general solution. Consider R2 as a base space with the symmetry group G = SO(2) and trivial representations ρin and ρout, i.e. maps between two scalar fields. The equivariance condition then becomes p(gx) = p(x), so the polynomial p has to be rotation invariant. This is the case if and only if it can be written as a function of |x|2, i.e. p(x) = q(|x|2) for some q : R≥0 → R. Since we want p to be a polynomial, q needs to be a polynomial in one variable. Because |x|2 = x21 + x22, we get the PDO D(p) = q(∆), where ∆ := ∂21 + ∂ 2 2 is the Laplace operator. So the SO(2)-equivariant PDOs between two scalar fields are exactly the polynomials in the Laplacian, such as 2∆2 + ∆ + 3.
As a second example, consider PDOs that map from a vector to a scalar field, still with SO(2) as the symmetry group. There are two such PDOs that often occur in the natural sciences, namely the divergence, div v := ∂1v1 + ∂2v2, and the 2D curl, curl2D v := ∂1v2 − ∂2v1. Both of these are SO(2)-equivariant (see Appendix A). We get additional equivariant PDOs by composing with equivariant scalar-to-scalar maps, i.e. polynomials in the Laplacian. Specifically, q(∆) ◦ div and q(∆) ◦ curl2D are also equivariant vector-to-scalar PDOs, for any polynomial q. We will omit the ◦ from now on.
We show in Appendix A that these PDOs already span the complete space of SO(2)-equivariant PDOs from vector to scalar fields. Explicitly, the equivariant PDOs in this setting are all of the form
q1(∆) div + q2(∆) curl2D (5)
for polynomials q1 and q2. One example of this is the PDO shown in Fig. 2.
In these examples, as well as in other simple cases (see Appendix A), the equivariant PDOs are all combinations of well-known operators such as the divergence and Laplacian. That the notion of equivariance can reproduce all the intuitively “natural” differential operators suggests that it captures the right concept.

heading: 3 BASES FOR SPACES OF STEERABLE PDOS
text: For the purposes of deep learning, we need to be able to learn steerable PDOs. To illustrate how to achieve this, consider again SO(2)-equivariant PDOs mapping from vector to scalar field. We have seen in Section 2.3 that they all have the form
q1(∆) div + q2(∆) curl2D . (6)
In practice, we will need to limit the order of the PDOs we consider, so that we can discretize them. For example, we could consider only polynomials q1 and q2 of up to first order, i.e. q1(z) = c1z + c2 and q2(z) = c3z + c4. This leads to PDOs of up to order three (since ∆ is a second order PDO and div and curl2D are first order). The space of such equivariant PDOs is then
c1∆ div + c2 div + c3∆ curl2D + c4 curl2D , ci ∈ R . (7)
We can now train the real-valued parameters c1, . . . , c4 and thereby learn arbitrary equivariant PDOs of up to order three.
The general principle is that we need to find a basis of the real vector space of steerable PDOs; then we can learn weights for a linear combination of the basis elements, yielding arbitrary equivariant PDOs.
Different group representations are popular in equivariant deep learning practice; for example PDOeConvs and group convolutions correspond to so-called regular representations (see Appendix L) but quotient representations have also been used very successfully (Weiler & Cesa, 2019). We therefore want to find bases of steerable PDOs for arbitrary representations ρin and ρout. To do so, we make use of the existing work on solving the closely related G-steerability constraint for kernels. In this section, we will first give a brief overview of this kernel steerability constraint and then describe how to transfer its solutions to ones of the PDO steerability constraint.

heading: 3.1 THE STEERABILITY CONSTRAINT FOR KERNELS
text: The kernel G-steerability constraint characterizes when convolution with a kernel κ : Rd → Rcout×cin is G-equivariant, like the PDO steerability constraint Eq. (4) does for PDOs. Namely,
κ(gx) = |det g|−1ρout(g)κ(x)ρin(g)−1 ∀g ∈ G (8)
has to hold. This constraint was proven in (Weiler et al., 2018a) for orthogonal G and later in (Weiler et al., 2021) for general G. It is very similar to the PDO steerability constraint: the only differences are the determinant that appears here and the κ(gx) term on the LHS where we had (g−1)T instead of g.
An explanation for this similarity is that both constraints can be seen as special cases of a more general steerability constraint for Schwartz distributions, which we derive in Appendix E. Schwartz distributions are a generalization of classical functions and convolutions with Schwartz distributions can represent both classical convolutions and PDOs. As we prove in Appendix E, such distributional convolutions are the most general translation equivariant, continuous, linear maps between feature spaces, strictly more general than either PDOs or classical kernels. We also show how steerable PDOs can be interpreted as the Fourier transform of steerable kernels, which we use to explain the remaining differences between the two steerability constraints.
But for the purposes of this section, we want to draw particular attention to the fact that forG ≤ O(d), the two constraints become exactly identical: the determinant then becomes 1, and (g−1)T = g. So for this case, which is by far the most practically important one, we will use existing solutions of the kernel steerability constraint to find a complete basis of equivariant PDOs.

heading: 3.2 TRANSFERRING KERNEL SOLUTIONS TO STEERABLE PDO BASES
text: Solutions of the kernel steerability constraint have been published for subgroups of O(2) (Lang & Weiler, 2021; Weiler & Cesa, 2019) and O(3) (Lang & Weiler, 2021; Weiler et al., 2018a), and they all use a basis of the form
{καβ(x) := ϕα(|x|)χβ (x/|x|) |α ∈ A, β ∈ B} , (9)
whereA and B are index sets for the radial and angular part and the ϕα span the entire space of radial functions. The reason is that the steerability constraint Eq. (8) constrains only the angular part if G ≤ O(d), because orthogonal groups preserve distances. The angular functions χβ are only defined on the sphere Sd−1. But crucially, we show in Appendix I that they can all be canonically extended to polynomials defined on Rd.2 Concretely, we define χ̃β(x) := |x|lβχβ(x/|x|), where lβ is chosen minimally such that χ̃β is a polynomial. What we prove in Appendix I is that such an lβ always exists.
The radial part ϕα in the kernel basis consists of unrestricted radial functions. To get a basis for the space of polynomial steerable kernels, it is enough to use only powers of |x|2 for the radial part. Specifically, we show in Appendix G that{
pkβ := |x|2kχ̃β ∣∣∣ k ∈ N0, β ∈ B} (10)
is a basis for the space of polynomial steerable kernels. As a final step, we interpret each polynomial as a PDO using the isomorphism D defined in Section 2.2, yielding a complete basis of the space of steerable PDOs. The |x|2k terms become Laplacian powers ∆k, while in the examples from Section 2.3, χ̃β corresponds to divergence and curl.
In Appendix H, we apply this procedure to subgroups of O(2) and O(3) to obtain concrete solutions.

heading: 4 EXPERIMENTS
text: Implementation We developed the theory of steerable PDOs in a continuous setting but for implementation the PDOs need to be discretized, just like steerable kernels. The method of discretization is completely independent of the computation of the steerable PDO basis, so steerable PDOs can be combined with any discretization procedure. For our experiments, we compare three methods, finite differences (FD), radial basis function finite differences (RBF-FD) and Gaussian derivatives. In all three cases, the PDOs are discretized as stencils, with which we then convolve.
Finite differences are a generalization of the usual central difference approximation and are the method used by PDO-eConvs (Shen et al., 2020). One of their appeals is that they guarantee an asymptotic upper bound on the discretization error, for the distance between grid points approaching zero. This in turn places the same asymptotic bound on the equivariance error, as pointed out by Shen et al. (2020).
RBF-FD finds stencil entries by demanding that the discretization should become exact when applied to radial basis functions placed on the stencil points. Its advantage over FD is that it can be applied to structureless point clouds rather than only to regular grids.
Gaussian derivative stencils work by placing a Gaussian on the target point and then evaluating its derivative on the stencil points; these are the stencil entries. Like RBF-FD, this also works on stencils, and the Gaussian also has a slight smoothing effect, which is why this discretization is often used in Computer Vision. Formal descriptions of all three discretization methods can be found in Appendix K.
In addition to discretization, the infinite basis of steerable PDOs or kernels needs to be restricted to a finite subspace. For kernels, we use the bandlimiting filters by Weiler & Cesa (2019). For PDOs,
2This holds for the cases we discuss here, i.e. subgroups of O(2) and O(3). In higher dimensions, the situation is more complicated, but for those the kernel solutions have not yet been worked out anyway.
3As in the original PDO-eConv paper (Shen et al., 2020). Note that their performance is better, which is simply caused by their different architecture and hyperparameters.
we restrict the total derivative order to two for 3× 3 stencils and to three for 5× 5 stencils (except for PDO-eConvs, where we replicate the original basis that restricts the maximum order of partial derivatives).
Table 1: MNIST-rot results. Test errors ± standard deviations are averaged over six runs. Vanilla CNN is a solely translation equivariant model (G = {e}) with the same general architecture. See main text for details on the models.
Representation Method Stencil Error [%] Params
– VanillaCNN 3× 3 2.001± 0.030 1.1M 5× 5 1.959± 0.055
regular (our basis)
Kernels 3× 3 0.741± 0.036 837K 5× 5 0.683± 0.021 1.1M
FD 3× 3 1.196± 0.062 837K 5× 5 1.54 ± 0.32 941K
RBF-FD 3× 3 1.313± 0.065 837K 5× 5 1.475± 0.020 941K
Gauss 3× 3 0.795± 0.030 837K 5× 5 0.750± 0.017 941K
regular (PDOeConv) FD3 5× 5 1.98 ± 0.11 982K Gauss 5× 5 0.831± 0.039
quotient (our basis)
Kernels 3× 3 0.717± 0.026 877K 5× 5 0.670± 0.011 1.1M
FD 3× 3 1.143± 0.063 877K 5× 5 1.347± 0.026 951K
RBF-FD 3× 3 1.303± 0.077 877K 5× 5 1.422± 0.040 951K
Gauss 3× 3 0.825± 0.053 877K 5× 5 0.744± 0.040 951K
Finally, steerable PDOs and steerable kernels only replace the convolutional layers in a classical CNN. To achieve a fully equivariant network, all the other layers, such as nonlinearities or Batchnorm also need to be equivariant. Weiler & Cesa (2019) discuss in details how this can be achieved for various types of layers and for different group representations. In our experiments, we use exactly the same implementation they do and only replace the kernel layers with steerable PDOs. Care also needs to be taken with biases in the PDO layers in order to achieve equivariance. Here, we again follow (Weiler & Cesa, 2019) by adding a bias only to the trivial irreducible representations that make up ρout.
Rotated MNIST We first benchmark steerable PDOs on rotated MNIST (Larochelle et al., 2007), which consists of MNIST images that have been rotated by different angles, with 12k train and 50k test images. Our results can be found in Table 1. The models with 5 × 5 stencils use an architecture that Weiler & Cesa (2019) used for steerable CNNs, with six C16-equivariant layers followed by two fully connected layers. The first column gives the representation under which the six equivariant layers transform (see Appendix B for their definitions). PDO-eConvs implicitly use regular representations (see Appendix L), but with a slightly different basis than the one we present, so we test both bases. We also tested models that are D16-equivariant in their first layers and C16-equivariant in their last one but did not find any improvements, see Appendix M. For the models with 3× 3 stencils, we use eight instead of six C16-equivariant layers, in order to compensate for the smaller receptive field and keep the parameter count comparable. The remaining differences between kernel and PDO parameter counts come from the fact that the basis restrictions necessarily work slightly differently (via bandlimiting filters or derivative order restriction respectively). All models were trained with 30 epochs and hyperparameters based on those by Weiler & Cesa (2019), though we changed the learning rate schedule and regularization slightly because this improved performance for all models, including kernel-based ones. The training data is augmented with random rotations. Precise descriptions of the architecture and hyperparameters can be found in Appendix N.
STL-10 The rotated MNIST dataset has global rotational symmetry by design, so it is unsurprising that equivariant models perform well. But interestingly, rotation equivariance can also help for natural images without global rotational symmetry (Weiler & Cesa, 2019; Shen et al., 2020). We therefore benchmark steerable PDOs on STL-10 (Coates et al., 2011), where we only use the labeled portion of 5000 training images. The results are shown in Table 2. The model architecture and hyperparameters are exactly the same as in (Weiler & Cesa, 2019), namely a Wide-ResNet-16-8 trained for 1000 epochs with random crops, horizontal flips and Cutout (DeVries & Taylor, 2017) as data augmentation. The group column describes the equivariance group in each of the three residual blocks. For example, D8D4D1 means that the first block is equivariant under reflections and 8 rotations, the second under 4 rotations and the last one only under reflections. All layers use regular
representations. The D8-equivariant layers use 5 × 5 filters to improve equivariance, whereas the other layers use 3× 3 filters.
Discussion While all equivariant models improve significantly over the non-equivariant CNN, the method of discretization plays an important role for PDOs. The reason that FD and RBF-FD underperform kernels is that they don’t make full use of the stencil, since PDOs are inherently local operators. When a 5 × 5 stencil is used, the outermost entries are all very small compared to the inner ones, and even in 3 × 3 kernels, the four corners tend to be closer to zero (see Appendix M for images of stencils to illustrate this). Gaussian discretization performs significantly better and almost as well as kernels because its smoothing effect alleviates these issues. This fits the observation that kernels and Gaussian methods profit from using 5× 5 kernels, whereas these do not help for FD
and RBF-FD (and in fact decrease performance because of the smaller number of layers). We also observe a very small but rather consistent advantage of quotient representations over regular ones, demonstrating the practical usefulness of non-regular representations.

heading: 5 CONCLUSION
text: We have described a general framework for equivariant PDOs acting on feature fields over Euclidean space. With this framework, we found strong similarities between equivariant PDOs and equivariant convolutions, even unifying the two using convolutions with Schwartz distributions. We exploited these similarities to find bases for equivariant PDOs based on existing solutions for steerable kernels.
Our experiments show that the locality of PDOs can be a disadvantage compared to convolutional kernels. However, our approach for equivariance can easily be combined with any discretization method, and we show that using Gaussian derivatives for discretization alleviates the issue. Equivariant PDOs could also be very useful in cases where their strong locality is a desideratum rather than a drawback.
The theory developed in this work provides the necessary foundation for applications where equivariant PDOs, rather than kernels, are needed. For example, Probabilistic Numerical CNNs (Finzi et al., 2021) use PDOs in order to parameterize convolutions on continuous input data. Finzi et al. also derive a constraint to make these PDOs equivariant, which is a special case of our PDO G-steerability constraint Eq. (4). The solutions to this constraint presented in this paper are the missing piece for implementing and empirically evaluating Probabilistic Numerical CNNs – neither of which Finzi et al. do.
Another promising application of PDOs is an extension to manifolds. Gauge CNNs (Cohen et al., 2019b; Kicanaoglu et al., 2019; Haan et al., 2021; Weiler et al., 2021) are a rather general framework for convolutions on manifolds; see Weiler et al. (2021) for a thorough treatment and literature review. As in the Euclidean case, Gauge CNNs use feature fields that transform according to some representation ρ. The kernels are defined on the tangent space and are still constrained by the G-steerability constraint. Because of that, our approach to equivariant Euclidean PDOs is very well-suited for generalization to manifolds and our steerable PDO solutions will still remain valid. One advantage of PDOs in this setting is that they require no Riemannian structure and can achieve equivariance with respect to arbitrary diffeomorphisms, as is common in physics, instead of mere isometry equivariance.

heading: REPRODUCIBILITY STATEMENT
text: All Propositions and Theorems in the main paper and in the appendix explicitly state all their assumptions. In cases where theoretical results are stated only informally in the main paper (such as
our results involving Schwartz distributions), the precise claims can be found in the appendix. The appendix also contains complete proofs for all of our theoretical claims.
The code necessary to reproduce all of our experiments is included in the supplementary material and will be published upon publication. The required datasets are downloaded and preprocessed automatically. The exact hyperparameters we used are available as pre-defined configuration options in our scripts. We also include a version lockfile for installing precisely the right versions of all required Python packages. Our implementation of steerable PDOs is easy to adapt to different use cases and fully documented, allowing other practitioners to test the method on different datasets or tasks.

heading: SUPPLEMENTARY MATERIAL

heading: A STEERABLE PDOS BETWEEN VECTOR AND SCALAR FIELDS
text: We give bases for the space of steerable PDOs for many important cases in Appendix J. However, the generality of the description there obscures the connection to well-known PDOs such as the gradient or divergence. So to complement the general solutions, we discuss a few simple cases in much more detail in this section. We will see that Laplacian, gradient, divergence, and curl are all rotation equivariant and, more interestingly, that all rotation equivariant PDOs can be constructed by combining these (in the simple settings we cover in this section). In many cases, we rederive the solutions even though all of them would follow immediately from the general case in Appendix J, in order to provide some intuition on why these are the only equivariant PDOs. Readers who are only interested in an overview of the results may wish to skip to the end of this section.

heading: A.1 SCALAR TO SCALAR PDOS
text: We have already argued in Section 2.3 that the SO(2)-equivariant PDOs between two scalar fields are precisely polynomials in the Laplacian, i.e. of the form q(∆) for an arbitrary real polynomial q ∈ R[x]. The derivation given there applies without changes to SO(d) for d > 2 and to O(d) as well, so the same holds in these cases.

heading: A.2 SCALAR TO VECTOR PDOS
text: We start by considering the case where ρin is trivial (with cin = 1) and ρout is the vector field representation (i.e. cout = d). We can then represent PDOs by a d× 1 matrix of polynomials, i.e. a column vector. The PDO steerability constraint becomes
P (gx) = gP (x) (11)
since ρin is trivial. A good mental model for this subsection is to think of P as a vector field P : Rd → Rd whose entries happen to be polynomials. The steerability constraint simply states that this vector field must “look the same” after a rotation, see Fig. 3 for examples.
We begin by discussing the case G = SO(2), which is somewhat different from O(2) and from d > 2. Any rotation equivariant vector field on R2 is fully determined by its values on the ray {(x1, 0) |x1 > 0}. Specifically, the rotation equivariant vector fields v are precisely those that can in polar coordinates be written as
v(r, ϕ) = ( cosϕ − sinϕ sinϕ cosϕ )( f1(r) f2(r) ) (12)
for arbitrary radial parts f1, f2. In Cartesian coordinates, we can write this as
v(x, y) = ( x −y y x )( f̃1(x 2 + y2)
f̃2(x 2 + y2)
) (13)
where f̃i(z) := fi( √ z)√ z
. We need v1 and v2 to be polynomials in x, y, which means that the f̃i need to be polynomials. If we then apply the D isomorphism, we get: Proposition 3. The SO(2)-equivariant differential operators from a scalar to a vector field are exactly those of the form(
∂1 −∂2 ∂2 ∂1 )( q1(∆) q2(∆) ) = q1(∆) ( ∂1 ∂2 ) + q2(∆) ( −∂2 ∂1 ) (14)
where q1 and q2 are arbitrary polynomials.
In words, the SO(2)-equivariant PDOs are all linear combinations of the gradient and the transpose of the 2D curl, (−∂2, ∂1)T , with coefficients being polynomials in the Laplacian (rather than just real numbers).
4Created using https://www.desmos.com/calculator/eijhparfmd
The gradient is also equivariant under reflections, i.e. O(2)-equivariant, and it easily generalizes to higher dimensions. However, the transpose 2D curl only appears in this particular setting: it is not reflection-equivariant, and it does not have an analogon in higher dimensions.5 We summarize this in the following result:
Proposition 4. Let G = O(d) for d ≥ 2 or G = SO(d) for d > 2. Then the G-equivariant differential operators from a scalar to a vector field are exactly those of the form
q(∆) grad (15)
for an arbitrary polynomial q ∈ R[x].
A possible intuition for why SO(2) is a special case is that SO(2) ∼= S1 whereas SO(d) 6∼= Sd−1 for d > 2 and O(d) 6∼= Sd−1 for d ≥ 2. Our construction above heavily makes use of the fact that rotations of R2 correspond one-to-one to angles, i.e. points on S1, and this construction thus does not generalize to any other cases.
Before we prove Proposition 4, we show a helpful lemma:
Lemma 5. Let d > 2. Then for any linearly independent vectors v, w ∈ Rd, there is a rotation g ∈ SO(d) such that gv = v but gw 6= w. For d = 2, there is such a g ∈ O(2).
Proof of lemma. Since v and w are linearly independent, we can write Rd = span(v)⊕W , where W is a linear subspace containing w. Pick an element g̃ ∈ SO(d− 1) (or O(d− 1) = O(1) = {±1} in the d = 2 case) such that g̃w 6= w, where g̃ acts on W (this always exists, note that w 6= 0). Then there is a g ∈ SO(d) (or O(d)) that restricts to g̃ on W and is the identity on span(v) (in terms of matrices with respect to a basis of the form {v, . . .}, g would be block-diagonal, with a 1× 1 identity block and a (d− 1)× (d− 1) block for g̃).
Proof of Proposition 4. One direction is clear: the gradient is induced by the matrix of polynomials P (x) = x, x ∈ Rd, which is clearly equivariant. We have also seen that polynomials of the Laplacian are equivariant, and we can compose these equivariant PDOs to get another equivariant PDO. So what remains to show is that no other equivariant PDOs exist.
5To avoid confusion, we remark that the 3D curl of course exists but maps between two vector fields. The 2D curl is in fact closely related to the 3D curl, and the fact that it does not have a higher-dimensional analogon corresponds to the fact that the 3D curl, as a vector to vector PDO cannot easily be generalized to higher dimensions.
So let P be G-equivariant, with G as in Proposition 4. Let furthermore x ∈ Rd be arbitrary and g ∈ G be an element in the stabilizer of x, i.e. gx = x. Then we have
gP (x) = P (gx) = P (x) . (16)
In other words, for any g ∈ G with gx = x, we also have gP (x) = P (x). By the lemma, x and P (x) are thus linearly dependent, i.e. P (x) = cx for this particular x and some c ∈ R. We can apply this reasoning to any x ∈ Rd, so there is a function c : Rd → R such that P (x) = c(x)x. Furthermore,
c(x)gx = gP (x) = P (gx) = c(gx)gx , (17)
so c has to be rotation invariant. As we have already argued, this implies that c(x) = q(x21 + . . .+x 2 d) for some function q, and since c needs to be a polynomial, so does q. Then
D(P ) = q(∆) grad (18)
as claimed.
If d = 1, the lemma also holds and the proof goes through – this is just the scalar case from the previous section, which is generalized here. But forG = SO(2), this argument does not work because SO(1) = {1}, so the decisive step in the proof of the lemma fails. That is what allows the additional equivariant PDOs.
A.3 VECTOR TO SCALAR PDOS
Equivariant PDOs mapping from vector to scalar fields are simply the transpose of those mapping from scalar to vector fields (for orthogonal groups G): P is now a 1× d matrix and the steerability constraint is
P (gx) = P (x)g−1 = P (x)gT . (19)
By transposing, we get PT (gx) = gPT (x) , (20)
which is the equivariance condition for scalar to vector layers. As solutions, we get the divergence (as the transpose of the gradient) and for G = SO(2) also the 2D curl. They are again combined linearly with polynomials in the Laplacian as coefficients.
A.4 VECTOR TO VECTOR PDOS
For PDOs mapping between two vector fields, the steerability constraint is
P (gx) = gP (x)g−1 . (21)
Since the solution space in this case is somewhat more complicated, we will only cover G = SO(2) in these examples; see Appendix J for more solutions. In principle, we could apply the same method that we already used before for SO(2), choosing the radial components of P freely and using the steerability constraint to determine the angular components. But since the computations in this case are more involved and don’t yield much additional insight, we will instead use the solutions from Appendix J. Vector fields correspond to frequency 1 irreps, and writing out the solutions for those explicitly, we get that the equivariant PDOs are exactly linear combinations of(
1 0 0 1
) , ( 0 −1 1 0 ) , ( ∂2x − ∂2y 2∂x∂y 2∂x∂y ∂ 2 y − ∂2x ) , ( −2∂x∂y ∂2x − ∂2y ∂2x − ∂2y 2∂x∂y ) (22)
with polynomials in the Laplacian as coefficients.
The first two operators are simply the identity and a π2 rotation matrix, both zeroth order PDOs. Note that the rotation matrix rotates the fibers of the vector field, it does not act on the base space. The other two operators are less interpretable, but we can replace them through a change of basis:
1 2 ∆ ( 1 0 0 1 ) + 1 2 ( ∂2x − ∂2y 2∂x∂y 2∂x∂y ∂ 2 y − ∂2x ) = ( ∂2x ∂x∂y ∂x∂y ∂ 2 y ) , (23)
which is the matrix describing the composition grad ◦ div. Similarly,
−1 2 ∆ ( 0 −1 1 0 ) + 1 2 ( −2∂x∂y ∂2x − ∂2y ∂2x − ∂2y 2∂x∂y ) = ( −∂x∂y ∂2x −∂2y ∂x∂y ) , (24)
which is the matrix describing grad ◦ curl2D. So if we write R for the π2 rotation matrix (interpreted as a PDO), then the equivariant PDOs mapping between vector fields are exactly linear combinations of id, R, grad ◦ div, grad ◦ curl2D , (25) as always with polynomials in the Laplacian as coefficients.
We can be even more economical and describe these PDOs with fewer building blocks. For two vectors P,Q of one-dimensional PDOs, e.g. P = ∂ = (∂1, ∂2)T , we write P ⊗ Q for the 2 × 2 PDO with entries (P ⊗ Q)ij = PiQj . Then we can write for example grad ◦ div = ∂ ⊗ ∂. We furthermore note that
R∂ = ( 0 −1 1 0 )( ∂x ∂y ) = ( −∂y ∂x ) = curl2D . (26)
This means we can write the basis from above as
id, R, ∂ ⊗ ∂, ∂ ⊗R∂ . (27)

heading: A.5 SUMMARY OF RESULTS
text: • SO(d)- or O(d)-equivariant PDOs between two scalar fields are exactly polynomials in the Laplacian.
• G-equivariant PDOs mapping from scalar to vector fields are exactly those of the form q(∆) grad for G = SO(d) with d > 2 or G = O(d). For PDOs from vector to scalar fields, we similarly get q(∆) div.
• SO(2)-equivariant PDOs from vector to scalar fields are exactly those of the form
q1(∆) div + q2(∆) curl2D . (28)
For PDOs from scalar to vector fields, we get grad instead of div and the transpose of the 2D curl instead.
• SO(2)-equivariant PDOs between two vector fields are linear combinations of the identity, a π 2 rotation, grad ◦ div and grad ◦ curl2D, with polynomials in the Laplacian as coefficients.

heading: B REPRESENTATION THEORY PRIMER
text: This section introduces the fundamental definitions of representation theory that we need.

heading: B.1 BASIC DEFINITIONS
text: Definition 1. A group representation of a group G, or representation for short, is a group homomorphism ρ : G→ GL(V ) for some vector space V . This means that ρ(gg′) = ρ(g)ρ(g′) for all g, g′, so multiplication of group elements in G is represented as matrix multiplication in GL(V ).
For this paper, we only need V = Rc, so we focus on this case. In particular, some of the following definitions make use of the fact that we only consider finite-dimensional representations.
Given multiple representations of the same group, we can “stack them together” using direct sums:
Definition 2. Let ρ1 : G → GL(Rc1) and ρ2 : G → GL(Rc2) be two representations of G. Then we define the direct sum representation ρ1 ⊕ ρ2 : G→ GL(Rc1+c2) by
(ρ1 ⊕ ρ2)(g) := ( ρ1(g)
ρ2(g)
) , (29)
which acts independently on the subspaces Rc1 and Rc2 of Rc1+c2 .
This corresponds exactly to stacking feature fields of different types. For example, we can stack a vector and a scalar field into one four-dimensional field, such that its vector and scalar part transform independently. The representation of this four-dimensional field will be the direct sum of the vector and scalar field representations.
Often, two representations are formally different but can be transformed into one another using a change of basis; they then behave the same in all relevant aspects. For example, SO(2) has a representation
ρ(ϕ) := ( cosϕ − sinϕ sinϕ cosϕ ) (30)
that represents a rotation angle ϕ by a counterclockwise rotation matrix (this representation is the one used for vector fields). However, we could just as well use
ρ(ϕ) := ( cosϕ sinϕ − sinϕ cosϕ ) , (31)
where the rotation is clockwise. Using one or the other is pure convention and we would like to treat them as “the same” representation. This is formalized as follows: Definition 3. Two representations ρ1, ρ2 : G → GL(Rc) are equivalent if there is a matrix Q ∈ GL(Rc) such that
ρ2(g) = Q −1ρ1(g)Q (32)
for all g ∈ G.
Intuitively, ρ1 and ρ2 differ only by a change of basis, which is given by Q.

heading: B.2 DECOMPOSITION INTO IRREDUCIBLE REPRESENTATIONS
text: We can now discuss irreducible representations, which play an important role for solving the kernel and PDO steerability constraints. Definition 4. A linear subspace W ⊆ Rc is called invariant under a representation ρ : G→ GL(Rc) if ρ(g)w ∈ W for all g ∈ G and w ∈ W . In this case, we can define the restriction ρ|W : G → GL(W ), called a subrepresentation of ρ. Definition 5. A representation ρ is called irreducible if all its subrepresentations are either ρ itself or representations G→ GL({0}), where {0} is the trivial vector space.
For example, ρ1 and ρ2 are both subrepresentations of ρ1 ⊕ ρ2, so direct sums are never irreducible. A natural question is the converse: if a representation is not equivalent to a direct sum, does that mean that it is irreducible? In other words, can all representations be split into a direct sum of irreducible ones? In general, this is false, but it holds for the cases that interest us: Definition 6. A topological group is a group G equipped with a topology, such that the group multiplication G×G→ G and the inverse map G→ G are continuous with respect to that topology. A compact group is a topological group that is compact as a topological space. Theorem 6. Let G be a compact group. Then every finite-dimensional representation of G is equivalent to a direct sum of irreducible representations.
Since we only consider compact subgroups of O(2) and O(3), this theorem applies to all the cases we solve.
So we can always write
ρin = Q −1 in ⊕ i∈Iin ψi Qin
ρout = Q −1 out ⊕ i∈Iout ψi Qout (33)
where the ψi are irreducible representations. It is then easy to show (Weiler & Cesa, 2019) that a kernel k solves the G-steerability constraint
k(gx) = ρout(g) −1 k(x) ρin(g) −1 (34)
if and only if κ := Qout k Q−1in solves a block-wise steerability constraint between irreducible representations. Concretely,
κij(gx) = ψi(g) −1κij(x)ψj(g) −1 ∀i, j (35)
where κij(x) is the submatrix of κ(x) that belongs to ψi and ψj .
The approach to solving the steerability constraint is thus to solve Eq. (35) for arbitrary irreducible representations ψi and ψj . For general (not necessarily irreducible) ρout and ρin, we then first find the decompositions into irreducible representations. Each basis element κij of the solution to Eq. (35) is then padded with zeros to the right size and finally transformed via k = Q−1out κQin to get the final basis elements. See Weiler & Cesa (2019) for a more detailed discussion and visualization.
Clearly, this procedure works just as well for PDOs as it does for kernels: for PDOs, we need to find the restriction of the kernel solution space to polynomials, and it does not matter whether we restrict on the level of irreducible representations and then combine them, or first combine irreducible representations and then restrict.

heading: B.3 SPECIFIC REPRESENTATIONS
text: We now define all the types of representations that occur in the main paper:
• As already mentioned in the paper, scalar fields are described by the trivial representation ρ(g) := 1 and vector fields by the representation ρ(g) := g (for G ≤ GL(Rd)).
• For a finite group G, the regular representation is ρ : G→ GL(R|G|), defined by
ρ(g)eg′ := egg′ , (36)
where (eg)g∈G is the canonical basis of R|G| (for some ordering ofG). So this representation associates one basis vector to each group element and then acts by permuting these basis vectors. ρ(g) is thus always a permutation matrix.
• Quotient representations are a generalization of regular representations, defined as follows: let G be a finite group and H ≤ G a subgroup. Then we define the quotient representation ρ G/H quot : G→ GL(R|G|/|H|) by
ρ G/H quot (g)eg′H := egg′H , (37)
where we now use a basis indexed by the cosets gH ∈ G/H for g ∈ G. For H = {e}, we recover regular representations. Appendix C by Weiler & Cesa (2019) provides some intuition for these quotient representations in the case where G and H are both cyclic groups CN and CM .
C INTUITION FOR THE GROUP ACTION ON POLYNOMIALS
Our work makes heavy use of terms of the form p(gx), where p is a polynomial p : Rd → R, g ∈ G is a group element, and x ∈ Rd. We would now like to provide a bit more intuition for this action of the group G on polynomials. Note that we will only cover scalar-valued polynomials in this appendix, i.e. using trivial representations. The general case is straight-forward: the group acts on each polynomial in a matrix of polynomials the way we describe here, while ρin(g) and ρout(g) act via matrix multiplication.
To prevent confusion, let us reiterate that we can think of polynomials in two different ways. The first is as a formal expression, where x is a placeholder for things to be plugged in. This is the approach we take when connecting polynomials to PDOs, where we plug in differential operators for x. The second is as a specific type of function on Rd—in which case x ∈ Rd is simply the argument of that function.
This second perspective is the one in which the group action on polynomials is easiest to understand. Specifically, the action of g on x is simply matrix multiplication. Similar to how x 7→ p(x) defines a function on Rd, x 7→ p(gx) defines a different function on Rd. We could think of it as composing the function p with the group action of g on Rd.
Crucially, this new function is still a polynomial in the components of x, just with different coefficients. For purposes of illustration, consider a very simple example with d = 2 and G = SO(2). We will use the polynomial p(x) = x22 + x1 (this is just meant to be a simple non-trivial polynomial, it does not have special equivariance properties).
We can represent g as a rotation matrix parameterized by an angle θ. The action on x is then given by
gx = ( cos θ − sin θ sin θ cos θ )( x1 x2 ) = ( x1 cos θ − x2 sin θ x1 sin θ + x2 cos θ ) . (38)
The left-hand side is what we now plug into our polynomial, which yields
p(gx) = (gx)22 + (gx)1 = (x1 sin θ + x2 cos θ) 2 + (x1 cos θ − x2 sin θ) . (39)
We can expand this and collect the coefficients for each power of x1 and x2, which yields
p(gx) = sin2(θ)x21 + 2 sin(θ) cos(θ)x1x2 + cos 2(θ)x22 + cos(θ)x1 − sin(θ)x2 . (40)
If desired, we can now switch to the first perspective, and interpret this polynomial a formal expression defined by its coefficients. In that perspective, g acts on p by modifying its coefficients, rather than by composition. Computing the new coefficients is straightforward analytically (given a matrix representation of g); we just gave a simple example for a polynomial of order two. The general case proceeds along the same lines, it just requires using the binomial theorem (or the multinomial theorem for d > 2).

heading: D BACKGROUND ON DISTRIBUTIONS
text: In Appendix E, we will describe our framework for equivariant maps using convolutions with Schwartz distributions. To facilitate that, we now give the necessary background on Schwartz distributions. We restrict ourselves to what is absolutely necessary for our purposes; for a much more thorough introduction and for proofs, see e.g. Treves (1967).

heading: D.1 BASIC DEFINITIONS
text: Definition 7. For U ⊂ Rd open,D(U) := C∞c (U) is called the space of test functions on U (C∞c (U) is the space of compactly supported smooth functions U → R). A sequence (ϕn) in D(U) is defined to converge to 0 iff
(i) there is a compact subset K ⊂ U such that the support of each ϕn is contained in K and
(ii) ∂αϕn → 0 uniformly for all multi-indices α.
One can define the so-called canonical LF topology on D(U). This topology induces the notion of convergence just given. However, constructing this topology explicitly is unnecessary for our purposes since knowing when sequences converge will be enough.
Definition 8. A linear functional T : D(U)→ R is defined as continuous if for every sequence ϕn that converges to 0 in D(U), Tϕn → 0 (in R). Such a continuous functional is called a distribution on U . The space of all distributions on U is written as D′(U).
This notion of continuity is also induced by the canonical LF topology and then D′(U) is the topological dual of D(U) as the notation suggests. So intuitively, a distribution is a “reasonable” way of linearly assigning a number in R to each test function in C∞c .
A function f : U → R induces a distribution Tf ∈ D′(U), defined by
Tfϕ := ∫ U fϕdλd . (41)
Since ϕ has compact support, we don’t need many restrictions on f (for example, any locally integrable function f ∈ L1,loc(U) works). We will also use the duality pairing 〈T, ϕ〉 := T (ϕ) (42) and under slight abuse of notation also write
〈f, ϕ〉 := Tf (ϕ) = ∫ U fϕdλd . (43)
Note that this coincides with the inner product on L2(U), but it allows functions f that are not in L2 while in exchange requiring ϕ to have compact support.

heading: D.2 CONVOLUTIONS
text: We define the translation τx : Rd → Rd by τx(y) = y + x and set τxf := f ◦ τ−x for functions f on Rd. This is exactly the same as the action of (Rd,+) we used in the main paper, just with more explicit notation. Furthermore, we write f̌(x) := f(−x). Then we can define the convolution between a distribution and a test function: for f ∈ D(U) and T ∈ D′(U), the convolution T ∗ f ∈ C∞(U) is defined by
(T ∗ f)(x) := T (τxf̌) . (44)
Explicitly, (τxf̌)(y) = f(x− y). This immediately shows that this notion of convolution extends the classical one, i.e.
Tg ∗ f = g ∗ f . (45)

heading: D.3 DERIVATIVES
text: Definition 9. For T ∈ D′(U), we define the distributional derivative ∂αT as the unique distribution on U for which
〈∂αT, ϕ〉 = (−1)|α|〈T, ∂αϕ〉 . (46)
The distributional derivative of all orders always exists, i.e. “distributions are infinitely differentiable”, but of course only in this distributional sense. At least for f ∈ C∞(U) (but also under much weaker assumptions), this definition also extends the definition of derivatives of functions, in the sense that T∂αf = ∂ αTf .

heading: D.4 COMPOSITION WITH DIFFEOMORPHISMS
text: LetF : U → U be a diffeomorphism and T ∈ D′(U). Then we define the composition T◦F ∈ D′(U) as
〈T ◦ F,ϕ〉 := 〈T, ∣∣detDF−1∣∣ϕ ◦ F−1〉 (47)
where detDF−1 is the inverse Jacobian.
As in the previous constructions, this extends the definition for classical functions in the sense that
Tf◦F = Tf ◦ F . (48)
This follows immediately from the transformation theorem for integrals: 〈Tf◦F , ϕ〉 = ∫ U (f ◦ F ) · ϕdλ (49)
= ∫ U f · (ϕ ◦ F−1) ∣∣detDF−1∣∣ dλ (50)
= 〈Tf , (ϕ ◦ F−1) ∣∣detDF−1∣∣〉 (51)
= 〈Tf ◦ F,ϕ〉 . (52)
Furthermore, this type of composition is associative: 〈(T ◦ F ) ◦G,ϕ〉 = 〈T ◦ F, ∣∣detDG−1∣∣ϕ ◦G−1〉 (53)
= 〈T, ∣∣detDF−1∣∣∣∣detDG−1 ◦ F−1∣∣ϕ ◦G−1 ◦ F−1〉 (54)
= 〈T, ∣∣det(DG−1 ◦ F−1)DF−1∣∣ϕ ◦ (F ◦G)−1〉 (55)
= 〈T, ∣∣detD(F ◦G)−1∣∣ϕ ◦ (F ◦G)−1〉 (56)
= 〈T ◦ (F ◦G), ϕ〉 . (57)
We are particularly interested in the case where F is a linear transformation, i.e. F ∈ GL(Rd). Then we have 〈T ◦ F,ϕ〉 = ∣∣detF−1∣∣〈T, ϕ ◦ F−1〉 (58)
(note that in this case, we can pull out the determinant because it is just a constant). For translations, we get
〈T ◦ τx, ϕ〉 = 〈T, ϕ ◦ τ−x〉 = 〈T, τxϕ〉 . (59)

heading: D.5 THE DIRAC DELTA DISTRIBUTION
text: A very simple but important distribution is the following: Definition 10. For any x ∈ Rd, the Dirac delta distribution δx ∈ D(Rd) is defined by
δx[ϕ] := ϕ(x) . (60)
It is clear from the definition that the distributional derivatives of the delta distribution are given by
(∂αδx)[ϕ] = (−1)|α|∂αϕ(x) . (61)

heading: D.6 CONVERGENCE OF DISTRIBUTIONS
text: Definition 11. We say that a sequence Tn of distributions converges to T , written Tn → T , if it converges pointwise, i.e. 〈Tn, ϕ〉 → 〈T, ϕ〉 for n→∞ (62) for all ϕ ∈ D(U).
Abusing notation a bit, we also write fn → T for functions fn if Tfn → T . One important type of function sequences are Dirac sequences, which converge to the Delta distribution: Lemma 7. Let fn be a sequence in L1(Rd) such that
(i) fn ≥ 0,
(ii) ‖fn‖1 = 1, and (iii) ∫ Rd\Bε(0) fndλ d → 0 for all ε > 0.

heading: Then fn → δ0 in the sense of distributions.
text: We also note that convergence plays together nicely with derivatives and with convolutions: Lemma 8. If Tn → T , then ∂αTn → ∂αT for all α.
Proof.
〈∂αTn, ϕ〉 = (−1)|α|〈Tn, ∂αϕ〉 → (−1)|α|〈T, ∂αϕ〉 = 〈∂αT, ϕ〉 .
(63)
Lemma 9. If Tn → T , then Tn ∗ f → T ∗ f pointwise for all f ∈ D(U).
Proof.
(Tn ∗ f)(x) = Tn(τxf̌) → T (τxf̌) = (T ∗ f)(x) .
(64)

heading: D.7 TEMPERED DISTRIBUTIONS AND THE FOURIER TRANSFORM
text: The Schwartz space S(Rd) is the space of functions for which all derivatives decay very quickly as |x| → ∞. More precisely, a smooth function f : Rd → R is in S(Rd) iff
sup x∈Rd ∥∥xβ∂αf(x)∥∥ <∞ . (65) Intuitively, all derivatives of f must decay more quickly than any polynomial. Examples are compactly supported functions or Gaussians.
Its dual space S ′(Rd) is called the space of tempered distributions and can be continuously embedded into D′(Rd). Tempered distributions are still very general; in particular the delta distribution, distributions induced by functions, and derivatives of tempered distributions are all tempered.
The reason we’re interested in tempered distributions is that there is a Fourier transform defined by
〈F {T} , ϕ〉 := 〈T,F {ϕ}〉 (66)
for tempered distributions T . This is an automorphism on S ′(Rd). We use the following convention for the Fourier transform on functions:
F {f} (ξ) := (2π)−d/2 ∫ Rd f(x) exp(−ix · ξ) dx , (67)
which means the inverse is given by F−1 {f} (x) := (2π)−d/2 ∫ Rd f(x) exp(ix · ξ) dξ . (68)
We will later need the Fourier transform of derivatives of the Dirac delta distribution:
〈F {∂αδ0} , ϕ〉 = 〈∂αδ0,F {ϕ}〉 = (−1)|α|∂αF {ϕ} ∣∣∣ 0
= (−i)|α|F {xαϕ} ∣∣∣ 0
= (−i)|α|(2π)−d/2 ∫ xαϕ(x) exp(−ix · 0)dx
= (−i)|α|(2π)−d/2 ∫ xαϕ(x)dx ,
(69)
so F {∂αδ0} ∝ xα (or more precisely, the distribution induced by the function x 7→ xα).

heading: D.8 MATRICES OF DISTRIBUTIONS
text: Analogously to how we defined matrices of PDOs or of polynomials, we will need matrices of distributions. We will write D′(U,Rcout×cin) for the space of cout × cin dimensional matrices with entries in D′(U). We get a pairing
〈·, ·〉 : D′(U,Rcout×cin)×D(U,Rcin)→ Rcout (70)
defined by 〈T, f〉i := ∑ j 〈Tij , fj〉 (71)
The convolution of T ∈ D′(U,Rcout×cin) and f ∈ D(U,Rcin) is defined analogously to the scalar case, i.e. (T ∗ f)(x) := 〈T, τxf̌〉, where we now use the more general duality pairing just defined. We can also define the multiplication of matrices of distributions by matrices over R. For a matrix A ∈ Rcout×cout and a matrix of distributions T ∈ D′(Rd,Rcout×cin), we write
(AT )ij := ∑ l AilTlj ∈ D′(Rd,Rcout×cout) (72)
(multiplying a distribution by a scalar obviously defines another distribution). Analogously, we define TB for B ∈ Rcin×cin . It immediately follows that
A〈T,Bϕ〉 = 〈ATB,ϕ〉 (73)
holds.
Composition with diffeomorphisms of Rd can simply be defined component-wise for matrices of distributions. This commutes with multiplication by constant matrices, meaning that
ATB ◦ F = A(T ◦ F )B . (74)
Finally, we also define the Fourier transform component-wise, and this commutes with multiplication by matrices in the same way.

heading: E DISTRIBUTIONAL FRAMEWORK FOR EQUIVARIANT MAPS
text: In this section, we develop two main results: first that all linear continuous translation equivariant maps between feature spaces are convolutions with some distribution, and then the equivariance constraint for such convolutions. See Appendix D for background on distributions. We equip D(U) with the canonical LF topology throughout this section.

heading: E.1 TRANSLATION EQUIVARIANT MAPS ARE CONVOLUTIONS WITH DISTRIBUTIONS
text: We begin by showing that the framework using convolutions with Schwartz distributions encompasses all translation equivariant continuous linear maps between feature spaces. As preparation, we prove a simple Lemma on the reflection map s:
Lemma 10. The map s : D(U)→ D(U) given by s(f) := f̌ is linear and continuous.
Proof. Linearity is clear:
s(af + bg)(x) = af(−x) + bg(−x) = as(f)(x) + bs(g)(x) . (75)
For continuity, we use the fact that a linear map from D(U) to itself is continuous if and only if it is sequentially continuous (Treves, 1967, Proposition 14.7). So take any sequence fn → 0 in D(U). This means that
(i) there is a compact subset K ⊂ U such that the support of each fn is contained in K
(ii) ∂αfn → 0 uniformly for all multi-indices α
We set −K := {−x |x ∈ K}, then the support of s(fn) (which is just the mirror of the support of fn) is contained in −K, and −K is compact. Additionally,
∂αs(fn) = (−1)αs (∂αfn) . (76)
Since ∂αfn → 0 uniformly, the same is true for s (∂αfn). It follows that s is sequentially continuous, and thus continuous, which concludes the proof.
We first prove our main generality result for the special case of one-dimensional fibers: Lemma 11. Let Φ : D(U)→ D(U) be a translation equivariant continuous linear map. Then there is a distribution T ∈ D′(U) such that Φ(f) = T ∗ f .
Proof. Let Φ : D(U)→ D(U) be any continuous linear map, where continuity is understood with respect to the canonical LF topology. We then define
T : D(U)→ R, T (f) := Φ(f̌)(0) . (77) Equivalently, we can write T = δ0 ◦ Φ ◦ s, understood as normal composition of functions. δ0 and Φ are linear and continuous, as is s by Lemma 10. Therefore, T is also linear and continuous, and thus a Schwartz distribution T ∈ D′(U). We will also need that (
~τ−xf ) (y) = (τ−xf)(−y)
= f(x− y) = f̌(y − x) = (τxf̌)(y) .
(78)
Now using the assumption that Φ is translation equivariant, i.e.
τx ◦ Φ = Φ ◦ τx , (79) it follows that convolution with T is given by
(T ∗ f)(x) = 〈T, τxf̌〉
= 〈T, ~τ−xf〉 = Φ(τ−xf)(0)
= τ−xΦ(f)(0)
= Φ(f)(x) .
(80)
This shows that convolution with T is Φ, which concludes the proof.
Finally, we generalize to multi-dimensional feature fields: Theorem 12. Let Φ : D(U,Rcin)→ D(U,Rcout) be a translation equivariant continuous linear map. Then there is a distribution T ∈ D′(U,Rcout×cin) such that Φ(f) = T ∗ f .
Proof. We write ej for the j-th canonical basis vector of Rcin . For f ∈ D(U,Rcin), we write fj for its j-th component. Then we have
Φ(f) = Φ  cin∑ j=1 fjej  = cin∑ j=1 Φ (fjej) . (81)
We now define the components Φij : D(U)→ D(U) by Φij(g) = (Φ(gej))i . (82) Here, gej ∈ D(U,Rcin) is the map with g in its j-th output component and 0 in the other ones. We can then write
Φ(f)i = cin∑ j=1 Φ(fjej)i = cin∑ j=1 Φij(fj) . (83)
It is clear that the components Φij of Φ are still translation equivariant, continuous and linear. Therefore, by the previous Lemma, there are distributions Tij ∈ D′(U) such that Φij(g) = Tij ∗ g. We then get
Φ(f)i = cin∑ j=1 Tij ∗ fj = T ∗ f (84)
where T ∈ D′(U,Rcout×cin) is the matrix of distributions with entries Tij .

heading: E.2 EQUIVARIANCE CONSTRAINT FOR DISTRIBUTIONS
text: We will now characterize the distributions T ∈ D′(Rd,Rcout×cin) for which the operator
T∗ : D(Rd,Rcin)→ C∞(Rd,Rcout) (85)
given by f 7→ T ∗ f is H-equivariant. The codomain of T∗ may contain functions that are not compactly supported for some T . We are mostly interested in distributions T for which the codomain can be restricted to D(Rd,Rcout) but for the discussion of equivariance this does not make any difference, so we keep the derivation general by not restricting the distribution T .
First, we can note that convolution with a distribution is always translation equivariant: Proposition 13. The map f 7→ T ∗ f is translation equivariant for any distribution T .
Proof. The general equivariance condition is
(h · (T ∗ f))(x) = (T ∗ (h · f))(x) (86)
for all f ∈ D(Rd,Rcin), h ∈ H,x ∈ Rd. More explicitly, this means
ρout(g)〈T, τh−1xf̌〉 = 〈T, τx(ρin(g) f ◦ h−1)〉 , (87)
where g is the linear component of h.
Let h = t ∈ Rd, i.e. a pure translation. The condition then becomes
〈T, τ(x−t)f̌〉 = 〈T, τx (f ◦ τ−t)〉 = 〈T, τxτ−tf̌〉 , (88)
which always holds (since τ(x−t) = τxτ−t). So T∗ is translation equivariant by construction.
Therefore, it suffices to consider pure linear transformations g ∈ GL(Rd). For those, we have the following result, which generalizes the steerability constraint for PDOs and the one for kernels: Theorem 14. Let T ∈ D′(Rd,Rcout×cin). Then the map f 7→ T ∗ f is G-equivariant if and only if
T ◦ g = |det g|−1ρout(g)Tρin(g)−1 . (89)
As is shown in Appendix D, the definitions of convolution with distributions and of composition with a diffeomorphism are compatible with those for classical functions. So if T is a classical kernel, then the constraint on T is given by the same equation, which already gives us the steerability constraint for kernels.
Proof. We start with the constraint Eq. (87) from above with h = g ∈ G and will transform this into the desired form Eq. (89). First, to simplify notation, we can remove all the reflections ·̌ in Eq. (87). This is now possible because f ◦ g−1 = f̌ ◦ g−1, which was not true for translations. Furthermore, we can pull the translations to the other side of the duality pairing:
ρout(g)〈T ◦ τg−1x, f〉 ! = 〈T ◦ τx, ρin(g)f ◦ g−1〉 = |det g|〈T ◦ τx ◦ g, ρin(g)f〉 .
(90)
Now we’ll use τx ◦ g = g ◦ τg−1x , (91)
since x+ gy = g(g−1x+ y) . (92)
Plugging this into Eq. (90), we get
ρout(g)〈T ◦ τg−1x, f〉 ! = |det g|〈T ◦ g ◦ τg−1x, ρin(g)f〉 . (93)
We want to have only f on the right side of the duality pairing, so we use the notation we introduced for multiplying distributions by matrices and get
〈ρout(g)T ◦ τg−1x, f〉 ! = |det g|〈T ◦ g ◦ τg−1xρin(g), f〉 . (94)
This holds for all f iff we have equality of distributions,
ρout(g)T ◦ τg−1x ! = |det g|T ◦ g ◦ τg−1xρin(g) . (95)
Finally, multiplication with matrices and composition with diffeomorphisms commute, so we can cancel the τg−1x. This means our final constraint on T is
ρout(g)T ! = |det g|(T ◦ g)ρin(g) , (96)
which we can slightly rewrite as
T ◦ g != |det g|−1ρout(g)Tρin(g)−1 . (97) This is exactly Eq. (89). All steps of the derivation work equally well in the other direction, which proves the “if and only if”.

heading: E.3 DIFFERENTIAL OPERATORS AS CONVOLUTIONS
text: It is clear that convolutions with classical kernels are a special case of convolutions with distributions (see Appendix D). But we have claimed that convolutions with distributions also cover PDOs, which is what we show now.
We have already seen in Appendix D that derivatives of the delta distribution are closely related to PDOs. So we calculate the convolution with such derivatives:
((∂αδ0) ∗ f)(x) = 〈∂αδ0, τxf̌〉 = (−1)|α|∂α(τxf̌)(0) = (−1)|α|τx(∂αf̌)(0)
= τx
( }∂αf ) (0)
= ( }∂αf ) (−x)
= ∂αf(x) .
(98)
Noting that the map T 7→ T∗ is R-linear, we see that D(p)f = p(∂)δ0 ∗ f for polynomials p. It then also immediately follows that this is true for matrices of polynomials P because
(D(P )f)i = cin∑ j=1 D(Pij)fj
= cin∑ j=1 Pij(∂)δ0 ∗ fj
= P (∂)δ0 ∗ f .
(99)
This shows that PDOs can be interpreted as convolutions with distributions, as claimed.

heading: E.4 THE FOURIER DUALITY BETWEEN KERNELS AND PDOS
text: Interpreting PDOs as convolutions with distributions allows us to relate them to classical convolutional kernels via a Fourier transform. We have already seen in Appendix D that
F {∂αδ0} = (2π)−d/2(−i)|α|xα , (100) which by linearity of the Fourier transform immediately implies
F {P (∂)δ0}ij = (2π) −d/2 ∑ α (−i)|α|cijα xα . (101)
for Pij = ∑ α c ij α x
α. So the Fourier transform of derivatives of the delta distribution are polynomials, and of course vice versa via the inverse Fourier transform. Since PDOs are convolutions with such delta distribution derivatives, we can also interpret them as convolutions with the (inverse) Fourier transform of polynomials. We will use this interpretation to give a derivation of the PDO steerability constraint that sheds some light on the similarities and differences to the kernel steerability constraint.
We begin by proving a basic fact about the Fourier transform of a composition of functions:
Lemma 15. Let g ∈ GL(Rd) and ϕ ∈ D(U). Then F {ϕ ◦ g} = ∣∣det g−1∣∣F {ϕ} ◦ g−T , (102) where we use the shorthand g−T := ( g−1 )T .
Proof. First we apply the transformation theorem for integrals: F {ϕ ◦ g} (ξ) = (2π)−d/2 ∫ ϕ(gx) exp(−ix · ξ)dx
= ∣∣det g−1∣∣(2π)−d/2 ∫ ϕ(x) exp(−ig−1x · ξ)dx . (103)
Now we just rewrite g−1x · ξ = x · g−T ξ, which gives the desired result.
This result holds more generally for tempered distributions (a subset of distributions for which the Fourier transform can be defined, see Appendix D):
Proposition 16. For a tempered distribution T and g ∈ GL(Rd), F {T ◦ g} = ∣∣det g−1∣∣F {T} ◦ g−T . (104)
Proof.
〈F {T ◦ g} , ϕ〉 = 〈T ◦ g,F {ϕ}〉 = ∣∣det g−1∣∣〈T,F {ϕ} ◦ g−1〉
(1) = 〈T,F { ϕ ◦ gT } 〉 = 〈F {T} , ϕ ◦ gT 〉 = ∣∣det g−1∣∣〈F {T} ◦ g−T , ϕ〉 ,
(105)
where we used Lemma 15 for (1).
Now note that the equivariance constraint for distributions, Eq. (89), is equivalent to the constraint we get when we take the Fourier transform on both sides. That’s because the Fourier transform is an automorphism on the space of tempered distributions. By applying the result we just proved, we then get the constraint ∣∣det g−1∣∣F {T} ◦ g−T = ∣∣det g−1∣∣ρout(g)F {T} ρin(g)−1 . (106) We can cancel the determinants, which gives the equivariance constraint in Fourier space:
F {T} ◦ g−T = ρout(g)F {T} ρin(g)−1 . (107)
Now let T = P (∂)δ0, so that T ∗ f = D(P )f for some matrix of polynomials P . We’ve seen above that in this case F {T} = P ◦m−i (up to a constant coefficient), where m−i is multiplication by the negative imaginary unit −i. So D(P ) is equivariant iff
P ◦m−i ◦ g−T = ρout(g)P ◦m−iρin(g)−1 , (108)
but since m−i is invertible and commutes with the other maps, we can cancel it and get the equivariance condition P ◦ g−T = ρout(g)Pρin(g)−1 , (109) which is precisely the PDO steerability constraint. The reason that it differs slightly from the kernel steerability constraint can now be traced back to Lemma 15. Intuitively speaking, since PDOs are in a sense the Fourier transform of convolutional kernels, they transform differently under GL(Rd), which leads to superficial differences in the steerability constraints. However, Fourier transforms commute with rotations and reflections (i.e. transformations from O(d)), which is why for G ≤ O(d), the two steerability constraints coincide.

heading: E.5 PDOS AS THE INFINITESIMAL LIMIT OF KERNELS
text: Let ψn be any sequence of functions such that ψn → δ0, for example a Dirac sequence. Then for a polynomial p = ∑ α cαx α, we also have
p(∂)ψn → p(∂)δ0 (110)
because of Lemma 8. Then Lemma 9 implies that
p(∂)ψn ∗ f → p(∂)δ0 ∗ f = D(p)f , (111)
where the convergence is understood pointwise. Therefore, any sequence of kernels that approximates the delta distribution (by becoming “increasingly narrow”) can be used to approximate arbitrary PDOs by convolving with the derivatives of the kernels.
One example is a sequence of Gaussians
ψε(x) := 1
(2πε)d/2 e−x 2/ε (112)
(now indexed by ε > 0 instead of natural numbers). For ε → 0, we have ψε → δ0, as is easy to check with Lemma 7. This naturally leads to the “derivative of Gaussian” discretization used in our experiments: we can approximate the PDO as convolution with a derivative of a Gaussian kernel and then simply discretize this Gaussian derivative by sampling it on the grid points.
The discussion in this subsection focused on 1× 1 PDOs and kernels, i.e. cin = cout = 1, but since convergence for multi-dimensional PDOs and kernels works component-wise, everything generalizes immediately to that setting.

heading: F PROOF OF THE PDO EQUIVARIANCE CONSTRAINT
text: In this section, we prove Proposition 1 and Theorem 2, the characterizations of translation- and G-equivariance for PDOs. Appendix E already contains a proof of Theorem 2 that is arguably more insightful, but in this section we present an elementary proof that does not rely on Schwartz distributions.

heading: F.1 TRANSLATION EQUIVARIANCE
text: We begin by proving that translation equivariance of a PDO is equivalent to spatially constant coefficients, Proposition 1.
Let Φij = ∑ α c ij α ∂
α, where cijα ∈ C∞(Rd) and i, j index the cout × cin matrix describing the PDO. The equivariance condition Eq. (2) for the special case of translations is
Φ(tBin f) = tBout Φ(f) ∀t ∈ Rd, f ∈ Fin . (113)
Using the definition of the B action in Eq. (1), this becomes∑ α cα(x)∂ αf(x− t) = ∑ α cα(x− t)∂αf(x− t) ∀x, t ∈ Rd, f ∈ Fin . (114)
Here, cα are matrix-valued, f is vector valued, and we have a matrix-vector product between the two. Explicitly, this means∑
j ∑ α cijα (x)∂ αfj(x− t) = ∑ j ∑ α cijα (x− t)∂αfj(x− t) (115)
for all indices i. If f is zero in all but one component, the sum over j reduces to one summand. We therefore get a simpler scalar constraint∑
α
cijα (x)∂ αg(x− t) = ∑ α cijα (x− t)∂αg(x− t) (116)
that has to hold for all indices i, j, where g ∈ C∞(Rd) is now scalar-valued.
Since this must hold for all functions g, we get an equality of differential operators,∑ α cijα (x)∂ α = ∑ α cijα (x− t)∂α ∀x, t ∈ Rd . (117)
This implies cijα (x− t) = cijα (x) for all x, t and α, which in turn means that cijα must be constants. From Eq. (114), it is also apparent that constant coefficients are sufficient for translation equivariance, which proves the converse direction.
F.2 G-EQUIVARIANCE
Secondly we prove Theorem 2, theG-steerability constraint for PDOs. Recall that theG-equivariance condition for D(P ) is
D(P )(g Bin f) = g Bout (D(P )f) , (118)
where f ∈ Fin and g ∈ G. If we write out the definition of the induced action B, this becomes
D(P )(ρin(g)f ◦ g−1) = ρout(g)(DP f) ◦ g−1 . (119) We write P = ∑ α cαx
α, where cα are constant matrix-valued coefficients cα ∈ Rcout×cin . The LHS of Eq. (119) is then
D(P )(ρin(g)f ◦ g−1) = ∑ α cαρin(g)∂ α(f ◦ g−1) . (120)
In Lemma 17, we will show that the last term in Eq. (120) is given by ∂α(f ◦ g−1) = (( g−T∂ )α f ) ◦ g−1 , (121) where g−T := ( g−1 )T and ∂ is understood as a vector that g−T acts on. Plugging this into Eq. (120), we get
D(P )(ρin(g)f ◦ g−1) = ∑ α cα (( g−T∂ )α ρin(g)f ) ◦ g−1 (122)
for the LHS of Eq. (119). For the matrix of polynomials P = ∑ α cαx α, we now define
g · P := ∑ α cα(g −1x)α , (123)
which is again a matrix of polynomials, each one rotated by g. Then Eq. (122) can be written more compactly as D(P )(ρin(g)f ◦ g−1) = (D(gT · P )ρin(g)f) ◦ g−1 . (124) Plugging this back into Eq. (119), canceling the g−1 and using the fact that this has to hold for all f , we get
D(gT · P )ρin(g) = ρout(g)D(P ) (125) as an equality of differential operators. We move the ρin(g) to the other side and use the fact that the map D is bijective, which yields our final constraint
P (g−Tx) = ρout(g)P (x)ρin(g) −1 . (126)
This concludes the proof of the G-steerability constraint for PDOs.
Finally, we prove the Lemma we just made use of, a higher-dimensional chain rule for the special case we need:
Lemma 17. Let f ∈ C∞(Rd,Rc) and g ∈ GL(Rd). Then for any multi-index α,
∂α(f ◦ g−1) = (( g−T∂ )α f ) ◦ g−1 , (127)
where g−T := ( g−1 )T .
Proof. In general, for a linear map A ∈ GL(Rd), we have
∂iAj(x) = Aji ∀x . (128)
Therefore,
∂i(f ◦ g−1) = ∑ j (∂jf) ◦ g−1∂i(g−1)j (129)
= ∑ j (∂jf) ◦ g−1 · (g−1)ji (130)
= ((( g−1 )T ∂ ) i f ) ◦ g−1 . (131)
We can apply this iteratively to show that ∂α(f ◦ g−1) = ((( g−1 )T ∂ )α f ) ◦ g−1 . (132)

heading: G TRANSFERRING STEERABLE KERNEL BASES TO STEERABLE PDO BASES
text: In this section, we develop the method presented in Section 3.2 in more detail and with proofs.
We fix a group G ≤ O(d) and representations ρin and ρout. Then we write K for the space of G-steerable kernels. Because the steerability constraints for kernels and PDOs are identical in this setting, the space of equivariant PDOs is the image under the isomorphism D of the intersection
Kpol := R[x1, . . . , xd] ∩ K . (133)
In words, the space of equivariant PDOs is isomorphic to the space of polynomial steerable kernels. Both spaces are infinite-dimensional real vector spaces. The question we tackle now is how we can find a basis of Kpol given a basis of K under certain conditions. It will vastly simplify our discussion to treat K and Kpol as modules over invariant kernels instead of as real vector spaces, at least for now. A module is a generalization of a vector space, where the scalars for scalar multiplication can form a ring instead of a field. Because the radial parts of steerable kernels are unrestricted, it makes sense to think of K as a module over the ring of radial functions. Formally: Lemma 18. K is a C∞(R≥0)-module, with scalar multiplication defined by
(fκ)(x) := f ( |x|2 ) κ(x) (134)
for κ ∈ K and f ∈ C∞(R≥0). Similarly, Kpol is an R[r2]-module, where r2 := x21 + . . .+ x2d and multiplication is simply multiplication of polynomials.
Proof. The steerability constraint
κ(gx) = ρout(g)κ(x)ρin(g) −1 (135)
is clearly R-linear and in particular, K is closed under addition and 0 ∈ K. Furthermore, for κ ∈ K,
(fκ)(gx) = f(|gx|2)κ(gx) = f(|x|2)ρout(g)κ(x)ρin(g)−1
= ρout(g)(fκ)(x)ρin(g) −1 ,
(136)
so K is also closed under the given scalar multiplication. The proof for Kpol is exactly analogous.
In Lemma 18 and in the following, we write r2 instead of |x|2 simply to emphasize its role as a polynomial; so when reading r2, think of it as a polynomial, and when reading |x|2 simply as a function of x.
A basis of a module is defined analogously to a basis of a vector space, as a set of linearly independent vectors that span the entire module. However, linear combinations now allow coefficients in the ring of radial functions, instead of only real numbers. This means that fewer vectors are needed to span the entire space, because the coefficients “do more work”.
In contrast to vector spaces, not every module has a basis. K and Kpol do have a basis in the cases we consider in the paper but for this section that doesn’t matter to us: we will simply assume that K has a basis with certain properties and then transfer this basis to Kpol. That the method developed here is indeed applicable to subgroups of O(2) and O(3) will be the topic of Appendix I.
We roughly proceed in two steps:
1. We show that a basis of K (over C∞(R≥0)) that consists only of polynomials (and fulfills a few other technical conditions) is also a basis of Kpol, but this time of course over R[r2].
2. We then show how to turn this module basis into a vector space basis of Kpol.
The first step is formalized as follows: Proposition 19. Let B ⊂ R[x1, . . . , xd]cout×cin be a basis of K such that no matrix of polynomials in B is divisible by r2 and each one is homogeneous. Then B is also a basis of Kpol as an R[r2]-module.
Here, we say that b ∈ R[x1, . . . , xd]cout×cin is divisible by r2 if every component is divisible by r2 (as a polynomial). We call it homogeneous if all its entries are homogeneous polynomials of the same degree.
Proof. Linear independence is obvious: R[r2] ⊂ C∞(R≥0), so if B is linearly independent over C∞(R≥0), then it is also linearly independent over R[r2].
To show that B generates Kpol, first let p ∈ Kpol be homogeneous of degree l. Because p ∈ K, there are fi ∈ C∞(R≥0) and κi ∈ B such that
p(x) = ∑ i fi(|x|2)κi(x) . (137)
We want to show that fi ∈ R[z], i.e. that fi is a polynomial. Since each κi is homogeneous of degree li, we get ∑
i
fi(λ 2|x|2)λliκi(x) = λl ∑ i fi(|x|2)κi(x) (138)
because of p(λx) = λlp(x) for any λ ≥ 0. Because the κi are linearly independent, we must have
fi(λ 2|x|2) = λl−lifi(|x|2) . (139)
Thus, fi is homogeneous of degree l−li2 and as we show in Lemma 20 below, this implies
fi(z) = cz (l−li)/2 , (140)
or alternatively fi(|x|2) = ci|x|l−li . (141) What remains to show is that l−li2 is a natural number, i.e. that l − li is even and non-negative. To prove this, we divide all the i into two groups: those for which l − li is even and those for which it is odd. Then we get an expression of the form
p(x) = |x| ∑ k akr 2mkκk + ∑ k bkr 2nkκk . (142)
Each summand is a rational function (the quotient of two polynomials), so both sums are rational functions. p(x) is also rational, so the first term on the RHS has to be rational. This is only possible if ak = 0 for all k, otherwise we could divide by the (rational) sum and should get a rational function, but |x| is not rational. This shows that l − li is even for all i. Now let lmax := maxi li. Then
p(x) =
∑ i ci|x| lmax−liκi(x)
|x|lmax−l . (143)
It is not possible to cancel any terms: because κi is not divisible by r2 and for one i, lmax − li = 0, the enumerator is not divisible by r2 (as a polynomial). Since the denominator is a power of r2, the fraction can’t be simplified. But we know that p(x) is a polynomial. Therefore, lmax ≤ l. In summary, we’ve shown that
fi(z) = ciz (l−li)/2 (144)
where l − li is even and non-negative. Therefore, all fi are polynomials. Now recall that we assumed p to be homogeneous. But this is no significant restriction: we can write any polynomial as a sum of homogeneous polynomials, each of which can be written as a linear combination of the κi, as we just showed. Adding those up leads to a linear combination of the κi for arbitrary polynomials. This complete the proof.
Lemma 20. Let f : R≥0 → R be homogeneous of degree l ∈ R, meaning that f(λx) = λlx for all λ ≥ 0. Then f(z) = czl for some c ∈ R (and if l ∈ N, then f is a polynomial).
Note that this Lemma does not hold in higher dimensions – in general there are many more homogeneous functions than polynomials!
Proof. For any z ≥ 0, f(z) = zlf(1) , (145)
which proves the claim by setting c := f(1).
In Appendix I, we will show that the construction of angular basis elements described in Section 3.2 leads to a basis B of K with the properties required by Proposition 19. It then follows that this also defines a basis of Kpol as a module. We now come to the second step, turning this module basis into a vector space basis.
Proposition 21. Let B be a basis of Kpol as an R[r2]-module. Then{ r2kb ∣∣ k ∈ N≥0, b ∈ B} (146) is a basis of Kpol as a real vector space.
Proof. Let v ∈ Kpol. By assumption, there are then basis vectors b1, . . . , bn ∈ B and coefficients p1, . . . , pn ∈ R[r2] such that
v = n∑ i=1 pibi . (147)
For each i = 1, . . . , n there are also real coefficients a(i)0 , . . . , a (i) mi such that
pi = mi∑ k=0 a (i) k r 2k . (148)
Combining these equations, we get
v = n∑ i=1 mi∑ k=0 a (i) k r 2kbi . (149)
This is a real linear combination of elements of the form r2kb for b ∈ B. So since v was arbitrary, the set of such elements does indeed span Kpol. To prove linear independence, let
n∑ i=0 air 2kibi = 0 (150)
for some choice of real coefficients ai. Now we only need to note that air2ki ∈ R[r2]. So we can interpret this expression as a linear combination over R[r2]. Because B is a basis, it follows that air 2ki = 0 for all i, and thus ai = 0. That proves linear independence.
As a final note, we show that the condition in Proposition 19 that the basis elements are not divisible by r2 is purely technical; any basis can easily be transformed into one that fulfills it in a canonical way, as formalized by the following lemma. We do not formally need this result anywhere but it may prove useful if this approach is extended to other groups G because it clarifies which parts of the conditions in Proposition 19 are actually important. Lemma 22. Let B ⊂ R[x1, . . . , xd]cout×cin be a basis of K. For b ∈ B, we write b = r2kb′, where k is chosen maximally such that there is an (automatically unique) polynomial matrix b′ ∈ R[x1, . . . , xd]cout×cin . Then B′ := {b′ | b ∈ B} is a basis of K and no b′ is divisible by r2. Furthermore, if b is homogeneous, then so is b′.
Proof. First, note that b′ is in fact well-defined: For k = 0, b′ = b works, while for k > deg(b)2 , no fitting b′ exists6. So there is some maximal k ≥ 0 with the desired property. b′ is clearly unique, namely b′ = |x|−2kb as functions on Rd. Furthermore, b′ is not divisible by r2 because k was chosen maximally.
B′ is also clearly a generating set: κ = ∑ i fibi = ∑ i(fir 2ki)b′i.
To prove linear independence, assume that∑ i fib ′ i = 0 . (151)
We can write k := maxi ki and then get 0 = r2k ∑ i fib ′ i = ∑ i fir 2(k−ki)bi . (152)
This means that fir
2(k−ki) = 0 ∀i, (153) which implies fi = 0 for all i, since r2(k−ki) is non-zero everywhere except in the origin and fi is continuous. Therefore, B′ is linearly independent and hence a basis of K.
Now we can formulate Proposition 19 more generally: for any basis B of K consisting only of homogeneous matrices of polynomials, the corresponding B′ is a basis of Kpol.

heading: H SOLUTIONS FOR IMPORTANT GROUPS
text: This and the next two appendices describe the solutions of the PDO equivariance constraint for subgroups of O(2) and O(3). In this appendix, we describe the general form of these solutions; we recommend readers who are not interested in all the details focus on this one. Appendix I contains proofs for some claims we make in this appendix where these proofs do not provide as much insight. Finally, Appendix J contains tables with concrete solutions, it is mainly relevant to the implementation of steerable PDOs.
H.1 SOLUTIONS FOR SUBGROUPS OF O(2)
To solve the steerability constraint for all (compact) subgroups G ≤ O(2) and for arbitrary representations ρin and ρout, Weiler & Cesa (2019) derive explicit bases only for irreducible representations, since the general case can then easily be computed (see Appendix B for details). This works exactly the same for PDOs as well.
The angular basis elements χβ that they describe for irreps are all matrices with entries of the form cos(kϕ) and sin(kϕ), where k differs between basis elements but is the same for all entries of one matrix. For example, for G = SO(2), ρin the frequency n irrep, i.e. ρin(g) = gn, and ρout trivial, the angular basis elements are
χ1 = (cos(nϕ), sin(nϕ)) and χ2 = (− sin(nϕ), cos(nϕ)) . (154) 6Here, the degree of a matrix of polynomials is the maximum of the degrees of all components (though in
this case, even if k is larger than the minimum degree, no b′ exists).
We will show how to find the steerable PDO basis using this example, but the method works exactly the same in all cases. Tables with all explicit solutions can be found in Appendix J.
As described in Section 3.2, we now need to multiply these matrices with the smallest power of |x| such that all entries become polynomials. We show in Appendix I that the necessary coefficient for cos(nϕ) and sin(nϕ) is |x|n, so in the example above, we get
χ̃1 = ( |x|n cos(nϕ), |x|n sin(nϕ) ) and χ̃2 = ( −|x|n sin(nϕ), |x|n cos(nϕ) ) .
(155) The entries are written in polar coordinates here, but they are in fact polynomials in the Cartesian coordinates x1 and x2. More precisely, we show in Appendix I that they are closely related to Chebyshev polynomials, based on which we derive the following explicit expressions:
|x|n cos(nϕ) = ∑
i≤n even
(−1) i2 ( n
i
) xn−i1 x i 2 and |x| n sin(nϕ) = ∑ i≤n odd (−1) i+1 2 ( n i ) xn−i1 x i 2 . (156)
This Cartesian form then allows us to interpret the polynomial as a differential operator, by applying the ring isomorphism D, i.e. plugging in ∂1 and ∂2 for x1 and x2.
If we set n = 1 in our example above, which corresponds to a vector field, we get simply |x| cos(ϕ) = x1 and |x| sin(ϕ) = x2, so we recover the angular PDO basis
D(χ̃1) = (∂1, ∂2) = div and D(χ̃2) = (−∂2, ∂1) = curl2D (157)
that we already derived in Section 2.3. To get a complete basis, we combine these PDOs with powers of the Laplacian.
H.2 SOLUTIONS FOR O(3) AND SO(3)
We now turn to the other cases for which steerable kernel solutions have been published, namely O(3) and SO(3). Like for O(2), we only need to consider pairs of irreducible representations. As described in (Weiler et al., 2018a; Lang & Weiler, 2021), we can build the corresponding angular parts χβ out of real spherical harmonics Ylm using Clebsch-Gordan coefficients. We show in Appendix I that we can apply this procedure to |x|lYlm instead of Ylm to obtain the corresponding χ̃β . We then build the basis by combining with powers of |x|2, as we described in Section 3.2 and already did for O(2). To find the corresponding differential operators D(χ̃β), we only need a Cartesian representation of the polynomials |x|lYlm, which is fortunately well-known (Varshalovich et al., 1988). This then again leads to a complete basis of the space of steerable PDOs, with the same general form containing powers of the Laplacian.
I PROOFS FOR SOLUTIONS FOR SUBGROUPS OF O(2) AND O(3)
In this section, we apply the method from Appendix G to subgroups of O(2) and O(3), by describing bases for the space of steerable kernels that satisfy the conditions of Proposition 19.
I.1 O(2)
As described in Appendix H.1, it is sufficient to consider irreducible representations ρin and ρout and for these representations, the angular part of the kernel has a basis consisting of matrices with entries of the form cos(nϕ) and sin(nϕ) (where n ∈ Z). Crucially, n is the same for all entries of one matrix (though it may differ between basis elements).
We claim that multiplying each angular basis element with r|n| gives a basis of the space of polynomial steerable kernels Kpol (as a R[r2]-module), as described in Appendix G. This follows from Proposition 19 if we can prove that r|n| cos(nϕ) and r|n| sin(nϕ) are
(i) polynomials (in x, y), (ii) homogeneous (of degree n),
(iii) and not divisible by r2.
To show that r|n| cos(nϕ) and r|n| sin(nϕ) are polynomials, we use Chebyshev polynomials, which can among other things be seen as a generalization of the addition theorems for sin(2ϕ) and cos(2ϕ). Specifically, they are families of polynomials Tn (first kind) and Un (second kind) defined for n ≥ 0 with the property that
cos(nϕ) = Tn(cosϕ) (158) sin(nϕ) = Un−1(cosϕ) sinϕ . (159)
We extend the definition to negative n by setting
T−n := Tn for n ≥ 0 (160) U−1 := 0 (161)
U−n−1 := −Un−1 for n ≥ 1 . (162)
Then Eq. (158) holds for all n ∈ Z, as follows immediately from the parity of sin and cos. Motivated by the close relation to Chebyshev polynomials, we then define the following notation for the matrix entries we are considering:
T̃n := r |n|Tn(cosϕ) = r |n| cos(nϕ) (163)
Ũn := r |n| sin(ϕ)Un−1(cosϕ) = r |n| sin(nϕ) . (164)
We now prove claim (i), that T̃n and Ũn are polynomials in the Cartesian coordinates x, y, by using a few well-known facts about Chebyshev polynomials: First, Tn has degree |n|, so the highest order term in T̃n is r|n|(cosϕ)|n| = x|n|. Similarly, Un−1 has degree |n| − 1 and the highest order term in Ũn is thus r|n|(cosϕ)|n|−1 sinϕ = x|n|−1y. Lower order terms have additional powers of r that aren’t “matched” by a cosine or sine. But the second fact about Chebyshev polynomials is that they are either even or odd, so all of these powers of r are even (i.e. powers of r2) and thus themselves polynomials.
That T̃n and Ũn are homogeneous of degree n – claim (ii) – is immediately clear from their definition.
It remains to show claim (iii), that they are not divisible by r2. For that, we use the following Lemma:
Lemma 23. Let p be a non-zero, harmonic, homogeneous polynomial. Then p is not divisible by r2.
Here, p is harmonic if its Laplacian vanishes.
Proof. For any homogeneous polynomial p, there is a unique decomposition
p = r2q + h , (165)
such that h is harmonic and homogeneous. Since q = 0, h = p is one such decomposition for p harmonic, there is no solution with h = 0 (unless p = 0). So non-zero homogeneous harmonic polynomials are not divisible by r2.
Now we only need to show that T̃n and Ũn are in fact harmonic. This can be done with a brief calculation in polar coordinates:
∆T̃n = ( ∂2r + 1
r ∂r +
1
r2 ∂2ϕ
) (rn cos(nϕ))
= n(n− 1)rn−2 cos(nϕ) + nrn−2 cos(nϕ)− n2rn−2 cos(nϕ) = 0 .
(166)
The same holds for Ũn:
∆Ũn = ( ∂2r + 1
r ∂r +
1
r2 ∂2ϕ
) (rn sin(nϕ))
= n(n− 1)rn−2 sin(nϕ) + nrn−2 sin(nϕ)− n2rn−2 sin(nϕ) = 0 .
(167)
In summary, we have shown that T̃n and Ũn are homogeneous polynomials not divisible by r2, which makes Proposition 19 applicable.
In order to make the basis practically applicable, we also give explicit Cartesian expressions for T̃n and Ũn. We restrict ourselves to non-negative n to keep the notation less cluttered. Equation (160) immediately gives the case for n < 0. An explicit formula for the Chebyshev polynomials of the first kind is
Tn(x) = bn2 c∑ k=0 ( n 2k )( x2 − 1 )k xn−2k . (168)
It follows that
T̃n = r nTn(cosϕ)
= √ x2 + y2 n Tn ( x√
x2 + y2
)
= bn2 c∑ k=0 ( n 2k )( −y2 )k xn−2k
= bn2 c∑ k=0 ( n 2k ) (−1)ky2kxn−2k .
(169)
Similarly, the Chebyshev polynomials of the second kind are given by
Un(x) = bn2 c∑ k=0 ( n+ 1 2k + 1 )( x2 − 1 )k xn−2k , (170)
which means that
Ũn = r n sinϕUn−1(cosϕ)
= y √ x2 + y2 n−1 Un−1 ( x√
x2 + y2
)
= y bn−12 c∑ k=0 ( n 2k + 1 )( −y2 )k xn−1−2k
= bn−12 c∑ k=0 ( n 2k + 1 ) (−1)ky2k+1xn−(2k+1) .
(171)
I.2 O(3)
As already mentioned in Appendix H.2, for O(3) and SO(3), the irreps angular basis between irreducible representation can be built by combining real spherical harmonics Ylm using ClebschGordan coefficients. We refer to Weiler et al. (2018a); Lang & Weiler (2021) for details on how this works since it is exactly the same procedure whether one uses kernels or PDOs. The only fact we need to know here is that each basis element is built from spherical harmonics with the same l.
The necessary steps are now very similar to those for O(2) and its subgroups. We will use rlYlm where we had T̃n and Ũn for O(2). Here it becomes important that inside one basis element, only one l appears, because this means we can multiply the entire matrix of polynomials by rl.
We would then need to show that these are homogeneous polynomials not divisible by r2 but for spherical harmonics, it is already very well known that rlYlm are homogeneous harmonic polynomials of degree l; Lemma 23 then implies that they are not divisible by r2. For explicit Cartesian expressions, we refer to e.g. Varshalovich et al. (1988).

heading: J SOLUTION TABLES
text: Using the results from Appendices G and I, we now very easily get complete bases for all subgroups of O(2) and for all irreducible representations by transferring the solutions by Weiler & Cesa (2019). The appendix of (Weiler & Cesa, 2019) contains tables with all kernel solutions; we simply replace every term of the form cos(kϕ) with T̃k and sin(kϕ) with Ũk to get the following tables.
For the sake of readability, we only write the polynomials T̃k and Ũk inside the tables, though the PDOs themselves should of course be D(T̃k) and D(Ũk). Explicit formulas for T̃k and Ũk were given in Eq. (156).
All bases described here are module bases for the space of steerable PDOs as an R[∆]-module. See Section 3.2 and Appendix G for details.
J.1 SPECIAL ORTHOGONAL GROUP SO(2)
The irreducible representations of SO(2) are the trivial representation ψ0 and those of the form
ψk : SO(2)→ GL(R2), g 7→ gk . (172)
The bases for all the combinations of these irreps are:
Out In ψ0 ψn, n ∈ N>0
ψ0 (1) ( T̃n Ũn ) , ( −Ũn T̃n )
ψm, m ∈ N>0
( T̃m Ũm ) , ( −Ũm T̃m
) ( T̃m−n −Ũm−n Ũm−n T̃m−n ) , ( −Ũm−n −T̃m−n T̃m−n −Ũm−n ) ,(
T̃m+n Ũm+n Ũm+n −T̃m+n
) , ( −Ũm+n T̃m+n T̃m+n Ũm+n )
J.2 ORTHOGONAL GROUP O(2)
Elements of O(2) can be written as a tuple (r, s) consisting of a rotation r ∈ SO(2) and a flip s ∈ {±1}. The irreducible representations with frequency k > 0 are similar to those for SO(2), only with the additional flip:
ψ1,k : SO(2)→ GL(R2), (r, s) 7→ rk ◦ s . (173)
Here, s is understood as a map acting on R2 either as the identity or by flipping along a certain axis. In contrast to SO(2), there are now two irreducible representations for k = 0, namely the trivial representation ψ0,0(r, s) = 1 and the representation ψ1,0(r, s) = s. The solutions for all possible combinations of these irreducible representations are as follows:
Out In ψ0,0 ψ1,0 ψ1,n, n ∈ N>0
ψ0,0 (1) ∅ ( −Ũn T̃n ) ψ1,0 ∅ (1) ( T̃n Ũn
) ψ1,m,
m ∈ N>0
( −Ũm T̃m ) ( T̃m Ũm ) ( T̃m−n −Ũm−n Ũm−n T̃m−n ) , ( T̃m+n Ũm+n Ũm+n −T̃m+n )

heading: J.3 REFLECTION GROUP ±1
text: The reflection group has only two irreducible representations: the trivial representation ψ0 and the representation ψ1(s) = s. Both are one-dimensional, so all the PDOs are only 1× 1 matrices:
Out In ψ0 ψ1
ψ0 ( T̃µ ) ( Ũµ ) ψ1 ( Ũµ ) ( T̃µ )
To get the full basis, µ needs to take on all natural numbers (in practice, we use all µ up to some maximum value). Note that we assume that the reflection is with respect to the x1-axis, both here and later for DN . To get other reflection axes, the PDOs simply need to be rotated.

heading: J.4 CYCLIC GROUP CN
text: The irreducible representations of CN are the same as those of SO(2) but only up to a frequency of k = bN−12 c. If N is even, there is an additional one-dimensional irreducible representation, namely
ψN/2(θ) = cos
( N
2 θ
) . (174)
We then get the following solutions, where t ranges over Z and t̂ over N:
Out In ψ0 ψN/2 (if N even) ψn, 1 ≤ n < N/2
ψ0
( T̃t̂N ) ,(
Ũt̂N )
( T̃(t̂+1/2)N ) ,(
Ũ(t̂+1/2)N
) (−Ũn+tN T̃n+tN) ,( T̃n+tN Ũn+tN )
ψ1 (N even)
( T̃(t̂+1/2)N ) ,(
Ũ(t̂+1/2)N
) (T̃t̂N) ,( Ũt̂N ) (−Ũn+(t+1/2)N T̃n+(t+1/2)N) ,( T̃n+(t+1/2)N Ũn+(t+1/2)N )
ψn 1 ≤ n < N/2
( −Ũm+tN T̃m+tN ) ,(
T̃m+tN Ũm+tN
) ( −Ũm+(t+1/2)N T̃m+(t+1/2)N ) ,(
T̃m+(t+1/2)N Ũm+(t+1/2)N
) ( T̃m−n+tN −Ũm−n+tN Ũm−n+tN T̃m−n+tN ) , ( −Ũm−n+tN −T̃m−n+tN T̃m−n+tN −Ũm−n+tN ) ,(
T̃m+n+tN Ũm+n+tN Ũm+n+tN −T̃m+n+tN
) , ( −Ũm+n+tN T̃m+n+tN T̃m+n+tN Ũm+n+tN )

heading: J.5 DIHEDRAL GROUP DN
text: Similarly toCN , the irreducible representations ofDN are the same as those of O(2) up to a frequency of k = bN−12 c. If N is even, there are two additional one-dimensional irreducible representations, namely
ψ0,N/2(θ, s) = cos
( N
2 θ
) ,
ψ1,N/2(θ, s) = s cos
( N
2 θ
) .
(175)
The solutions are (again with t ∈ Z and t̂ ∈ N):
Out In ψ0,0 ψ1,0 ψ0,N/2 (N even ) ψ1,N/2 (N even) ψ1,n, 1 ≤ n < N/2
ψ0,0 ( T̃t̂N ) ( Ũt̂N ) ( T̃(t̂+1/2)N ) ( Ũ(t̂+1/2)N ) ( −Ũn+tN T̃n+tN ) ψ1,0 ( Ũt̂N ) ( T̃t̂N ) ( Ũ(t̂+1/2)N ) ( T̃(t̂+1/2)N ) ( T̃n+tN Ũn+tN
) ψ0,N/2,
(N even)
( T̃(t̂+1/2)N ) ( Ũ(t̂+1/2)N ) ( T̃t̂N ) ( Ũt̂N ) ( −Ũn+(t+1/2)N T̃n+(t+1/2)N ) ψ1,N/2,
(N even)
( Ũ(t̂+1/2)N ) ( T̃(t̂+1/2)N ) ( Ũt̂N ) ( T̃t̂N ) ( T̃n+(t+1/2)N Ũn+(t+1/2)N )
ψ1,m, 1 ≤ m ≤ N/2 ( −Ũm+tN T̃m+tN ) ( T̃m+tN Ũm+tN ) ( −Ũm+(t+1/2)N T̃m+(t+1/2)N ) ( T̃m+(t+1/2)N Ũm+(t+1/2)N
) (T̃m−n+tN −Ũm−n+tN Ũm−n+tN T̃m−n+tN ) ,(
T̃m+n+tN Ũm+n+tN Ũm+n+tN −T̃m+n+tN
)

heading: K DISCRETIZATION METHODS FOR PDOS

heading: K.1 FINITE DIFFERENCES
text: Finite difference methods are common in machine learning; for example, the discretization of ddx as [−1 1] or [−1 0 1], or of d 2
dx2 as [1 −2 1] all use finite difference methods. To understand where these filters come from, we need the following well-known result: Proposition 24. Let x1, . . . , xN ∈ R be arbitrary but distinct grid points. Then for m ≤ N − 1, there are unique coefficients w(m)n such that the approximation
f (m)(0) ≈ 1 hm N∑ n=1 w(m)n f(hxn) (176)
has an error O(hN−m) for any f ∈ CN (R).
The coefficients w(m)n are called finite difference coefficients and approximating derivatives using Eq. (176) is the finite difference method. We will soon describe how to generalize this to higher dimensions as well.
We remark that O(hN−m) is an asymptotic upper bound on the error, and it can sometimes be lower, even for all f . For example, the central difference discretization of d 2
dx2 as [1 −2 1] uses N = 3 grid points but still achieves an error of O(h2), rather than O(h). For details on when such a “boosted” order of accuracy occurs, see (Sadiq & Viswanath, 2014).
Note that Eq. (176) can be generalized to
f (m)(x) ≈ 1 hm N∑ n=1 w(m)n f(hxn + x) . (177)
This follow immediately because the coefficients w(m)n don’t depend on the function f , so we can apply Proposition 24 to τ−xf .
Particularly interesting for us is the case of a regular grid. We can use infinitely many grid points xn ∈ Z as long as we demand that w(m)n is zero for almost all n. Then we get
f (m)(x) ≈ 1 hm ∑ n w(m)n f(hn+ x) . (178)
Fixing h = 1, this is exactly the cross correlation
f (m) ≈ w(m) ? f (179)
if we interpret w(m) as a function n 7→ w(m)n . This is why, in the end, we discretize a derivative by convolving with some stencil, such as [1 −2 1], at least on a regular 1D grid.
Generalizing Proposition 24 to higher dimensions does not work in a straightforward way. However, if we restrict ourselves to regular grids, then finite difference methods can be easily applied to PDOs. The idea is very simple: a PDO such as ∂x∂2y can be interpreted as first applying ∂ 2 y and then ∂x (or the other way around). So we discretize each of these with the one-dimensional finite difference method described before, and then we convolve with both filters one after the other.7 We can also combine the two one-dimensional filter into one two-dimensional filter, the outer product of the two.8 The asymptotic error of this discretization will simply be the highest asymptotic error along all the dimensions, so we get similar guarantees.

heading: K.2 RBF-FD
text: As mentioned, Proposition 24 does not directly generalize to higher dimensions. So to discretize a PDO on arbitrary point clouds in higher dimensions, a somewhat different approach is needed.
RBF-FD is one such method and works as follows: we still want to approximate a derivative using
∂αf(0) ≈ N∑ n=1 wαnf(xn) , (180)
similar to finite difference methods. Here, xn ∈ Rd are arbitrary (but again distinct) points. The idea is now that we require this approximation to be exact if f(x) = ϕ(‖x− xn‖), where ϕ is an arbitrary but fixed radial basis function. In words, the approximation should become exact for a certain radial basis function centered on any of the points xn. This leads to a linear system, which is solved for the coefficients wαn .
9 For more details on both finite differences and RBF-FD, see for example Fornberg & Flyer (2015).

heading: K.3 GAUSSIAN DERIVATIVES
text: Discretizing PDOs using derivatives of Gaussians is very simple to describe: given grid points xn ∈ Rd, we approximate using
∂αf(0) ≈ N∑ n=1 (∂αG(xn;σ)) f(xn) (181)
whereG(x;σ) is a Gaussian kernel with standard deviation σ centered around 0. σ is a free parameter; larger σ will lead to a stronger denoising effect.
On regular grids, this again turns into a cross-correlation, with the filter being the derivative ∂αG(xn;σ) evaluated on the grid coordinates.
In Appendix E we briefly touch on a possible interpretation of this discretization method using the distributional framework for PDOs.

heading: L RELATION TO PDO-ECONVS
text: In this section, we describe how PDO-eConvs (Shen et al., 2020) fit into the framework of steerable PDOs. We mostly follow the original notation from Shen et al. (2020) when describing PDO-eConvs but do make some minor changes to avoid clashes and confusion with our own notation.
As in our presentation, Shen et al. (2020) use polynomials to describe PDOs. One difference is that they never explicitly use matrices of polynomials, because they model the feature space somewhat differently (which we will discuss in a moment). They write H for the polynomial describing a PDO
7As we have seen, finite difference methods can most immediately be seen as performing a cross-correlation rather than a convolution. However, we can easily switch to convolutions by flipping the filter.
8For a simple PDO such as ∂x∂2y , this may be undesirable for computational reasons. But in practice, we have PDOs that are sums of such pure terms and thus don’t factorize.
9In practice, one often solves an extended linear system containing additional low-order polynomials, but we won’t discuss that here.
(where we would write e.g. p) and write χ(A) for the corresponding PDO transformed by A ∈ O(d). In our notation,
χ(A) := D(A ·H) = D(H ◦A−1) . (182)
PDO-eConvs use two types of PDO layers. The first one, Ψ, can be interpreted as a steerable PDO with ρin trivial and ρout regular. It maps the scalar input to the network to the internally used regular representation. The second layer type, Φ maps between regular representations and is used for hidden layers. At the end, pooling is performed to obtain a scalar output again.
The first layer type is defined as
Ψ : C∞(Rd)→ C∞(Ẽ(d)), Ψ(f)(x,A) := (χ(A)f)(x) . (183)
Here, Ẽ(d) := Rd o S with S ≤ O(d); in practice, S needs to be a finite subgroup, i.e. CN or DN . Elements of Ẽ(d) can be uniquely written as (x,A) with x ∈ Rd and A ∈ S.
There is an obvious bijection C∞(Ẽ(d)) ∼= C∞(Rd,Rc), where c := |S| is the order of S, i.e. the number of group elements. Concretely, we define
Θ : C∞ ( Ẽ(d) ) → C∞(Rd,Rc), f 7→ ( x 7→ ( f(x,A1), . . . , f(x,Ac) )) (184)
where A1, . . . , Ac is an enumeration of the group elements of S. We will therefore interpret C∞(Ẽ(d)) as a c-dimensional feature field over Rd, and we will show that using regular representations for this field (and trivial ones for the input) makes the PDO-eConv layers equivariant and thus steerable PDOs.
First, note that under the Θ bijection, the first PDO-eConv layer type Ψ becomes a c× 1 matrix of PDOs, namely
D(HΨ) := D  A1 ·H... Ac ·H   = χ (A1) ... χ(Ac)  . (185) What we mean by this is that the diagram
C∞(Rd) C∞(Ẽ(d))
C∞(Rd) C∞(Rd,Rc)
Ψ
Θ
D(HΨ)
commutes. Concretely, we have
Ψ(f)(x,Ai) = (χ (Ai)f)(x) = (D(HΨ)f)(x)i . (186)
So we need to check whether HΨ satisfies the PDO steerability constraint for trivial to regular PDOs:
HΨ(Ax) = ρregular(A)HΨ(x) . (187)
Using the definition of ρregular(A), we can rewrite the RHS as
ρregular(A)HΨ(x) = c∑ k=1 ρregular(A)Hψ(x)keAk
= c∑ k=1 (Ak ·H)(x)eAAk
= c∑ k=1 H(A−1k x)eAAk
= c∑ l=1 H(A−1l Ax)eAl .
(188)
Here, we use basis vectors eA1 , . . . , eAc for Rc, with the same enumeration A1, . . . , Ac of S used to define Θ. For the final step, we reparameterized the sum with Al := AAk.
The LHS of Eq. (187) can be written as
HΨ(Ax) = (A1 ·H(Ax), . . . , Ac ·H(Ax)) T = ( H(A−11 Ax), . . . , H(A −1 c Ax) )T .
(189)
This is the same as the final term in Eq. (188), which proves that the PDO steerability constraint is satisfied.
The second PDO-eConv layer type, mapping between regular representations, is defined as
Φ : C∞(Ẽ(d))→ C∞(Ẽ(d)), Φ(e)(x,A) := c∑ j=1 ( χ (A) Aj e ) (x,AAj) . (190)
χ (A) Aj are c different PDOs and they act on e ∈ C∞(Ẽ(d)) by acting on each of the c components separately. Under the Θ bijection, Φ becomes
Θ(Φ(e))(x)i = c∑ j=1 ( χ (Ai) Aj e ) (x,AiAj)
= c∑ j=1 ( χ (Ai) Aj Θ(e)AiAj ) (x)
= c∑ j=1 ( χ (Ai) A−1i Aj Θ(e)Aj ) (x) .
(191)
We can therefore represent it as a c× c PDO D(HΦ) with (HΦ)ij (x) = HA−1i Aj (A−1i x) =: Hij(x) , (192)
where HA−1i Aj is the polynomial that induces χ (I) A−1i Aj . This makes the diagram
C∞(Ẽ(d)) C∞(Ẽ(d))
C∞(Rd,Rc) C∞(Rd,Rc)
Φ
Θ Θ
D(HΦ)
commute, similar to the case discussed above. So again, we need to check that HΦ satisfies the PDO steerability constraint, this time for ρin and ρout both regular:
HΦ(Ax) = ρregular(A)HΦ(x)ρregular(A −1) . (193)
Writing out HΨ in its components, this becomes HΨ(Ax) ! = ∑ i,j ρregular(A)Hij(x)eAie T Ajρregular(A) −1
(1) = ∑ i,j Hij(x)eAAie T AAj
= ∑ i,j HA−1i Aj (A−1i x)eAAie T AAj
(2) = ∑ i,j HA−1i AA−1Aj (A−1i Ax)eAie T Aj
= ∑ i,j HA−1i Aj (A−1i Ax)eAie T Aj
= ∑ i,j Hij(Ax)eAie T Aj .
(194)
The first and the last term are the same, just written out in components on the RHS, so the steerability constraint is again satisfied. For (1), we used that ρregular(A) is orthogonal, and thus
eTj ρregular(A) −1 = eTj ρregular(A) T = (ρregular(A)ej) T . (195)
(2) was again a reparameterization of the sum, with AAi 7→ Ai and AAj 7→ Aj . The other steps are only simplifications and plugging in definitions.
In conclusion, we have shown that there is a simple bijection between the feature spaces used for PDO-eConv hidden layers and the feature fields we use, and that under this bijection, PDO-eConvs correspond to steerable PDOs with regular representations (and trivial representations for the input).
It is relatively easy to adapt the argument we present for the converse direction: every equivariant PDO between two regular feature fields (or from a scalar to a regular one) can be interpreted as a PDO-eConv layer.

heading: M ADDITIONAL EXPERIMENTAL RESULTS

heading: M.1 RESTRICTION EXPERIMENTS
text: Table 3 shows additional results on MNIST-rot. The general architecture and hyperparameters are the same as in the experiments in Section 4 with regular representations, using our basis. However, in the experiments in this section, the first five layers areD16-equivariant, while the final PDO/convolutional layer is C16-equivariant. The motivation for this is that while the input images do not have global reflectional symmetry, such symmetry occurs on smaller scales, so that stronger equivariance in earlier layers might be helpful.
However, we don’t observe clear improvements over pure C16-equivariance. A reason could be that even the C16-equivariant networks are already very parameter efficient compared to classical CNNs, so that parameter efficiency and equivariance are not a bottleneck anymore. It is also possible that the architecture would need to be adapted slightly to profit from the D16-equivariant layers.

heading: M.2 STENCIL IMAGES
text: Figure 4 contains examples of stencils used during the MNIST-rot experiments. For the 3× 3 stencil, the different methods yield qualitatively similar results (though FD and RBF-FD have fewer stencils that make use of the four corners). But for 5 × 5 stencils, kernels and Gauss discretization make significantly more use of the outer stencil points than FD and RBF-FD.

heading: N DETAILS ON EXPERIMENTS

heading: N.1 MNIST-ROT EXPERIMENTS
text: For the MNIST-rot experiments, we use an architecture similar to one from (Weiler & Cesa, 2019). Table 4 contains a listing of all the layers. Each conv block consists of a 3× 3 or 5× 5 steerable layer, either convolutional or a PDO, followed by batch-normalization and an ELU nonlinearity. The output fields are the number of C16-regular feature fields that is used; in the case of the Vanilla CNN and for quotient representations, the number of fields is adjusted so that the parameter count is approximately preserved.
For the quotient experiments, we use 5ρregular⊕ 2ρC16/C2quot ⊕ 2ρ C16/C4 quot ⊕ 4ρtrivial as the representation, where the numbers are scaled to reach the same parameter count as the model with only regular representations. This combination of representations is the same one used by Weiler & Cesa (2019) and we refer to their appendix for motivation on why we need to combine different representations when using quotients.
We trained all MNIST-rot models for 30 epochs with Adam (Kingma & Ba, 2015) and a batch size of 64. The training data was normalized and augmented using random rotations. For the final training runs, we used the entire set of 12k training plus validation images, as is common practice on MNIST-rot. The initial learning rate was 0.05, which was decayed exponentially after a burn-in of 5 epochs at a rate of 0.7 per epoch. We used a dropout of 0.5 after the fully connected layer, and a weight decay of 1e-7.
These hyperparameters are based on those used in (Weiler & Cesa, 2019); the main difference is that we use another learning rate schedule, which works better than the original one for all models.
For the Gaussian discretization models, we use a standard deviation of σ = 1 for 3× 3 stencils and σ = 1.3 for 5× 5 stencils; we chose these values by visual inspection of the stencils, with the aim that full use is made of the stencil (see Appendix M). The RBF-FD discretization uses third-order polyharmonic basis functions, i.e. ϕ(r) = r3.

heading: N.2 STL-10 EXPERIMENTS
text: For the STL-10 experiments, we used exactly the same architecture and hyperparameters as Weiler & Cesa (2019), which in turn are essentially those of DeVries & Taylor (2017). This means we train a Wide-ResNet-16-8 (Zagoruyko & Komodakis, 2016) for 1000 epochs, with SGD and Nesterov momentum of 0.9, a batch size of 128 and weight decay of 5e-4. We begin with a learning rate of 0.1 and divide it by 5 after 300, 400, 600 and 800 epochs. For data augmentation, we pad the image by 12 pixels, then randomly crop a 96× 96 pixel patch, randomly flip horizontally, and apply Cutout (DeVries & Taylor, 2017) with a cutout size of 60× 60 pixels.

heading: N.3 COMPUTATIONAL REQUIREMENTS
text: We performed our experiments on an internal cluster with a GeForce RTX 2080 Ti and 6 CPU cores. A single run of an MNIST-rot model took about 12 minutes and a run of the STL-10 model about 5.5h. Multiplying this by the number of experiments we did and by six runs with different seeds, the MNIST-rot results took about 34 hours and the STL-10 results about 264 hours to produce. The initial tests we did to debug and find a good learning rate schedule took much less time than that. So we estimate that producing this paper took a bit over 300 GPU-hours on the GeForce RTX 2080 Ti.

Review summaries:
The paper develops the theory of equivariant partial differential operators -- a steerable PDO which is equivariant under any given symmetry. Interestingly, this work reveals the relation between convolution and PDO by unifying them into a single framework. The work also provides the rigorous theoretical analysis and proof. The paper also provides experimental results to validate its theoretical analysis. 

This paper introduces equivariant Partial Differential Operators as a drop in replacement for equivariant kernels in $G$-steerable networks. The paper provides a thorough treatment on the subject matter and gives the necessary equivariant constraints for PDO's in a similar flavor to $G$-steerable kernels. In addition, explicit solutions to these constraints are provided for subgroups of $O(2)$, $O(3)$ and $SO(3)$. Finally, the paper also empirically validates equivariant PDO's on rotation MNIST and STL-10 highlighting important discussion points on the discretization error.

The authors developed the general theory of equivariant partial differential operators (PDOs) between feature fields on the Euclidean space. Given an arbitrary group G, the G-steerability constraint is derived that fully characterized when a PDO is equivariant between vector fields for given representations. Experiments on the rotated-MNIST and STL-10 datasets have been conducted to compare the performance of equivariant-PDOs and equivariant steerable kernels.