Описание задачи

Рынок электроэнергии имеет различные цены в разное время дня. Если у вас есть генераторы, предоставляющие электричество, можно использовать в своих интересах эту переменную оценку путем планирования генераторов, чтобы действовать, когда цены высоки. Предположим, что вы управляете двумя генераторами. Каждый генератор имеет три уровня мощности (прочь, низко, и высоко). Каждый генератор имеет заданный уровень расхода топлива и выработки энергии на каждом уровне мощности. Расход топлива 0, когда генератор выключен.

Можно присвоить уровень мощности каждому генератору для каждого получасового временного интервала в течение дня (24 часа, таким образом, 48 интервалов). На основе хронологических записей примите, что вы знаете доход на мегаватт-час (МВт·ч), который вы получаете в каждом временном интервале. Данными для этого примера является от австралийского Оператора Энергетического рынка https://www.nemweb.com.au/REPORTS/CURRENT/ в середине 2013, и используется в их термины https://www.aemo.com.au/privacy-and-legal-notices/copyright-permissions.

load dispatchPrice; % Get poolPrice, which is the revenue per MWh
bar(poolPrice,.5)
xlim([.5,48.5])
xlabel('Price per MWh at each period')

Существует стоимость, чтобы запустить генератор после того, как это было выключено. Кроме того, существует ограничение на максимальное топливное использование в течение дня. Это ограничение существует, потому что вы покупаете свое топливо день загодя, таким образом, можно использовать только, что вы только купили.

Проблемное обозначение и параметры

Можно сформулировать проблему планирования как бинарную проблему целочисленного программирования. Задайте индексирует iJ, и k, и двоичный файл, планируя векторный y, можно следующим образом:

Определите, когда генератор запустится будучи прочь. Для этого задайте вспомогательную бинарную переменную z(i,j) это указывает, взимать ли за включение генератора j в период i.

Вам нужен способ установить z автоматически на основе настроек y. Линейное ограничение ниже обрабатывает эту установку.

Вам также нужны параметры проблемы для затрат, уровней генерации для каждого генератора, уровней потребления генераторов и доступного топлива.

Вы получили poolPrice когда вы выполнили load dispatchPrice;. Установите другие параметры можно следующим образом.

fuelPrice = 3;
totalFuel = 3.95e4;
nPeriods = length(poolPrice); % 48 periods
nGens = 2; % Two generators
gen = [61,152;50,150]; % Generator 1 low = 61 MW, high = 152 MW
fuel = [427,806;325,765]; % Fuel consumption for generator 2 is low = 325, high = 765
startCost = 1e4; % Cost to start a generator after it has been off

КПД генератора

Исследуйте КПД этих двух генераторов в их двух рабочих точках.

efficiency = gen./fuel; % Calculate electricity per unit fuel use
rr = efficiency'; % for plotting
h = bar(rr);
h(1).FaceColor = 'g';
h(2).FaceColor = 'c';
legend(h,'Generator 1','Generator 2','Location','NorthEastOutside')
ax = gca;
ax.XTick = [1,2];
ax.XTickLabel = {'Low','High'};
ylim([.1,.2])
ylabel('Efficiency')

Заметьте, что генератор 2 немного более эффективен, чем генератор 1 в его соответствующих рабочих точках (низко и высоко), но генератор 1 в его высокой рабочей точке более эффективен, чем генератор 2 в его низкой рабочей точке.

Переменные для решения

Чтобы настроить проблему, необходимо закодировать все проблемные данные и ограничения в проблемной форме. Переменные y(i,j,k) представляйте решение проблемы и вспомогательные переменные z(i,j) укажите, взимать ли за включение генератора. y nPeriods-by-nGens-by-2 массив и z nPeriods-by-nGens массив. Все переменные являются двоичным файлом.

y = optimvar('y',nPeriods,nGens,{'Low','High'},'Type','integer','LowerBound',0,...
    'UpperBound',1);
z = optimvar('z',nPeriods,nGens,'Type','integer','LowerBound',0,...
    'UpperBound',1);

