CC BY
59Преобразование математической модели процесса сушки для управления сушильной установкой барабанного типа
Ю. В. Янюк1, Е. А. Питухин Петрозаводский государственный университет
АННОТАЦИЯ
Статья посвящена проблеме адекватного преобразования сложной нелинейной математической модели процесса сушки в сушильной установке барабанного типа в линейную модель, имеющую аналитическое решение. Рассматриваются существующие модели. Предлагается новый подход к описанию процесса сушки. Приводится аналитическое решение и пример идентификации коэффициентов новой модели.
Ключевые слова: процесс сушки, математическая модель, идентификация, системы управления.
SUMMARY
The paper is devoted to the problem of an adequate transformation of complicated nonlinear mathematical model of drying process in rotary drum dryer to the linear model, which has an analytical solution. Existing models are considered. The new approach to description of drying process is offered. The analytical solution and an example of new model's factors identification are given.
Keywords: drying process, mathematical model, identification, control systems.
Изучение процесса сушки в установке барабанного типа проводилось рядом авторов [1-4] в течение многих лет. Результатом исследований стал вывод, что общая модель сушильной барабанной установки составлена из двух вспомогательных моделей, из которых одна является подробным описанием поведения твердого материала, а другая описывает барабан [2].
Первая модель включает в себя характеристики твердого материала, такие как, например, динамика сушки. Вторая модель, в свою очередь, является описанием оборудования и определяет время распределения материала в барабане и тепловую передачу. При комбинировании этих двух моделей получается несколько математических формул, решение которых дает представление о сушке твердого вещества в барабанной сушильной установке [1].
Общая модель, включающая в себя параметры как твердого материала, так и барабана, обычно состоит из некоторого набора дифференциальных уравнений,
1 Авторы - соответственно преподаватель кафедры технологии металлов и ремонта и инженер ПФ © Ю. В. Янюк, Е. А. Питухин, 2003
описывающих передачу массы и тепла между газообразной и твердой фазами. Такая модель является динамической и наиболее сложной для решения. Обычно она упрощается с целью получения удобных для решения линейных дифференциальных уравнений. Существуют также и статические модели, которые могут быть использованы для определения кривых содержания влаги и температуры для твердого материала и сушащего воздуха в осевом направлении. Однако их адекватность доказать достаточно трудно, из-за того, что внутри барабана сложно измерить содержание влаги и температуру.
Процесс сушки может быть описан кривой, которая показывает изменение значения условного коэффициента скорости сушки в зависимости от уменьшения содержания влаги X в материале [2]. Типичная кривая скорости сушки для конвекционного процесса сушки представлена на рисунке 1.
Рис. 1. Типичная кривая скорости сушки
Как видно из рисунка 1, процесс сушки может быть разделен на три части: начальный период сушки, период сушки с постоянной скоростью и период снижения скорости сушки, зависящие от характеристик твердого материала. Кривая А - В иллюстрирует изменение влажности в материале во время начального периода сушки. В течение этого периода температура твердого материала и его влажной поверхности ниже, чем температура сушки, в результате скорость сушки в области А - В будет возрастать до тех пор, пока температура поверхности твердого материала не достигнет температуры, соответствующей линии В - С. Если температура влажного материала выше, чем температура сушки, начальный период будет соответствовать линии А' - В. Начальный период обычно очень короткий, и поэтому чаще всего им пренебрегают. Скорость сушки за период В - С постоянна и равна прямому участку В - С кривой влажности. Содержание влаги в точке перехода от периода постоянной скорости к пе-
риоду падения скорости называется критическим содержанием влаги X сг . В этой критической точке содержание влаги в материале линейно снижается, и прямая линия становится кривой, которая асимптотически приближается к равновесному содержанию влаги в материале Xeq .
В качестве объекта исследования рассматривается наиболее адекватно описывающая процесс сушки динамическая модель, представленная в виде дифференциальных уравнений в частных производных с распределенными параметрами, как для влажностей, так и для температур газа и твердого материала [4, 5]. Модель охватывает все периоды процесса сушки.
Рассмотрим данную модель [4]. Она получена на основе приемлемых допущений, связанных с физической стороной процесса сушки:
- коэффициенты передачи тепла и массы постоянны;
- передача тепла из-за проводимости материала и сушащего газа незначительна;
- диффузия водяного пара в осевом направлении не принимается во внимание;
- передача тепла из-за радиации незначительна;
- скорость проникновения сушащего газа в осевом направлении постоянна;
- распределение размеров гранул материала постоянно;
- в материале не происходят химические реакции во время сушки;
- температуры материала и сушащего воздуха и влажность материала являются лишь функциями времени и осевой координаты.
дХ д1
■ + v„
дХ д1
= ~Я:
дУ дУ
д(с,м д(сгтг) «к
дг
■ + уг
д!
0 ^т)'
г' —
/й1/в-' , ,, ^^Рт^ш) _ /гр _гр ч _
\ 2 т)
д(РтТт)
■Шу.
+ у„
д!
Используются следующие обозначения: X - влажность материла, кг (H2O)/ кг (материала); Y - сухость воздуха, кг (H2O)/ кг (материала); Tg -температура сушащего воздуха, К; Тш - температура
материала, К; Vm - скорость материала в осевом направлении, м / с ; Vg - скорость сушащего газа в
осевом направлении, м / с; Cm - удельная теплоемкость материала, Дж / кг • К ; - удельная теплоемкость газа, Дж / кг • К ; От - линейная плотность материала, кг / м ; Gg - линейная плотность газа,
кг / м ; Vv - удельный объем барабана, м / м ; av -удельный коэффициент передачи тепла, Дж/(м -К-с); Я- теплота парообразования, Дж / кг ; Rv - скорость сушки, 1 / с .
Назовем данную общую модель М_Ш\ Эта модель не имеет аналитического решения из-за наличия частных производных. Кроме того, в правых частях всех уравнений содержится величина Rv , так называемая скорость сушки, которая может быть представлена в сле-дуюгцемвиде [5]:
ш
К^, =
ная теплоемкость воды, Дж / кг ■ К ; сы сушащего газа, К.
где ^ - удель-
т,
¿р - точка ро-
Это выражение содержит переменную Тш в минус первой степени, что делает уравнения общей модели нелинейными и еще более осложняет решение систе-
Поиск аналитического решения обусловлен задачей получения оптимальных алгоритмов управления сушильной установкой. Преимуществом базирующихся на аналитических зависимостях алгоритмов является высокая скорость выработки управляющих воздействий в условиях управления в реальном времени по сравнению с алгоритмами, созданными на основе численных решений.
Для получения аналитического решения необходимо провести ряд модификаций общей модели, которые вели бы к ее упрощению. При этом на каждом шаге нужно проводить сравнение результатов расчетов упрощенных моделей с результатами расчетов общей модели, полученных численными методами.
Ранее [1, 3, 5, 6] преобразование общей модели сводилось как к преобразованию левой части уравнений, так и к представлению коэффициента скорости сушки
где
ввидеЬЦ = ^Х + к2Тт + к3Тё [8]
коэффициенты к^к^кз определяются экспериментально и зависят от физических свойств материалов, подвергаемых сушке. Представление коэффици-
ента К у в виде К У или К у зависит от выбора модели.
Решение общей системы М_ЬТ было получено нами численным методом, путем разложения системы уравнений в частных производных в систему разностных уравнений, при представлении искомых функций дискретным множеством своих значений. Вводилась подходящая сетка значений координат по времени 1 и длине барабана Ь . При этом с учетом длины барабана Ь =3 ми времени прохождения материала через барабан (время задержки), равном 627 с, разбиение (соответствующее скорости прохождения материала
Ут =0,00478 м/с) вдоль оси Ь составляло 30 шагов,
а по оси 1 - 6270 шагов. Такое соотношение количества узлов сетки выбиралось из соображений пригодности решения, выявленного при исследовании устойчивости решения системы, полученной разностной аппроксимацией.
Обозначим общую модель с коэффициентом К У как а с Я - М_LT_RV2.
Результаты численного решения на примере общей модели М_LT_RV1 представлены на рисунках 2, 3 и 4.
Из анализа поверхностей, представляющих решения, видно, что в качестве начальных условий интегрирования системы М_ЬТ выбираются нулевые начальные условия по всей длине барабана за исключением точки входа. Физически это означает, что на момент начала подачи материала и газа барабан пуст. Такое допущение принимается во всех нижерассматриваемых моделях для обеспечения одинаковых начальных условий моделирования.
Рис. 3. Зависимость температуры сушащего газа Tg от времени и длины барабана
Рис. 4. Зависимость температуры материала Tm от времени и длины барабана
Для сравнения результатов решения общих моделей 3-мерные графики решений преобразованы в 2-мерные во временной точке начала выхода материала из барабана. Результаты численного решения общих моделей М_LT_RV2 и М_LT_RV1 на конечный момент времени представлены соответственно на рисунках 5 и 6.
Рис. 2. Зависимость влажности материала X от времени и длины барабана
Рис. 5. Графики решения системы М_LT_RV2
Рис. 6. Графики решения системы М_ЬТ_ЯУ1
Совместное изображение двух данных решений представлено на рисунке 7. Из рисунка видно, что поведение решений качественно похоже. Решение системы
М_ЬТ_ЯУ2 со "сложным" Яу^ теоретически более
верно отражает суть процесса сушки, но практическая пригодность общей модели появляется с коэффициентом Я у , который можно настроить под скорость сушки конкретного материала.
600
300
100
-
* * -А. А.,
.1 .■ ■ .и«-1 •' .Гь^л-1 .1 '.В .1
V
0,03 -- 0,02 0,01 0
-0,02 -0,03
^ -0,06 0,07
0,5
1,5
2,5 Ц т 3
-Тд2 ■
■Тд1
-Х2 ■
Рис. 7. Графики решения систем М_ЬТ_ЯУ1 (XI, Тд 1, Тт1) и системы М_ LT_RV2 (Х2, Тд2, Тт2)
Модели М_LT_RV1 и М_LT_RV2 ведут себя асимптотически одинаково. Это предположение подтверждается тем, что можно добиться полного совпадения их
решений, подбирая значения коэффициентов кц, к 2 и кз так, чтобы коэффициент Я™ заменил коэффициент Я у в правой части модели [7].
Следует подчеркнуть, что сведение решения системы М_LT_RV1 к решению системы М_LT_RV2 имеет теоретическое значение в качестве проверки асимптотического равенства. Учитывая прикладной характер задачи, при решении системы М_LT_RV1 должны отражаться характеристики процесса сушки каждого конкретного материала. Это означает, что коэффициенты кц, к2 и кз должны быть различными для
различных материалов. Определение этих коэффициентов есть задача идентификации и настройки параметров для каждой анализируемой математической модели. Таким образом, в качестве коэффициента
сушки Я у выбирается Я у , а в качестве базовой
модели для сравнения с ней всех последующих модификаций выбирается модель М_LT_RV1.
Как было отмечено выше, ряд авторов занимались поиском решения общей системы М_ЬТ, при этом они проводили модификацию модели, основываясь на различных допущениях. Так, в работе [6] автор, допуская,
что скорость прохождения материала У т вдоль барабана постоянна, избавляется в уравнениях от производной по длине барабана и, следовательно, превращает уравнения в линейные дифференциальные. Результаты решения данной модели, обозначенной М_TP_RV1, представлены на рисунке 8.
600 500 400 300 200 100
-
X, ------ ------
ч А.
—
4
0,025
0,02
0,015
0,01
0,005
0,5
1,5
2,5 ^ т 3
-Тд ---♦--■Тд1
-Тт ---■--■Тт1
-Х---А---Х1
Рис. 8. Графики решения систем М_LT_RV1 (Х1, Тд1, Тт1) и системы М_TP_RV1 (X, Тд, Тт)
В работе [1] автор предлагает провести линеаризацию модели вокруг рабочей точки, упрощая таким образом модель. Избавляется от параметра 1 - длины барабана и записывает зависимости от параметра 1 - времени сушки. Уравнения превращаются в линейные дифференциальные уравнения. Найдя решение системы в нескольких рабочих точках, автор получает кривые влажности и температур сушащего газа и материала в зависимости от длины барабана. Результаты представлены на рисунке 9.
Решения систем М_TP_RV1 и М_TY_RV1 близко подходят к решению базовой системы для переменных
Тт и Tg, но расходятся для переменной X . При
этом модели систем М_TP_RV1 и М_TY_RV1 основаны на дополнительных допущениях и решались численно, что не позволяло ввести в решение функциональные зависимости от параметров управления
( Т , Vg и Ут) и настройки ( кц, к2 и кз ), где
Т
де в барабан
начальная температура сушащего газа на вхо-
500
400
0
0
0
2
200
0
0
2
Х1
rJK
500 400 300 200
X, kg/kg
к
» « 4-1-
* А.
----------
0,025
0,02
0,015
0,01
0,005
2,5 L, m 3
-Tg -------Tg1
-Tm ■
■Tm1
-X-
X1
Рис. 9. Графики решения систем М_ТУ_ЯУ1 (X, Тд, Тт) и системы М^Т_ЯУ1 (XI, Тд1, Тт1)
Для решения проблемы управления сушильной установкой в реальном времени необходимо стремиться к получению аналитического решения. Поиск математической модели, дающей адекватное аналитическое решение в виде функциональных зависимостей, является задачей авторов статьи.
Преобразование математической модели, предлагаемое в данной работе, основано на разделении процесса сушки на два этапа - стационарный и нестационарный. Из экспериментальных данных видно [1], что 1-й этап сушки представляет собой процесс, где скорость сушки практически не зависит от времени. Этот этап, как показывают эксперименты, является самым продолжительным по времени и представляет собой собственно сушку, т. е. испарение влаги из материала.
Основываясь на допущении о неизменности скорости сушки во времени и применении асимптотического подхода к решению общей системы М_ЬТ [8], из уравнений общей модели исключаются члены, содержащие параметр времени. На основании этого модель преобразуется в систему линейных дифференциальных уравнений М_Ь_УР с производной по длине барабана 1.
V .
т д1 V
Такая модифицированная система линейных дифференциальных уравнений с коэффициентом решается строгим аналитическим способом. Общее аналитическое решение этой системы было найдено интегрированием дифференциальных уравнений в квадратурах. Частное решение (задача Коши) получено с использованием программных средств с символьными вычислениями и в силу громоздкости не приводится в работе.
Чтобы оценить, какая из моделей ряда М_Ь_УР является наиболее адекватной общей модели М_ЬТ, проведем ряд сравнений. Для этого сравним решения модифицированной системы с различными коэффициентами Rv с соответствующими базовыми моделями,
то есть модель М_Ь_УР_Ю/1 - с М_ЬТ_Ю/1, а модель М_Ь_УР_Ю/2 с - М_ЬТ_Ю/2.
Результаты сравнения аналитического решения модели М_Ь_УР_ЯУ1 с численным решением базовой модели М_ЬТ_ЯУ1 представлены на рисунке 10.
Рис. 10. Графики решения систем М_Ь_УР_ЯУ1 (X, Тд, Тт) и системы М_ЬТ_ЯУ1 (XI, Тд1, Тт1)
При выборе модели М_Ь_УР_ЯУ2 предлагается два способа решения для сравнения с базовой моделью: численный и аналитический. Сравнение с численным решением представлено на рисунке 11.
ÔY
д(СтТт) avVv
m
m m dl
G
m
d(Ç T•) av
V -2-2- = —}LJL(T -T )-
S dl G, & m
g
; m t? - m ~F< V
500
300
-
0,03 -■ 0,02
-0,02
-0,06 -0,07
0,5
1,5
2,5 L, m 3
-Tg2 —....Tg -
- Tm2 ---■--■ Tm ■
-X2 .......X
Рис. 11. Графики решения систем М_Ь_УР_ЯУ2 (X, Тд, Тт) и системы М_ЬТ_ЯУ2 (X2, Тд2, Тт2)
00
0
0
0
0,5
,5
2
600
400
200
00
0
0
2
Часть аналитического решения системы М_L_YP_RV2, для которой возможно было решить задачу Коши, представлена на рисунке 12.
600
400 300 200 100
T, K
X, kg/kg
0,03 -- 0,02 0,01 0
-0,01 -0,02 -0,03 -- -0,04 -0,05 ^ -0,06 0,07
2,5 L, m 3
-Tg2 ■■■
■Tg ■
-Tm2--
■Tm
-X2 ■■■
■X
Рис. 12. Графики решения систем М_L_YP_RV2 (X, Тд, Тт) и системы М_LT_RV2 (Х2, Тд2, Тт2)
Из рисунков видно, что наибольшее визуальное совпадение модифицированных моделей с базовыми наблюдается на рисунке 10. Количественно это подтверждается относительными оценками точности совпадения графиков, приведенными в таблице. Оценки вычислялись по формуле
1 п П! = 1
х(1)_х(0)
X
(0)
• 100% , где Х-0^ - точка ба-
-(1)
Л/Ч1/
зовой кривой; X - точка сравниваемой кривой; п - число точек кривой.
Таблица
Процентное отклонение для графиков переменных на рисунках 10, 11, 12
Переменная Рисунок 10 Рисунок 11 Рисунок 12
X 5,449% 14,545% 42,316%
T 0,139% 0,345% 2,861%
T 1,703% 5,562% 17,985%
Данное сравнение показывает, что наименьшее отклонение в процентном отношении наблюдается у модели М_L_YP_RV1 (рис. 10). Это позволяет использовать ее в дальнейших задачах идентификации и оптимального управления сушильной установкой.
Проверка возможности решения модели М_L_YP_RV1 соответствовать некоторым экспериментальным кривым сушки или теоретическим законам сушки осуществлялась путем настройки коэффициента Я У . С этой целью решалась задача идентификации параметров к1, к 2 и кз . При этом использовались вероятностно-статистические методы
оптимизации (метод Монте-Карло по наилучшей пробе и метод Монте-Карло с сужением интервала поиска). Значения коэффициентов, например, сводящие решения модели М_L_YP_RV1 к М_LT_RV2, равны:
к1 =1,059*10-3; к2 =0; к3 =2,117*10-7.
Попытки применить классическую теорию управления для барабанной сушильной установки предпринимались в работе [1]. Все алгоритмы управления строились на базе моделей, не имеющих аналитического решения. Это не позволяет найти прямую зависимость между параметрами системы и не дает возможности получить выражение для плотности распределения выходных переменных.
Предлагаемая модель М_L_YP_RV1 имеет устойчивое аналитическое решение с выражением параметров настройки и управления в виде функциональной зависимости. Таким образом, данная модель наилучшим образом подходит для создания на ее базе оптимального алгоритма управления барабанной сушильной установкой [9].
СПИСОК ЛИТЕРАТУРЫ
1. Yliniemi L. Advanced Control of a Rotary Dryer / Oulun Yliopisto. Oulu, 1999.
2. Balchen J. G., Kenneth I. M. Process control: structures and applications. New York: Van Nostrand Reinhold company, cop. 1988. 540 p.
3. Anilkumar S. M., Arun S. Mujumdar. Drying of Solids: Principles, Classification, and Selection of Dryers // Drying 94: proceedings of the 9th International Drying Symposium (IDS '94). Gold Coast, Australia, 1994. August 1-4.
4. Brambilla A., Nardini, G., Stabert Z. Dynamic model of rotary dryer // 2 nd symposium on the use of computers in chemical engineering. Tchekhoslovakia, 1973. Vol. II. Р. 218-231.
5. Yliniemi L., Arola A., Jutila E., Leiviska K. Modelling and Simulation of a Pilot-Plant Rotary Drier // Proceedings of Process Engineering Department / University of Oulu. Oulu, 1980.
6. Питухин Е. А. Математическая модель управления качеством работы сушильной установки барабанного типа // Труды Петрозаводского государственного университета. Сер. Прикладная математика и информатика. Вып. 6. 1997. С. 71-76.
7. Маслов В. П. Асимптотические методы и теория возмущений. М.: Наука, 1988. 312 с.
8. Арнольд В. И. Обыкновенные дифференциальные уравнения. М.: Наука, 1984. 272 с.
9. Pitukhin E. A. Optimal Design of the Control System Elements Admissions // Proceedings of TOOLMET'97 - Tool Environment and Development Methods for Intelligent Systems / University of Oulu. Control Engineering Laboratory, Report A. No 6, April 1997.
500
0
0
0,5
,5
2
