The Algorithm of Calculation of Density of a Stream of Radon(Rn) from the Surface of the Ground .pdf В последние десятилетия уделяется повышенное внимание к поиску предвестников землетрясений в сейсмоэманационных эффектах. В работах [1 - 4] показана перспективность таких исследований и предложен ряд новых подходов, в которых рассматриваются методы исследования массопереноса почвенного радона, с целью поиска предвестников сильных землетрясений, причем расчеты проводятся с помощью диффузионно-конвективной модели [5]. На Петропавловск-Камчатском геодинамическом полигоне регистрация Rn в каждом пункте осуществляется в двух точках на различных глубинах в рыхлых отложениях, которые рассматриваются как относительно однородная и пористая среда [6]. Как показано в работе [4], для такой среды в динамике Rn будет слабо выражена реакция на изменение напряженно-деформированного состояния геосреды. Это связано с тем, что изменение скорости конвекции не может привести к росту объемной активности (ОА) Rn выше ее фонового значения [7].Авторы работы [7] предлагают для повышения эффективности выделения предвестниковых аномалий использовать параметр, связанный с ростом скорости конвекции - плотность потока радона (ППР) с поверхности земли. Увеличение скорости конвекции приводит к увеличению конвективного потока и как следствие к увеличению ППР на земной поверхности. Расчеты показали, что увеличение скорости конвекции в два раза приводят к росту ППР в 4 раза, в тоже время рост концентрации Rn увеличивается всего в два раза, причем в неоднородной среде ППР ещё более чувствительна к вариациям скорости конвекции. В работе [8] была реализована эта идея, где расчет ППР с поверхности земли проводился с помощью эмпирической формулы [9].Естественно, что при поисках предвестниковых аномалий необходимо учитывать влияние вариаций метеорологических величин, в частности атмосферного давления. Исходные данные, полученные на сети станций, подвергались бароком-пенсации по методике А.А. Любушина.В настоящей работе рассмотрен алгоритм расчета ППР с поверхности земли, который представляет собой решение некорректной обратной задачи массопере-носа радона в приземный слой атмосферы. Данный алгоритм реализован в программе для обработки геофизических данных РЭКСЭМ [11], где также реализован метод расчета ППР рассмотренный выше. Программа проводит обработку экспериментальных данных и визуализацию расчетов ППР, выделяет аномальные значения в зависимости от превышения над уровнем фона. Данная методика прошла тестирование на данных сети станций мониторинга подпочвенного радона на Камчатке [10].1. Вычисление плотности потока на земной поверхности как решение обратной задачи массопереноса радона1.1. Постановка задачиЗадача нахождения ППР на земной поверхности по измеренным значениям концентраций радона Nex на различных глубинах с краевыми условиями имеет следующий вид [5] := q (t),Здесь N = Q'n?s Бк/м3 - фоновая концентрация радона на заданной глубине для данного района мониторинга; Q, Бк/(м3-с) - скорость образования (эманирования) радона в грунте; n- пористость грунта; ^, с-1 - постоянная полураспада радона; D, м2/с - коэффициент диффузии радона в грунте; De=D/n, м2/с - эффективный коэффициент диффузии; v, м/с - скорость конвективного потока радона; qc, Бк/(м2-с) - постоянная единичная ППР на земной поверхности; N(z,t), Бк/м3 - объемная концентрация радона. В данной модели предполагается, что процесс массо-переноса радона на земную поверхность в некоторый момент времени является установившимся. Поэтому параметры модели являются константами и также предполагаются известными.Для задачи (1) оценка ППР на поверхности по экспериментальным данным является обратной некорректной задачей геофизики, так как не задано ни одного граничного условия при z = 0.1.2. Алгоритм решения обратной задачи массопереноса радонаРассмотрим алгоритм решения задачи по определению ППР, основываясь на приближенном соответствии между измеренной концентрацией Rn для нескольких разноглубинных датчиков и рассчитанной концентрации по диффузионно-конвективной модели [5,12].Алгоритм решения можно разделить на следующие этапы:1))аппроксимация искомой ППР постоянной единичной ППР и решение соответствующей прямой задачи по определению концентрации Rn в точках наблюдения;2)представление решения прямой задачи массопереноса Rn как суперпозиции его элементарных блоков с помощью интегральной формулы Дюамеля;3)минимизация невязки между концентрацией Rn, полученной экспериментальным путем и рассчитанной по модели методом наименьших квадратов.1)1.3. Аппроксимация ППР ступенчатой единичной функциейДля решения прямой задачи массопереноса радона аппроксимируем ППР на земной поверхности постоянной единичной функцией:(0, t < 0,q (t} = qc (1, t > 0.На поверхности ППР аппроксимируется ступенчатой единичной плотностью. На интервалах времени т1/2, т3/2,...., тМ-1/2 используются значения плотности потока qbq2,..., qMв моменты времени tM: от 0 до т1, от т1 до т2, от тМ-1 до тМ.1.4. Решение прямой задачи массопереноса радонаПрямая задача массопереноса радона на земную поверхность представляет собой дифференциальное уравнение в частных производных с начальными и граничными условиями [12]:а^(М = De- v-XN (z, t) + Q,dtdz dzIN(z,0) = No, t = 0, (2) -De ndN^ + vnN (0, t ) = qc.zРешение задачи (2) можно получить интегральным преобразованием: N(z,t) = (N„ - qc) m2DeH/(Xn), H = Hi + H2- H3, Hi = (m + n)exp((m - n)z) erfc(z/(4Det)1/2 - n(Det)1/2),H2 = (m - n) exp((m + n)z) erfc(z/(4Det)1/2 + n(Det)1/2), (3) H3 = 2exp(-zm - Xt) erfc(z/(4Det)1/2 - m(Det)1/2), m = v/(2De), n = (v2 + 4DeX)1/2/(2De). Это решение описывает процесс изменения концентрации Rn в рыхлых отложениях под воздействием постоянной единичной плотности потока qc в моменты времени t. Дифференцируя решение (3) по qc, получаем выражение для коэффициента чувствительности cp(z,t) = SN(z,t)/öqc = - m2DeH/(Xn), причем если коэффициенты чувствительности малы или коррелированны, то задача оценивания становится трудной и чувствительной к погрешностям измерений.1.5. Интегральная формула ДюамеляСогласно аппроксимации ППР в виде единичной функции и теореме Дюамеля о представлении решения в виде суперпозиции элементарных блоков его решений по соответствующим интервалам времени имеемN(z,tM) = N,x,+q1[cp(z,tM-To) - q>(z,tM- T1)]+q2[q>(z,tM- T1) - cp(z,tM - T2XI+....+ + q^[?(z,tM - Тм-1) - ?(z,tM - Tm)], q>(z, tM - Tm) = ?(z,0)], M = 0, 1, 2,..., /.(4) Здесь cp(z, tM) - элементарное решение в tM момент времени. Выражение (4) можно упростить, записав его в следующей дискретной форме:MN(z,Im) = Nco + Z 9„Аф(z,?m-„), Аф(z,*м_„) = Ф(z,?m-„+i)"Ф(z>*м-икоторое называется интегральной формулой Дюамеля. Опуская индексы его можно переписать в следующей простой форме:МNM = Nco + Z qnАфМ-n (6) n=1При решении обратных задач массопереноса радона на земную поверхность, выражение (6) является важным, так как устанавливает удобную связь между концентрацией Rn и ППР на поверхности. Так как ППР изменяется по времени, то выражение (6) дает численный результат.1.6. Метод наименьших квадратов для обратной задачи массопереноса радонаЕсли у нас имеется J датчиков измерения концентрации Rn, расположенных на различных глубинах, то тогда для каждой точки наблюдения получаем:N1m = N1m|?m =о + p11qm >N2m = N2M\qM =0 +p21qM '(7)njm = njm|qm =0 +pJ 1qM ■Между измеренной и рассчитанной концентрациями Rn довольно сложно установить точное согласование, так как единственную составляющую ППР трудно выбрать чтобы выполнялось равенство NjM = Nex. Поэтому согласование должно выполняться в некотором усредненном приближении. Оценим qM с учетом (7) методом наименьших квадратов:S = X (( - NjM [m =0 - 9ylqM) . (8)Минимизируем функционал невязки (8) по составляющей плотности потока qM, что сводится к нахождению нулей производных по параметру qM:J ГM-1 ЛZ I Nex - Z qßАф jm-i - N„ 1ф jЧм = ^*=S- .(9)J 2k=1Замечание. Предлагаемый алгоритм решения обратной задачи по нахождению ППР имеет один недостаток - чувствительность к погрешностям измерения при малых шагах по времени (секунды) [12]. Однако регистрация экспериментальных данных сети мониторинга радона осуществляется с шагом дискретизации 6,0 ч-1, поэтому при реализации данного алгоритма чувствительность к погрешностям измерения является не высокой.2. Экспериментальные данные и результаты расчетовВ пункте «Институт» (ИНС), располагающемся рядом со зданием Института вулканологии и сейсмологии ДВО РАН, осуществляется регистрация объемной активности (ОА) Rn в подпочвенном воздухе на глубине один метр от пола подземного бункера глубиной 2,5 м (зона аэрации) и с поверхности пола. Над бункером на поверхности земли установлен металлический контейнер, причем подземный бункер через контейнер и систему труб вентилируется за счет естественной конвекции воздуха. В наземном контейнере и подземном бункере установлены датчики температуры и атмосферного давления. Регистрация всех параметров осуществляется с помощью двух измерительных приборов ALMEMO 2590-9 с частотой дискретизации 6,0 ч-1 [8].Как показало сравнение метеорологических величин в подземном бункере и наземном контейнере, в подземном бункере в летнее время, за счет больших вариации температуры в наземном контейнере, четко прослеживается суточные колебания атмосферного давления, которые хорошо видны на рис. 1, а, где показана динамика ОА Rn за период 20 июля - 29 августа 2006 г. в зоне аэрации и с поверхности. На рис. 1, а видно, что с 13 по 22 августа наблюдается увеличение амплитуды объемной активности Rn на обоих каналах, причем в зоне аэрации превышение составляет 115%, а на поверхности 130% от фона. На рис. 1, б приведена расчетная кривая ППР, вычисленная по описанной выше методике с помощью программы «РЭКСЭМ» со следующими значениями параметров: h1 = 0,2 м, hi = 1,2 м, D = 5,510-7 м2/с, v = 5-10-5 м/с, I = 2,110-6 с-1. В расчетных кривых ППР превышение амплитуды над фоном для аномального периода составляет более 160%, что подтверждает более высокую чувствительность ППР по сравнению с динамическими параметрами ОА Rn.3 яО0,6В дальнейшем планируется использовать предложенную методику расчета ППР с поверхности земли в реальном времени с целью мониторинга изменений напряженно-деформированного состояния геосреды в районе Петропавловск-Камчатского геодинамического полигона для прогноза сильных землетрясений южной Камчатки.ЛИТЕРАТУРА

Скачать электронную версию публикации