Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Reorganises code and removes unused routines #178

Merged
merged 7 commits into from
Nov 28, 2023
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
34 changes: 5 additions & 29 deletions minimisers/DE/deopt.m
Original file line number Diff line number Diff line change
Expand Up @@ -77,24 +77,16 @@
% General Public License can be obtained from the
% Free Software Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [FVr_bestmem,problem] = deopt(fname,problem,problemDefLimits,problemDefCells,controls,S_struct)
function [FVr_bestmem,problem] = deopt(fname,problem,problemDefCells,controls,S_struct)


%function FVr_bestmem = rascal_deopt(fname,problem,PlotIt,controls,S_struct)

%[FVr_bestmem,S_bestval,I_nfeval]
str = struct('I_nc',0,'FVr_ca',0,'I_no',0,'FVr_oa',0);
S_val = repmat(str,S_struct.I_NP,1);
%coder.varsize(S_val(:),[Inf 1],[1 0]);

%-----This is just for notational convenience and to keep the code uncluttered.--------


coder.varsize('problemDef.resample',[Inf,1],[1 0]);
coder.varsize('FVr_bestmem',[1 Inf],[0 1]);
coder.varsize('FVr_bestmemit',[1 Inf],[0 1]);
%coder.varsize('FM_pop',[S_struct.I_NP,2],[1 0]);


stopflag = 0;
I_best_index = 1;
Expand All @@ -111,9 +103,6 @@
I_strategy = S_struct.I_strategy;
I_refresh = S_struct.I_refresh;
I_plotting = S_struct.I_plotting;
%coder.varsize('FM_pop',[20,2],[0 0]);
%FM_pop = zeros(I_NP,2);


%-----Check input variables---------------------------------------------
if (I_NP < 5)
Expand Down Expand Up @@ -141,42 +130,29 @@
FM_pop(k,:) = FVr_minbound + rand(1,I_D).*(FVr_maxbound - FVr_minbound);
end

%FM_popold = zeros(size(FM_pop)); % toggle population
%FVr_bestmemit = zeros(1,2);% best population member in iteration
I_nfeval = 0; % number of function evaluations

%------Evaluate the best member after initialization----------------------
%str = struct('I_nc',0,'FVr_ca',0,'I_no',0,'FVr_oa',0);
% S_MSE.FVr_ca = [];
% S_MSE.I_no = [];
% S_MSE.FVr_oa(1) = [];

str = struct('I_nc',0,'FVr_ca',0,'I_no',0,'FVr_oa',0);
S_val = repmat(str,I_NP,1);


%intrafun(p,problemDef,controls,problemDefCells,problemDefLimits);

coder.varsize('I_best_index',[1 1],[0 0]);
I_best_index = 1; % start with first population member
S_val(1) = fname(FM_pop(I_best_index,:),problem,controls,problemDefCells,problemDefLimits);
S_val(1) = fname(FM_pop(I_best_index,:),problem,controls,problemDefCells);
S_bestval = S_val(1); % best objective function value so far
I_nfeval = I_nfeval + 1;
for k=2:I_NP % check the remaining members
S_val(k) = fname(FM_pop(k,:),problem,controls,problemDefCells,problemDefLimits);
S_val(k) = fname(FM_pop(k,:),problem,controls,problemDefCells);
I_nfeval = I_nfeval + 1;
if (leftWin(S_val(k),S_bestval) == 1)
I_best_index = k; % save its location
S_bestval = S_val(k);
end
end
%val = [0 0];
val = FM_pop(I_best_index,:);
FVr_bestmemit = val; % best member of current iteration

% FVr_bestmemit = FM_pop(I_best_index,:);
% S_bestvalit = S_bestval; % best value of current iteration

FVr_bestmem = FVr_bestmemit; % best member ever

%------DE-Minimization---------------------------------------------
Expand Down Expand Up @@ -281,7 +257,7 @@
FM_origin = FM_pm3;
FM_ui = FM_popold.*FM_mpo + FM_ui.*FM_mui; % crossover
else % either-or-algorithm
if (rand < 0.5); % Pmu = 0.5
if (rand < 0.5) % Pmu = 0.5
FM_ui = FM_pm3 + fWeight*(FM_pm1 - FM_pm2);% differential variation
FM_origin = FM_pm3;
else % use F-K-Rule: K = 0.5(F+1)
Expand All @@ -308,7 +284,7 @@
end
%=====End boundary constraints==========================================

