Module:Temperament data: Difference between revisions
Jump to navigation
Jump to search
CompactStar (talk | contribs) No edit summary |
CompactStar (talk | contribs) No edit summary |
||
Line 104: | Line 104: | ||
function p.pseudoinv(a) | function p.pseudoinv(a) | ||
return | return matmul(transpose(a), matinv(matmul(a, transpose(a)))) | ||
end | end | ||
return p | return p |
Revision as of 08:10, 15 October 2023
Note: Do not invoke this module directly; use the corresponding template instead: Template:Temperament data.
local rat = require('Module:Rational')
local p = {}
local function matadd(a, b)
local result = {}
for i = 1, #a do
result[i] = {}
for j = 1, #(b[1]) do
result[i][j] = a[i][j] + b[i][j]
end
end
return result
end
local function matsub(a, b)
local result = {}
for i = 1, #a do
result[i] = {}
for j = 1, #(b[1]) do
result[i][j] = a[i][j] - b[i][j]
end
end
return result
end
local function matmul(a, b)
local result = {}
for i = 1, #a do
result[i] = {}
for j = 1, #(b[1]) do
result[i][j] = 0
for k = 1, #(a[1]) do
result[i][j] = result[i][j] + (a[i][k] * b[k][j])
end
end
end
return result
end
local function scalarmatmul(a, b)
local result = {}
for i = 1, #a do
result[i] = {}
for j = 1, #(a[1]) do
result[i][j] = a[i][j] * b
end
end
return result
end
function p.matinv(a)
dbl_identity = {}
for i = 1, #a do
dbl_identity[i] = {}
for j = 1, #a do
if i == j then
dbl_identity[i][j] = 2
else
dbl_identity[i][j] = 0
end
end
end
xn = scalarmatmul(a, 0.0001)
for i = 1, 30 do
xn = matmul(xn, matsub(dbl_identity, matmul(a, xn)))
end
return xn
end
local function matinv(a)
local dbl_identity = {}
for i = 1, #a do
dbl_identity[i] = {}
for j = 1, #a do
if i == j then
dbl_identity[i][j] = 2
else
dbl_identity[i][j] = 0
end
end
end
local xn = scalarmatmul(a, 0.0001)
for i = 1, 30 do
xn = matmul(xn, matsub(dbl_identity, matmul(a, xn)))
end
return xn
end
local function transpose(a)
local result = {}
for i = 1, #a[1] do
result[i] = {}
for j = 1, #a do
result[i][j] = a[j][i]
end
end
return result
end
function p.pseudoinv(a)
return matmul(transpose(a), matinv(matmul(a, transpose(a))))
end
return p