Алгоритмы и программное обеспечение для локации источников акустической эмисси в системе геомеханического мониторинга Текст научной статьи по специальности «Электротехника, электронная техника, информационные технологии»

------------------------------------ © И.Ю. Рассказов, А.Ю. Искра,

К.А. Кянно, 2007

УДК 622.831:519.688

И.Ю. Рассказов, А.Ю. Искра, К.А. Кянно

АЛГОРИТМЫ И ПРОГРАММНОЕ ОБЕСПЕЧЕНИЕ ДЛЯ ЛОКАЦИИ ИСТОЧНИКОВ АКУСТИЧЕСКОЙ ЭМИССИ В СИСТЕМЕ ГЕОМЕХАНИЧЕСКОГО МОНИТОРИНГА *

Одной из наиболее актуальных проблем при ведении подземных горных работ на больших глубинах и в сложных горно-геологических условиях является предотвращение опасных проявлений горного давления, для прогнозирования которых в последние годы широко применяются автоматизированные системы контроля горного давления (АСКГД), в том числе основанные на непрерывной регистрации сейсмоакустических событий в массиве горных пород.

На основе выявленных современных мировых тенденций построения измерительно-вычислительных комплексов в ИГД ДВО РАН разработана цифровая автоматизированная сейсмоакустиче-ская система геомеханического мониторинга «Prognoz ADS» , состоящая из подземной и поверхностной частей и включающая в себя цифровые приемные преобразователи-акселерометры (геофоны), подсистему питания/синхронизации геофонов и центр приёма и обработки потока акустических импульсов, управления датчиками и контроля всех узлов и трактов системы на базе персонального компьютера [1, 2].

Работу АСКГД обеспечивают алгоритмы и программное обеспечение, включая разработанный программный комплекс «GeoAcoustics ADS», предназначенный для автоматизации процесса регистрации и обработки данных по источникам акустической эмиссии (ИАЭ), хранения и представления результатов геомехани-ческого мониторинга. Структура «GeoAcoustics ADS» приведена на

* Работа выполнена при финансовой поддержке Фонда развития МП в НТС (гос-контракт № 2749р/5112).

рис. 1. Комплекс условно можно разделить на следующие основные функциональные модули:

Модуль формирования файлов данных о сигналах АЭ, зарегистрированных системой контроля “Прогноз”

Модуль импорта исходных данных

в формат БД

Модуль

селекции

исходных

данных

Модуль расчета параметров АЭ событий (координат и энергии)

Модуль оценки геомехани-ческого состояния массива пород

Г-*С ' аблица

База данных

Архив первичных сигналограмм

К

Таблица зарегистрированных локационных серий

Таблица параметров АЭ событий (обработанных данных)

>

критических ситуации

Таблицы настроек параметров системы контроля

>

Таблица координат геофонов

Таблица конфигураций антенн

Таблица параметров среды

Таблица условий селекции

Таблица пространственных зон

Таблица критериев удароопасности

ти

Модуль настройки параметров системы

Модуль просмотра и анализа данных

Модуль формирования отчетов

Структура программного комплекса «GeoAcoustics ADS» АСКГД

- модуль конфигурации системы и инициализации программного комплекса;

- модуль управления и приёма импульсов АЭ, регистрируемых аппаратурой;

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

- модуль расчета координат и оценки энергии ИАЭ;

- модуль анализа группировки ИАЭ в очаги микроразрушений;

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

Исходная база данных (БД) содержит сведения о дате и времени каждого АЭ-события, временах прихода акустических импульсов на геофоны и амплитудах сигналов, а также ряд других макропараметров АЭ-сигналов (длительность, крутизна переднего фронта, форма сигналограммы и др.).

Эффективность геомеханического мониторинга во многом зависит от точности локации естественных источников АЭ, что в условиях высокой структурной и техногенной нарушенности массива горных пород требует применения специальных алгоритмов, учитывающих существующую анизотропию условий распространения акустических волн, различия в сопряжении геофонов с массивом. В разработанной системе геомеханического мониторинга используются три различных алгоритма локации (два из которых являются итерационными), основанные на определении местоположения источников АЭ по значениям разности времен прихода (РВП) импульсной упругой волны, зарегистрированной не менее, чем 4-5 геофонами из большой группы.

Для удобства дальнейшего анализа следует рассматривать индексацию геофонов порядком прохода через них упругой волны от источника АЭ в модели сферического фронта. Ближайший к источнику АЭ датчик принимается как i=0, для наименьшего РВП i=0 и так далее до самого удаленного датчика с i=n-1 (n>4). Тогда входной набор РВП переменной длины n по результатам регистрации импульсов от любого ИАЭ может быть записан как вектор:

