УДК 621.78.68.3
АНАЛИЗ РАСЧЕТНЫХ РЕЖИМОВ ОПТИМАЛЬНОГО УПРАВЛЕНИЯ ВАКУУМНОЙ ЦЕМЕНТАЦИЕЙ1
М.Ю. Деревянов, М.Ю. Лившиц, Д.М. Федорченко2
Самарский государственный технический университет 443100, г. Самара, ул. Молодогвардейская, 244
E-mail: entcom@samgtu.ru
Разработана математическая модель вакуумной цементации. Произведена реструктуризация модели для использования её в оптимизационных процедурах. Проанализированы зависимости расчетных режимов оптимального управления от технологических параметров процесса.
Ключевые слова: вакуумная цементация, математическая модель, коэффициент мас-сопереноса, коэффициент диффузии, управляющее воздействие, оптимизационные процедуры, альтернансный метод оптимизации, метод функций Грина.
Введение
Поверхностное упрочнение деталей путем их химико-термической обработки (ХТО) широко распространено в машиностроении в силу универсальности и экономичности этого метода [1]. Среди различных видов ХТО наиболее широко применяется цементация, причем в последнее время вместо традиционной газовой цементации используют вакуумную, обеспечивающую ряд преимуществ перед другими методами цементации:
• не происходит внутреннего окисления деталей из-за отсутствия кислородсодержащих компонентов в атмосфере;
• поверхность деталей после цементации оказывается светлой;
• отсутствует необходимость в газоприготовительных установках и приборах контроля углеродного потенциала;
• сокращается длительность процесса в результате его проведения при более высокой температуре;
• уменьшается удельный расход электроэнергии и технологического газа.
Математическая модель процесса вакуумной цементации
В ходе вакуумной цементации происходит диффузия углерода из атмосферы печи через поверхность детали. При этом в силу незначительной глубины диффузионного слоя по сравнению с толщиной детали для большинства обрабатываемых изделий, поверхность которых не имеет малых радиусов скруглений, кромок и т.п., процесс диффузии может быть описан в соответствии со вторым законом Фика, краевой задачей параболического типа для полубесконечной пластины [2]:
1 Работа выполнена в рамках реализации федеральной целевой программы «Научные и научно-педагогические кадры инновационной России» на 2009 - 2013 годы (ГК №П1509 от 03.09.2009 г.).
2 Деревянов Максим Юрьевич - к.т.н., преподаватель кафедры «Управление и системный анализ в теплоэнергетике».
Лившиц Михаил Юрьевич - д.т.н., профессор кафедры «Управление и системный анализ в теплоэнергетике».
Федорченко Дмитрий Михайлович - аспирант.
dC(x,t) „ d2C(x,t)
-= D--VA te 0,oo , xe 0,oo . (1)
dt 8x2
Граничными условиями, наиболее адекватно отражающими физику переноса углерода из газовой фазы к поверхности детали, являются граничные условия третьего рода [2] :
дС x,t . .
-D-—;-х=о = Р' Ф * ~с Uo ' (2)
ох
dC(x,t) i „ „ i „
—--,Cx,t\x^aD^C0=const, (3)
uX
C(x,0) = /(x). (4)
Здесь C x,t - концентрация углерода, %; t - время, с.; x - глубина слоя, м.; D- коэффициент диффузии, м2/с; (3 - коэффициент массопереноса, м/с; <р / -углеродный потенциал атмосферы, %; f x - начальное распределение углерода, %.
Постановка задачи оптимального управления
В результате процесса вакуумной цементации необходимо получить распределение концентрации углерода C( x, t)| (i ) с наименьшем отклонением от заданного по
эксплуатационным требованиям распределения C * x . Профиль концентрации C * x определяет распределение по глубине детали твердости, предел прочности, износостойкость [1]. Превышение локального отклонения С(х,т-^) от С* х сверх допустимого уровня р х приводит к повышенному трещинообразованию, к снижению твердости, браку [3]. Поэтому в качестве критерия оптимальности принят минимаксный критерий, обеспечивающий, в отличие от среднеквадратичного критерия, абсолютное отклонение результирующего профиля концентрации углерода от заданного.
В ходе диффузионного насыщения нет необходимости и технической возможности обеспечить в конце процесса точную реализацию заданного профиля C * x , так как в производственных условиях всегда на процесс воздействует ряд неконтролируемых возмущений: вариация начального содержания углерода в стали C0, нестабильность газового состава среды, неравномерность потока атмосферы и т.д. Кроме того, заданный профиль C* x может не принадлежать решениям краевой
задачи (1)-(4), что, вообще, свидетельствует о его принципиальной недостижимости. Поэтому в реальных условиях допустимых диапазонов изменения параметров и состояния модели (1)-(4) требуемое результирующее состояние процесса цементации
трансформируется из заданного распределения концентрации C* x в некоторую
область Q - «трубку» допустимых отклонений С* x ± р х , которая характеризуется Чебышевской мерой [4]:
C* x -р x , (5)
где h^ - глубина диффузионного слоя и |р х | > 0 почти везде на отрезке 188
s = max
хе 0,h
хе 0;ИСЛ . Для получения максимально высокой износостойкости величина £
должна быть минимальной. Удовлетворение условия
тт
0<Ф 1 <фшах
шах
хе 0;/г„„
С х -С хЛ
(6)
где тг - время окончания процесса цементации, - не единственное технологическое требование. Для цементации изделий массового производства актуальна проблема достижения максимальной производительности печи для вакуумной цементации при сохранении удовлетворительного качества, т.е. минимума времени цементации
т['' при условии в = сп'1Ш . При этом сокращение времени цементации благоприятно влияет на структуру металла [1].
В случае превышения предельного значения концентрации углерода Стах на цементируемой поверхности возможно образование карбидной сетки, которая характеризуется высокой хрупкостью и снижает износостойкость детали, поэтому вводится ограничение на максимальный уровень концентрации углерода [1,2]:
С хЛ <С7
0,т,
О, да
(7)
Конструкция и ресурсные возможности печи для вакуумной цементации обуславливают ограничения на максимальный расход ацетилена, диссоциация которого определяет максимальный углеродный потенциал сртах атмосферы печи [2,5]:
0<ср Г <Фшах, (8)
рассматриваемый в качестве управляющего воздействия.
Таким образом, для объекта управления (1)-(4) в условиях ограничений (7), (8), технологически обоснованы задачи: быстродействия:
Г,
= тт
0<ф ? ¿ф^
(9)
еП
максимальной точности:
./;;"'" = тт
0<Ф 1 <сртах
тах
ХЕ 0,Й„„
С х,х, -С х
t = х,
(10)
Поставленные задачи представляют собой задачи оптимального управления с подвижным правым концом траектории в бесконечномерной негладкой области
□ НС х,/ : тах
хе 0,к..
I *
С Х,Х; -С X
<8
(11)
допустимых результирующих состоянии для заданной в = е3 или предельно достиг
жимой 8 = с|11П1 точности в области допустимых управлений i -того класса [3,4].
Для решения поставленных задач оптимизации используется альтернансный метод оптимизации (АМО), выгодно отличающийся своей эффективностью при параметрической оптимизации [4].
Применение АМО требует параметризации задач (9) и (10). В работе [5] с помощью принципа максимума Понтрягина установлено, что в условиях ограничений
х
(8) на предельный уровень углеродного потенциала решение задачи быстродействия
(9), как и решение задачи максимальной точности (10), сводится к поиску количества / и длительности Ап' , « = 1,2,...,/ интервалов постоянства управления
ф ' =ф (рис.1).
Р и с. 1. Общий вид управляющего воздействия
Доказательная часть теории АМО рассмотрена в опубликованных работах [4, 6]. Метод основан на использовании специфических свойств результирующих состояний оптимального процесса. Эти свойства позволяют получить достаточное количество уравнений вида
' = ±е,„
(12)
= 0
для определения параметров Д;гг , smin поставленных оптимальных задач (9), (10).
■-Xi
0 > Z3J -
Здесь 9 г,т/ = С х,т/ -С * х у/ ФШаХ - С0 ; г = х\Н; £т=е/ 9ШаХ - С точки экстремума функции 9 г,х- ; г. = 1,/?С7 - граничные точки функции 0 ; / = 1,2,...,./;5 = 7 при в^1 > в = в* > в^п ; 5' = / + 1 при в = в^п;
I
гк] = гЭ] и 2г. ' = ^ Ля' . При этом наименьшая точность в = стш достигается ./=1.2...' ,2=1
на подмножестве одноинтервальных управлений / = 1, а наибольшая точность в = вп^п = !;П1|- > 0 обычно достигается на подмножестве достаточно большого количества интервалов 1=3. Таким образом, решением последовательности задач максимальной точности (9) формируется ряд неравенств:
4 <8„^п <-<8тш <-<8тт = (13)
Очевидно, для решения определяющей системы (12) необходимо получить прямое решение краевой задачи (1)-(4) относительно С х, ? для управления в форме
(см. рис.1):
I
Ф , п = 1,3,5,...
ф / =ф А' = ,п = 1,2,...,/. (14)
п 1 '0, п = 2,4,6,...
Приведенная математическая модель (1)-(4) удовлетворительно описывает на качественном уровне распределение концентрации углерода в цементованном слое. Однако коэффициенты массопереноса (3 и диффузии £) изменяются в достаточно широких пределах в зависимости от конструкции печи, свойств обрабатываемого материала и других факторов. Поэтому для адекватного математического описания конкретного технологического процесса необходимо идентифицировать указанные параметры [3, 4, 5].
Не касаясь здесь процедуры идентификации, рассмотренной в работе [5], отметим, что использование в математической модели постоянных коэффициентов массопереноса (3 и диффузии £) приводит к неудовлетворительному отклонению расчетного профиля распределения углерода от экспериментального. Это объясняется тем, что механизм массопереноса углерода с поверхности детали в глубину на стадиях насыщения и диффузии различен [5, 7] из-за диссоциации ацетилена на поверхности металла в циклах поступления ацетилена и разуглероживания поверхности в циклах пауз. Поэтому в математической модели вакуумной цементации необходимо использовать переменные коэффициенты массопереноса на стадиях насыщения и диффузии (15).
Р ' =
тг-1'т/ , при г =1,3,5,...
(15)
, при г =2,4,6,...
где \ - количество полных циклов в программах цементации; п - номер стадии, нечетный для стадии насыщения и четный для стадии диффузии.
Таким образом, структура управления (14) определяет краевую задачу (1)-(4) как задачу с переменными во времени коэффициентами (3 в соответствии с зависимостью (15).
Решение краевой задачи с помощью метода функций Грина
При многоинтервальном управлении (/>1) (см. рис. 1) получить решение операционным методом в компактном аналитическом виде не представляется возможным из-за сложного начального распределения углерода для каждого из интервалов
управления Аиг .
В этом случае целесообразно воспользоваться методом функций Грина. Решение линейной неоднородной краевой задачи (1)-(4) известно [8] и имеет вид
С(х,0 = ]/©С(3)(х,$,0^+ |р(/)ф(/)С(3)(х,0,/-т)Л. (16)
о о
В работе [7] показано, что использование решения (16) для произвольных начальных условий / \ , выраженного через функцию Грина (¡'31 (х. <;. I) для полубесконечной краевой задачи с граничными условиями третьего рода, в проблемно-ориентированной на оптимизационные процедуры математической модели затруднительно в силу сложности вычисления интегралов свертки для многоинтервального
управления с учетом (14) и (15). Эта сложность во многом обусловлена зависимостью функции Грина в (16) от коэффициента массопереноса р(0 [10,11]:
1
2 -ч/ттБ?
(е
41) Г
+ е
41) Г
-21|в
2) -1
о
4 1>г I)
(17)
который, в свою очередь, непостоянен в соответствии с (15). Это делает выражение (17) малопригодным для многократного обращения к модели в ходе оптимизации в режиме реального времени.
Для упрощения вычисления преобразуем модель объекта (1)-(4) с использованием функции Грина для полубесконечного тела с граничными условиями второго рода:
дС хЛ
дх
= д ( ,
ас(х,о
х=0
дх
= 0, С(х,?) = / х .
X—>00
В этом случае решение краевой задачи примет форму
С(х,0 = ¡ДуО12)(х,^Щ + 1)[д(0О12)(х, 0,/ - т)Л,
(18)
(19)
где функция Грина О ~ х, имеет вид [8, 9]
О2 = 1
4Dt , 4Dt е +е
(20)
Для этого поток д ? сформируем согласно граничным условиям третьего рода
(2):
/ =Р / • ср ( -С I,/
д-0
(21)
В этом случае решения (16) и (19) будут идентичными, что позволяет использовать полученную математическую модель с граничными условиями второго рода, избегая вычислительных трудностей при итерационном процессе моделирования. Структура преобразованной математической модели представлена на рис. 2.
Р и с. 2. Функционально-ориентированная структура математической модели Анализ результатов
С помощью преобразованной математической модели (19), подставленной в определяющую систему уравнений (12), решены задачи максимального быстродействия (9) и максимальной точности (10). 192
Проанализируем зависимость расчетных режимов оптимального управления от режимных и технологических параметров процесса.
На рис. 3 приведены зависимости погрешности науглероживания е
1 2 Р Р
ПШ1' 111П1
л 2 .2
от длительностей интервалов постоянства Д1 , Д2
управляющего воздействия
Ф ? в форме (14) и общего времени процесса т2 при различных коэффициентах
_7 _1
диффузии В. Расчеты проведены для исходных данных: Р1 =2.22-10 м-с ,
р2 =6.9 -10"%/ -с
= 4.1%С.
о.б
о £
и 0.5
!04 I- 0.3
в 0.2
е
о
о
В 0.1
■ I
А 1 2 1 0.1152 0.0338 ППЙ47
V
ч* V Ч,. V ч *. ч 1 /
ч * V К- V »
г б
5200
С 0 1040 2030 3120 4160
Длительность интервала постоянства Л*; ,с.
Рис. 3. Зависимость погрешности науглероживания 8 от длительностей интервалов постоянства Д^ (а), Д22 (б) и от общего времени процесса т22 (в) при различных коэффициентах диффузии В :
1 - £> = 4.17-10"
2 - £> = 6.19-10
-11 2 -1
Л1 ■ с ,
-11 2 -1.
—1 1 9 —1 3 - £> = 8.3-10 Л/ -с
времяпроцес
Анализ полученных зависимостей (см. рис. 3) показывает, что с увеличением коэффициента диффузии В длительности интервалов постоянства Д|2 , Д22 и общее время процесса т22 уменьшаются. Это объясняется тем, что с увеличением коэффициента диффузии В увеличивается скорость проникновения атомов углерода с поверхности в глубину стали, что уменьшает общее время процесса. Величина рассогласования между кривыми 1, 2 и 3 при увеличении коэффициента диффузии В постепенно уменьшается. Для металлов с большим коэффициентом диффузии В время достижения заданной трубки погрешностей (5) становится мало зависимым от величины В.
На рис. 4 приведены зависимости погрешности науглероживания е <
1 2 - р
" ПШ1' 111П1
управляющего воздействия
от длительностей интервалов постоянства Д1 , I
2
Ф I (14) и общего времени процесса т2 при различных максимальных значениях углеродного потенциала Фтах • Расчеты проведены для исходных данных:
в =2.22-10м-с~х, р2 =6.9-10"
м ■ с
£> = 6.19-10 "пл/2-с"1.
О.бг
0.5
0.4
0.3
£ 0-2
0.1
■
1 2 1 0.0829 0.0838 0.0843
Хх /
2' 6
0
1040 2080 3120 4160 5200
Длительность интерьала постоянства Д.
Р и с. 4. Зависимость погрешности науглероживания 8 от длительностей интервалов постоянства Д^ (а), Д22 (б) и от общего вре-
2 / ч
мени процесса х2 (в) при различных максимальных значениях углеродного потенциала Фтах : 1 " Фтах = 3 ■ 5%С ; 2" Фтах =4.1%С ; 3- Фтах = 4.7%С
Из полученных решением оптимальных задач зависимостей (см. рис. 4) видно, что увеличение максимального значения углеродного потенциала фтах приводит к
уменьшению длительности интервала постоянства Д^ . Это объясняется тем, что с увеличением фтах за счёт увеличения расхода ацетилена на поверхности изделия образуется больше атомов углерода, следовательно, для насыщения поверхности детали до максимально допустимого уровня требуется меньше времени, а значит,
Д^ уменьшается. Для диффузии образовавшихся атомов углерода с поверхности в глубь изделия требуется больше времени с увеличением максимального значения углеродного потенциала фтах, так как при одинаковых коэффициентах диффузии
В увеличивается крутизна спада концентрации в глубь детали в конце нечётного
А 2 А 2
интервала Д1 , поэтому длительность интервала постоянства Д2 увеличивается.
БИБЛИОГРАФИЧЕСКИЙ СПИСОК
1. Лахтин Ю.М., Арзамасов Б.Н. Химико-термичекая обработка металлов. - М.: Металлургия, 1985. -216 с.
2. Деревянов М.Ю., Лившиц М.Ю., Липкинд В.Я. Системная оптимизация упрочнения поверхности контактирующих деталей методами ХТО // Вестник Самар. гос. техн. ун-та. Сер. Технические науки. - 2005. - №33. - С. 28-34.
3. Бутковский А.Г. Теория оптимального управления системами с распределенными параметрами. -М.: Наука, 1965. - 474 с.
4. Рапопорт Э.Я. Альтернансный метод в прикладных задачах оптимизации. - М.: Наука, 2000. -336 с.
5. Деревянов М.Ю. Оптимальное управление процессом вакуумной цементации деталей буровых долот: Дис. ... канд. техн. наук. - Самара: Самар. гос. техн. ун-т, 2007. - 184 с.
6. Рапопорт Э.Я. Задача оптимального по быстродействию управления нестационарным процессом теплопроводности // Изв. вуз. Математика. - 1976. - №11. - С. 112.
7. Головской А.Л., Деревянов М.Ю., Ищук А.Г., Лившиц М.Ю., Муратов В.С. Оптимизация вакуумной цементации как объекта управления с распределенными параметрами // Вестник Самар. гос. техн. ун-та. Сер. Технические науки. - 2007. - №1(19). - С. 152-158.
8. Дилигенский Н.В., Темников А.В., Девяткин А.Б., Слесаренко А.П. Современные методы математического моделирования теплопроводности в теплоэнергетике и машиностроении. - Самара: Сам-ГТУ, 1995. - С. 45-54.
9. Бутковский А.Г. Характеристики систем с распределенными параметрами. - М.: Наука, 1979. -224 с.
Статья поступила в редакцию 2 февраля 2010 г.
UDC 621.78.68.3
RATED CONDITIONS ANALYSIS OF OPTIMAL CONTROL OF VACUUM CEMENTATION PROCCESS
M. Yu. Derevyanov, M. Yu. Livshitz, D.M. Fedorchenko
Samara State Technical University
244, Molodogvardeyskaya str., Samara, 443100
The mathematical model of vacuum cementation in is worked out. Conversion of model for using it in optimization problem is made. Relation of rated conditions of optimal control from technological characteristics of process is analyzed.
Keywords: vacuum carburizing, a mathematical model, the coefficient of mass transfer, diffusion coefficient, control action, optimization procedures, alternansny optimization method, the Green's functions method.
Mihail Yu. Derevyanov - Candidate of Technical Sciences, Senior Lecture. Mihail Yu.Livshitz - Doctor of Technical Sciences, Professor. Dmitriy M. Fedorchenko - Postgraduate student.