Some Tensor Code

Riemannian Geometry on the Computer

Pick a Riemannian manifold, perhaps a Minkowsy manifold (with signature <− + + +>). Pick a general point on the manifold and then a coordinate system for which the metric tensor is η_ij at the chosen point. (η_ij = δ_ij except for η₀₀ which is −1.) In a neighborhood of the point expand the metric tensor in a taylor series, keeping two or three terms. Each element of this tensor is a Multi Variate Polynomial (mvp). We must compute the matrix inverse of g_ij. Even for matrices we have
(1 − x)⁻¹ = Σ{j = 0, ...} x^j. We want (η − x)⁻¹ in a similar form. Note that x and η don’t commute.
(η − x)⁻¹ = (1 − η⁻¹x)⁻¹η⁻¹ = Σ{j = 0, ...} (η⁻¹x)^jη⁻¹
In the case at hand η⁻¹ = η and so we get:
(η − x)⁻¹ = Σ{j = 0, ...} (ηx)^jη.
Terms beyond the third vanish, because we retain only the terms of degree less than 4. There is no question of convergence. This means that we can compute the contravariant g^ij exactly without even a divide! I find this bizarre! The Scheme routine inv performs the above inversion.

I will use the following definition of the Riemann curvature tensor:
R^μ_ναβ = ∂_αΓ^μ_νβ − ∂_βΓ^μ_να + Γ^μ_σαΓ^σ_νβ − Γ^μ_σβΓ^σ_να

Behavior of some tensor code

g1, g2, g3, g4
gn takes a procedure of n arguments and returns a procedure of n arguments that stores dimⁿ values.
if 0 <= i, j, k, l < dim and the free variables of zot are among i, j, k, l then ((g4 (lambda (i j k l) zot)) i j k l) yields zot. (gN fun) is subsequently better than fun merely for the economy of not recomputing the values.

Operators for multivariate polynomials

((padd n k) p1 p2) adds two polynomials, p1 & p2, in n variables. The inputs must be of at least degree k−1 and the yield will be of degree k−1. ((pmul n k) p1 p2) multiplies polynomials discarding terms of order higher than k−1. (pgen n k ng z) generates a polynomial in n variables of degree k−1. It uses the yield of the 0 argument procedure ng for an element, except for the constant term, which is z.

Tentative embodiment: A degree n polynomial is a vector of the n+1 coefficients where element 0 is the constant coefficient. A polynomial in k variables has coefficients which are polynomials in k−1 variables. By the “degree” of a multivariate polynomial, we mean the maximum, over the terms, of the sums of the exponents of the variables.

Example

(pgen 3 5 (let ((x 0))(lambda()(set! x (+ x 1)) x)) −9) yields
#(#(#(−9 34 33 32 31) #(27 30 29 28) #(24 26 25) #(22 23) #(21))
#(#(11 20 19 18) #(15 17 16) #(13 14) #(12)) #(#(5 10 9) #(7 8) #(6))
#(#(2 4) #(3)) #(#(1)))

See the code.

((iter n)(lambda (k) (... k ...))) performs (... k ...) for k descending from n−1 thru 0.
(ip e a b) returns the inner product of two 1D mvps; e terms.
(mm a b) is the matrix product of 2 2D matrices.
(ftl g) turns a 1D vector into a list.
(curry g) curries the function g. If f = (curry g) then ((f x) y z) = (g x y z).
((ml n) g) returns an n deep list of values from a n-ary function g.

Alas the program reports that Bianchi’s identity is false. I have not had time to find the bug.

A classic paper by Riemann