fmincon - penalty function

9 次查看（过去 30 天）

显示更早的评论

Matthias K 2022-3-17

0
链接

此问题的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1673624-fmincon-penalty-function

评论： Matthias K 2022-3-31

采纳的回答： Torsten

在 MATLAB Online 中打开

Hi,

I have the following problem:

I have a data set with half-hourly data for gas production and electricity demand. These are now to be brought together via a CHP (with a certain capacity). The produced gas is to be converted into electricity to cover the electricity demand as well as possible. So I am looking for a power P at each point in time.

To do this, I have calculated the deviation between demand and power in an objective function e and then summed it up.

q = fmincon(@(x) obj_func(x,P,demand,z,P_max,V_level,V_initial,V_size,V_proz,production,BHKW_mode,d),x0,[],[],[],[],lb,ub);
function e = obj_func(x,P,demand,z,P_max,V_level,V_initial,V_size,V_proz,production,BHKW_mode,d)
e = (demand - P).^2;
e = sum(e) + sum(d);

This objective function is then to be minimized with fmincon. This is done and at any time I get : P = demand (why should it not be so).

But now my constraint comes into play.

The gas comes first into a storage and if gas is converted into electricity, the storage content becomes smaller by this quantity.

V_level(1) = V_initial + production(1) - P(1);
for i = 2:z
    V_level(i) = V_level(i-1) + production(i) - P(i);
end
V_proz = 100*V_level/V_size;

And this storage may never be fuller than 100 % and never emptier than 0 %. I tried to realize this with a penalty function, which is added to e.

for i = 1:z
    if V_proz(i) > 100
        d(i) = 100000;
    elseif V_proz(i) < 0
        d(i) = 100000;
    else
        d(i) = 0;
    end
end

But no matter how big I choose the penalty for non-compliance with this condition: Nothing changes in the result, it is not considered by the optimizer at all....

Does anyone have any idea where my error lies or how I can solve it differently?

Thank you for your help, Matthias

0 个评论
显示 -2更早的评论隐藏 -2更早的评论

请先登录，再进行评论。

请先登录，再回答此问题。

采纳的回答

Torsten 2022-3-17

0
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1673624-fmincon-penalty-function#answer_920264

编辑：Torsten 2022-3-17

在 MATLAB Online 中打开

As far as I can see, you could formulate your problem as

min: sum_i (abs(demand(i)-P(i)))

s.c.

Pmin <= P(i) <= Pmax
0 <= V_level(i) <= V_size
V_level(i) = V_level(i-1) + production(i) - P(i)

If the difference to your former objective

min: sum_i (demand(i)-P(i))^2

is acceptable, this can be formulated as a linear optimization problem in the unknown P and V_level.

Use linprog to solve.

9 个评论
显示 7更早的评论隐藏 7更早的评论

Torsten 2022-3-31

编辑：Torsten 2022-3-31

And which of the inputs to patternsearch did you fill ?

A, b, Aeq, beq, lb, ub, nonlcon ?

From the description, patternsearch should turn out to be quite safe, but extremely slow.

Matthias K 2022-3-31

Just lb and ub.

I went back to just determine the values for P. so: x=(P(1),...,P(n))

Patternsearch didn't get along well with my many linear equality constraints, because then it doesn't find an initial point for its search which complies with the l.e. constraints.

请先登录，再进行评论。

更多回答（1 个）

Matt J 2022-3-17

1
链接

此回答的直接链接

https://ww2.mathworks.cn/matlabcentral/answers/1673624-fmincon-penalty-function#answer_919894

Problems

(1) The objective function code you've shown appears to depend on everything except the unknown variables x

(2) Assuming V_proz was supposed to be one of the unknowns, your penalty d is a discontinuous function of it.

10 个评论
显示 8更早的评论隐藏 8更早的评论

Matthias K 2022-3-17

编辑：Matthias K 2022-3-17

在 MATLAB Online 中打开

All variables apart from P, V_level and V_proz, which are calculated from x, are known problem data.

Here is my full code, if you want to look into it...

clear all
clc
 file_data = 'Daten_neu2.xlsx';
%% Read data
% Read data
opts = detectImportOptions(file_data,'Sheet','Results');
opts.VariableUnitsRange = 'A2';
data = readtable(file_data,opts,'Sheet','Daten');
% Assign variables
demand = data{:,2};     %Residual load
t = datenum(data{:,1});
t = t-t(1,1);           %time vektor
z = numel(t);
production = data{:,4}; %gas production
P_max = 300;           % Max capacity [kW]
P_min = P_max * 0.25;  % Min capacity [kW]
V_initial = 50;     % start gas storage level  [%]
V_size = 6;         % Gas storage size [h]
V_size = V_size * P_max/2;  %Gas storage size [kWh]
V_initial = V_initial/100*V_size;
P = zeros(size(t,1),1); % power
V_level = zeros(size(t,1),1);   %gas storage level [kWh]
V_proz = zeros(z,1);        %gas storage level [%]
d = zeros(z,1);
x0 = ones(z,1)*100; % Starting vector for optimization
%% Constraints
lb = zeros(z,1);
ub = ones(z,1);
ub = ub*P_max;
%% Solver
 
q = fmincon(@(x) obj_func(x,P,demand,z,P_max,V_level,V_initial,V_size,V_proz,production,d),x0,[],[],[],[],lb,ub);
%% calculation
P = q;
idx = find(P < P_min);
P(idx)=0;
V_level(1) = V_initial + production(1) - P(1);
for i = 2:z
    V_level(i) = V_level(i-1) + production(i) - P(i);
end
V_proz = 100*V_level/V_size;
for i = 1:z
    if V_proz(i) > 100
        d(i) = 100000;
    elseif V_proz(i) < 0
        d(i) = 100000;
    else
        d(i) = 0;
    end
end
Diff = (demand - P).^2;
SSQ = sum(Diff) + d;
%% Plot results
tiledlayout(2,1)
ax1 = nexttile;
plot(ax1,t,demand,'r');
hold on
plot(t,production,'Color','k');
hold on
plot(t,P,'Color','b');
hold off
ax2 = nexttile;
plot(ax2,t,V_proz,'r');
%% Objective function
function e = obj_func(x,P,demand,z,P_max,V_level,V_initial,V_size,V_proz,production,d)   
P_min = P_max * 0.25;
idx = find(x < P_min);
x(idx)=0;
P = x;
V_level(1) = V_initial + production(1) - P(1);
for i = 2:z
    V_level(i) = V_level(i-1) + production(i) - P(i);
end
V_proz = 100*V_level/V_size;
d = zeros(z,1);
for i = 1:z
    if V_proz(i) > 100
        d(i) = 100000;
    elseif V_proz(i) < 0
        d(i) = 100000;
    else
        d(i) = 0;
    end
end
e = (demand - P).^2;
e = sum(e) + sum(d);
end