Hi, I am facing a problem by trying to optimize a thermodynamic problem with parallel pattersearch. It happens that, by changing the starting point of the optimization, i get a different optimum value of the objective function, and, being just 2 the number of variables, the optimum is always got with a value of the second variable equal to the starting one. Why does this happen? I'd be very grateful for every idea/suggestion.
The code is pretty long, basically it is divided into 3 file: objective funtion, constraint and finally the caller of the optimizer.
Thanks in advance GC
OBJECTIVE FUNCTION
function [ cashflow ] = modello( X )
%MODELLO è lo stesso del main ma per ottimizzare
m_cane=430; %ton/h
p_1t=65.7; %[bar]
t_1t=520; %[C]
p_cond=0.1;
p_reg=X(1); %40
m_rig=X(2); %5
%prepration
[ Pm_desfibr ] = preparaz( m_cane );
%milling
[ P_MOL1T2,P_MOL2T2,P_MOL3T2,P_MOL4T2,P_MOL5T2,P_MOL6T2, P_tandem ] = mill( m_cane );
%electric need
[ P_el_need ] = electric( m_cane, P_tandem, Pm_desfibr );
%fuel heat power
[ Q_bag ] = power_bag( m_cane );
%process demand
Q25_need=103800; %[KW]
%cogeneration bottom
heat_fract_bot=1; %in realtà si assesta sull 80%
[m25_bot, Q25_bot heat_fract_bot, heat_fract_top, m_vap_bot, m_att_bot, P_el_gross_bot, P_el_net_bot, Q_bag_left ] = cog_bot_repower1( m_cane, P_el_need, Q_bag, Q25_need, heat_fract_bot )
%cogeneration topping
[t_8t, t_7t,t_7tsat, m_vap_process_top, m25_top, n_ad_top, x_4t, Q25_top, m_rig, m_vap_top, P_el_gross_top, P_el_net_top, m_boil_new ] = cog_top_repower1( m_cane, P_el_need, Q_bag, Q25_need,p_1t, t_1t, p_reg, p_cond, m_rig, heat_fract_top, Q_bag_left )
%PERFORMANCES
[ P_el_surplus, n_thermal, n_el, n_global, n_II, cashflow ] = performances( Q_bag, P_el_need, P_el_net_bot, P_el_net_top, Q25_need )
cashflow=-cashflow;
end
THIS FUNCTION COMPUTES THE CONSTRAINTS
function [ c, ceq ] = constraint_prova( X )
%CONSTRAINT computes the constraints. here both bottom and topping because
%bottom has to compute things that are required by bottom
Q25_needC=103800;
heat_fract_botC=1;
Q_bagC=2.9563e+05;
p_1tC=65.7;
t_1tC=520;
p_condC=0.1;
p_regC=X(1);
m_rigC=X(2);
%%BOTTOM
%%PARAMETERS
n_mech_turbC=0.98;
n_el_turbC=0.95;
n_hyd_pumpC=0.85;
n_mech_pumpC=0.95;
boiler_t_lossC=0.36;
boiler_p_lossC=0.01;
n_elC=0.95;
n_ad_turb_TGC=0.68;
P_el_nom_botC=8500+3500; %[KW]
%%CYCLE VARIABLE, che però in sto caso non si muovono
p_1C=28.6;
t_1C=410;
t_in_turbC=370;
p_VEC=2.72;
%%CYCLE
%check if we are in saturation condition to calculate properties
t_sat_at_p1C=XSteam('Tsat_p', p_1C); %[C], t sat at 1
if t_1C==t_sat_at_p1C %check to be in vapour phase to calculate properties
h_1C=XSteam('hV_p', p_1C);
s_1C=XSteam('sV_p', p_1C);
else
h_1C=XSteam('h_pT', p_1C, t_1C); %[kj/kg]
s_1C=XSteam('s_pT', p_1C, t_1C); %[kj/kgK] specific entropy
end
% INLET TURBINES
p_2C=p_1C;
t_2C=t_in_turbC;
%check not to be in saturation
t_sat_at_p2C=XSteam('Tsat_p', p_2C); %[C], t sat at 2
if t_2C==t_sat_at_p2C %check to be in vapour phase to calculate properties
h_2C=XSteam('hV_p', p_2C);
s_2C=XSteam('sV_p', p_2C);
else
h_2C=XSteam('h_pT', p_2C, t_2C); %[kj/kg]
s_2C=XSteam('s_pT', p_2C, t_2C); %[kj/kgK] specific entropy
end
%OUTLET TURBINES
p_3_TGC= p_VEC; %[bar]
s_3_iso_TGC=s_2C; %[kj/kgK] specific entropy of the 3 isoentropic
h_3_iso_TGC=XSteam('h_ps', p_3_TGC, s_3_iso_TGC); %[kj/kg] specific enthalpy of point 3 isoentropic
h_3_TGC=h_2C-n_ad_turb_TGC*(h_2C-h_3_iso_TGC); %[kj/kg] specific enthalpy of point 3
t_3_TGC=XSteam('T_ph', p_3_TGC, h_3_TGC); %[C] temperature of point 3
%check, if it is in double phase, i need vapour fraction to be able to
%calculate the entropy
x_3_TGC=XSteam('x_ph', p_3_TGC, h_3_TGC); %steam fraction
if x_3_TGC<1
s_3_TGC=x_3_TGC* XSteam('sV_p', p_3_TGC) + (1-x_3_TGC)*XSteam('sL_p', p_3_TGC); %[kj/kgK] specific entropy, combination of the saturated values
else
s_3_TGC=XSteam('s_ph', p_3_TGC, h_3_TGC); %[kj/kgK] specific entropy
end
w_turb_TG_topC=(h_2C-h_3_TGC); %[kj/kg]
% PUMP
p_4C=p_1C/(1-boiler_p_lossC);
w_pumpC=(0.001*(p_4C-p_3_TGC)*100000/1000)/n_hyd_pumpC; %[kj/kg]
h_4C=XSteam('hL_p', 2.72)+w_pumpC;
%%HEAT, sistema lineare
%calcolo dal Q frct la massa totale di vapore+att, poi ho solo due
%incognite e faccio il sistema lineare
Q25_botC=Q25_needC*heat_fract_botC; %amount of heat provided by bot cycle
m25_botC=Q25_botC/(XSteam('hV_p', p_VEC)-XSteam('hL_p', p_VEC)); %[Kg/s] sia di vapor che di att insieme
h_watC=XSteam('h_pT', p_VEC, 30);
AC=[1 1; h_3_TGC h_watC];
BC=[m25_botC; m25_botC*XSteam('hV_p', p_VEC)];
XXXC=AC\BC;
m_vap_botC=XXXC(1);
m_att_botC=XXXC(2);
%%POTENZE
P_el_gross_botC=(m_vap_botC)*w_turb_TG_topC*n_elC*n_mech_turbC; %[KW]
heat_fract_topC=1-heat_fract_botC;
% CONTROLLO CHE STIA NEI LIMITI DI PORTATA POSSIBILE
if P_el_gross_botC>P_el_nom_botC
P_el_gross_botC=P_el_nom_botC; %[KW]
m_vap_botC=P_el_gross_botC/(n_elC*w_turb_TG_topC*n_mech_turbC); %[kg/s]
m_att_botC=m_vap_botC*(XSteam('hV_p', p_VEC)-h_3_TGC)/(h_watC-XSteam('hV_p', p_VEC));
Q25_botC=(m_vap_botC+m_att_botC)*(XSteam('hV_p', p_VEC)-XSteam('hL_p', p_VEC));
heat_fract_botC=Q25_botC/Q25_needC;
heat_fract_topC=1-heat_fract_botC;
end
P_el_pumpC=m_vap_botC*w_pumpC/(n_elC*n_mech_pumpC); %[KW]
P_el_net_botC=P_el_gross_botC-P_el_pumpC; %[KW]
%%Q BOILER
Q_in_rankine_botC=m_vap_botC*(h_1C-h_4C);
Q_fuel_botC=Q_in_rankine_botC/(1-boiler_t_lossC); %[kW] taken by bottom cycle
Q_bag_leftC=Q_bagC-Q_fuel_botC; %[KW] of fuel still available for topping
%%TOPPING
%%PARAMETRI VARI
n_mech_turbC=0.98;
n_el_turbC=0.95;
n_hyd_pumpC=0.85;
n_mech_pumpC=0.95;
boiler_t_lossC=0.32;
boiler_p_lossC=0.01;
n_elC=0.95;
p_VEC=2.72;
%%THERMODYNAMIC CYCLE
% POINT 1
%check if we are in saturation condition to calculate properties
t_sat_at_p1tC=XSteam('Tsat_p', p_1tC); %[C], t sat at 1
if t_1tC==t_sat_at_p1tC %check to be in vapour phase to calculate properties
h_1tC=XSteam('hV_p', p_1tC);
s_1tC=XSteam('sV_p', p_1tC);
else
h_1tC=XSteam('h_pT', p_1tC, t_1tC); %[kj/kg]
s_1tC=XSteam('s_pT', p_1tC, t_1tC); %[kj/kgK] specific entropy
end
% POINT 3
p_3tC=p_VEC; %[bar]
t_3tC=XSteam('Tsat_p', p_3tC)+5; %[C] 5 degrees of supersat
h_3tC=XSteam('h_pT', p_3tC, t_3tC);
s_3t_isoC=s_1tC;
h_3t_isoC=XSteam('h_ps', p_3tC, s_3t_isoC);
n_ad_topC=(h_1tC-h_3tC)/(h_1tC-h_3t_isoC);
% POINT 2
p_2tC=p_regC;
s_2t_isoC=s_1tC;
h_2t_isoC=XSteam('h_ps', p_2tC, s_2t_isoC);
h_2tC=h_1tC-n_ad_topC*(h_1tC-h_2t_isoC);
t_2tC=XSteam('T_ph', p_2tC, h_2tC);
% POINT 4
p_4tC=p_condC;
s_4t_isoC=s_1tC;
h_4t_isoC=XSteam('h_ps', p_4tC, s_4t_isoC);
h_4tC=h_1tC-n_ad_topC*(h_1tC-h_4t_isoC);
x_4tC=XSteam('x_ph', p_4tC, h_4tC);
% POINT 5
p_5tC=p_4tC;
h_5tC=XSteam('hL_p', p_5tC);
% POINT 5B process output
p_5BtC=p_3tC;
h_5BtC=XSteam('hL_p', p_5BtC);
% POINT 5C rigenerator
p_5CtC=p_2tC;
h_5CtC=h_2tC;
% POINT 6
p_6tC=p_3tC;
% POINT 7
p_7tC=p_2tC;
% POINT 8
p_8tC=p_1tC/(1-boiler_p_lossC);
h_8tC=100;
error_h8C=10;
while error_h8C>1
h_8tC=h_8tC+0.1;
%%LAVORI SPECIFICI
w_pump1C=(0.001*(p_6tC-p_5tC)*100000/1000)/n_hyd_pumpC; %[KJ/kg] pumps
w_pump2C=(0.001*(p_7tC-p_6tC)*100000/1000)/n_hyd_pumpC;
w_pump3C=(0.001*(p_8tC-p_7tC)*100000/1000)/n_hyd_pumpC;
w_turb_hpC=h_1tC-h_2tC; %[Kj/kg]
w_turb_mpC=h_2tC-h_3tC;
w_turb_lpC=h_3tC-h_4tC;
%%MASS BALANCES
Q25_topC=Q25_needC*heat_fract_topC; %[KW]
m25_topC=Q25_topC/(XSteam('hV_p', p_VEC)-XSteam('hL_p', p_VEC)); %[Kg/s] sia di vapor che di att insieme
h_watC=XSteam('h_pT', p_VEC, 30);
Related Question
Best Answer