ФЕДЕРАЛЬНОЕ ГОСУДАРСТВЕННОЕ АВТОНОМНОЕ ОБРАЗОВАТЕЛЬНОЕ УЧРЕЖДЕНИЕ ВЫСШЕГО ОБРАЗОВАНИЯ
«БЕЛГОРОДСКИЙ ГОСУДАРСТВЕННЫЙ НАЦИОНАЛЬНЫЙ
ИССЛЕДОВАТЕЛЬСКИЙ УНИВЕРСИТЕТ»
(НИУ «БелГУ»)
ИНСТИТУТ ИНЖЕНЕРНЫХ ТЕХНОЛОГИЙ И ЕСТЕСТВЕННЫХ НАУК
КАФЕДРА ОБЩЕЙ МАТЕМАТИКИ
МАТЕМАТИЧЕСКОЕ МОДЕЛИРОВАНИЕ И ОПТИМИЗАЦИЯ
ОБРАБОТКИ КЕРАМИЧЕСКИХ ИЗДЕЛИЙ В КАМЕРНОЙ СУШИЛКЕ
Магистерская диссертация
обучающегося по направлению подготовки
01.04.01 Математика
очной формы обучения,
группы 07001534
Ярославцева Сергея Викторовича
Научный руководитель
к.п.н. доцент Ерина Т.А.
Рецензент
к.т.н., доцент Красовская Л.В.
БЕЛГОРОД 2017
СОДЕРЖАНИЕ
УCЛОВНЫЕ ОБОЗНАЧЕНИЯ ........................................................................ 4
ВВЕДЕНИЕ........................................................................................................ 6
1 АНАЛИЗ ОБЪЕКТА МОДЕЛИРОВАНИЯ И ВОЗМОЖНОСТИ ЕГО
ОПТИМИЗАЦИИ ...................................................................................................... 10
1.1 Математические методы исследования реальных объектов ............. 10
1.1.1. Классификация математических моделей .............................................. 10
1.1.2 Методы экспериментальных исследований ............................................ 15
1.2 Моделирование физических процессов .................................................. 17
2
ПОСТРОЕНИЕ
МАТЕМАТИЧЕСКОЙ
МОДЕЛИ
ОБРАБОТКИ
КЕРАМИЧЕСКИХ ИЗДЕЛИЙ ................................................................................ 28
2.1 Разработка математической модели распределения теплоносителя в
камерной сушилке периодического действия ........................................................ 28
2.2 Модификация математической модели термовлажностной обработки
керамической плитки ................................................................................................ 47
3
РАЗРАБОТКА
И
РАСЧЕТ
ОПТИМАЛЬНОГО
ПО
БЫСТРОДЕЙСТВИЮ РЕЖИМА СУШКИ КЕРАМИЧЕСКИХ ИЗДЕЛИЙ ..... 60
ЗАКЛЮЧЕНИЕ ............................................................................................... 77
СПИСОК ИСПОЛЬЗОВАННЫХ ИСТОЧНИКОВ ..................................... 79
ПРИЛОЖЕНИЕ А ........................................................................................... 86
Программа моделирования распределения тепловых потоков в пустой
камере ......................................................................................................................... 86
ПРИЛОЖЕНИЕ Б............................................................................................ 92
Результаты моделирования процесса распределения тепловых потоков в
пустой камере ............................................................................................................ 92
ПРИЛОЖЕНИЕ В ........................................................................................... 94
2
Программа
моделирования
распределения
тепловых
потоков
в
заполненной изделиями камере ............................................................................... 94
3
УCЛОВНЫЕ ОБОЗНАЧЕНИЯ
Обозначение Наименование
Кн
Коэффициент, характеризующий неравномерность сушки
V
Вектор, показывающий скорость движения теплоносителя в
ρ
сушильной камере
Плотность,
характерная для воздушной среды камеры
T
Температура, характерная для теплоносителя в сушильной
µ
камере
Динамический
коэффициент, характеризующий вязкость
среды
Функция тока
Функция, характеризующая завихренность
α
Коэффициент, характеризующий температуропроводность
λ
(влагопроводность)
Коэффициент,
характеризующий теплопроводность
Cр
(влагопроводность)
Показатель теплоемкости воздуха
х, y
Пространственные координаты
S
Площадь входных отверстий
Q
Производительность подающего вентилятора
H
Высота камеры
D
Длина камеры
U
Влагосодержание изделий
K
Коэффициент, характеризующий потенциалопроводность
γ0
Объемный вес абсолютно сухого материала
δ
Коэффициент, характеризующий термоградиентность
с
Показатель теплоемкости керамической плитки
γ
Удельный вес материала
d
Толщина материала
h
Коэффициент, характеризующий теплоотдачу (влагоотдачу)
l
Линейный размер тела
ψ
Коэффициент, характеризующий чувствительность сырья
4
Е
Модуль упругости
σ
Напряжения в материале
Ku
Коэффициент инерционности воздушного потока
L
Длина поперечного сечения канала
F
Площадь поперечного сечения канала
R
Аэродинамическое сопротивление
H
Депрессия
Нижний и верхний индексы
min
Минимальная величина
max
Максимальная величина
BX
Величина на входе
ВЫХ
Величина на выходе
X
Величина проекции вектора на ось х
У
Величина проекции вектора на ось у
i,j
Величины для решения сеточных задач
В
Величины параметров теплоносителя
ТП
Величины температурных полей
ВП
Величины полей влагосодержаний
ср
Средняя величина
5
ВВЕДЕНИЕ
Актуальность работы. При проектировании различных технологий
наряду с опытно-экспериментальными исследованиями часто становится
необходимым
создание
имитационных
позволяющих
осуществлять
комплексное
и
математических
исследование
моделей,
процесса
для
определения его технологических возможностей и получения оптимального
результата. Известно, что керамические изделия являются экологически
чистыми и надежными в эксплуатации материалами. Современное развитие
строительной индустрии вызывает необходимость увеличения производства
изделий из керамики, которое относится к разряду материалоемких и
энергоемких производств [25].
Технологический процесс производства керамических изделий состоит из
нескольких последовательных этапов, среди которых этап сушки является
наиболее важным: и с точки зрения его энергоемкости (около 30 %), и с точки
зрения длительности (до 90 %), и с точки зрения обеспечения качества
выпускаемой
продукции.
Поэтому
показатели,
соответствующие
технологическому этапу сушки, способствуют определению экономической
эффективности
производства
в
целом.
Используемые
режимы
сушки
керамических изделий зачастую оказываются нерациональными для улучшения
показателей,
Трудностью
характеризующих
при
организации
эффективность
процесса
сушки
процесса.
Затруднения
является
существенная
неравномерность протеканий теплофизических процессов, происходящих как
внутри отдельных сушильных камер, так и в сушильной установке в целом.
Кроме того, проведение стадии сушки изделий нужно осуществлять с учетом
климатических условий добычи и предварительной обработки исходного
материала, а также его физико-химических характеристик [16, 21, 26, 40, 41, 44,
45, 46, 55, 78].
Осуществление технологического процесса сушки керамических изделий
необходимо реализовывать с учетом широкого диапазона изменения входных
6
параметров
[6,
42-44,
Поэтому
57].
для
повышения
эффективности
функционирования камерных сушилок, используемых для сушки изделий из
керамики, нужно определение оптимальных технологических регламентов,
действующих в широком диапазоне изменений входных параметров.
Анализ научной литературы показал, что используемые практически
методы
определения
технологических
регламентов
при
керамическом
производстве недостаточно учитывают взаимосвязь разных теплофизических
процессов, которые протекают в камерных сушилках. К ним относятся:
распространение
теплоносителя
термовлажностная
теплоносителя,
обработка
в
объеме
керамических
поступающего
из
общего
сушильной
изделий;
канала,
камеры;
распределение
между
камерами;
распределенность сушильной установки [15, 18, 22, 24, 34, 35, 37, 76, 84]. В
работах А.Д. Альтшуля, В.И. Бодрова, А.В. Золотарского, И.С. Кашкаева, В.
Каста, Р.Б. Кея, А.В Лыкова, П.Д Лебедева, М.В. Месрова, В.Т. Морозовского,
В.В. Перегудова, К.А. Похратяна, В.В. Полякова, О.С. Соболева, А.Ф.
Чижекого, А.Д. Цепина, Е.Ш. Шейнмана, Л.А. Щукина, Т.К. Шервуда, Р.
Френкса, и др. отдельно исследованы теплофизические и аэродинамические
процессы, которые применяются в условиях сушки керамических изделий,
однако
не
решена
задача
определения
оптимальных
технологических
регламентов на конкретном технологическом оборудовании. Решить данную
задачу можно только применяя метод математического моделирования. Для
этого необходима разработка математического описания процессов, которые
протекают, как в камере сушилки, так и в межкамерном пространстве [29, 69].
Объект исследования – камерная сушилка периодического действия.
Предмет исследования – процесс распределения тепловых потоков в
камерной сушилке периодического действия.
Целью работы является математическое моделирование и оптимизация
процесса сушки керамических изделий в камерной сушилке.
Задачи магистерской диссертации:
7
1.
Изучить
теоретико-методические
основы
моделирования
распределения тепловых потоков.
2. Обосновать выбор средств для моделирования распределения тепловых
потоков.
3. Разработать математическую модель распределения теплоносителя в
камерной сушилке периодического действия.
4. Произвести расчет параметров оптимального управления режимом
обработки керамических изделий.
Научная новизна работы заключается в следующем:
Предложена математическая модель распределения теплоносителя в
камерной сушилке периодического действия в процессе термовлажностной
обработки керамических изделий.
Разработан
алгоритм
поиска
оптимального
управления
технологическим режимом сушки керамических изделий в сушильной камере
периодического действия, позволяющий сократить длительность обработки при
соблюдении требуемого качества продукции.
Практичная значимость работы заключается в том, что предложенный
алгоритм поиска оптимального управления технологическим режимом сушки
керамических изделий, разработанный на основе математической модели,
позволяет учесть взаимодействие камеры и физическую распределенность
параметров объекта, тем самым сократить длительность обработки.
Методы исследования. В работе использовали методы интегрального и
дифференциального
исчисления,
математической
физики,
гидро-
и
аэродинамики, теории многосвязных однотипных и оптимальных систем,
математического моделирования, линейного программирования, а также
методы проектирования аналоговых и цифровых систем автоматического
управления.
Структура и объем работы. Диссертационная работа состоит из
введения, трех глав, заключения и 3 приложений, списка литературы из 88
8
наименований, и содержит 105 страниц, в том числе 85 страниц основного
текста, 29 рисунков и 6 таблиц.
9
1 АНАЛИЗ ОБЪЕКТА МОДЕЛИРОВАНИЯ И ВОЗМОЖНОСТИ
ЕГО ОПТИМИЗАЦИИ
1.1 Математические методы исследования реальных объектов
1.1.1. Классификация математических моделей
Интенсификация
необходимость
высокопроизводительных
использовать
специальное
процессов
программное
определяет
обеспечение
и
вычислительные комплексы, которые выполняют функции автоматического
управления и оптимизации режимов работы оборудования. Для этого
необходим анализ технологического процесса и работы оборудования с
применением
современных
методов
их
исследований
(моделирование
процесса).
Моделирование
является процессом исследований реальной системы,
который включает построение модели, исследование ее свойств и перенос
полученных сведений на моделируемую реальную систему. Общим для
функций моделирования является описание, объяснение и прогнозирование
поведения реальной системы.
В науке, технике и производстве при решении разнообразных задач
применяются новые эффективные методы исследования. Эффективный анализ
механизмов явлений и способов управления процессами производств возможен
при выявлении взаимосвязей между факторами, которые определяют ход
процессов, и представлении их в количественном виде – в качестве
математической модели (рисунок 1.1), отображающей самые значимые стороны
процесса.
Использование математических методов при планировании производства
продукции позволяет определить оптимальное соотношение затрат на материал
и оборудование, необходимых для достижения поставленной цели, и получения
наибольшей прибыли, обеспечивая при этом выпуск качественного сырья.
Благодаря реализации полученных математических моделей и результатов при
10
помощи использования ЭВМ, происходит не только развитие отдельно взятого
предприятия, но и экономики в целом.
Объект
(исследуемый процесс)
Знания об
аналогичных
Наблюдения за
объектом и
эксперименты
объектах
Общие принципы,
законы сохранения
и пр.
Абстрагирование,
идеализация, анализ и
синтез
Математическая модель
Рисунок 1.1 – Источники информации, используемые при построении математической
модели
Модель является совокупностью уравнений, условий и алгоритмических
правил, позволяющих:
– получать информацию о процессах, которые протекают в объектах;
– рассчитать, проанализировать и спроектировать системы;
– получать информацию, которая может быть использована для
управления моделируемым объектом.
Различают модели следующих типов:
1. По типу представления:
а)
статистические
(эмпирические)
модели
можно
получить
при
статистических обработках экспериментальных данных, которые собраны на
изучаемом объекте. Они имеют простую структуру, представляются в виде
полиномов;
11
б) физико-химические модели (аналитические или теоретические)
получают в результате обширных теоретических исследований в виде сложных
систем уравнений (системы алгебраических, обыкновенных дифференциальных
уравнений или уравнений в частных производных), достаточно точно
описывающих процессы.
2. По поведению моделей во времени:
а) статические (стационарные) модели описывают поведение объекта,
процесса или системы в какой-либо момент времени (позволяют произвести
расчет статических характеристик и анализ установившихся процессов в
системах);
б) динамические модели (модели динамики) описывают переходные
процессы, нестационарные состояния, имеют вид передаточных функций или
систем дифференциальных уравнений (используются как при анализе
переходных
процессов
в
системе
и
синтезе
систем
с
требуемыми
динамическими показателями, так и при разработке систем управления
объектом регулирования).
3. По виду закона:
а) линейные;
б) нелинейные.
4. По характеру отображаемых свойств объекта:
а) функциональные модели, в которых величины, характеризующие
явление или объект, выражаются количественно. Такие модели предназначены
для отображения физических или информационных процессов, протекающих в
объекте;
б)
структурные
−
характеризуют
структуру
сложного
объекта,
представленного в виде отдельных частей, имеющих определенные связи,
которые не поддаются количественному измерению. В этом случае наиболее
целесообразным является использование теории графов.
5. В зависимости от характера исследуемых реальных процессов и
систем:
12
а)
детерминированные
−
предполагают
отсутствие
случайных
воздействий, данные модели имеют точно установленные элементы, а также
достаточно точно определяется и поведение системы. Важно отметить, что при
изменении
любого
параметра
системы
значения
выходных
величин
определяется однозначно. Для построения таких моделей в основном
используются функциональные зависимости (алгебраические и интегральные
уравнения);
б) стохастические − носят случайный характер процессов в исследуемых
объектах и системах, результат определяется заданными параметрами с
определенной степенью достоверности. Для описания моделей этого типа
используют методы теории вероятности и математической статистики [18, 35,
60, 71].
Можно выделить 2 класса моделей:
1 класс – вещественные модели:
•
натурные модели, включающие в себя реальные объекты, процессы
и системы, которые используются при выполнении научных, технических и
производственных экспериментов;
•
физические – макеты или муляжи, воспроизводящие физические
свойства оригиналов;
•
математические
–
аналоговые,
структурные,
геометрические,
графические, цифровые и кибернетические модели;
2 класс −
•
идеальные модели:
наглядные (схемы, карты, чертежи, графики, графы, аналоги,
структурные и геометрические модели);
•
знаковые
(символы,
алфавит,
языки
программирования,
упорядоченная запись, топологическая запись, сетевое представление);
•
математические (аналитические, функциональные, имитационные,
комбинированные модели) [5, 7, 31, 56, 65, 71].
Описывая один и тот же объект или процесс, используя различный набор
параметров,
можно
разработать
несколько
13
моделей,
соответствующих
различным целям изучения и исследования. При этом необходимо учитывать,
что модель должна быть наглядной, удобной для проведения экспериментов и
затраты на ее создание не должны превышать затрат на создание оригинала.
Существует ряд
требований,
по которым можно оценить правильно
построенную модель:
1.
требование адекватности – соответствия модели реальному объекту
или технологическому процессу, результатам исследований, сопоставимым
уровню сложности и организации, множеству свойств;
2.
соответствие полученной модели решаемой задаче
– модель
должна строиться для решения определенного класса задач или конкретной
задачи исследования системы;
3.
упрощение при сохранении существенных свойств системы –
модель должна быть абстрагирована от второстепенных деталей реальной
системы;
4.
соответствие
между
требуемой
точностью
результатов
моделирования и сложностью модели;
5.
баланс
погрешностей
различных
видов
(например,
баланс
систематической погрешности моделирования и погрешности исходных
данных);
6.
многовариантная реализация элементов модели;
7.
блочное строение – разделение модели по этапам и режимам
функционирования системы.
Подходы к построению моделей:
– анализ функционирования системы;
– проведение ограниченного эксперимента на самой системе;
– использование аналога;
– анализ исходных данных.
Модели можно также классифицировать по степени адекватности
описания поведения реальной системы. Уровень абстрагирования в описании
объекта определяет иерархический уровень. Моделирование большинства
14
технических объектов можно выполнить на микро-, макро- и метауровнях [47,
74].
Математические модели по принципам построения бывают:
а) имитационные − вычислительные эксперименты, имитирующие
поведение реальных объектов, процессов или систем.
б) аналитические − представлены в виде явных функциональных
зависимостей, которые в зависимости от математической проблемы можно
представить в виде:
•
алгебраических,
трансцендентных,
дифференциальных
и
интегральных уравнений,
•
аппроксимационных задач, таких как интерполяция, экстраполяция,
численное интегрирование и дифференцирование,
•
задач оптимизации,
•
стохастических проблем.
По виду входной информации:
•
непрерывные,
характеризуются
наличием
устойчивых
математических связей, а информация и параметры – непрерывны,
•
дискретные – модели, противоположные непрерывным.
1.1.2 Методы экспериментальных исследований
Научные исследования реальных процессов возможно проводить с
использованием
теоретических
или
экспериментальных
методов.
Для
современных условий развития науки и техники характерно проведение
комплексного исследования объекта на основе новых методологий и
технологий научных исследований.
Массовое
использование
электронно-вычислительных
машин
при
математическом моделировании является очень мощной теоретической и
экспериментальной базой, позволяющей считать вычислительный эксперимент
15
новой технологией и методологией при выполнении научных и прикладных
исследований.
Для вычислительного эксперимента характерно то, что он является
экспериментом над математической моделью объекта, выполненным на
электронно-вычислительных машинах. Суть его заключается в том, что одни
параметры модели вычисляются благодаря другим параметрам. В результате
делаются
выводы
о
свойствах
явлений,
описываемых
с
помощью
математической модели [68, 77].
Необходимость
проведения вычислительного эксперимента возникает
при условии, если невозможен натурный эксперимент. Использование
вычислительного эксперимента, по сравнению с натуральным, обладает
несколькими существенными преимуществами: он более дешевый и более
доступный, его гораздо легче изменить, он даѐт более подробную информацию,
на его подготовку и проведение требуется значительно
меньше времени.
Кроме того, можно отметить, что при выполнении вычислительного
эксперимента можно выявить границы применения математической модели,
которые позволяют спрогнозировать эксперимент в естественных условиях. В
этой
связи,
применение
выше
названного
эксперимента
ограничено
математическими моделями, которые участвуют в выполнении научных
исследований. Это не даѐт возможности для полной замены натурального
эксперимента вычислительным. Выход из этой ситуации – проводить сложные
эксперименты, сочетающие вычислительные и натурные эксперименты, при
этом использующие широкий спектр математических моделей: прямые задачи,
обратные задачи, оптимизированные задачи, задачи идентификации.
Применение
на
основных
этапах
решения
задач
электронно-
вычислительных машин можно считать одним технологическим циклом
выполняемого эксперимента. Наблюдается тесная взаимосвязь на всех этапах
технологического
цикла
вычислительного
эксперимента.
В
результате
получаем продукт заданной точности в короткий период с адекватным
16
количественным описанием свойств изучаемых реальных объектов при тех или
иных условиях.
Для проведения натурного эксперимента существует множество методик
[48-55], таких как адаптационная оптимизация дискретных и непрерывных
технологических процессов, метод крутого восхождения, ортогональное
планирование,
регрессионный
композиционное
метод
прогнозирования
анализ
при
производственных
пассивном
ротатабельное
процессов,
эксперименте,
униформное
центральное
планирование,
симплекс-
планирование, метод случайного баланса, эволюционное планирование,
экстремальный эксперимент и др. [13, 85].
1.2 Моделирование физических процессов
Для математического моделирования характерным является то, что
решить любую физическую задачу возможно с использованием теоретического
пути. При этом использование теоретического решения задачи зависит от
степени сложности ее математической модели. Сложность математической
модели
взаимосвязана
со
сложностью
описываемого
с
ее
помощью
физического процесса. Чем он более сложен, тем более проблемным является
применение такой модели в качестве расчетов.
В
большинстве
математической
случаев
сложности
сложно
используемой
решить
модели.
задачу
Здесь
вследствие
чаще
всего
используют численные методы решения задачи, эффективную реализацию
которых возможно осуществить, используя компьютер.
Таким образом, для физических исследований, построенных на основе
сложных математических моделей, характерным является использование
компьютерного математического моделирования.
Математическое моделирование в физике – чрезвычайно важный метод
исследования. Вычислительная физика (computational physics) является одним
17
из фундаментальных разделов физики. Причиной этого является максимальное
проникновение в физику математических методов, так как реальное решение
математических задач с помощью только традиционных методов очень
ограничено. Из множества причин можно выделить 2 наиболее часто
встречающиеся:
это
нелинейность
многих
физических
процессов
и
необходимость исследований совместных движений многих тел, для которых
нужно решить системы большого числа уравнений [10].
Исследования с помощью компьютерной техники физических процессов
называется вычислительным экспериментом. Для вычислительной физики
характерным является осуществление взаимосвязи между теоретической
физикой, которая поставляет математические модели, и экспериментальной
физикой,
реализующей
компьютерах.
С
виртуальные
использованием
физические
компьютерной
эксперименты
графики,
на
обработка
полученных результатов исследований становится наглядной, что не только
улучшает их восприятие, но и помогает интерпретировать [9].
Практическое использование метода математического моделирования
возможно для определения оптимальных технологических регламентов сушки
керамических изделий.
Существующие сушилки могут быть подразделены на группы по ряду
признаков [32-34]. Наиболее характерными из них являются:
цикличность: непрерывного и периодического действия;
способ
передачи
тепла
к
высушиваемому
материалу:
конвективные, контактные, радиационные и высокочастотные;
характер движения сушильного агента: с рециркуляцией и без
рециркуляции сушильного агента;
вид теплоносителя: горячий воздух, дымовые газы, пар и
электрический ток высокой частоты;
технологическое назначение и конструктивные отличия: в основу
кладется форма рабочего пространства и характер перемещения в нем
18
материала.
Соответственно
различают
сушилки
камерные,
подовые,
туннельные, конвейерные и пневматические.
Рассмотрим подробнее основные типы сушильных агрегатов и установим
степень их автоматизации в современных условиях производства [20, 50, 83].
Камерные сушилки используются в виде отдельных камер или блоков,
состоящих из нескольких камер (см. рисунок 1.2, а). Размеры камеры: длина от
8 до 13,5 м., ширина от 1.2 до 1.5 м. и высота от 2.3 до 3 м. Процесс сушки
осуществляется циклично, загрузка и выгрузка материала происходит
периодически. Сушки могут работать как на дымовых газах, так и на
подогретом воздухе. Конструктивно камерные сушилки снабжены каналами
принудительной циркуляции сушильного агента (см. рисунок 1.2, б). Несмотря
на
зональную
циркуляцию,
в
сушилках
наблюдается
значительная
неравномерность термовлажностной обработки по объему камеры. Наиболее
низкую температуру и малую интенсивность сушки изделия имеют в середине
камеры [17, 23].
Рисунок 1.2 - Камерная сушилка
Основным
достоинством
камерных
сушилок
является
простота
конструкции и эксплуатации. К недостаткам следует отнести [27].
большие потери времени на загрузку и выгрузку изделий;
большие потери тепла в период выгрузки и загрузки изделий;
19
периодичность
работы,
затрудняющая
создание
поточного
производства изделий;
неравномерность сушки изделий.
В
туннельных
сушилках
изделия,
установленные
на
вагонетках,
специальными толкателями перемещаются по рельсам по длине туннеля с
одной позиции на другую, при этом меняется их влажность. Туннельные
сушилки представляют собой обычно блок от 3 до 8 туннелей. Длина каждого
туннеля от 24 до 38 м, ширина от 1.1 до 1.65 до 1.75 м. (см. рисунок 1.3). Когда
изделия достигают последней позиции туннеля, их влажность соответствует
требуемой и вагонетки выходят из сушилки. Интервал изделий зависит от рода
занятий и интенсивности процесса.
Рисунок 1.3 - Туннельная сушилка
К простому типу относятся противоточные туннельные сушилки с
однократным использованием агента сушки в туннеле. В них сушильный агент
подают с одного конца туннеля и отводят с другого, а вагонетки движутся
навстречу потоку газов [37].
Распределенный подвод сушильного агента, применяемый в некоторых
конструкциях сушилок нельзя признать удачным, так как он приводит к
расслоению потока газов по высоте туннеля и неравномерности сушки
вследствие малых скоростей газов на значительном участке туннеля, что
является основным недостатком туннельных сушилок [67]. К конструктивным
недостаткам
сушилок
некоторых
старых
типов
следует
отнести
и
многопутность туннелей. Отсутствие разделительных стенок между путями
20
приводит к плохой организации потока газов и их нерациональному
использованию [70].
Конвейерные сушилки относятся к сушилкам непрерывного действия –
сырые изделия непрерывно загружаются с одного конца сушилки, а готовые
непрерывно выгружаются с другого конца (см. рисунок 1.4).
Рисунок 1.4 - Конвейерная сушилка
Увеличивающийся удельный вес конвейерных сушилок в сушильных
хозяйствах
керамических предприятий объясняется в
первую очередь
небольшими размерами таких сушилок и непрерывностью их действия, что
позволяет создать поточное производство [30, 39].
К недостаткам конвейерных сушилок следует отнести большие удельные
потери тепла через ограждения, интенсивный воздухообмен между сушилкой и
внешней средой, сравнительно большой удельный расход тепла на сушку [75].
Вне
зависимости
от
того,
какая
сушилка
используется
в
производственном процессе, можно сформулировать общие требования,
которым она должна удовлетворять. Соответственно этими положениями
следует
руководствоваться
при
проектировании
новых
производств
строительной керамики.
Сушилка
должна
обеспечивать
максимальную
скорость
сушки
высушиваемого материала при соблюдении высокого качества, минимальный
расход тепла и электроэнергии на 1 кг испаряемой влаги, равномерность сушки
по всему объему сушилки, обладать, возможно, большей напряжѐнностью
объема по влаге, т.е. количеством испаряемой влаги на 1
объема, легкостью
регулирования параметров сушильного агента. Она должна быть оснащена
механизмами для загрузки, выгрузки и перемещения материала, снабжена
21
приборами теплового контроля и автоматикой и удовлетворять санитарным
нормам.
Одним из основных требований, предъявляемых к сушилкам, является
равномерность сушки изделий по всему объему сушильного пространства.
Она определяется коэффициентом неравномерности сушки Кн, т.е.
отношением конечной влажности двух или нескольких высушенных изделий,
расположенных в различных местах сушилки: изделий с наибольшей конечной
влажностью umax к изделиям с наименьшей конечной влажностью umin, при этом
начальная влажность этих изделий принимается одинаковой: Кн =
.
Коэффициент неравномерности всегда Кн > 1и с увеличением неравномерности
сушки возрастает, при теоретически равномерной сушке Кн = 1.
Несмотря на увеличение доли туннельных и конвейерных сушилок в
общем количестве используемых сушильных агрегатов, камерные сушилки еще
достаточно
широко
используется
на
предприятиях
по
производству
керамических изделий, при этом автоматизация процесса зачастую минимальна
или вообще отсутствует. В связи с этим появляется необходимость в
модернизации существующих производств, основанных на использовании
камерных сушилок в качестве сушильных аппаратов.
Сушильная установка, как было отмечено, состоит из блока камер,
работающих
параллельно.
конструктивные
данные,
При
этом
все
характеризуются
камеры
имеют
близкими
одинаковые
динамическими
свойствами. Выполняя одну и ту же технологическую операцию, камеры
оказываются взаимосвязанными через общие каналы подачи свежего и отвода
отработавшего теплоносителя, что позволяет сделать вывод о принадлежности
объекта к классу многосвязных однотипных систем.
Основным признаком, отличающим однотипные связанные системы от
многосвязных других типов, является то, что основные каналы передачи
воздействий в объекте одинаковы [33, 61, 62, 73].
При организации процесса сушки керамических изделий можно выделить
22
общий источник энергии для всех камер – теплогенератор, вырабатывающий
свежий (горячий) теплоноситель, который нагревается в подающий канал.
Кроме того, имеется еще один общий канал отвода теплоносителя,
предназначенный для удаления холодного влажного воздуха из камер.
Для
построения
математической
модели
достаточно
разобрать
математическое описание для одной камеры и впоследствии применять его для
любой из камер. Особенности описания определяются действующими между
сепаратными системами перекрестными связями, их физической природой. Для
рассматриваемого объекта природа перекрестных связей представляет собой
распределение потока теплоносителя из общего канала по отдельным камерам.
При этом необходим учет их индивидуального потребления в каждый момент
времени.
В каждом агрегате регулируется одна переменная технологического
процесса – количество поступающего в камеру теплоносителя, число агрегатов
может варьироваться, но в общем случае составляет n.
Перемещение любого из регулирующих органов, например органа с
координатой u1, приводит к изменению давления в трубопроводе и как
следствие к изменению расходов g1,g2,…,gn воздушного агента на все камеры.
Таким образом, значение каждого из расходов определяется положением не
только своего, но и соседних шиберов, являющихся регулирующими органами.
Поскольку
процессы
в
пневматической
безынерционны по сравнению с процессами
сети
обычно
практически
в сушильных камерах, то
перераспределение расходов происходит мгновенно вслед за изменениями
положения регулирующего органа [48]. В линейном приближении зависимость
расходов от положений шиберов описывается алгебраическим уравнением
(1.1):
(1.1)
где
– коэффициент передачи j-го регулирующего органа к расходу через k-й
агрегат.
23
Коэффициенты передачи прямых каналов
могут различаться из-за
индивидуальных особенностей характеристик шиберов, а перекрестных связей
– из-за особенностей пневматической сети.
Изменение каждой из регулируемых переменных
происходит под
действием изменения расхода на камеру gk и описывается уравнением (1.2):
(1.2)
где
- возмущение;
– передаточные функции агрегата по каналам
управляющего расхода и возмущения.
Подставив данные уравнения (1.1) в уравнение (1.2), получим уравнение
для объекта-блока камер (1.3):
(1.3)
Поскольку процессы являются безынерционными, то перераспределение
расходов по каналам от каждого из регулирующих органов ко всем
регулируемым переменным различаются только коэффициенты передачи
.
В соответствии с типовой структурой системы автоматизированного
управления [42, 43], для каждой камеры должен быть установлен регулятор с
передаточной
функцией
,
обеспечивающий
организацию
технологического регламента производства, так что уравнение i-го регулятора
имеет вид (1.4):
(1.4)
где
- задание.
В уравнениях (1.3), (1.4) передаточные функции динамических звеньев не
отмечены
индексами,
соответствующих
звеньев
поскольку
динамические
сепаратных
систем,
как
характеристики
уже
отмечалось,
приближенно могут считаться одинаковыми.
Идентичность сепаратных каналов передачи воздействий позволяет
вместо громоздкой полной структурной схемы системы ограничиться
изображением звеньев и связей, соответствующих только одной сепаратной
системе.
24
Такое представление структурной схемы соответствует ее уравнениям
(1.3), (1.4), относящимся к одной сепаратной системе, при использовании
бегущего индекса для распространения уравнений на все n камер сушильной
установки.
Переведем описание системы к более компактной матричной форме
записи уравнений. Для этого введем векторы переменных:
;
;
;
,
где y – вектор регулируемых переменных; x – вектор заданий; u – вектор
управлений; f – вектор возмущений.
А также передаточные матрицы объекта: по основным каналам передачи
воздействий (1.5), по каналам возмущений (1.6) и матрицу регулятора (1.7):
=
=
=
(1.5)
=
(1.6)
=
(1.7)
Одинаковые передаточные функции звеньев системы в передаточных
матрицах вынесены в качестве скалярных сомножителей при матрицах
коэффициентов передачи: некоторой числовой А и единичной Е.
Каждая из передаточных матриц описывает всю совокупность связей
25
между переменными, входящими в состав соответствующей пары векторов;
передаточные функции каналов, связывающих между собой компоненты
векторов, играют роль элементов матриц и занимают в них определенные
места.
Вводные матрицы позволяют перейти в (1.3) и (1.4) к матричной записи:
(1.8)
)
Дальнейший
анализ
(1.9)
системы
будет
направлен
на
определение
передаточной функции сепаратной системы, структуры и параметров ее
регулятора, а также расчет коэффициентов перекрестных связей.
В зависимости от особенностей матрицы коэффициентов передачи,
однотипные связанные системы делятся на [51, 53]:
1.
симметричные;
2.
антисимметричные;
3.
циркулянтные.
Технологический процесс сушки керамических изделий (рисунок 1.5)
конструктивно не обладает симметрией расположения камер как относительно
подающего
теплогенератора, так и системы отводящих вентиляторов,
связанных с общими каналами транспортировки теплоносителя.
Рисунок 1.5 - Схема технологического процесса сушки керамических изделий.
Поэтому рассматриваемую систему можно отнести к циркулянтной
однотипной
связанной
Циркулянтные
матрицы
системе
автоматического
коэффициентов
26
передачи
регулирования
взаимосвязей
(САР).
(1.10)
формируются по следующему правилу: каждая последующая строка повторяет
предыдущую при сдвиге всех элементов на одно место вправо:
(1.10)
Таким
образом,
априорно
сушильная
установка
обладает
всеми
критериями, присущими многосвязным однотипным системам автоматического
регулирования:
1.
уравнения сепаратных САР совпадают;
2.
передаточные функции перекрестных связей, берущих начало и
приложенных в идентичных точках сепаратных САР, совпадают с точностью
до коэффициентов;
3.
матрицы коэффициентов передачи различных звеньев перекрестных
связей системы перестановочны.
Следовательно, для дальнейшего исследования сушильной установки и
синтеза законов управления можно применять методы, предусмотренные
теорией многосвязных однотипных систем автоматического управления.
Выводы
Таким образом, для физических исследований, построенных на основе
сложных математических моделей, характерным является использование
компьютерного математического моделирования.
Проведенный
анализ
возможности
создания
автоматизированной
системы управления технологическим процессом сушки керамических изделий
в камерных сушилках периодического действия устанавливает узкие места в
эффективной реализации указанной системы, которые настоятельно требуют
своего решения. Использование подобной системы на стадии сушки
керамических
изделий
позволит
количественно
определить
теплотехнологические характеристики изделий и среды, а также применить
27
полученную информацию для оптимизации режима работы сушильной камеры
периодического действия с целью уменьшения удельных энергозатрат на сушку
керамических изделий требуемого качества и сокращения времени сушки. Из
этого
вытекает
следующая
цель
научного
исследования:
построить
математическую модель, разработать алгоритмы и систему управления
технологическим процессом сушки керамических изделий, которые позволят
повысить эффективность функционирования камерной сушилки.
2 ПОСТРОЕНИЕ МАТЕМАТИЧЕСКОЙ МОДЕЛИ ОБРАБОТКИ
КЕРАМИЧЕСКИХ ИЗДЕЛИЙ
2.1 Разработка математической модели распределения теплоносителя
в камерной сушилке периодического действия
При проведении анализа объекта было установлено, что построение
общей математической модели для сушки керамических изделий, вначале
нужно иметь описание внутренних теплофизических процессов, которые
протекают в камерной сушилке. Здесь необходимой является модель, которая
отражает как количественные, так и качественные аэродинамические свойства
теплового агента с момента поступления его в камеру до подачи в общий
канал. Такую модель можно быть получить при решении уравнений НавьеСтокса. Известны несколько подходов решения данных уравнений для разных
сред с различными химическими и теплофизическими свойствами агента.
Представим более подходящие, с нашей точки зрения, для решения
поставленной задачи [72].
Используя уравнение Навье-Стокса для несжимаемой жидкости получим
их аналог для сжимаемой жидкости, рассматривая жидкость сжимаемой. Таким
образом, используя данные уравнения для несжимаемой жидкости, у нас будет
частный случай уравнений Навье-Стокса для сжимаемой жидкости. Данные
28
уравнения для несжимаемой жидкости при наличии постоянных свойств и
отсутствии массовых сил и подвода тепла извне, можно записать следующим
образом:
Уравнение непрерывности:
(2.1)
Уравнение движения:
(2.2)
Уравнение энергии:
(2.3)
Здесь V – вектор скорости, м/с;
– плотность среды, кг/м3;
статистическое давление, Па; Т – температура,
о
–
С; µ - динамический
коэффициент вязкости, мПа*с; Ф – диссипативная функция.
Для
данных
уравнений
характерна
смешанная
эллиптически-
параболическая система относительно неизвестных (V, p, T). При этом для
температуры характерным является то, что данный параметр рассматривается
только
в
уравнении
энергии,
поэтому
данное
уравнение
возможно
рассматривать в отдельности [2].
Для
камерной
сушилки характерна
симметричность по
ширине,
следовательно, здесь нами могут быть использованы двумерные уравнения
Навье-Стокса в декартовой системе координат:
Уравнение неразрывности:
(2.4)
Уравнение движения по координате x:
(2.5)
Уравнение движения по координате y:
(2.6)
где u, v – компоненты вектора скорости.
29
Для данных уравнений характерна запись относительно так называемых
примитивных переменных р, µ, v. Один из распространенных методов решения
уравнений
Навье-Стокса
для
несжимаемой
жидкости
примитивные переменные меняются на функцию тока
вышеназванные
и завихренность
.
Определение функции тока и завихренности в декартовых двумерных
координатах следующее:
(2.7)
(2.8)
С использованием новых независимых переменных и комбинируя выше
обозначенные уравнения движения, исключив давление, получаем:
(2.9)
Данное параболическое уравнение в частных производных называют
уравнением переноса завихренности. Чтобы решить его, возможно применить
численные методы решения нелинейного уравнения Бюргерса [2].
Подставив выражения для функций тока, получаем дополнительное
уравнение для независимых переменных
и :
(2.10)
Данное эллиптическое уравнение в частных производных - уравнение
Пуассона. Известен ряд методы его решения [1, 68].
Результат этой замены переменных − возможность разделить смешанную
эллиптически-параболическую
несжимаемой
жидкости
на
систему
одно
уравнений
параболическое
Навье-Стокса
уравнение
и
для
одно
эллиптическое уравнение. На наш взгляд, данный вариант является самым
удачным для решения рассматриваемой задачи, так как цель моделирования получить информацию движения потока теплоносителя внутри камеры в
процессе сушки и рассчитать температурное поле, которое создает этот агент.
Чтобы найти значение температуры нужно решить уравнение энергии. ,
Исходные данные для него − поля скоростей воздушного потока в сушильной
30
камере. Для нахождения значений u и v вектора скорости V оптимально
применить подход, в котором используются завихренность и функция тока. В
данном случае нам не нужно будет вычислять значения давления p,
являющиеся промежуточными для решения нашей задачи с помощью метода с
использованием примитивных переменных и предназначенных для нахождения
компонента вектора скорости в узловых точках сетки. В данном случае
давление - ни управляющий, ни управляемый параметр объекта. Кроме этого,
дифференциальные
уравнения,
содержащие
примитивные
переменные,
являются более сложными и громоздкими, по сравнению с уравнениями в
функциях вихря и тока. Здесь также необходимо учитывать, что для решения
уравнения энергии и нахождения значения температуры, используют как
компоненты вектора скорости, так и значения вспомогательных функций вихря
и тока.
Температуру можно найти из уравнения энергии [68]:
(2.11)
Данное уравнение может также иметь следующий вид [17]:
(2.12)
где Е – число Эккерта, Pe – число Пекле, Re – число Рейнольдса.
При замене в уравнении (2.12) Т на
, Pe на Re и исключении
диссипативного члена Ф, в результате будем иметь уравнение переноса вихря.
Исходя из того, что для газов число Прандтля Pr=0(I) [68], то Pe=0(Re), таким
образом,
для
решения
уравнения
энергии
можно
применить
методы
аналогичные методам для уравнения переноса вихря.
В нашем случае в процессе сушки керамических изделий удобнее всего
применить уравнение энергии в виде уравнения теплопроводности:
(2.13)
31
– коэффициент температуропроводности, м2/с; λ – коэффициент
где
теплопроводности, Вт/м*оС; ρ – плотность среды (воздуха), кг/м3; Сp –
теплоемкость воздуха, Дж/кг*оС.
С целью определения граничных условий нам будет нужно уравнение
потока тепла:
(2.14)
Данная модель - это система дифференциальных уравнений (2.9), (2.10)
и (2.13). Такой подход к анализу процесса распространения агента в сушильной
камере ранее не применяли, а для оценки количественных характеристик агента
использовали
усредненные
значения,
которые
получали
с
датчиков
температуры в камерах.
Решение системы дифференциальных уравнений в частных производных
такого
типа
можно
проводить
различными
методами,
наиболее
распространѐнные из них – метод конечных элементов и метод конечных
разностей.
На наш взгляд, нам
целесообразнее использовать метод конечных
разностей, так как заданная область с достаточно простой конфигурацией и
данный метод является более простым в реализации, но одинаково точным.
Решение полученных уравнений проще и удобнее осуществляется в
стационарном виде. Такой вид уравнений решается с помощью обычного
метода установления. Исходя из того, что на процесс сушки требуется
сравнительно очень много времени (более 100 ч), можно ограничиться
установившимся режимом, образующимся в камере по окончании переходных
процессов, которые занимают небольшой промежуток времени в сравнении со
временем цикла сушки.
Общеизвестно, что для разработки математической модели настоящего
технического объекта или процесса требуются установления ограничений и
допущений
при
моделировании,
которые
32
являются
несущественными
характеристиками. Поставим следующие допущения и ограничения, часть из
которых будет изъята на ряде этапов математического моделирования:
решается плоская задача;
влажность поступающего воздуха не берется во внимание и
рассматривается как постоянная;
изучается
стационарный
процесс
скорости
воздуха
и
распространения температуры;
регенерация
вихрей
в
условиях
отсутствия
возмущающих
воздействий в камере является маловероятной;
взаимодействие частиц и их влияние на поле скоростей не
учитывается;
рассматривается
несжимаемый
газ,
изменение
плотности
отсутствует (p = const);
имеется только конвективный теплообмен;
действие градиента избыточного давления воздуха очень мало;
в камерной сушилке идеальные теплотехнические характеристики;
физические параметры и коэффициенты материала и сушильного
агента в процессе сушки являются постоянными;
физические параметры и коэффициенты по всему объему изделий
считаются однородными и постоянными.
Цель моделирования – изучение процесса распространения воздушных
потоков теплоносителя внутри сушилки, который аналогичен процессу
распространения воздуха в помещении [11]. Важным является получение полей
скоростей
для
расчета
температурного
поля.
Этот
процесс
является
необходимым для внесения правок и определения начальных и граничных
условий при математическом описании сушки керамических изделий по всей
камере.
33
С этой целью мы приводим решение системы уравнений Навье-Стокса и
уравнения неразрывности. При определенных ограничениях и допущениях
данная система будет иметь следующий вид [82]:
(2.15)
Для решения задачи нами будут использованы новые переменные:
завихренность
и функция тока :
(2.16)
Для плоской задачи:
(2.17)
Функция тока определяется из уравнений:
,
(2.18)
Тогда исходная система запишется в следующем виде:
(2.19)
или:
(2.20)
Можно предположить, что в сушилке отсутствуют вихри, т.к. воздух
является маловязкой средой. В камере отсутствуют различные сопротивления
и возмущающие потоки воздуха, поэтому здесь регенерация вихрей является
маловероятной [82], таким образом, движение является ламинарным и
=0.
Исходя из этого, поле скоростей можно описать одним уравнением Лапласа для
функции тока:
(2.21)
34
Данное уравнение можно решить методом конечных разностей. В
качестве плоскости для решения этой задачи принимаем сечение, которое
проходит вдоль длины камеры, так как она превосходит ширину и самым
интересным является распространение теплоносителя в направлении высотадлина. Исходя из этого, значения скорости и температуры по ширине можно
принять
постоянными
и
равными
соответствующим
значениям
для
продольного сечения. Линии тока, которые получены в разрезе сушилки,
можно считать постоянными, при этом картина не изменяется от разреза к
разрезу. Подача воздуха будет осуществляться через отверстия, которое
расположено в полу сушилки, тогда как отвод отработавшего агента будет
проходить через отверстие, которое расположено на потолке сушилки.
Скорость потока на входе представляем как V1, на выходе - V2.
Для определения значения скорости воздуха, который поступает в
камеру, нужно знать производительность подающего вентилятора в площадь
входного отверстия. Производительность подающего вентилятора, как правило,
составляет 7740-13500 м3/ч. Площадь входного отверстия обозначим S, она
равна 0.0314 м2.
Можем найти скорость воздуха, который поступает в камеру:
Vmin = Qmin/2*Sподачи = 2.15/0.0628 = 34.24 [м/с],
Vmax = Qmax/2*Sподачи = 3.75/0.0628 = 59.71 [м/с].
Внутри сушильной камеры необходимо существование баланса воздуха,
то есть Qвх = Qвых:
Q = V*S
Q1 = V1*H1
Qвх = Q1
Qвых = Q2 = V2*H2
V1*H1 = V2*H2
Vвых = V2 = V1*H1/H2
(2.22)
35
Необходимо задать краевые условия для скорости и функции тока. В
нашем случае Vy = 0 на всех стенках, Vx = 0 на стенках, а в проемах Vx = Vвх и
Vx = Vвых соответственно. Функция тока
непрерывна и постоянна:
(2.23)
Можем построить расчетную сетку, на основе габаритных размеров
сушильной камеры (см. рисунок 2.1). Задаем число узлов сетки по осям,
учитывая то, что на проемы должно приходиться по 2 и более узлов. По оси
OX: N=40 точек, шаг hx = H/M = 0.027; по оси OY: M = 26 точек, шаг hy = D,N =
0.027. Тогда текущие координаты x и y будут пересчитываться по формулам:
xi = i*hx,i = 0,1,…,M; yj = j*hy, j = 0,1,…,N.
Значения H1, H2 (2.22) находятся в узлах.
Рисунок 2.1 – Размеры элементов сушильной камеры и ориентация осей
Теперь можно описать сеточные функции:
36
(2.24)
(2.25)
Получаем сеточные аналоги для уравнений:
(2.26)
для внутренних точек сетки iЄ[1;N-1],jЄ[1;M-1].
На границах условия в сеточных функциях запишутся так:
1)
левая стенка: j=0, 0≤ i ≤ M,
2)
правая стенка: j=N, 0≤ i ≤ M,
3)
потолок вне проема: i=0,
в проеме: i=0,
4)
;
пол: вне проема: i = M,
,
в проеме i = 0,
Начальные условия для скорости внутри сетки
функции тока:
или
.
Признак окончания пересчета значений функции тока:
Чтобы не хранить весь массив значений
следует считать не
; для
< ξ.
, на новом временном слое
, а другую функцию, например, Q, сравнивать результат
< ξ, а потом присваивать значение Q:
= Q.
После достижения точности подсчитываем значения скоростей во
внутренних точках:
(2.27)
37
В результате вычислений получаем значения скоростей и функций тока в
узлах сетки, которые необходимы для дальнейших расчетов и получения
температурного поля.
Уравнение теплопроводности первоначально имеет вид (2.13) или [44]:
(2.28)
где T(x, y, t) – значение температуры внутри расчетной области, зависящее от
геометрических координат и времени.
Так как рассматривается стационарный режим, то
поэтому в
стационарном случае уравнение теплопроводности примет вид:
(2.29)
При этом значение функции T теперь будет зависеть только от
геометрических координат: T(x, y).
Сеточная функция температуры:
T(x, y) T(
(2.30)
Сеточный аналог уравнения теплопроводности:
(2.31)
где h = hx = hy = 0.027 м.
Начальные условия для функции температуры внутри сетки и на
границах
= 20 .
Начальными условиями во входных отверстиях являются:
, где
– максимальное значение температуры нагретого воздуха, подаваемого в
камеру.
Признак окончания пересчета значений функции температуры идентичен
признаку перерасчета значений функции тока.
Граничные условия температуры определяю используя уравнение потока
тепла (2.14) и теплопроводности (2.29). Для твердой стенки поток тепла равен
38
нулю: П = 0. Для потока тепла, идущего от стены вовнутрь камеры или наружу:
П = ρ Cρ*T Vn, если Vn <0 или П = ρ Cρ*Tmax Vn, если Vn >0.
Рассмотрим вывод граничных условий для потолка камеры: i =0, j =0..N.
Для потолка уравнение теплопроводности примет вид:
i,j
+
i,j
-
i,j =
0,
(2.32)
так как на потолке Vy = 0.
Для потолка получается:
Ti+1,j = Ti,j +
i,j
+
i,j ,
(2.33)
откуда
i,j
=
Ti+1,j -
Ti,j -
i,j
.
(2.34)
Подставив выражение (2.34) в уравнение теплопроводности (2.32),
получим:
Ti+1,j -
Ti,j -
i,j
+
i,j
-
i,j
=0
(2.35)
Из конечно-разностного представления [51] следует, что:
i,j
=
(2.36)
Подставив (2.36) в (2.35) и сгруппировав подобные челны, получим
следующее выражение:
-
Ti,j +
Ti+1,j +
Ti,j-1 +
Ti,j+1 – (
i,j
= 0.
(2.37)
Для потолка при потоке тепла, который поступает внутрь камеры,
получаем выражение П = ρ Cρ*Imax Vn, при Vn <0.
Уравнение потока тепла для потолка будет следующим [44]:
= -λ
+ Vx ρ CρT.
(2.38)
Подставив выражение для потока тепла (2.38) в уравнение для потолка
(2.33) и выразив
i,j
, получим следующее выражение для
= Ti,j*Vx - Tmax*Vx.
:
(2.39)
39
Подставив выражение (2.39) для
в уравнение (2.37), после некоторых
простых преобразований и выражения Ti,j получим выражение для граничных
условий на потолке в проеме:
Ti,j =
(2.40)
Или в координатах сетки:
T0,j =
.
(2.41)
В случае отсутствия потока тепла на потолке, то есть на твердой
поверхности Пx = 0. Тогда выражение (2.39) примет вид:
i,j
= Ti,j*Vx
(2.42)
Выражение граничных условий на потолке примет вид:
Ti,j =
(2.43)
Или в координатах сетки:
T0,j =
(2.44)
Аналогично можно вывести граничные условия для стен и пола.
Граничные условия для правой стены:
Ti,N =
(2.45)
Граничные условия для левой стены:
Ti,1 =
(2.46)
Граничные условия для пола вне проема:
TM,j =
(2.47)
Граничные условия для пола в проеме:
TM,j =
(2.48)
40
Чтобы представить результаты и показать анализ составленной нами
математической модели, необходимо ввести нужные параметры сушильной
установки: высота камеры H = 0.7 м, длина D = 1 м, длины проемов H1 = H2 =
0.2 м (см. рисунок 2.1).
Моделирование будем производить в системе инженерных и научных
расчетов SolidWorks и MATLAB (приложение А). Для программирования
математической модели данные значения можно представить в виде матриц
размерностью Mˣ N, элементы которых V(i, j) и T(i, j) - значения скорости и
температуры в определенной точке. С помощью MATLAB можно наглядно
представить результаты моделирования в виде линий уровня и трехмерного
изображения результирующих матриц в заданной системе координат.
Анализ
результатов
моделирования
возможно
посмотреть
по
получившимся полям при разных начальных параметрах (приложение Б). Для
примера можно посмотреть некоторые из рассчитанных зависимостей (см.
рисунки 2.3 – 2.5), которые соответствуют принятой расчетной сетке и
ориентации координатных осей (см. рисунок 2.2).
Рисунок 2.2 – Изображение линий уровней скорости в условиях максимальной
производительности подающего вентилятора
41
В случае увеличения скорости входного потока показатели скоростей
внутри камеры увеличиваются, при этом перемещение воздуха от входа к
выходу происходит более интенсивно. Около стенок показатели скорости
интенсивно уменьшается, а на стенках – равны нулю, что описывается
граничными условиями. На линии, которая соответствует
показатели скорости также равны нулю, что обуславливается
симметрии,
взаимным
гашением потоков воздуха, которые направлены друг другу от входного
отверстия. Самый интенсивный поток воздуха и наибольшие показатели
скорости регистрируются у входного и выходного отверстий.
Поле скоростей – вспомогательное поле, которое является необходимым
для расчета температурного поля, представляющего самый высокий интерес.
В камере наименьшая температура воздуха около входного отверстия и в
направлении самого интенсивного распространения воздуха, как и в случае со
значениями
скорости.
Воздух
и
поток
тепла
интенсивнее
всего
распространяются вдоль оси х от отверстия входа к потолку камеры.
Показатели, характеризующие температуру воздуха, уменьшаются, воздух
нагревается
приближаясь
к
центру
камеры,
что
обуславливается
незначительной длиной сушилки. Остывание наблюдается по оси x, y пола
температура незначительно больше, чем у потолка [64].
С увеличением значений температуры воздуха, который поступает в
камеру, средняя температура в камере увеличивается, при этом для поля
характерен большой диапазон значений. Одновременно наблюдаем охлаждение
воздуха при движении к центру камеры и в направлении оси x (см. рисунки 2.3,
2.4).
42
Рисунок 2.3 – Изображения линий уровня температуры потоков при наименьшей
производительности подающего вентилятора и температуре входящего воздуха 500С
Рисунок 2.4 - Температурное поле при наименьшей производительности подающего
вентилятора и температуре входящего воздуха 500С
Средние значения температур в сушильной камере, которая не заполнена
изделиями (сопротивлениями), при разных значениях температуры входящего
воздуха представлены в таблице 2.1.
43
Таблица 2.1. Показатели температуры в пустой камере
Температура поступающего
Средняя температура в камере
воздуха, оС
о
С
30
21,85
50
28,11
100
71,44
Таким образом, используя математическую модель, мы получили поля
скоростей
и
температур
внутри
сушильной
камеры
без
внутренних
сопротивлений, описывающих поведение теплоносителя в сушилке. С их
помощью представляется возможным получить нужные значения параметров
процесса сушки в любой точке расчетной сетки.
С целью определения количественных показателей высушенных изделий
– внутреннего напряжения и остаточного влагосодержания, нужно знать
значения температуры на гранях образцов в течение всего процесса сушки.
Если будут внесены сопротивления высушиваемых образцов внутрь сушилки –
поле скоростей и температур значительно изменится. Рассчитывая показатели
высушенных образцов,
нужно брать параметры температуры и скорости с
учетом сопротивлений внутри камеры. С этой целью нужно незначительно
поменять методику расчета, в которой не пересчитываются значения функции
тока и температуры в точках, где имеются сопротивления. На границах, где
находятся блоки с высушиваемыми образцами, также существуют граничные
условия.
Во внутренних точках блоков задаются следующие условия:
44
На границах блоков условия следующего вида:
Значение функций температур
на границе блоков рассчитывают по
формулам граничных условий для температуры. Здесь верхнюю грань блока
рассматривают как пол камеры, нижнюю – как потолок, левую грань – как
правую стенку, а правую – как левую стенку.
Необходимо дополнить описание сушильной камеры значениями и
показателями внутренних сопротивлений для непосредственных изделий,
которые подвергаются термовлажностной обработке, и произвести их
относительно к сеточной модели. В сушильной камере имеется цилиндр с
отверстиями для хорошего прохождения воздуха. В саму камеру помещают
керамическую плитку в семь рядов, каждый из которых является полкой, на
которую кладут 15 плиток. При полной загрузке камера насчитывает 105 шт.
керамической плитки. Расстояние между отдельными плитками равно 20 мм, от
крайних образцов до границ камеры примерно 120 мм. От верхних плиток до
потолка и от нижней полки до пола около 100 мм.
Для нахождения поля ряда скоростей и температур в камерной сушилке,
учитывая сопротивления (образцы), нужно расположить их в расчетной сетке.
Для шага сетки характерно расстояние, составляющее 12 мм. Это находится в
соответствии с условиями сходимости и устойчивости решений дифуравнений
[55]. При этом значение шага сетки непосредственно оказывает влияние на
показатель времени для решения задачи на ЭВМ. Исходя из этого, на расчетной
сетке сопротивления будут иными, по сравнению с сопротивлениями, которые
происходят в реальной сушильной камере. Это нужно учитывать для расчета
значений скоростей и температур между испытуемыми объектами и полками.
С целью учета изделий, которые располагаются внутри камеры, вносим
корректировки в программу расчета теплотехнологических параметров (см.
45
приложение
В).
Некоторые
результаты
уточненного
математического
моделирования представлены на рисунке 2.5.
Для внутреннего сопротивления характерно препятствие распространения
воздуха в камерной сушилке. Самым интенсивным распространением воздуха
является распространение, которое можно наблюдать у входного и выходного
отверстия, а также между стенами и изделиями, между сопротивлениями и у
полом.
Рисунок 2.5 – Изображения линий, характеризующих уровень температуры потоков при
наименьшей производительности подающего вентилятора и температуре входящего воздуха
1000С при наличии внутренних сопротивлений
Значения температурного поля изменяются, по сравнению с теми
значениями, которые рассчитаны без учета внутренних сопротивлений.
Вследствие изменения полей скорости и изменения самих сопротивлений
распространение тепла по камере происходит менее интенсивно. Наблюдаем
более интенсивным остывание по направлению оси х и к центру камеры, по
сравнению с показателями, которые регистрируем в условиях отсутствия
46
внутренних сопротивлений. Средняя температура в камере снижается (таблица
2.2).
Таблица 2.2. Показатели температуры заполненной камеры
Температура поступающего
Средняя температура в камере
воздуха, оС
о
С
30
19,13
50
25,92
100
61
Используя
математическую
модель,
которая
показывает,
как
распространяется теплоноситель по камерной сушилке, можем определить ряд
значений температуры на гранях экспериментальных образцов, чтобы
рассчитать количественные характеристики продукции. Это следующая часть
построения
математической
модели
термовлажностной
обработки
керамической плитки.
2.2 Модификация математической
обработки керамической плитки
модели
термовлажностной
Движение влаги внутри объекта в процессе сушки происходит под
воздействием разных сил, которые обуславливаются разностью концентраций
влаги в разных точках тела, градиентом температур и избыточного давления
пара в изделии. Поверхность изделия испаряется, влага удаляется в
окружающую среду. Для внутреннего массообмена характерна неразрывная
связь с внешним [56].
Механизм и интенсивность движения влаги внутри изделий находится в
прямой зависимости от ряда факторов. К ним относятся: форма связи влаги с
исследуемым объектом, строение материалов, его температуры, от уровня
влажности, пористости и др. параметров. Перемещение влаги в объекте может
происходить как в виде жидкости, так и в виде пара.
47
Закон переноса влаги в материалах можно выразить уравнением [19, 32]:
(2.49)
где m – плотность потока или удельный поток влаги, т.е. количество влаги,
переместившейся в единицу времени через единицу поверхности внутри тела,
U – влагосодержание изделия; Т – температура изделия, оС; P –
кг/с*
давление водяного пара, Па; К – коэффициент потенциаловлагопроводности,
учитывающий перенос пара и жидкости, м2/с;
– объемный вес абсолютно
сухого материала, кг/м3; δ – термоградиентный коэффициент; D – коэффициент
молярного переноса пара.
Первому слагаемому уравнения (2.49) соответствует перемещение влаги в
виде пара и жидкости в условиях влажностного градиента. Для переноса влаги
в виде жидкости характерно присутствие диффузионно-осмотических и
капиллярных сил. Согласно
А.В. Лыкову, данное явление является
влагопроводностью. Соответственно, движение влаги в виде пара называют
паропроводностью.
Второму слагаемому уравнения (2.49) соответствует перемещение влаги
под действием температурного фактора. Влага перемещается из места более
высокой
температуры
к
месту
более
низкой
температуры.
Данная
закономерность перемещения влаги в глиняном объекте по направлению
теплового потока впервые была исследована А.В. Лыковым, он назвал еѐ
термовлагопроводностью. Данное явление наблюдал и использовал
при
выполнении исследований В.Н.Зимин, занимающийся сушкой стекловаренных
шамотных горшков [3].
Движение влаги в виде жидкости по направлению потока тепла
осуществляется как термодиффузией, так и уменьшением поверхностного
натяжения жидкости в капиллярах при увеличении температуры и благодаря
пузырькам воздуха, которые находятся внутри капилляров.
Третьим слагаемым уравнения (2.49) является перенос влаги в материале,
который обусловлен градиентом избыточного давления пара. Данный фактор
48
является решающим для переноса влаги при глубинном нагреве объекта
исследования от 100оС до 150оС и выше.
Получение поля температур в исследуемом объекте может быть
достигнуто, благодаря использованию уравнения Фурье для теплопередачи
[60]:
(2.50)
где Т – температура, оС; х – текущая координата в направлении потока, м; с –
теплоемкость тела, Дж/кг*оС; λ – коэффициент теплопроводности, Вт/м*оС; t –
время, с; γ – удельный вес, кг/м3.
Для начальных и граничных условий уравнения теплопроводности
характерен следующий вид [81]:
(2.51)
где
- начальная температура тела, оС; d – толщина материала, м;
коэффициент теплоотдачи, Вт/м*оС;
–
– температура воздуха, оС.
Чтобы получить поля влагосодержания в объекте, можно использовать
следующее уравнение влагопроводности [43]:
(2.52)
где U – влагосодержание изделия; T – температура изделия, оС; х – текущая
координата в направлении потока, м; t – время, с; К – коэффициент
потенциалопроводности, учитывающий перенос пара и жидкости, м 2/с; δ –
термоградиентный коэффициент; γ – объемный вес абсолютно сухого
материала, кг/м3.
Уравнению
влагопроводности
начальных
соответствует вид [43]:
49
и
граничных
условий
(2.53)
где
– начальное влагосодержание тела;
Вт/м*оС;
– коэффициент влагоотдачи,
– коэффициент влагопроводности, Вт/м*оС; UВ – влагосодержание
воздуха.
На сегодня разработанными являются одномерные модели произведения
процессов сушки, которые не способны обеспечивать достаточную точность
самого
процесса.
Благодаря
им,
возможно
получить
только
общее
представление о физико-химических изменениях, которые происходят в период
термовлажностной обработки. Переходя к двумерной модели, можно будет
прийти к более качественным результатам процесса сушки.
Для
математической
модели
использовали
ниже
обозначенные
допущения и ограничения:
используется только конвективный теплообмен;
действие градиента избыточного давления пара незначительно;
для камерной сушилки характерны идеальные гидродинамические и
теплотехнические характеристики;
распределение воздуха по всей камере определяется результатами
математической модели распределения теплоносителя по камерной сушилке;
для физических параметров и коэффициентов материала характерно
постоянство;
для физических параметров и коэффициентов по всему объему
плитки характерна однородность и постоянство.
Согласно уравнениям (2.49) – (2.53), можно получить математическое
описание процесса сушки (2.54).
Данной модели соответствует система дифференциальных уравнений
второго порядка в частных производных с граничными условиями Дирихле [13,
50
47, 53], еѐ решение можно выполнить с помощью двух методов: конечных
элементов и конечных разностей.
Мы воспользуемся методом конечных разностей, исходя из того, что для
заданной
области
характерна
простая
конфигурация,
аналогичная
конфигурации, которая используется для расчета полей распределения
теплоносителя внутри сушильной камеры.
(2.54)
Решению системы уравнений (2.54) способствует расчетная схема (см.
рисунок 2.6), предполагающая вычисление значений функций F(I, J, t), через
значения F(I-1, J, t-1), F(I, J, t-1), F(I+1, J, t-1), F(I, J-l, t-1), F(I, J+l, t-1).
51
Рисунок 2.6 – Схема для расчета
В работе было использовано программное средство реализации
−
система компьютерной математики MATLAB, являющаяся лидером в области
матричных операций и имеющая множество программных приложений,
которые позволяют провести расчеты разные по сложности и длительности [52,
63], а также имеющая широкие возможности графического отображения
полученных результатов.
Последовательным решением уравнений (2.54) методом сеток, можно
получить поле влагосодержания и температуры (см. рисунки 2.7 – 2.8) в
изделии в процессе сушки.
Рисунок 2.7 – Температурный режим сушки
52
Рисунок 2.8 – Влажностный режим сушки
Чтобы учесть неравномерность распределения температуры по камерной
сушилке и осуществить моделирование самого
процесса, необходимо
распределить результаты модели распространения теплоносителя по камере на
модель сушки одного образца нужное число раз. Первую модель используем
для нахождения средних значений температуры окружающего воздуха в
отдельности для каждого объекта исследования. Далее для остальных объектов
исследования находим показатели характеристик процесса сушки.
В камерной сушилке можно выделять образцы, сушка которых
происходит
при
максимальной
и
минимальной
температурах,
Они
соответствуют пороговым значениям показателей готовой продукции. С целью
экономии ресурсов, возможно использовать данные, при которых получаются
образцы, высушиваемые при данных значениях температуры. В итоге получим
пороговые значения остаточного влагосодержания и внутренних напряжений.
Параметры данных характеристик для других образцов можем расположить в
пределах этих пороговых значений. Исходя из этого, качество продукции будет
соответствовать
пороговым
значениям
остаточных
влагосодержаний
и
внутреннего напряжения. Моделирование можно произвести для двух
пограничных по температуре объектов. Исходя из того, что для сушильной
камеры характерна симметричность относительно центрального отводящего
53
отверстия, мы можем рассмотреть только половину сушилки и объекты,
которые здесь находятся, чтобы сузить рассматриваемые блоки.
Максимальная температура сушки по периметру соответствует образцу,
находящемуся в первом ряду крайнего верхнего слоя, т.е. у подводящего тепло
отверстия.
Минимальная
температура
сушки
характерна
для
образца,
находящегося практически напротив отводящего тепло в отверстии. Для
нижних образцов температурный режим практически идентичен, в связи с этим
для определенности используем данные, характерные для нижнего объекта.
Результаты моделирования процесса для 2-х образцов в условиях
используемого в производстве технологического режима получены следующие
данные (см. рисунки 2.9 – 2.12):
Рисунок 2.9 - Поля температур при окончании процессов сушки
а) поле верхнего образца; б) поле нижнего образца
54
Рисунок 2.10 - Поля остаточных влагосодержаний на 120 часов
а) поле верхнего образца; б) поле нижнего образца
Рисунок 2.11 – Показатель средней температуры изделий в процессе сушки
Рисунок 2.12 – Среднее содержание влаги в изделиях при сушке
55
Неравномерное распределение влаги по объекту при сушке способствует
неравномерной усадке разных его слоев. Высушенные наружные слои
уменьшаются в размерах, а внутренние препятствуют усадке поверхностных
слоев. При этом происходит их усадка, которая способствует поверхностному
растягиванию. При этом внутри тела наблюдаем сжатия, способные вызвать
деформации, коробления изделия, а также трещины. В связи с этим каждому
материалу и изделию необходимо установить свой рациональный режим
сушки, включающую определенную скорость, температуру и параметры,
обеспечивающие качественную сушку [8, 28].
Объекты, которые высушиваются в различных местах камерной сушилки,
обладают разными значениями средних остаточных емкостей влаги, это
свидетельствует о неравномерном распределении тепла по сушильной камере.
Крайнему
правому
верхнему
образцу
соответствует
остаточное
влагосодержание иср = 5.0316 %, а нижнему центральному - иср = 5.0220 %.
Для большинства объектов при сушке характерно уменьшение в
размерах, т.е. усадка. Для глины и глиняных изделий, а также для других
материалов усадка осуществляется во время постоянной скорости сушки и
окончание ее регистрируется по достижению ими критической влажности [12].
Для процесса усадки глины характерен следующий вид: имеется
предположение, что влага, которая заполняет поры в исследуемом образце,
способствует образованию на границе «изделие - воздух» вогнутых менисков
во внутренних капиллярах. При сушке влага испаряется, а ее поверхностное
натяжение в капиллярах повышается, приводя к сужению изделие. Параметры
величин изделий уменьшаются до тех пор, пока глиняные частицы не будут
взаимно соприкасаться, и между ними будет происходить трение. При
достижении
величины
трения,
которая
будет
превосходить
силу
поверхностного натяжения влаги в капиллярах, уменьшение размеров изделия
будет прекращаться. Это момент, когда средняя влажность изделия достигнет
значения критической влажности. Далее наблюдаем прекращение усадки, при
этом дальнейшее удаление влаги будет вести только к увеличению пористости
56
изделий, т.к. капилляры будут освобождаться от влаги, которая в них находится
[21, 30].
На
основе
полученного
в
результате
моделирования
поля
влагосодержания можем рассчитать напряжения, которые возникают в изделии
в процессе сушки (см. рисунки 2.13-2.15), с этой целью можно использовать
следующие соотношения [39, 81]:
(2.55)
где
- линейный размер тела при влажности U, м;
сухой вес);
- влажность тела (на
- линейный размер абсолютно сухого материала (когда усадка
происходит на протяжении всего периода сушки или в период падающей
скорости), м; ls - средний размер, м; ψ - коэффициент чувствительности сырья;
Е – модуль упругости, кг/м2; σ – напряжения в материале, Па.
Рисунок 2.13 – Параметры остаточных напряжений в изделиях
а) верхний образец; б) нижний образец
57
Рисунок 2.14 – Показатели максимальных напряжений в изделиях
Рисунок 2.15 – Показатели минимальных напряжений в изделиях
Согласно А.Д. Цепину [81], основная причина, влияющая на размер,
прочность,
увеличение
пористости,
качества
и
свойств
изделия,
это
пластическая деформация объекта исследования, происходящая при начальной
стадии сушки, то есть в самые первые моменты, когда происходит усадка
объекта.
Повышению
показателей
исследуемого
58
объекта
соответствует
небольшая скорость сушки с удалением от 4 до 6 % начальной влажности,
после чего следует осуществить процесс с большими скоростями. Качество
изделий должно соответствовать ГОСТ.
Выводы
Использование методов математической физики для решения уравнений
Навье-Стокса, которые описывают состояния объектов, и методов решения
уравнений тепло- и влагопроводности, которые составляют математическую
основу
термовлажностной
обработки
материала,
а
также
методов
математического моделирования в ходе разработки математической модели
теплового объекта, теории однотипных многосвязных систем управления,
современных методов теории однотипных многосвязных систем управления
способствовало решению следующих задач:
проведение анализа функционирования теплового объекта и
разработка его математической модели;
модификация
математической
модели
для
процесса
термовлажностной обработки отдельного образца;
совмещение модели обработки изделия с моделью теплового
объекта в условиях соблюдения технологического режима сушки.
Для решения таких задач как:
оптимизация технологического регламента сушки керамических
изделий для многосвязного теплового объекта управления и
разработка
технологических
алгоритма
регламентов
управления
сушки
керамических
для
оптимальных
изделий
имеется
необходимость разработать и рассчитать оптимальный по быстродействию
режим сушки керамических изделий.
59
3
РАЗРАБОТКА
И
РАСЧЕТ
ОПТИМАЛЬНОГО
ПО
БЫСТРОДЕЙСТВИЮ РЕЖИМА СУШКИ КЕРАМИЧЕСКИХ ИЗДЕЛИЙ
Для заводов, производящих керамическую плитку, характерен выпуск
изделий определенной номенклатуры. Каждому виду продукции соответствуют
стандартизированные габаритные размеры, поэтому для процесса удаления
влаги из керамических изделий характерна индивидуальность. Наряду с этим,
габаритные размеры образцов влияют на загрузку сушильной камеры. Данное
обстоятельство ведет к различной плотности расположений изделий в камерной
сушилке и, как результат, к периодическому изменению аэродинамических
характеристик установки.
Исходя их того, что для процесса производства керамических изделий
характерна круглогодичность, проведение сушки керамических изделий нужно
осуществлять с учетом климатических условий добычи и предварительной
обработки исходного сырья, а также его физико-химических характеристик.
Качеству глины соответствует место ее добывания. Исходя из этого, при
обработке
глины
показатели
коэффициентов
теплоемкости,
тепло-
и
влагопроводности сформированных образцов, которые подвергаются сушке,
будут меняться. Помимо этого, для климатических условий, в которых
добывают
сырье,
необходимо
определить
начальную
температуру
и
влагосодержание изделий [14, 58, 80, 87, 88].
Исходя из этого, для технологического процесса сушки керамических
изделий существует
широкий диапазон изменения входных параметров.
Однако регламенту термовлажностной обработки изделий на действующих
предприятиях, как правило, присуща неизменность. Результатом является
появление бракованных изделий, число которых составляет до 30% от общего
объема выпуска [36, 38, 54, 86].
Для снижения доли бракованных изделий и повышения эффективности
использования
камерных
сушилок
нужно
определить
оптимальные
технологические регламенты, которые могут действовать в широком диапазоне
60
изменений входных параметров. Решение этой задачи вероятно только в
условиях применения метода математического моделирования.
С целью нахождения режима
оптимального управления для сушки
керамических изделий нужно рассмотреть технологический процесс его
термовлажностной обработки в камерной сушилке.
В ходе расчета необходимо устанавливать ограничения, которые связаны
непосредственно с процессом обработки, а также с показателями характеристик
изделий с целью корректной постановки задачи оптимального управления.
В качестве состояния Хi рассмотрим поля температур Ti(t,xi,yi) и
влагосодержаний
Ui(t,xi,yi)
i-го
изделия
термовлажностной обработке – Хi
(i=
),
которое
подвергается
= f(Ti(t,xi,yi), Ui(t,xi,yi)). Фазовым
координатам присуще описание непосредственных высушиваемых изделий.
Каждому виду продукции prodi, соответствуют определенные габаритные
размеры (lxi, lyi), находящиеся в рассматриваемом двумерном пространстве,
которые установлены требованиями ГОСТа 6141. Ограничения на координаты
изделия можно представить в виде (3.1):
(3.1)
Вместе с тем порядковому номеру изделий i соответствует его
пространственное местоположение в камере (x, y), которое ограничено
размерами сушилки (3.2):
(3.2)
Т.о., ограничениям на координаты, согласно (3.1) и (3.2), соответствуют
системе неравенств (3.3), показывающей широкий диапазон изменений
фазовых координат, существенно зависящий от входных параметров PROD:
(3.3)
61
Для анализа камерной сушилки как объекта управления характерны 3
группы ограничений:
1.
Ограничения
производительности
оборудования,
которое
обеспечивает подготовку теплоносителя и его подачу в камерную сушилку
(3.4):
Ограничения производительности теплогенератора QТГ = 70000
м3/ч;
Ограничения производительности вентилятора, который подает в
камеру теплоноситель Qmin
Q ≤ Qmax, где Qmin = 7400 м3/ч, Qmax = 13500 м3/ч;
Ограничения номинальной производительности вентилятора Qном =
10000 м3/ч;
Ограничения скорости поступающего в камеру теплоносителя Vmin
≤ V ≤ Vmax, где Vmin = 34.24 м/с, Vmax = 59.71 м/с,
(3.4)
2.
Ограничения
параметров используемого теплоносителя в ходе
термовлажностной обработки объектов (3.5):
Ограничения
Ограничения
подается в общий канал
параметров
начальной температуры объектов
параметров температуры теплоносителя, который
= 150 ;
Ограничения параметров температуры i – го образца (
) в
ходе сушки
(3.5)
3.
Ограничения, которые обуславливают качество образцов в конце
стадии сушки (3.6):
Ограничения остаточного влагосодержания i – го изделия (
где
= 0.04,
= 0.06,
62
)
(3.6)
Рассмотрев выражения (3.4) – (3.6), получаем векторы равенств и
неравенств, которые характеризуют ограничения параметров технологического
процесса керамической сушки объектов в сушильных камерах (3.7):
(3.7)
Критерий
поиска
оптимальных
технологических
режимов
сушки
керамических изделий в сушильной камере с периодическим действием можно
представить следующим образом [4]:
(3.8)
где
– весовой коэффициент.
Кроме быстрого осуществления процесса, нужно постараться увеличить
качество готовой продукции за счет снижения показателя, характеризующего
остаточное влагосодержание изделий и удержание напряжений, которые
возникают в объекте при термовлажностной обработке, на допустимом уровне
[70].
Ограничения, которые обуславливают уровень качества изделий в ходе и
по окончании стадии сушки (3.9) включают:
Ограничения остаточного влагосодержания i – го изделия (
где
= 0.04,
)
= 0.06,
Ограничения допустимых напряжений растяжений в плитке в
процессе сушки
МПа;
Ограничения допустимых напряжений сжатий в плитке в течение
процесса сушки
МПа;
Ограничения уровня допустимых текущих напряжений i – го
изделия (
)
(3.9)
63
Ограничения, которые характеризуют качество готовой продукции с
соблюдением требований ГОСТа, можно представить в следующем виде:
(3.10)
где
– весовой коэффициент.
Для оптимального управления технологического режима сушки нужно
использовать одновременно два критерия (3.8) и (3.10), так как для них
характерна
противоречивость
Целью
[79].
является
максимальное
сокращение времени для сушки керамической плитки. Это можно достичь при
мгновенном увеличении температуры сушильной камеры и далее при условии
еѐ поддержания на максимальном уровне. Наряду с этим, в изделиях во время
прогрева вероятнее всего будут возникать разрушающие напряжения и в
результате сушки продукция будет бракованной. Целью критерия
является
обеспечение максимального качества изделий, которое может быть обеспечено
благодаря «мягкому» постепенному увеличению температуры в камере, что
будет способствовать увеличению процесса сушки. Т.о., покажем функционал
качества (3.11) оптимального управления технологическим режимом сушки
керамической плитки, который представляет свертку цели и ограничений по
формуле Беллмана-Заде «решение есть слияние цели и ограничений»:
(3.11)
где
– нормирующий коэффициент.
Также необходимо осуществить свертку ограничений (3.7) и (3.9):
(3.12)
Чтобы вычислить собственные значения критериев
,
и задания
весовых коэффициентов при учете технологических ограничений (3.12), можно
воспользоваться подготовленной моделью распределения теплоносителя в
камерной сушилке (2.9), (2.10), (2.13) и двумерной моделью, которая
показывает
обеспечивание
процесса
керамической плитки (2.54).
64
термовлажностной
обработки
Т.о., мы можем представить формулировку задачи оптимального
управления
технологическим
режимом
сушки
керамической
плитки:
необходимо осуществить управление, которое обеспечит минимум критерия
(3.11) используя связи в виде уравнений математической модели (2.54) и
ограничений (3.3), (3.12).
Сложностью
при
решении
названной
задачи
является
то,
что
математическая модель термовлажностной обработки керамической плитки
представляется как нелинейные уравнения в частных производных, критерии
оптимизаций являются составными. В этой связи является невозможным
использование
таких
приемов
как
принцип
максимума
Портнягина,
динамического программирования, методов вариационного исчисления [49, 59,
66].
При этом выполненные исследования показывают, что управление можно
представить как зависимость температуры обработки от времени (см. рисунок
3.1), характеризующуюся 2-мя моментами переключения:
переключение продолжительности периода подготовки плитки в
начале процесса сушки (
;
переключение длительности временного интервала при увеличении
температуры теплоносителя от начальной до максимальной (
.
Таким образом, поиск оптимального управления можно осуществлять
методом математического программирования [66].
65
Рисунок 3.1 – Показатели переменных параметров технологического процесса при
оптимизации режима сушки
Используемые нами параметры
и
технологического режима можно
считать зависимыми, так как в течение времени
в камерной сушилке
происходит «естественная» сушка, способствующая начальной выдержке
изделия с целью выравниваний полей содержания влаги внутри объекта. Затем
идет начало процесса искусственной сушки материала, при котором
регистрируется появление значительных градиентов температуры и содержания
влаги в изделии и, как результат, можем наблюдать опасные для целостности
плитки напряжения. Скорости увеличения температуры в керамической
сушилке принадлежит очень большая роль. Для нее характерна взаимосвязь с
величиной времени предварительной выдержки изделия [26]. Величина
параметра
определяет время выхода процесса на постоянную скорость
сушки, и, следовательно, величину градиентов полей температуры и
содержания влаги в материале.
Качество изделий и его прочность определяются интенсивностью
процесса термовлажностной обработки в начале сушки [79]. Исходя из этого,
целесообразным будет проведение внешнего цикла алгоритма оптимизации по
параметру
, а вложенного цикла – по параметру
Для параметра
параметра
.
характерно изменение в пределах от 2 до 20 ч, для
– от 1 до 10 ч. Выбору данных граничных значений
и
способствует то, что решение поставленной задачи оптимизации нужно
проводить
в
условиях
ужесточения
режима
сушки
относительно
существующего. Существующими параметрами являются следующие:
= 24 ч,
= 10 ч. Расчетные модели показали следующие результаты (см. рисунки 3.2. –
3.3):
66
а)
б)
Рисунок 3.2 – Показатели зависимости от значений параметров технологического режима
остаточного влагосодержания
а) верхнее изделие; б) нижнее изделие
67
Рисунок 3.3 - Показатели зависимости от значений параметров технологического режима
напряжений сжатия для верхнего и нижнего изделий
Напряжения растяжения на нижнем образце принимают более высокие
значения в зависимости от значений параметров
и
(при их малых
значениях), по сравнению с напряжениями растяжения для верхнем изделии.
Напряжения сжатия верхнего и нижнего образцов имеют одинаковые
параметры в зависимости от значений
.
Чтобы сделать удобным анализ и определение границ изменений
параметра нужно спроецировать построенные зависимости на плоскость
параметров
-
и рассмотреть линии уровня для разных значений
напряжений.
Согласно приведенным зависимостям, более жесткие ограничения на
параметры технологического процесса накладывают требования к напряжениям
растяжения, особенно это касается нижнего изделия. То есть параметры
и
определим для самого плохого изделия с точки зрения наибольших напряжений
растяжения. Выберем
и
так, чтобы при их комбинировании значения
68
напряжений растяжения не превысили 1.62 МПа (при учете запаса прочности в
10 %).
Исходя
из
выше
изложенного,
технологических параметров
и
получаем
следующие
значения
(см. таблицу 3.1).
Таблица 3.1 - Допустимые значения параметров
№ n/n
T1 , ч
T2, ч
1
8
1
2
7
2
3
6
3
4
5
4
5
4
5
6
3
6
7
4
6
8
3
7
9
2
8
Ниже приведены изменения среднего влагосодержания изделия во
времени при указанных в таблице 3.1 допустимых значениях параметров
(см. рисунок 3.4).
69
и
Рисунок 3.4 - Изменение среднего влагосодержания в процессе сушки
1 – исходный режим; 2 – новые режимы
Вычисляем
по
математической
влагосодержаний
изделий
параметров
(см. таблицу 3.2).
и
при
модели
выбранных
значения
величинах
технологических
Таблица 3.2 - Остаточные влагосодержания изделий
№ n/n
T1, ч
T2, ч
Uост %
1
8
1
5.0125
2
7
2
5.0125
3
6
3
5.0124
4
5
4
5.0123
5
4
5
5.0122
6
3
6
5.0121
7
4
6
5.0124
8
3
7
5.0124
9
2
8
5.0123
70
остаточных
Определяем время выхода системы по значению влагосодержания на
уровень остаточного влагосодержания исходного технологического процесса
(5.0220 %). Это время (таблица 3.3) можно считать общей продолжительностью
процесса
термовлажностной
обработки
при
заданных
значениях
технологических параметров (см. рисунок 3.5).
Рисунок 3.5 - Определение времени выхода на уровень влагосодержания исходного процесса
Таблица 3.3 - Определение общей продолжительности процесса сушки
№ n/n
T1, ч
T2, ч
Tобщ %
1
8
1
98.833
2
7
2
98.666
3
6
3
98.333
4
5
4
98.166
5
4
5
97.833
6
3
6
97.666
7
4
6
98.500
8
3
7
98.166
9
2
8
98.000
71
С целью проектирования системы управления очень важным является
увеличение быстродействия, т.е. уменьшение времени, затрачиваемого на
процесс сушки, так как данный процесс длится по времени до 90 % от времени
всего технологического цикла производства керамической плитки [29]. Таким
образом, сокращая период сушки керамических изделий, возможно сократить
время всего процесса и повысить производительность как процесса в
отдельности, так и производительность предприятия в целом.
Для равноценного учета границ изменения критериев
и
нужно
рассчитать нормируемый коэффициент, который определяется следующим
образом:
,
(3.13)
где n = 10 – количество рассматриваемых вариантов режимных параметров
сушки керамической плитки (таблица 3.1, включая исходный режим);
С
целью
определения
самого
лучшего
соотношения
.
весовых
коэффициентов, нужно вычислить значения функционала качества (см. рисунок
3.6) для допустимых вариантов параметров сушки керамической плитки,
приведенных в таблице 3.3.
Рисунок 3.6 - Значения функционала качества в зависимости от соотношения весовых
коэффициентов
Согласно данной таблице, общий минимум критерия наблюдаем при к 1 =
60% и к 2 = 40% для всех изучаемых наборов режимных параметров, что
72
свидетельствует о большой важности критерия быстродействия, который
определяет энергозатраты на сушку. С учетом выражения (3.13) критерий (3.11)
системы имеет вид:
(3.14)
С целью подбора параметров оптимального управления режимом
обобщим полученные результаты, которые мы получили с использованием
модели в таблицу, и для каждого отдельного случая проведем расчет значений
критерия, таблица 3.4.
Таблица 3.4 – Показатели определения оптимального режима
№ n/n
T1, ч
T2, ч
Tобщ %
Uоста, %
J
1
8
1
98.833
5.0125
98.9747
2
7
2
98.666
5.0125
98.8244
3
6
3
98.333
5.0124
98.5245
4
5
4
98.166
5.0123
98.3740
5
4
5
97.833
5.0122
98.0741
6
3
6
97.666
5.0121
97.9236
7
4
6
98.500
5.0124
98.6748
8
3
7
98.166
5.0124
98.3742
9
2
8
98.000
5.0123
98.2246
Минимизируя функционал, стало возможно вычислить условия самого
лучшего режима для термовлажностной обработки: T1 = 3ч, T2 = 6ч.
Использование
выбранных
основными
моделирования
показателях
показателями
режимных
качества
технологического
процесса
при
параметров,
получили
результаты,
которых
являются:
остаточное
влагосодержание (см. рисунок 3.7) и напряжения в них (см. рисунок 3.8), а
также
линии
уровня
температуры
73
потоков
при
наименьшей
производительности подающего вентилятора и температуре входящего воздуха
1000С при наличии внутренних сопротивлений (см. рисунок 3.9).
Рисунок 3.7 – Показатели остаточного влагосодержания в изделиях (120 ч)
а) верхний образец; б) нижний образец
Рисунок 3.8 – Показатели напряжения в керамических изделиях в ходе сушки
а) верхний образец; б) нижний образец
74
Рисунок 3.9 – Линии, показывающие уровень температуры потоков в условиях наименьшей
производительности подающего вентилятора и температуре входящего воздуха 100°С с
наличием внутренних сопротивлений
Выводы
В качестве итогового технологического режима термовлажностной
обработки
керамической
плитки
принимаем
режим
со
следующими
характеристиками:
•
общая продолжительность процесса: 98 часов;
•
параметр T1: 3 часа;
•
параметр Т2: 6 часов;
•
остаточная влажность изделий на 120 часов: 5.0177 % и 5.0121 %;
75
•
остаточная влажность изделий на 98 часов: 5.0314 % и 5.0218 %.
Данный режим удовлетворяет лишь одному набору исходных данных:
камера полностью загружена плиткой с начальной влажностью 20% и
температурой окружающей среды 20°С. Поэтому при изменении вида
керамической плитки, ее параметров или условий окружающей среды
необходимо повторно использовать алгоритм оптимизации.
76
ЗАКЛЮЧЕНИЕ
При выполнении научно-исследовательской работы нами получены
следующие результаты:
1.
Произведен анализ возможности использования методов теории
многосвязных однотипных систем для исследования камерных сушилок
периодического действия на стадии сушки керамической плитки.
2.
сушильной
Результаты моделирования распространения теплоносителя в
камере
показали,
что
теплоноситель
имеет
существенно
неравномерное распределение по объему как пустой, так и заполненной
изделиями камеры.
3.
Расчеты по математической модели термовлажностной обработки
керамической плитки позволило определить, что для анализа качества процесса
сушки в камере достаточно проанализировать характеристики небольшого
количества отдельно взятых образцов.
4.
Исследованы возможные варианты режимных параметров процесса
сушки керамической плитки, и с учетом оптимальности по быстродействию
при поддержании требуемого качества изделий найден их наилучший набор.
5.
Разработан
алгоритм
поиска
оптимального
управления
технологическим режимом сушки керамических изделий в сушильной камере
периодического действия, позволяющий сократить длительность обработки при
соблюдении требуемого качества продукции.
6.
Разработан алгоритм оптимизации регламента работы связанной
однотипной системы управления сушильных камер, позволяющей учесть
взаимодействие камеры и физическую распределенность параметров объекта.
Разработанная система управления технологическим процессом сушки
керамической плитки в камерных сушилках периодического действия
на
основе математической модели, при внедрении в производственный процесс
будет способствовать повышению эффективности работы действующего
оборудования на предприятиях, а также увеличению объема выпуска
77
продукции в условиях соблюдения требований государственного основного
стандарта.
78
СПИСОК ИСПОЛЬЗОВАННЫХ ИСТОЧНИКОВ
1.
Абрамович Г.Н. Прикладная газовая динамика: часть 1 - 5-е издание
переработанное и дополненное. - М.: Наука, 1991. - 597 с.
2.
Андерсон
Д., Таннехилл Дж., Плетчер
Р. Вычислительная
гидромеханика и теплообмен: В 2-х томах. Том 2: Перевод с английского. - М.:
Мир, 1990. – 728 - 392 с.
3.
Арутюнов В.А. и др. Математическое моделирование тепловой
работы промышленных печей. - М.: Высшая школа, 1985. - 236 с.
4.
Атанс М., Фабл П. Оптимальное управление. - М.: Наука, 1967. -
5.
Афифи
560 с.
А.,
Эйзен
С.
Статистический
анализ:
Подход
с
использованием ЭВМ. Пер. с англ.. - М.: Мир, 1982. - 488 с.
6.
Баскаков С.В. Сушка кирпича. - М.: Изд-во литературы по
строительству, 1966. - 322 с.
7.
Брандт 3. Анализ данных. Статистические и вычислительные
методы для научных работников и инженеров: Пер. с англ. - М.: Мир, ООО
«Издательство АСТ», 2003. - 686 с.
8.
Буров Ю.С. Технология строительных материалов и изделий.
Учебник для ВТУЗов. - М.: Высшая школа, 1972. - 464 с.
9.
Вентцель Е.С., Овчаров Л.А. Теория вероятностей и ее инженерные
приложения. - М.: Академия, 2003. - 464 с.
10.
Вентцель Е.С. Теория вероятностей: Учебник для студентов вузов
Изд. 9¬е. - М.: Академия, 2003. - 576 с.
11.
Гримитлин М.И.'Распределение воздуха в помещении. - М.:
Стройиздат, 1982. - 170 с.
12.
Демиденко
Н.Д.
Моделирование
и
оптимизация
тепломассобменных процессов в химической технологии. - М.: Наука, 1991. 239 с.
13.
Демидович Б.П., Марон И.А. Основы вычислительной математики -
М.: Наука, 1970. - 664 с.
79
14.
Духовный М.Л., Коен Г.Н., Копп В.Г., Орделли М.А., Юцис М.Л.
Сушка строительной керамики. - М.: Стройиздат, 1967. - 164 с.
15.
Дытнерский Ю.И. Процессы и аппараты химической технологии. -
Ч. II. - М.: Химия, 1995. - 560 с.
16.
Дьяконов В.П. MATLAB 6: учебный курс. - СПб.: Питер, 2001. - 592
17.
Жуков Д.В. Основы теории и техники сушки теплоизоляционных
с.
изделий. - М.: Стройиздат, 1974. - 206 с.
18.
Жуков Д.В. Скоростная сушка кирпича-сырца. - М.: Стройиздат,
1959. - 324 с.
19.
Закгейм А.Ю. Введение в моделирование химико-технологических
процессов. - М.: Химия, 1989. - 182 с.
20.
Зеличенок Г.Г. Автоматизация технологических процессов и учета
на предприятиях строительной индустрии. - М.: Высшая школа, 1975. - 351 с.
21.
Золотарский А.В., Шейнман Е.Ш. Производство керамического
кирпича. - М.: Стройиздат, 1975. - 386 с.
22.
Золотарский А.В., Шейнман Е.Ш. Производство керамического
кирпича. - М.: Высшая школа, 1989. - 264 с.
23.
Зорохович В.С. Системы управления машинами для производства
стеновой керамики. - Л.: Стройиздат, 1984. - 132 с.
24.
Зорохович В.С., Шукуров Э.Д. Производство кирпича. - Л.:
Стройиздат, 1988. - 232 с.
25.
Калиткин Н.Н. Численные методы. - М.: Наука, 1989. - 516 с.
26.
Канаев В.К. Новая технология строительной керамики. - М.:
Стройиздат, 1990. - 264 с.
27.
Касаткин А.Г. Основные процессы и аппараты химической
технологии. - М.: Химия, - 1971. - 482 с.
28.
Кафаров В.В., Глебов М.Б. Математическое моделирование
основных процессов химических производств: Учеб. Пособие для вузов. - М.:
Высшая школа, 1991. - 400 с.
80
29.
Кашкаев И.С., Никитин И.А., Володин Н.Н. Производство лицевых
керамических изделий. - М.: Стройиздат, 1977. - 176 с.
30.
Кашкаев И.С., Шейнман Е.Ш. Производство керамического
кирпича. - М.: Высшая школа, 1983. - 223 с.
31.
Кокрен У. Методы выборочного исследования. - М.: Статистика,
1976. - 440 с.
32.
Конвективный тепло- и массоперенос: Единое описание для
течения в каналах и внеш. обтекания тел любой формы и расположения/ В.
Каст, О. Кришер, Г. Райнике, К. Винтермантель. - М.: Энергия, 1980. - 46 с.
33.
Крамарухин Ю.Е. Приборы для измерения температуры. - М.:
Машиностроение, 1990. - 208 с.
34.
Кришер О. Научные основы техники сушки. - М.: изд-во инностр.
литературы, 1961. - 462 с.
35.
Лебедев П.Д. Расчет и проектирование сушильных установок. - М.-
Л.: Госэнергоиздат, 1963. - 320 с.
36.
Лебедев
П.Д.
Теплообменные,
сушильные
и
холодильные
установки. Тепломассообменные и холодильные установки. - М.: Энергия,
1972. - 320 с.
37.
Лебедев П.Д., Щукин А.А. Промышленная теплотехника. - М.-Л.:
Госэнергоиздат, 1956. - 384 с.
38.
Лебедев П.Д., Щукин А.А. Теплоиспользующие
установки
промышленных предприятий. - М.: Энергия, 1970. - 408 с.
39.
Луцык
Р.В.,
Ментковский
Ю.Л.,
Холод
В.П.
Взаимосвязь
деформационно-релаксационных и тепломассобменных процессов. - Киев:
Вища школа, 1992. - 184 с.
40.
Лыков А.В., Берковский Б.М. Конвекция и тепловые волны. - М.:
Энергия, 1974. - 335 с.
41.
Лыков А.В., Михайлов Ю.А. Теория тепло- и массопереноса. - М.:
Госэнергоиздат, 1963. - 535 с.
42.
Лыков А.В. Теория сушки. - М.: Энергия, 1968. - 472 с.
81
43.
Лыков А.В. Теория тепло - массопереноса. - М.: Высшая школа,
1963. – 486 с.
44.
Лыков А.В. Теория теплопроводности - М.: Высшая школа, 1967. -
45.
Лыков А.В. Тепломассообмен: справочник. - М.: Энергия, 1978. -
46.
Магергут B.3., Вент Д.П., Кацер И.А. Инженерные методы выбора и
600 с.
480 с.
расчета оптимальных настроек промышленных
регуляторов.
-
Новомосковск.: НФ РХТУ им. Д.И. Менделеева, 1994. - 158 с.
47.
Маликов В.Т., Кветный Р.Н. Вычислительные методы и применение
ЭВМ: учеб, пособие. - Киев: Выща шк., 1989. - 600 с.
48.
Мееров М.В., Ахметзянов А.В., Берщанский Я.М., Кулибанов В.Н.
Многосвязные системы управления. - М.: Наука, 1990. - 264 с.
49.
Методы классической и современной теории автоматического
управления. В 5-ти т. Т.4. Теория оптимизации систем автоматического
управления / Под ред. К.А. Пупкова, Н.Д. Егупова. - М.: Издательство МГТУ
им. Н.Э.Баумана, 2004. - 744 с.
50.
Мороз И.И. Автоматизация производства строительной керамики. -
Киев: Литература по строительству и архитектуре, 1961. - 208 с.
51.
Морозовский В.Т.
Многосвязные
системы
автоматического
регулирования. - М.: Энергия, 1970. - 320 с.
52.
Мэтьюз Джон Г., Финк Д., Куртис Д. Численные методы,
использование MATLAB. 3-е издание; пер. с англ. - М.: изд. дом «Вильямс»,
2001. - 720 с.
53.
Никитенко Н.И. Исследование процессов тепло- и массо- обмена
методом сеток. - Киев: Наук, думка, 1978. - 212 с.
54.
Никитенко Н.И. Теория тепломассопереноса. - Киев: Наук. думка,
1983. - 351 с.
55.
Нохратян К.А. Сушка и обжиг в промышленности строительной
керамики. - М.: Госстройиздат, 1962. - 603 с.
82
56.
Паниотто В.И. Качество социологической информации. - Киев:
Наукова думка, 1986. - 208 с.
57.
Перегудов В.В., Роговой М.И. Тепловые процессы и установки в
технологии строительных изделий и деталей. - М.: Стройиздат, 1983. – 416 с.
58.
Перегудов В.В. Теплотехника и теплотехническое оборудование. -
М.: Стройиздат, 1990. - 336 с.
59.
Понтрягин
Л.С.,
Болтянский
В.Г.,
Гамкрелидзе
Р.В.
Математическая теория оптимальных процессов. - М.: Наука, 1976. - 392 с.
60.
Потемкин В.Г. Вычисления в среде Matlab. - М.: Диалог-МИФИ,
2004. - 720 с.
61.
Преображенский В.П. Теплотехнические измерения и приборы:
Учебник для ВУЗов по специальности «Автоматизация теплоэнергетических
процессов» - 3-е издание, переработанное. - М.: «Энергия», 1978. - 704 с.
62.
Прокопенко М.Н. Верификация математической модели процесса
распределения
теплоносителя
в
камерной
сушилке//
Информационные
технологии в управлении и моделировании: сб. докл. Международной науч.технич. Интернет-конференции. - Белгород: Изд-во БГТУ им. В.Г. Шухова,
2005. - С. 164-167.
63.
Прокопенко
М.Н.
Математическая
модель
технологического
процесса термовлажностной обработки керамического кирпича// «Вестник
БГТУ им. В.Г. Шухова» Научно-теоретический журнал. - Ч. 111, 6/2003. -
С.
194-198.
64.
Прокопенко М.Н., Окунева Г.Л., Прасол Д.А. Моделирование
процесса распределения теплоносителя в камерной сушилке // Проектирование
инженерных и научных приложений в среде MATLAB: Сб. тр. Второй
Всероссийской научной конференции. - М.: ИПУ РАН, 2004. –
65.
сушки
С. 682-687.
Прокопенко М.Н. Оценка адекватности математической модели
керамических
изделий
//
Наука
Материалы VI Международной научной конференции.
Беловский полиграфист, 2006. - С. 520-525.
83
и
образование:
- 4.1. - Белово:
66.
Прокопенко М.Н. Постановка и решение задачи ' оптимизации
технологического
регламента
сушки
керамического
кирпича//
Вестник
Воронежского государственного технического университета. - т.2. - 2006. - №8.
- С. 146-153.
67.
Роговой
М.И. Теплотехническое оборудование керамических
заводов: учебник для техникумов. - М.: Стройиздат, 1983. - 367 с.
68.
Роуч П. Вычислительная гидродинамика. - М.: Мир, 1980. - 616 с.
69.
Рубанов В.Г., Прокопенко М.Н. Определение интервальное™
параметров
сушки
керамического
кирпича
в
камерных
сушилах//
Перспективные задачи инженерной науки: Сб. науч. тр. III Международной
научной конференции, Выпуск 3. - Днепропетровск: GAUDEAMUS, 2002. - С.
253-259.
70.
Сажин В.С. Основы техники сушки. - М.: Химия, 1984. - 320 с.
71.
Самарский А.А. Теория разностных схем. - М.: Наука, 1977. – 656 с.
72.
Самойлович Г.С. Гидрогазодинамика: Учебник для ВУЗОВ. - 2-е
издание, переработанное и дополненное. - М.: Машиностроение, 1990. – 384 с.
73.
Соболев О.С. Однотипные связанные системы регулирования. - М.:
Энергия, 1973. - 136 с.
74.
Строительные машины. Справочник. В 2т. Т.2: Оборудование для
производства строительных материалов и изделий / В.Н. Лямин, М.Н.
Горбовец, И.И. Быховский и др. - М.: Машиностроение, 1991. - 496 с.
75.
Сушильное оборудование для химических производств / под ред.
А.А. Корягина. - М.: НИИхиммаш, 1988. - 119 с.
76.
Сушка керамических стройматериалов пластического формования /
И.М. Пиевский, В.В. Гречина, Г.Д. Назаренко, А.И. Степанова. – Киев.:
Наукова думка, 1985. - 144 с.
77.
Тепловые процессы в технологии силикатов/ А.В. Ралко, А.А.
Крупа, Н.Н. Племянников. - Киев: Вища школа, 1986. - 232 с.
78.
Хигерович М.И., Байер В.Е. Производство глиняного кирпича. - М.:
Стройиздат, 1984. - 96 с.
84
79.
Черноруцкий И.Г. Методы оптимизации в теории управления:
Учеб, пособие. - СПб.: Питер, 2004. - 256 с.
80.
Чернявский E.B. Производство кирпича. - M.: Издательство
литературы по строительству, 1966. - 176 с.
81.
Чижский А.Ф. Сушка керамических материалов и изделий. - М.:
Стройиздат, 1971. - 208 с.
82.
Шаптала
В.Г.,
Окунева
Г.Л.
Численное
моделирование
воздухообмена производственных помещений на основе уравнений НавьеСтокса//
Математическое
моделирование
в
технологии
строительных
материалов: Сб. науч. тр. БТИСМ. - Белгород, 1992. - С. 49-54.
83.
Шувалов В.В., Огаджанов Г.А., Голубятников В.А. Автоматизация
производственных процессов в химической промышленности. - М.: Химия,
1991. - 478 с.
84.
Юшкевич М.О., Роговой М.И. Технология керамики. - М.:
Стройиздат, 1969. - 320 с.
85.
Abbas Emami-Naeini Feedback Control of Dynamic Systems. - USA:
Prentice Hall, 2002. - 680 p.
86.
Franks Roger G. Е. Mathematical modeling in chemical engineering. -
New York: John Wiley & Sons INC, 1971. - 272 p.
87.
Keey R.B. Introduction to industrial drying operations. - Oxford:
Pergamob Press, 1980. - 262 p.
88.
Sherwood Thomas K., Pigford Robert L., Wilke Charles R. Mass
transfer. - New York: McGraw-Hill, 1975. - 695 p.
85
ПРИЛОЖЕНИЕ А
Программа моделирования распределения тепловых потоков в
пустой камере
clear all;
close all;
clc;
D = 1;
%длина камеры
H = 0.7; %высота камеры
H1 = 0.2, % ширина входного отверстия
Н2 = 0.25; % ширина выходного отверстия
V1 = 34.24; % скорость воздуха, поступающего ч/з входное отверстие
V2 = (V1*H1) / Н2;% скорость воздуха, выходящего из камеры
Е = 0.01;% точность расчета скорости
Et = 0.01,% точность расчета температуры
N = 40;%кол-во точек по оси у (по длине)
М = 26; %кол-во точек по оси х (по высоте)
hy = D/N;% шаг по оси у
hx = Н/М;%шаг по оси х
F = zeros (M+1, N+1); % матрица значений ф-ии тока
Ux = zeros (M+1, N+1); % матрица значений х-составляющей скорости Uy zeros (M+1, N+1), % матрица значении у-составляющей скорости
U = zeros (M+1, N+1); % матрица значений скорости
Т = zeros (M+1, N+1), % матрица значений температуры
Тmах = 150, % значение максимальной температуры поступающего воздуха
lambda = 0 026, %коэффициент теплопроводности
ro = 1 2; % плотность воздуха
Ср = 1005, % теплоемкость воздуха
а = lambda/ (rо*Ср),% коэффициент температуропроводности
% Рассчет скорости
%левая стенка
for i=1:М+1
Vx (i, 1) = 0;
Vy (i, 1) = 0;
F (i, 1) = V1*D/2;
end
%правая стенка
for i=1:М+1
Vx (i, N+1) = 0;
Vy (i, N+1) = 0;
F (i, N+1) = V2*D/2;
end
%потолок
for j = 1:6
Vx (1, j) = 0;
Vy (1, j) = 0;
F(1, j) = 0;
86
end
for j = 7:11
Vx(1, j) = V1;
Vy (1, j)= 0;
F (1, j) = V1*( (N / 2+1) - j) *hy;
end
for j = 12 :120
Vx (1, j) = 0;
Uy (1, j) = 0;
F (1, j) = 0;
end
for j = 121:125
Vx (1, j) = V2;
Vy(1, j) = 0;
F (1, j) = V2* (j-(N/2 + 1))*hy;
end
for j = 126:N+1
Vx (1,j) = 0;
Vy (1, j) = 0;
F (1, j) = 0;
end
% пол
for j = 1:62
Vx (M+1, j) = 0;
Vy (M+1, j) = 0;
F (M+1, j) = V1*D/2;
end
for j = 63:68
Vx (M+1, j) = V3;
Vy (M+1, j) = 0;
F (M+1, j) = V3*(abs(N/2+1-j))*hy;
end
for j = 69:N+1
Vx (M+1, j) = 0;
Vy (M+1, j) = 0;
F(M+1, j) = U2*D/2;
end
for i = 2:M
for j = 2:N
f(1, j) = (F (1, j) + F(M+1, j ))/M;
end
end
% линия симметрии
for i = 1:M
Uy (i,(N/2+1)) = 0;
F (i,(N/2+1)) = 0;
end
i = 0,
R0 = E + 1,
whi1e R0 > E
87
R0 = 0;
i = i + 1;
for i = 2:M
for j = 2:N/2
Q = (hy*hy*(F(i+1, j) + f(i-1, j)) + hx*hx*(F(i, j+i) + f(i, j-1)))/(2*(hx*hx+hy*hy));
R = abs (Q - F(i, j) );
if R > R0
R0 = R;
end
F(i, j) = Q,
end
end
end
i = 0;
R0 = Е + 1;
whi1e R0 > E
R0 = 0;
i=i+1;
for 1 = 2 : M
for j = 67-N
Q=(hy*hy*(F(i+1, j)+F(i-1, j))+hx*hx* (F (i, j+1) + F(i, j-1))) / (2*(hx*hx+hy*hy));
R = abs (Q - F(1, j));
if R > R0
R0 = R;
end
F(i, j) = Q,
end
end
end
for i = 2:M
for j = 2: N
Uy (i, j) = (F(i+1, j) - F(i-1, j)) / (-2*hx);
Ux (i, j) = (F(i, J+1) - F(i, j-1))/ (2*hy);
U (i, j) = sqrt (Ux(i, j) *Ux (i, j) + Uy(i, j)*Uy(i, j));
end
end
for j = 1:N+1
U(i, j) = sqrt (Ux(1, j)*Ux(1, j) + Uy(1, j)*Uy(1, j));
U(M+1, j) = sqrt (Ux(M+1, j)*Ux (M+1, j) + Uy (M+1, j)*Uy(M+1, j));
end
X = 1:1:N+1;
У = 1:1:M+1;
figure(1)
C = U,
surf (X, Y, U, C)
%colormap jet
%colorbar
Axis ([1;N+1;1;M+1;0;150])
xlabel ('y')
ylabel ('x')
zlabel ('V, M/c')
Tit1e (‘Speed Fie1d’)
V = [0.1 0 2 0 3 0 4 0 5 0.6 0 7 0 8 0 9 1 2 3 4 5 6 7 8 9 10 11 12]
% 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41
% 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60,
Figyre (2)
88
[C, h] = contour (X, Y, U, V);
c1abe1(C, h, V)
%Рассчет температуры
for i=1:M+1
for j=1:N+1
T(i, j) = 20,
end
end
%потолок
for j= 7:11
T (1, j) = Tmax;
end
for j = 121:125
T (1, j) = Tmax;
end
i = 0,
R0 = Et + 1;
whi1e R0 > Et
R0 = 0,
i=i+1;
%потолок
for j = 2:6
T(1,j) = (2*T(2,j)+T(1,j-1)+T(1,j+1))/(4+(Ux(1,j)*Ux(1,j)*hx*hx)/(a*a) + (2*Ux(1,j)*hx)/a);
end
for j = 7:11
T(1,j)
=
(2*T(2,j)+T(1,j-1)+T(,j+1)
+
(Ux(1,j)*Ux(1,j)*Tmax*hx*hx)/(a*a)
+
(2*Ux(1,j)*Tmax*hx)/a)/(4+(Ux(1,j)* Ux(1,j) *hx*hx)/(a*a) + (2*Ux(1, j)*hx)/a);
end
for j = 12:120
T(1,j) = (2*T(2,j)+T(1,j-1)+T(1,j+1))/(4+(Ux(1,j)*Ux(1, j)*hx*hx)/(a*a) + (2*Ux(1, j) *hx) /a);
end
for j = 121:125
T(1, j) =(2*T(2,j)+T(1,j1)+T(1,j+1)+(Ux(1,j)*Ux(1,j)*Tmax*hx*hx)/(a*a)+(2*Ux(1,j)*Tmax*
hx)/a)/(4+(Ux(1,j)* Ux(1,j)*hx*hx)/(a*a)+(2*Ux(1,j)*hx)/a);
end
for : = 126:N
T (1,j) = (2 *T(2,j)+T(1,j-1)+T(1,j+1))/(4+(Ux(1,j)*Ux(1, j)*hx*hx)/(a*a) + (2*Ux(1,j)*hx)/a);
end
%пол
for j = 2:62
T(M+1, j) = (2*T(M,j)+T(M+1,j-1)+T(M+1,j + 1))/ (4+(Ux(M+1,j)*Ux(M+1,j)*hx*hx)/(a*a)
-(2*Ux(M+1,j)*hx)/a);
end
for j = 63:68
T(M+1, j) = (2*Т(М,:)+Т(М+1,j-1)+Т(М+1,j+1))/4;
end
for j = 69:N
89
T(M+1, j) = (2*T(M,j)+T(M+1,j-1)+T(M+1,j + 1))/(4+(Ux(M+1,j)*Ux(M+1,j)*hx*hx)/(a*a)
- (2*Ux (M+1,j)*hx)/a);
end
%левая стенка
for 1 = 2:M
T(i, 1) = (2 *T(i,2)+T(i-1,1)+T(i + 1,1))/(4+(Uy(i,1)*Uy(i,1)*hy*hy)/ (a*a) + (2*Uy (i, 1)*hy)
/a);
end
%правая стенка
for 1 = 2:M
T(i, N+1) = (2*T(i,N)+T(i-1,N+1) +T(i + 1,N+1))/(4+(Uy(i,N+i) *Uy(i,N+1)*hy*hy)/(a*a)(2*Uy(i,N+1)*hy)/a);
end
%Углы
T(1, 1) =(2*T(2,1)+T(2,1)+T(1,2))/(4+(Ux(1,1)*Ux(1,1)*hx*hx)/(a*a)+(2*Ux(1,1)*hx)/a);
T(1, N+1)=(2*T(2,N+1)+T(1,N)+T(2,N+1))/(4+(Ux(1,N+1)*Ux(1,N+1)*hx*hx)/(a*a)+(2*Ux
(1,N+1)*hx) /a);
T(M+1, 1) = (2*T(M,1)+T(M,1)+T(M+1,2))/(4+(Ux(M+1,1)*Ux(M+1,1)*hx*hx)/(a*a)-(2*Ux
(M+1, 1) *hx) /a);
T(M+1, N+1) =(2 *Т(М,N+1)+T(M+1,N)+T(M,N+1))/(4+(Ux(M+1,N+1)*Ux(M+1,N+1)*hx*
hx)/(a*a)-(2*Ux(M+1,N+1)*hx)/a);
for i = 2:M
for j= 2:N
Q = (T(i-1, j)+T(i+1,j) +T(i, j-1)+T(,j+1)/4 - (T(i+1, j)-T(i-1,j))*(hx*Ux(i,j)/8*a) (T(i,
j+1)-T(i, j-1))*(hy*Uy(i,j)/8*a);
R = abs(Q - T(i, j));
if R > R0
R0 = R;
end
T(i, j) = Q;
end
end
clc
end
X = 1:1:N+1;
Y = 1:1:M+1;
figure(3)
C = T;
surf(X, Y, T, C)
colormap jet
colorbar
axis([1;N+1;1;M+1;0;150])
xlabel('у')
ylabel('x')
zlabel ('T’, C')
title('Temperature field')
W = [20 20.1 20.2 20.3 20.4 20.5 20.6 20.7 20.8 20.9 21 22 22.5 23 24 25 27.5 30 32.5 35 37.5 40
90
45 50 55 60 65 70 75 80 85 90 95 100 110 120 130 140 150 155];
figure (4)
[C, h] = contour (X, Y, T, W);
clabel(C, h, W)
Tsr = 0;
for i=1:M+1
for j=1:N+1
Tsr = Tsr + T(i, j);
end
end
Tsr = Tsr/(M*N);
Tsr
Tmin=min(min(T( , )));
Tmin
z=Tmax-Tmin;
z
91
ПРИЛОЖЕНИЕ Б
Результаты
моделирования
процесса
распределения
тепловых
потоков в пустой камере
а)
б)
Рисунок Б.1 – Линии уровня скорости при различной производительности подающего
вентилятора а) минимальной б) максимальной
Рисунок Б.2 – Теплотехнологические характеристики агента в камере при минимальной
производительности подающего вентилятора и исходной температуре 300С
а) линии уровня температуры б) температурное поле
92
Рисунок Б.3 – Теплотехнологические характеристики агента в камере при минимальной
производительности подающего вентилятора и исходной температуре 500С
а) линии уровня температуры б) температурное поле
Рисунок Б.4 – Теплотехнологические характеристики агента в камере при минимальной
производительности подающего вентилятора и исходной температуре 1000С
а) линии уровня температуры б) температурное поле
93
ПРИЛОЖЕНИЕ В
Программа моделирования распределения тепловых потоков в
заполненной изделиями камере
clear all;
close all;
c1c;
D = 1;
%длина камеры
H = 0.7; %высота камеры
H1 = 0.2, % ширина входного отверстия
Н2 = 0.25; % ширина выходного отверстия
V1 = 34.24; % скорость воздуха, поступающего ч/з входное отверстие
V2 = (V1*H1) / Н2;% скорость воздуха, выходящего из камеры
Е = 0.01;% точность расчета скорости
Et = 0.01,% точность расчета температуры
N = 40;%кол-во точек по оси у (по длине)
М = 26; %кол-во точек по оси х (по высоте)
hy = D/N;% шаг по оси у
hx = Н/М;%шаг по оси х
F = zeros (M+1, N+1); % матрица значений ф-ии тока
Ux = zeros (M+1, N+1); % матрица значений х-составляющей скорости
Uy - zeros (M+1, N+1), % матрица значении у-составляющей скорости
U = zeros (M+1, N+1); % матрица значений скорости
Т = zeros (M+1, N+1), % матрица значений температуры
Тmах = 150, % значение максимальной температуры поступающего воздуха
lambda = 0 026, %коэффициент теплопроводности
ro = 1 2; % плотность воздуха
Ср = 1005, % теплоемкость воздуха
а = lambda/ (rо*Ср),% коэффициент температуропроводности
% Рассчет скорости
%левая стенка
for i = 1:M+1
Ux(i, 1) = 0;
Uy(i, 1) = 0;
F(i, 1) = V1*D/2;
end
%правая стенка
for i = 1:M+1
Ux(i, N+1) = 0;
Uy(i, N+1) = 0;
F(i, N+1) = V2*D/2;
end
%потолок
for j = 1:6
Ux(1, j) = 0;
Uy(1, j) = 0;
F(1, j) = 0;
94
end
for j = 7:11
Ux(1, j) = V1;
Uy(1, j) = 0;
F(1, j) = V1*((N/2+1) - j)*hy;
end
for j = 12:120
Ux(1, j) = 0;
Uy(1, j) = 0;
F(1, j) = 0;
end
for j = 121:125
Ux(1, j) = V2;
Uy(1, j) = 0;
F(1, j) = V2*(j-(N/2+1))*hy;
end
for j = 126:N+1
Ux(1, j) = 0;
Uy(1, j) = 0;
F(1, j) = 0;
end
%пол
for j = 1:62
Ux(M+1, j) = 0;
Uy(M+1, j) = 0;
F(M+1, j) = V1*D/2;
end
for j = 63:68
Ux(M+1, j) = V3;
Uy(M+1, j) = 0;
F(M+1, j) = V3*(abs(N/2+1-j))*hy;
end
for j = 69:N+1
Ux(M+1, j) = 0;
Uy(M+1, j) = 0;
F(M+1, j) = V2*D/2;
end
for i = 2:M
for j = 2:N
F(i, j) = (F(1, j) + F(M+1, j))/M;
end
end
%линия симметрии
for i = 1:M
Uy(i, (N/2+1)) = 0;
F(i, (N/2+1)) = 0;
end
%Сопротивления
95
for n = 0:11
for i=5:21
for j = 7:15
Ux(i, j+10*n) = 0;
Ux(i+1, j+10*n) = 0;
Ux(i+2, j+10*n) = 0;
Uy(i, j+10*n) = 0;
Uy(i+1, j+10*n) = 0;
Uy(i+2, j+10*n)=0;
F(i, j+10*n) = 0;
F(i+1, j+10*n) = 0;
F(i+2, j+10*n) = 0;
end
end
end
i=0
R0 = E+1;
while R0 >E
R0 = 0;
i = i+1;
for i = 2:M
for j = 2:6
Q = (hy*hy*(F(i+1, j) + F(i-1, j)) + hx*hx*(F(i, j+1) + F(i, j-1))) / (2(hx*hx+hy*hy));
R = abs(Q - F(i, j));
if R > R0
R0 = R;
end
F(i , j) = Q;
end
end
for i = 2:4
for j = 7:N/2
Q = (hy*hy**(F(i+1, j) + F(i-1, j)) + hx*hx(F(i, j+1) + F(i, j-1))) / (2*(hx*hx+hy*hy));
R = abs(Q - F(i, j));
if R >R0
R0 = R;
F(i, j) = Q;
end
end
for i = 24:25
for j = 7:N/2
Q = (hy*hy**(F(i+1, j) + F(i-1, j)) + hx*hx(F(i, j+1) + F(i, j-1))) / (2*(hx*hx+hy*hy));
R = abs(Q - F(i, j));
if R >R0
R0 = R;
F(i, j) = Q;
end
end
96
for m = 0:5
for i = 8:20
for j = 7:15
Q = (hy*hy*(F(i+1, j+10*m) + F(i-1, j+10*m)) + hx*hx*(F(i, j+1+10*m) + F(i, j1+10*m))) / (2*(hx*hx+hy*hy));
R = abs(Q - F(i, j+10*m));
if R > R0
R0 = R;
end
F(i, j+10*m) = Q;
end
end
end
for i = 5:23
for j = 16:56
Q = (hy*hy*(F(i+1, j) + F(i-1, j)) + hx*hx*(F(i, j+1) + F(i, j-1))) / (2*(hx*hx+hy*hy));
R = abs(Q - F(i, j);
if R > R0
R0 = R;
end
F( i, j) = Q;
end
end
end
i =0
R0 = E +1;
while R0 > E
R0 = 0,
i = i +1;
for i = 2:4
for j = 67:N
Q = (hy*hy*(F(i+1, j) + F(i-1, j)) + hx*hx*(F(i, j+1) + F(i, j-1)))/
2*(hx*hx+hy*hy));
R = abs(Q - F(i, j),
if R > R0
R0 = R;
end
F(i, j) = Q;
end
end
for i = 24:25
for j = 67:N
Q = (hy*hy*(F(i+1, j) + F(i-1, j)) + hx*hx*(F(i, j+1) + F(i, j-1)))/
(2*(hx*hx+hy*hy));
R = abs(Q - F(i, j));
if R > R0
R0 = R;
97
end
F(i, j) = Q;
end
end
for m = 0:5
for i = 8:20
for j = 67:75
Q = (hy*hy*(F(i+1, j+10*m) + F(i-1, j+10*m)) + hx*hx*(F(i, j+1+10*m) +
F(i, j-1+10*m)))/(2*(hx*hx+hy*hy));
R = abs(Q - F(i, j+10*m));
if R > R0
R0 = R;
end
F(i, j+10*m) = Q;
end
end
end
for i = 5:23
for j = 76:116
Q = (hy*hy*(F(i+1, j) + F(i-1, j)) + hx*hx*(F(i, j+1) + F(i, j-1)))/
(2*(hx*hx+hy*hy));
R = abs(Q - F(i, j));
if R > R0
R0 = R;
end
F(i, j) = Q;
end
end
for i = 2:M
for j = 126:N
Q = (hy*hy*(F(i+1, j) + F(i-1, j)) + hx*hx*(F(i, j+1) + F(i, j-1)))/
(2*(hx*hx*hy*hy));
R = abs(Q - F(i, j));
if R > R0
R0 = R;
end
F(i, j) = Q;
end
end
end
end
for i = 2:M
for j = 2:6
Uy(i, j) = (F(i+1, j) - F(i-1, j)) / (-2*hx);
Ux(i, j) = (F(i, j+1) - F(i, j-1)) / (2*hy);
U(i, j) = sqrt(Ux(i, j)*Ux(i, j) + Uy(i, j)*Uy(i, j));
end
end
98
for i = 2:M
for j = 126:N
Uy(i, j) = (F(i+1, j) - F(i-1, j)) / (-2*hx);
Ux(i, j) = (F(i, j+1) - F(i, j-1)) / (2*hy);
U(i, j) = sqrt(Ux(i, j)*Ux(i, j) + Uy(i, j)*Uy(i, j));
end
end
for i = 2:4
for j = 7:125
Uy(i, j) = (F(i+1, j) - F(i-1, j)) / (-2*hx);
Ux(i, j) = (F(i, j+1) - F(i, j-1)) / (2*hy);
U(i, j) = sqrt(Ux(i, j)*Ux(i, j) + Uy(i, j)*Uy(i, j));
end
end
for i = 24:25
for j = 7:125
Uy(i, j) = (F(i+1, j) - F(i-1, j)) / (-2*hx);
Ux(i, j) = (F(i, j+1) - F(i, j-1)) / (2*hy);
U(i, j) = sqrt(Ux(i, j)*Ux(i, j) + Uy(i, j)*Uy(i, j));
end
end
for i = 5:3
for j = 16:56
Uy(i, j) = (F(i+1, j) - F(i-1, j)) / (-2*hx);
Ux(i, j) = (F(i, j+1) - F(i, j-1)) / (2*hy);
U(i, j) = sqrt(Ux(i, j)*Ux(i, j) + Uy(i, j)*Uy(i, j));
end
end
for i = 5:23
for j = 76:116
Uy(i, j) = (F(i+1, j) - F(i-1, j)) / (-2*hx);
Ux(i, j) = (F(i, j+1) - F(i, j-1)) / (2*hy);
U(i, j) = sqrt(Ux(i, j)*Ux(i, j) + Uy(i, j)*Uy(i, j));
end
end
for m = 0:11
for i = 8:20
for j = 7:15
Uy(i, j+10*m) = (F(i+1, j+10*m) - F(i-1, j+10*m)) / (-2*hx);
Ux(i, j+10*m) = (F(i, j+1+10*m) - F(i, j-1+10*m)) / (2*hy);
U(i, j+10*m) = sqrt(Ux(i, j+1+10*m)*Ux(i, j+10*m) + Uy(i , j+10*m)*Uy(i, j+10*m));
end
end
end
for j = 1:N+1
99
U(1, j) = sqrt (Ux(1, j) *Ux(1, j) + Uy(1, j)*Uy(1, j));
U(M+1, j) = sqrt(Ux(M+1, j)*Ux(M+1, j) + Uy(M+1, j)*Uy(M+1, j));
end
X = 1:1:N+1;
Y =1:1:M+1;
figure(1)
C=U
surf (X, Y, U, C)
%colormap jet
%colorbar
axis([1;N+1;1;M+1;0;100])
xlabel('у')
ylabel('x')
zlabel('U, М/с')
title('Speed field')
V= [1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25];
%V = [0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 2 3 4 5 6 7 8 9 10 11 12 13];
%21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 ...
%42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60];
figure (2)
[C, h] = contour(X, Y, U, V);
clabel(C, h, V)
xlabel('y, points')
ylabel('x, points')
%Рассчет температуры
for i =1:M+1
for j = 1:N+1
T(i, j) = 20,
end
end
%потолок
for j = 7:11
T (1, j) = Tmax;
end
for j = 121:125
T (1, j) = Tmax;
end
i = 0,
R0 = Et + 1,
while R0 > Et
R0 = 0,
i = i + 1,
%потолок
for j =2:6
T(1, j) = (2*T(2, j)+T(1, j-1)+T(1, j+1))/(4+(Ux(1,j)*Ux(1,j)*hx*hx)/(a*a)+(2*Ux(1,
j)*hx)/a);
end
for j = 7:11
T(1, j) = (2*T(2,j)+T(1,j-1)+T(1,j+1) + (Ux(1,j)*Ux(1,j)*Tmax*hx*hx)/(a*a) +
(2*Ux(1,j)*Tmax*hx)/a)/(4+(Ux(1, j)* Ux (1, j) *hx*hx) / (a*a) + (2*Ux(1, j) *hx) /a);
end
for j = 12:120
100
T(1, j) = (2 *T(2,j)+T(1,j-1)+T (1,j+1)) / (4+ (Ux(1,j)*Ux(1,j)*hx*hx)/(a*a) +
(2*Ux(1,j)*hx)/a);
end
for j = 121:125
T(1, j) = (2*T(2,j)+Т(1,j-1)+Т(1,j +1) + (Ux(1,j)*Ux(1,j)*Tmax*hx*hx)/(а*а) +
(2*Ux(1,j)*Tmax*hx)/а)/(4+(Ux(1,j)* Ux(1,j)*hx*hx)/(а*а)+(2*Ux(1,j)*hx)/а);
end
for j = 126:N
T(1, j) = (2*T(2,j)+T(1,j-1)+T(1,j+1)) / (4 + (Uх(1, j) *Ux(1, j) *hx*hx)/ (a*a) +
(2*Ux(1,j) *hx) /a);
end
%пол
for j = 2:62
T(M+1,j)=(2*Т(М,j)+Т(М+1,j1)+T(M+1,j+1))/(4+(Ux(M+1,j)*Ux(M+1,j)*hx*hx)/(a*a)(2*Ux(M+1,j)*hx)/a);
end
for j = 63:68
T(M+1, j) = (2*Т(М,j)+Т(М+1,j-1)+T(M+1,j+1)) /4,
end
for j = 69:N
T(M+1,j)=(2*Т(М,j)+Т(М+1,j-1)+T(M+1,j+1))/(4+(Ux(M+1,j)*Ux(M+1,j)*hx*hx)
/(а*а)-(2*Uх(М+1,j)*hx)/a);
end
%левая стенка
for i = 2:M
T((i,1)=(2*T(i,2)+T(i1,1)+T(i+1,1))/(4+(Uy(i,1)*Uy(i,1)*hy*hy)/(a*a)+(2*Uy(i,1)*hy);
end
%правая стенка
for i = 2:M
T(i,N+1)=(2*T(i,N)+T(i-1,N+1)+T(i+1,N+1))/(4+(Uy(i,N+1)*Uy(i,N+1)*hy*hy)/(a*a)(2*Uy (i,N+1)*hy) /a);
end
%Сопротивления
%верхняя грань
for n=0:11
for i=5:21
for j=8:14
T(i, j+10*n) = (2*T(i-1,j+10*n)+T(i,j-1+10*n)+T(i,j+1+10*n))/(4+(Ux(i,j+10*n)
*Ux(i,j+10*n)*hx*hx)/(a*a)-(2*Ux(i,j+10*п)*hx)/a);
end
end
end
%нижняя грань
for n=0:11
for i=7:23
for j=8:14
T(i, j+10*n) = (2*T(i+1,j+10*n)+T(i,j-1+10*n)+T(i,j+1+10*n))/(4+(Ux
(i,j+10*n)*Ux(i,j+10*n)*hx*hx)/(a*a)+(2*Ux(i,j+10*n)* hx)/a);
end
end
end
101
%левая грань
for n=0:4
for j=7:117
for i=5:7
T (i + 4*n, j) = (2*T(i+4*n,j-1)+T(i-1+4*п,j)+Т(i+1+4*п,j))/(4+(Uy(i+4*
n,3)*Uy(i+4*n,j)*hy*hy)/(a*a)(2*Uy(i+4*n,j);
end
end
end
%правая грань
for n=0:4
for j=15:125
for i=5:7
T (i + 4 *n, j) = (2*T(i+4*n,j+1)+T(i1+4*n,j)+T(i+1+4*n,j))/(4+(Uy(i+4*n,j)*Uy(i+4*n,j)*hy*hy)/(a*a)+(2*Uy(i+4*n,j)*hy)/a);
end
end
end
%Углы
T(1,1)=(2*T(2,1)+T(2,1)+T(1,2))/(4+(Ux(1,1)*Ux(1,1)*hx*hx)/(a*a)+(2*Ux(1,1)*hx)/a)
T(1, N+1) =(2*T(2,N+1)+T(1,N)+T(2,N+1))/(4+(Ux(1,N+1)*Ux(1,N+1)*hx*hx)/
(a*a)+(2*Ux(1,N+1)*hx)/a);
T(M+1, 1) = (2*T(M,1)+T(M,1)+T(M+1,2))/(4+(Ux(M+1,1)*Ux(M+1,1)*hx*hx)/(a*a)- (2
*Ux(M+1,1)*hx)/a);
T(M+1,N+l)=(2*T(M,N+1)+T(M+1,N)+T(M,N+1))/(4+(Ux(M+1,N+1)*Ux(M+1,N+1)*hx
*hx)/(a*a)- (2*Ux(M+1,N+1)*hx)/a);
for i = 2:M
for j =2:6
Q = (T(i-1, j)+T(i+1, j)+T(i, j-1)+T(i, j+1))/4 – (T(i+1, j) – T(i-1, j))*(hx*Ux(i ,j)/8*a) (T(i, j+1) - T(i, j-1))*(hy*Uy(i, j)/8*a);
R = abs(Q - T(i, j));
if R > R0
R0 = R;
end
T(i, 3) = Q;
end
end
for i = 2 : 4
for j=7:125
Q = (T(i-1, j)+T(i+1, j)+T(i, j-1)+T(i, j+1))/4 - (T(i+1, j) - T(i-1, j))*(hx*Ux(i ,j)/8*a) - (T(i,
j+1) - T(i, j-1))*(hy*Uy(i, j)/8*a);
R = abs(Q - T(i, j ) ) ;
if R > R0
R0 = R;
end
T(i, j ) = Q ;
end
end
for i = 24:25
for j=7:125
Q = (T(i-1, j)+T(i+1, j)+T(i, j-1)+T(i, j+1))/4 - (T(i+1, j) - T(i-1, j))*(hx*Ux(i ,j)/8*a) - (T(i,
102
j+1) - T(i, j-1))*(hy*Uy(i, j)/8*a);
R = a b s (Q - T(i, j));
if R > R0
R0 = R;
end
T(i, j) = Q;
end
end
for i = 2 M
for j = 126:N
Q = (T(i-1, j)+T(i+1, j)+T(i, j-1)+T(i, j+1))/4 - (T(i+1, j) - T(i-1, j))*(hx*Ux(i ,j)/8*a) - (T(i,
j+1) - T(i, j-1))*(hy*Uy(i, j)/8*a);
R = abs(Q - T(i, j));
if R > R0
R0 = R;
end
T(i, j) = Q;
end
end
for i = 5:23
for j = 16:10:116
Q = (T(i-1, j)+T(i+1, j)+T(i, j-1)+T(i, j+1))/4 - (T(i+1, j) - T(i-1, j))*(hx*Ux(i ,j)/8*a) - (T(i,
j+1) - T(i, j-1))*(hy*Uy(i, j)/8*a);
R = abs(Q - T(i, j));
if R > R0
R0 = R;
end
T(i, j) = Q;
end
end
for m = 0:11
for i = 8.4:20
for j = 7:15
Q = (T(i-1, j+10*т)+T(i+1, j+10*m)+T(i, j-1+10*m)+T(i, j+1+10*m))/4 - (T(i+1,
j+10*m)-T(i-1,j+10*m))*(hx*Ux(i,j+10*m)/8*a) - (T(i,j+1+10*m)-T(i,j-1+10*m))*(hy*Uy
(i,j+10*m)/8*a);
R = abs(Q - T(i, j+10*m));
if R > R0
R0 = R;
end
T(i, j+10*m) = Q;
end
end
end
clc
1
end
1
X = 1;1;N+1;
Y = 1;1;M+1;
figure (3)
C = T,
103
surf (X, Y, T, C)
colormap jet
colorbar
axis([1;N+1;1;M+1;0;150])
xiabei('у')
yiabei('x')
ziabei (' T, C')
titie('Поле температур')
%W = [20 22 24 26 28 30 32 34 36 38 39 40 41 42 43 44 45 46 47 48 50];
W = [45 50 55 60 65 70 75 80 82 85 87 90 93 95 97 100],
figure (4)
[C, h] = contour (X, Y, T, W),
ciabei(C, h, W)
xiabei('у')
yiabei('x')
grid on
Tvsr=zeR0s(10, 6),
for n=0:5
for i=5:21
for j=8:14
T(i+1,j+10*n)=(T(i, j+10*n)+T(i+2,j+10*n))/2;
end
end
end
for n=0:5
for i=1:5
for j=7:15
Tvsr(i*2-1,n+1)=Tvsr(i*2-1,n+1)+T(1+4*i,j+10*n)+T(1+4*i+1,j+10*n);
Tvsr(i*2,n+1)=Tvsr(i*2,n+1)+T(1+4*i+1,j+10*n)+T(1+4*i+2,j+10*n);
end
end
end.
for i=1:10
for j=1:6
Tvsr(i, j)=Tvsr(i, j)/18,
end
end
Tvsr
Tsredn = 0,
for i = 2:M
for j = 2 : 6
Tsredn = Tsredn + T(i, j ),
end
end
for i = 2:4
for j = 7:125
Tsredn = Tsredn + T(i, j ),
end
end
for i = 24:25
for j = 7:25
Tsredn = Tsredn + T(i, j );
104
end
end
for i = 2 :M
for j = 126:N
Tsredn = Tsredn + T(i, j );
end
end
for i = 5:23
for j = 16:116
Tsredn = Tsredn + T(i, j );
end
end
for m = 0 : 1 1
for i = 8 . 4 :20
for j = 7:15
Tsredn = Tsredn + T(i, j + 1 0 * m );
end
end
end
Tsredn = Tsredn/((M*N)-(3*9*60))
Tsredn
%T(3,8 )
%T(2 0 ,8 )
%Т(12,33)
%Т(24,64)
%Т(24,20)
%Т(12,40)
er=zeros(М+1, N+1),
for i =2 :М
for j=2 : N
er(i ,j)=Т(i j)-Tsredn,
end
end
Min=abs(er(2,2)),
for i=2:M
for j = 2 : 6 6
if abs(er(i, j)) < Min
Min=abs(er(i, j)),
verh=i,
prav=j;
end
end
end
Min
verh
prav
T(verh,prav)
T(verh,34)
105
Отзывы:
Авторизуйтесь, чтобы оставить отзыв