#!/usr/bin/env luajit
local target = assert(tonumber((...)),'Need target value!')
local consts = {
{'pi',math.pi},
{'e',math.exp(1)},
{'phi',(math.sqrt(5)+1)/2},
{'1/2',1/2}}
for i = 1,9 do table.insert(consts, {tostring(i),i}) end
local unary = {
{'ln',math.log,inv='exp'},
{'exp',math.exp,inv='ln'},
{'sqrt',math.sqrt,inv='sqr'},
{'sqr',function(x)return x^2 end,inv='sqrt'},
{'sinpi',function(x)return math.sin(math.pi*x)end,noint=true},
{'cospi',function(x)return math.cos(math.pi*x)end,noint=true},
{'tanpi',function(x)return math.tan(math.pi*x)end,noint=true},
{'sin',math.sin},
{'cos',math.cos},
{'tan',math.tan},
{'-',function(x)return -x end,inv='-'},
{'1/',function(x)return 1/x end,inv='1/'}}
local binary = {
{'logab',math.log,norat=true},
{'atan2',math.atan2,norat=true},
-- {'root',function(a,b)return b^(1/a) end,norat=true},
{'+',function(a,b) return a+b end,comm=true},
{'*',function(a,b) return a*b end,comm=true},
{'^',function(a,b) return a^b end},
{'-',function(a,b) return a-b end,norat=true},
{'/',function(a,b) return a/b end,norat=true}}
local rhs = {}
for _,const in ipairs(consts) do table.insert(rhs,{val=const[2],str=const[1],is='c',c=1,side='r'}) end
local lhs = {{val=target,str='x',is='x',c=1,side='l'}}
local C = 1
local function dist(a,b) return math.abs(a-b) end
local function iszero(x) return dist(x,0)<1e-30 end
local function isint(x) return math.min(dist(x,math.floor(x)),dist(x,math.ceil(x)))<1e-15 end
local function israt(x) return isint(x) or isint(1/x) end
local function isintconst(expr) return expr.is=='c' and isint(expr.val) end
local function form1(expr,op)
local opn,v=op[1],expr.val
local nc = 1+expr.c if nc ~= C then return end
if op.noint and isint(expr.val) then return end
if op.inv == expr.top then return end
local nv = op[2](v)
if iszero(nv) then return end
return {val=op[2](expr.val),str=op[1]..'('..expr.str..')',c=nc,side=expr.side,top=opn} end
local function form2(expr,expr2,op)
local opn,v,v2=op[1],expr.val,expr2.val
local nc = 1+expr.c+expr2.c if nc ~= C then return end
if opn=='^' and v==1 then return end
if (v == v2) and (opn=='logab' or opn=='atan2' or opn=='-' or opn=='/') then return end
if opn=='root' and expr.is=='x' then return end
if op.norat and israt(v/v2) then return end
local nv = op[2](v,v2) if iszero(nv) or (opn=='logab' and israt(nv)) then return end
return {val=nv,c=nc,side=expr.side,
str=#opn==1 and '('..expr.str..')'..opn..'('..expr2.str..')'
or opn..'('..expr.str..','..expr2.str..')'} end
local function cinsert(list,expr)
if expr then local x = expr.val if x==x and x~=x+1 then
table.insert(list,expr) end end end
local function nextstep(exprlist)
local n = #exprlist
for ix=1,n do
local expr=exprlist[ix]
for _,op in ipairs(unary) do
cinsert(exprlist,form1(expr,op))end
for ix2=1,n do local expr2=exprlist[ix2]
for _,op in ipairs(binary) do if ix<=ix2 or not op.comm then
cinsert(exprlist,form2(expr,expr2,op))end end end end end
local function put(lx,rx,d)
local s = ('%s {%d} = {%d} %s'):format(lx.str,lx.c,rx.c,rx.str)
local function sp(n) return (' '):rep(n) end
local sx = sp(30-#lx.str)..s
local sd = ('d = %.2e'):format(d)
print(sx..sp(70-#sx)..sd) end
local threshold=0.0001
for c=1,6 do
C = c nextstep(lhs) nextstep(rhs)
local all = {}
table.move(lhs,1,#lhs,1,all) table.move(rhs,1,#rhs,#lhs+1,all)
table.sort(all,function(a,b) return a.val<b.val end)
for i=1,#all-1 do local ex,ex2=all[i],all[i+1] local d=math.abs(ex.val-ex2.val)
if ex.side ~= ex2.side and d < threshold then
if d>0 then threshold=d end
if ex.side=='l' then put(ex,ex2,d) else put(ex2,ex,d) end
if 0<d and d<1e-15 then goto done end end end end ::done::
post a comment