Планирование систем Распределения водой с использованием обучения с подкреплением

В этом примере показано, как узнать оптимальную политику планирования насоса для системы распределения воды с помощью обучения с подкреплением (RL).

Система распределения

На следующем рисунке показана система распределения.

Здесь:

  • QSupply количество воды, подаваемой в бак с водой из резервуара.

  • QDemand количество воды, вытекающей из бака для удовлетворения потребности в использовании.

Цель агента обучения с подкреплением состоит в том, чтобы запланировать количество работающих насосов, чтобы минимизировать энергопотребление системы и удовлетворить потребность в использовании (h>0). Диафрагма системы бака определяется следующим уравнением.

Adhdt=QSupply(t)-QDemand(t)

Вот, A=40m2 и hmax=7m. Потребность в течение 24-часового периода является функцией времени, заданного как

QDemand(t)=μ(t)+η(t)

где μ(t)является ожидаемым спросом и η(t) представляет неопределенность спроса, которая выбирается из равномерного случайного распределения.

Подача определяется количеством работающих насосов, a{0,1,2,3}согласно следующему отображению.

QSupply(t)=Q(a)={0a=0164a=1279a=2344a=3  cmh

Чтобы упростить задачу, потребление степени определяется как количество работающих насосов, a.

Следующая функция является вознаграждением для этого окружения. Чтобы избежать переполнения или опорожнения бака, добавляются дополнительные затраты, если высота воды близка к максимальному или минимальному уровню воды. hmaxили hmin, соответственно.

r(h,a)=-10(h(hmax-0.1))-10(h0.1)-a

Сгенерируйте профиль потребности

Чтобы сгенерировать и профиль потребности в воде, основанный на количестве рассматриваемых дней, используйте generateWaterDemand функция, заданная в конце этого примера.

num_days = 4; % Number of days
[WaterDemand,T_max] = generateWaterDemand(num_days);

Просмотрите профиль потребности.

plot(WaterDemand)

Figure contains an axes. The axes with title Time Series Plot:Water Demand contains an object of type line.

Откройте и сконфигурируйте модель

Откройте распределение систему Simulink модели.

mdl = "watertankscheduling";
open_system(mdl)

В дополнение к агенту обучения с подкреплением в блоке Control law MATLAB Function задан простой контроллер базовой линии. Этот контроллер активирует определенное количество насосов в зависимости от уровня воды.

Задайте начальный уровень воды.

h0 = 3; % m

Задайте параметры модели.

SampleTime = 0.2;
H_max = 7; % Max tank height (m)
A_tank = 40; % Area of tank (m^2)

Создайте интерфейс окружения для агента RL

Чтобы создать интерфейс окружения для модели Simulink, сначала задайте спецификации действий и наблюдений, actInfo и obsInfo, соответственно. Действие агента является выбранным количеством насосов. Наблюдение агента - этот уровень воды, которая эмазируется как сигнал непрерывного времени.

actInfo = rlFiniteSetSpec([0,1,2,3]);
obsInfo = rlNumericSpec([1,1]);

Создайте интерфейс окружения.

env = rlSimulinkEnv(mdl,mdl+"/RL Agent",obsInfo,actInfo);

Задайте пользовательскую функцию сброса, которая задана в конце этого примера, которая рандомизирует начальный уровень воды и потребность в воде. Это позволяет обучать агента на различных начальных уровнях воды и функциях спроса на воду для каждого эпизода.

env.ResetFcn = @(in)localResetFcn(in);

Создание агента DQN

Агент DQN аппроксимирует долгосрочное вознаграждение, заданные наблюдения и действия, используя представление функции Q-значения критика. Чтобы создать критика, сначала создайте глубокую нейронную сеть. Для получения дополнительной информации о создании представления функции ценности глубокой нейронной сети, смотрите, Создают политику и представления функции ценности.

% Fix the random generator seed for reproducibility.
rng(0);

Создайте глубокую нейронную сеть для критика. В данном примере используйте непериодическую нейронную сеть. Чтобы использовать рецидивирующую нейронную сеть, задайте useLSTM на true.

useLSTM = false;
if useLSTM
    layers = [
        sequenceInputLayer(obsInfo.Dimension(1),"Name","state","Normalization","none")
        fullyConnectedLayer(32,"Name","fc_1")
        reluLayer("Name","relu_body1")
        lstmLayer(32,"Name","lstm")
        fullyConnectedLayer(32,"Name","fc_3")
        reluLayer("Name","relu_body3")
        fullyConnectedLayer(numel(actInfo.Elements),"Name","output")];
else
    layers = [
        featureInputLayer(obsInfo.Dimension(1),"Name","state","Normalization","none")
        fullyConnectedLayer(32,"Name","fc_1")
        reluLayer("Name","relu_body1")
        fullyConnectedLayer(32,"Name","fc_2")
        reluLayer("Name","relu_body2")
        fullyConnectedLayer(32,"Name","fc_3")
        reluLayer("Name","relu_body3")
        fullyConnectedLayer(numel(actInfo.Elements),"Name","output")];
end

Задайте опции для создания представления критика.

criticOpts = rlRepresentationOptions('LearnRate',0.001,'GradientThreshold',1);

Создайте представление критика с помощью заданных глубокой нейронной сети и опций.

critic = rlQValueRepresentation(layerGraph(layers),obsInfo,actInfo,...
    'Observation',{'state'},criticOpts);

Создание агента DQN

Чтобы создать агента, сначала задайте опции агента. Если вы используете сеть LSTM, установите длину последовательности равной 20.

opt = rlDQNAgentOptions('SampleTime',SampleTime); 
if useLSTM
    opt.SequenceLength = 20;
else
    opt.SequenceLength = 1;
end
opt.DiscountFactor = 0.995;
opt.ExperienceBufferLength = 1e6;
opt.EpsilonGreedyExploration.EpsilonDecay = 1e-5;
opt.EpsilonGreedyExploration.EpsilonMin = .02;

Создайте агента с помощью заданных опций и представления критика.

agent = rlDQNAgent(critic,opt);    

Обучите агента

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

  • Запустите обучение для 1000 эпизодов с каждым эпизодом, длящимся в ceil(T_max/Ts) временные шаги.

  • Отображение процесса обучения в диалоговом окне Диспетчер эпизодов (установите Plots опция)

Задайте опции для обучения с использованием rlTrainingOptions объект.

trainOpts = rlTrainingOptions(...
    'MaxEpisodes',1000, ...
    'MaxStepsPerEpisode',ceil(T_max/SampleTime), ...
    'Verbose',false, ...
    'Plots','training-progress',...
    'StopTrainingCriteria','EpisodeCount',...
    'StopTrainingValue',1000,...
    'ScoreAveragingWindowLength',100);

В то время как вы не делаете это для этого примера, можно сохранить агентов в процессе обучения. Например, следующие опции сохраняют каждый агент со значением вознаграждения, большим или равным -42.

Save agents using SaveAgentCriteria if necessary
trainOpts.SaveAgentCriteria = 'EpisodeReward';
trainOpts.SaveAgentValue = -42;

Обучите агента с помощью train функция. Обучение этого агента является интенсивным в вычислительном отношении процессом, который занимает несколько часов. Чтобы сэкономить время при запуске этого примера, загрузите предварительно обученного агента путем установки doTraining на false. Чтобы обучить агента самостоятельно, установите doTraining на true.

doTraining = false;
if doTraining
    % Train the agent.
    trainingStats = train(agent,env,trainOpts);
else
    % Load the pretrained agent for the example.
    load('SimulinkWaterDistributionDQN.mat','agent')
end

Следующий рисунок показывает процесс обучения.

Симулируйте агент DQN

Чтобы подтвердить производительность обученного агента, моделируйте его в среде бака с водой. Для получения дополнительной информации о симуляции агента смотрите rlSimulationOptions и sim.

Чтобы симулировать эффективность агента, соедините блок RL Agent путем переключения ручного блока switch.

set_param(mdl+"/Manual Switch",'sw','0');

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

NumSimulations = 30;
simOptions = rlSimulationOptions('MaxSteps',T_max/SampleTime,...
    'NumSimulations', NumSimulations);

Чтобы сравнить агент с контроллером базовой линии при тех же условиях, сбросьте входный случайный seed, используемое в функции сброса окружения.

env.ResetFcn("Reset seed");

Симулируйте агент против окружения.

experienceDQN = sim(env,agent,simOptions);

Симулируйте базовый контроллер

Чтобы сравнить агента DQN с базовым контроллером, вы должны симулировать базовый контроллер с помощью тех же опций симуляции и начального случайного seed для функции сброса.

Включите базовый контроллер.

set_param(mdl+"/Manual Switch",'sw','1');

Чтобы сравнить агент с контроллером базовой линии при тех же условиях, сбросьте случайный seed, используемое в функции сброса окружения.

env.ResetFcn("Reset seed");

Симулируйте контроллер базовой линии относительно окружения.

experienceBaseline = sim(env,agent,simOptions);

Сравнение агента DQN с базовым контроллером

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

resultVectorDQN = zeros(NumSimulations, 1);
resultVectorBaseline = zeros(NumSimulations,1);

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

for ct = 1:NumSimulations
    resultVectorDQN(ct) = sum(experienceDQN(ct).Reward);
    resultVectorBaseline(ct) = sum(experienceBaseline(ct).Reward);
end

Постройте график совокупных вознаграждений.

plot([resultVectorDQN resultVectorBaseline],'o')
set(gca,'xtick',1:NumSimulations)
xlabel("Simulation number")
ylabel('Cumulative Reward')
legend('DQN','Baseline','Location','NorthEastOutside')

Figure contains an axes. The axes contains 2 objects of type line. These objects represent DQN, Baseline.

Совокупное вознаграждение, получаемое агентом, постоянно составляет около -40. Это значение намного больше, чем среднее вознаграждение, полученное базовым контроллером. Поэтому агент DQN последовательно превосходит контроллер базового уровня с точки зрения экономии энергии.

Локальные функции

Функция потребности в воде

function [WaterDemand,T_max] = generateWaterDemand(num_days)

    t = 0:(num_days*24)-1; % hr
    T_max = t(end);

    Demand_mean = [28, 28, 28, 45, 55, 110, 280, 450, 310, 170, 160, 145, 130, ...
        150, 165, 155, 170, 265, 360, 240, 120, 83, 45, 28]'; % m^3/hr

    Demand = repmat(Demand_mean,1,num_days);
    Demand = Demand(:);

    % Add noise to demand
    a = -25; % m^3/hr
    b = 25; % m^3/hr
    Demand_noise = a + (b-a).*rand(numel(Demand),1);

    WaterDemand = timeseries(Demand + Demand_noise,t);
    WaterDemand.Name = "Water Demand";
end

Функция сброса

function in = localResetFcn(in)

    % Use a persistent random seed value to evaluate the agent and the baseline
    % controller under the same conditions.
    persistent randomSeed
    if isempty(randomSeed)
        randomSeed = 0;
    end
    if strcmp(in,"Reset seed")
        randomSeed = 0;
        return
    end    
    randomSeed = randomSeed + 1;
    rng(randomSeed)
    
    % Randomize water demand.
    num_days = 4;
    H_max = 7;
    [WaterDemand,~] = generateWaterDemand(num_days);
    assignin('base','WaterDemand',WaterDemand)

    % Randomize initial height.
    h0 = 3*randn;
    while h0 <= 0 || h0 >= H_max
        h0 = 3*randn;
    end
    blk = 'watertankscheduling/Water Tank System/Initial Water Height';

    in = setBlockParameter(in,blk,'Value',num2str(h0));

end