Lattices with Pari/GP

I was using Sage before, which is a huge application that
draws together many free and, if you have them, proprietary
and expensive mathematical packages. I noticed that most
of the clever functions I was using came from Pari. I'm
now using GP, which is a front end for Pari and nothing
else.

Here's some code for the parametric extension of scalar
badness that I'm now calling Cangwu badness:

Cangwu metric(plimit,Ek) = {
local(size, H, W, scale);
size = primepi(plimit);
H = log(primes(size))/log(2);
W = matrix(size,size,i,j,(i==j)/H[i]);
scale = H*W^2*H~;
W^2*(1+Ek^2)/scale-W^2*H~*H*W^2/scale^2
}

Cangwu badness(mapping, Ek) = {
local(plimit, operator);
plimit = prime(matsize(mapping)[1]);
operator = Cangwu metric(plimit, Ek);
sqrt(matdet(mapping~ * operator * mapping))
}

See http://x31eq.com/badness.pdf for my best explanation of
what it's for. The mapping is by columns instead of rows
because Pari prefers things that way. As an example, try

Cangwu badness([31;49;72;87;107;115], 1/1200)*1200

which duplicates an example from page 5 of the PDF. Or

Cangwu badness([31,41;49,65;72,95;87,115;107,142;115,152],
1/1200)*1200

(all on one line) for a Miracle from page 9.

One nice thing about Cangwu badness is that it constitutes
a quadratic form, and qualifies as a lattice norm. It
happens that Pari knows about such things, and has a
function for finding equal temperaments:

G = Cangwu metric(11, 1/12e2);
qfminim(G, 0.005, 10, 2)[3]

This gives

[-31 -12 41 72 15 46 27 58 -14 22]

[-49 -19 65 114 24 73 43 92 -22 35]

[-72 -28 95 167 35 107 63 135 -32 51]

[-87 -34 115 202 42 129 76 163 -39 62]

[-107 -42 142 249 52 159 94 201 -48 76]

They're as columns, so the first row is the number of steps
to the octave.

A different algorithm only works with integers, but doesn't
require the cutoff to be specified.

qfminim(round(G*1e10), 0, 10, 0)[3]

gives only the 31 note mapping.

You can also get the best equal temperaments for a given
rank 2 mapping:

meantone = [31, 12; 49, 19; 72, 28; 87, 34; 107, 42]
meantone*qfminim(meantone~*G*meantone, 0.006, 10, 2)[3]

gives

[31 12 -19]

[49 19 -30]

[72 28 -44]

[87 34 -53]

[107 42 -65]

For another trick, find the simplest ETs according to
scalar complexity:

complexity metric(plimit) = {
local(size, H, W);
size = primepi(plimit);
H = log(primes(size))/log(2);
W = matrix(size,size,i,j,(i==j)/H[i]);
W^2 / (H*W^2*H~)
}

complexity(mapping) = {
/* scalar complexity */
local(plimit, operator);
plimit = prime(matsize(mapping)[1]);
operator = complexity metric(plimit);
sqrt(matdet(mapping~ * operator * mapping))
}

meantone*qfminim(meantone~*K*meantone, 3, 10, 2)[3]

gives

[1 2]

[1 3]

[0 4]

[-3 4]

[-7 4]

We can also look at stupidly high prime limits. It turns
out that the error parameter has to be stupidly low to get
what pass for sensible results. The simplest equal
temperaments it'll vouch for over all 2 digit ratios have
5,882,736 and 3,386,123 steps to the octave. They probably
aren't consistent. You can get the mappings from

qfminim(Cangwu metric(100, 1e-9), 0.0008, 10, 2)[3]

Graham