У = (At0, Мг, М2,...Мп_1), (1)

где ti — время регистрации фронта АЭ-импульса датчиком г (г > 0, Д^, > А: ); Д^ = t — tГ, — РВП по отношению к ближай-

^ ’ г+1 г г г 0

шему от ИАЭ геофону с t0; п — количество точек регистрации АЭ-

Н

импульса (“длина” вектора ? ).

Все Ati могут содержать разнородные ошибки (погрешности)

разного веса и различной природы в векторе I измеренных значений РВП. Важнейшим качеством алгоритма локации является способность выявить уровень погрешностей и их распределение по

компонентам t , поскольку они по своей величине и месту (г) разделяются на существенные и малые. Уровень последних таков, что их совокупное влияние в данном конкретном алгоритме приводит к малой пространственной ошибке решения Аг<е (для зоны контроля размером 150x150 м высокая точность обеспечивается при е = 1 м).

Наиболее сложным аспектом обработки АЭ-импульсов, регистрируемых на сети из множества геофонов, является выделение однозначного подмножества импульсов, относящихся к одному источнику АЭ, из общей совокупности регистрируемых в определенный интервал времени акустических сигналов, включая помеховые. При априорной равномерности распределения источников АЭ по зоне контроля и за ее пределами и ненаправленной диаграмме чувствительности геофонов, периферийная часть объема массива горных пород, окружающая зону контроля, продуцирует в десятки раз больше сигналов АЭ, чем компактный центр зоны контроля [3]. Дистанция до таких источников АЭ от ближайших геофонов может достигать 200 м и более. При этом упругий импульс в столь удаленной точке регистрации утрачивает крутизну переднего фронта, приобретает дополнительную изрезанность огибающей и несколько увеличивает собственную длительность [4, 5], а время его прохождения через “антенное поле” геофонов имеет максимальное значение, как правило, превышающее 60 мс. Все они часто играют роль помех.

Поэтому в современных системах геоакустического контроля необходимо реализовывать многоточечные алгоритмы локации, которые обеспечивают устойчивый расчет координат источников АЭ и позволяют обнаруживать ложные импульсы и ошибочные отсчеты РВП во входных наборах данных. В реализованном в системе «Prognoz АDS» программном обеспечении предусмотрена обра-

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

Рассмотрим апробированные в шахтных условиях удароопасных рудников Дальнего Востока алгоритмы локации источников АЭ, дающие устойчивые решения и реализованные в программном комплексе «GeoAcoustics ADS».

Алгоритм локации с использованием метода наименьших квадратов.

Координаты источника АЭ в МГП определяются решением системы уравнений для расстояний от центра источника до каждого из n датчиков:

1^-= R + di, (i>0) (2)

Г - КI = R;

I 0 u\ ’

где n — количество датчиков, допускаемых в расчет координат ИАЭ; i — индексы геофонов, не превышающие n; Г — координаты i-го датчика (ra = |Га|, где а — любой индекс); К — координаты центра ИАЭ (xu, yu, zu); Г0 — координаты датчика, первым принявшим импульс АЭ; R — расстояние от ИАЭ до ближайшего к нему датчика; v — скорость распространения (продольной) упругой волны от ИАЭ; Ati — измеренная РВП на i-ом датчике,

Ati = tt —10, (At0 = 0); di — дистанция запаздывания волны от источника до i-го датчика, dt = v • Ati.

Следует особо подчеркнуть, что в условиях анизотропии свойств реального горного массива (даже в пределах зон 150-200 м в поперечнике), равенство расстояний между двумя различными парами точек не означает равенства соответствующих времен пробега упругих волн между ними. После несложных стандартных преобразований в (2), получаем следующую систему линейных уравнений, записанных в матрично-векторной форме:

А •К b, (3)

где А — матрица размерности 3х(п-1), т.е. столбец вида Ai = (р-р), г=1^(п-1); Ь — вектор-столбец размерности п-1 с компонентами:

(г2 -г2 -и2)

Ь = ^---------------------0-!— - R • И, ¡=Щп-\). (4)

Такая система линейных уравнений (3) при п>4 избыточна, а для измеренных (дискретных) приближенных РВП в строгом смысле несовместна. Используя метод наименьших квадратов (МНК) для наилучшего приближенного решения системы, получим, что такое решение (3) должно удовлетворять следующей системе линейных уравнений:

