Синтез прогнозного регулятора для процесса выращивания объемных монокристаллов кремния для солнечных ФЭП Текст научной статьи по специальности «Математика»

УДК 004.421: 548.55

А.П. ОКСАНИЧ, В.Р. ПЕТРЕНКО, С.А. ВОЛОХОВ

СИНТЕЗ ПРОГНОЗНОГО РЕГУЛЯТОРА ДЛЯ ПРОЦЕССА ВЫРАЩИВАНИЯ ОБЪЕМНЫХ МОНОКРИСТАЛЛОВ КРЕМНИЯ ДЛЯ СОЛНЕЧНЫХ ФЭП

На основе использования подхода Бокса-Дженкинса к синтезу моделей стохастических линейных динамических процессов разрабатывается АРМАХ-модель процесса вытягивания монокристаллического слитка, связывающая вариации скорости вытягивания с вариациями диаметра слитка. Полученная модель используется для синтеза оптимального прогнозного регулятора процесса выращивания на этапе вытягивания цилиндрической части слитка. Результаты моделирования работы регулятора подтвердили его работоспособность,

1. Введение

Высокочистый монокристаллический кремний является основным полупроводниковым материалом, используемым в современной электронной промышленности.

Дополнительный интерес к монокристаллическому кремнию появился в связи с решением проблем альтернативных источников энергии. Научные разработки возобновляемых источников энергии в последние десятилетия достигли практической реализации и стали активно внедряться в мировую экономику. Среди таких разработок особое место принадлежит солнечной энергетике, базирующейся на использовании полупроводниковых фотоэлектрических преобразователей (ФЭП).

По прогнозам аналитиков предусматривается, что общая мощность кремниевых фотоэлектрических станций в странах ЕС, в США и в Японии достигнет в 2010 году 9-11 ГВт (из них 7,5-9,2 ГВт будет получено с использованием ФЕП на основе кремния, причем около 40% придется на слитки монокристаллов), а до 2020 года увеличится еще в 20-30 раз, достигнув 15-17% от общих объемов энергопотребления в этих странах [1,2].

Таким образом, сегодня актуальной проблемой является создание новой, стратегически важной для Украины, отрасли производства - фотоэнергетического приборостроения, способной обеспечить потребности страны в альтернативных воспроизводимых источниках энергии и реализовать полный технологический цикл от производства кремния через производство солнечных элементов к конечному продукту солнечной энергетики - фотоэлектрическим модулям, батареям и автономным солнечным электростанциям различного назначения.

Первым этапом при реализации упомянутого технологического цикла в Украине является производство «солнечного» кремния диаметром 200мм. В ООО «Силикон» для этой цели была разработана и внедрена модернизированная установка «РЕДМЕТ-60А», которая базируется на наиболее распространённой технологии получения монокристаллических слитков кремния больших диаметров, реализующей метод Чохральского [3] . Ее эксплуатация показала, что для достижения более высоких значений показателей эффективности процесса выращивания монокристаллических слитков необходимо усовершенствование системы управления.

В работах [4,5] отмечено, что для повышения эффективности контроля и управления процессами промышленного производства С2-8ьмонокристаллов целесообразно использовать стохастический подход к описанию объекта управления (ОУ), который позволяет учитывать естественную неопределенность исходных данных и повышать степень адекватности моделей, используемых при синтезе системы управления. В связи с этим появилась необходимость в разработке более эффективных регуляторов (по сравнению с используемыми ПИД-регуляторами), в частности регуляторов, базирующихся на концепции прогнозного управления. Синтез такого регулятора для процесса выращивания объемных Се^ монокристаллов на установке «РЕДМЕТ-60А» и является целью данной работы.

Метод прогнозного управления был разработан в конце семидесятых годов прошлого века [6]. На сегодня известно значительное количество модификаций этого метода: EPSAC, GPC, MUSMAR, MAC, PFC, QDMC, SOLO [7-13]. В этих алгоритмах используется одна и та же концепция прогнозного управления: наличие внутренней модели, метод отступающего горизонта и вычисление последовательности прогнозных оптимальных сигналов управления. Алгоритмы отличаются применением различных моделей реальной системы, возмущающих воздействий и критериев оптимальности управления.

2. Постановка задачи

