# The first Weyl algebra

Consider the vector space $${\mathcal A}$$ of linear operators on univariate functions; $${\mathcal A}$$ can be made into an algebra where multiplication (denoted by juxtaposition) is defined as operator composition:

$({\mathcal O}_1{\mathcal O}_2)f={\mathcal O}_1({\mathcal O}_2f)$

We consider the algebra generated by by $$\partial\colon f\longrightarrow f'$$ [that is, $$(\partial f)(x) = f'(x)$$] and $$x\colon f\longrightarrow xf$$ [that is, $$(xf)(x) = xf(x)$$]. This is known as the (first) Weyl algebra. We observe that the Weyl algebra is not commutative: $$\partial xf=(xf)'=f+xf'$$ but $$x\partial f=xf'$$, so $$\partial x=x\partial+1$$. The algebra generated by $$\left\lbrace x,\partial\right\rbrace$$ will include elements such as $$7\partial + 4\partial x\partial^3 x$$, which would map $$f$$ to $$7f' + 4\left(x\left(xf\right)'''\right)'$$.

It can be shown that any element of the Weyl algebra can be expressed in the standard form

$\sum_i a_i \partial^{p_i}x^{q_i}$

for real $$a_i$$ and nonnegative integers $$p_i,q_i$$. Converting a general word to standard form is not straightforward but we have

$\partial x^n = x^n\partial + nx^{n-1}$

and

$\partial^n x = x\partial^n + n\partial^{n-1}.$

We can apply these rules recursively to find standard form for products $$(\partial^i x^j)(\partial^l x^m)$$. Alternatively we may follow Wolf and use the fact that

$(\partial^i x^j)(\partial^lx^m)= \sum_{r=0}^j{j\choose r}{l\choose r}\partial^{i+l-r}x^{j+m-r}.$

All this is done in the package. Consider $$d_1=\partial x + 2\partial^3$$ and $$d_2=3+7\partial -5x^2\partial^2$$. Observing that $$d_1$$ and $$d_2$$ are in standard form, package idiom to create these operators would be:

library(weyl)
d1 <- weyl(spray(rbind(c(1,1),c(0,3)),c(1,2)))
d2 <- weyl(spray(rbind(c(0,0),c(1,0),c(2,2)),c(3,7,-5)))
d1
## A member of the Weyl algebra:
##  x d     val
##  0 3  =    2
##  1 1  =    1
d2
## A member of the Weyl algebra:
##  x d     val
##  2 2  =   -5
##  1 0  =    7
##  0 0  =    3

Observe that, like the spray package, the order of the matrix rows is not defined. We may apply the usual rules of arithmetic to these objects:

d1*d2
## A member of the Weyl algebra:
##  x d     val
##  2 5  =  -10
##  1 4  =  -60
##  2 2  =  -10
##  0 3  =  -54
##  1 3  =   14
##  1 0  =    7
##  0 2  =   42
##  3 3  =   -5
##  2 1  =    7
##  1 1  =    3

We may choose to display results in symbolic form:

options(polyform=TRUE)
(d1^2 + d2) * (d2 - 3*d1)
## A member of the Weyl algebra:
## 9 -63*x^3*d^2 +168*d^5 +20*x^4*d^4 +87*x^3*d^3 -70*x^2*d -24*d^9
## +126*d^2 +49*x^2 -9*x*d -20*x^2*d^8 +112*x*d^3 -234*d^3 +49*x
## -178*x^2*d^5 -480*x*d^4 +24*x^2*d^2 -696*d^6 +28*x*d^6 -276*x*d^7
## -20*x^3*d^6 +28*x^2*d^4

### Comparison with mathematica

Mathematica can deal with operators and we may compare the two systems’ results for $$\partial^2x\partial x^2$$:

In[1] := D[D[x*D[x^2*f[x],x],x],x] // Expand

Out[1] := 4 f[x] + 14 x f'[x] + 8 x^2 f''[x] + x^3f'''[x]

x <- weyl(cbind(0,1))
D <- weyl(cbind(1,0))
x^2*D*x*D^2
## A member of the Weyl algebra:
## 4 +x^3*d^3 +8*x^2*d^2 +14*x*d

# Further Weyl algebras

The package supports multivariate Weyl algebras. Consider:

options(polyform=FALSE)  # revert to default print method
set.seed(0)
x <- rweyl()
x
## A member of the Weyl algebra:
##   x  y  z dx dy dz     val
##   2  0  1  2  0  1  =    3
##   0  1  2  2  0  1  =    2
##   1  0  2  1  0  1  =    1

Above, object x is a member of the operator algebra generated by $$\left\lbrace\partial_x,\partial_y,\partial_z,x,y,z\right\rbrace$$. We may verify associativity of multiplication:

x <- rweyl(n=1,d=2)
y <- rweyl(n=2,d=2)
z <- rweyl(n=3,d=2)
options(polyform=TRUE)
x*(y*z)
## A member of the Weyl algebra:
## +6*x*y^2*dx^2*dy +36*x*y^3*dx*dy^3 +2*x^2*y^4*dx*dy^4
## +4*x^2*y^2*dx*dy^5 +36*x*y^2*dx*dy^2 +x*y^4*dx*dy^3 +4*x*y^3*dx*dy^2
## +x^2*y^4*dx^2*dy^3 +12*x^2*y^2*dx*dy^2 +2*x^2*y^2*dx^2*dy^4
## +2*x*y^2*dx*dy^4 +2*x*y^2*dx*dy +2*x^2*y^2*dx^2*dy +4*x^2*y^3*dx^2*dy^2
## +3*x*y^4*dx^2*dy^3 +12*x^2*y^3*dx*dy^3 +6*x*y^4*dx*dy^4
## +12*x*y^3*dx^2*dy^2
(x*y)*z
## A member of the Weyl algebra:
## +2*x*y^2*dx*dy^4 +2*x^2*y^2*dx^2*dy^4 +3*x*y^4*dx^2*dy^3
## +x^2*y^4*dx^2*dy^3 +12*x^2*y^2*dx*dy^2 +x*y^4*dx*dy^3
## +2*x^2*y^4*dx*dy^4 +4*x^2*y^3*dx^2*dy^2 +2*x^2*y^2*dx^2*dy
## +2*x*y^2*dx*dy +4*x*y^3*dx*dy^2 +12*x^2*y^3*dx*dy^3 +6*x*y^4*dx*dy^4
## +12*x*y^3*dx^2*dy^2 +6*x*y^2*dx^2*dy +36*x*y^3*dx*dy^3
## +36*x*y^2*dx*dy^2 +4*x^2*y^2*dx*dy^5

Comparing the two results above, we see that they apparently differ. But the apparent difference is due to the fact that the terms appear in a different order, a feature that is not algebraically meaningful. We may verify that the expressions are indeed algebraically identical:

x*(y*z) - (x*y)*z
## A member of the Weyl algebra:
## the NULL multinomial of arity 4

# Derivations

A derivation $$D$$ of an algebra $${\mathcal A}$$ is a linear operator that satisfies $$D(d_1d_2)=d_1D(d_2) + D(d_1)d_2$$, for every $$d_1,d_2\in{\mathcal A}$$. If a derivation is of the form $$D(d) = [d,f] = df-fd$$ for some fixed $$f\in{\mathcal A}$$, we say that $$D$$ is an inner derivation:

$D(d_1d_2) = d_1d_2f-fd_1d_2 = d_1d_2f-d_1fd_2 + d_1fd_2-fd_1d_2 = d_1(d_2f-fd_2) + (d_1f-fd_1)d_2 = d_1D(d_2) + D(d_1)d_2$

Dirac showed that all derivations are inner derivations for some $$f\in{\mathcal A}$$. The package supports derivations:

f <- rweyl()
D <- as.der(f)  # D(x) = xf-fx

Then

d1 <- rweyl()
d2 <- rweyl()
D(d1*d2) == d1*D(d2) + D(d1)*d2
## [1] TRUE

# Low-level considerations and generalizations

In the package, the product is customisable. In general, product a*b [where a and b are weyl objects] is dispatched to the following sequence of functions:

• weyl_prod_multivariate_nrow_allcolumns()
• weyl_prod_multivariate_onerow_allcolumns()
• weyl_prod_multivariate_onerow_singlecolumn()
• weyl_prod_univariate_onerow()
• weyl_prod_helper3() [default]

In the above, “univariate” means “generated by $$\left\lbrace x,\partial_x\right\rbrace$$” [so the corresponding spray object has two columns]; and “multivariate” means that the algebra is generated by more than one variable, typically something like $$\left\lbrace x,y,z,\partial_x,\partial_y,\partial_z\right\rbrace$$.

The penultimate function weyl_prod_univariate_onerow() is sensitive to option prodfunc which specifies the recurrence relation used. This defaults to weyl_prod_helper3():

weyl_prod_helper3
## function (a, b, c, d)
## {
##     f <- function(r) {
##         factorial(r) * choose(b, r) * choose(c, r)
##     }
##     ind <- numeric(0)
##     val <- numeric(0)
##     for (r in 0:b) {
##         ind <- rbind(ind, c(a + c - r, b + d - r))
##         val <- c(val, f(r))
##     }
##     spray(ind, val, addrepeats = TRUE)
## }
## <bytecode: 0x561fc10f09a8>
## <environment: namespace:weyl>

Function weyl_prod_helper3() follows Wolf. This gives the univariate concatenation product $$(\partial^a x^b)(\partial^c x^d)$$ in terms of standard generators:

$\partial^a x^b \partial^c x^d=\sum_{r=0}^b r!{b\choose r}{c\choose r} \partial^{a+c-r}x^{b+d-r}$

The package also includes lower-level function weyl_prod_helper1() implementing $$\partial^a x^b \partial^c x^d=\partial^ax^{b-1}\partial^cx^{d+1} + c\partial^ax^{b-1}\partial^{c-1}x^d$$ (together with suitable bottoming-out). I expected function weyl_prod_helper3() to be much faster than weyl_prod_helper1() but there doesn’t seem to be much difference between the two.

# Generalized commutator relations

We can exploit this package customisability by considering, instead of $$\left\lbrace x,\partial\right\rbrace$$, the algebra generated by $$\left\lbrace e,\partial\right\rbrace$$, where $$e$$ maps $$f$$ to $$e^xf$$: if $$f$$ maps $$x$$ to $$f(x)$$, then $$ef$$ maps $$x$$ to $$e^xf(x)$$. We see that $$\partial e-e\partial=e$$. With this, we can prove that $$\partial^ne=e(1+\partial)^n$$ and $$e^n\partial=e^n\partial+ne^n$$ and, thus

$(e^a\partial^b)(e^c\partial^d) =e^{a+1}(1+\partial)^be^{c-1}\partial^d =e^{a}\partial^{b-1}e^{c}\partial^{d+1}+ce^{a}\partial^{b-1}e^{c}\partial^d$

We may implement this set in package idiom as follows:

weyl_e_prod <- function(a,b,c,d){
if(c==0){return(spray(cbind(a,b+d)))}
if(b==0){return(spray(cbind(a+c,d)))}
return(
Recall(a,b-1,c,d+1) +
c*Recall(a,b-1,c,d)  # cf: c*Recall(a,b-1,c-1,d)) for regular Weyl algebra
)
}

Then, for example, to calculate $$\partial^2e=e(1+2\partial+\partial^2)$$:

options(prodfunc = weyl_e_prod)
options(weylvars = "e")  # changes print method
d <- weyl(spray(cbind(0,1)))
e <- weyl(spray(cbind(1,0)))
d*d*e
## A member of the Weyl algebra:
## +e +2*e*d +e*d^2
d^2*e
## A member of the Weyl algebra:
## +e +2*e*d +e*d^2

By way of verification:

d^5*e == e*(1+d)^5
## [1] TRUE

which verifies that indeed $$\partial^5e=e(1+\partial)^5$$. Another verification would be to cross-check with Mathematica, here working with $$\partial e\partial^2e$$:

In[1] := D[Exp[x]*D[D[Exp[x]*f[x],x],x],x]

Out[1] := 2E^2x f[x] + 5E^2x f'[x] + 4E^2xf''[x] + E^2x f'''[x]

options(polyform = TRUE)
d*e*d^2*e
## A member of the Weyl algebra:
## +2*e^2 +5*e^2*d +4*e^2*d^2 +e^2*d^3

We can manipulate more complicated expressions too. Suppose we want to evaluate $$(1+e^2\partial)(1-5e^3\partial^3)$$:

o1 <- weyl(spray(cbind(2,1)))
o2 <- weyl(spray(cbind(3,3)))
options(polyform = FALSE)
(1+o1)*(1-5*o2)
## A member of the Weyl algebra:
##  e d     val
##  5 3  =  -15
##  3 3  =   -5
##  5 4  =   -5
##  2 1  =    1
##  0 0  =    1

And of course we can display the result in symbolic form:

options(polyform = TRUE)
(1+o1)*(1-5*o2)
## A member of the Weyl algebra:
## 1 -15*e^5*d^3 -5*e^3*d^3 -5*e^5*d^4 +e^2*d

## References

• S. C. Coutinho 1997. The many avatars of a simple algebra. The American Mathematical Monthly, 104(7):593-604. DOI https://doi.org/10.1080/00029890.1997.11990687.

• K. B. Wolf 1975. The Heisenberg-Weyl ring in quantum mechanics. In Group theory and its applications, Volume III, Academic Press, Inc