Exponenta Banner

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

Открыть скрипт

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

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

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

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

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

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

$$\frac{d C_A}{dt} = \frac{F}{A * h}(C_{feed} - C_A)-r*C_A$$

$$\frac{d T}{dt} = \frac{F}{A * h}(T_{feed} - T) - \frac{H}{c_p \rho} r
-\frac{U}{c_p*\rho*h}(T - T_{cool})$$

$$r = k_0*e^{\frac{-E}{R*T}}$$

  • $C_A$, и$C_{feed}$ - Концентрации в CSTR и в канале [kgmol/m^3]

  • $T$$T_{feed}$, и$T_{cool}$ - CSTR, канал и температуры хладагента [K]

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

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

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

  • $E$ и$H$ - Энергия активации и тепло реакции для$A\rightarrow B$ [kcal/kgmol]

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

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

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

open_system('sdoCSTR');

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

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

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

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

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

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

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

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

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

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

Ограничьте площадь поперечного сечения областью значений [0.2 2] м^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')

Похожие темы

Документация Simulink Design Optimization
Поддержка