В общем виде алгоритм прогнозного управления может быть представлен последовательностью следующих шагов :

1. Предсказание в момент времени t (на основе модели объекта управления) значений

выходной переменной системы yt (t + k), где k=1,...,N1. (При этом выход зависит от будущих управляющих воздействий xt (t+k), k = 0,1,..., N2).

2. Выбор целевой функции управления и оптимизация с ее помощью xt(t+k),k = 0,1,...,N2 .

3. Реализация управления xt = xt (t).

4. Переход в момент времени (t+1) к шагу 1 и повторение шагов 1-4 до достижения цели управления.

Критерий оптимизации, используемый на шаге 2, может быть выбран из достаточно широкого класса критериев. Мы будем применять критерий вида:

N1 N2-1

J(t) = £ [yt (t + k) - rt(t + k)]2q1(k) + £ [Xt(t + k)]2q2(k).

k=1 k=0

Значения штрафов на управление и его ошибку можно изменять с помощью весовых коэффициентов q2(k) и q1(k). Предполагается, что всегда N2<N1 и что Axt(t + k) = 0 для k>N2.

Для обеспечения возможности реализации принципа прогнозного управления возникает необходимость в разработке адекватной модели объекта управления.

3. Разработка математической модели

Для построения математической модели предлагается использовать класс стохастических разностных моделей с дискретным временем [14]. Предположим, что модель передаточной функции

Yt =v(B)Xt + Nt (1)

может быть экономично параметризована в виде

Yt =5-1(B)ra(B)Xt_b +Nt, (2)

где b - чистое запаздывание; B - оператор сдвига назад, т.е. BYt = Yt_1;

5(B) = 1 - 51B - S2B2 -... - 5rBr; ro(B) = ю0 - ro1B - ra2B2 -... - rosBs; Yj, i = t,t- 1,....,t-r и Xj, j = t-b,t-b- 1,...,t-b-s - соответственно отклонения выхода и входа от равновесных состояний; N t - шум, генерируемый некоторым процессом авторег-рессии-проинтегрированного скользящего среднего (АРПСС) и статистически независимый от Xt .

Таким образом, для построения модели необходимо по доступным наблюдениям (X1,Y1), (X2,Y2),..., ( Xn,Yn ) определить оценки параметров r, s, b и начальные оценки параметров 5j, i = s, r и юj, j = 0,s , а также идентифицировать и оценить модель шума.

Заметим, что путем взятия конечных разностей над процессами Xt,Yt модель (1) можно привести к виду

yt = v0xt + v1xt-1 + v2xt-2 +... + (3)

где yt =VdYt,xt =VdXt,nt =VdNt - стационарные процессы с нулевыми средними значениями, d-порядок разности.

Процедура построения модели передаточной функции в соответствии с подходом Бок-са-Дженкинса сводится к выполнению следующих основных этапов:

- получение грубых оценок ^ импульсного отклика в (3) с помощью алгоритма, основанного на выравнивании спектра входа;

- определение оценок г,ё,Ь параметров г, 8, Ь на основе анализа поведения последовательности Vj ;

- вычисление начальных оценок , 1 = 1, г и coj, j = 0,8 на основании оценок V, ГДЬ ;

- определение структуры и начальных оценок параметров модели шума;

- уточнение оценок параметров комбинированной модели;

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

Для определения оценок ^ параметров V модели (3) использовалось соотношение

гар (.¡)

, j=0Д,2,..., (4)

а

где гар (j) - выборочная взаимная корреляционная функция процессов аг и рг, а яр и яа -выборочные стандартные отклонения для этих процессов.

Процесс аг определяется путем подгонки АРСС-модели к процессу хг , т е.

фх (Б)е-1х (Б)х1 = аг , (5)

а процесс рг - результат применения преобразования фх (Б)е-1х (Б) к процессу уг, т.е.

Рг =Фх (Б)е-1х (Б)у|. (6)

Модель (3) при этом может быть представлена в виде

Рг = V (Б)а г + ^ (7)

где ^ = фх (Б)е-1х (Б)пг.

