Например, Бобцов

Разработка математической модели пекарной камеры как объекта с сосредоточенными параметрами

УДК-681.513.5
Разработка математической модели пекарной камеры как объекта с сосредоточенными параметрами
В.Б. Данин, А.Ю. Кириков
Санкт-Петербургский государственный университет низкотемпературных и пищевых технологий
В данной работе рассмотрены вопросы построения математической модели пекарной камеры тоннельной типа. Авторами предлагается математическая модель пекарной камеры, представляющей собой объект с распределенными параметрами как объект со сосредоточенными параметрами. На базе разработанной эквивалентной модели с сосредоточенными параметрами можно сформулировать и дать способ решения задачи оптимальной стабилизации температурного режима пекарной камеры.
Ключевые слова: пекарная камера, система уравнений, математическая модель.
В качестве объекта взята тоннельная печь для выпечки хлеба в объединении ОАО “Каравай”.
Температурное поле тоннельной печи условно разделено на 3 последовательных температурных участка для выпечки движущихся хлебобулочных заготовок с температурами: 120 – 140°C; 270 - 290°C; 180 - 220°C соответственно. Такое деление на участки определяет и конфигурацию математической модели пекарной камеры, которая рассматривается в виде моделей отдельных участков. Схема математической модели показана на рис.1, которая состоит из их последовательного соединения с указанием тепловых потоков между ними и управлений, соответствующих изменению подачи топлива.
1

Рис. 1. Схема модели температурных участков. GП,GГ,GВ – расход продукта, газа и воздуха для каждого из участков; ТП,ТД – температура продукта и дымовых газов для каждого из участков; QК,QР – коли-чество тепла , излучаемого
кладкой и роликами для каждого из участков.
На основании геометрических размеров пекарной камеры (печи) и материалов кладки (футеровки стенок), собственно хлебной заготовки и конвейерного транспорта (роликов) необходимо определить численные значения коэффициентов математической модели для последующих расчетов. В этой связи, целесообразно показать расчет коэффициентов излучения для участков, полагая, что в режиме стабилизации основной поток теплообмена осуществляется излучением. При расчете коэффициентов излучения, важно заранее знать от каких поверхностей передается тепло излучением, а какие поверхности его поглощают. Направленность тепловых потоков будет зависеть от абсолютной температуры поверхностей и объемов. Наибольшей температурой будет обладать продукты сгоревшего топлива (назовем их дымовыми газами) и от них тепловой поток будет направлен к выпекаемому продукту, стенками в виде кладки (или футеровки) и к роликовому транспорту.
При определении математической модели процесса выпекания сделаны следующие предположения:
Температура излучающих поверхностей принимается средней на каждом участке;
Коэффициенты излучения и теплопередачи являются независимыми от температуры в пределах каждого участка;
Теплоемкость тел и воздушной массы также не зависят от температуры, в пределах каждого участка;
Потери в окружающую среду постоянны.
2

Все эти допущения справедливы при стабилизации температурного поля

внутри пекарной камеры, так как отклонения температуры от расчетных режимов

не могут быть большими.

Известно что пекарная камера имеет четыре аккумулятора тепла:

- объем нагретой газовоздушной массы внутри;

- пекарный продукт (полуфабрикат);

- транспортная конвейерная лента (ролики);

- кладка внутренней поверхности камеры.

Поэтому каждый участок в динамике будет описываться системой из четырех

дифференциальных уравнений.

Уравнение теплового баланса для объема нагретой газовоздушной массы

где:

М Г

⋅ GГ

⋅ d ⋅TГ



⋅σТ

⋅d

⋅t

+ σ Г1 ⋅ CГ1 ⋅ TГ1 ⋅ dt

+σВ

⋅ CВ

⋅ TВ

⋅ dt

+

( ) ( )GП1 ⋅ CП1 ⋅ TП1 ⋅ dt − −FК ⋅σ ГК ⋅ TГ 4 − TК 4 ⋅ dt − Fр ⋅σ Гр ⋅ TГ 4 − Tр4 ⋅ dt −

