和声搜索 HarmonySearch

和声搜索 HarmonySearch matlab

function [BestGenBestFitnessgx]=HarmonySearch

% This code has been written with Matlab 7.0
% You can modify the simple constraint handlening method using more efficient
% methods. A good review of these methods can be found in:
% Carlos A. Coello Coello“Theoretical and numerical constraint-handling techniques used with evolutionary algorithms: a survey of the state of the art“

clc
clear;
global NVAR NG NH MaxItr HMS HMCR PARmin PARmax bwmin bwmax;
global HM NCHV fitness PVB BW gx;
global BestIndex WorstIndex BestFit WorstFit BestGen currentIteration;


NVAR=4;         %number of variables
NG=6;           %number of ineguality constraints
NH=0;           %number of eguality constraints
MaxItr=5000;    % maximum number of iterations
HMS=6;          % harmony memory size
HMCR=0.9;       % harmony consideration rate  0< HMCR <1
PARmin=0.4;      % minumum pitch adjusting rate
PARmax=0.9;      % maximum pitch adjusting rate
bwmin=0.0001;    % minumum bandwidth
bwmax=1.0;      % maxiumum bandwidth
PVB=[1.0 4;0.6 2;40 80;20 60];   % range of variables

% /**** Initiate Matrix ****/
HM=zeros（HMSNVAR）;
NCHV=zeros（1NVAR）;
BestGen=zeros（1NVAR）;
fitness=zeros（1HMS）;
BW=zeros（1NVAR）;
gx=zeros（1NG）;
% warning off MATLAB:m_warning_end_without_block

MainHarmony;

% /**********************************************/

    function sum =Fitness（sol）
        
        sum = 0.6224*sol（1）*sol（3）*sol（4）+1.7781*sol（2）*sol（3）^2+3.1661*sol（1）^2*sol（4）+19.84*sol（1）^2*sol（3）+ eg（sol）;  %F（x） = f（x） + penalty  
        
    end

% /*********************************************/
    function sum=eg（sol）
        
        % constraints g（x） > 0
        gx（1）=sol（1）-0.0193*sol（3）;     %  x1 - 0.0193 x3 > 0
        gx（2）=sol（2）-0.00954*sol（3）;
        gx（3）=3.14*sol（3）^2*sol（4）+（4/3）*3.14*sol（3）^3 - 1296000;
        gx（4）=-sol（4）+240;
        gx（5）=sol（1） - 1.1;
        gx（6）=sol（2） - 0.6;
        
        % we use static penalty function to handle constraints
        sum = 0;
        for i=1:NG
            if（gx（i）<0）
                sum = sum - 1000 * gx（i）;
            end
        end
    end

% /*********************************************/

    function initialize
        % randomly initialize the HM
        for i=1:HMS
            for j=1:NVAR
                HM（ij）=randval（PVB（j1）PVB（j2））;
            end
            fitness（i） = Fitness（HM（i:））;
        end
    end

%/*******************************************/

    function MainHarmony
        % global NVAR NG NH MaxItr HMS HMCR PARmin PARmax bwmin bwmax;
        % global HM NCHV fitness PVB BW gx currentIteration;
        
        initialize;
        currentIteration  = 0;
        
        while（StopCondition（currentIteration））
            
            PAR=（PARmax-PARmin）/（MaxItr）*currentIteration+PARmin;
            coef=log（bwmin/bwmax）/MaxItr;
            for pp =1:NVAR
                BW（pp）=bwmax*exp（coef*curr

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----
     文件        6214  2010-09-29 03:56  HarmonySearch.m
     文件        1333  2010-09-29 03:56  license.txt