При известных значения параметров г , я , Ь модели (2) можно оценить, используя следующие известные факты [1-8]: для модели вида (2) веса Vj импульсного отклика состоят из Ь нулевых значений У0, V!,..., уь_ , последующих я-г+1 значений уь , уь+1,. , уь+8_г с произвольным поведением (таких значений нет, если я<г) и значений V j при j > Ь + я _ г +1, поведение которых определяется разностным уравнением г-го порядка с г начальными значениями Уь+8 , УЬ+8_1 ^. ^ УЬ+8_Г+1.

Определение начальных оценок параметров 5;, 1 = 1, г и юj, j = осуществляется путем использования следующей системы уравнений:

у j = 0, ] < Ь;

V =51vj_1 +52 vj_2 +••• + 5г vj_г + ю0, ) = Ь;

V =5lVj_l +52Vj_2 +... + 5гVj_г + ю j_ь , j = Ь + 1,Ь + 2,...,Ь + я ; (8)

V =51vj_1 +52 vj_2 +... + 5гvj_г,j > Ь + я .

При известных оценках параметров передаточной функции можно найти оценки пг с помощью соотношения

пг = уг -5_1(Б)й(Б)х1 _ь . (9)

Далее с помощью известных методов идентификации АРСС-процессов определяется структура модели для пг и начальные оценки ее параметров.

На этапе оценивания модели решается задача одновременного эффективного оценивания параметров Ь, 5, та,ф, е ранее идентифицированной модели вида

Уг = 5-1 (Б)ю(Б)хг_ь + ф-1 (Б)е(Б)аг. (10)

Эта задача решается путем минимизации условной суммы квадратов

Б0(Ъ,8,ю,ф,0) = ^ а2(Ъ,8,ю,ф,0 |х0,у0,ао) щ)

г=и+р+1

где и - большее из г и 8+Ъ. Для поиска оценок параметров, минимизирующих (11), использовался хорошо известный метод Марквардта [15]. При этом ковариационная матрица оценок V определяется формулой

V = (XX)-1 а2, (12)

здесь X - регрессионная матрица в линеаризованной модели, вычисленная на последней итерации в процедуре Марквардта, а с2 - оценка остаточной дисперсии.

Диагностическая проверка соответствия комбинированной модели анализируемым данным основана на исследовании поведения остаточных ошибок аг, которые можно определить с помощью следующего соотношения:

аг =0-1(Б)ф(Б) (Уг -8-1(Б)ю(Б)хг-ъ ), (13)

где г > и + 1,и = тах {г, 8+Ъ} и все а^ ) = 1,и + р равны нулю. При этом вычисляется статистика

Р = т Егаа(к) , (14)

к=1

т - количество используемых в расчетах значений аг (обычно т = п - и - р, п - объем

выборки); К - задержка, для которой справедливо, что при к > К автокорреляции пренебрежительно малы.

В [14] указано, что величина Р распределена примерно как %2 с К - р - q степенями свободы. Если р меньше табличного значения для заданного уровня значимости, то принимается гипотеза об адекватности разработанной комбинированной модели.

Ниже приведены результаты основных этапов синтеза модели рассматриваемого класса, связывающей отклонения диаметра растущего кристалла с отклонениями скорости вытягивания. В качестве исходных данных для решения задачи синтеза модели использовались два временных ряда, образованных наблюдаемыми значениями диаметра кристалла (Уг) и скорости вытягивания (Хг). Данные снимались на этапе выращивания цилиндрической части слитка диаметром 200 мм с интервалом 1 мин. на установке „РЕДМЕТ-60А". Фрагменты рядов приведены на рис.1.

204.500 204 000 203.500 203.000 202.500 202.000 201.500 201.000 200.500

Диаметр слитка

1 5 9 1 3 1 7 2 1 25 2 9 3 3 3 7 4 1 45 4 9 5 3 5 7 6 1 6 5 6 9 7 3 77 81 85 89

Время, мин

н и 1.400

2 1.200

?

м 1.000

,ь 0.800

т

с о 0.600

р ок 0.400

кС 0.200

Скорость вытягивания

1 5 9 1 3 1 7 2 1 25 2 9 3 3 3 7 4 1 4 5 49 53 57 61 6 5 6 9 7 3 7 7 8 1 8 5 89

Время, мин

.600

0.000

Рис. 1. Фрагменты исследуемых рядов

