Module:Temperament data

Revision as of 10:43, 4 December 2023 by CompactStar (talk | contribs)
Module documentation[view] [edit] [history] [purge]
Note: Do not invoke this module directly; use the corresponding template instead: Template:Temperament data.
Todo: add documentation

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

local function 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 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

local function pseudoinv(a)
	return matmul(transpose(a), matinv(matmul(a, transpose(a))))
end

local function get_te_generator(subgroup, mapping)
	local w = {}
	for i = 1, #subgroup do
		w[i] = {}
		for j = 1, #subgroup do
			if i == j then
				w[i][j] = math.log(2)/math.log(subgroup[i])
			else
				w[i][j] = 0
			end
		end
	end
	
	local jw = {{}}
	for i = 1, #subgroup do
		jw[1][i] = 1
	end
	local vw = matmul(mapping, w)
	local g = matmul(jw, pseudoinv(vw))
	return g
end

local function gcd(a,b)
  if type(a) == "number" and type(b) == "number" and 
        a == math.floor(a) and b == math.floor(b) then
    if b == 0 then
      return a
    else
      return gcd(b, a % b) -- tail recursion
    end
  else
    error("Invalid argument to gcd (" .. tostring(a) .. "," .. 
          tostring(b) .. ")", 2)
  end
end



function p.get_comma_list(subgroup, mapping)
	local te_generator = get_te_generator(subgroup, mapping)
	-- mw.logObject(te_generator)
	local te_tuning_map = {}
	for i = 1, #subgroup do
		te_tuning_map[i] = te_generator[1][1] * mapping[1][i] + te_generator[1][2] * mapping[2][i]
	end
	
	-- mw.logObject(te_tuning_map)
	
	local commas = {}
	if #subgroup == 3 then
		for i=-10, 10 do
			for j=-10, 10 do
				for k=-10, 10 do
					local comma = {i, j, k}
					local comma_num = 0
					for i = 1,#subgroup do
						comma_num = comma_num + (math.log(subgroup[i])/math.log(2))*comma[i]
					end
					local mapped_comma_num = 0
					for i = 1,#subgroup do
						mapped_comma_num = mapped_comma_num + te_tuning_map[i]*comma[i]
					end
				
					if comma_num >= 0 and mapped_comma_num >= 0 and mapped_comma_num < 0.000000001 then
						local comma_gcd = gcd(gcd(comma[1], comma[2]), comma[3])
						if comma_gcd == 1 or comma_gcd == -1 then
							commas[#commas + 1] = comma
						end
					end
				end
			end
		end
	end
	
	return commas
end

return p