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

Введение. Для решения задач автоматизации технологических процессов дальней транспортировки природного газа характерна ситуация, когда используемая для управления или диагностики математическая модель процесса течения газа не является до конца параметрически определённой. Это касается ряда трудноопределимых констант, связанных с удельным сопротивлением трубопроводов, осреднением скорости потока по сечению, теплообменом их с окружающей средой, осреднением температуры по длине трубопровода при использовании изотермической модели и др. Поэтому аналитически построенную и, следовательно, структурно адекватную процессу математическую модель необходимо идентифицировать параметрически. Это особенно важно для описания динамики магистральных газопроводов (МГП), т. к. посредством исследований установлено, что структура связи переменных состояния этого распределённого объекта такова, что его модель очень чувствительна к некоторым параметрам и имеет физический смысл только в определённых и довольно узких диапазонах [1—3]. Исследования по определению этих диапазонов позволяют уточнить (параметрически идентифицировать) ряд трудноопределимых параметров, а также, опираясь на найденные диапазоны их возможных значений, решить задачу робастного выбора конструктивных и технологических параметров технологических участков МГП [3—6].

Формулировка задачи исследования. Для решения комплекса обозначенных здесь задач необходимо исследовать те аналитические возможности, которые предоставляют получаемые на стадии математического моделирования динамические характеристики распределённого объекта. При этом задачу целесообразно привязать к конкретному технологическому объекту — участку магистрального газопровода. Кроме того, её достаточно локализовать на оценке установившихся значений переменных процесса транспортировки газа, а также на отдельных конструктивных параметрах трубопровода. Необходимо привлечь и исследовать те взаимозависимости переменных состояния и параметров объекта, которые позволяет выделить математическая модель его установившегося состояния, т. е. модель стационарного процесса течения газа по участку МГП [3—6]. Исходная математическая модель задачи исследования. Решение сформулированной выше задачи опирается на математическую модель, характеризующую изотермический процесс1 пе- редачи потоком газа массы вещества и преобразования его механической энергии, описывается математической системой следующего вида:

d p ( x, t )      R . T d q ( x , t )_

.

—

.

d t

Ц

d x

;

( a )

d q ( x, t ) _

—

d t

R . T q ( x , t ) d q ( x , t )

.

—

Ц

p ( x,t )

.

d x

+