( ) ( )FК ⋅ dК ⋅ TГ − TК ⋅ dt − Fр ⋅ d р TГ 4 − TК 4 ⋅ dt + GП2 ⋅ СП2 ⋅ TП2 ⋅ dt − β ⋅ GП ⋅ dt −

GГ 2

⋅С Г2

⋅ TГ 2

(1)

Mг - масса дымовых газов в объеме участка; Gг - теплоемкость дымовых газов; Tг - температура дымовых газов; γ - удельная теплота сгорания топлива;

Gв - расход воздуха на сгорание топлива; Gт - расход топлива; Gг - расход дымовых газов с предыдущего участка; Tг1- температура дымовых газов из предыдущего участка; Gп - расход выпекаемого продукта; Cп1 - теплоемкость продукта на входе в участок; Cп2 - теплоемкость продукта на выходе из участка; Tп2 -температура продукта на выходе из участка; Tк, Tр - температура поверхности камеры и роликов соответственно;

σгк,

σ гр

-

коэффициенты

излучения;

αк,αр - коэффициенты теплопередачи конвекцией;

β - коэффициент химических реакций;

3

Gт2,Cт2,Tг2 - расход, теплоемкость, температура дымовых газов на выходе из участка соответственно.

Уравнение (1) показывает, что в каждый момент времени нагрев газовой сме-

си (дымовых газов) происходит за счет сгорания топлива, поступления тепла с

дымовыми газами с предыдущего участка, поступления тепла с воздухом, поступ-

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

отдачи тепла роликами и кладкой камеры за счет поглощения тепла на химиче-

ские реакции, уноса тепла с продуктом и с дымовыми газами, уходящими на по-

следующий участок.

В зависимости от характера теплообмена уравнение (1) будет разным для раз-

ных участков. На участке 1 следует учитывать как передачу тепла излучением, так

и конвекцией.

На участках 2 и 3 преобладает передача тепла излучением, а передачей тепла

конвекцией можно пренебречь. Проведя линеаризацию уравнения (1) и сгруппи-

ровав подобные члены выражение примет вид

[ ]( )M

Г

⋅СГ



d∆TГ dt

+ 4 ⋅ TГ 3 FК ⋅ σ ГК + FР ⋅ σ ГР

+ FК ⋅ d К + Fр ⋅ d р + GГ 2 ⋅ СГ 2

⋅ ∆TГ =

( )γ ⋅ ∆GТ + GГ1 ⋅ СГ1 ⋅ ∆TГ1 + GВ ⋅ CВ ⋅ ∆TВ + 4 ⋅ TК 3 ⋅ FК ⋅ GГК + FК ⋅ α К ⋅ ∆TК +

.

( ) [ ]4 ⋅ TР3 ⋅ FР ⋅σ ГР + FР ⋅α Р ⋅ ∆TР − 2GП (TП1 − TП ) − β ⋅ ∆GП + 2 ⋅ GП ⋅ СП ⋅ ∆TП1

Введем обозначение

( )A = 4 ⋅ TГ 3 ⋅ FК ⋅ σ ГК + FРσ ГР + FК ⋅ α К + FР ⋅ α Р + GГ 2 ⋅ СГ 2 .
и приведем уравнение (1) к виду

M Г ⋅СГ A

⋅ d∆TГ dt

+ ∆TГ

=

γ A





+

GГ1



С Г1

A

⋅ ∆TГ1

+ GВ ⋅ CВ A

⋅ ∆TВ

+

4 ⋅ TК

⋅ FК

⋅ σ ГК A

+



⋅αК

⋅ ∆TК

+

4 ⋅ TР

⋅ FР

⋅ σ ГР A

+



⋅αР

⋅ ∆TР

+

2 ⋅ GП

⋅ (TП1 − TП ) −
A

β

+

2 ⋅ GП ⋅ CП A

⋅ ∆TП1

(2)

обозначим

M Ã ⋅CÃ A

= T1

- постоянная времени газового объема,

γ - коэффициент передачи объекта по управлению, за которое примем расA

ход топлива;

4

Введем обозначения:

4 ⋅ TК FК ⋅ GГК A