Выборочная автокорреляционная и частная автокорреляционная функции входа хг представлены на рис. 2 и 3 соответственно.

Рис. 2. Автокорреляционная функция

Рис. 3. Частная автокорреляционная функция Анализ автокорреляционной и частной автокорреляционной функций позволяет предположить, что вход системы может быть представлен моделью авторегрессии первого порядка. Для оценивания параметров модели входа х1 использовался метод, основанный на применении уравнений Юла - Уоккера для оценивания параметров оператора авторегрессии. Результаты приведены на рис. 4.

+Г Анализ динамических характеристик объектов управления

Данные Преобразование Корреляционный анализ Синтез АРПСС-моделей Синтез передаточных функций Помощь

N1 Н2 Преобразое Разностный АЕТокоЕари Автокорреляц Частная автоко АР

0 0.0055

1 1.345 гоз.за- 1.345 0.0355 0,0037 0 68 -0.34 06838

2 1.342 У 1.342 0.0325 0.0016 0 28 ПК

3 1.347 203.532; 1;.;З'45' 0.0375 ■ШЙН --0.01 -о. Два Щ -0Р6

1 '1 .ЗЖ шаш- 0Д1В5 даг

5 1:35 203.25 1.35 0.0405 -О.0О1;3- •о:гз 0.11

6 1.374 202.702 1..374 0.0645 -0.0003 ■0.1-5- 0.08

7 1Ж 203.11? 1.365 0.0555 а 0 -0.04

е 1.325 203.-324- 1.328 0.0185 0.0005 01 0.04

9 ¡ЙШ в -Д:£$45 адсю |р -014

10 1'.;265 203.1 Г? 1.265 -РР445 0.0003 0 05 Л. 03

11 р.зЗг 203.532; 0.932 -р;з?75 -О.ООРЗ -0.06 -0.13

1?. 1.055 1.224 203.702 1.055 -0.2545 -0.0'01 0'

13 204.001 "024 -0.0355 -0.fl01.J- •о: 24 008

14 кш Л-Ш (Й025 •о.те ИВ

15 1.344 202.909 1.344. 0.0345 0 0 -0.06

16 ли 1.535 202.11:7- 1.395 0.0655 0.0005 0 08 0,07 2}

Процесс АРПСС пхолэ

р с| ц

ШГ пу 1Г-

Оценка дисперсии белого шума

Стандартные отклонения

10.054

входа выхода

10.293

Оценка остаточной дисперсии

10.0736

^•статистика По автокоррел. |Т1

По взаим.корр И 7-2

Число степеней свободы

