Очевидно, последний интеграл в (12) равен нулю, а первый не превышает величины
m T ^ _
Е Iа W(T)j-Д f (X (T-x))dx j=1 0
<
m
T
J (AWj)2
0
T _
J(Дf(x))2dx .
V о
Рассматривая периодические воздействия To получаем
A i =
a < EAi;
i=1
iT0 _
J Д W (x)-Дf(X(T-x))dx (i 1)To
T k ’
m
^E
j=1
iT0 d - -
J d- Wj(C ij)-A fj(X(T-x) .[X-(i - 1)To]dx)
(i-1)To dT
Аналогично (9) и (10) каждое слагаемое ограничено
величиной
2k2C'’C2'T2-
Отсюда
T2
Л - ^Т'm • max{C1j C2j}. 2k j
Заключение
(13)
Показано, что как и в рассмотренном раннее случае линейной модели [2], главным членом нелинейного функционала биоэффекта является дозовый функцио-
нал, однако не обязательно позитивный. Формула (13) обобщает (11) на случай m действующих факторов и определяет погрешность от замены при нормировании функционала биоэффекта дозовым функционалом. Отметим, что соответствующий дозовый функционал не является чисто «объективным» показателем внешней среды, поскольку учитывает также некоторые средние характеристики биосистем W0 .
Научная новизна состоит в том, что предлагаемое обобщение линейного функционала биоэффекта на нелинейный многомерный случай является актуальным и важным результатом. Предлагаемые модели более точно отражают реальную ситуацию, связанную с воздействием нескольких факторов на систему (организм человека). Благодаря этому появилась возможность учитывать нелинейность и многофакторность в моделях биовоздействия внешних условий на организм человека, в чем и состоит практическая значимость предложенных моделей типа Гаммерш-тейна.
Литература: 1. Зависимость биоэффектов микроволнового облучения от интенсивности и длительности воздействия/ НИИ гигиены труда и профзаболеваний АМН СССР. М., 1976. 477с. 2. Основы безопасности эргатичес-ких систем / Б.В.Дзюндзюк. К.: УМК ВО. 1990. 56 с.
Поступила в редколлегию 15.02.2006
Рецензент: д-р техн. наук, проф. Самойленко Н.И.
Сердюк Наталья Николаевна, аспирантка кафедры охраны труда ХНУРЭ. Научные интересы: управление условиями труда на рабочем месте. Адрес: Украина, 61166, Харьков, пр.Ленина, 14, тел. 7021-360.
УДК389+517.958:532.5
ПАРАМЕТРИЧЕСКАЯ ИДЕНТИФИКАЦИЯ И МЕТРОЛОГИЧЕСКАЯ АТТЕСТАЦИЯ МАТЕМАТИЧЕСКОЙ МОДЕЛИ ГАЗОПЕРЕКАЧИВАЮЩЕГО АГРЕГАТА
СЕНДЕРОВ О.А.___________________________
Предлагается способ параметрической идентификации и метрологической аттестации математической модели газоперекачивающего агрегата при использовании ее для нахождения оценок степени сжатия, объемной производительности и потребляемой мощности центробежных нагнетателей.
Актуальность. В настоящее время для оценивания параметров режима работы газоперекачивающего агрегата (ГПА) используется множество математических моделей (ММ) различной степени сложности, но все они перед непосредственным применением в задачах нахождения параметров работы технологического оборудования должны быть идентифициро-
ваны, чтобы максимально соответствовать реальности, т.е. быть адекватными. И хотя идентифицированные модели обладают меньшей погрешностью, все равно необходимо знать, с какой точностью получены оценки тех или иных параметров режима работы ГПА с учетом погрешностей прямых измерений и погрешности идентифицированных параметров ММ ГПА.
Цель: разработать эффективный способ проведения параметрической идентификации и метрологической аттестации математической модели газоперекачивающего агрегата.
Задачи. Провести анализ предлагаемой ММ ГПА на предмет возможности максимально точной параметрической идентификации, обоснованно выбрать, какие параметры ММ ГПА и каким способом будут идентифицированы, оценить точность проведенной идентификации. Разработать способ метрологической аттестации ММ ГПА при использовании её в задаче нахождения оценок степени сжатия, объемной производительности и потребляемой мощности центробежных нагнетателей и апробировать предложенный способ на реальных данных.
РИ, 2006, № 1
113
1. Математическая модель ГПА
Рассматриваемая математическая модель ГПА, состоящая из математических моделей центробежного нагнетателя (ЦБН) и газотурбинной установки (ГТУ), применима к стационарным режимам работы ГПА. Математическая модель ЦБН основывается на его паспортных безразмерных характеристиках.
К безразмерным характеристикам ЦБН относятся: коэффициент подачи Kv, коэффициент напора Кн, степень сжатия є и КПД "Л. К размерным характеристикам ЦБН относятся: угловая скорость u или число оборотов n, объемная производительность V, мощность N цбн .
Математическая модель ГПА состоит из следующих уравнений [5,6]:
Р є = Рвых Рвх , (1)
Kv = fKV(KH), (2)
Рпол = fЛпол (KH) , (3)
V KV-Л-d2 • u2 V 4 ’ (4)
m V 'P вх 'H пол N ЦБН - , Рпол ■ Рм (5)
Рвх О = вх ^вх Z RT , вх вх (6)
m-1 KH •U2 = , • Zвх •R • Твх •(s m m -1 -1) =
= -E- • R • (Твых - Твх ), m -1 (7)
m lg(e) m _ 1 ^(Твых /Твх ) , (8)
n-D2 • n U2 = 2— 2 60 ’ (9)
Zвх = fZ(Pвх,Твх,Pc,Ха,Xу) , (10)
fNгтупр (Рза КВД 'P) NГТУ = „ , а-p (11)
а =
T
* пр
^ ТпередКНД ’
(12)
Р
пр
P
(13)
NГТУ _ NЦБН + Nмех.потерь , (14)
где Рвх ,Рвых - давление на входе и выходе ЦБН; Твх, Твых - температура на входе и выходе ЦБН; fKH (Kv ) - функциональная зависимость коэффициента напора от коэффициента подачи из паспорта ЦБН; fлпол (KV) - функциональная зависимость политропического КПД от коэффициента подачи из паспорта ЦБН; Hпол - политропический напор; n - частота вращения ЦБН; D2,u2 - диаметр колеса нагнетателя и угловая скорость его вращения; лм - механический КПД ЦБН; R - постоянная газа; Z^ - коэффициент сжимаемости газа; ^(Рвх ,Твх > Pc, X а , X у) - функциональная зависимость коэффициента сжимаемости газа от параметров газа [10]; Ха ,Xу - молярные доли азота и углерода в природном газе; рс - плотность газа при стандартных условиях; m - показатель политропы; Рза квд - давление газа за КВД ГТУ; fNj-j^ (РзаКВД ’Р) - функциональная зависимость эффективной приведенной мощности ГТУ от приведенного давления за КВД из паспорта ГТУ; N гту -эффективная мощность ГТУ; ^еХЛотерь - механические потери мощности на валу; а - коэффициент приведения по температуре перед КНД; Тпр - температура перед КНД приведения; ТШредКНД - температура перед КНД; Р - коэффициент приведения по атмосферному давлению; Рпр - атмосферное давление приведения; Ратм - атмосферное давление.
Уравнения (1)-(10) относятся к ММ ЦБН, (11)-(13) -к ММ ГТУ, а выражение (14) является уравнением-связкой, которое увязывает ММ ЦБН и ММ ГТУ.
Почему именно так была составлена система, будет понятно из целей, которые мы хотим достичь, проводя идентификацию в п. 3.
Рассмотрим основополагающие утверждения о свойствах безразмерных характеристик.
Утверждение 1. По форме безразмерные характеристики ничем не отличаются от размерных и при соответствующем выборе масштабов совпадают с последними.
Утверждение 2. Безразмерные характеристики не зависят от геометрически размеров машины, числа оборотов и начальных параметров газа (в чем и заключается их преимущество перед размерными и приведенными).
Утверждение 3. В области автомодельности по числу Re характеристики машин однопараметрические. В качестве определяющего параметра был выбран коэффициент напора Kh (более удобен, чем другие в данном случае).
Из последних двух утверждений следует, что
KV = fKV (Kh) и Лпол. = fГ|пол. (Kh) . (15)
114
РИ, 2006, № 1
Утверждение 4. Число маха м и показатель изоэн-тропы kизоэн являются критериями подобия течений, и их изменение влечет изменение характеристик компрессора.
Это утверждение следует из общей теории подобия в аэродинамике. Поэтому на характеристиках компрессора обязательно должны быть указаны начальные параметры газа Рн и Тн, чтобы этими характеристиками можно было пользоваться.
Число маха равно
M = -U2, (16)
a
где a - скорость звука в данном газе, определяемая по формуле
Аппроксимируем данные из таблицы 1 в виде функциональных зависимостей Ку = fKy(KH)
и "Л пол. = fЛпол. (KH) сплайнами [1] вида
q-1 . n
f(x) = jxj + Ёd;
j=0 i=1
(x - Xj)j:q 1 (2q -1)! ,
(18)
где q - 1 - степень полинома, на основе которого строится сплайн; n - количество точек, через которые должен „пройти” сплайн;{(хі,Уі)}і=!, ,п - сами точки;
I x при x > 0 ,
(х) + = [о при x < о ; {аj }j=0,.. ,q—15{di}i=1,..,n - коэффициенты, которые находятся из решения следующей системы n + q линейных уравнений:
a = Vk • z • R • T . (17)
Поскольку окружная скорость u 2 пропорциональна числу маха м, то может использоваться вместо него как определяющий параметр.
Экспериментально проверено, что окружная скорость для центро бежных нагнетателей незначительно влияет на безразмерные характеристики.
2. Построение функциональных зависимостей математической модели ГПА по паспортным данным
Исходными данными для построения модели ЦБН в этом случае являются паспортные газодинамические безразмерные характеристики, построенные по результатам натурных испытаний ЦБН на стенде, проведенных заводом-производителем СМНПО им. Фрунзе. Эти результаты представлены в виде таблицы значений для ЦБН типа ГПА-Ц-16С (табл. 1). Про
Таблица 1
Массив коэффициентов расхода К у
0.0347 0.0422 0.0515 0.0549 0.0594 0.0620 0.0672 0.0749
Массив политропических КПД ц
0.722 0.772 0.813 0.819 0.809 0.799 0.750 0.556
Массив коэффициентов напора Кн
1.037 1.034 0.985 0.945 0.877 0.829 0.731 0.449
данный тип ЦБН известно, что диаметр колеса D2 = 0.83 (м), число маха M = 0.557, номинальная частота вращения n0 = 5200 (об / мин), а характеристики получены при температуре на входе Т = 288 (К), степени сжатия е = 1.44, давлении на выходе
Дж
Рвых = 7.455, газовой постоянной R пр = 506.84-—,
кг • К
коэффициенте сжимаемости = 0 9254 .
n 1
Ёdi(~i)k = 0, k = 0,...,q-1
‘ i=1
f(yi) = уь і = 1, .,n •
(19)
Коэффициенты сплайн-аппроксимации, найденные из решения системы (19), для функциональных зависимостей Ку = fKy (KH) , Лпол. = fипол. (KH) и
NГТУпр = fNГТУпр (Рза КВДпр ) выглядят следующим
образом: {ai }м,...,п+q ={0.0749, 2.5, -2.3102, -0.104796, -0.0188235, -0.0120098, -0.00110544, -0.0257563,
-0.027305}, {bi>i=u.,n+q ={0.556, 16.6667, -15.8299, -0.686735, -0.297059, -0.0612745, -0.291667, -0.187943,0.687943}, {ФЫ,...,р2+q ={3750.96,
842.496, 4.05628, 61.7485, -113.317, 47.5121}.
3. Параметрическая идентификация математической модели ГПА
Цель проведения идентификации в нашем случае заключается в том, чтобы решить две проблемы:
1) получить фактические характеристики ЦБН и ГТУ путем адаптации паспортных к реальным данным;
2) согласовать полезную мощность, выдаваемую ГТУ, с потребляемой мощностью ЦБН.
Поскольку паспорт на ЦБН завод-изготовитель (СМНПО им. Фрунзе) предоставляет нам типовой, а на ГТУ завод-изготовитель (НПО “ Заря”) предоставляет индивидуальный, т. е. для каждой ГТУ свой, то паспортные характеристики ГТУ более точны, чем характеристики ЦБН. Кроме того, мощность ГТУ замерялась экспериментально с помощью электрогенератора, который раскручивала турбина, а это достаточно точный метод измерения мощности.
Под фактическими характеристиками в данном случае будем понимать следующие скорректированные паспортные характеристики: fKy (Kh), ^пол.(Кн), ^гТУпр (Рза КВДпр), полученные сплайн-аппроксимацией паспортов ЦБН и ГТУ, заданных в табличном виде.
РИ, 2006, № 1
115
Предположим, что мы обладаем набором различных реальных стационарных режимов работы конкретного ГПА, каждый из которых характеризуется следующей совокупностью значений измерений:
X _ {Рвх,Рвых,РзаКВД, nЦБН,Твх,Твых}i, І _ 1,..., Nреж с известной погрешностью в виде среднеквадратичных отклонений (СКО) {стр, ст-р, ст n }.
Если мы рассчитаем через модель (1)-(14) параметры
{Рвх,Рвых,РзаКВД,nЦБН,Твх,Твых}i, i _ 1,...,Nреж
при нормальных условиях, то обнаружим, что невязки
чин (СВ) сводится к методу наименьших квадратов [3], следующую целевую функцию:
^еж Nx ~ n 2
F(c) = £ £((fj(X)-Xij)/aj)2 = i=1 j=1
^реж n n 2 n n 2
= £ _P5^i)/CTP) + С^ВыхС^О _ївЬЩіУCTP +
i=1
+ ((РзаКВД (c,Xi) _РзаКВД,i)/CTP)2 +
+ ((nЦБН (c,Xi) _ nЦБН,і)/ ctn ) +
^X _ {Pвх _ рвх ,_Рвых, Рза КВД _ Рза КВД,
ЦБН _ n ЦБН, 'Гвх _ Твх, 'Гвых _ Твых }i,
где i = 1,---, , изменяются в зависимости от режи-
ма, т. е. они не линейно зависят от режима работы ГПА. Обычно нелинейность упомянутой зависимости может быть незначительной, тогда ею можно пренебречь и скорректировать паспортные характеристики fKV (KH), ^пол. (KV),fNppy (РзаКВД) только мультипликативной и аддитивной поправками. В таком случае фактические характеристики будут иметь вид:
^факт (x) = c1 + с2 'fпасп (x), (21)
где С1 - аддитивная поправка, С2 - мультипликативная поправка.
Если же нелинейность зависимости невязок AX от режима X будет значительной, то необходимо подобрать функцию, которая описывала бы эту нелинейность. В общем виде фактические характеристики будут выглядеть:
^факт (x) = g(f пасп (x)) , (22)
здесь g - корректирующая функция.
В нашем случае фактические характеристики будут иметь вид:
f факт fKV f факт
Л пол.
(Kh) = gKV(fK7 (Kh)), (23)
(kh) = g Лпол.^ Oh)), (24)
fNA6o(P?^®^A) gNА66(fNдоб(Р^еаА )) , (25)
где gKV,gЛпол.,gNГТУ - корректирующие функции с параметрами сkv , сг|пол ,сNГТУ соответственно, причем количество параметров сkv , сЛп л , сNrтy может быть разное. Обозначим
С = {С Kv, сЛпол.’6 Nn-У }. (26)
Таким образом, необходимо оценить параметры корректирующих функций С . Для этого составим по методу максимального правдоподобия, который в случае нормально распределенных случайных вели-
+ ((Твх(C,Xi) -Твх,і)/ар)2 +
+ ((Твых(С,Xi) _ТвыхдУ СТТ) +
+ ((ТпередКНД(С, Xi) _ ТпередКНД;і)/ стт)2]-> ™in ,(27)
где
f(c, Xi) = {Рвх (СДО, Рвых (c,Xi), n цбн (С, Xi),
Твх (c Xi), Твых (С Xi), ТпередКНД (С Xi)}
- функциональные зависимости, полученные из модели (1)-(14) путем подстановки выражений для фактических характеристик (23), (24), (25); Nx - количество измеренных данных X; N - количество оцениваемых параметров модели С .
Замечание. Хотя описанный способ идентификации предоставляет возможность ее проведения на нескольких режимах, это не означает, что чем больше режимов мы будем использовать при идентификации, тем более достоверные оценки параметров модели мы получим, потому что ГПА реально работает большую часть времени на одном режиме, а остальные режимы кратковременны и возникают при каких-то определенных обстоятельствах, которые могут больше и не повторяться. В то же время, чем больше различных режимов участвует в идентификации, тем менее адекватной будет модель основному номинальному режиму. Правда, бывает, что изначально в проекте закладывается возможность работы ГПА в двух режимах: то на один газопровод (50 кгс/см2), то на другой газопровод (75 кгс /см2). В последнем случае мы имеем два основных режима, поэтому обязаны использовать при идентификации оба режима.
Проведем идентификацию ГПА-Ц-16С по реальным данным, взятым с КС „Тарутино”, где установлены такого типа ГПА, которые эксплуатируются в одном номинальном режиме. Были выбраны два режима с такими параметрами:
1) ~ = 5127.38(об/мин), РзаКВД = 15.45 (кгс/см2) , Рвх = 44.32 (кгс/см2) , Рвых = 64.43(кгс/см2) ,
Твх = 25.37 (° С), Твых = 57.39 (° С);
116
РИ, 2006, № 1
2) n = 4921.13(об/мин), Рзаквд = 14.59 (кгс/см2), Рвх = 46.22 (кгс/см2) , Рвых = 65.06(кгс/см2) ,
Твх = 26.69 (°С), Твых = 56.03 (°С).
В результате минимизации функционала (27) на одном режиме получаем следующие оценки:
С = (ску, с w сКгТУ } = {0.000221086,1.00637,
- 0.00471175,0.993474, - 65.7814,0.994026} .
Оценим точность применения МНК. Асимптотически несмещенной оценкой СКО [3] является
^ МНК -
F(5)
^ Nреж ■ NX _ Nс
0.517979
(28)
где Nx = 7 - количество измеренных данных X; Nc = 6 - количество оцениваемых параметров модели с .
Определим дисперсии полученных оценок и ковариации. Оценку матрицы асимптотической ковариации [8] будем определять через линеаризацию функциональных зависимостей модели f (c,X;):
к =стМнк • (GT X G)Л (29)
где G - матрица производных размерности
(N реж ■ Nx) х Nc;
етДх;)
G = gCk , 1 = 1,...,N реж, j = 1,...,Nx, k = 1,...,Nc.
(30)
Матрица K размерности Nc х Nc представлена в табл. 2.
Таблица 2
7.275e-12 2.525e-13 -3.588e-13 -2.591e-13 -2.383e-17 -2.624e-13
2.525e-13 8.760e-15 -1.245e-14 -8.989e-15 -8.269e-19 -9.105e-15
-3.588e-13 -1.245e-14 1.769e-14 1.278e-14 1.175e-18 1.294e-14
-2.591e-13 -8.989e-15 1.278e-14 9.224e-15 8.485e-19 9.343e-15
-2.383e-17 -8.269e-19 1.175e-18 8.485e-19 7.804e-23 8.594e-19
-2.624e-13 -9.105e-15 1.294e-14 9.343e-15 8.594e-19 9.463e-15
Известно, что на главной диагонали матрицы K располагаются оценки дисперсий параметров c .
На рис. 1-3 изображены графики паспортных функциональных зависимостей fKV(KH), ^пол (KH) , fNГГУпр (Рза КВДпр) до проведения идентификации (сплошная линия) и после (пунктирная линия), а „жирные” точки - это те паспортные точки, по которым строились (см. п. 2) сплайны.
Значение целевой функции (27) до идентификации составляло 17330817.37, а после - 0.268302. Соответственно рассогласованность Niry _(БТцбн+^ехпотері) до идентификации составляла 404.651 (кВт), а после - -0.0488385 (кВт).
РИ, 2006, № 1
"Л пол
Рис. 1. График функциональной зависимости Л пол = крпол (KH)
Kv
KH
Рис. 2. График функциональной зависимости Kv = fKV(KH)
N ГТУпр
Рис. 3. График функциональной зависимости N ГТУпр = %гтупр (Р2 пр )
4. Метрологическая аттестация математической модели ГПА
Будем проводить метрологическую аттестацию предложенной ММ ГПА (1)-(14) на примере расчета основных параметров режима работы ГПА: степени сжатия є, объемной производительности V и потребляемой мощности нагнетателя Nцбн . Под метрологической аттестацией подразумевается определение статистических хар актеристик оценок параметров режима роботы ГПА.
Нам известны „классы точности” 5 х приборов, используемых для проведения прямых измерений. В
117
нашем случае прямые измерения производятся по температуре, давлению на входе и выходе ЦБН, частоте вращения ТН, давлению за КВД, температуре перед КНД. Класс точности датчиков следующий: по давлению 5р = 0.1%, по температуре 5Т = 0.1%, по оборотам 5n = 0.05%. Максимальная ошибка датчика определяется [2]
д х =
Xmax ' 8Х 100%
(31)
В предположении, что наши результатні прямых измерений представляют собой СВ, распределенные по нормальному закону, СКО выражается
A X
CTX «-f- (32)
Тогда, объединяя формулы (31) и (32), получаем
АX _ Xmax ' 8X
3 _ 3 -100%
(33)
Рассчитав по (33) СКО для датчиков оборотов, давления и температуры, получим
стn = 1.667 об/мин ,
ctp = 0.00333(ММП;, ctt = 0.0333(°C). (34)
4.1. Определение статистических свойств и характеристик степени сжатия
Расчет значения степени сжатия есть косвенное измерение, и величина степени сжатия, равная по определению отношению двух независимых нормальных случайных величин, также является случайной со своим распределением, схожим с распределением Коши:
Рвых (ю)
рвх (ю)
pE (m Е , СТЕ ),
где Рвых(ю) ~ ^Рвых , СТР),РвХ(ю) ~ N(mPbx ,СТ Р^ p Е (m Е, стЕ)- плотность р аспределения с параметр ами mE (центр распределения) и стЕ (масштаб распределения).
Построим функцию распределения по определению Fs(x) = р| рьіх < х| (рис. 4):
Подставив в (36) плотности Рр (zx) и Рр (z), которые заданы как нормальные распределения, получим аналитическое выражение для pЕ (x), которое довольно громоздкое, поэтому здесь не приводится, но по своим свойствам похоже на распределение отношения двух независимых центрированных нормальных случайных величин q(ro) = ^(ю) [4]:
Г|(ю)
ps (x) =---2----
л(ст2 +Ст^x2) ■
У
(37)
є(ю) = x
Рис. 4. Область интегрирования функции распределения Fs (x)
Распределение (37) является частным случаем распределения Коши.
Распределение (36), как и распределение Коши, не имеет конечных моментов, т. е. не сходятся интегралы моментов. Но, учитывая симметричность обоих распределений, для оценки центра распределения mЕ можно использовать оценки выборочного математического ожидания и медианы. Кроме того, мы можем оценить доверительный интервал, задавшись уровнем значимости а и решив уравнение
P(“x1-a/2 ^ x ^ x1-a/2) =
x1-a/2 a
= 2 J p E (x)dx = 1 - — . 0 2
(38)
I p I да zx
Fg(x) = Pj-p^ <xUj Jp(y)p(z)dydz +
І Рвх J 0 -да
0 да
+ J Jp(y)p(z)dydz. (35)
-да zx
Продифференцировав (35), находим
да
pE (x) = Jz • Ррвыx (zx) • pp^ (z) dz -0
0
-I z • РPвыx (zx) • ^вх (z)dz. (36)
—да
Например: если a = 0.005, Рвых = 75 кгс / см2, Рвх = 50 кгс / см2, Стр = 0.03333; то x1_a/2 = 0.22, 0.0033
є = 1.5 ±0.0033 , т.е. 8 = -100% = 0.22% .
В отличие от распределения Коши распределение (36) зависит от оценок математических ожиданий гпРвых и mP случайных величин Рвых (ю) и Рвх (ю), поскольку в общем случае гпРвых — 111 Рвх . Таким образом, оценка медианы £ и её доверительный интервал будет изменяться вместе с режимом работы ГПА.
118
РИ, 2006, № 1
4.2. Определение статистических свойств и характеристик объемной производительности ЦБН
Расчет объемной производительности ЦБН с точки зрения метрологии относится к совокупным измерениям [9], поскольку находится вместе с еще несколькими параметрами режима ГПА из решения системы нелинейных уравнений (2)-(4), (7)—(10).
Поскольку последняя система уравнений довольно громоздкая, чтобы аналитически или хотя бы численно (путем линеаризации) определить закон распределения в предположении, что измеренные данные имеют нормальное распределение, применим метод имитационного моделирования.
Генерируем выборку объемом Nген = 5000 системы независимых случайных величин
X ген (ю) = (n(ro),PBX (ю)
, Рвых (ю)Двх (“) , Твых (ю)Ь распределенных по нормальному закону с математическими ожиданиями в виде значений параметров режима, выбранного для идентификации (см. п.3), и СКО (34).
Вычисляем значения объемной производительности, подставляя наборы СВ в систему (2)-(4), (7)-(10), и
получаем выборку СВ (У(ю)} объемом N ген. Оценим математическое ожидание выборки Ы(У(ю)} = 254.214 (м3/ммин) и сравним с объемной продуктивностью, рассчитанной при подстановке в модель математических ожиданий измеренных величин V(M{X ген (ю)}) = У(Х) = 254.213 (м3/ммин) Видно, что модель значимого смещения не вносит.
Теперь построим гистограмму и проверим с помощью критерия Шапиро-Уилка SW гипотезу о нормальности распределения. Критерий Шапиро-Уилка SW является более мощным [7], чем критерий
согласия X2 или критерий Колмогорова, применительно к проверке гипотезы о нормальном распределении, и его результат не зависит от количества интервалов разбиения выборки для построения гистограммы, как зависит, например, результат критерия согла-2
сия X или критерия Колмогорова.
тем более достоверный результат, чем больше выборка.
Таким образом, мы принимаем гипотезу о нормальном распределении СВ {У(ю)}, а поскольку нормальное распределение симметричное и характеризуется двумя параметрами: математическим ожиданием и
дисперсией стУ, то при уровне значимости а = 0.005 доверительный интервал оценки МО будет симметричным +3сту . Оценка СКО сту будет изменяться
вместе с режимом, поэтому найдем функциональную зависимость СКО от параметров режима (прямых
измерений) ay = f (X, стх). Аналитически определить функцию fCy (х,стх) не представляется возможным из-за громоздкости системы (2)—(4), (7)—(10), но численно, путем линеаризации, это возможно.
Выражаем явно из системы (2)-(4), (7)-(10) объемную производительность как функцию от измеряемых данных и поправочных коэффициентов, которые были найдены при проведении идентификации
V = fV(X,6). (39)
Замечание. Все прямые измерения есть независимые СВ, а поправочные коэффициенты c как минимум коррелированны, что видно из матрицы ковариации (см. табл. 2), и это необходимо учитывать при построении функции fCy (х, стх).
Формула дисперсии от линеаризации функции у = fy (X, с) в общем виде [11] имеет вид:
Dy(X,ctx) (3fy(X,C)/5X;)2-ctX +
i=1 Xl
Nc i ~ Л ~ Л
+ 2EE (Sfy (X, C)/ dCi )(9fy (X, C) / 3cj )Kc c i=1j=1 " j
. (40)
Применив (40) к (39) и подставив данные из примера для идентификации и оценки поправочных коэффициентов, полученных при идентификации, получим
В нашем случае на выборке из 5000 реализаций статистика Шапиро-Уилка SW = 0.999754771, что соответствует вероятности принятия гипотезы о нормальном распределении p = 0.8583 , а менее мощный критерий Шапиро -Фр ансия, но тоже инвариантный к разбиению на интервалы, дает следующие результаты: статистика SF = 0.99978, вероятность принятия гипотезы о нормальном распределении выборки p = 0.99999 . Важно отметить, что оба критерия дают
стУ1 =д/Dv(5~,^5) = 0.0849077 (м3/ммин.) (41)
Для того чтобы узнать, „снизу” или „сверху” мы получили оценку для СКО через линеаризацию, вычисляем оценку дисперсии выборки по известной формуле
1
ТГЧ2
1 N -1 ■ ,
1,исп * i=1
Z (Уі - у)
ctV2 _
0.083(м 3/ммин.)
РИ, 2006, № 1
119
Поскольку сту2 <avl, следовательно, оценка СКО cjv1 сверху, т.е. с небольшим запасом, что нас больше устраивает, нежели оценка снизу.
Таким образом, мы получили оценку с доверительным интервалом
V(X) + 3ст Vl = 254.213 + 0.254723 (м3/ммин),
относительная ошибка составила 8 = 0.100201 %.
4.3. Определение статистических свойств и характеристик потребляемой мощности ЦБН
Расчет потребляемой мощности ЦБН с точки зрения метрологии, также как и расчет объемной производительности, относится к совокупным измерениям, поскольку находится вместе с еще несколькими параметрами режима ГПА из решения системы нелинейных уравнений (2)—(10).
Как уже было сказано в предыдущем пункте, система уравнений (2)-(10) довольно громоздкая, чтобы аналитически или хотя бы численно (путем линеаризации) определить закон распределения в предположении, что измеренные данные имеют нормальное распределение, поэтому таким же образом применим метод имитационного моделирования.
Используем уже сгенерированную выборку объемом N ген = 5000 системы независимых случайных величин :
СВ {Nцбн(ю)} мы принимаем, а поскольку нормальное распределение симметричное и характеризуется двумя параметрами: математическим ожиданием и
дисперсией ст
ЦБН
то при уровне значимости
а
= 0.005 доверительный интервал оценки МО будет
симметричным +3стцбн . Оценка СКО СТцбн будет изменяться вместе с режимом, поэтому найдем функциональную зависимость СКО от параметров режи-
ма (прямых измерений) стn = f (X,ctx) . Ана-
ЦБН °Nцбн
литически определить функцию
стN = f (X,ctx) не представляется возмож-
ЦБН ^ЦБН
ным из-за громоздкости системы (2)-(l 0), но численно, путем линеаризации, это возможно.
Линеаризацию производим аналогично, как и в предыдущем пункте. Применяя (40) и подставляя данные из пример а для идентификации и оценки попр авочных коэффициентов, полученных при идентификации, имеем:
хБНД =<! °n4BH (X>5)=20-5961 (ккВт). (42)
Для того чтобы узнать, „снизу” или „сверху” мы получили оценку для СКО через линеаризацию, вычисляем оценку дисперсии выборки по известной формуле
X ген (ю) = {п(ю) ,рвх (ю)
, Рвых (ю)
,Твх (ю)
Двых (Ю)},
распределенных по нормальному закону с математическими ожиданиями в виде значений параметров режима, выбранного для идентификации (см. п.3), и
СКО (34).
Вычисляем значения потребляемой мощности, подставляя наборы СВ в систему (2)—(10), и получаем выборку СВ {N цбн (ю)} объемом N ген.
Оценим математическое ожидание выборки M{Nцбн(ю)} = 10725.5 (ккВт) и сравним с потребляемой мощностью ЦБН, рассчитанной при подстановке в модель математических ожиданий измеренных величин:
N цбн (M{X ген (ffl)}) = N цбн (X) = 10725.7 (ккВт ). Видно, что модель значимого смещения не вносит.
Рассчитаем статистику Шапиро-Уилка
SW = 0.9996982 , что соответствует вероятности принятия гипотезы о нормальном распределении p = 0.7013, и рассчитаем для контроля статистику Шапиро-Франсия SF = 0.99967, вероятность принятия гипотезы о нормальном распределении выборки p = 0.99628 . По сравнению с применением названных критериев к выборке объемной производительности, вероятность принять гипотезу о нормальном распределении меньше, но достаточно велика и поэтому выдвинутую гипотезу о нормальном распределении
I l ^сп IZ 2
%Щ2 = J ^—!Г £ (NЦБНi -NЦБН) = 20-57(ккВт).
4 Ї ^исп 1 i=l
Поскольку ctn < ctn , следовательно, оцен-
N ЦБН,2 МЦБН,2
ка СКО стм сверху, т.е. с небольшим запасом, что
nHBH,1
нас больше устраивает, нежели оценка снизу.
Таким образом, мы получили оценку с доверительным интервалом
N ЦБН (X) + 3ст NЦБН,l = 10725.7 + 61.7882 (ккВт ), относительная ошибка составила 8 = 0.576076% . Выводы
Разработан и опробован на реальных данных способ параметрической идентификации и метрологической аттестации математической модели газоперекачивающего агрегата при использовании её в задаче нахождения степени сжатия є, объемной производительности ЦБН Q пр и потребляемой мощности ЦБН N цбн .
Практическая ценность предложенного способа заключается в том, что он может быть использован в проведении метрологических аттестаций аналогичных моделей ГПА.
Литература: 1. Лоран П.Ж. Аппроксимация и оптимизация: Пер. с французского Завьялова Ю.С., Звягиной Р.А., Квасова Б.И., Шмырева В.И. /Под ред. Рубинштейна Г.Ш., Яненко Н.Н. М.: Мир, 1975. 405 с. 2. Свешников А.А.
120
РИ, 2006, № l
Основы теории ошибок. Санкт-Петербург: Изд-во Санкт-Петербургского ун-та, 1972. 122 с. 3. ЛинникЮ.В. Метод наименьших квадратов и основы математико-статистической теории обработки наблюдений. М.: Физматгиз, 1962. 344с. 4. Гнеденко Б.В., Курс теории вероятностей. Л.: Физматгиз, 1961. 406 с. 5. Апанасенко А.И., Крившич Н.Г., Федоренко Н.Д. Монтаж, испытания и эксплуатация газоперекачивающих агрегатов в блочно-контейнерном исполнении. Л.: Недра, 1991. 361 с. 6.ШерстюкА.Н. Насосы, вентиляторы и компрессоры. М.: Высшая школа, 1972. 344 с. 7. Хан Г., Шапиро С. Статистические модели в инженерных задачах /Пер. с англ. М.: Мир, 1969. 297 с. 8. Fox J. Nonlinear regression and nonlinear least squares. Appendix to An R and S-PLUS Companion to Applied Regression. January 2002. 9. РМГ 29-99. Государственная
система обеспечения единства измерений. Метрология. Основные термины и определения. М: Изд-во стандартов, 1999. 124 с. 10. Газ природный. Методы расчета физических свойств. Определение коэффициента сжимаемости. ГОСТ 30319.2-96. Минск: ИПК Изд-во стандартов, 1997. 21с. 11. Тевяшев А.Д., Коток В.Б., Сендеров О.А. Метрологическая аттестация математически моделей газоперекачивающего агрегата. // АСУ и приборы автомати-ки.2004. Вип.2.
Поступила в редколлегию 01.03.2006
Рецензент: д-р техн. наук, проф. Гинзбург М.Д.
Сендеров Олег Александрович, аспирант кафедры прикладной математики ХНУРЭ. Научные интересы: системный анализ. Адрес: Украина, 61171, Харьков, Салтов-ское шоссе, 240, тел. 711-27-17.
УДК 621.396: 510.62
РАСПОЗНАВАНИЕ
РАДИОЛОКАЦИОННЫХ ОТМЕТОК ПО СПЕКТРАЛЬНОМУ ИЗОБРАЖЕНИЮ С АДАПТИВНЫМИ ВЕСОВЫМИ КОЭФФИЦИЕНТАМИ
ЖИРНОВ В.В., СОЛОНСКАЯ С.В.__________
Обосновывается возможность использования алгебры предикатов для распознавания радиолокационного спектрального изображения в целях выделения отметок подвижных объектов на фоне дискретных мешающих отражений типа «ангел-эхо». При этом эффективность узнаваемости спектральной картины достигается с помощью адаптивной весовой обработки спектра сигнала, где весовые коэффициенты определяются и постоянно корректируются на основе многообзорного анализа и накопления информации о расположении и поведении радиолокационного спектрального изображения мешающих отражений.
1. Введение
Низкая эффективность систем радиолокационного распознавания и выделения малозаметных целей (МЗЦ) на фоне дискретных мешающих отражений типа « ангел-эхо» (АО) объясняется малой величиной отношения сигнал/помеха [1]. В результате этого ухудшается узнаваемость спектральной картины из-за уменьшения контрастности ее изображения, снижается вероятность выделения полезных сигналов подвижных объектов и увеличивается вероятность их пропуска системами автоматического обнаружения и слежения. Для того чтобы улучшить узнаваемость спектральной картины, предлагается использовать адаптивную весовую обработку спектра сигнала, где весовые коэффициенты определяются и постоянно корректируются на основе многообзорного анализа и накопления информации о расположении и поведении радиолокационного спектрального изображения мешающих отражений.
2. Цель и задачи исследования
Цель - разработка системы распознавания радиолокационного спектрального изображения с адаптивными весовыми коэффициентами, основанная на совмещении сигнального (энергетического) и логического спектрального анализа с адаптацией параметров распознавания к статистике и к типу спектра помех в окрестности анализируемого элемента обработки.
Задачи. Сформулировать особенности использования алгебры предикатов для распознавания радиолокационного спектрального изображения в алгоритмах выделения отметок подвижных объектов на фоне дискретных мешающих отражений типа « ангел-эхо». Оценить возможности улучшения узнаваемости спектральной картины при адаптивной весовой обработке спектра сигнала и определить механизм формирования вектора предикатов Aj,A2,...,Ar с учетом адаптивных весовых коэффициентов спектральных каналов в виде логических уравнений, связывающих предикатные переменные Xj,X2,...,Xr. Обосновать алгоритм интеллектуальной системы идентификации спектральных типов Sj с применением системы предикатных признаков.
3. Предикатное представление спектрального изображения с учетом адаптивных весовых коэффициентов спектральных каналов
Для идентификации спектральной картины радиолокационного объекта по признакам ее изображения используется компараторная схема, которая реализует предикат P(y 1,У2, — ,Ук) = t, соответствующий отношению P к типу изображения спектральной картины (рис. 1). К входам компаратора (элемента распознавания) подключаются своими выходами идентифицируемые информационные процессы fj,f2,... ,fk (входные преобразователи сигнальной информации в спектральных каналах в логические предикатные функции). Изображение спектральной картины преобразуется с помощью адаптивных весовых коэффициентов спектральных каналов, а затем формируется предикатная функция с использованием
РИ, 2006, № 1
121