+ FК ⋅ α К

= a12 ,

2 ⋅ GП ⋅ СП A

=

g2 ,

4 ⋅ Tр Fр ⋅ GГр A

+ Fр ⋅α р

= a13 ,

2 ⋅ GП ⋅ (TП1 − TП ) = g1,
A

2



GГ1



С Г1

A

=

g3 ,

2







С В

A

=

g4 .

С учетом введенных обозначений уравнение (1) примет вид для объема газо-

вой смеси внутри камеры

T1

d

⋅ ∆T dt

+ ∆TГ

= b ⋅ ∆GГ

+ a12

⋅ ∆TК

+ a13 ⋅ ∆TР

+

g1 ⋅ ∆GП

+

g2

⋅ ∆TП1

+,

g3 ⋅ ∆TГ1 + g4 ⋅ ∆TB

(3)

здесь ∆Gп - расход выпекаемого продукта; ∆Tг1 - температура входящих газовой смеси; ∆Tп1 - температура входящего в печь продукта; ∆Tв - температура окружающего воздуха, которую как правило, считают по-
стоянной.

Изменение подачи топлива (сгораемого газа) считают управляющим воздей-

ствием.

Если все возмущения принять равным нулю, а приращение управления не

равно нулю, то уравнение (3) вырождается в уравнение вида

T1

d∆TГ dt

+ ∆TГ

= b ⋅ ∆GГ

+ a12 ⋅ ∆TК

+ a13 ⋅ ∆TР .

(4)

Если все возмущения, кроме ∆Gп, в уравнении равны нулю, то уравнение (4) приобретает вид

T1

d∆TГ dt

+ ∆TГ

= a12 ⋅ ∆TК

+ a13 ⋅ ∆TР

+ g1∆GП .

(5)

Аналогично следует поступить и при определении уравнения теплового ба-

ланса для выпекаемого продукта.
( )( )M П ⋅ СП ⋅ dTП = FП ⋅ GГП TГ 4 − TП 4 ⋅ dt + FП ⋅ GКП (TК 4 − TП 4 ) ⋅ dt + ,
FП ⋅ α П ⋅ TГ − TП

(6)

где Mп и Cп - масса и теплоемкость выпекаемого продукта; Fп - площадь продукта (заготовки);

σ гп

и

σ кп

-

коэффициенты

излучения;

αп - коэффициент теплопередачи;

5

Tп,Tг,Tк - температуры продукта, газовоздушной смеси, кладки внутри камеры.

Уравнение (6) показывает, что выпекаемый продукт нагревается за счет из-

лучения от газовоздушного объема (от сгорания газа) и от стенок камеры, и за

счет теплопередачи конвекцией от сгоревшего газа.

Проводя линеаризацию уравнения (6) и перегруппировав члены, получим

M П ⋅ СП B

⋅ d∆TП dt

+ ∆TП

=

4 ⋅ TГ 3 ⋅ FП

⋅ GГП B

+ FП

⋅αП

⋅ ∆TГ

+

4 ⋅ TК

⋅ FК B

⋅ GКП

⋅ ∆TК ,

где B = 4 ⋅ TП ⋅ (FП ⋅ GГП + FП ⋅ GКП ) + FП ⋅α П .

(7)

С учетом обозначений:

M Ï ⋅ CÏ B

= T2 - постоянная времени выпекаемого продукта,

4TГ ⋅ FП ⋅ GГП B

+ FП ⋅ α П

= a21

;

4TК

⋅ FК ⋅ GК B

= a22 ;

уравнение (7) примет вид

T2

d∆TП dt

+ ∆TП

= a21 ⋅ ∆TГ

+ a22 ⋅ ∆TК .

(8)

Уравнение теплового баланса роликового конвейерного транспорта
( ) ( )( )M Р ⋅ CР ⋅ dTР = FР ⋅σ Р TГ 4 − TР 4 ⋅ dt + FР ⋅σ КР TК 4 − TР4 ⋅ dt + ,
FР ⋅α Р TГ − TР ⋅ dt − Qр ⋅ dt