Ат • А • р = Ат • Ь , (5)

где Ат — транспонированная матрица А.

Представим векторы искомого решения С , а также Ь каждый в виде подобных линейных комбинаций двух различных векторов:

Ь= Ь+ R •Ь; Ь = Ь+ Я • И. (6)

и е ’ V '

Тогда, подставляя из (6), систему уравнений (5) можно переписать так:

Ат • А • р+ R • Ат • А • С = Ат • р + R • Ат • И, (7)

где С — вектор-столбец с компонентами вида с= Ь

= (г2 - г02 - И2)/2; И — вектор-столбец с компонентами (-И).

В результате из (7) получаем две подсистемы линейных уравнений:

Ат• А• р= Ат • С , (8)

Ат• А • С= Ат • И. (9)

После решения этих систем и определения векторов Ги¥е

необходимо определить параметр R, который определяется из условия:

R = \Р + R • р - р| = |р + R • р|, (10)

где Г = Г - Г0.

Применив очевидные преобразования, из (10) получаем квадратное уравнение относительно R>0:

(г2 -1)• R2 + 2R• (Р,Р) + г2 = 0. (11)

Решив данное уравнение по R и подставив его значение в Ь, получаем искомое наилучшее решение Г системы уравнений (3). Единственной мерой отклонения этого решения от неизвестной ло-цируемой точки (центра ИАЭ) является, в частности, удельная временная квадратичная невязка по РВП, которую необходимо вычислить для координат и Ги 2 , полученных от каждого из корней (11) - RI и R2, в целях исключения неоднозначности решения, по формуле:

I (*. - )2

' -, (12)

ЕГ = 1

п -1

где 8їі — вычисленная (модельная) РВП для геофона с измеренным Аіі по известной V от источника с координатами Г.

Увеличение расстояния между сопоставляемыми по (12) точками от 1 м до 10 м в направлении от ИАЭ всегда приводит к многократному увеличению значения Етг при достаточно большой антенне по составу (более 7 датчиков). Следовательно, практически приемлемое решение задачи локации Г лежит вблизи глобального

минимума значения Етг для полученного от источника АЭ вектора Р

і (при пространственной ошибке в пределах 0,5-3 м). При применении в АСКГД малых единиц измерения времени для целочисленного представления величин РВП с высокой точностью (микросекунды), сумма в (12) служит быстрой и простой мерой пространственного отклонения координат любой точки в массиве горных пород от лоцируемого ИАЭ.

Если полученные оба координатных решения системы уравнений (5) по значению невязки Етг близки, то локационная неоднозначность сохраняется (в особенности при сопоставимости величин

R1 и R2) и для выбора окончательных координат фактического ИАЭ остается только детальный анализ амплитуд и форм импульсных сигналограмм.

Достоинствами реализации данного алгоритма является высокая скорость вычислений координат (при п< 15) и возможность наращивания размерности исходной системы уравнений для улучшения решения за счет добавления компонент РВП, измеренных по временам вступления волн 5-типа. Учет обоих р- и 5-типов волн в алгоритме эквивалентен 2-кратному измерению каждого РВП и приводит к повышению достоверности локации источника АЭ.

Довольно детальное рассмотрение влияния скоростной анизотропии МГП (которая, в принципе, сводится к вариациям значений РВП) для 5-точечных антенн на высокоточные решения систем уравнений, подобных (5), содержится в работе [6].

Применительно к погрешностям измерения РВП во входном векторе (1), следует отметить, что, при регистрации фронта упругого импульса во множестве точек массива (свыше 7-8), часть полученных РВП имеет наибольшие погрешности, искажающие координатное решение системы (5) в максимальной степени, а их исключение (когда это возможно) из состава системы улучшает качество локационного решения. Обычно наибольшую погрешность в измерения вносят самые отдаленные от ИАЭ датчики вследствие ослабления амплитуды импульса АЭ и искажения формы его огибающей, хотя импульсные локальные помехи сказываются в любой точке приема АЭ-сигнала. Их устранение из входного набора данных требует анализа и сопоставления формы и спектров всех импульсов, вошедших в локационный пакет данных.

Алгоритм локации методом семи точек применяется в сейсмоакустике более 20 лет в силу его логической простоты и нетребовательности к вычислительным возможностям компьютера (в нем задействованы только простые арифметические операции). Алгоритм локации 7 точек по существу заключается в следующем итерационном поиске решения Г, соответствующего измеренному

С и

РВП-набору I . В качестве исходной (центральной) точки может