Линейные ограничения

Чтобы гарантировать, что уровень мощности имеет не больше, чем один компонент, равный 1, устанавливает линейное ограничение неравенства.

powercons = y(:,:,'Low') + y(:,:,'High') <= 1;

Производственные затраты на период являются стоимостью топлива в течение того периода. Для генератора j работа с уровнем k, стоимостью является fuelPrice * fuel(j,k).

Создайте выражение fuelUsed это составляет все используемое топливо.

yFuel = zeros(nPeriods,nGens,2);
yFuel(:,1,1) = fuel(1,1); % Fuel use of generator 1 in low setting
yFuel(:,1,2) = fuel(1,2); % Fuel use of generator 1 in high setting
yFuel(:,2,1) = fuel(2,1); % Fuel use of generator 2 in low setting
yFuel(:,2,2) = fuel(2,2); % Fuel use of generator 2 in high setting

fuelUsed = sum(sum(sum(y.*yFuel)));

Ограничение состоит в том, что используемое топливо является не больше, чем доступным топливом.

fuelcons = fuelUsed <= totalFuel;

Установите переменные индикатора запуска генератора

Как можно заставить решатель устанавливать z переменные автоматически, чтобы совпадать с активными/от периодами y переменные? Вспомните, что условием удовлетворить является z(i,j) = 1 точно, когда sum_k y(i,j,k) = 0 и sum_k y(i+1,j,k) = 1.

Заметьте, что sum_k ( - y(i,j,k) + y(i+1,j,k) ) > 0 точно, когда это необходимо, z(i,j) = 1.

Поэтому включайте эти линейные ограничения неравенства в формулировку задачи.

sum_k ( - y(i,j,k) + y(i+1,j,k) ) - z(i,j) < = 0.

Кроме того, включайте z переменные в целевой функции стоятся. С z переменные в целевой функции, решатель пытается понизить их значения, означая, что это пытается установить их всех, равняются 0. Но для тех интервалов, когда генератор включает, линейное неравенство обеспечивает z(i,j) равняться 1.

Создайте вспомогательную переменную w это представляет y(i+1,j,k) - y(i,j,k). Представляйте неравенство запуска генератора в терминах w.

w = optimexpr(nPeriods,nGens); % Allocate w
idx = 1:(nPeriods-1);
w(idx,:) = y(idx+1,:,'Low') - y(idx,:,'Low') + y(idx+1,:,'High') - y(idx,:,'High');
w(nPeriods,:) = y(1,:,'Low') - y(nPeriods,:,'Low') + y(1,:,'High') - y(nPeriods,:,'High');
switchcons = w - z <= 0;

Задайте цель

Целевая функция включает затраты на топливо для выполнения генераторов, дохода от выполнения генераторов, и стоит за запуск генераторов.

generatorlevel  = zeros(size(yFuel));
generatorlevel(:,1,1) = gen(1,1); % Fill in the levels
generatorlevel(:,1,2) = gen(1,2);
generatorlevel(:,2,1) = gen(2,1);
generatorlevel(:,2,2) = gen(2,2); 

Входящий доход = y.*generatorlevel.*poolPrice.

revenue = optimexpr(size(y));
for ii = 1:nPeriods
    revenue(ii,:,:) = poolPrice(ii)*y(ii,:,:).*generatorlevel(ii,:,:);
end

Общая стоимость топлива = fuelUsed*fuelPrice.

fuelCost = fuelUsed*fuelPrice;

Стоимость запуска генератора = z*startCost.

startingCost = z*startCost;

Прибыль = входящий доход - общая стоимость топлива - запуск стоится.

profit = sum(sum(sum(revenue))) - fuelCost - sum(sum(startingCost));

Решите задачу

Создайте задачу оптимизации и включайте цель и ограничения.

dispatch = optimproblem('ObjectiveSense','maximize');
dispatch.Objective = profit;
dispatch.Constraints.switchcons = switchcons;
dispatch.Constraints.fuelcons = fuelcons;
dispatch.Constraints.powercons = powercons;