(20 ГГэ

Рис.4. Исходные данные и результаты оценивания модели входа

Была получена следующая оценка параметра ф1: 0,684. Таким образом, в нашем конкретном случае модель (5) приняла вид

(1 - 0.684Б)хг = аг,

(15)

с эа = 0.29, а модель (6) и в модели (7) - соответственно

Рг = (1 - 0.63В)уг, (16)

^ = (1 - 0.684Б)пг. (17)

Полученные оценки отклика на единичный импульс приведены на рис.5. Анализируя поведение функции отклика на единичный импульс и следуя рекомендациям [14], получаем следующие оценки структурных параметров модели передаточной функции: г = 2, э = 1, Ъ = 1.

Таким образом, модель передаточной функции принимает вид

(1 -5ХБ-52Б2)у( = (Ю0 -ю1Б)х(-1 + п(. (18)

Начальные оценки «левосторонних» параметров (5;) можно получить путем решения системы уравнений

А5 = Ь, (19)

Гиъ+8+;-г 8 +1 > ^

-1п[ х|

где aij =

=иъ+8+;; У = 1,2,..г.

Рис. 5. Оценки отклика на единичный импульс Для получения начальных оценок «правосторонних» параметров (с^) использовались

следующие соотношения:

ю0 =иъ;

(20)

если г > э, то со! = £ 5;иЪ+- и

1=1

ъ+!

если

г<э, то с^ = £5; -иоъ+! для ! = 1,2,...г;

1=1

с^ = £5;ииъ+-иъ+! для ! = г + 1,...э. 1=1

Используя (19), (20), а также полученные ранее оценки структурных параметров модели и оценки отклика на единичный импульс, получаем следующие начальные оценки параметров модели передаточной функции:

ô1 = 0.63; S2 = -0.06; â0 = -1.28; â1 = 0.9.

Построение модели шума основывается на восстановлении последовательности nt путем использования (3) и полученных оценок отклика на единичный импульс, т.е.

nt = Yt - U0 • xt - U1 • xt-i -...- Uk • xt-k,

где значение k должно удовлетворять условию uk+i = 0, i = 1,2,... В нашем случае, как видно из рис. 5, можно взять k = 14 . К полученной последовательности nt осуществляется подгонка АРСС-модели тем же самым способом, с помощью которого мы строили АРСС-модель для последовательности x t .

Оценки автокорреляций и частных автокорреляций шума представлены соответственно на рис. 6 и 7. Их анализ позволил определить структуру модели шума в виде (2, 0, 0). После предварительного оценивания модели шума она приняла следующий вид:

(1 - 1.11B + 0.39B2)nt = at, (21)

sa2 = 0.078.

Уточнение оценок параметров комбинированной модели было выполнено по хорошо известному методу Марквардта [15], позволяющему получить эффективные оценки параметров b , ô , â, Ф , 0 комбинированной модели путем минимизации условной суммы квадратов

S0 (b, ô,â,ф,0)= i a2t (b, ô,â,ф,0 |x0,Y0,a0).

u+p+1

-1П1 Xj

Рис. 6. Выборочная автокорреляционная функция шума

Рис. 7. Выборочная частная автокорреляционная функция шума

После этого комбинированная модель приняла следующий вид:.

2 1

(1 - 1.3В + 0.76В2)у, = (-0.85 + 0.21В)х, , +-- а,

1-1 1 - 0.76В - 0.08В2 1

или после несложных преобразований

(1 -2.06В + 1.67В2 -0.48В3 -0.06В4)у, = -(0.85-0.86В + 0.09В2 + 0.02В3)х(-1 + а,; (22)

8- = 0.073.

В целях диагностической проверки модели (22) с использованием (14) было определено значение статистики Q. Оно оказалось равным 11 при 20 степенях свободы. Критическое значение этой статистики при указанном числе степеней свободы и уровне значимости 0.05 равно 31.41. Таким образом, можно считать справедливой гипотезу об адекватности разработанной модели.

4. Синтез прогнозного регулятора

Сначала опишем процедуру синтеза регулятора в общем виде для класса моделей

Лу(0 = Вх^ -1) + а(0, (23)

где Л = 1 + а1Б + а2Б2 + а3Б3 +••• + а„а8"а ;В = Ь0 + Ь1Б + Ь2Б2 + Ь3Б3 +••• + Ьпь8"ь , Б - оператор сдвига назад, т.е. 8у(1;)=у(1;-1), а затем применим полученные результаты к модели (22)(легко заметить, что она принадлежит классу моделей (23)). Будем также считать, что горизонты прогнозирования и управления равны между собой (N^N2) . Представим (23) в виде

у(,) = Л-х(! -1) + Л-а(1). (24)

Тогда уравнение прогноза с упреждением к будет иметь вид

В 1

уО- + к) = —х0 + к -1) + — а0 + к). (25)

Л Л

Используя диофантово тождество [12], имеющее для модели (24) вид

ЕкЛ = 1 - Б^, (26)

где Ек = е0 + е^ + е2Б2 +••• + е^^; Бк = # + ^ + £2^ +••• + ;

к > 1; пе = к -1; п = па -1; па - порядок полинома А, перепишем уравнение (25) в виде

у0- + к) = ВЕкх(; + к -1) + Бку(,) + Ека(г + к) . (27)

При этом прогнозное значение выхода с упреждением к можно получить с помощью соотношения

у0- + к / г) = ВЕкх(, + к -1) + Бку(,) . (28)

Реальный выход системы может быть записан в виде

у(г + к) = у^ + к/г) + Ека(г + к). (29)

Полагая, что нужно получить прогнозы для некоторого диапазона значений к (от к=1 до к=К ), запишем уравнение (27) в векторной форме, используя подход, предложенный в [16]:

У = ОХ + £ + а, (30)

где

Ут = [у(, + 1),у(, + 2), •••,у(, + К)];

Хт = [х(0, х(г + 1), •• •, х(г + N -1)];

£т = [Цуа) + ^, Б2у(1) + ё2, •••, БкУ(,) + ёк];

ат = [Еха(1 +1), Е2а(, + 2), •• •, + К)];

