CC BY
93
УДК 550.42:546
ДИФФУЗИОННО-АДВЕКТИВНЫЙ ПЕРЕНОС РАДОНА В МНОГОСЛОЙНЫХ ГЕОЛОГИЧЕСКИХ СРЕДАХ
В.С. Яковлева
Томский политехнический университет E-mail: [email protected]
Представлены решения стационарного диффузионно-адвективного уравнения переноса радона в многослойных геологических средах численными методами. Для построения однородных консервативных разностных схем использовали интегро-интерпо-ляционный метод (метод баланса). Полученные явную и неявную конечноразностные схемы решали методом правой и левой прогонки. Произведен выбор оптимального шага сетки из сопоставления данных численных расчетов и аналитического решения, полученного в пакете Mathematica 5 для стационарного случая. С помощью полученной численной модели диффузионно-адвективного переноса радона можно рассчитывать распределения поровой активности радона и плотности потока радона по глубине и во внешнюю атмосферу. Приведены результаты численных расчетов функции поровой активности радона по глубине для случая трехслойной среды, представленной разными типами грунтов и пород с сильно различающимися физико-геологическими характеристиками.
Ключевые слова:
Радон, перенос, диффузия, адвекция, численное моделирование, слой, геологическая среда. Key words:
Radon, transfer, diffusion, advection, numerical modeling, layer, geological medium.
Актуальность радоновой тематики в различных областях науки и практики до сих пор продолжает расти. Все еще остаются не до конца изученными вопросы, связанные с выявлением и описанием процессов и механизмов переноса радона в различных средах, факторов, обуславливающих временную и пространственную динамику радонового поля. Это способствует активному развитию методов математического моделирования процессов переноса радона и его дочерних продуктов распада в различных средах. На основе моделей производят прогнозные оценки, находящие свое приложение в различных областях, например, в строительстве для введения поправок к существующим нормам и правилам, в геологоразведке при поиске урансо-держащих руд, в радиоэкологии при оценках радо-ноопасности территорий и зданий, в геофизике при изучении литосферно-атмосферных связей, поиске предвестников землетрясений.
Для моделирования переноса радона широко используют диффузионно-адвективную (диффузионно-конвективную) модель, которая имеет простое аналитическое решение в некоторых частных случаях, при многих ограничениях и допущениях (стационарный случай, однородная среда и пр.). Термин «адвекция», активно используемый за рубежом [1-10], был введен для обозначения перемещения радона под воздействием внешних сил, в основном, градиента давления, а также включая процессы конвекции, фильтрации и другие, имеющие место в рассматриваемой ситуации, кроме молекулярной диффузии. В России авторами работ [11-15] был введен условный термин «конвекция», для обозначения перемещения радона под действием тех же внешних воздействий, который на некоторое время прижился в России и за рубежом среди ученых, занимающихся радоновой тематикой, и продолжает использоваться [16-21] в на-
стоящее время. Однако, чтобы избежать путаницы в подразумеваемых под условным термином «конвекция» процессах переноса радона, в данной работе используется термин «адвекция».
Моделирование переноса радона в геологических средах, максимально приближенных при их описании к реальным, является сложной задачей, поскольку геологическая среда является неоднородной, слоистой с заметно различающимися физико-геологическими характеристиками. При рассмотрении неоднородных сред, особенно, когда коэффициенты уравнения переноса являются функциями, а не постоянными, аналитическое решение невозможно. В данном случае для решения уравнений переноса радона используют численные методы, чему и посвящена данная работа. Далее последовательно рассмотрены численные решения стационарного уравнения переноса радона в однородной и слоистой геологических средах.
Модель переноса радона в стационарном случае в однородной пористой среде с постоянными коэффициентами
В стационарном случае одномерное диффузионно-адвективное уравнение переноса радона в пористой среде имеет вид (ось г направлена вниз от поверхности земли, г>0):
Ое+ и-АА(г) + АА„ = 0, 0<г<ю, (1)
5z2
5z
с граничными условиями
A(0) = 0, A(x) = A
(2)
где А(г) - поровая активность радона, приходящаяся на единицу объема порового пространства, Бк-м-3; и - скорость адвективного переноса радона, м-с-1; Бе - эффективный (объемный) коэффи-
циент диффузии радона, м2-с-1; X - постоянная распада радона, с-1; A - поровая активность радона, находящегося в радиоактивном равновесии с
226Ra, равная Ая = КетАаPs (1 ~пП, где Km - коэф-
П
фициент эманирования радона, отн. ед.; АRa -удельная активность 226Ra, Бк-кг-1; ps - плотность твердых частиц грунта, кг-м-3; п - пористость грунта, отн. ед.
Для решения такой задачи пригоден метод конечных элементов. Рассмотрим равномерную сетку с шагом h на полубесконечном интервале 0<z<<». Раскладывая производные по формуле Тейлора в окрестности произвольного внутреннего узла сетки z=zi с шагом h и используя обозначения, опуская аргументы A;=A(z), получим выражение для левой разностной производной (опуская индекс i)
Ax s = dZ-(z) -2^(z) + °(h2)' (3)
h dz 2 dz и для правой разностной производной
A.,, A, dA , „ h d A , Л„ 2, ✓ л\ Ax = '+1, ' =—(zi) + -—(zi) + °(h2). (4) h dz 2 dz
Решая систему из уравнений (3) и (4), находим выражения для первой и второй производных функции A(z), которые подставляем в исходное уравнение (1)
D .A'+1 + A-1 -24 +u4+1 -A-1 _xa =-XAa. (5)
h2
2h
Поскольку в поставленной задаче определены граничные условия (2), для его решения можно использовать метод левой или правой прогонки [22]. Метод прогонки заключается в решении трехто-
чечного разностного уравнения (6), к которому приводим уравнение (5), обозначив Ai+l=ym, Аи=УИ, Ai=y,
АУ- - СУ + Вум I = 1,2,..., N-1, (6)
где коэффициенты уравнения определены как
= D '= h 2h :
4 = D ' h2 2h
С' = + X, F =X4
(7)
Алгоритм решения уравнения (6) с коэффициентами (7) по методу правой прогонки выглядит следующим образом
В.
- С,-а4 а = о, д = о, у = а+1 у ,+1+в
yN = A •
в+1 =
A в+ F С-аА
' = 1,2,..., N-1,
' = N -1,N - 2....0,
(8)
Аналогично строится алгоритм левой прогонки [22].
Расчеты функции распределения поровой активности радона по глубине были произведены по алгоритмам левой и правой прогонки для территории, геологическая структура которой относительно однородна, т. е. поверхностный слой мощностью до 6 м представлен суглинками со следующими параметрами: эффективный (объемный) коэффициент диффузии радона De=0,03 см2-с-1; коэффициент эманирования радона ^=0,2; удельная активность 22!^а в грунте А&=30 Бк-кг-1; плотность твердых частиц грунта р^2,7 г-см-3; пористость грунта п=0,45; скорость адвекции для расчетов
глубина, см
Рис. 1. Распределение поровой активности радона по глубине грунта 68
принята и=10-4 см-с-1. Радоновые аномалии, т. е. залежи пород с высоким содержанием урана, крупные разломы и трещины в земной коре отсутствуют. Результаты численного моделирования с разным шагом пространственной сетки были сопоставлены с результатом аналитического решения [23], в итоге для дальнейших расчетов был выбран метод правой прогонки и шаг й=1 см, при котором наблюдали максимальное совпадение результатов. Метод левой прогонки давал существенное расхождение численных и аналитических расчетов на правой границе сетки (рис. 1).
Решение стационарного уравнения переноса радона в слоистой геологической среде с непостоянными коэффициентами
Процесс переноса радона посредством диффузии и адвекции характеризуются интегральным законом сохранения массы. В случае неоднородной, слоистой среды, когда коэффициенты уравнения, характеризующие каждый слой грунта (породы), на границе двух сред резко изменяются, необходимо рассматривать уравнение переноса как уравнение с разрывными коэффициентами, которые зависят от пространственной координаты -глубины грунта. Согласно этому замечанию уравнение переноса радона (1) будет иметь следующий вид
А Гц („ )+„, „ -
-ХЛ( I) + ХЛоо (I) = 0.
(9)
Для решения задачи (2,9) необходимо использовать интегро-интерполяционный метод. Согласно этому методу, переход от дифференциального уравнения (9) к разностному уравнению осуществляется с помощью некоторого интегрального соотношения (уравнения баланса), выражающего закон сохранения для малого объема [24, 25]. Входящие в уравнение баланса интегралы и производные следует заменить приближенными разностными выражениями.
Запишем для уравнения (9) уравнение баланса на отрезк 11-цг<1<11+цг
I 1Ц(2]—Г + 1 ^^ -
дх
-X I Л(2)йх = -Х I Лю (I)йх.
(10)
#/+1/2 4-1-1/2
= И Ц+1/2У,+1 - (Ц+1/2 + Ц-1/2 )У/ + Ц-1/2 У/-1 ). (11)
И
Произведем далее аппроксимацию второго интеграла, отражающего адвективный поток радона, с использованием квадратурной формулы трапеций
2 ЖI)
I "(I)-
дх
-¿Ь !
«-^(У, - У/-1) У/+1 - у ). (12)
Произведем аппроксимацию других членов уравнения (распад и образование радона) функциями, значения которых будем искать в целых узлах сетки
1+1/2 2,+1/2 X I Л(I)(к « ХПу,, X I Ла (I)сЪ « ХИЛ(13)
Выражения (11), (12) и (13) подставим в уравнение (10), и сразу приведем к трехточечному разностному уравнению вида (6) для определения его коэффициентов
—Ц
-¿-Л- 1 /о
„2 /-1/2
1
(Ц-1/2 + Цг+1/2) +
2И и
У,-1
/+1/2 2И
/-1/2 2И
X
х У,
2- Ц .
И2 Ц+1/2
и
/+1/2 2И
У+1 = ~ХЛт.. (14)
(14)
Введем следующие обозначения для выражения
Л =1Ц
^ 1/2
2И
С, = И2( Ц+1/2 + Ц /—1 / 2 ) +
+ 2-("+1/2 +и-1/2) + Х
2И
В.=4 Ц
и
/+1/2
И2 г+1/2 2И
Разностное уравнение (14) с коэффициентами (15), полученное при аппроксимации краевой задачи переноса радона в многослойной геологической среде, решаем методом правой прогонки по алгоритму (8) с теми же граничными условиями.
Рассмотрим практическую задачу, для решения которой может быть использована описанная численная модель. Большой интерес представляет оценка влияния на величину объемной активности приповерхностного радона (на глубине измерения, обычно до 1 м) ниже залегающей породы с повышенным содержанием урана. Произведем расчеты функции распределения поровой активности радо-
Fi = ХЛ„..
(15)
Балансное соотношение (10) отражает закон сохранения для отрезка £;-1/2<г<г;+1/2. Для получения разностного уравнения из балансного соотношения (10) необходимо использовать те или иные восполнения сеточных функций. Функцию решения будем искать в целых узлах (у(т), 1=1), а диффузионный и адвективный потоки - в полуцелых. Представим первый интеграл разностью диффузионных потоков (д(г), 1=1±\/2) в полуцелых узловых точках и запишем их аппроксимацию в соответствии с [25]
г
1-1/ 2
на по глубине для территории, геологическая структура которой является слоистой.
Рассмотрим среду, состоящую из 3-х слоев мощностью 4, 2 и 4 м, параметры которых для упрощения анализа взяты одинаковыми, как для суглинков из описанного выше примера, кроме удельной активности радия. В первом и третьем слое Аа=30 Бк-кг-1. Во втором слое активность радия, продукта распада урана, была выбрана для расчетов из характерного для урансодержащих пород диапазона АКа=1000 Бк-кг-1. Численные расчеты были произведены при разных значениях скорости адвекции, взятые из диапазона значений, полученного в работах [17, 23]. При положительных значениях скорости и адвективный поток направлен к поверхности земли и складывается с диффузионным, увеличивая суммарный поток радона в атмосферу. При отрицательных значениях и адвективный поток радона направлен вглубь земной поверхности снижая суммарный поток радона в атмосферу и, соответственно, значение поровой активности радона у поверхности земли.
На рис. 2 показаны функции зависимости поровой активности радона от глубины при положительной и отрицательной скорости адвекции |и|=10-4 см-с-1. Там же для сравнения показана функция распределения поровой активности по глубине в случае однородной среды, т. е. все три слоя представлены суглинком с одинаковым содержание радия АКа=30 Бк-кг-1, скорость адвективного потока взята раной 10-4 см-с-1.
Для проверки правильности модели, было произведено сравнение данных численного моделиро-
вания и аналитических расчетов, полученных в пакете МаШешаНеа при решении системы 3-х уравнений типа (1), записанных для каждого слоя с учетом условий непрерывности потоков на границе слоев и равенства объемных активностей. Результаты аналитических и численных расчетов имеют хорошее согласие при шаге пространственной сетки 1 см.
Анализ функций распределения поровой активности (рис. 2) показал, что влияние второго высокоактивного слоя, залегающего даже на небольшой глубине от 4 до 6 м, может, при определенных условиях, не проявляться на результатах измерения по-ровой активности радона, производимых, например, в радиоэкологических или геологоразведочных целях, при оценках радоноопасности территорий. При отрицательной скорости конвекции по-ровая активность радона до глубины 2 м не превышает значений, которые наблюдаются в относительно однородной среде.
Таким образом, если измерения производить в те периоды времени, когда адвективный поток направлен вглубь земли, можно получить ошибочные результаты по оценке радоноопасности обследуемой территории. И, наоборот, при больших адвективных потоках в атмосферу, поровая активность радона на глубине измерения (обычно 1 м) может увеличиваться до 5 раз.
Анализ результатов численного моделирования, а также экспериментальные данные [17, 23, 26, 27], показывающие, что адвективные потоки радона и, соответственно, скорость адвекции, могут изменяться в широком диапазоне значений даже в течение суток, доказывают необходимость
0,35
0,3
п
о 0,25 ■2 ш
.о
о 0,2
0
X
Ой
1 0,15
к ш
Ой
о
§■ 0,1 с
0,05
1 слой 2 слой V 3 слой
\ ---однородная среда
/ д * \ / \ 1 \ / \ % \ -скорость адвекции положительная * \ — - скорость адвекции отрицательная \
/ \ 1 1 % \ %
/ 1 \ \
/ 1 . \ 1 Ч
/ щ / / \ \ 1 Ч \ \
/ / _ / # / / / \ 4 ч \ % \ ^—. _
( __. - " '
0
0 83 166 249 332 415 498 581 664 747 830 913 996 1079 1162 1245 1328 1411 1494
глубина, см
Рис. 2. Распределение поровой активности радона по глубине грунта
проведения интегральных измерений активности радона в почвенном воздухе для получения достоверной оценки радонового потенциала обследуемой территории.
На основе описанной численной модели переноса радона в неоднородных слоистых геологических средах был разработан алгоритм и написана программа «SimRaTran» [28], позволяющая моделировать перенос радона в средах, представленных несколькими (до 20-ти) эманирующими слоями.
Выводы
Получены разностные уравнения переноса радона с учетом процессов диффузии и адвекции в относительно однородной и неоднородной (слоистой) средах, определены коэффициенты уравнения. В случае однородной среды достаточно использовать метод конечных элементов в сочетании с методом правой прогонки, в случае неоднород-
СПИСОК ЛИТЕРАТУРЫ
1. Etiope G., Martinelli G. Migration of carrier and trace gases in the geosphere: an overview // Physics of The Earth and Planetary Interiors. - 2002. - V. 129. - P. 185-204.
2. Lerman A. Geochemical Processes-Water and Sediment Environments. - New York: Wiley-Interscience, 1979. - 480 р.
3. Gascoyne M., Wuschke D.M., Durrance E.M. Fracture detection and groundwater flow characterization using He and Rn in soil gases, Manitoba, Canada // Applied Geochemistry. - 1993. - V. 8. -P. 223-233.
4. Klusman R.W Soil-Gas and Related Methods for Natural Resource Exploration. - N.Y.: Wiley, 1993. - 483 р.
5. Hinkle M. Environmental conditions affecting concentrations of He, CO2, O2 and N2 in soil-gases // Appl. Geochem. - 1994. - V. 9.
- P. 53-63.
6. Clements WE., Wilkening M.H. Atmospheric pressure effects on ra-don-222 transport across the Earth-air interface // J. Geophys. Res.
- 1974. - V. 79. - P. 5025-5029.
7. Wilkening M.H. Radon transport mechanism below the Earth's surface // The Natural Radiation Environment III: Technical Information Center, US Department of Environment, Springfield, 1980. -P. 90-104.
8. Volckaert G., Put M., Ortiz L., De Canniere P., Horseman S.T., Harrington J.F., Fioravante V., Impey M., Worgan K. MEGAS -Modelling and experiments on gas migration in repository host rocks // Proc. of the Meeting of Koln on PEGASUS Project. 3-4 June 1993. CEC EUR Series, 1993.
9. Mosley R. B., Menetrez M.Y., Snoddy R., Brubaker S.A. The influence of diffusion and advective flow on the distribution of radon activity within USEPA's soil chamber // Environmental International.
- 1996. - V. 22. (Suppl. 1) - P. 521-533.
10. Van der Spoel W. H., Van der Graaf E.R., De Meijer R.J. Diffusive transport of radon in a column of moisturized sand // Health Phys.
- 1999. - V. 77 (2). - P. 163-177.
11. Граммаков А.Г., Никонов А.И., Тарфеев Г.П. Радиометрические методы поисков и разведки урановых руд. - М.: Госгеол-техиздат, 1957. - 610 с.
12. Булашевич Ю.П, Хайретдинов Р.К. К теории диффузии эманации в пористых средах // Известия АН СССР. Серия геофизическая. - 1959. - № 12. - С. 1787-1792.
ных слоистых сред, параметры которых скачкообразно меняются на границе двух сред (слоев), хорошие результаты дает интегро-интерполяционный метод в сочетании с методом правой прогонки. Из сравнения результатов численных и аналитических расчетов был определен оптимальный для задач переноса радона шаг разностной сетки равный 1 см.
Разработанная численная модель имеет следующие преимущества по сравнению с аналитическим решением:
• коэффициенты уравнения переноса радона могут задаваться функциями, отдельно для каждого слоя;
• можно проводить расчеты для большого количества слоев грунта с различными физико-геологическими характеристиками.
Работа выполнена при поддержке Аналитической ведомственной целевой программы «Развитие научного потенциала высшей школы» № 2.1.1/544.
13. Булашевич Ю.П, Карташов Н.П. Определение коэффициента диффузии радона в горных породах методом мгновенного источника // Известия АН СССР. Физика Земли. - 1967. - № 10.
- С. 71-76.
14. Новиков Г.Ф., Капков Ю.Н. Радиоактивные методы разведки.
- Л.: Недра, 1965. - 759 с.
15. Новиков Г.Ф. Радиометрическая разведка. - Л.: Недра, 1989. -407 с.
16. Хайкович И.М. Математическое моделирование процессов миграции радона // Аппаратура и новости радиационных измерений. - 1996/97. - № 3 (9). - С. 99-107.
17. Iakovleva V.S., Ryzhakova N.K. A method for estimating the con-vective radon transport velocity in soils // Radiation Measurements.
- 2003. - V. 36. - № 1. - P. 389-391.
18. Паровик Р.И., Ильин И.А., Фирстов П.П. Обобщенная одномерная модель массопереноса радона и его эксхаляция в приземный слой атмосферы // Математическое моделирование. -2007. - № 11. - Т. 19. - С. 43-50.
19. Козлова И.А. Определение диффузионных параметров среды методом мгновенного источника в конвективном режиме измерений // VIII Уральская молодежная научная школа по геофизике: Сб. матер. - Пермь: ГИ УрО РАН, 2007. - С. 126-129.
20. Jiranek M., Svoboda Z. The verification of radon protective measures by means of a computer model // Proceedings of Building Simulation. - 1997. - V. 2. - P. 165-171.
21. Leung J.K.C., Tso M.Y.W., Hung L.C. Radon migration by depres-surization of concrete walls and slabs // Health Phys. - 1999. -V. 77. - № 4. - P. 420-426.
22. Самарский А.А. Введение в теорию разностных схем. - М.: Наука, 1971. - 552 с.
23. Яковлева В.С. Процессы переноса радона в неравновесных средах: Автореф. дис. ... канд. физ.-мат. наук. - Томск: ТПУ, 2002. - 23 с.
24. Тихонов А.Н. О разностных схемах для уравнений с разрывными коэффициентами // Доклады АН СССР. - 1956. - Т. 108. -№ 3. - С. 393-396.
25. Самарский А.А., Вабишевич П.Н. Численные методы решения задач конвекции-диффузии. - M.: УРСС, 2004. - 248 с.
26. Firstov P.P., Yakovleva V.S., Shirokov V.A., Rulenko O.P., Filip-pov Yu.A., Malysheva O.P. The nexus of soil radon and hydrogen dy-
namics and seismicity of the northern flank of the Kuril-Kamchatka subduction zone // Annals of Geophysics. - 2007. - V. 50. - № 4. - P. 547-556.
27. Фирстов П.П., Чернева Н.В., Пономарев Е.А., Бузевич А.В. Подпочвенный радон и напряженность электрического поля атмосферы в районе Петропавловск-Камчатского геодинами-чечского полигона // Вестник КРАУНЦ. Науки о Земле. -2006. - № 1. - Вып. 7. - С. 102-109.
28. Программа для численного моделирования переноса радона в пористых многослойных средах «SimRaTran» / В.С. Яковлева. Свидетельство от отраслевой регистрации разработки № 11729 от 10.11.2008 г. Федеральное агентство по образованию. ФГНУ «Государственный координационный центр информационных технологий».
Поступила 06.04.2009г.
УДК 543.427.2
ОЦЕНКА ПОРОГА ОБНАРУЖЕНИЯ РЕНТГЕНОТРАНСМИССИОННОГО МЕТОДА ОПРЕДЕЛЕНИЯ СЕРЫ В НЕФТИ И НЕФТЕПРОДУКТАХ
Т.Н. Стрежнева, А.А. Лобова, Н.А. Антропов, Ю.Ю. Крючков
Томский политехнический университет E-mail: [email protected]
Описан рентгенотрансмиссионный метод определения содержания серы в нефти и нефтепродуктах в технологических трубопроводах. Проведена оценка порога обнаружения серы. Показано, что метод позволяет уменьшить погрешность определения серы в таких важных объектах, как нефть и нефтепродукты.
Ключевые слова:
Рентгенотрансмиссионный метод, нефть, концентрация серы, количественный анализ, порог обнаружения. Key words:
X-raytransmission method, oil, sulfur content, quantitative analysis, detecting threshold.
Рентгеновская трансмиссия - прохождение рентгеновского излучения через вещество - находит применение в дефектоскопии, рентгеновской томографии, толщинометрии, контроле параметров. Ослабление интенсивности излучения зависит от элементного состава, плотности вещества, толщины анализируемой пробы. Как аналитический, рентгенотрансмиссоный метод не нашел широкого применения, поскольку для веществ, имеющих сложный химический состав, практически невозможно выделить вклад отдельного элемента в общее ослабление первичного излучения. В то же время в частных случаях этот метод обладает рядом преимуществ, ставящих его вне конкуренции по сравнению с другими аналитическими методами. Одним из таких случаев является определение содержания серы в нефти и нефтепродуктах в технологических трубопроводах под давлением до 50 атм в автоматическом режиме. В настоящее время рентгенотрансмиссионный метод является единственным методом, который может быть использован для решения данной задачи. Наиболее близким конкурентом этому методу является рентгенофлуо-ресцентный анализ, однако для измерения в технологических трубопроводах он не годится, т. к. энергия характеристического рентгеновского излучения серы очень мала, ~2,3 кэВ, и она полностью поглощается любым окном, разделяющим пробу и детектор.
Важной характеристикой аналитического метода является порог обнаружения измеряемой величины, в нашем случае серы, который определяет нижний край диапазона измеряемых концентраций серы. Хотя теория взаимодействия излучения с веществом достаточно хорошо разработана, в литературе не встречаются методики оценки параметров рентгеновской трансмиссии как аналитического метода для многокомпонентных сред переменного состава, к которым можно отнести и нефть. В силу названной причины использование этого метода для анализа в различных технологических процессах является недооцененным. В данной работе проведена оценка порога обнаружения для вышеназванных сред, учтены факторы, влияющие на порог обнаружения.
В дальнейших расчетах за основу взяты параметры поточного анализатора серы производства американской фирмы «Asoma», модель 682Т-НР, у которого заявленный порог обнаружения серы -0,1 % (в массовых долях %), время анализа 200 с, толщина анализируемого слоя х=0,022 м. Этот прибор является единственным представителем данного класса средств измерений на мировом рынке.
В общем случае, интенсивность параллельного моноэнергетического пучка после прохождения однородного слоя образца определяется формулой [1]:
