УДК 556.3
ОПЕРАТИВНОЕ УПРАВЛЕНИЕ РЕЖИМАМИ ЭКСПЛУАТАЦИИ ВОДОЗАБОРОВ ПОДЗЕМНЫХ ВОД
© 2008 г. А.В. Малков
ООО «Нарзан-гидроресурсы», 357743, Ставропольский край, г. Кисловодск, ул. Кирова, 43
«Narzan-hydroresursy», 357743, Stavropolsky region, Kislovodsk, Kirov St., 43
Рассмотрена схема управления, основанная на принципах отрицательной обратной связи, которая широко используется в теории автоматического управления. Схема включает в себя регулятор (пропорциональный или ПИД-регулятор) и модель объекта (гидравлическую или дискретную математическую), обратная связь формируется через сеть наблюдательных скважин. Оптимальность режимов обосновывается критериальной функцией, которая определяется конкретными геолого-гидрогеологическими условиями. Такие схемы могут быть очень полезны в сложных условиях с ограниченными запасами подземных вод.
Ключевые слова: подземные воды, оперативное управление, режим эксплуатации, отрицательная обратная связь.
In the work is considered the conducting plan which is based on the principles of the denied opposite connection. It includes the regulator and the model of the object, the opposite connection is formed through the system of the observing pores. The optimum of regimes is substantiated by the criterial function which is formed by the concrete condition.
Keywords: underground waters, effective control, conditions of exploitation, negative opposite connection.
Эксплуатация водозаборов подземных вод осуществляется в соответствии с технологическими схемами разработки, в которых в обязательном порядке предусматривается обоснование рациональных режимов эксплуатации. Рассчитываются такие режимы на основании данных детальной разведки месторождения и результатов оценки эксплуатационных запасов. Это так называемые долгосрочные прогнозы, которые строятся на принципах программного управления и охватывают периоды упреждения в несколько десятков лет. Основное их назначение - обоснование целесообразности капитальных вложений и возможности покрытия потребности в воде заданного качества за счет данных водоисточников. Оптимальность же рассчитанных режимов довольно условна, поскольку степень изученности месторождения на этом этапе не отличается высокой достоверностью в силу ограниченности сроков изучения реакции объекта на возмущение. Многие факторы, формирующие гидродинамический и гидрогеохимический режимы водоносных систем на ранних стадиях исследований, никак не проявляются и естественно не могут быть учтены в расчетах. Это можно сделать на стадии эксплуатации, которая характеризуется значительно более высокой активностью и длительностью возмущения. Полученные решения по долгосрочному прогнозу систематически корректируются и выдаются в виде рекомендаций текущей (оперативной) эксплуатации каптажных сооружений, охватывающих перспективу до 1 - 3 лет (краткосрочное прогнозирование). В них учитываются фактические и ретроспективные режимы эксплуатации водоносных горизонтов, тенденции в динамике изменения контролируемых параметров, даются рекомендации по дальнейшей эксплуатации каптажей с учетом изменения сезонных нагрузок.
Общая оптимизационная задача формулируется следующим образом [1, 2]. Положим, что состояние системы зависит от N параметров х1, х2,...хп, на которые накладываются ограничения а < х, </З^. Рассматривается некоторая функция Е (критериальная), зависящая от этих параметров F = х2, ...,хп). Требуется найти точку (х00,) = (х1°, х20... хп0} в Ж-мерном пространстве, принадлежащую области V¿, в которой значение критерия оптимальности экстремально:
Е = Е(х0,х0...х0 ) ^ ехГ;
12 п '
а < х. < в.;
< о 0 0, ТТ
Оптимальный режим эксплуатации объекта не означает, что функционирование всех его элементов должно быть оптимальным. Для получения решения вполне достаточно, чтобы один из критериев, наиболее важный для данной задачи и принятый за целевую функцию, был оптимален. Остальные могут представляться в виде системы ограничений или не учитываться вообще.
Если дебит является фиксированным и определяется потребностью в водных ресурсах, то с математической позиции задача выглядит следующим образом:
£ Q = Q ;
^ сум
i = 1
n
F = £ i=1
(
К -
Hi®
H-A id /
\2
^ min;
0 < Q < Q ;
^i ^max'
m . < m. < m , min i max'
n
где ш„
- соответственно минимальное и мак-
симально допустимое значение минерализации подземных вод; тс - среднее значение минерализации; И$) - текущее положение динамического уровня; И^ -допустимое положение динамического уровня в эксплуатационных скважинах; к0 - некоторый коэффициент пропорциональности.
При ^ = 0:
IQг = Q,
сум '
H
id
^ min;
(1)
1=1
п
F = 2 ' = 1 .
0 < Я; < Ятах;
тт,п < < ттах ■
Фактически (1) означает, что оптимальным режимом является такой, который обеспечивает на любой момент времени одинаковое соотношение динамического уровня к предельно допустимому во всех скважинах, а предельное положение в них будет достигнуто одновременно [3].
Решение задач подобного класса известно. Здесь могут применяться методы линейного программирования, например, симплекс-метод [4] или др. Целесообразность использования линейного программирования определяется тем, что оно позволяет найти максимальное или минимальное значение некоторой линейной функции при заданной системе ограничений, представляющих собой систему алгебраических неравенств. Однако для оперативного управления такая постановка малопригодна. Во-первых, неизвестно, какое должно быть положение динамического уровня на текущий момент времени. Ясно, что в диапазоне времени 0 < / < 4 оно будет меньше предельно допустимого, однако какое конкретно - неизвестно. Во-вторых, предельное положение уровня привязывается к расчетным срокам эксплуатации (1,), которое согласно нормативам принимается равным, как правило, 25 или 50 лет, однако соответствует ли это фактическому положению дел, вопрос открытый.
То же самое можно сказать и о предельном уровне. Величина его устанавливается исходя из конкретного геолого-гидрогеологического строения объекта. Ограничения могут быть вызваны как экологическими, так и технологическими причинами: подтягивание некондиционных вод из более глубоких горизонтов, ухудшение микробиологических свойств за счет латерального во-допритока или подпитки из вышележащих обводненных отложений, интенсивное пескование скважин, снижение устойчивости кровли водоносных горизонтов. Установить на ранних стадиях, что более существенно, не всегда возможно. Следует иметь в виду, что с повышением степени изученности объекта величина предельного снижения уровня может меняться.
Рассмотрим задачу в следующей постановке. Положим, что известно предельное положение уровня в любом водозаборном сооружении, никак не привязанное к нормативным срокам эксплуатации объекта. Известна заявленная потребность в воде, а следовательно, и суммарный отбор воды. Независимо от конечных сроков эксплуатации, обеспечим такое распределение общего дебита между скважинами, чтобы
на любой текущий момент времени (() выполнялось условие (1).
Расчеты прогнозного положения уровня могут производиться гидравлическим методом или методом моделирования. В зависимости от вида используемых моделей схемы управления несколько различаются. В гидравлических моделях динамические характеристики объекта не учитываются, и для построения системы управления используется пропорциональный регулятор.
Динамика изменения уровня в любой рассматриваемой скважине определяется выражением [3, 5]:
Hi (t) — Hi0 - Ci x Qi -I (Cij x Qj)
(2)
j—i
где НА) - положение уровня в рассматриваемой скважине на момент времени ; Н,(0) - начальное (статическое) положение уровня; С, - удельное гидравлическое сопротивление (удельное понижение) скважины; С- коэффициент гидравлического взаимодействия }-х скважин с рассматриваемой; Qi Qi -дебит 1-й и}-й скважин соответственно.
Разделим обе части уравнения (1) на предельно допустимое положение уровня:
Hi (t)_, _Но
• — К* — -
Ci x Qi I (C,j X Qj )
H
H
H
id
H
Для определения коэффициента к0 положим на первом этапе вычислений последний член в правой части (2) равным нулю. Это допущение справедливо в том случае, если матрица гидравлических коэффициентов С,- будет диагонально доминантной [6, 7], т.е. коэффициенты, расположенные на диагонали матрицы, будут больше суммы коэффициентов данной строки. Это условие, как показывают расчеты, выполняется при расположении скважин на расстоянии не менее 0,3^0,5 радиуса влияния. Тогда Qв выразится следующим образом:
Qi —
Hi0 - К0 x Hid C,
(3)
Просуммировав дебиты всех скважин, будем
n H,„ n H
иметь QCyM — I HL - Ко xI
i—1 Ci i—1 Cj
ти kn:
k0 —
I—-Q I g Q1
THid
сум
откуда можно наи-
(4)
Далее можно определить оптимальные нагрузки на скважины. Это предварительное распределение, поскольку не учтено взаимодействие скважин. По этой причине приходится прибегать к итерационной процедуре, суть которой заключается в следующем. Определив предварительные нагрузки, производят учет взаимодействия скважин. Полученные значения уровней и дебитов принимаются за начальные, и снова рассчитывается коэффициент пропорциональности к0, но уже при QCуM= 0.
В процессе вычислений по некоторым скважинам могут быть получены дебиты, имеющие отрицательный знак или превышающие возможности скважин ^тах). Таким скважинам присваиваются крайние зна-
ш
шах
чения (0 либо Qmax). Полученный дисбаланс представляется в дальнейших расчетах как А((сум, и процедура решения продолжается по общей схеме. Счет прекращается после выполнения условия:
W(p) =
im im ko — k 0
km+1
< i,
(5)
km.
дИ Л
ду J+
Hi(t)s
Блок сравнения
к„хНш
ь Регулятор Модель
t
Hi(tf
hM
Hi,(t)'
Режимная сеть
, Hi(t)
Рис. 1. Структурная схема управления
Параметры регулятора и закон регулирования выбираются с учетом особенностей объекта управления. Динамика развития понижения уровня в скважине при ее возмущении единичным дебитом формируется по схеме, довольно точно описываемой апериодическим звеном, передаточная функция которого W(p) имеет вид
Tp +1
-pr
(7)
где к0т, к0т+1 - коэффициент пропорциональности соответственно на т и т+1 шаге итерационного процесса; £ - заданная точность.
С переходом к новому шагу прогнозирования задача решается вновь по указанному алгоритму.
Если описание объекта дается в форме дискретной математической модели, то алгоритм принципиально не меняется, но имеются некоторые различия.
Математические модели строятся на принципах сохранения баланса и описываются дифференциальным уравнением в частных производных второго порядка [5, 8]:
_дИ д (, дИ Л д
¡л *-= — I ктх-I +--
Л дх ^ дх ) ду ^
+Ък (Ик - И)+Ьп (Ип - И), (6)
где ц* - водоотдача пласта; ктх, кту - водопроводи-мость разреза по плановым координатам; Н, Нк, Нп -положение уровня в рассматриваемом водоносном горизонте и в смежных, залегающих в кровле и подошве соответственно; Ьк, Ьп - параметры перетекания разделяющих пластов в кровле и подошве; х, у - плановые координаты; / - время от начала эксплуатации объекта. Уравнение дополняется условиями однозначности (начальные и граничные условия).
Решение (6) выполняется в конечных разностях, дискретный аналог дифференциала Л во много раз меньше, чем период прогнозирования, и позволяет учитывать взаимодействие скважин параллельно с решением оптимизационной задачи. Поэтому необходимость в итерационной процедуре отпадает.
Рассогласование желаемого и фактического уровня АБА) = Нфз - Н$) через блок сравнения поступает на регулятор (рис. 1), где вырабатывается управляющее воздействие на модель в виде изменения дебита скважин А(2(). Воздействие измененной нагрузки приводит к перераспределению напоров в водоносном пласте и формированию новой депрессионной воронки, размеры которой контролируется через режимную сеть скважин, образующих обратную связь.
где р - оператор Лапласа; к,Т,т- соответственно коэффициент усиления объекта, постоянная разгона и время задержки. Константы к, Т, т определяются по графикам переходного процесса. Зная их, можно построить амплитудно-фазово-частотные характеристики объекта и выбрать конструктивные параметры регулятора. Для апериодического звена амплитуда А(о) и фаза р(со) рассчитываются из соотношений [6, 7, 9]:
Л(о)= , k
( ) 41+T2 о (8)
р(о ) = -arctg( охт) — охт,
где о - круговая частота.
Гидродинамические процессы могут эффективно управляться пропорциональными (П-регуляторами) или комбинированными высокоточными пропорционально-интегрально-дифференциальными (ПИД-регу-ляторами). Закон регулирования ПИД-регулятора имеет вид [6]
Qi(t) = kp xAS(t) + ±-\hSl(t)dt+ТЯ«*Ш
Tu 0 dt
или в дискретной форме
Qi(t) = kp xASi(t) + -1 £ASi(t)xAt + Tn A(ASl(t)) , (9)
Tu At
где Qi(t) - управляющее воздействие; kp - коэффициент усиления регулятора; AAS(t) - рассогласование управляемого параметра; Ти - время изодрома; Тп -время предварения; At - временной шаг; Ти, Тп, kp -параметры настройки, определяемые по частотным характеристикам.
Рассмотрим пример построения схемы управления на водозаборе артезианских вод «Куюлус» (Р. Казахстан). Водозабор представляет собой линейную систему скважин, вытянутую в меридиональном направлении с шагом расположения скважин 1,5 км. Общее количество скважин - 30. Для простоты расчетов скважины сгруппированы в восемь блоков по 3-4 скважины в каждом. Эксплуатация месторождения сопровождается истощением эксплуатационных запасов и проявляется в систематическом снижении пьезометрических уровней. С начала эксплуатации месторождения понижение достигло более чем 100 м и продолжает нарастать со средней скоростью 0,51,5 м/год. Это основная проблема, вызванная интенсивным техногенным воздействием на гидролитосферу. Требуется построить схему управления гидродинамическим режимом эксплуатации. Исходная информация представлена в табл. 1.
В схеме управления использовалась дискретная математическая модель, описываемая дифференциальным уравнением (6). Фильтрационные параметры и граничные условия пласта установлены по результатам геологоразведочных и опытно-фильтрационных работ, адаптация модели выполнена по ретроспективным данным наблюдений за режимом эксплуатации за период 1964 -1984 гг. Дискретизация по плановым координатам принята Ах = Ау = 5000 м. Всего выделено 30х30 блоков. Расчеты проводились в программном комплексе PLAST.
k
e
Таблица 1
Исходные данные по водозабору «Куюлус»
0,007
№ блока С„ сут/м2 АХ, м АУ, м ß* Нi0, м Ни Q, м3/сут
1 0,0061 5000 5000 0,0002 128 250 5500
2 0,0061 5000 5000 0,0002 139 250 5500
3 0,0061 5000 5000 0,0002 143 250 5500
4 0,0061 5000 5000 0,0002 142 250 5500
5 0,0061 5000 5000 0,0002 145 250 5500
6 0,0061 5000 5000 0,0002 145 250 5500
7 0,0061 5000 5000 0,0002 166 250 5500
8 0,0061 5000 5000 0,0002 166 250 5500
X 44000
Для получения амплитудно-частотных характеристик на модели в одном из блоков был задан постоянный дебит и на различные моменты времени сняты соответствующие понижения уровня. График переходного процесса изображен на рис. 2. Получены следующие параметры:
- коэффициент усиления к = = 0,007;
- время разгона Т = 26 сут;
- время задержки г» 6 сут.
8/0, сут/м2
0,0050
0,0025
0,00
4 /Т=76_!
25
50
75
t, сут
т « 6,0
По Рис- 2.График переходного процесса итудно-фазовово-частотные характеристики (рис. 3), которые принимались в расчетах ПИД-регулятора.
Здесь необходимо соблюдать некоторые условия, связанные с особенностями геологических объектов. Дело в том, что, желая получить максимальное быстродействие регулятора, амплитуду выбирают из расчета, что запас устойчивости по фазе будет не больше (0,2^0,3)л, что соответствует точке (Е) на графике частот и точке (Р) на графике амплитуд (рис. 3). Это хорошо для технических объектов, но для геологических, характеризующихся высокой инерционностью, принимать крайнее значение не следует. В противном случае могут вырабатываться управляющие воздействия, лишенные физического смысла, например, отрицательные дебиты. Для геологических объектов запас устойчивости по фазе следует принимать гораздо больше.
В данном случае параметры ПИД-регулятора рассчитывались для 0,8 я; что соответствует частоте Lg(ю0) = -1,6 (ю0 = 0,025 сут-1). Для этой частоты рассчитывались амплитуда Л(а0) = - к
-¡\ + T2 хю,2
Vi + 262 х 0,0252
= 0,0061 и коэффициент усиления
регулятора kD =1 = —1— = 164.
p k 0,0061
-2,6 -2,4 -2,2 -2 -1,8 -1,6 -1,4 -1,2 -1 L3(°>)
А(ю) -2,0 -2,5 -3,0
^ю) -0,2 -0,4 -0,6 -0,8 —я
-2,6 -2,4 -2,2 -2 -(,8 -(,6 -(,4 -(,2 -( Ш®)
Рис. 3. Амплитудно-фазово-частотные характеристики блока скважин
Прямые 1 и 2 на рис. 3, являющиеся асимптотами интегрального и дифференциального регуляторов, имеют угловой коэффициент, равный единице. Пространственное расположение их симметрично относительно вертикальной плоскости, проходящей через точку А(ю0). Выбрав произвольно положение одной, вторую рассчитывают, в соответствии с законами симметрии. Пересечение их в точках С1 и С2 с горизонтальной линией, проходящей через точку С0, дает проекции на ось частот: Lg(au) = Lg(1/Тu) = -1,8; Lg(юn) = Lg(1/Тn) = -1,4, что позволяет рассчитать время изодрома и предварения: Ти = 63 сут и Тп = 25 сут. С учетом рассчитанных характеристик, ПИД-регулятор имеет вид (при = 30 сут)
0 = 164х АБ(г) + 0,48х 2 АБ((г) + 0,83х А(АБгО.
Оперативное управление начато с момента t = =8000 сут от ввода в эксплуатацию водозабора. Шаг дискретизации времени принят 30 сут. С момента
1 = 8060 сут предусматривается увеличение общего отбора с 44000 до 49000 м3/сут. Расчеты (в сокращении) сведены в табл. 2. Как видно, после каждого шага решения задачи расчетный водоотбор незначительно отличается от заданного. Дисбаланс вызван тем, что ПИД-регулятор выфабагавает воздействие, стремящееся максимально быстро достичь поставленной цели. Учитывается он на следующем шаге в виде Qсум при расчете к0. Этого не наблюдается при использовании пропорционального регулятора, но и выход на заданный режим требует большего времени.
На рис. 4. изображены хронологические графики рекомендуемых нагрузок на скважины, из которых видно, что наращивание суммарного дебита следует осуществлять преимущественно за счет блока скважин № 1. Объясняется это более близким расположением его к области выхода отложений на дневную поверхность, где пласт получает инфильтрационное питание.
Рис. 4. Динамика рекомендуемых дебитов блоков скважин
Таким образом, из вышеизложенного можно сделать следующие выводы.
Под оптимальным (гидродинамическим) режимом эксплуатации водозаборов понимается режим, обеспечивающий на любой момент времени одинаковое соотношение текущего динамического уровня к предельно возможному положению его во всех эксплуатационных скважинах. При создании системы управления целесообразно применение методов теории автоматического управления, математический аппарат которой разработан наиболее полно. При этом в качестве моделей объекта могут использоваться как гидравлические, так и дискретные математические.
Литература
1. Гавич И.К., Теория и практика применения моделирования в гидрогеологии. М., 1980.
2. Методы охраны подземных вод от загрязнения и истощения / Под ред. И.К. Гавич. М., 1985.
3. Бочевер Ф.М., Веригин Н.И. Методическое пособие по расчетам эксплуатационный запасов подземных вод для водоснабжения. М., 1961.
4. Плотников Н.И. Эксплуатационная разведка подземных вод. М., 1973.
5. Шестаков В.М. Динамика подземных вод. М., 1979.
6. Мамсуров А.Х. Киптелая Л.В. Основы автоматики и автоматизации производственных процессов в общественном питании. М., 1980.
7. Першин И.М. Анализ и синтез систем с распределенными параметрами. Пятигорск, 2007.
Поступила в редакцию_
Таблица 2
Расчет режимов эксплуатации водозабора «Куюлус»
№ блока Q, м2/сут Н, м Шаг 1 (t = 8030 сут, Q = 44000 м2/сут)
Kg AS, м XAS AQ Q
1 5500 128 0,587 18,75 18,75 3100 8600
2 5500 139 0,587 7,75 7,75 1281 6781
3 5500 143 0,587 3,75 3,75 620 6120
4 5500 142 0,587 4,75 4,75 785 6285
5 5500 145 0,587 1,75 1,75 289 5789
6 5500 145 0,587 1,75 1,75 289 5789
7 5500 166 0,587 -19,25 -19,25 -3182 2318
8 5500 166 0,587 -19,25 -19,25 -3182 2318
X 44000 0 44000
№ блока Q, м2/сут Н, м Шаг 2 (t = 8060 сут, Q = 49000 м2/сут)
Ко AS, м XAS AQ Q
1 8600 138 0,594 10,44 29,19 1719 10318
2 6781 146 0,594 2,44 10,19 400 7181
3 6120 147 0,594 1,44 5,19 236 6356
4 6285 145 0,594 3,44 8,19 566 6852
5 5789 146 0,594 2,44 4,19 402 6191
6 5789 145 0,594 3,44 5,19 567 6357
7 2318 146 0,594 2,44 -16,81 409 2727
8 2318 144 0,594 2,44 -14,81 740 3058
X 44000 5040 49040
№ блока Q, м2/сут Н, м Шаг 3 (t = 8090 сут, Q = 49000 м2/сут)
Ко AS, м XAS AQ Q
1 10318 145 0,588 1,97 31,16 331 10649
2 7181 152 0,588 -5,03 5,16 -829 6352
3 6356 151 0,588 -4,03 1,16 -665 5691
4 6852 149 0,588 -2,03 6,16 -335 6517
5 6191 149 0,588 -2,03 2,16 -336 5856
6 6357 147 0,588 -0,03 5,16 -5 6351
7 2727 142 0,588 4,97 -11,84 811 3539
8 3058 141 0,588 5,97 -8,84 976 4034
X 49040 -51 48989
№ блока Q, м2/сут Н, м Шаг 7 (t = 8210 сут, Q = 49000 м2/сут)
Ко AS, м XAS AQ Q
1 11020 150 0,600 0,11 33,26 34 11054
2 5535 151 0,600 -0,89 -0,74 -146 5389
3 5198 150 0,600 0,11 -1,74 17 5215
4 6525 150 0,600 0,11 6,26 21 6546
5 5693 150 0,600 0,11 1,26 20 5712
6 6852 150 0,600 0,11 8,26 22 6874
7 3681 150 0,600 0,11 -10,74 15 3695
8 4510 150 0,600 0,11 -5,74 16 4526
X 49015 -1 49014
8. Гидрогеологические расчеты на ЭВМ / Под ред. Р.С. Штенгелова М., 1994.
9. Бутковский А.Г. Характеристика систем с распределенными параметрами. М., 1979.
_6 марта 2008 г.