ОПРЕДЕЛЕНИЕ НЕКОТОРЫХ ПАРАМЕТРОВ ОБЪЕМНЫХ ТЕПЛОВЫХ ИСТОЧНИКОВ В МАССИВЕ ГОРНЫХ ПОРОД ПО РАСПРЕДЕЛЕНИЮ ГРАДИЕНТА ЭЛЕКТРИЧЕСКОГО ПОТЕНЦИАЛА НА ПОВЕРХНОСТИ
УДК 622.8222: 622.02: 519.6
Н.В. Трушникова
Процессы тепло и массообмена в угольных средах (например, на разрезах, в шахтах, местах хранения угля и т. д.) могут привести к образованию источников тепла за счет самонагревания углей. В дальнейшем становится возможным процесс самовозгорания угля, что может привести как к экономическим потерям, так и к экологическим проблемам [1-3].
Физико-химические механизмы возникновения очагов тепловыделения зависят от большого количества факторов: геологических условий залегания пластов и их мощности, строения и химического состава углей, диффузии кислорода из атмосферы и т.д. [1,2,4].
Следует отметить, что единой теории и соответствующей математической модели, описывающей физико-химические процессы при самонагревании и самовозгорании углей, не существует. Поэтому актуальна задача локации объемных источников тепла в угольных средах, решение которой может быть найдено с помощью методов обратных задач с использованием экспериментальных данных. В частности, из решения обратных задач можно получить информацию о геометрических размерах и форме очага тепловыделения [6], что может быть использовано при проведении профилактических работ по тушению очагов самовозгорания на начальных стадиях про цесса [1,2].
Существуют различные экспериментальные методы регистрации очагов тепловыделения. Например, измерение потенциала электрического поля (или его градиента), создаваемого тепловым источником [3,5,7,8], или же определение распре деления температуры на дневной поверхности тепловизорами. Переход от температуры к электрическому потенциалу осуществляется по известной формуле связи между вариациями потенциала электрического поля и полем температуры [3,5,8]. Перечисленные способы имеют свои преимущества и недостатки, что обусловлено физическими процессами и конструктивными особенностями аппаратуры. Поэтому, желательно их "совместное" применение при проведении экспериментов. Это позволит повысить достоверность получаемых экспериментальных данных необхо-
димых при решении обратных задач.
В [10] приведены результаты по восстановлению границы теплового источника для плоской задачи на основе экспериментально измеренного электрического потенциала из работы [5].
В данной работе рассматривается пространственный случай - восстановление поверхности объемного теплового источника в массиве.
На рис. 1. представлена геометрическая иллюстрация решаемой пространственной задачи: источник тепловыделения расположен в области £> с границей 1((р,в) , где (р, 0 - значения соответствующих углов в сферической системе координат {(р, 9, I) , привязанной к декартовой системе (Т, Т]) ; Л(р(х,у) - градиент потенциала (р, измеряемый на дневной поверхности 2=0 в
в пространстве Я3.
области [-М,М]х[-М,М] .
Решение обратной задачи сводится к восстановлению поверхности области по измеряемому градиенту потенциала электрического поля Л(р на дневной поверхности:
\rnir V п)-д(Р(хР'Уо^о)
о=-
где Хо^о^о — кооРДинаты точек измерений.
Предполагается, что область £) - выпуклая и известны координаты, по крайней мере, одной точки, находящейся внутри £), занимаемой тепловым источником. Используя сферическую систему
координат {(р, в, Я) , область О может быть задана в виде:
\( (р;0;Я ) :0<Я<£((р,в),} 0<(р<2л, 0<в<к \
где функция /= 1((р,в) - уравнение поверхности области D.
Искомая функция поверхности очага тепловыделения /= 1((р, в) определяется из следующего нелинейного интегрального уравнения первого рода (см. [6]):
2/г к
¥\( Хо>Уо)= \dcp\Ki х0,у0,(р,в,£( (р,в))йв, о о
(1)
где
¥\ (хо >Уо) = - А ^(х0,у0,0) / к,
к - эффективная плотность электрических зарядов на поверхности
К(*0 >Уо'<Р> Ц(р,0)) = БЫ 6 *
(а2 -Ь2 )^£2((р,в)-2Ы((р,0) + а21
х(а2[НЬ + с(2а2 -ЗЬ2 )] + [Н(а2 -2Ь2 )--cb(5a2 -6Ь2 )]Ц(р,в)+ + с(а2 -Ь2 )£2 ((р,О )) + (Н -ЗЬс)х
х1п
Г(<р,в)-2М((р,е) + а2 + Ц(р,0)-Ъ
а[НЬ + с(2а2-ЗЬ2 )] а2-Ъ2
+
а2 = х02+у02+Н2, Ь = (х0 соз(р + уъ 8т(р)зт6+Н созв, с-созв, а>|£|.
Так как входные данные задачи А(р(хо,уо,0) определяются из эксперимента с погрешностью 5 > 0 и задача (1) некорректна (неустойчива), то для численного решения следует применять ре-гуляризующие алгоритмы.
Численный алгоритм решения уравнения (1), реализованный в виде т- функций в среде МАТЬАВ 6.Х , состоит из следующих шагов.
1) В области
£>! = {О < ср < 2п, 0 < в < я)
вводится сетка:
срк = кк(р, к(р=2л/Ых;
в]=]кв,кв=я/Ы2; к = \,2,...,Ых; у = 1,2,...,
N1 и N2 - число точек разбиения.
2) Двойной интеграл в (1) заменяется ку-батурной формулой
Уг(хоя,Уот) =
N.
N.
(2)
= ^(р^К(х0»,У™,<рк,ОГ1к.)Ив к = 1 у = 1
где 1к]=1((рев^ Х0",Лт-узлы в се-точной области
<х0" <МГ
-М2<у0т <М2};
М] и М2 - характерные размеры области, в которой измеряется градиент потенциала электрического поля
Аср = А <р(х0п, у™ ,0); п,т = -Ы.. N.
3) Система нелинейных уравнений (2) переписывается в виде:
= _ (3)
где F(...) - вектор-столбец из (2); £ = } -
вектор размерности N1 х N2 (искомое решение), а 2 - {у/х (х0", Уот)} — задаваемая из эксперимента правая часть (вектор размерности (2Ы+1) ). Ясно, что должно быть выполнено условие совместности системы: (2Ы+1)2 > N¡N2 , т.е. число уравнений больше числа неизвестных.
4) Система (3) решается с помощью регу-ляризованной модификации итеративного метода Левенберга-Маркварда (см. [9], где приведены также условия сходимости и устойчивости данного метода).
-10 -10
Рис. 2. Распределение градиента потенциала на дневной поверхности г=0
6—
1|-
у.м
-6 -4 •2 0 2 4 6 8
1.М 1 | 1 I
к
_1_ Х-, » ' У У.»
15
•« •6 -4 -2 0 2 4 6 8
Рис. 3. Распределение градиента потенциала Л(р(0,у) (верхний рисунок) и результаты восстановления формы теплового источника при начальном приближении - сфере радиуса Яо=0.5 (нижний рисунок)
Л*0 .у) 11111 11111 у.м 1
а -6 -4 -2 0 2 4 6 е
х.и _1_ О _1_1_1_1_1_ у. и
-а -6 -4 -2 0 2 4 $ 8
Рис. 5. Распределение градиента потенциала А(р(0,у) (верхний рисунок) и результаты восстановления формы теплового источника при начальном приближении - сфере радиуса Яо=2 (нижний рисунок)
Была проведена серия тестовых и модельных расчетов для случая, когда очаг тепловыделения, распределенный в области £>, представлен в виде эллипсоида произвольной формы при различных уровнях среднеквадратической погрешности
входных данных 5.
На рис. 2 дано распределение градиента
электрического потенциала А(р(х0П ,у0"г) (в
относительных единицах) для тестового примера, где тепловой источник имел форму эллипсоида с центром в точке (0, 0, -Н) (£0 =т0=7]0=0),
полуосями ао=3, Ьо=1, Со=2 и углом поворота £0=45° относительно оси т.
Так как решаемое уравнение (1) нелинейное, то существенным фактором, при использовании
Рис. 4. Графики точного (сплошная линия) и восстановленного (пунктирная линия) решений при начальном приближении - сфере радиуса Я0=0.5
Рис. 6. Графики точного (сплошная линия) и восстановленного (пунктирная линия) решений при начальном приближении - сфере радиуса Яо=2
итеративных методов, является выбор начального приближения (нулевой итерации). За начальное приближение выбиралась сфера радиуса Яо .
На рис. 3-6 представлены результаты сравнения точного решения (граница очага источника тепловыделения) и приближенно восстановленной границы в сечениях эллипсоида полуплоскостями (р= 90° и (р=210° (как в декартовой (рис. 3,5), так и в сферической системе координат (рис. 4,6)).
Результаты получены при начальных приближениях - сферах радиусов Яо=0.5 и Яо=2 в прямоугольной области [-7,7] х[-7,7] в тестовом примере.
В верхней части рис. 3,5 даны распределения
градиента электрического потенциала Л(р(х,у) истинного (сплошная линия) и "создаваемого" восстановленным источником тепловыделения
.............. (*:=<>»)
- (Л: =»)
1§г0(к) ю'
Рис. 7. Поведение невязки = - 2:
в логарифмическом масштабе. Начальное и конечное значения функционала невязки: г0(нач)=124.9693, г0(кон)=0.4163-Ю 3 (Я0=0.5), г0(нач)=12.9912, г0(кон)=0.7616-1()-4 (Я0=2)
(пунктирная линия).
В нижней части рис. 3,5 представлены "истинная" граница источника (сплошная линия) и
Лр(0 .у) г 1 1 I I I
/ / / / \\
/ / \ \
/ / \ \
/у \ 1
1 * \\ \\ V V
А л
л // Л
■ // // \ \ » \ V \ -
\\
✓ / \ \ \ \ \ 1
л/
_1_ _1_ .1_1_ _1_1_ _1_
Рис. 8. Распределения градиентов потенциала:
"точного" (-) и "возмущенного" (----) при
3*10%.
г.м 1 1 . 1
чГ\ V
1 « 1 г" !
1 1 1 1 1 1 1_
Рис. 9. Графики начального (. ....), точного
(--) и "восстановленного" (----) решений
при 3*10%.
восстановленная (пунктирная линия) в декартовых координатах при начальных приближениях Яо=0.5 и Яо=2 соответственно.
Наряду с выбором числа итераций, проводился "контроль" счета по невязке
где £ - к-я итерация.
Поведение логарифма невязки Го(к) в зависимости от числа итераций к и начальных приближений Яо=0.5 и Яо=2 дано на рис. 7. Удовлетворительные результаты в большинстве тестов были получены в среднем за 40 итераций.
Выполненный численный анализ позволяет сделать следующий практический вывод: начальное приближение следует выбирать в виде сферы с внутренним источником тепловыделения.
Далее на модельных примерах рассмотрена обратная задача, когда относительная среднеквад-ратическая погрешность измерительной аппаратуры ЗтЮ. На основе тестовых расчетов оценивается уровень погрешности д\ при котором можно получить удовлетворительные результаты по восстановлению формы очага тепловыделения. Ясно, что полученная оценка величины 3>0 , носит приближенный характер (из теории некорректных
задач ддолжна стремиться к нулю). Тем не менее, эта оценка может быть использована при выборе параметров измерительной аппаратуры и при обработке данных в конкретных задачах определения формы очагов тепловыделения.
В модельных задачах в качестве точного решения (так называемого квазирешения по В.К. Иванову [11]) используется эллипсоид с центром в точке (0, 0, -Н) = т0 = Г/0 = 0) , полуосями ао=3, Ьо=1, со=2 и углом поворота £0=45° относительно оси Т.
Подставляя задаваемое таким образом решение в (1) и вычисляя правую часть уравнения > У о )' будем вносить случайные ошибки
в значения (Х()п, у0т ) по формуле:
Ых0я,у0т)=у1(х0я,у0т)-(1 + е-£тп).
п, т=-И... И, где N - число узлов разбиения в сеточной области Зд/ , ~ псевдо-
случайные числа с равномерным законом распределения, £>0 - задаваемая в модельных расчетах величина погрешности измерений у/х.
Среднеквадратическая погрешность ^
(хо" ,Уот )-Ыхоп,Уот)Р =
п т
- с2
У12(хоя.Уот)-
п т пт
Таким образом, в расчетах, задавая различные значения £, можно моделировать разный уровень относительной погрешности входных данных:
ё =----100%
( ИТ^2(х0п,у0т))1/2
п т
Из (4) и очевидного неравенства \%пт\<1
(для всех п,т ), имеем оценку: 3/100<£ .
С помощью разработанных компьютерных программ пеИп.т и ed3.ni были проведены модельные расчеты при различных "размерах" эллипсоида и уровнях погрешности 6.
На рис. 8 изображены распределения градиента потенциала (сплошная линия), соответствующее точному решению и А(р (пунктирная
линия), на которое накладывалось псевдослучайное возмущение с относительной погрешностью
&~10% .Форма источника восстанов-
ленного по возмущенному потенциалу, а также точное решение 1((р, в) (в сечении эллипсоида
полуплоскостями ф1=90° и (р1=210° ) представлены на рис. 9. За начальное приближение бралась сфера радиуса Яо=2 .
На рис. 10 приведены графики поведения невязки Го(к) для соответствующих значений б. Анализ полученных результатов показывает
'О-Г------ т -т--- г -г----------• —г—
Рис. 10. Поведение невязки Уо(к) при разных значениях погрешности 6. удовлетворительное согласование "истинной" формы (границы) очага тепловыделения 1((р, 6) и "восстановленной" 1((р,в).
В заключение, отметим, что применение регу-ляризующих итеративных алгоритмов может быть эффективным средством как для локации очагов тепловыделения, так и при решении других обратных задач: сейсмики, теплопроводности, медицины и т.д.
Автор благодарит профессора Д.В. Алексеева за обсуждение работы и ценные замечания.
СПИСОК ЛИТЕРАТУРЫ
1. Линденау Н.И., Маевская В.М., Крылов В.Ф. Происхождение, профилактика и тушение эндогенных пожаров в угольных шахтах - М.: Недра, 1977, с. 319
2. Глузберг Е.И. Теоретические основы прогноза и профилактики шахтных эндогенных пожаров -М: Недра, 1986, с. 160
3. Тарасов Б.Г., Иванов В.В., Дырдин В.В., Фокин А.Н. Физический контроль массивов горных пород - М.: Недра, 1994, с. 240
4. Агроскин A.A. Физические свойства углей - М.: Металлургиздат, 1961, с. 308
5. Кроль Г. В. Разработка электрометрического способа контроля самонагревания и самовозгорания каменного угля на разрезах. Канд. дисс., ВостНИИ, Кемерово, 1983, с. 195
6. Трушникова Н.В. Об определении формы объемного теплового источника в массиве горных пород по измерениям электрического потенциала на поверхности // Геодинамика и напряженное состояние недр Земли. Труды международной конференции. -Новосибирск, ИГД СО РАН, 2001, с. 165-167
7. Хямяляйнен В.А., Простое С.М., Сыркин П.С. Геоэлектрический контроль разрушения и инъекционного упрочнения горных пород - М.: Недра, 1996, с. 287
8. Алексеев Д.В., Егоров П.В. Механизм формирования квазистационарного электрического поля в массиве горных пород при наличии естественных и техногенных тепловых источников // ФТПРПИ, 1994, № 5, с. 3-7
9. Бакушинский А.Б., Кокурин М.Ю. Итеративные методы решения некорректных операторных уравнений с гладкими операторами - М.: УРСС, 2002, с. 190
10. Трушникова Н.В. Восстановление формы теплового источника в массиве по измеренному электрическому потенциалу на поверхности// Вестн. КузГТУ, № 2, 2005, с. 91 - 93
11. Васин В.В., Агеев АЛ. Некорректные задачи с априорной информацией. -Екатеринбург: УИФ, «Наука», 1993, с. 263
□Автор статьи:
Трушникова Надежда Васильевна - ст. преп. каф. высшей математики