diff options
Diffstat (limited to 'cmaes.lua')
| -rw-r--r-- | cmaes.lua | 62 |
1 files changed, 31 insertions, 31 deletions
@@ -1,16 +1,16 @@ require 'torch' require 'math' -local BestSolution = {} ---[[ An implementation of `CMAES` (Covariance Matrix Adaptation Evolution Strategy), +local BestSolution = {} +--[[ An implementation of `CMAES` (Covariance Matrix Adaptation Evolution Strategy), ported from https://www.lri.fr/~hansen/barecmaes2.html. - + Parameters ---------- ARGS: -- `opfunc` : a function that takes a single input (X), the point of - evaluation, and returns f(X) and df/dX. Note that df/dX is not used +- `opfunc` : a function that takes a single input (X), the point of + evaluation, and returns f(X) and df/dX. Note that df/dX is not used - `x` : the initial point - `state.sigma` float, initial step-size (standard deviation in each @@ -20,16 +20,16 @@ ARGS: - `state.ftarget` float, target function value - `state.popsize` - population size. If this is left empty, + population size. If this is left empty, 4 + int(3 * log(|x|)) will be used -- `state.ftarget` +- `state.ftarget` stop if fitness < ftarget - `state.verb_disp` int, display on console every verb_disp iteration, 0 for never RETURN: - `x*` : the new `x` vector, at the optimal point -- `f` : a table of all function values: +- `f` : a table of all function values: `f[1]` is the value of the function before any optimization and `f[#f]` is the final fully optimized value, at `x*` --]] @@ -50,13 +50,13 @@ function optim.cmaes(opfunc, x, config, state) local min_iterations = state.min_iterations or 1 local lambda = state.popsize -- population size, offspring number - -- Strategy parameter setting: Selection + -- Strategy parameter setting: Selection if state.popsize == nil then lambda = 4 + math.floor(3 * math.log(N)) end local mu = lambda / 2 -- number of parents/points for recombination - local weights = torch.range(0,mu-1):apply(function(i) + local weights = torch.range(0,mu-1):apply(function(i) return math.log(mu+0.5) - math.log(i+1) end) -- recombination weights weights:div(weights:sum()) -- normalize recombination weights array local mueff = weights:sum()^2 / torch.pow(weights,2):sum() -- variance-effectiveness of sum w_i x_i @@ -69,18 +69,18 @@ function optim.cmaes(opfunc, x, config, state) local cmu = math.min(1 - c1, 2 * (mueff - 2 + 1/mueff) / ((N + 2)^2 + mueff)) -- and for rank-mu update local damps = 2 * mueff/lambda + 0.3 + cs -- damping for sigma, usually close to 1 - -- Initialize dynamic (internal) state variables + -- Initialize dynamic (internal) state variables local pc = torch.Tensor(N):zero():typeAs(x) -- evolution paths for C local ps = torch.Tensor(N):zero():typeAs(x) -- evolution paths for sigma - local B = torch.eye(N):typeAs(x) -- B defines the coordinate system + local B = torch.eye(N):typeAs(x) -- B defines the coordinate system local D = torch.Tensor(N):fill(1):typeAs(x) -- diagonal D defines the scaling - local C = torch.eye(N):typeAs(x) -- covariance matrix + local C = torch.eye(N):typeAs(x) -- covariance matrix if not pcall(function () torch.symeig(C,'V') end) then -- if error occurs trying to use symeig - error('torch.symeig not available for ' .. x:type() .. + error('torch.symeig not available for ' .. x:type() .. " please use Float- or DoubleTensor for x") end local candidates = torch.Tensor(lambda,N):typeAs(x) - local invsqrtC = torch.eye(N):typeAs(x) -- C^-1/2 + local invsqrtC = torch.eye(N):typeAs(x) -- C^-1/2 local eigeneval = 0 -- tracking the update of B and D local counteval = 0 local f_hist = {[1]=opfunc(x)} -- for bookkeeping output and termination @@ -90,7 +90,7 @@ function optim.cmaes(opfunc, x, config, state) local function ask() - --[[return a list of lambda candidate solutions according to + --[[return a list of lambda candidate solutions according to m + sig * Normal(0,C) = m + sig * B * D * Normal(0,I) --]] -- Eigendecomposition: first update B, D and invsqrtC from C @@ -117,9 +117,9 @@ function optim.cmaes(opfunc, x, config, state) Parameters ---------- - `arx` + `arx` a list of solutions, presumably from `ask()` - `fitvals` + `fitvals` the corresponding objective function values --]] -- bookkeeping, preparation counteval = counteval + lambda -- slightly artificial to do here @@ -142,7 +142,7 @@ function optim.cmaes(opfunc, x, config, state) local c = (cs * (2-cs) * mueff)^0.5 / sigma ps = ps - ps * cs + z * c -- exponential decay on ps - local hsig = (torch.sum(torch.pow(ps,2)) / + local hsig = (torch.sum(torch.pow(ps,2)) / (1-(1-cs)^(2*counteval/lambda)) / N < 2 + 4./(N+1)) hsig = hsig and 1.0 or 0.0 --use binary numbers @@ -155,23 +155,23 @@ function optim.cmaes(opfunc, x, config, state) for i=1,N do for j=1,N do local r = torch.range(1,mu) - r:apply(function(k) + r:apply(function(k) return weights[k] * (arx[k][i]-xold[i]) * (arx[k][j]-xold[j]) end) local Cmuij = torch.sum(r) / sigma^2 -- rank-mu update - C[i][j] = C[i][j] + ((-c1a - cmu) * C[i][j] + + C[i][j] = C[i][j] + ((-c1a - cmu) * C[i][j] + c1 * pc[i]*pc[j] + cmu * Cmuij) end end -- Adapt step-size sigma with factor <= exp(0.6) \approx 1.82 - sigma = sigma * math.exp(math.min(0.6, + sigma = sigma * math.exp(math.min(0.6, (cs / damps) * (torch.sum(torch.pow(ps,2))/N - 1)/2)) end - local function stop() - --[[return satisfied termination conditions in a table like - {'termination reason':value, ...}, for example {'tolfun':1e-12}, - or the empty table {}--]] + local function stop() + --[[return satisfied termination conditions in a table like + {'termination reason':value, ...}, for example {'tolfun':1e-12}, + or the empty table {}--]] local res = {} if counteval > 0 then if counteval >= maxEval then @@ -184,7 +184,7 @@ function optim.cmaes(opfunc, x, config, state) res['condition'] = 1e7 end if fitvals:nElement() > 1 and fitvals[fitvals:nElement()] - fitvals[1] < 1e-12 then - res['tolfun'] = 1e-12 + res['tolfun'] = 1e-12 end if sigma * torch.max(D) < 1e-11 then -- remark: max(D) >= max(diag(C))^0.5 @@ -206,8 +206,8 @@ function optim.cmaes(opfunc, x, config, state) end if iteration <= 2 or iteration % verb_modulo == 0 then local max_std = math.sqrt(torch.max(torch.diag(C))) - print(tostring(counteval).. ': ' .. - string.format(' %6.1f %8.1e ', torch.max(D) / torch.min(D), sigma * max_std) + print(tostring(counteval).. ': ' .. + string.format(' %6.1f %8.1e ', torch.max(D) / torch.min(D), sigma * max_std) .. tostring(fitvals[1])) end @@ -224,7 +224,7 @@ function optim.cmaes(opfunc, x, config, state) fitvals[i] = objfunc(candidate) end - tell(X) + tell(X) disp(verb_disp) end @@ -245,7 +245,7 @@ end -BestSolution.__index = BestSolution +BestSolution.__index = BestSolution function BestSolution.new(x, f, evals) local self = setmetatable({}, BestSolution) self.x = x |