S_tempval = fname(FM_ui(k,:),problem, controls,problemDefCells,problemDefLimits); % check cost of competitor
S_tempval = fname(FM_ui(k,:),problem, controls,problemDefCells); % check cost of competitor
I_nfeval = I_nfeval + 1;
if (leftWin(S_tempval,S_val(k)) == 1)
FM_pop(k,:) = FM_ui(k,:); % replace old vector with new one (for new iteration)
Expand Down
283 changes: 122 additions & 161 deletions minimisers/DE/runDE.m
Original file line number Diff line number Diff line change
@@ -1,169 +1,130 @@
function [problemDef,problem,result] = runDE(problemDef,problemDefCells,problemDefLimits,controls)

[problemDef,~] = fitsetup(problemDef,problemDefCells,problemDefLimits,controls);
F_VTR = controls.targetValue; %Value to reach
I_D = length(problemDef.fitpars);

FVr_minbound = problemDef.fitconstr(:,1)';
FVr_maxbound = problemDef.fitconstr(:,2)';
I_bnd_constr = 1; %1: use bounds as bound constraints, 0: no bound constraints

% I_NP number of population members
I_NP = controls.populationSize;

% I_itermax maximum number of iterations (generations)
I_itermax = controls.numGenerations;

% fWeight DE-stepsize fWeight ex [0, 2]
fWeight = controls.fWeight;

% F_CR crossover probability constant ex [0, 1]
F_CR = controls.crossoverProbability;

% I_strategy 1 --> DE/rand/1:
% the classical version of DE.
% 2 --> DE/local-to-best/1:
% a version which has been used by quite a number
% of scientists. Attempts a balance between robustness
% and fast convergence.
% 3 --> DE/best/1 with jitter:
% taylored for small population sizes and fast convergence.
% Dimensionality should not be too high.
% 4 --> DE/rand/1 with per-vector-dither:
% Classical DE with dither to become even more robust.
% 5 --> DE/rand/1 with per-generation-dither:
% Classical DE with dither to become even more robust.
% Choosing fWeight = 0.3 is a good start here.
% 6 --> DE/rand/1 either-or-algorithm:
% Alternates between differential mutation and three-point-
% recombination.

I_strategy = 5;

% I_refresh intermediate output will be produced after "I_refresh"
% iterations. No intermediate output will be produced
% if I_refresh is < 1
I_refresh = 1;

% I_plotting Will use plotting if set to 1. Will skip plotting otherwise.
I_plotting = 0;

%-----Definition of tolerance scheme--------------------------------------
%-----The scheme is sampled at I_lentol points----------------------------
I_lentol = 50;
FVr_x = linspace(-1,1,I_lentol); %ordinate running from -1 to +1
FVr_lim_up = ones(1,I_lentol); %upper limit is 1
FVr_lim_lo = -ones(1,I_lentol); %lower limit is -1

%Tell compiler abut variable sizes
coder.varsize('S_struct.I_lentol', [Inf 1],[1 0]);
coder.varsize('S_struct.FVr_x', [1 Inf],[0 1]);
coder.varsize('S_struct.FVr_lim_up', [1 Inf],[0 1]);
coder.varsize('S_struct.FVr_lim_lo', [1 Inf],[0 1]);

coder.varsize('S_struct.I_NP', [1 1],[0 0]);
coder.varsize('S_struct.fWeight', [1 1],[0 0]);
coder.varsize('S_struct.F_CR', [1 1],[0 0]);
coder.varsize('S_struct.I_D', [1 1],[0 0]);
coder.varsize('S_struct.FVr_minbound', [1 Inf],[0 1]);
coder.varsize('S_struct.FVr_maxbound', [1 Inf],[0 1]);
coder.varsize('S_struct.I_bnd_constr', [1 1],[0 0]);
coder.varsize('S_struct.I_itermax', [1 1],[0 0]);
coder.varsize('S_struct.F_VTR', [1 1],[0 0]);
coder.varsize('S_struct.I_strategy', [1 1],[0 0]);
coder.varsize('S_struct.I_refresh', [1 1],[0 0]);
coder.varsize('S_struct.I_plotting', [1 1],[0 0]);
coder.varsize('S_struct.FM_pop',[Inf 2],[1 0]);
coder.varsize('S_struct.FVr_bestmem',[1 Inf],[0 1]);
coder.varsize('FVr_bestmem',[1 Inf],[0 1]);