g0 0 0 0

g1 g0 0 0

G = g2 g1 g0 ■ 0

0

gN-1 gN-2 gN-3 • • g0

gi = hi,i = 0,N-1; H = BEn = ho + h^1 + h2S2 + • + hN_1S(N-1) + •••.

dk = (BEk - (h0 + h1S1 + h2S2 + • + hk-1S(k-1)))x(t + k- 1),k = 1N. Используя функцию цели и векторную форму модели системы (30), запишем критерий оптимальности управления в виде

E{J(t)} = E{(GX + f + a - R)T Q1 (GX + f + a - R) + (XTQ2X)}>min, (31)

где RT = [r(t +1), r(t + 2), •••,r(t + N)] - заданное движение системы; E{x}- математическое ожидание х; Q1 и Q2 - диагональные матрицы размерности NxN с элементами на главных диагоналях соответственно Q1 (1), Q1 (2), • • • ,Q1 (N) и Q2 (1), Q2 (2), • ■ ■ ,Q2 (N) .

Дифференцируя (31) по Х и приравнивая производную нулю, находим оптимальный вектор Х:

X = (GTQ1G + Q2)-1GtQ1(R - f). (32)

Выражение (32) определяет оптимальный прогнозный регулятор.

Приведем некоторые ранее полученные нами [17] полезные формулы, упрощающие расчеты, связанные с синтезом оптимального регулятора.

Коэффициенты полинома Ek (его порядок равен k-1 )можно получить с помощью следующего соотношения:

ek =-1^^ ,k = 1,N-1,eo =1,e: =0Vj<0, (33)

i=1

где na - порядок полинома А; а ai - его коэффициенты.

Коэффициенты полинома Fk, порядок которого nf = na -1, можно вычислить по формуле

fjk = -IIek-iaj+i, j = 0,na -1,k = 1,N,an =0Vn >na . (34)

i=1

Параметры матрицы G определяются с помощью следующего соотношения:

gk-1 = hk-1 = b0ek-1 + b1ek-2 + L + bnb ek-nb-1 ,k = 1,N,ej = 0Vj < (35)

Для вычисления значений величин dk целесообразно использовать следующее соотношение:

dk = (b1ek-1 + b2ek-2 + L + bnb ek-nb )x(t - 1) +

+ (b2ek-1 + Vk-2 + L + bnb ek-nb +1)X(t - 2) +

• (36)

+ bnbek-1x(t - nb), k = 1N, ej = 0Vj < 0. Рассмотрим синтез прогнозного регулятора для объекта, описываемого моделью (22), при N=3. В этом случае имеем:

na = 4; a0 = 1; a1 = -2.06; a2 = 1.67; a3 = -0.48; a4 = -0.06; nb = 3;

b0 = -0.85; b1 = 0.86; b2 = -0.09; b3 = -0.02; nf = 3; ne =2. .

С помощью (33) находим e0 = 1; e1 = 2.06; e2 = 2.574. Используя (34), находим fjk :

= 2.06; Г/ = -1.67; ^ = 0.48; $ = 0.06;

г,! = 2.574; Г2 = -2.96; Г22 = 1.049; Г32 = 0.124;

$3 = 2.341; Г/3 =-3.249;Г3 = 1.359; = 0.154. Далее находим элементы матрицы О, используя (35):

В0 =-0.85;в1 =-0.891; В2 =-0.506. С помощью (36) вычисляем величины

ё1 = 0.86x0 -1) - 0.09x0 - 2) - 0.02x0 - 3); ё2 = 1.682x0 -1) - 0.205x0 - 2) - 0.041x0 - 3); ё3 = 2.008x0 -1) - 0.273x0 - 2) - 0.052x0 - 3).

Формируем компоненты вектора Г Г(1 +1) = ^у« + ^уО -1) + Г!У(1 - 2) + Г3У(1 - 3) + ^ =

