GEATbx: | Main page Tutorial Algorithms M-functions Parameter/Options Example functions www.geatbx.com |
When using the toolbox the implementation of an objective function consumes most of the work. Inside this function the calculation of the objective values depending on the variables is performed. Here you must implement your problem! The kind of implementation determines how good the evolutionary algorithm can work on and solve the problem.
Included in the distributed version of the toolbox are many examples of objective functions. All included objective functions follow the naming convention obj*.m (see also Naming Convention, Section 5.1). These functions implement a broad class of parameter optimization problems (and a few ordering/scheduling problems). When using these functions as a template it will be easier to implement your own functions/problems. For an overview of the mathematical description see Examples of Objective Functions.
Consider the following tasks:
3.1 Parametric optimization functions |
Let's finish the theory. Here is an example.
Consider implementing the simple quadratic function (sum of quadrate or bowl function; known as De Jong's function 1 - objfun1). This functions works on real value variables.
function objval = objfun1(x) objval=sum((x.*x),2); % vectorized, thus fast % for i=1:size(x,1), objval(i)=sum(x(i,:).^2); end
objfun1 returns the objective values of all individuals and because of the vectorization it will be fast. The second line (commented) shows the implementation of this function unvectorized. Every individual will be computed separately. The result is the same, however, Matlab takes a considerably longer time.
A fully documented version of the above objective function is implemented in objfun1.
Other examples of objective functions (using real value variables, vectorized implementation and a definable number of dimensions) are:
These functions are very often used as a standard set of test functions for evaluation of the performance of different evolutionary algorithms.
3.2 Defining default values of the objective function |
The next step is defining default values for boundaries (domain of variables), best objective value and a description for the function. During the work on the toolbox the following style developed and is implemented inside all provided objective functions: Call the objective function with one special individual only (first variable is NaN, the optional second defines the type of info requested, for example x = [NaN, 1]).
function objval = objfun1(x, option) [Nind,Nvar]=size(x); % if Chrom is [NaN xxx] define size of boundary-matrix and others if all([Nind == 1, isnan(Chrom(1))]), % If only NaN is provided if length(Chrom) == 1, option = 1; else option = Chrom(2); end % Default dimension of objective function Dim = 20; % return text of title for graphic output if option == 2, ObjVal = ['DE JONGs function 1']; % return value of global minimum elseif option == 3, ObjVal = 0; % define size of boundary-matrix and values else % lower and upper bound, identical for all n variables ObjVal = repmat([-512; 512], [1 Dim]); end else % Compute objective function end
The line objfun1([NaN,1]) or objfun1([NaN]) will return the default boundaries (includes implicitly the dimension) of the objective function. This is used inside a number of functions to retrieve the default boundaries and implicitly the number of variables for a given objective function.
objfun='objfun1'; VLUB = feval(objfun, [NaN, 1]); geamain2(objfun, [], VLUB);
Getting the descriptive name of the objective function (objfun1([NaN,2])) or the minimal objective value (objfun1([NaN,3])) is not used in the distributed version of the toolbox (because I can not guaranty that you write all your objective functions according to this scheme). However, when comparing different algorithms and defining a termination criterion against the minimal objective value it is very useful (and I use it quite often for my own work).
A fully documented version of the above objective function is implemented in objfun1.
3.3 Optimization of dynamic systems |
Often the solution of a problem involves the simulation of a system or the call of other functions. Consider the optimization of the control vector of a double integrator (push cart system). An overview of this system is given in Examples of Objective Functions.
function objval=objdopi(Chrom,option) [Nind,Nvar]=size(Chrom); XINIT=[0;-1];XEND=[0;0]; TSTART=0; TEND=1;TIMEVEC=linspace(TSTART,TEND,Nvar)'; if Nind==0, % see above example or objdopi else STEP=abs((TEND-TSTART)/(Nvar-1))); for indrun=1:Nind control=[TIMEVEC [Chrom(indrun,:)]']; [t x]=rk23('simdopiv', [TSTART TEND], XINIT, ... [1e-3,STEP,STEP],control); % Calculate objective function objval(indrun)=sum(abs(x(size(x,1),:)'-XEND))+ ... trapz(t,Chrom(indrun,:).^2)); end end
At the beginning of the calculation of the objective values problem specific parameters are defined (XINIT, XEND, TSTART, TEND, TIMEVEC). The direct calculation is done separately for every individual (rk23 is written for only one system at one time). Every individual is converted in the form required by rk23 (vector of time values in the first column, next column(s) contain control values at specified time). The simulation function is called with appropriate parameters, the state values are returned. (rk23 is a Matlab4 specific function. The calling and setup of the integration routine changed under Matlab5 - see help for sim and simset.)
The objective value is a combination of two terms:
The objective function is finished. However, there is still another function needed - the s-function for the simulation routine simdopiv. (for an introduction to writing s-functions see the Simulink reference guide or look at the provided s-functions of the GEATbx sim*.m).
function [sys,x0]=simdopiv(t,x,u,flag); if abs(flag) == 1 sys(1,:) = u(1,:); sys(2,:) = x(1,:); elseif abs(flag)==0 sys=[2,0,0,1,0,0]; x0 = [0; -1]; end
Let's go one step further. The toolbox provides the possibility to pass up to 10 parameters to the objective function. In the above example it could be useful to have the chance of changing TSTART, TEND, XINIT, XEND from outside and thus optimizing different situations of the double integrator system. The beginning of the function would be changed to:
function objval=objdopi(Chrom,option,TSTART,TEND,XINIT,XEND) [Nind,Nvar]=size(Chrom); if nargin<3, TSTART=0; end if nargin<4, TEND=1; end if nargin<5, XINIT=[0;-1]; end if nargin<6, XEND=[0;0]; end
The problem specific parameters are checked and if not provided set to default values. Thus, if necessary the parameters can be passed to the function or default values will be used. For optimizing objdopi the script would be (implemented in demodopi):
objfun='objdopi'; TSTART=0; TEND=2; XINIT=[0;-2]; XEND=[0;0]; tbxmpga(objfun, [], [], 0, TSTART, TEND, XINIT, XEND)
An even more advanced parameter checking could be done with (example):
if nargin<3, TSTART=[]; end if isempty(TSTART), TSTART=0; end if nargin<4, TEND=[]; end if isempty(TEND), TEND=1; end
The parameter passing mechanism offers the possibility of getting multiple solutions at once without editing the objective function:
for i=1:10, XINIT=[0;-i]; geamain2(objfun,GeaOpt,VLUB,[],0,[],[],XINIT) end
Undefined parameters (TSTART, TEND, XEND) will be set to default values - really useful in day to day work.
One problem remains - vectorization. Many simulation problems are not vectorizable and thus slow. The evolutionary algorithm spends most of the time calculating the objective values. For simulation functions there are actually two problems: The Matlab-provided integration functions (rk23, rk45 or sim) are not vectorized. And it's often difficult to vectorize the s-function (see sim*.m). However, for this example both of these problems are solved. The s-function simdopiv is vectorized and together with the GEATbx a vectorized integration routine, intrk4, is provided. For more information see the documentation of intrk4 and method=12 inside objdopi.
If the objective function computes quite a few temporary results that are not part of the objective value it is good style to return them as additional output parameters.
function [objval, t, x]=objdopi(...)
What is the benefit? This opens the possibility of writing problem specific result plotting or special result computing routines. An example is implemented for the double integrator, which serves at the same time as an example for these features.
For using these features the appropriate GeaOpt parameters (Special.InitDo, Special.InitFunction, Output.StatePlotInterval, Output.StatePlotFunction) must be defined. When calling the evolutionary algorithm Special.InitDo should be set to 1 and Output.StatePlotInterval should be set to 1 or greater. The corresponding *Function options contain the manes of the special initialization or state plot functions. The concept of setting the parameter in the options structure for using (or not) the special feature and defining the name of the special function is very flexible. For every objective function a special function for initialization or state plot can be defined (every problem needs it's own initialization or state plot function). Via parameter the use can be switched on or off. An example of all this is provided with the demo function demodopi, the objective function objdopi, the initialization function initdopi, and the state plot function plotdopi.
3.4 Remark |
Two things should be stated at the end:
GEATbx: | Main page Tutorial Algorithms M-functions Parameter/Options Example functions www.geatbx.com |