Commit 3e51ed32 authored by gonciarz's avatar gonciarz

Stuff from archive

parent a4dda14d

Too many changes to show.

To preserve performance only 1000 of 1000+ files are displayed.

File added
function xmin=EAforStartPoint % (mu/mu_w, lambda)-CMA-ES
% -------------------- Initialization --------------------------------
% User defined input parameters (need to be edited)
strfitnessfct = 'evaluateRegions'; % name of objective/fitness function
N = 3; % number of objective variables/problem dimension
xmean = rand(N,1); % objective variables initial point
sigma = 0.3; % coordinate wise standard deviation (step size)
stopfitness = 0; % stop if fitness < stopfitness (minimization)
stopeval = 1e3*N^2; % stop after stopeval number of function evaluations
% Strategy parameter setting: Selection
lambda = 4+floor(3*log(N)); % population size, offspring number
mu = lambda/2; % number of parents/points for recombination
weights = log(mu+1/2)-log(1:mu)'; % muXone array for weighted recombination
mu = floor(mu);
weights = weights/sum(weights); % normalize recombination weights array
mueff=sum(weights)^2/sum(weights.^2); % variance-effectiveness of sum w_i x_i
% Strategy parameter setting: Adaptation
cc = (4+mueff/N) / (N+4 + 2*mueff/N); % time constant for cumulation for C
cs = (mueff+2) / (N+mueff+5); % t-const for cumulation for sigma control
c1 = 2 / ((N+1.3)^2+mueff); % learning rate for rank-one update of C
cmu = min(1-c1, 2 * (mueff-2+1/mueff) / ((N+2)^2+mueff)); % and for rank-mu update
damps = 1 + 2*max(0, sqrt((mueff-1)/(N+1))-1) + cs; % damping for sigma
% usually close to 1
% Initialize dynamic (internal) strategy parameters and constants
pc = zeros(N,1); ps = zeros(N,1); % evolution paths for C and sigma
B = eye(N,N); % B defines the coordinate system
D = ones(N,1); % diagonal D defines the scaling
C = B * diag(D.^2) * B'; % covariance matrix C
invsqrtC = B * diag(D.^-1) * B'; % C^-1/2
eigeneval = 0; % track update of B and D
chiN=N^0.5*(1-1/(4*N)+1/(21*N^2)); % expectation of
% ||N(0,I)|| == norm(randn(N,1))
% -------------------- Generation Loop --------------------------------
counteval = 0; % the next 40 lines contain the 20 lines of interesting code
while counteval < stopeval
disp('pitsi');
disp(counteval);
% Generate and evaluate lambda offspring
for k=1:lambda
arx(:,k) = xmean + sigma * B * (D .* randn(N,1)); % m + sig * Normal(0,C) for i=1:N
if(arx(1,k)<0)
arx(1,k)=0;
end
if(arx(1,k)>0.3)
arx(1,k)=0.3;
end
if(arx(2,k)<0)
arx(2,k)=0;
end
if(arx(2,k)>1)
arx(2,k)=1;
end
if(arx(3,k)<0)
arx(3,k)=0;
end
arfitness(k) = feval(strfitnessfct, arx(:,k)); % objective function call
counteval = counteval+1;
end
% Sort by fitness and compute weighted mean into xmean
[arfitness, arindex] = sort(arfitness); % minimization
xold = xmean;
xmean = arx(:,arindex(1:mu))*weights; % recombination, new mean value
% Cumulation: Update evolution paths
ps = (1-cs)*ps ...
+ sqrt(cs*(2-cs)*mueff) * invsqrtC * (xmean-xold) / sigma;
hsig = norm(ps)/sqrt(1-(1-cs)^(2*counteval/lambda))/chiN < 1.4 + 2/(N+1);
pc = (1-cc)*pc ...
+ hsig * sqrt(cc*(2-cc)*mueff) * (xmean-xold) / sigma;
% Adapt covariance matrix C
artmp = (1/sigma) * (arx(:,arindex(1:mu))-repmat(xold,1,mu));
C = (1-c1-cmu) * C ... % regard old matrix
+ c1 * (pc*pc' ... % plus rank one update
+ (1-hsig) * cc*(2-cc) * C) ... % minor correction if hsig==0
+ cmu * artmp * diag(weights) * artmp'; % plus rank mu update
% Adapt step size sigma
sigma = sigma * exp((cs/damps)*(norm(ps)/chiN - 1));
% Decomposition of C into B*diag(D.^2)*B' (diagonalization)
if counteval - eigeneval > lambda/(c1+cmu)/N/10 % to achieve O(N^2)
eigeneval = counteval;
C = triu(C) + triu(C,1)'; % enforce symmetry
[B,D] = eig(C); % eigen decomposition, B==normalized eigenvectors
D = sqrt(diag(D)); % D is a vector of standard deviations now
invsqrtC = B * diag(D.^-1) * B';
end
% Break, if fitness is good enough or condition exceeds 1e14, better termination methods are advisable
if arfitness(1) <= stopfitness || max(D) > 1e7 * min(D)
break;
end
end % while, end generation loop
xmin = arx(:, arindex(1)); % Return best point of last iteration.
% Notice that xmean is expected to be even
% better.
% ---------------------------------------------------------------
function f=frosenbrock(x)
if size(x,1) < 2 ,error('dimension must be greater one'); end
f = 100*sum((x(1:end-1).^2 - x(2:end)).^2) + sum((x(1:end-1)-1).^2);
\ No newline at end of file
This diff is collapsed.
function [Vol] = Vol_lp(N,r,p)
% N - dimension
% r - radius
% p - p-norm
if p > 100
Vol = (2*r).^N;
else
Vol = (2*gamma(1/p+1).*r).^N/(gamma(N/p+1));
end
\ No newline at end of file
function [out_averageCov,warningCell]= averageCov(out,numLast)
p_test = out.opts.pn;
dim = size(out.xAcc,2);
% save estimated volume from MC at each iteration (of the last numLast iter)
Vol_test_all = nan(numLast,1);
% save how close eigenvectors are to first set of eigenvectors
Vec_similar_first = nan(numLast-1,dim);
Vec_similar_previous = nan(numLast-1,dim);
% save average eigenvalues and eigenvectors
eigVal_m_previous = nan(dim,1);
eigVec_m_previous = nan(dim);
eigVal_m_first = nan(dim,1);
eigVec_m_first = nan(dim);
eigValAll_first = nan(numLast,dim);
eigValAll_previous = nan(numLast,dim);
warningCell=cell(numLast,1);
warningCnt = 1;
for k=1:numLast
%numLast
%size(out.QCell)
%get MC Covariance in each step
Q = out.QCell{end-numLast+k};
r = out.rVec(end-numLast+k);
C = r^2*(Q*Q');
%estimated Volume from MC
Vol_test_all(k) = Vol_lp(dim,r,p_test)*out.P_empVecAll(end-numLast+k);
%eigen decomposition
[eigVec,eigVals] = eig(C);
if k==1
%to save average eigenvalues and eigenvectors
eigVal_m_previous=diag(eigVals);
eigVec_m_previous=eigVec;
eigVal_m_first = diag(eigVals);
eigVec_m_first = eigVec;