Проектируйте оптимизацию с неопределенными переменными (код)

Этот пример показов, как оптимизировать проект, когда существуют неопределенные переменные. Вы оптимизируете размерности постоянно перемешиваемого Бака реактора (CSTR), чтобы минимизировать изменение концентрации продукта и производственные затраты в случае изменения, или неопределенности, сырья.

Бак с непрерывным перемешиванием (CSTR)

Баки с непрерывным перемешиванием (CSTR) распространены в технологической промышленности. Модель 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}$ - Концентрации A в CSTR и в сырье [кгмоль/м ^ 3]

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

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

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

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

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

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

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

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

open_system('sdoCSTR');

Модель включает каскадный ПИД-регулятор в Controller подсистема. Контроллер регулирует температуру реактора $T$и остаточную концентрацию реактора.$C_A$

Проблема проекта 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'});

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

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

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

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

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

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

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.

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

rng('default'); %For reproducibility
uSmpl = sdo.sample(uSpace,60);
sdo.scatterPlot(uSmpl)

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

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

Задайте требования проект

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

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

  • Концентрация 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];

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

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

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

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

The evalDesign функция:

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

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

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

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

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

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

Оценка первоначального проекта

Вызовите evalDesign функция с начальными размерностями CSTR.

dInit = evalDesign(p)
dInit = 

  struct with fields:

              F: 11.3356
       costConc: 6.4390
    costCoolant: 4.8965

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

sdoCSTR_plotModelResponse(p,simulator,pUnc,uSmpl);

The sdoCSTR_plotModelResponse графики функций модели отклика. Чтобы просмотреть эту функцию, введите edit sdoCSTR_plotModelResponse.

Оптимизируйте проект

Передайте целевую функцию и начальные размерности CSTR, чтобы sdo.optimize.

pOpt = sdo.optimize(evalDesign,p)
 Optimization started 27-Jan-2021 14:55:02

                               max                     First-order 
 Iter F-count        f(x)   constraint    Step-size    optimality
    0      4      5.17535            0
    1      8      3.81522            0         2.01         7.83
    2     12      2.63963            0         0.57         3.03
    3     16      2.53426            0        0.159        0.292
    4     20      2.50995            0        0.017        0.432
    5     24      2.50266            0        0.125         1.44
    6     28      2.50097            0       0.0966        0.234
    7     32      2.47473            0       0.0135       0.0719
    8     37      2.47473            0     0.000955       0.0694
Local minimum possible. Constraints satisfied.

fmincon stopped because the size of the current step is less than
the value of the step size tolerance and constraints are 
satisfied to within the value of the constraint tolerance.
 
pOpt(1,1) =
 
       Name: 'A'
      Value: 2
    Minimum: 1
    Maximum: 2
       Free: 1
      Scale: 0.5000
       Info: [1x1 struct]

 
pOpt(2,1) =
 
       Name: 'h'
      Value: 2.2850
    Minimum: 1
    Maximum: 3
       Free: 1
      Scale: 2
       Info: [1x1 struct]

 
2x1 param.Continuous
 

Оценка оптимизированного проекта

Вызовите evalDesign функция с оптимизированными размерностями CSTR.

dFinal = evalDesign(pOpt)
dFinal = 

  struct with fields:

              F: 2.4747
       costConc: 1.4505
    costCoolant: 1.0242

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

sdoCSTR_plotModelResponse(pOpt,simulator,pUnc,uSmpl);

Похожие примеры

Чтобы узнать, как использовать анализ чувствительности для исследования проектного пространства CSTR и выбора начального проекта для оптимизации, см. Исследование проекта с использованием выборки параметров (код).

Ссылки

[1] Бекет, B.W. Динамика процесса: моделирование, анализ и симуляция. 1-я ред. Верхняя Седловая река, Нью-Джерси: Prentice Hall, 1998.

% Close the model
bdclose('sdoCSTR')