参考: パソコンで見る天体の動き, pp.112-119
function staticVals(funcs, inits, h) --local function emptyMat(numOfRows) local sute = {} for r = 1, numOfRows do sute[r] = {} end return sute end local dim = #funcs local w = {} local f = {{},{},{},{}} local temp = {} for i = 1, 4 do for j = 1, dim do f[i][j] = funcs[j](unpack(inits[i])) end end for j = 1, dim do w[j] = 2 * inits[4][j] - inits[3][j] + h*h * (204 * f[4][j] + 14 * f[3][j] + 4 * f[2][j] - f[1][j]) / 240 temp[j] = 2 * f[4][j] - f[3][j] end return w, temp -- {変化しない部分}, {假定値} end -- 確かめる。 do -- この聯立二階微分方程式で確かめる。 local function xDotDot(x, y) return -2.95912208e-4 * x / (math.sqrt(x^2 + y^2)^3) end local function yDotDot(x, y) return -2.95912208e-4 * y / (math.sqrt(x^2 + y^2)^3) end -- 4 組の出発値にはこの値を使う。 {{x1, y1}, {x2, y2}, {x3, y3}, {x4, y4}} local inits = {{1.0509143, -0.4038388}, {1.0989720, 0.0 }, {1.0509143, 0.4038388}, {0.9168453, 0.7754147}} local w, temp = staticVals({xDotDot, yDotDot}, inits, 20) print("{"..table.concat(w ,",").."}") -- 変化しない部分 print("{"..table.concat(temp,",").."}") -- 假定値 end