быть взят ближайший к ИАЭ геофон, либо точка в его окрестностях. Вокруг центральной точки хс, ус, zc (ее индекс 7=0) анализируются еще 6 пробных точек со смещенными коор-динатами

(/'=1,2,..,6): xc-h, уа zc; х+h, уа zc; ха yc-h, zc; ха у+h, zc; ха уа zc-h;

хс, ус, zc+h, где h — величина смещения, первоначально равная 10 %

поперечника антенны.

Показателем большего соответствия пробной /-точки к измеренным значениям РВП является её “фактор качества”, которым служит квадратичная временная невязка по РВП этой точки q, вычисляемая по формуле (15) ниже. В используемой изотропной модели массива горных пород заметно меньшей невязке соответствует более близкое к координатам ИАЭ положение пробной точки (источника) при n>4 и оптимальной конфигурации антенны. Формулы системы измерений и расчётов используемых параметров:

AtiJ = tij - t0 j , tij = dij /V , (13)

dij = д/(xi - xj )2 + (yi - yj )2 + (zi - zj )2 , (14)

qj =g (A t- St,)2, (15)

i=0

где j =0,..,6 — индексы пробных точек, индекс i соответствует РВП У t

из t ; li/ — время пробега упругой волны от j-го “источника” до датчика i; t0 — время пробега волны от j-го “источника” до ближайшего датчика; d/ — рассчитанное расстояние между датчиком i и j-м “источником” с координатами х, yj, z/; Ati — измеренная РВП от АЭ-источника до датчика i; StjJ — рассчитанная РВП между приходом волны к i-му датчику (i<n) по отношению к ближайшему датчику (t0 ) для “источника” j.

Суть алгоритма заключается в подборе наилучшей точки из семи пробных по критерию минимизации qj. В ходе итерационного процесса решения локационной задачи в качестве центральной на следующем шаге выбирается та из пробных точек, для которой qj оказался минимальным. Если наименьший qj оказался у центральной точки, то пространственный шаг смещения h уменьшается по выбранному правилу (обычно — вдвое). Критерий назначения следующей центральной точки сохраняется в ходе всех итераций. Финальное координатное решение (XU,YU,ZU) получают, останавливая итерации на стадии, когда достигается h<0,5 м, т.е. обеспечивается

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

В усовершенствованной нами последней версии алгоритм отличается от известной “классической” реализации выбором начальных условий. В качестве исходной пробной точки - к данному 1

t - выбирается точка (в пределах зоны контроля), грубо минимизирующая невязку Егг из (12), а не ближайший к источнику АЭ датчик, как прежде. После такой модификации алгоритму требуется, как правило, вдвое меньше итераций (па-раметр к) для получения решения. Наш подход часто позволяет сократить время поиска (или улучшить точность расчета) коор-динат ИАЭ и одновременно сделать вывод о качестве входного РВП-набора, в частности — о

его физической реализуемости. Следовательно, помеховые (и лож-

1

ные) импульсы могут быть от-селектированы из t перед заключительным этапом получения решения, снижая его координатную погрешность (отклонение) до 2-3 раз.

Поскольку не все датчики равноценны по качеству измерений принятого АЭ-импульса, в целях учета этого фактора, нами в формулу (15) введены весовые коэффициенты, рассчитываемые по значениям РВП, удаленности от пробной точки и другим макропараметрам, а окончательное решение локационной задачи, таким образом, минимизирует функционал вида:

где В, — весовой коэффициент для 7-го геофона (0<В7<50 - оценочно); £ £ — радиус-векторы 7-го геофона и “источника” соответственно; гк —I £ - £ — длина радиус-вектора от 7-го геофона до

“источника”; рі — \Лґі - бґ7| — ошибка сличения РВП с модель-

ной для 7-го геофона; є7 0 — оценка акустической энергии от АЭ-

источника в 7-м геофоне.

При выполнении итераций от к = 1 до нескольких десятков “центральная точка” проходит пространственную траекторию, схо-

(17)

7—0

дящуюся к решению Г, при этом Bi является монотонноубывающей функцией от rik и R (при возрастании аргумента) и, например, линейно-возрастающей функцией от pi. Вид функции должен быть достаточно простым для расчета, поскольку в ходе решения значения Bi вычисляются не менее 100 раз (до 200-300).

Введение системы весов для ошибок по РВП предполагает ослабление влияния на локационный результат датчиков, не обеспечивших качественные измерения параметров АЭ (и даже позволяет реализовать при k>10 их полное исключение).

Более детальный анализ макропараметров (РВП, площадь и

т.д.) каждого АЭ-импульса и корреляции их с макропараметрами

первоначального импульса АЭ в данном АЭ-событии позволяет

1

отбраковывать из обрабатываемого вектора t заведомо ложные (помеховые) импульсы до локационного расчета.

В программе «GeoAcoustics ADS» также реализован и апробирован алгоритм локации методом “сферических итераций”. Его особенность заключается в следующем. Четыре локационных геофона вне плоскости располагаются в точках с координатами Xt, Y, Zt, где i = 0,1,2,3 (следовательно, здесь уже отсутствуют избыточные РВП).

Данный алгоритм основан на приближенном расчете расстояния R0 от ИАЭ до датчика, принявшего в такой антенне АЭ-сигнал первым (с его минимальной оценкой) и уточнении этой величины в сторону увеличения за 3-5 итерации. Избегая изложения довольно громоздкого вывода методики решения для Xu,Yu,Zu, приведем конечные соображения об эффективности такой локации.

На практике, для 5-точечных антенн, доля наиболее успешных решений по критерию минимальности значения Ert у данного алгоритма при сопоставлении с двумя вышеописанными не превышает 3 % от общего числа решений при локации естественных АЭ-событий последовательно каждым из них. (Для сравнения: первый алгоритм обеспечивает лучшее решение в 60-75 % локаций ИАЭ, в зависимости от помеховой обстановки и конфигурации антенн.).

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

Для них далее необходимо провести усреднение или выборку наилучшего. В результате неоправданно увеличивается совокупное время расчета искомых координат фактического источника АЭ, превышая временные затраты первых двух алгоритмов.

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

1. В процессе геоакустического контроля установлено, что наилучшее качество локации при сложной морфологической структуре значительных пространственных зон требует многофакторного анализа входных наборов данных по макропараметрам и итерационной оптимизации состава используемых в расчетах подмножеств регистрируемых импульсов АЭ.

2. Реализованные в прикладном ПО «GEOACUSTICS ADS» алгоритмы локации ИАЭ в массиве горных пород минимизируют квадратичную временную невязку точек-решений, что обеспечивает допустимую погрешность определения координат ИАЭ (не более 3 % от базы контроля).

3. Для неизбежных погрешностей определения параметров АЭ-импульсов в реальном массиве горных пород наилучшую устойчивость локационных решений имеет алгоритм с применением МНК, а наилучшую адаптивность к условиям локации, - алгоритм локации “методом 7 точек”, одновременно упрощающий дискриминацию неприемлемых для локации импульсов АЭ (при наличии избыточности данных).

4. Локация при параллельном отборе данных по макропараметрам позволяет уточнять координатные решения и энергетические оценки АЭ-источников.

---------------------------------------- СПИСОК ЛИТЕРАТУРЫ

1. Совершенствование систем геоакустического мониторинга при ведении подземных горных работ / А.Ю. Искра, И.Ю. Рассказов, Г.А. Калинов, Ю.И. Болотин // Горный журнал. - 2006. - № 6. - С. 72-77.

2. Акустический измерительно-вычислительный комплекс для геомеханического мониторинга массива пород при ведении горных работ / Г.А. Калинов, И.Ю. Рассказов, А.Ю. Искра, Д.А. Куликов, К.О. Харитонов // Физическая акустика. Распространение и дифракция волн. Геоакустика: Сборник трудов XVI сессии Российского акустического общества. Т. 1. - М.: ГЕОС, 2005. - С. 351-354.

3. Куксенко В.С. Возможности акустической эмиссии в прогнозировании разрушения горных пород // Системы контроля горного давления. - М.: ИПКОН

АН СССР, 1989. - С. 5-22.

4. Ямщиков В.С. Волновые процессы в массиве горных пород. - М.: Недра,

5. Искра А.Ю., Рассказов И.Ю. Развитие средств акустического контроля проявлений горного давления // Горный журнал. - 2007. - № 1. - С. 85-88.

6. Локация источников акустической эмиссии / К.Н. Трубецкой, Д.М. Бронников, С.В. Кузнецов, В.А. Трофимов // Физика Земли. - 1994. - № 7-8. - С. 78-83.

— Коротко об авторах --------------------------------------------------

Рассказов И. Ю. - доктор технических наук, заместитель директора по научным вопросам,

Искра А. Ю. - старший научный сотрудник,

Кянно К. А. - ведущий инженер,

Институт горного дела ДВО РАН, г. Хабаровск.

1984.