%-----tie all important values to a structure that can be passed along----
S_struct.I_lentol = I_lentol;
S_struct.FVr_x = FVr_x;
S_struct.FVr_lim_up = FVr_lim_up;
S_struct.FVr_lim_lo = FVr_lim_lo;

S_struct.I_NP = I_NP;
S_struct.fWeight = fWeight;
S_struct.F_CR = F_CR;
S_struct.I_D = I_D;
S_struct.FVr_minbound = FVr_minbound;
S_struct.FVr_maxbound = FVr_maxbound;
S_struct.I_bnd_constr = I_bnd_constr;
S_struct.I_itermax = I_itermax;
S_struct.F_VTR = F_VTR;
S_struct.I_strategy = I_strategy;
S_struct.I_refresh = I_refresh;
S_struct.I_plotting = I_plotting;
S_struct.FM_pop = zeros(I_NP,2);
S_struct.FVr_bestmem = [0 0];

[res,problemDef] = deopt(@intrafun,problemDef,problemDefCells,controls,S_struct);
problemDef.fitpars = res;
problemDef = unpackparams(problemDef,controls);
[problem,result] = reflectivityCalculation(problemDef,problemDefCells,controls);

if ~strcmpi(controls.display,'off')
fprintf('Final chi squared is %g\n',problem.calculations.sum_chi);
end

[problemDef,~] = fitsetup(problemDef,problemDefCells,problemDefLimits,controls);
F_VTR = controls.targetValue; %Value to reach
I_D = length(problemDef.fitpars);

FVr_minbound = problemDef.fitconstr(:,1)';
FVr_maxbound = problemDef.fitconstr(:,2)';
I_bnd_constr = 1; %1: use bounds as bound constraints, 0: no bound constraints

% I_NP number of population members
I_NP = controls.populationSize;

% I_itermax maximum number of iterations (generations)
I_itermax = controls.numGenerations;

% fWeight DE-stepsize fWeight ex [0, 2]
fWeight = controls.fWeight;

% F_CR crossover probability constant ex [0, 1]
F_CR = controls.crossoverProbability;

% I_strategy 1 --> DE/rand/1:
% the classical version of DE.
% 2 --> DE/local-to-best/1:
% a version which has been used by quite a number
% of scientists. Attempts a balance between robustness
% and fast convergence.
% 3 --> DE/best/1 with jitter:
% taylored for small population sizes and fast convergence.
% Dimensionality should not be too high.
% 4 --> DE/rand/1 with per-vector-dither:
% Classical DE with dither to become even more robust.
% 5 --> DE/rand/1 with per-generation-dither:
% Classical DE with dither to become even more robust.
% Choosing fWeight = 0.3 is a good start here.
% 6 --> DE/rand/1 either-or-algorithm:
% Alternates between differential mutation and three-point-
% recombination.

I_strategy = 5;

% I_refresh intermediate output will be produced after "I_refresh"
% iterations. No intermediate output will be produced
% if I_refresh is < 1
I_refresh = 1;

% I_plotting Will use plotting if set to 1. Will skip plotting otherwise.
I_plotting = 0;

%-----Definition of tolerance scheme--------------------------------------
%-----The scheme is sampled at I_lentol points----------------------------
I_lentol = 50;
FVr_x = linspace(-1,1,I_lentol); %ordinate running from -1 to +1
FVr_lim_up = ones(1,I_lentol); %upper limit is 1
FVr_lim_lo = -ones(1,I_lentol); %lower limit is -1

%Tell compiler abut variable sizes
coder.varsize('S_struct.I_lentol', [Inf 1],[1 0]);
coder.varsize('S_struct.FVr_x', [1 Inf],[0 1]);
coder.varsize('S_struct.FVr_lim_up', [1 Inf],[0 1]);
coder.varsize('S_struct.FVr_lim_lo', [1 Inf],[0 1]);

coder.varsize('S_struct.I_NP', [1 1],[0 0]);
coder.varsize('S_struct.fWeight', [1 1],[0 0]);
coder.varsize('S_struct.F_CR', [1 1],[0 0]);
coder.varsize('S_struct.I_D', [1 1],[0 0]);
coder.varsize('S_struct.FVr_minbound', [1 Inf],[0 1]);
coder.varsize('S_struct.FVr_maxbound', [1 Inf],[0 1]);
coder.varsize('S_struct.I_bnd_constr', [1 1],[0 0]);
coder.varsize('S_struct.I_itermax', [1 1],[0 0]);
coder.varsize('S_struct.F_VTR', [1 1],[0 0]);
coder.varsize('S_struct.I_strategy', [1 1],[0 0]);
coder.varsize('S_struct.I_refresh', [1 1],[0 0]);
coder.varsize('S_struct.I_plotting', [1 1],[0 0]);
coder.varsize('S_struct.FM_pop',[Inf 2],[1 0]);
coder.varsize('S_struct.FVr_bestmem',[1 Inf],[0 1]);
coder.varsize('FVr_bestmem',[1 Inf],[0 1]);

