eval_update eval字符串调用函数
function [idealpoint, subproblems]= eval_update(idealpoint, subproblems, inds)
%EvaluationUpdate Summary of this function goes here
% Detailed explanation goes here
objs = [inds.objective];
weights = [subproblems.weight];
idealpoint = min(idealpoint, min(objs,[],2));
for i=1:length(inds)
subobjs = subobjective(weights, objs(:,i), idealpoint, 'ws');
%update the values.
C = subobjs<[subproblems.optimal];
if any(C)
ncell = num2cell(subobjs(C));
[subproblems(C).optimal] = ncell{:};
[subproblems(C).optpoint] = deal(inds(i));
end
end
end
evaluate
function [v, x] = evaluate( prob, x )
%EVALUATE function evaluate an individual structure of a vector point with
%the given multiobjective problem.
% Detailed explanation goes here
% prob: is the multiobjective problem.
% x: is a vector point, or a individual structure.
% v: is the result objectives evaluated by the mop.
% x: if x is a individual structure, then x's objective field is modified
% with the evaluated value and pass back.
% TODO, need to refine it to operate on a vector of points.
if isstruct(x)
v = prob.func(x.parameter);
x.objective=v;
else
v = prob.func(x);
end
gaussian_mutate
function ind = gaussian_mutate( ind, prob, domain)
%GAUSSIAN_MUTATE Summary of this function goes here
% Detailed explanation goes here
if isstruct(ind)
x = ind.parameter;
else
x = ind;
end
parDim = length(x);
lowend = domain(:,1);
highend =domain(:,2);
sigma = (highend-lowend)./20;
newparam = min(max(normrnd(x, sigma), lowend), highend);
C = rand(parDim, 1)<prob;
x(C) = newparam(C);
if isstruct(ind)
ind.parameter = x;
else
ind = x;
end
genetic_op
function ind=genetic_op(subproblems, index, domain, params)
%GENETICOP function implemented the DE operation to generate a new
%individual from a subproblems and its neighbours.
% subproblems: is all the subproblems.
% index: the index of the subproblem need to handle.
% domain: the domain of the origional multiobjective problem.
% ind: is an individual structure.
neighbourindex = subproblems(index).neighbour;
%The random draw from the neighbours.
nsize = length(neighbourindex);
si = ones(1,3)*index;
si(1)=neighbourindex(ceil(rand*nsize));
while si(1)==index
si(1)=neighbourindex(ceil(rand*nsize));
end
si(2)=neighbourindex(ceil(rand*nsize));
while si(2)==index || si(2)==si(1)
si(2)=neighbourindex(ceil(rand*nsize));
end
si(3)=neighbourindex(ceil(rand*nsize));
while si(3)==index || si(3)==si(2) || si(3)==si(1)
si(3)=neighbourindex(ceil(rand*nsize));
end
%retrieve the individuals.
points = [subproblems(si).curpoint];
selectpoints = [points.parameter];
oldpoint = subproblems(index).curpoint.parameter;
parDim = size(domain, 1);
jrandom = ceil(rand*parDim);
randomarray = rand(parDim, 1);
deselect = randomarray<params.CR;
deselect(jrandom)=true;
newpoint = selectpoints(:,1)+params.F*(selectpoints(:,2)-selectpoints(:,3));
newpoint(~deselect)=oldpoint(~deselect);
%repair the new value.
newpoint=max(newpoint, domain(:,1));
newpoint=min(newpoint, domain(:,2));
ind = struct('parameter',newpoint,'objective',[], 'estimation',[]);
%ind.parameter = newpoint;
%ind = realmutate(ind, domain, 1/parDim);
ind = gaussian_mutate(ind, 1/parDim, domain);
%clear points selectpoints oldpoint randomarray deselect newpoint neighbourindex si;
end
get_structure
function str = get_structure( name )
%STRUCTURE Summary of this function goes here
%
% Structure used in this toolbox.
%
% individual structure:
% parameter: the parameter space point of the individual. it's a column-wise
% vector.
% objective: the objective space point of the individual. it's column-wise
% vector. It only have value after evaluate function is called upon the
% individual.
% estimation: Also a structure array of the individual. It's not used in
% MOEA/D but used in MOEA/D/GP. For every objective, the field contains the
% estimation from the GP model.
%
% estimation structure:
% obj: the estimated mean.
% std: the estimated standard deviation for the mean.
%
% subproblem structure:
% weight: the decomposition weight for the subproblem.
% optimal: the current optimal value of the current structure.
% curpoiont: the current individual of the subproblem.
% optpoint: the point that gain the optimal on the subproblem.
%
switch name
case 'individual'
str = struct('parameter',[],'objective'[],'estimation'[]);
case 'subproblem'
str = struct('weight',[],'optimal',[],'curpoint',[],'optpoint',[]);
case 'estimation'
str = struct();
otherwise
end
init_weights
function subp=init_weights(popsize, niche, objDim)
% init_weights function initialize a pupulation of subproblems structure
% with the generated decomposition weight and the neighbourhood
% relationship.
subp=[];
for i=0:popsize
if objDim==2
p=struct('weight',[],'neighbour',[],'optimal', Inf, 'optpoint',[], 'curpoint', []);
weight=zeros(2,1);
weight(1)=i/popsize;
weight(2)=(popsize-i)/popsize;
p.weight=weight;
subp=[subp p];
elseif objDim==3
%TODO
end
end
% weight = lhsdesign(popsize, objDim, 'criterion','maximin', 'iterations', 1000)';
% p=struct('weight',[],'neighbour',[],'optimal', Inf, 'optpoint',[], 'curpoint', []);
% subp = repmat(p, popsize, 1);
% cells = num2cell(weight);
% [subp.weight]=cells{:};
%Set up the neighbourhood.
leng=length(subp);
distanceMatrix=zeros(leng, leng);
for i=1:leng
for j=i+1:leng
A=subp(i).weight;B=subp(j).weight;
distanceMatrix(i,j)=(A-B)'*(A-B);
distanceMatrix(j,i)=distanceMatrix(i,j);
end
[s,sindex]=sort(distanceMatrix(i,:));
subp(i).neighbour=sindex(1:niche)';
end
end
moead
function pareto = moead( mop, varargin)
%MOEAD runs moea/d algorithms for the given mop.
% Detailed explanation goes here
% The mop must to be minimizing.
% The parameters of the algorithms can be set through varargin. including
% popsize: The subproblem's size.
% niche: the neighboursize, must less then the popsize.
% iteration: the total iteration of the moead algorithms before finish.
% method: the decomposition method, the value can be 'ws' or 'ts'.
starttime = clock;
%global variable definition.
global params idealpoint objDim parDim itrCounter;
%set the random generator.
rand('state',10);
%Set the algorithms parameters.
paramIn = varargin;
[objDim, parDim, idealpoint, params, subproblems]=init(mop, paramIn);
itrCounter=1;
while ~terminate(itrCounter)
tic;
subproblems = evolve(subproblems, mop, params);
disp(sprintf('iteration %u finished, time used: %u', itrCounter, toc));
itrCounter=itrCounter+1;
end
%display the result.
pareto=[subproblems.curpoint];
pp=[pareto.objective];
scatter(pp(1,:), pp(2,:));
disp(sprintf('total time used %u', etime(clock, starttime)));
end
function [objDim, parDim, idealp, params, subproblems]=init(mop, propertyArgIn)
%Set up the initial setting for the MOEA/D.
objDim=mop.od;
parDim=mop.pd;
idealp=ones(objDim,1)*inf;
%the default values for the parameters.
params.popsize=100;params.niche=30;params.iteration=100;
params.dmethod='ts';
params.F = 0.5;
params.CR = 0.5;
%handle the parameters, mainly about the popsize
while length(propertyArgIn)>=2
prop = propertyArgIn{1};
val=propertyArgIn{2};
propertyArgIn=propertyArgIn(3:end);
switch prop
case 'popsize'
params.popsize=val;
case 'niche'
params.niche=val;
case 'iteration'
params.iteration=val;
case 'method'
params.dmethod=val;
otherwise
warning('moea doesnot support the given parameters name');
end
end
subproblems = init_weights(params.popsize, params.niche, objDim);
params.popsize = length(subproblems);
%initial the subproblem's initital state.
inds = randompoint(mop, params.popsize);
[V, INDS] = arrayfun(@evaluate, repmat(mop, size(inds)), inds, 'UniformOutput', 0);
v = cell2mat(V);
idealp = min(idealp, min(v,[],2));
%indcells = mat2cell(INDS, 1, ones(1,params.popsize));
[subproblems.curpoint] = INDS{:};
clear inds INDS V indcells;
end
function subproblems = evolve(subproblems, mop, params)
global idealpoint;
for i=1:length(subproblems)
%new point generation using genetic operations, and evaluate it.
ind = genetic_op(subproblems, i, mop.domain, params);
[obj,ind] = evaluate(mop, ind);
%update the idealpoint.
idealpoint = min(idealpoint, obj);
%update the neighbours.
neighbourindex = subproblems(i).neighbour;
subproblems(neighbourindex)=update(subproblems(neighbourindex),ind, idealpoint);
%clear ind obj neighbourindex neighbours;
clear ind obj neighbourindex;
end
end
function subp =update(subp, ind, idealpoint)
newobj=subobjective([subp.weight], ind.objective, idealpoint, 'te');
oops = [subp.curpoint];
oldobj=subobjective([subp.weight], [oops.objective], idealpoint, 'te' );
C = newobj < oldobj;
[subp(C).curpoint]= deal(ind);
clear C newobj oops oldobj;
end
function y =terminate(itrcounter)
global params;
y = itrcounter>params.iteration;
end
randompoint
function ind = randompoint(prob, n)
%RANDOMNEW to generate n new point randomly from the mop problem given.
if (nargin==1)
n=1;
end
randarray = rand(prob.pd, n);
lowend = prob.domain(:,1);
span = prob.domain(:,2)-lowend;
point = randarray.*(span(:,ones(1, n)))+ lowend(:,ones(1,n));
cellpoints = num2cell(point, 1);
indiv = struct('parameter',[],'objective',[], 'estimation', []);
ind = repmat(indiv, 1, n);
[ind.parameter] = cellpoints{:};
% estimation = struct('obj', NaN ,'std', NaN);
% [ind.estimation] = deal(repmat(estimation, prob.od, 1));
end
realmutate
function ind = realmutate(ind, domains, rate)
%REALMUTATE Summary of this function goes here
% Detailed explanation goes here
% double rnd, delta1, delta2, mut_pow, deltaq;
% double y, yl, yu, val, xy;
% double eta_m = id_mu;
eta_m=20;
numVariables = size(domains,1);
if (isstruct(ind))
a = ind.parameter;
else
a = ind;
end
for j = 1:numVariables
if (rand() <= rate)
y = a(j);
yl = domains(j,1);
yu = domains(j,2);
delta1 = (y - yl) / (yu - yl);
delta2 = (yu - y) / (yu - yl);
rnd = rand();
mut_pow = 1.0 / (eta_m + 1.0);
if (rnd <= 0.5)
xy = 1.0 - delta1;
val = 2.0 * rnd + (1.0 - 2.0 * rnd) * (xy^(eta_m + 1.0));
deltaq = (val^mut_pow) - 1.0;
else
xy = 1.0 - delta2;
val = 2.0 * (1.0 - rnd) + 2.0 * (rnd - 0.5) * (xy^ (eta_m + 1.0));
deltaq = 1.0 - (val^mut_pow);
end
y = y + deltaq * (yu - yl);
if (y < yl)
y = yl;
end
if (y > yu)
y = yu;
end
a(j) = y;
end
end
if isstruct(ind)
ind.parameter = a;
else
ind = a;
end
end
subobjective
function obj = subobjective(weight, ind, idealpoint, method)
%SUBOBJECTIVE function evaluate a point's objective with a given method of
%decomposition.
% Two method are implemented by far is Weighted-Sum and Tchebesheff.
% weight: is the decomposition weight.(column wise vector).
% ind: is the individual point(column wise vector).
% idealpoint: the idealpoint for Tchebesheff decomposition.
% method: is the decomposition method, the default is 'te' when is
% omitted.
%
% weight and ind can also be matrix. in which have two scenairos:
% When weight is a matrix, then it's treated as a column wise set of
% weights. in that case, if ind is a size 1 column vector, then the
% subobjective is computed with every weight and the ind; if ind is also
% a matrix of the same size as weight, then the subobjective is computed
% in a column-to-column, with each column of weight computed against the
% corresponding column of ind.
% A row vector of subobjective is return in both case.
if (nargin==2)
obj = ws(weight, ind);
elseif (nargin==3)
obj = te(weight, ind, idealpoint);
else
if strcmp(method, 'ws')
obj=ws(weight, ind);
elseif strcmp(method, 'te')
obj=te(weight, ind, idealpoint);
else
obj= te(weight, ind, idealpoint);
end
end
end
function obj = ws(weight, ind)
obj = (weight'*ind)';
end
function obj = te(weight, ind, idealpoint)
s = size(weight, 2);
indsize = size(ind,2);
weight((weight == 0))=0.00001;
if indsize==s
part2 = abs(ind-idealpoint(:,ones(1, indsize)));
obj = max(weight.*part2);
elseif indsize ==1
part2 = abs(ind-idealpoint);
obj = max(weight.*part2(:,ones(1, s)));
else
error('individual size must be same as weight size, or equals 1');
end
end
testmop
function mop = testmop( testname, dimension )
%Get test multi-objective problems from a given name.
% The method get testing or benchmark problems for Multi-Objective
% Optimization. The implemented problems included ZDT, OKA, KNO.
% User get the corresponding test problem, which is an instance of class
% mop, by passing the problem name and optional dimension parameters.
mop=struct('name',[],'od',[],'pd',[],'domain',[],'func',[]);
switch lower(testname)
case 'kno1'
mop=kno1(mop);
case 'zdt1'
mop=zdt1(mop, dimension);
otherwise
error('Undefined test problem name');
end
end
%KNO1 function generator
function p=kno1(p)
p.name='KNO1';
p.od = 2;
p.pd = 2;
p.domain= [0 3;0 3];
p.func = @evaluate;
%KNO1 evaluation function.
function y = evaluate(x)
y=zeros(2,1);
c = x(1)+x(2);
f = 9-(3*sin(2.5*c^0.5) + 3*sin(4*c) + 5 *sin(2*c+2));
g = (pi/2.0)*(x(1)-x(2)+3.0)/6.0;
y(1)= 20-(f*cos(g));
y(2)= 20-(f*sin(g));
end
end
%ZDT1 function generator
function p=zdt1(p,dim)
p.name='ZDT1';
p.pd=dim;
p.od=2;
p.domain=[zeros(dim,1) ones(dim,1)];
p.func=@evaluate;
%KNO1 evaluation function.
function y=evaluate(x)
y=zeros(2,1);
y(1) = x(1);
su = sum(x)-x(1);
g = 1 + 9 * su / (dim - 1);
y(2) =g*(1 - sqrt(x(1) / g));
end
end
147
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