2.060у0) - 1.67у0 -1) + 0.48у0 - 2) + 0.06у0 - 3) + 0.86x(t -1) - 0.09x0 - 2) - 0.02x0 - 3); Г0 + 2) = fo2y(t) + ^уО -1) + Г22у0 - 2) + fз2y(t - 3) + а! =

2.574у0) - 2.96y(t -1) + 1.049у0 - 2) + 0.124у0 - 3) + 1.682x(t -1) - 0.205x0 - 2) - 0.041x0 - 3); Г0+3) = Г3у0)+Г3у0 -1)+f23y(t - 2)+fз3y(t - 3)+а3 =

2.341у0) - 3.249y(t -1) + 1.359y(t - 2) + 0.154у0 - 3) + 2.008x(t -1) - 0.273x0 - 2) - 0.052x0 - 3).

В качестве матриц Q1 и Q2 возмем диагональные матрицы, все диагональные элементы которых равны д1=0.5 и д2=0.1 соответственно .

Подставив все полученные величины в (32) и выполнив необходимые матричные операции, получим следующее представление оптимального регулятора в векторной форме:

" x(t) " " -0.77 -0.206 0.037 " " г0 +1) - ^+1)

x(t +1) = 0.602 -0.616 -0.206 r(t + 2) - ^ + 2)

x(t + 2) -0.135 0.602 -0.77 r(t + 3) - ^ + 3)

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

Для x(t) мы имеем

x(t) = -0.77(гО + 1) - ^ +1)) - 0.206(г0 + 2) - ^ + 2)) + 0.03700 + 3) - ^ + 3)) = = -0.77(гО + 1) - (2 060у0) - 1.67у0 -1) + 0.48у0 - 2) + 0.06у0 - 3) + + 0.86x0 -1) - 0.09x0 - 2) - 0.02x0 - 3))) -- 0.026(гО + 2) - (2.574у0) - 2.96у0 -1) +1.049у0 - 2) + 0.124у0 - 3) + + 1.682x0 -1) - 0.205x0 - 2) - 0.041x0 - 3))) + + 0.037(гО + 3) - (2.341у0) - 3.249у0 -1) +1.359у0 - 2) + 0.154у0 - 3) + + 2.008x0 -1) - 0.273x0 - 2) - 0.052x0 - 3))).

Результаты экспериментального исследования работоспособности разработанного прогнозного регулятора представлены на рис.8 и 9.

2.5 2 1.5 1

0.5

>

7 0

~ -0.5 -1 -1.5 -2 -2.5

0 20 40 60 80 100 120 140

t

Рис.8. Работа регулятора при отсутствии шума на выходе Я1=0.5 я2=0.1 ва=0.01 N=3

3 2.5 2 1.5 1

>

7 0.5 Т 0 "" -0.5 -1 -1.5 -2

0 20 40 60 80 100 120 140

t

2

Рис. 9. Работа регулятора при наличии «белого» шума (Ба =0.01) на выходе Выводы

1. Синтезирована адекватная АРМАХ-модель процесса вытягивания монокристаллических слитков в условиях промышленной установки «РЕДМЕТ-60А», которая может быть использована для регулирования диаметра слитка по каналу скорости вытягивания.

2. Усовершенствована методика синтеза оптимального прогнозного регулятора на основе применения соотношений, позволяющих легко алгоритмизировать и автоматизировать процедуру синтеза оптимальных прогнозных регуляторов для установок выращивания слитков «солнечного» кремния.

1 1 1 Л 1 1 .1 -•- ' 1 1 1

1 1 1 1 1 1

1 1 1 1 1 Ъ; 1 1

1 Ту Ц ¿5 Г: V 'I 1 1 1 1 1 4 ___!_ 1 ---1? 1 1 1 " 1 1 1 1 1 1 : | ,,• "

1 1 1 1 1 1 1

1 1 1 1 1 1 1 1

3. На основе использования разработанной ARMAX-модели процесса выращивания выполнен синтез оптимального прогнозного регулятора для установки «РЕДМЕТ-60А», позволяющего учитывать прогнозируемые состояния регулируемого процесса с заданным упреждением, и получены результаты, подтверждающие его необходимое влияние на диаметр выращиваемого слитка кремния.

4. Дальнейшие работы следует посвятить экспериментальному исследованию синтезированного регулятора.