(9)

где Mр,Ср - масса и теплоемкость роликов транспорта соответственно; Fp – площадь роликов;

σ σгр, кр- коэффициенты излучения;
αр - коэффициент теплопередачи; Tр,Tк,Tг - температура роликов, стенок печи и газовоздушной смеси внутри печи соответственно; QР - потери в окружающую среду через площадь роликового транспорта. Уравнение (9) показывает, что ролики нагреваются за счет излучения от нагретого газом воздушного объема и стенок печи и за счет теплопередачи от сгоревшего газа, а охлаждаются за счет потерь через концы в окружающую среду.

6

Проведя линеаризацию уравнения (9) и группируя соответствующим образом члены, получим уравнение вида

M Р ⋅CР D



d∆Tр dt

+

∆TР

=

4 ⋅ TС 3

⋅ FР

⋅ σ ГР D

+



⋅αР

⋅ ∆TГ

+

4 ⋅ TК 3

⋅ σ КР D

⋅ FР

⋅ ∆TК

− ∆Qр D

(10)

где обозначены

M Р ⋅ CР D

= T3

-

постоянная

времени

роликов;

4 ⋅T 3Г

⋅ FР ⋅ σ ГР D

+ FР ⋅ α Р

= a31;

4 ⋅ TК 3 ⋅ FР ⋅ σ КР D

= a32

;

D = 4 ⋅ TР ⋅ (FР ⋅σ ГР + FР ⋅σ КР ) + FР ⋅α Р .

1 D

=



;

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

порта

T3



d∆TР dt

+ ∆TР

=

a31 ⋅ ∆TГ

+ a32

⋅ TК .

Уравнение теплового баланса для стенок (кладки) печи примет вид

( ) ( )M

К

⋅ CК



dTК dt

= FК ⋅ σ ГК ⋅ TГ 4 − TК 4

⋅ dt + FК ⋅ σ ПК ⋅ TП 4 − TК 4

⋅ dt +

( )FК ⋅σ КР ⋅ TР4 − TК 4 ⋅ dt + FК ⋅α К ⋅ (TГ − TК )⋅ dt − QК ⋅ dt

(11)

где Mк,Ск - масса и теплоемкость кладки стенок соответственно; Fк – площадь внутренней поверхности кладки;

σ σ σгк, кр, пр - коэффициенты излучения,
αк - коэффициент теплопередачи; Tк,Tп,Tг,Tр - температура поверхности кладки, продукта, газовой смеси внутри камеры, роликового транспорта соответственно.

Уравнение (11) показывает, что кладка печи нагревается за счет излучения от газового объема и теплопередачи от газового объема и охлаждается за счет потерь

в окружающую среду, нагревается и охлаждается за счет излучения на транспорт

и продукт.

7

Проведя линеаризацию уравнения (11) и группируя соответствующим обра-

зом члены уравнения, получим

M К ⋅CК E

⋅ d∆TК dt

+ ∆TК

=

4 ⋅ TГ 3

⋅ FК

⋅ GГК E

+



⋅αК

⋅ ∆TГ

+
,

4 ⋅ TП 3 ⋅ GК E

⋅ FКП

⋅ ∆TП

+

4 ⋅ TР3 ⋅ FК E

⋅ GКР

∆TР



∆QК E

E = 4 ⋅ TR ⋅ (FR ⋅σ ГК + FК ⋅σ КП ) + FК ⋅α К .

Обозначим:

M К ⋅CК E

= T4

-

постоянная времени кладки;

(12)

4 ⋅ T 3Г ⋅ FК ⋅ GГК E

+ FК ⋅ α К

= a41 ;

4 ⋅ TП 3 ⋅ FК E

⋅ GПР

= a44 ;

4 ⋅ TР3 ⋅ FР ⋅ GЕР E

= a43 .

С учетом перечисленных обозначений уравнение (12) примет вид

T4



d∆TК dt

+ ∆TК

= a41 ⋅ ∆TГ

+ a44 ⋅ ∆TП

+ a43 ⋅ ∆TР .

(13)