Чтобы оставить свободное место, подавите итеративное отображение.

options = optimoptions('intlinprog','Display','final');

Решите задачу.

[dispatchsol,fval,exitflag,output] = solve(dispatch,'options',options);

Solving problem using intlinprog.

Optimal solution found.

Intlinprog stopped because the objective value is within a gap tolerance of the
optimal value, options.AbsoluteGapTolerance = 0 (the default value). The intcon
variables are integer within tolerance, options.IntegerTolerance = 1e-05 (the
default value).

Исследуйте решение

Постройте решение в зависимости от времени.

subplot(3,1,1)
bar(dispatchsol.y(:,1,1)*gen(1,1)+dispatchsol.y(:,1,2)*gen(1,2),.5,'g')
xlim([.5,48.5])
ylabel('MWh')
title('Generator 1 Optimal Schedule','FontWeight','bold')
subplot(3,1,2)
bar(dispatchsol.y(:,2,1)*gen(1,1)+dispatchsol.y(:,2,2)*gen(1,2),.5,'c')
title('Generator 2 Optimal Schedule','FontWeight','bold')
xlim([.5,48.5])
ylabel('MWh')
subplot(3,1,3)
bar(poolPrice,.5)
xlim([.5,48.5])
title('Energy Price','FontWeight','bold')
xlabel('Period')
ylabel('$ / MWh')

Генератор 2 запуска дольше, чем генератор 1, который вы ожидали бы, потому что это более эффективно. Генератор 2 запуска на его мощном уровне каждый раз, когда это включено. Генератор 1 запуск в основном на его мощном уровне, но падениях вниз к малой мощности для одной единицы измерения времени. Каждый генератор запуски для одного непрерывного набора периодов ежедневно, и, поэтому, несет только одни расходы запуска каждый день.

Проверяйте что z переменная 1 в течение периодов, когда генераторы запускаются.

starttimes = find(round(dispatchsol.z) == 1); % Use round for noninteger results
[theperiod,thegenerator] = ind2sub(size(dispatchsol.z),starttimes)

theperiod = 2×1

    23
    16

thegenerator = 2×1

     1
     2

Периоды, когда генераторы начинают матч графики.

Сравните с более низким штрафом за запуск

Если вы задаете нижнее значение для startCost, решение включает несколько периодов генерации.

startCost = 500; % Choose a lower penalty for starting the generators
startingCost = z*startCost;
profit = sum(sum(sum(revenue))) - fuelCost - sum(sum(startingCost));
dispatch.Objective = profit;
[dispatchsolnew,fvalnew,exitflagnew,outputnew] = solve(dispatch,'options',options);

Solving problem using intlinprog.

Optimal solution found.

Intlinprog stopped because the objective value is within a gap tolerance of the
optimal value, options.AbsoluteGapTolerance = 0 (the default value). The intcon
variables are integer within tolerance, options.IntegerTolerance = 1e-05 (the
default value).

subplot(3,1,1)
bar(dispatchsolnew.y(:,1,1)*gen(1,1)+dispatchsolnew.y(:,1,2)*gen(1,2),.5,'g')
xlim([.5,48.5])
ylabel('MWh')
title('Generator 1 Optimal Schedule','FontWeight','bold')
subplot(3,1,2)
bar(dispatchsolnew.y(:,2,1)*gen(1,1)+dispatchsolnew.y(:,2,2)*gen(1,2),.5,'c')
title('Generator 2 Optimal Schedule','FontWeight','bold')
xlim([.5,48.5])
ylabel('MWh')
subplot(3,1,3)
bar(poolPrice,.5)
xlim([.5,48.5])
title('Energy Price','FontWeight','bold')
xlabel('Period')
ylabel('$ / MWh')

starttimes = find(round(dispatchsolnew.z) == 1); % Use round for noninteger results
[theperiod,thegenerator] = ind2sub(size(dispatchsolnew.z),starttimes)

theperiod = 3×1

    22
    16
    45

thegenerator = 3×1

     1
     2
     2