Exponenta Banner
exponenta event banner

Оптимизация конструкции с неопределенными переменными (код)

Открыть сценарий

В этом примере показано, как оптимизировать конструкцию при наличии неопределенных переменных. Оптимизируйте размеры реактора с непрерывно перемешиваемым резервуаром (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}$ - Концентрации А в 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');

Модель включает каскадный PID-контроллер в 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];

Создание функции цели/ограничения

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

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

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

evalDesign функция:

  • Имеет один входной аргумент, указывающий размеры CSTR

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

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);

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] Becette, B.W. Process Dynamics: моделирование, анализ и моделирование. 1-я ред. река Верхнее Седло, Нью-Джерси: Прентис Холл, 1998.

% Close the model
bdclose('sdoCSTR')

Документация по оптимизации проектирования Simulink

Поддержка