Запишем систему уравнений для участка печи, приняв все возмущения рав-

ным нулю, кроме подачи в печь пекарной заготовки (продукта):

T1



d∆TГ dt

+ ∆TГ

=b

⋅ ∆GГ

+ a12 ⋅ ∆TК

+ a13 ⋅ ∆TР + g1 ⋅ ∆GП ;

T2

⋅ d∆TП dt

+ ∆TП

= a21 ⋅ ∆TГ

+ a22 ⋅ ∆TК ;

T3

⋅ d∆TР dt

+ ∆TР

=

a31 ⋅ ∆TГ

+ a32

⋅ ∆TК ;

(14)

T4



d∆TК dt

+ ∆TК

= a41 ⋅ ∆TГ

+ a43 ⋅ ∆TР

+ a44∆TП .

Если поставить задачу расчета изменения температур при изменении управ-

ляющего воздействия на одинаковую величину во всех участках, и принять воз-

мущения равными нулю, то уравнение статики примет вид

∆TГ = b ⋅ ∆GГ + a12 ⋅ ∆TК + a13 ⋅ ∆TР ;

∆TП = a21 ⋅ ∆TГ + a22 ⋅ ∆TК ;

(15)

∆TР = a31 ⋅ ∆TГ + a32 ⋅ ∆TК ;

∆TК = a41 ⋅ ∆TГ + a43 ⋅ ∆TР + a44∆TП .

Дав приращение в подаче топлива ∆Gт получим приращение температур на каждом участке.

8

Заключение
Полученная модель пекарной камеры, состоящая из разбитых (условно) на три участка газового объема, состоящая из 12-ти дифференцированных уравнений и имеющая три управляющих воздействия является весьма информативной.
С ее помощью можно имитировать воздействие тепловых потоков и излучающих объемов и поверхностей, что оказывает неоценимую помощь при проектировании систем управления.
Однако при определении коэффициентов дифференциальных уравнений большую погрешность вносят усреднения температур излучающих поверхностей. Устранить этот недостаток в какой-то мере можно, если участки разбить на секции, но тогда возрастет и существенно количество уравнений и их решение вызовет значительные трудности.
Можно избрать и иной путь – составить математическую модель объекта как объекта с распределенными параметрами. Но тогда решение проблем сводится к решению нелинейных дифференциальных уравнений в частных производных с несколькими управлениями и сложными краевыми и начальными условиями, что является практически необозримой задачей, даже при современном уровне вычислительных средств.
Однако можно пойти на компромисс и получить сравнительно грубую, но линейную модель, позволяющую сравнительно легко анализировать тепловые процессы в объекте и проектировать для нее систему управления. В дополнение следует отметить тот факт, что все температурные поля (по участкам) предлагались равномерными. Однако это довольно грубое приближение и допущение. Температурное поле на каждом участке неравномерно и зависит от расположения источников тепла (газовых горелок) смешения продуктов сгорания, местных подсосов и потерь тепла через стенки и транспорт. Поэтому реальные температуры могут отличаться от расчетных. Это также вносит определенные трудности при реализации алгоритмов управления. Для практических нужд эта модель может быть использована в нескольких направлениях. Если ориентироваться на использование средств вычислительной техники, то представленная модель даст возможность построить температурное поле всего теплового объекта в статическом режиме, исходя из этого рассчитать расходы топлива на каждом участке, компенсирующие отклонения температуры от заданной.
9

Сложным является вопрос о стабилизации температуры на каждом из участ-

ков из-за многосвязности объекта, его переменных коэффициентов и многочис-

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

ского конструирования регуляторов (следуя Шаталову или Калману) с неполной

информацией о векторе состояния, применив наблюдателя на основе полученной

математической модели.

Можно пойти и традиционным путем конструирования локальных регулято-

ров с типовыми алгоритмами. Для этих целей модель может быть трансформиро-

вана в передаточные функции (по управлению, по возмущению).

Для реализации поставленной цели в общем виде введем следующие обозна-

чения: x1 – температура дымовых газов, x2 – температура кладки печи, x3 – температура продукта, x4 – температура роликового транспорта. С .учетом введенных обозначений запишем в общем виде систему уравнений для отдельного участка

