Разработайте исследование Используя выборку параметра (код)

Этот пример показывает, как выбрать и исследовать пробел проекта. Вы исследуете проект Постоянно Реактор Смесителя, чтобы минимизировать изменение концентрации продукта и стоимость изготовления. Проект включает неуверенность запаса канала.

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

Можно также использовать выбранное пространство проекта и оценку Монте-Карло вывод, чтобы анализировать влияние параметров проекта на проекте, видеть "Анализ чувствительности для Оценки Параметра (Код)"

Модель Постоянно реактора смесителя (CSTR)

Постоянно Реакторы Смесителя (CSTRs) распространены в перерабатывающей промышленности. Модель Simulink, sdoCSTR, моделирует покрытое кожухом диабатическое (т.е. неадиабатическая) реактор корпуса, описанный в [1]. CSTR принят, чтобы быть отлично смешанным с одной экзотермической и необратимой реакцией первого порядка., реагент, преобразован в, продукт.

В этом примере вы используете следующую модель CSTR с двумя состояниями, которая использует основные принципы учета и энергосбережения:

  • , и - Концентрации в CSTR и в канале [kgmol/m^3]

  • , и - CSTR, канал и температуры хладагента [K]

  • и - Объемный расход [m^3/h] и плотность материала в CSTR [1/m^3]

  • и - Высота [m] и нагретая площадь поперечного сечения [m^2] CSTR.

  • - Предэкспоненциальный нетепловой фактор для [1/h] реакции

  • и - Энергия активации и тепло реакции для [kcal/kgmol]

  • - Газовая константа Больцманна [kcal / (kgmol * K)]

  • и - Теплоемкость [kcal/K] и коэффициенты теплопередачи [kcal / (m^2 * K * h)]

Откройте модель Simulink.

open_system('sdoCSTR');

Проблема проектирования CSTR

Примите, что CSTR является цилиндрическим с хладагентом, к которому применяются основа цилиндра. Настройте площадь поперечного сечения CSTR и высоту CSTR, чтобы удовлетворить следующим целям проекта:

  • Минимизируйте изменение в остаточной концентрации. Изменения в остаточной концентрации негативно влияют на качество продукта CSTR. Минимизация изменений также улучшает прибыль CSTR.

  • Минимизируйте среднюю температуру хладагента. Нагревание или охлаждение температуры хладагента конверта являются дорогими. Минимизация средней температуры хладагента улучшает прибыль CSTR.

Проект должен допускать изменения в качестве концентрации канала предоставления и питать температуру. CSTR питается каналом от различных поставщиков. Качество канала отличается от поставщика поставщику и также отличается в каждом пакете предоставления.

Задайте переменные проекта

Выберите следующие параметры модели как переменные проекта:

  • Цилиндрическая площадь поперечного сечения

  • Цилиндрическая высота

p = sdo.getParameterFromModel('sdoCSTR',{'A','h'});

Ограничьте площадь поперечного сечения областью значений [0.2 2] m^2.

p(1).Minimum = 0.2;
p(1).Maximum = 2;

Ограничьте высоту областью значений [0.5 3] m.

p(2).Minimum = 0.5;
p(2).Maximum = 3;

Выберите пробел проекта

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

pSpace = sdo.ParameterSpace(p);

Пространство параметров использует равномерные распределения по умолчанию для переменных проекта. Распределение нижние и верхние границы установлено в переменную проекта минимальное и максимальное значение соответственно.

Используйте функцию sdo.sample, чтобы сгенерировать выборки от пространства параметров. Вы используете выборки, чтобы оценить модель и исследовать пробел проекта.

rng default; % For reproducibility
pSmpl = sdo.sample(pSpace,30);

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

figure, sdo.scatterPlot(pSmpl)

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

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

pUnc = sdo.getParameterFromModel('sdoCSTR',{'FeedCon0','FeedTemp0'});

Создайте пространство параметров для неопределенных переменных. Используйте нормальные распределения для обеих переменных. Задайте среднее значение как текущее значение параметров. Задайте отклонение 5% среднего значения для концентрации канала и 1% среднего значения для температуры.

uSpace = sdo.ParameterSpace(pUnc);
uSpace = setDistribution(uSpace,'FeedCon0',makedist('normal',pUnc(1).Value,0.05*pUnc(1).Value));
uSpace = setDistribution(uSpace,'FeedTemp0',makedist('normal',pUnc(2).Value,0.01*pUnc(2).Value));

Концентрация канала обратно пропорционально коррелируется с температурой канала. Добавьте эту информацию в пространство параметров.

uSpace.RankCorrelation = [1 -0.6; -0.6 1];

Матрица порядковой корреляции имеет строку и столбец для каждого параметра с (i, j) запись, задающая корреляцию между i и j параметрами.

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

uSmpl = sdo.sample(uSpace,60);
sdo.scatterPlot(uSmpl)