Список литературы: 1. Hirshman W.P., SchmelaM. Silicon Shortage - so what! Photon International, 3, 2006. P. 100-125. 2. Proc. of the 3rd Solar Silicon Conference, April 3, 2006; Munich, Germany. 3. Шашков Ю.М. Выращивание монокристаллов методом вытягивания. М.: Металлургия, 1982. 312 с. 4. Разработка стохастических моделей передаточных функций для системы управления процессом выращивания монокристаллов кремния большого диаметра / А.П. Оксанич, В.Р. Петренко // Вестник Херсонского государственного технического университета. 2002. № 2(15). С. 360-363. 5. Оценивание адекватности стохастических моделей передаточных функций системы управления процессом выращивания монокристаллов кремния / А.П. Оксанич, В.Р. Петренко // Новi технологи. Науковий вюник 1нституту економши та нових технологш. 2004. № 3(6). С. 12-14. 6. Model predictive heuristic control / Richalet J., Rault A., Testud L.J., Papon J. // Automatica. 1978. Vol. 14. P. 413-428. 7. Generalized predictive control. Part 1 and 2. / Clarke D.W., Mohtadi C., Tuffs P.S. // Automatica. 1987. Vol. 23. P. 137-160. 8. Peterka V. Predictor-Based Self-Tuning Control / V. Peterka // Automatica. 1984. Vol. 20, № 1. P. 39-50. 9. Properties of Generalized Predictive Control / D.W. Clarke, C. Mohtadi // Automatica. 1989. № 25. P. 859-875. 10. Analysis and Tuning of Adaptive Generalized Predictive Control / Mcintosh A.R., S.L. Shah, D.G. Fisher // The Canadian Journal of Chemical Engineering. 1991. Vol. 69. P. 97-110. 11. Красовский А.А. Универсальные алгоритмы оптимального управления непрерывными процессами / А.А. Красовский, В.Н. Буко, В.С. Шендрик. М.: Наука, 1977. 272 с. 12. Generalized predictive control / D. Clarke, C. Mohtadi, P. Tu's // Automatica. 1987. Vol. 23. P. 137-160. 13. Allg'owerF, BadgwellTA., Qin J.S., Rawlings J.B., Wright S.J. Advances in Control/ Highlights of ECC'99 // Chapt. 12. Nonlinear Predictive Controls and Moving Horizon Estimation. Springer, London. 1999. P. 391-449. 14. Анализ временных рядов. Прогноз и управление / Бокс Дж., Дженкинс Г. / Под ред. В.Ф. Писаренко. М.: Мир, 1974. 197 с. 15. An Algorithm for least squares estimation of non-linear parameters / D.W. Marqvardt // J. Int. Appl. Math. 1963. № 11. Р. 431-440. 16. Adaptive general predictive controller for nonlinear systems / O.M. Zhu, K. Warwick, J.L. Douce // IEEE Proceedings-D. 1991. Vol. 138, № 1. P. 33-40. 17. Петренко В.Р. Синтез оптимального регулятора с предсказанием для процесса выращивания объемных Cz-Si монокристаллов // Склада системи i процеси. 2008. №2(14). С.64-76.

Поступила в редколлегию 28.08.2009

Оксанич Анатолий Петрович, д-р техн. наук, профессор, ректор, заведующий кафедрой компьютеризированных систем автоматики Кременчугского университета экономики, информационных технологий и управления. Научные интересы: методы и аппаратура контроля структурно-совершенных полупроводниковых монокристаллов. Адрес: Украина, 39600, Кременчуг, ул. Пролетарская, 24/37, тел. (05366) 3-11-44.

Петренко Василий Радиславович, д-р техн. наук, доцент, проректор по научной работе, заведующий кафедрой информатики Кременчугского университета экономики, информационных технологий и управления. Научные интересы: автоматизация процессов управления производством полупроводниковых материалов, информационные технологии. Адрес: Украина, 39600, Кременчуг, ул. Пролетарская, 24/37, тел. (05366) 3-11-44, E-mail: pvr@ient.net.

Волохов Сергей Александрович, генеральный директор ООО «Силикон». Научные интересы: технологии производства полупроводниковых материалов. Адрес: Украина, 36000, Свет-ловодск, ул. Заводская, 3.