www.gusucode.com > 几种多目标优化算法集合，包括MOEAD，MOPSO，NNIA，NSGA2等 > MOEAD/moead.m

    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)
    global params
    
    newobj=subobjective([subp.weight], ind.objective,  idealpoint, params.dmethod);
    oops = [subp.curpoint]; 
    oldobj=subobjective([subp.weight], [oops.objective], idealpoint, params.dmethod );
    
    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