R . t ( q ( x , t )

A

.

V

Ц

V

p ( x, t)

—

1 .

dp (x)

—

4 .X. R . T

ц. D

. q ^{( x} ) •

q(x)

p (x).

d x

( b )

Модель (1) построена в базисе (p,q), где p (x,t) — давление в произвольном сечении потока в произвольный момент времени (Па); q (x, t) — удельный (на единицу площади сечения трубопровода) массовый расход газа (кг / с • м2) ; R — универсальная газовая постоянная

( Дж / моль • К ) ; T — абсолютная его температура ( К ); ц — молярная масса газа ( кг / моль ) ; D — диаметр трубопровода; X — удельный (на единицу площади боковой поверхности трубопровода) коэффициент трения газа. Константы модели заданы следующими значениями: D = 1,22 м ; R = 8,3144; р = 0,01604 моль, а X и Т подлежат идентификации.

В связи с поставленной выше задачей решения проблем идентификации в рамках стационарного состояния объекта исследования модель (1) необходимо представить структурно в несколько другом виде, обусловленном допущением о затухании всех переходных процессов. Оно связано с равенством нулю производных по времени всех переменных моделируемого процесса. Таким образом, необходима модель стационарного течения газа по трубопроводу. В базисе модели (1) — ( p , q ) она имеет следующий вид:

dp ( x ) _ 4 .X. R . T . p ( x ) . q ( x ) ² ; (a ) ^dx D . ( R . T . q ( x ) 2 —ц. p ( x ) 2 )'

dqx _ ₀, dx

( b )

Интегрирование уравнения (2 b ) выявляет свойство инвариантности массового расхода к пространственной координате, которое моделируется условием V x е[ 0, L ]^ q ( x ) _ q ₀ _ const. Тогда, интегрируя (2 а ) при q ( x ) _ q ₀, можно получить систему алгебраических уравнений, эквивалентную системе дифференциальных уравнений (2)

q( x ) _ q 0 ;

^{( a} ) ■

_ 0. (b)

-----^ц----_г. [ p ₀² — p ² ) + ln [ — 1 — ^- ^ . x 2 . R . T . q 0^{2 ( У0 P} ) V p 0 J D

Математическая модель решения задачи идентификации. Структура зависимости (3b) позволяет получить явные выражения для температуры потока газа и коэффициента трения потока о стенки трубопровода. При необходимости их также можно преобразовать, для того чтобы в расчётах основываться не на фиксированном базисе, а использовать в качестве базисных переменных те переменные состояния потока, значения которых известны и которые удобны для исследователя.

Для демонстрации возможности такого преобразования и корректности расчёта, основанного на нём, приведём вариант преобразования базиса относительно скорости потока газа. В этом случае выражение (3b) примет вид (4) и будет содержать только одну переменную состояния по-

тока — скорость:

R - Т   fi v 1 । f v 1 4 - X      n -------у • 1 - -0- + ln —--- x = 0.                              (4) 2 -Ц- v 02 (    v 2 J     ( v 0 J    D                                       ()

Из (4) можно выразить X (5) и T (6):

x= D -Г -^Ц -f t - v 1 U f v- 11 ,                   (5) 4 - x L 2 -ц- v o2 ( v 2 J ( v о J J                                  () _ 2 -ц Г 4 -X      . f v )1 v 02 - vv T =    -- x - ln — - —0---7 R [ D        ( v 0 JJ v 2 - v o2                                ()

Объединим (5) и (6) в систему. Чтобы избежать линейной зависимости, используем два набора опытных данных при различных начальных условиях ( v ₀₁ , v ₁ и v ₀₂ , v ₂ ).

X =	D	" R	Т f	v 012 1 - vi?	+ In 1 v ^ 1 J      I v 01 J	,
	4 - x	2 -ц	v 2 1 v 01    \
1	г					(7)
T =	2m-	4 -X	x - In	v 2 1	v 022 - v 22
	R L	D		I v 02 J J	v 22 - v 022

Полученная ММ УМГП, разрешённая относительно температуры потока газа и коэффициента трения потока о стенки трубопровода, позволяет произвести идентификацию коэффициента трения потока газа о стенки трубопровода и средней температуры потока в ситуации, когда значения этих констант либо не замеряются на всём протяжении участка, либо количество замеров не позволяет корректно рассчитать их усреднённые оценки.

Хотя система (7) легко разрешима относительно неизвестных X и T , возникает проблема точности результата, поскольку погрешности в измерениях исходных величин могут привести к неточности итогового результата. Чтобы повысить точность идентификации, целесообразно использовать не два набора опытных данных, а n , т. е. v 0 i и v i при i = 1, n . Построим всевозможные сочетания этих измерений из n по 2 без повторений, разрешим для каждого из сочетаний систему (7) и вычислим усреднённую оценку X и T . С помощью этого, в соответствии с законом больших чисел, удастся получить искомый результат с высокой достоверностью.

Построение комбинаций опытных данных для решения системы (7) позволяет увеличить выборку без увеличения числа опытов. В табл. 1 указаны соответствия числа измерений ( n ) числу возможных сочетаний по два ( C_n ² ).

Таблица 1

Соответствие числа измерений числу сочетаний
n	2	3	4	5	6	7	8	9	10
C n 2	1	3	6	10	15	21	28	36	45

Пример решения задачи идентификации. Рассмотрим пример расчёта, произведённого согласно описанному методу. Исходные данные представлены в табл. 2.

Простейшим вариантом будет являться расчёт на основе лишь двух наборов опытных данных, необходимых и достаточных для однократного решения системы (7) и нахождения весьма приближённой оценки исследуемых констант. В более сложном варианте, т. е. при увеличении количества измерений, точность результатов возрастает. Поскольку этот случай представляет больший интерес, будем производить расчёты на основе 4 наборов опытных данных (измерений начальной и конечной скорости; табл. 3).

Таблица 2

Исходные данные для идентификации

L , м	D , м	кг н , моль
137000	1,195	0,016

Таблица 3

Начальная и конечная скорости

№	v 0 , м	v , м
1	13	81,0026
2	12	29,2682
3	11	20,066
4	10	15,3994

В рассматриваемом случае мы ограничились 4 измерениями. Как видно из табл. 1, эти измерения могут быть использованы в 6 сочетаниях.

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

Таблица 4

Расчёт комбинаций и средних значений

№	Изм. 1	Изм. 2	X	T , м
1	1	2	0,00096086	292,96
2	1	3	0,00096148	293,15
3	1	4	0,00096131	293,091
4	2	3	0,00096349	293,76
5	2	4	0,00096221	293,37
6	3	4	0,00096003	292,701
Средн.	—	—	0,00096156	293,17

Результаты исследований. Таким образом, коэффициент трения идентифицируется значением 0,00096156, температура — значением 293,17. Заметим, что реальное значение коэффициента трения при заданных параметрах 0,000961, а температуры — 293 К. Таким образом, точность идентификации констант при 4 наборах опытных данных оказалась достаточно высокой.

Таблица 5

Квадраты отклонения от среднего значения

№	( x-x ) 2	( t - T ) ‘
1	0,491717 · 10-12	0,0459939
2	0,00687775 · 10-12	0,000633233
3	0,066689 · 10-12	0,00621534
4	3,70436 · 10-12	0,3454766
5	0,421776 · 10-12	0,03927496
6	2,34615 · 10-12	0,2185451

Исследуем полученные результаты с точки зрения математической статистики, чтобы доказать, что использование комбинаций позволяет на основе относительно небольшого числа измерений получить статистически достоверный результат. Вычислим квадраты отклонений X и T от средних значений для каждой из комбинаций (табл. 5).

Рассчитаем среднеквадратическое отклонение с , дисперсию с ² , а также стандартное отклонение s (оценку среднеквадратичного отклонения относительно математического ожидания на основе несмещённой оценки дисперсии) для X и T (табл. 6).

Таблица 6

Статистические и вероятностные показатели

Показатель	X	T
f ( X i -X ) 2 i = 1	7,037567 · 10-12	0,656139
с	1,083018 · 10-6	0,330691
с 2	1,172928 · 10-12	0,109357
s	1,407513 · 10-12	0,362254
Среднее значение	0,00096156	293,17
Относительная ошибка	0,0011263	0,0011280

Заключение. Таким образом, теоретические расчёты и экспериментальная проверка показали, что предложенный метод позволяет достичь значительно более высокой точности результата, чем прямые методы идентификации, т. к. относительная ошибка оценки параметра немногим более 0,1 %.