печи

T1

dx1 dt

+

x1

=

b1

⋅u

+ a12



x2

+ a14

⋅ x4

T2

dx2 dt

+

x2

=

a21

⋅ x1

+ a23

⋅ x2

+ a24

⋅ x4

(16)

T3

dx 3
dt

+

x3

=

a31



x1

+

a32



x2

− b2

⋅G

T 4

dx 4 dt

+

x 4

=

a 41



x 1

+

a 42



x 2

Разделим все уравнения (16) на постоянную времени Ti (i=1..4) и получим систему уравнений в нормальной форме

x1

=

1 T1

⋅ x1

+

a 12
T1

⋅ x2

+

a 13
T1

⋅ x3

+

a 14
T1

+

b 1
T1

⋅u

+ 0⋅G ;

x2

=

a21` T2



x1

+

a12 T2

⋅ x2

+

a24 T2



x3

+

a24 T2



x4

+ 0⋅u

+ 0⋅G;

x3

=

a31 T3



x1

+

a32 T3



x



1 T1



x3

+

0⋅

x4

+

0⋅u

+

b2 T3

⋅G

;

(17)

x4

=

a41 T



x1

+

a42 T



x

+

0⋅

x



1 T



x

+

0⋅u

+

0⋅G.

44

4

10

В векторно-матричной форме уравнение (17) будет

 − a11 a12 a13 a14 

 T1 T1 T1 T 

A

=

      

a21 T2 a31 T3

− a 22 T2
a32 − T3

a23 T2 a33 T3

a24 T2 a34 T3

   ,   

 b1 



T 1



B

=

 

0

 

,

0 

0 

0 

0 

C

=

 

b2

 

.

 

T3 0

 



a 41
T4

a 42
T4

a 43
T4



a 44
T4



Составим характеристический полином в виде

A(s)

=

∆0 (s)

=

a 0



s4

+

a 1



s3

+

a 2



s2

+

a 3



s

+

a 4

.

Определим полином числителя при C=0

B(s)

=

b 0



s3

+

b 1



s2

+

b 2



s

+

b 3

.

полином числителя при B=0

Ñ(s)

=

c 0



s3

+

c 1



s2

+

c 2



s

+

c. 3

Тогда передаточная функция по управляющему воздействию U ≠ 0, G = 0

W 14

(s

)

=

x1 (s ) u(s)

=

B(s) A(s)

=

a0

b0 ⋅ s3 + b1 ⋅ s2 + b2 ⋅ s + b ⋅ s4 + a1 ⋅ s3 + a2 ⋅ s2 + a3 ⋅ s

+

a4

.

Передаточная функция по возмущению G ≠ 0,U = 0

( )W10

s

=

c 0



s3

+

c 1



s2

+

c 2



s

+

c 3

.

a0 ⋅ s4 + a1 ⋅ s3 + a2 ⋅ s2 + a3 ⋅ s + a4

Полученная математическая модель пекарной камеры (печи) как объекта с

сосредоточенными параметрами является весьма познавательной с ее помощью

можно понять взаимодействие тепловых потоков и излучающих объемов и по-

верхностей, что бывает весьма полезно при проектировании систем управления

тепловыми объектами.

Список литературы:
1. Михелев А.А. Справочник по хлебопекарному производству. е1. Оборудование и тепловое хозяйство. – 2-е изд., перераб. и доп. – М.;Пищевая промышленность, 1977 – 368 с.

11

A mathematical model of a baking chamber as an object with lumped parameters
Danin V.B., Kirikov A.U.
Saint-Petersburg State University of Refrigeration and Food Engineering
The paper considers the development of a mathematical model for a tunnel-type baking chamber. The workers propose to represent the mathematical model of a baking chamber (an object) with distributed constants as that with lumped parameters. On the basis of the developed equivalent model with lumped parameters it is possible to formulate and give a solution for the problem of optimal stabilization of temperature conditions in a baking chamber.
Keywords: baking chamber, a set of equations, mathematical model.
12