%-----tie all important values to a structure that can be passed along----
S_struct.I_lentol = I_lentol;
S_struct.FVr_x = FVr_x;
S_struct.FVr_lim_up = FVr_lim_up;
S_struct.FVr_lim_lo = FVr_lim_lo;

S_struct.I_NP = I_NP;
S_struct.fWeight = fWeight;
S_struct.F_CR = F_CR;
S_struct.I_D = I_D;
S_struct.FVr_minbound = FVr_minbound;
S_struct.FVr_maxbound = FVr_maxbound;
S_struct.I_bnd_constr = I_bnd_constr;
S_struct.I_itermax = I_itermax;
S_struct.F_VTR = F_VTR;
S_struct.I_strategy = I_strategy;
S_struct.I_refresh = I_refresh;
S_struct.I_plotting = I_plotting;
S_struct.FM_pop = zeros(I_NP,2);
S_struct.FVr_bestmem = [0 0];

%res = deopt(@intrafun,problemDef,controls,S_struct);

[res,problemDef] = deopt(@intrafun,problemDef,problemDefLimits,problemDefCells,controls,S_struct);
problemDef.fitpars = res;
problemDef = unpackparams(problemDef,controls);
[problem,result] = reflectivityCalculation(problemDef,problemDefCells,controls);

if ~strcmpi(controls.display,'off')
fprintf('Final chi squared is %g\n',problem.calculations.sum_chi);
end

end


function S_MSE = intrafun(p,problemDef,controls,problemDefCells,problemDefLimits)

% S_MSE.I_nc = [];
% S_MSE.FVr_ca = [];
% S_MSE.I_no = [];
% S_MSE.FVr_oa(1) = [];

coder.varsize('S_MSE.I_nc',[1 1],[0 0]);
coder.varsize('S_MSE.FVr_ca',[1 1],[0 0]);
coder.varsize('S_MSE.I_no',[1 1],[0 0]);
coder.varsize('S_MSE.FVr_oa',[1 1],[0 0]);

% data = problem.data;
% x = data(:,1);
% y = data(:,2);
% e = data(:,3);
%
% line = (p(1)*x) + p(2);
%
% fval = sum(((y-line).^2)./e);

problemDef.fitpars = p;
problemDef = unpackparams(problemDef,controls);
[problemDef,~] = reflectivityCalculation(problemDef,problemDefCells,controls);
fval = problemDef.calculations.sum_chi;

S_MSE.I_nc = 0;%no constraints THESE FIRST FEW VALS MAY BE WRONG
S_MSE.FVr_ca = 0;%no constraint array
S_MSE.I_no = 1;%number of objectives (costs)
S_MSE.FVr_oa = fval;

end


function PlotIt(FVr_bestmem,problem)

% problem.fitpars = FVr_bestmem;
% setappdata(0,'problem',problem);

p = FVr_bestmem;

data = problem.data;
x = data(:,1);
y = data(:,2);
e = data(:,3);

line = (p(1)*x) + p(2);

figure(1)
clf;hold on
%errorbar(x,y,e,'bo');
plot(x,line);
function S_MSE = intrafun(p,problemDef,controls,problemDefCells)

coder.varsize('S_MSE.I_nc',[1 1],[0 0]);
coder.varsize('S_MSE.FVr_ca',[1 1],[0 0]);
coder.varsize('S_MSE.I_no',[1 1],[0 0]);
coder.varsize('S_MSE.FVr_oa',[1 1],[0 0]);

problemDef.fitpars = p;
problemDef = unpackparams(problemDef,controls);
[problemDef,~] = reflectivityCalculation(problemDef,problemDefCells,controls);
fval = problemDef.calculations.sum_chi;

S_MSE.I_nc = 0; %no constraints THESE FIRST FEW VALS MAY BE WRONG
S_MSE.FVr_ca = 0; %no constraint array
S_MSE.I_no = 1; %number of objectives (costs)
S_MSE.FVr_oa = fval;

end
Loading