Идеально вы хотите оценить проект для каждой комбинации точек в проекте и неопределенных пробелах, который подразумевает 30*60 = 1 800 симуляций. Каждая симуляция занимает приблизительно 0,5 секунды. Можно использовать параллельные вычисления, чтобы ускорить оценку. Для этого примера вы вместо этого только используете выборки, которые имеют максимальную & минимальную концентрацию и температурные значения, уменьшая время оценки приблизительно до 1 min.

[~,iminC] = min(uSmpl.FeedCon0);
[~,imaxC] = max(uSmpl.FeedCon0);
[~,iminT] = min(uSmpl.FeedTemp0);
[~,imaxT] = max(uSmpl.FeedTemp0);
uSmpl = uSmpl(unique([iminC,imaxC,iminT,imaxT]) ,:)
uSmpl =

  4x2 table

    FeedCon0    FeedTemp0
    ________    _________

     9.4555      303.58  
     11.175      288.13  
     11.293      290.73  
     8.9308      294.16  

Создайте функцию оценки

Создайте функцию, которая оценивает проект для данной точки выборки на пробеле проекта. Проект оценен о том, как хорошо он минимизирует изменение в остаточной концентрации и средней температуре хладагента.

Задайте конструктивные требования

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

Регистрируйте следующие сигналы:

  • Концентрация CSTR, доступная во втором выходном порту блока sdoCSTR/CSTR

Conc = Simulink.SimulationData.SignalLoggingInfo;
Conc.BlockPath               = 'sdoCSTR/CSTR';
Conc.OutputPortIndex         = 2;
Conc.LoggingInfo.NameMode    = 1;
Conc.LoggingInfo.LoggingName = 'Concentration';
  • Температура хладагента, доступная при первом выводе блока sdoCSTR/Controller

Coolant = Simulink.SimulationData.SignalLoggingInfo;
Coolant.BlockPath               = 'sdoCSTR/Controller';
Coolant.OutputPortIndex         = 1;
Coolant.LoggingInfo.NameMode    = 1;
Coolant.LoggingInfo.LoggingName = 'Coolant';

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

simulator = sdo.SimulationTest('sdoCSTR');
simulator.LoggingInfo.Signals = [Conc,Coolant];

Функция оценки

Используйте анонимную функцию с одним аргументом, который вызывает функцию sdoCSTR_design.

evalDesign = @(p) sdoCSTR_design(p,simulator,pUnc,uSmpl);

Функция evalDesign:

  • Имеет один входной параметр, который задает размерности CSTR

  • Возвращает значение цели оптимизации

Функция sdoCSTR_design использует цикл for, который выполняет итерации через демонстрационные значения, заданные для концентрации канала и температуры. В цикле, функции:

  • Моделирует модель с помощью текущей точки проекта, концентрации канала и значений температуры канала

  • Вычисляет остаточное изменение концентрации и затраты температуры хладагента

Чтобы просмотреть целевую функцию, введите edit sdoCSTR_design.

type sdoCSTR_design
function design = sdoCSTR_design(p,simulator,pUnc,smplUnc)
%SDOCSTR_DESIGN
%
% The sdoCSTR_design function is used to evaluate a CSTR  design.
%
% The |p| input argument is the vector of CSTR dimensions.
%
% The |simulator| input argument is a sdo.SimulinkTest object used to
% simulate the |sdoCSTR| model and log simulation signals.
%
% The |pUnc| input argument is a vector of parameters to specify the CSTR
% input feed concentration and feed temperature. The |smplUnc| argument is
% a table of different feed concentration and temperature values.
%
% The |design| return argument contains information about the design
% evaluation that can be used by the |sdo.optimize| function to optimize
% the design.
%
% see also sdo.optimize, sdoExampleCostFunction
%

% Copyright 2012-2013 The MathWorks, Inc.

%% Model Simulations and Evaluations
%
% For each value in |smplUnc|, configure and simulate the model. Use
% the logged concentration and coolant signals to compute the design cost.
%
costConc    = 0;
costCoolant = 0;
for ct=1:size(smplUnc,1)
    %Set the feed concentration and temperature values
    pUnc(1).Value = smplUnc{ct,1};
    pUnc(2).Value = smplUnc{ct,2};
    
    %Simulate model
    simulator.Parameters = [p; pUnc];
    simulator = sim(simulator);
    logName   = get_param('sdoCSTR','SignalLoggingName');
    simLog    = get(simulator.LoggedData,logName);
    
    %Compute Concentration cost based on the standard deviation of the
    %concentration from a nominal value.
    Sig = find(simLog,'Concentration');
    costConc = costConc+10*std(Sig.Values-2);
    
    %Compute coolant cost based on the mean deviation from room
    %temperature.
    Sig = find(simLog,'Coolant');
    costCoolant = costCoolant+abs(mean(Sig.Values - 294))/30;
