projeto_ver2.m

%% fmincon: goal to maximize the power generated by the turbine flow
clear
clc

rng(1);
Vi = 50e6;
Q_in = randi([100,1000],1,50); % inflow
t = 3600;

for i = 1:50
    x0 = [0;0];
    A = [];
    B = [];
    Aeq = [];
    Beq = [];
    LB = [0;6];
    UB = [3;20];
    


    [x,fval,exitflag,output] = fmincon(@(x)func_ver2_1obj(x),x0,A,B,Aeq,Beq,LB,UB,@(x)rest_ver2_1obj(x,Vi));
    
     Vi = Vi + t*(Q_in(i) - pi*x(2)*x(1)^2);
     x
     fval
     exitflag
     
%     if Vi >= 55e6
%         fprintf('The dam just collapsed! RUN!!!')
%     end
end


% n = 0.75; %efficiency
% p = 998; %density of water
% g = 9.81;
% h = 110; %head
% fplot(@(x) (-(n * p * g * h * pi * 6 * x^2)));

%% multiobjective genetic algorithm: goal to maximize the power generated by the turbine flow and minimize the spill flow
clear
clc

rng(1);
Vi = 50e6;
Q_in = randi([100,1000],1,50); % inflow
t = 3600;
j = 0;
for i = 1:1
    
    x0 = [0;0;0];
    A = [];
    B = [];
    Aeq = [];
    Beq = [];
    LB = [0;6;0];
    UB = [3;20;5];

    op = optimoptions(@gamultiobj);
    op = optimoptions(op,'PopulationSize',100,'MaxStallGenerations',500,'MaxGenerations',1000);
    [x,fval,exitflag,output,population] = gamultiobj(@(x)func_ver2(x),3,A,B,Aeq,Beq,LB,UB,@(x)rest_ver2(x,Vi),op);
    
    Vi = Vi + t*(Q_in(i) - pi*x(2)*x(1)^2 -pi*x(2)*x(3)^2);
    x
    exitflag
    if Vi >= 55e6
        fprintf('The dam just collapsed! RUN!!!')
        disp(Vi)
    end
    
    if exitflag ~= 1
        j = j + 1;
    end
   
    fprintf('The number of points on the Pareto front was: %d\n', size(x,1));
    fprintf('The number of generations was : %d\n', output.generations);
    scatter3(x(:,1),x(:,2),x(:,3),'filled')
   
end



%% fminsearch: goal to maximize the power generated by the turbine flow. no constraints

%Remark: no finite solution!!!
clear
clc

x0=[-10,-10];
rng(1);
Vi = 50e6;
Q_in = randi([100,1000],1,50); % inflow
t = 3600;
op = optimset('MaxFunEvals',100000,'MaxIter',100000);%,'PlotFcns',@optimplotfval);

for i = 1:50
    %for j=80:50:2000
    for j=1:100:2000
        xk=fminsearch('func_ver2_pen',x0,op,j);
        fprintf('Penalty->%3i x->(%.5f,%.5f)\n',j,xk(1),xk(2));
        x0=xk;
    end
    
    %[x,fval,exitflag,output] = fminsearch(@(x)func_ver2_1obj(x),x0,op);
    
    Vi = Vi + t*(Q_in(i) - pi*xk(2)*xk(1)^2);
     
    if Vi >= 55e6
        fprintf('The dam just collapsed! RUN!!!')
    end
end

%% patternsearch: goal to maximize the power generated by the turbine flow. 

%Remark: for this algortihm to find feasible points in more runs, the max
%reservoir storage constraint must be removed
clear
clc

rng(1);
Vi = 50e6;
Q_in = randi([10,1000],1,50); % inflow
t = 3600;
j = 0;

x0 = [0;0];
A = [];
B = [];
Aeq = [];
Beq = [];
LB = [0;6];
UB = [3;20];

for i = 1:50
    
    [x,fval,exitflag,output] = patternsearch(@(x)func_ver2_1obj(x),x0,A,B,Aeq,Beq,LB,UB,@(x)rest_ver2_1obj(x,Vi));
    
    Vi = Vi + t*(Q_in(i) - pi*x(2)*x(1)^2);
    x
    fval
    exitflag
    
    if Vi >= 55e6
        fprintf('The dam just collapsed! RUN!!!')
    end
    if exitflag ~= 1
        j = j+1;
    end
end

%% PSwarm
clear
clc

rng(1);
Vi = 50e6;
Q_in = randi([100,1000],1,50); % inflow
t = 3600;

LB = [0;6];
UB = [3;20];

for i = 1:50
    op = optimoptions('particleswarm','SwarmSize',300,'HybridFcn',@fmincon,'MaxStallIterations',500,'MaxIterations',1000);
    [x,fval,exitflag,output] = particleswarm(@(x)func_ver2_1obj(x),2,LB,UB,op);
 
    Vi = Vi + t*(Q_in(i) - pi*x(2)*x(1)^2);
    x
    fval
    exitflag

    
end


%% PSwarm prof

clear
clc

rng(1);
Vi = 50e6;
Q_in = randi([100,1000],1,50); % inflow
t = 3600;

Problem.ObjFunction='func_ver2_1obj';
Problem.LB = [0;6];
Problem.UB = [3;20];
InitPop(1).x=[0;0];
Options.Size = 30;

for i = 1:50
    %op = optimoptions('particleswarm','SwarmSize',300,'HybridFcn',@fmincon,'MaxStallIterations',500,'MaxIterations',1000);
    [x,fx,RunData] = PSwarm(Problem,InitPop,Options);
 
    Vi = Vi + t*(Q_in(i) - pi*x(2)*x(1)^2);
    x

    
end

%% genetic algorithm: goal to maximize the power generated by the turbine flow

%Remark: for this algortihm to find feasible points in more runs, the max
%reservoir storage constraint must be removed
clear
clc

rng(1);
Vi = 50e6;
Q_in = randi([100,1000],1,50); % inflow
t = 3600;

for i = 1:50
    
    x0 = [0;0];
    A = [];
    B = [];
    Aeq = [];
    Beq = [];
    LB = [0;6];
    UB = [3;20];

    op = optimoptions(@ga);
    op = optimoptions(op,'PopulationSize',100,'PlotFcn', {@gaplotbestf,@gaplotrange});
    [x,fval,exitflag,output,population] = ga(@(x)func_ver2_1obj(x),2,A,B,Aeq,Beq,LB,UB,@(x)rest_ver2_1obj(x,Vi));
    
    x
    exitflag
    Vi = Vi + t*(Q_in(i) - pi*x(2)*x(1)^2);
    
    if Vi >= 55e6
        fprintf('The dam just collapsed! RUN!!!')
    end
   
end

%% Evolutionary Algorithm: single objective function
clear
clc

Problem.ObjFunction=@func_ver2_1obj;
Problem.Variables=2;
Problem.LB = [0 6];
Problem.UB = [3 20];
InitPop(1).x=[0 0];
Options.Rho=2;
Options.Mu=10;
Options.Verbosity=2;


[x,fval] = ES(Problem,[],Options)

%%

n = 0.75; %efficiency
p = 998; %density of water
g = 9.81;
h = 110; %head
[x1,x2] = meshgrid(0:3,6:20);
f = (n .* p .* g .* h .* pi .* x2 .* x1.^2);
figure
subplot(1,2,1)
surf(x1,x2,f)
[x2,x3] = meshgrid(6:20,0:5);
f2 = pi .* x2 .* x3.^2;
subplot(1,2,2)
surf(x2,x3,f2)