end

%% Return Total Cost 
%
% Compute the total cost as a sum of the concentration and coolant costs.
%
design.F = costConc + costCoolant;

%%
% Add the individual cost terms to the return argument. These are not used
% by the optimizer, but included for convenience.
design.costConc    = costConc;               
design.costCoolant = costCoolant;
end

Оценить

Используйте команду sdo.evaluate, чтобы оценить модель в демонстрационных точках проекта.

y = sdo.evaluate(evalDesign,p,pSmpl);
Model evaluated at 30 samples.

Просмотрите результаты оценки с помощью графика рассеивания. График рассеивания показывает попарные графики для каждой переменной проекта (A, h) и стоимость проекта. График включает общую стоимость, F, а также хладагент и затраты на концентрацию, costCoolant и costConc соответственно.

sdo.scatterPlot(pSmpl,y);

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

Создайте сетчатый график, показывающий общую стоимость функцией A и h.

Ftotal = scatteredInterpolant(pSmpl.A,pSmpl.h,y.F);
xR = linspace(min(pSmpl.A),max(pSmpl.A),60);
yR = linspace(min(pSmpl.h),max(pSmpl.h),60);
[xx,yy] = meshgrid(xR,yR);
zz = Ftotal(xx,yy);
mesh(xx,yy,zz)
view(56,30)
title('Total cost as function of A and h')
zlabel('Ftotal')
xlabel(p(1).Name), ylabel(p(2).Name);

График показывает, что высокие значения A и h приводят к более низким ценам. Лучший проект на выбранном пробеле соответствует проекту с самым дешевым значением.

[~,idx] = min(y.F);
pBest = [y(idx,:), pSmpl(idx,:)]
pBest =

  1x5 table

      F      costConc    costCoolant      A         h   
    _____    ________    ___________    ______    ______

    2.106     1.5505       0.55552      1.9271    2.3867

Совершенствуйте пробел проекта

Объемная поверхностная диаграмма общей стоимости показывает, что недорогие проекты являются проектами с в области значений [1.5 2] и h в области значений [2 3]. Измените дистрибутивы пространства параметров для A и h и передискретизируйте пробел проекта, чтобы фокусироваться на этой области.

pSpace = setDistribution(pSpace,'A',makedist('uniform',1.5,2));
pSpace = setDistribution(pSpace,'h',makedist('uniform',2,3));
pSmpl = sdo.sample(pSpace,30);

Добавьте pBest, найденный ранее к новым выборкам так, чтобы это была часть усовершенствованного пробела проекта.

pSmpl = [pSmpl;pBest(:,4:5)];
sdo.scatterPlot(pSmpl)

Оцените использующее Усовершенствованное Пространство Проекта

y = sdo.evaluate(evalDesign,p,pSmpl);
Model evaluated at 31 samples.

Создайте сетчатый график для этого раздела пробела проекта. Поверхность указывает, что лучшие проекты около = 1.9, h = 2,1 точки.

Ftotal = scatteredInterpolant(pSmpl.A,pSmpl.h,y.F);
xR = linspace(min(pSmpl.A),max(pSmpl.A),60);
yR = linspace(min(pSmpl.h),max(pSmpl.h),60);
[xx,yy] = meshgrid(xR,yR);
zz = Ftotal(xx,yy);
mesh(xx,yy,zz)
view(56,30)
title('Total cost as function of A and h')
zlabel('Ftotal')
xlabel(p(1).Name), ylabel(p(2).Name);

Найдите лучший проект от нового пробела проекта и сравните с лучшей точкой проекта, найденной ранее.

[~,idx] = min(y.F);
pBest = [pBest; [y(idx,:), pSmpl(idx,:)]]
pBest =

  2x5 table

      F       costConc    costCoolant      A         h   
    ______    ________    ___________    ______    ______

     2.106     1.5505       0.55552      1.9271    2.3867
    1.9754     1.4824       0.49295      1.9695    2.1174

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

Связанные примеры

Чтобы изучить, как исследовать пробел проекта CSTR с помощью Sensitivity Analysis Tool, см., "Что Исследование проекта использует Выборку Параметра (графический интерфейс пользователя)".

Чтобы изучить, как оптимизировать проект CSTR с помощью команды sdo.optimize, см. "Оптимизацию проекта с Неопределенными Переменными (Код)".

Чтобы изучить, как анализировать влияние параметров проекта на проекте с помощью команды sdo.analyze, см. "Анализ чувствительности для Оценки Параметра (Код)"

Ссылки

[1] Bequette, Б.В. Просесс Динэмикс: Моделирование, Анализ и Симуляция. 1-й редактор Верхний Сэддл-Ривер, NJ: Prentice Hall, 1998.

Закройте модель

bdclose('sdoCSTR')