Геоэкология. Инженерная геология, гидрогеология, геокриология, 2020, № 3, стр. 3-22

ПРОГНОЗ ВОЗДЕЙСТВИЯ АЭС НА РАДИОАКТИВНОСТЬ ПОВЕРХНОСТНЫХ И ПОДЗЕМНЫХ ВОД

В. Г. Румынин 12*, Л. Н. Синдаловский 12, А. А. Шварц 12, А. М. Никуленков 12, В. А. Ерзова 3**, Д. В. Бутырин 12

1 Санкт-Петербургское отделение Института геоэкологии им. Е.М. Сергеева РАН
199004 Санкт-Петербург, Средний пр., 41, Россия

2 Санкт-Петербургский государственный университет, Институт наук о Земле
199004 Санкт-Петербург, Средний пр., 41, Россия

3 Санкт-Петербургский горный университет
199121 Санкт-Петербург, 21 Линия, 2, Россия

* E-mail: office@hgepro.ru
** E-mail: valentina.valya-06@yandex.ru

Поступила в редакцию 15.11.2019
После доработки 01.12.2019
Принята к публикации 25.12.2019

Полный текст (PDF)

Аннотация

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

Ключевые слова: АЭС, продукты деления (ПД), выбросы, нормальная эксплуатация, аварийные выбросы, гидролого-гидрогеологическая модель, прогноз миграции радионуклидов

ВВЕДЕНИЕ

В районе расположения АЭС источником радионуклидов в окружающей среде являются газоаэрозольные выбросы продуктов деления (ПД) в атмосферу. При нормальной эксплуатации (НЭ) энергоблока выброс формируется за счет удаления вытяжного воздуха систем вентиляции зоны контролируемого доступа и технологических сдувок с оборудования из систем спецгазоочистки. В условиях НЭ АЭС дозовая нагрузка для критических групп населения от всех факторов радиационного воздействия газоаэрозольных выбросов энергоблоков АЭС за пределами промплощадки составляет менее 1% от основного предела годовой эффективной дозы облучения населения (1 мЗв/год), установленного в проектах российских АЭС для газоаэрозольных выбросов.

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

Тяжелые аварии, произошедшие на атомных электростанциях в США (Три-Майл-Айленд, 1979 г.), СССР (Чернобыль, 1986 г.) и Японии (Фукусима-1, 2011 г.) во многом изменили мнение о безопасности АЭС, вплоть до отказа некоторыми странами от программ развития ядерной энергетики.

В целом, представление о масштабе потенциального воздействия АЭС на окружающую среду дает табл. 1, где приведены проектные параметры различных выбросов в условиях эксплуатации реактора ВВЭР-1200 (модельные расчеты). Как следует из табл. 1, все аварии характеризуются выбросом радиоактивного йода (I-131) и двух изотопов цезия (Cs-134 и Cs-137). При ТЗПА в выбросе присутствует Sr-90.

Таблица 1.

Основные параметры высвобождения ПД (атмосферный выброс) в различных условиях эксплуатации АЭС с типом ядерного реактора ВВЭР-1200 (на примере проекта Нововоронежской АЭС-2, 2016 г.)

Условия эксплуатации НЭ ПА ЗПА ТЗПА
Класс DBC1,2 DBC4 DEC1 DEC2
Высота выброса, м 99.5 99.5 77.7 99.5 77.7 99.5 0
Продолжительность 60 лет 1, 4, 10 или 30 сут 56 ч
Обозн., ед. изм. $M{\text{'}}$, Бк/с M, Бк
Cs-134 0.4 1.7Е+06 1.7Е+08 1.5Е+07 1.5Е+09 6.4E+10 6.4E+12
Cs-137 0.6 6.3Е+05 6.2Е+07 6.8Е+06 6.8Е+08 4.0E+10 4.0E+12
I-131 140 5.3Е+07 5.3Е+08 1.5Е+09 1.5Е+10 3.5E+11 3.2E+13
Sr-90 1.1E+10 1.1E+12

Примечание: НЭ – нормальная эксплуатация; ПА – проектная авария; ЗПА – запроектная авария (DEC1 – без плавления активной зоны); ТЗПА – тяжелая запроектная авария (DEC2 – с плавлением активной зоны); ПД – продукты деления. DBC – Design Basis Conditions; DEC – Design Extension Conditions. Высота 99.5 м отвечает высотной трубе, 77.7 м – фильтру пассивной очистки, 0 м – байпас. $M{\text{'}}$ – интенсивность выброса в атмосферу (среднегодовая), M – суммарная активность выброса.

Класс аварий по интенсивности выброса образует в первом приближении ряд:

ТЗПА ≈ 10 000 · ЗПА, ЗПА ≈ 10 · ПА, ПА ≈ ≈ 10 · НЭ (в пересчете на 1 год эксплуатации АЭС).

Наиболее опасной является ТЗПА, характеризующаяся большой течью теплоносителя первого контура в сочетании с отказом системы аварийного охлаждения активной зоны при обесточивании станции на 56 ч. В данном режиме происходит расплавление и обрушение активной зоны с последующим разрушением днища реактора и выходом расплава в устройство локализации расплава.

Такая ТЗПА относится к аварии примерно уровня 5 по международной шкале INES (accidents with wider consequences) [12]. Например, уровень 5 по шкале INES присвоен аварии 1979 г. на АЭС Три-Майл-Айленд; Чернобыльская и Фукусимская аварии оцениваются по шкале INES уровнем 7 (катастрофическая авария – major/catastrophic accidents).

1. ОБОСНОВАНИЕ МАСШТАБА И КОНФИГУРАЦИИ МОДЕЛЬНОЙ ОБЛАСТИ

Существующие дисперсионные модели переноса газоаэрозольных выбросов в атмосфере позволяют довольно точно рассчитать поток осаждения радионуклидов на земную поверхность в районе размещения АЭС (условия НЭ) или конфигурацию следа и плотность осаждения аварийных радионуклидов при заданном сценарии выброса ПД (см., например, табл. 1) для различных погодных условий. Обычно постулируется размер области влияния, ограниченный “чернобыльским” радиусом 30 км.

Прогнозирование последствий НЭ предполагает использование статистически осредненных данных радиоактивных выпадений в районе АЭС. Для практических расчетов используется роза ветров, характеризующая повторяемость ветров различных направлений. Для прогнозов выбираются ландшафты и связанные с ними водосборы (см. ниже), наиболее значимые в экологическом отношении, для которых на основе моделирования рассчитывается параметр $N{\text{'}}$ – среднегодовой поток радионуклидов, осаждающихся на земную поверхность [Бк/м2 год].

Моделирование атмосферных выпадений при ПА и ЗПА выполняется для нескольких статистически обоснованных (“типовых”) условий, учитывающих время года, изменчивость направления и силы ветра, некоторые другие параметры состояния атмосферы. Для этого используются данные сайтов “погоды”. Модельный расчет дает значение характеристики ${{N}_{0}}$ – средняя плотность выпадения [Бк/м2].

Учитывая площадной характер источника загрязнения в форме радиоактивных следов на поверхности почвы, построение гидрологических и гидрогеологических моделей должно основываться на бассейновом принципе выделения границ модельных областей. Различные гидрологические бассейны в пределах 30-километровой зоны, в зависимости от конфигурации радиоактивного следа, могут иметь различную плотность загрязнения. Анализ гидрографов на замыкающих створах речных бассейнов, наряду с данными об атмосферных осадках, позволяет определить поверхностный склоновый сток, контролирующий поступление радионуклидов в поверхностные воды.

Альтернативные экспертные оценки выполняются посредством сочетания секторального и бассейнового принципов выделения границ модельных областей. Характеристика $N{\text{'}}$ рассчитывается по среднегодовой интенсивности выброса того или иного радионуклида $M{\text{'}}$ (см., например, табл. 1), которая распределяется по площади 30-километровой зоны в соответствии с розой ветров. Характеристика ${{N}_{0}}$ определяется по суммарной активности радионуклидов M, высвободившейся в процессе аварии (см. табл. 1) и равномерно распределенной в одном из секторов 30-километровой зоны с центральным углом 22.5°–45°.

2. ОСНОВНЫЕ МЕХАНИЗМЫ И МОДЕЛЬНЫЕ ПОДХОДЫ

Существуют два основных концептуальных подхода к моделированию гидрологических процессов на водосборных площадях [14]: модели с распределенными параметрами и модели с сосредоточенными параметрами.

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

− образование слоя воды на поверхности в результате избытка дождевых осадков;

− поток воды по земной поверхности в виде склонового и речного стока;

− перенос радионуклидов из почвы в поверхностный сток;

− эрозию почв, приводящую к появлению в воде дисперсных частиц – потенциальных транспортеров радионуклидов;

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

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

Рис. 1.

Схематическое представление процессов массообмена между почвой и склоновым потоком воды.

Рис. 2.

Схематическое представление двухкомпонентной модели для поверхностного-подземного стока.

Массообмен на поверхности раздела почва–вода определяется: процессами диффузии и сорбции–десорбции; турбулентностью потока, способствующей перемешиванию свободной и поровой воды, а также отделению частиц почвы с поверхности за счет механического воздействия на нее. При этом некоторые гидрогеологические модели основаны на концепции о существовании в верхней части почвенного разреза узкого активного слоя (рис. 3), в котором происходит взаимодействие (смешение) дождевой воды с почвенными водами [3, 9]. Толщина слоя смешения de принимается постоянной, и при самой общей постановке задачи допускается попадание химических веществ в поверхностный сток из подстилающих слоев почвы. Значение ${{d}_{e}}$, как правило, не превышает первые сантиметры [17].

Рис. 3.

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

Использование моделей с сосредоточенными параметрами оправдано в случаях, когда содержание растворенных компонентов в склоновом стоке или замыкающем створе в большей степени определяется кинетикой попадания радионуклидов с загрязненной поверхности почвы в поверхностный сток, нежели динамикой поверхностного потока воды (изменениями его скорости и мощности). Такие модели позволяют описать формирование слоя воды на поверхности при дождевых осадках или снеготаянии, инфильтрацию (поступление дождевых вод в почву и грунтовые воды), а также миграцию радионуклидов в поверхностных и подземных водах. Инерция водного потока, приводящая к временнóй задержке разгрузки склонового (водосборного) потока по отношению к избыточному дождевому питанию, может быть учтена через эффективные параметры в линейных или нелинейных кинетических уравнениях, или же может полностью не приниматься во внимание, если разгрузка водосборного потока принимается мгновенной, т.е. происходящей одновременно с дождевыми выпадениями. Временем миграции вещества в системе поверхностных водотоков можно пренебречь по сравнению с временем нахождения частиц в поверхностном стоке внутри водосборного суб-бассейна.

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

3. БАЛАНС РАДИОНУКЛИДОВ В ПОЧВЕННОМ СЛОЕ

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

(1)
${{r}_{e}} = r - i,$

– так называемая Хортоновская модель стокообразования [14]; где r – эффективные осадки; re – избыточные стокообразующие осадки, i – интенсивность инфильтрации; все потоковые компоненты имеют размерность скорости [м/сут].

Предположим, что расход воды, стекающей с единицы площади поверхности водосборного бассейна (после выпадения осадков), q1 [м/сут], является линейной функцией толщины слоя стока ${{S}_{1}} = {{S}_{1}}(t)$ [м]:

(2)
${{q}_{1}} = {{\mu }_{1}}{{S}_{1}},$
где μ1 – коэффициент скорости водообмена [сут–1], характеризующий общую задерживающую способность водосборного бассейна или его готовность отдавать воду. Коэффициент μ1 характеризует среднее время дренирования (MRT) водосборного бассейна:

${{\mu }_{1}} \approx 1{\text{/}}MRT.$

Уравнение водного баланса, описывающее изменение слоя стока ${{S}_{1}}$, может быть записано таким образом:

(3)
$\frac{{d{{S}_{1}}}}{{dt}} = {{r}_{e}} - {{q}_{1}}.$

Уравнение материального баланса с учетом радиоактивного распада в жидкой (слой воды на поверхности) и твердой (почва) фазах имеет вид:

(4)
$\begin{gathered} \frac{{d{{S}_{1}}{{C}_{1}}}}{{dt}} + \frac{d}{{dt}}(\theta {{d}_{e}}{{C}_{1}} + {{N}_{1}}) + \\ \, + \theta {{d}_{e}}\lambda {{C}_{1}} + \lambda {{N}_{1}} = N{\text{'}} - i{{C}_{1}} - q{{C}_{1}}, \\ \end{gathered} $
где ${{C}_{1}}$ – удельная активность радионуклидов в поверхностном стоке, принимаемая равной разбавленной концентрации в слое смешения (см. рис. 3) [Бк/м3]; $\lambda $ – константа радиоактивного распада [сут–1]; ${{N}_{1}}$ – удельная активность радионуклидов в сорбированном состоянии в слое смешения на единицу площади [Бк/м2]; $\theta $ – полная влагоемкость почвы; $N{\text{'}}$ – удельная плотность выпадения радионуклидов на поверхность [Бк/м2 сут].

Пусть обратимая сорбция подчиняется линейной изотерме ${{N}_{1}} = K_{{d1}}^{{}}{{\rho }_{b}}{{d}_{e}}{{C}_{1}}$ ($K_{{d1}}^{{}}$ – коэффициент сорбционного распределения в почве [см3/г]; ${{\rho }_{b}}$ – плотность почвы в сухом состоянии [г/см3]), тогда уравнение (4) может быть переписано в форме:

(5)
$\left( {\theta {{d}_{e}}{{R}_{1}} + \frac{{{{q}_{1}}}}{{{{\mu }_{1}}}}} \right)\frac{{d{{C}_{1}}}}{{dt}} + \lambda \left( {\theta {{d}_{e}}{{R}_{1}} + \frac{{{{q}_{1}}}}{{{{\mu }_{1}}}}} \right){{C}_{1}} = N{\text{'}} - r{{C}_{1}},$
(${{R}_{1}} = 1 + K_{{d1}}^{{}}{{\rho }_{b}}{\text{/}}\theta $ – фактор сорбционной задержки).

Частное решение системы уравнений (5), (3) для условий миграции стабильного компонента ($\lambda = 0$) при аппроксимации дождевого гидрографа r(t) ступенчатой функцией можно найти в работе [14]. Это решение учитывает растущий и падающий гидрографы поверхностного стока в периоды выпадения дождевых осадков различной интенсивности, или в периоды их отсутствия. При этом учитывается возможность как полного осушения поверхности в засушливые периоды, так и наличие на ней слоя воды, сформировавшегося в течение n-го периода, предшествующего расчетному n+1 периоду выпадения порции дождевых осадков.

На практике такая детализация процесса, обусловленная необходимостью учета конечной скорости реакции потока q1 на изменение интенсивности атмосферных осадков, может оказаться избыточной. При малых временах пребывания и удержания воды на поверхности, исследование поведения гидрологической системы можно ограничить анализом, основанным на принципе последовательной смены ее равновесных состояний. В соответствии с этим принципом считается, что отклик расхода разгрузки qn (цифровой индекс опускаем) на дождевые осадки является мгновенным, т.е. qn = ren. Такое состояние сохраняется до следующего дождевого пакета, когда наступает новое равновесие (qn + 1 = ren + 1), или в случае сухого периода qn + 1 = ren + 1 = 0. С математических позиций такое поведение системы достигается при весьма больших значениях кинетического коэффициента μ1, теоретически при ${{\mu }_{1}} \to \infty $. В этом случае уравнение (5) упрощается и принимает вид:

(6)
$\frac{{d{{C}_{1}}}}{{dt}} = \frac{{N{\text{'}}}}{{\theta {{d}_{e}}{{R}_{1}}}} - \left( {\lambda + \frac{r}{{\theta {{d}_{e}}{{R}_{1}}}}} \right){{C}_{1}},$
удобный для дальнейшего физико-математического анализа.

4. БАЗОВЫЕ РЕШЕНИЯ ДЛЯ КОНЦЕНТРАЦИОННОЙ ФУНКЦИИ В СЛОЕ СМЕШЕНИЯ

В контексте настоящей статьи, основной интерес представляют две типовые ситуации. Первая отвечает условиям НЭ АЭС, при которой основной характеристикой воздействия является удельная плотность выпадения радионуклидов на поверхность – источник активности, действующий в течение всего периода эксплуатации АЭС (десятки лет). Вторая – условиям проектных или запроектных аварий, при которых происходит мгновенный подъем активности в слое смешения de до значения, определяемого плотностью загрязнения поверхности N0 (активность радионуклидов на единицу площади сразу после аварии).

Решение уравнения (6) при $r = {\text{const}}$ и условии ${{C}_{1}}(t = 0) = {{C}_{{10}}}$ (начальная активность/концентрация) имеет вид:

(7)
$\begin{gathered} {{С}_{1}} = {{C}_{{10}}} + \left( {\frac{{N{\text{'/}}\theta {{d}_{e}}{{R}_{1}}}}{{\lambda + r{\text{/}}\theta {{d}_{e}}{{R}_{1}}}} - {{C}_{{10}}}} \right) \times \\ \, \times \left[ {1 - \exp ( - \lambda t)\exp \left( { - \frac{{rt}}{{\theta {{d}_{e}}{{R}_{1}}}}} \right)} \right]. \\ \end{gathered} $

Для условий аварийного осаждения, когда “фоновым” потоком радионуклидов можно пренебречь ($N{\text{'}} = 0$), интегрирование уравнения (6) выполняется для начальной концентрации:

(8)
${{C}_{1}}(t = 0) = {{C}_{{10}}} = \frac{{{{N}_{0}}}}{{(\theta + K_{{d1}}^{{}}{{\rho }_{b}}){{d}_{e}}}} = \frac{{{{N}_{0}}}}{{\theta {{d}_{e}}{{R}_{1}}}},$
что дает:

(9)
${{С}_{1}} = {{C}_{{10}}}\exp ( - \lambda t)\exp \left( { - \frac{{rt}}{{\theta {{d}_{e}}{{R}_{1}}}}} \right).$

Формально говоря, решения (7) и (9) могут быть обобщены на случай переменных во времени дождевых осадков. Если ${{r}_{n}} - {{r}_{{n - 1}}} = \pm \Delta {{r}_{n}}$, где $\Delta {{r}_{n}} = {\text{const}}$ в интервале ${{t}_{n}} - {{t}_{{n - 1}}}$, то концентрация в слое смешения ${{С}_{{1n}}}$ в момент ${{t}_{n}}$ может быть рассчитана, исходя из концентрации ${{С}_{{1n - 1}}}$ в предшествующий период ${{t}_{{n - 1}}}$:

– для первого сценария (НЭ):

(10)
$\begin{gathered} {{С}_{{1n}}} = {{C}_{{1n - 1}}} + \left( {\frac{{N{\text{'/}}\theta {{d}_{e}}{{R}_{1}}}}{{\lambda + {{r}_{n}}{\text{/}}\theta {{d}_{e}}{{R}_{1}}}} - {{C}_{{1n - 1}}}} \right) \times \\ \, \times \left[ {1 - \exp \left( { - \left( {\lambda + \frac{{{{r}_{n}}}}{{\theta {{d}_{e}}{{R}_{1}}}}} \right)({{t}_{n}} - {{t}_{{n - 1}}})} \right)} \right]; \\ \end{gathered} $
в случае переменных осадков в периоды (${{r}_{n}} = 0$) в расчетах должно учитываться накопление активности за счет “сухого” осаждения радионуклидов ${{С}_{n}} = {{C}_{{n - 1}}} + N{\text{'}}({{t}_{n}} - {{t}_{{n - 1}}}){\text{/}}\theta {{d}_{e}}{{R}_{1}}$;

– для второго сценария (ПА и ЗПА):

(11)
$\begin{gathered} {{С}_{{1n}}} = \exp ( - \lambda t){{C}_{{1n - 1}}}\exp \left( { - \frac{{{{r}_{n}}({{t}_{n}} - {{t}_{{n - 1}}})}}{{\theta {{d}_{e}}{{R}_{1}}}}} \right) = \\ \, = {{C}_{{1n - 1}}}\exp \left( { - \left( {\lambda + \frac{{{{r}_{n}}}}{{\theta {{d}_{e}}{{R}_{1}}}}} \right)({{t}_{n}} - {{t}_{{n - 1}}})} \right). \\ \end{gathered} $

Из структуры решений видно, что в засушливый период (${{r}_{n}}$ = 0) изменение активности обусловлено радиоактивным распадом.

В свете оценки воздействия площадного загрязнения на качество поверхностных и подземных вод необходима дифференциация потока r на компоненты поверхностного стока q1 и инфильтрации i1, о чем пойдет речь ниже в разделе 6.

5. МИГРАЦИЯ РАДИОНУКЛИДОВ С ПОДЗЕМНЫМ СТОКОМ

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

Соотношение между расходом разгрузки подземных вод q2 и расходом инфильтрационного питания f горизонта может быть представлено в виде линеаризованного уравнения фильтрации подземных вод [2].

(12)
$\frac{{d{{q}_{2}}}}{{dt}} = {{\mu }_{2}}(f - {{q}_{2}})\quad ({{q}_{2}} = {{\mu }_{2}}{{S}_{2}}),$
где ${{\mu }_{2}}$ – коэффициент скорости дренирования грунтового водоносного горизонта [сут–1]; для линеаризованного уравнения фильтрации подземных вод [2, 6] ${{\mu }_{2}} \approx (3{\kern 1pt} - {\kern 1pt} 10)T{\text{/}}\mu {\text{*}}L_{0}^{2}$ ($T$ – проводимость [м2/сут], $\mu {\text{*}}$ – водоотдача [безразмерная величина], $L_{0}^{{}}$ – средний размер бассейна в направлении основных траекторий движения поверхностных и подземных вод от водораздельной линии [м]); ${{S}_{2}} = {{S}_{2}}(t)$ – аккумулированный слой воды в горизонте.

Для описания миграции радионуклидов необходимо записать уравнение массового баланса:

(13)
${{\phi }_{2}}\frac{{d{{S}_{2}}{{C}_{2}}}}{{dt}} + \frac{{d{{N}_{2}}}}{{dt}} + \lambda {{S}_{2}}{{\phi }_{2}}{{C}_{2}} + \lambda {{N}_{2}} = f{{C}_{1}} - {{q}_{2}}{{C}_{2}},$
где ${{\phi }_{2}}$ – активная (дренируемая) пористость; концентрация C1 определяется из решений (7) или (8); ${{N}_{2}}$ – текущая активность радионуклида в сорбированном состоянии в грунтовой воде на единицу площади [Бк/м2].

Используя линейную сорбционную модель, ${{N}_{2}} \equiv K_{{d2}}^{{}}{{\rho }_{b}}{{S}_{2}}{{C}_{2}}$, можно преобразовать уравнение (13) к виду:

(14)
${{\phi }_{2}}{{R}_{2}}\frac{{d{{S}_{2}}{{C}_{2}}}}{{dt}} + \lambda {{\phi }_{2}}{{R}_{2}}{{S}_{2}}{{C}_{2}} = f{{C}_{1}} - {{q}_{2}}{{C}_{2}},$
где ${{R}_{2}} = 1 + K_{d}^{{}}{{\rho }_{b}}{\text{/}}{{\phi }_{2}}$; $K_{d}^{{}}$ – коэффициент сорбционного распределения для водоносных пород.

В упрощенной постановке, при больших ${{\mu }_{2}}$ (${{\mu }_{2}} \to \infty $), будет выполняться равенство ${{q}_{2}} = f$, и уравнение (13) принимает вид:

(15)
${{\phi }_{2}}{{R}_{2}}S_{2}^{0}\frac{{d{{C}_{2}}}}{{dt}} + \lambda {{\phi }_{2}}{{R}_{2}}S_{2}^{0}{{C}_{2}} = f({{C}_{1}} - {{C}_{2}}),$
где $S_{2}^{0}$ – средняя мощность водоносного горизонта.

Обыкновенное дифференциальное уравнение (15) может быть решено численно относительно ${{C}_{2}}(t)$ при ${{С}_{1}}(t)$ и $f(t)$, определяемых из решения задачи формирования поверхностного стока. В этом уравнении $f{{C}_{1}}$ играет роль функции-источника, имеющей кусочно-однородное представление (10) или (11).

6. МОДЕЛЬ ФРАГМЕНТАЦИИ АТМОСФЕРНЫХ ОСАДКОВ

Для того, чтобы воспользоваться результатами предыдущих разделов, где были получены концентрации ${{C}_{1}}$ и ${{C}_{2}}$, в прогнозных оценках необходимо определить потоковые компоненты q и i, для известной (заданной) функции r, т.е. воспользоваться одной из гидрологических моделей генерации поверхностного стока. В этом свете может оказаться полезным опыт анализа гидрологических данных на основе “модели номерных кривых” (curve number (CN) model), разработанной Службой охраны почв США (Soil Conservation Service (SCS), USA). Эта модель, основанная на обобщении значительного объема полевых данных, относится к классу эмпирических моделей бассейнов стока и используется для прогнозирования прямого стока или инфильтрации при избыточных дождевых осадках [16].

6.1. Стандартная модель SCS-CN

В модели используется две кумулятивные функции: слой осадков, P [мм], когда дождь рассматривается как событие, продолжительность которого и интенсивность в явном виде не учитываются; мощность избыточного слоя осадков или прямого (поверхностного, без подземной компоненты) стока, Q [мм], который является откликом на дождевые осадки [16, 5, 4, 10, 13 ].

Авторами метода вводится понятие потенциал поверхностного стока, определяемый как разница двух характеристик $P - {{I}_{a}}$, где ${{I}_{a}}$ – начальная доля осадков (initial abstraction), не участвующая в формировании поверхностного стока [мм], – параметр, отражающий тот факт, что пока $P \leqslant {{I}_{a}}$, $Q = 0$.

Характеристика ${{I}_{a}}$ включает кратковременные эффекты, а именно задержку осадков растительным покровом, испарение, поверхностную аккумуляцию влаги в локальных депрессиях, инфильтрационное впитывание влаги. Мощность избыточного слоя осадков или прямой поверхностный сток Q всегда меньше или равна мощности слоя осадков P. После того, как в водосборном бассейне формируется поверхностный сток, количество удерживаемой воды характеризуется величиной куммулятивной инфильтрации F [мм]. Значение F меньше или равно максимальной толщине слоя воды, удерживаемой водосборным бассейном S [мм] (максимальная величина впитывания).

В соответствии с исходным постулатом, отношения фактических значений Q и F к потенциально возможным $P - {{I}_{a}}$ и S должны быть равны, т.е.:

$\frac{Q}{{P - {{I}_{a}}}} = \frac{F}{S}$.

Решая это уравнение совместно с балансовым уравнением относительно Q,

(17)
$P = Q + {{I}_{a}} + F$
получаем

$Q = \frac{{{{{(P - {{I}_{a}})}}^{2}}}}{{P - {{I}_{a}} + S}},\quad {\text{при}}\quad P > {{I}_{a}},$
$Q = 0,\quad {\text{в}}\;{\text{других}}\;{\text{случаях}}\quad ({\text{при}}\;P \leqslant {{I}_{a}}).$

Аналогично для куммулятивной функции (при $P > {{I}_{a}}$):

(19)
$F = \frac{{S(P - {{I}_{a}})}}{{P - {{I}_{a}} + S}}.$

В изначальной формулировке метод SCS-CN предполагает, что уравнение (19) выполняется в конце дождевого события (ливневых осадков), т.е. выбранный временной масштаб соответствует длительности события, обычно – длительности ливневых осадков. Сумма $F + {{I}_{a}}$ – это часть выпавших осадков, не перешедшая в поверхностный сток. Обобщение значительного объема полевых данных (в основном по бассейнам стока в США) позволило прийти к выводу, что ${{I}_{a}}$ может быть определено как доля потенциально заполненной водой емкости почвы:

(20)
${{I}_{a}} = {{\lambda }_{a}}{\kern 1pt} S,$
где обычно принимается ${{\lambda }_{a}} = 0.2$ [13].

Эмпирическая формула, полученная для оценки способности водосборного бассейна удерживать воду, имеет вид [16]:

(21)
$S = 25.4\left( {\frac{{1000}}{{{\text{CN}}}} - 10} \right),$
где CN – номер кривой поверхностного стока, характеризующей тип почвы, растительный покров и условия землепользования, а также некоторые другие параметры, определяющие способность водосборного бассейна удерживать воду [8].

Существуют различные методики, позволяющие определить номер кривой стока для различных типов почв, покрытых растительностью, и условий использования земель. CN принимает значения от 0 (почва поглощает все выпадающие на ее поверхность осадки) до 100 (поглощение отсутствует).

6.2. Расширение дискретной SCS-CN модели

В вышеописанной классической SCS-CN модели P и Q – дискретные функции, характеризующие единичное событие, с которым в данном контексте ассоциируется дождевой период. В то же время нет формальных ограничений для рассмотрения функций P, Q и F в качестве непрерывных и дифференцируемых (в области $P > {{I}_{a}}$) [11, 10 ]. Их временные производные являются характеристиками удельных расходов: $r = dP{\text{/}}dt$ (интенсивность осадков), $q = dQ{\text{/}}dt$ (поверхностный сток) и $i = dF{\text{/}}dt$ (скорость инфильтрации). Если это так, то справедливы дифференциальные тождества:

(22)
$\frac{{dQ}}{{dP}} = \frac{{dQ{\text{/}}dt}}{{dP{\text{/}}dt}} = \frac{q}{r},\quad \frac{{dF}}{{dP}} = \frac{{dF{\text{/}}dt}}{{dP{\text{/}}dt}} = \frac{i}{r},$
(в двухкомпонентной модели подразумевается $q \equiv {{q}_{1}}$). Отсюда получаем:

(23)
$\begin{gathered} q = \frac{{(P - {{I}_{a}})(P - {{I}_{a}} + 2S)}}{{{{{(P - {{I}_{a}} + S)}}^{2}}}}r, \\ i = \frac{{{{S}^{2}}}}{{{{{(P - {{I}_{a}} + S)}}^{2}}}}r,\quad P > {{I}_{a}}. \\ \end{gathered} $

Здесь для упрощения записи опущен цифровой индекс у расхода потока ($q \equiv {{q}_{1}}$).

В формулe (23) входит характеристика $P - {{I}_{a}}$ (эффективные осадки), которая может быть определена из (19):

(24)
$P - {{I}_{a}} = \frac{{FS}}{{S - F}}.$

Подставляя (24) в (23), находим:

(25)
$q = \bar {V}(2 - \bar {V})r,\quad i = {{(\bar {V} - 1)}^{2}}r,\quad \bar {V} > {{I}_{a}}{\text{/}}S = {{\lambda }_{a}},$
где $\bar {V} = \frac{F}{S}$. Напомним, что $S = f({\text{CN}})$, более того, часто рекомендуемое значение ${{\lambda }_{a}}$ равно 0.2. Основываясь на (25), легко показать, что $q(\bar {V})$ и $i(\bar {V})$ связаны дифференциальным тождеством:

(26)
$\frac{{dq}}{{d\bar {V}}} + \frac{{di}}{{d\bar {V}}} = 0.$

Из балансового уравнения (17) при постоянных ${{I}_{a}}$ и $S$ следует уравнение неразрывности:

(27)
$\frac{{d\bar {V}}}{{dt}} = \frac{1}{S}(r - q).$

Подставляя в это уравнение выражение для $q(\bar {V})$ из (25), приходим к обыкновенному дифференциальному уравнению:

(28)
$\frac{{d\bar {V}}}{{dt}} = \frac{r}{S}{{(1 - \bar {V})}^{2}}.$

Разделяя переменные и интегрируя (28) при r = const (нулевые начальные условия), получаем решение для нахождения функции $\bar {V}$:

(29)
$\bar {V} = \frac{{\alpha t}}{{1 + \alpha t}},\quad \alpha = \frac{r}{S}.$

Нулевое начальное условие $\bar {V}(t = 0) = 0$ подразумевает отсутствие поверхностного стока в начале дождевого события. В более общем случае, $\bar {V}(t = 0) = {{\bar {V}}_{0}}$, т.е. когда выпадение осадков происходит при наличии слоя воды на поверхности, интегрирование (28) дает:

$\bar {V} = \frac{{{{{\bar {V}}}_{0}} + (1 - {{{\bar {V}}}_{0}})\alpha t}}{{1 + (1 - {{{\bar {V}}}_{0}})\alpha t}}.$

Подстановка (29) и (29а) в (25) позволяет определить удельный расход поверхностного стока и инфильтрацию в любой расчетный момент времени. Как следует из (25), в начальные моменты при ${{\bar {V}}_{0}} = 0$, $q \to 0$, $i \to r$; для больших моментов времени – картина обратная $q \to r$, $i \to 0$, что отвечает физике процесса. Особенностью решения (25) является обнуление функций q и i в момент прекращения дождя (когда r = 0), а также скачкообразный их рост или падение при дискретном приращении или убыли интенсивности дождевых осадков.

Вследствие эвапотранспирации воды из верхней части почвенного разреза существует различие между величиной инфильтрации i и величиной питания грунтовых вод f. Данный процесс может быть учтен в балансовой формуле (28) посредством включения дополнительного члена-источника [11]: $E = g(\bar {V},\;PET)$, где PET – потенциальная эвапотранспирация [м]. Очевидно, что рост $\bar {V}$ делает влагу почвы для процесса эвапотранспирации более доступной, причем при приближении $\bar {V}$ к предельному значению, равному 1, интенсивность эвапотранспирации $E\, = \,\Delta (ET){\text{/}}\Delta t$ равна величине ${{E}_{p}} = \Delta (PET){\text{/}}\Delta t$ [м/сут]. Этому условию, например, удовлетворяет простейшая линейная зависимость [1]:

(30)
$E = {{E}_{p}}\bar {V},$
так что расширенное уравнение массового баланса имеет вид:
(31)
$\frac{{d\bar {V}}}{{dt}} = \frac{r}{S}{{(1 - \bar {V})}^{2}} - \frac{{{{E}_{p}}}}{S}\bar {V}$
или

$\frac{{d\bar {V}}}{{dt}} = \alpha {{(1 - \bar {V})}^{2}} - \beta \bar {V},\quad \alpha = \frac{r}{S},\quad \beta = \frac{{{{E}_{p}}}}{S}.$

Уравнение (31а) также может быть решено посредством разделения переменных, что дает:

(32)
$\bar {V} = \frac{{1 - \bar {Q}}}{{{{g}_{1}} - {{g}_{2}}\bar {Q}}},$
где

$\begin{gathered} \bar {Q} = \frac{{1 - {{g}_{1}}{{V}_{0}}}}{{1 - {{g}_{2}}{{V}_{0}}}}\exp ( - \sqrt A t);\quad {{g}_{1}} = \frac{g}{{1 - P}}; \\ {{g}_{2}} = \frac{g}{{1 + P}};\quad g = \frac{{2\alpha }}{{2\alpha + \beta }};\quad P = \frac{{\sqrt A }}{{2\alpha + \beta }}; \\ A = \beta (4\alpha + \beta ). \\ \end{gathered} $

Таким образом могут корректироваться значения функций $q \equiv {{q}_{1}}$ и $i$ в (25) с поправкой на эвапотранспирацию. В этом случае инфильтрация может рассматриваться как величина питания подстилающего водоносного горизонта (f = i).

Решение (32) описывает частный случай эвапотранспирации при отсутствии атмосферных осадков (r = 0). В этом случае оно принимает вид решения уравнения кинетики первого порядка (распада):

$\bar {V} = {{\bar {V}}_{0}}\exp \left( { - \frac{{{{E}_{p}}}}{S}t} \right).$

Если r аппраксимируется ступенчатой функцией, то компоненты q и i в (25) на конец n-го расчетного периода определяются из представленных выше выражений для функции $\bar {V}$, в которых

(33)
$a \to {{a}_{n}} = \frac{{{{r}_{n}}}}{S}.$

7. ЧИСЛЕННО-АНАЛИТИЧЕСКИЕ РАСЧЕТЫ ДЛЯ ОДИНОЧНОГО БАССЕЙНА

7.1. Концентрационные функции

Математическая модель предполагает, что во время периодических дождевых событий, приводящих к образованию слоя воды на поверхности, растворенные радионуклиды (продукты выбросов с АЭС) из почвы поступают как в склоновый сток, так и резервуар грунтовых вод, причем в консервативных оценках можно допустить, что процесс такого вертикального массообмена протекает весьма быстро. Быстрый транспорт влаги – достаточно распространенное явление в структурированных (макропористых, растрескавшихся) почвах.

Суперпозиция двух массовых потоков позволяет определить среднюю (общую) концентрацию ${{C}_{q}}$ в выходном створе гидрологического бассейна в момент выпадения осадков и формирования расходов потоков (q1 и q2), а также суммарный массовых поток ${{j}^{{tot}}}$:

(34)
${{C}_{q}} = \frac{{{{q}_{1}}{{C}_{1}} + {{q}_{2}}{{C}_{2}}}}{{{{q}_{1}} + {{q}_{2}}}},\quad {{j}^{{tot}}} = ({{q}_{1}} + {{q}_{2}}){{C}_{q}}.$

Алгоритм расчетной модели предполагает аппроксимацию графика выпадения дождевых осадков P = P(t) ступенчатой функцией с постоянным шагом $\Delta t$ = 1 сут. В решениях для концентрационной функции в слое смешения C1 (10, 11) эффективные осадки определяются из выражения $r = [P - {{I}_{a}}(S)]{\text{/}}\Delta t$. Массив значений P(t) трансформируется в массивы ${{q}_{n}}$ и ${{i}_{n}}$ в соответствии с SCS-CN концепцией (см. разд. 6). Решение задачи производится последовательно для каждого временного шага, причем значения ${{q}_{{1n}}}$ и ${{q}_{{2n}}}$ (расходы поверхностного и подземного потока), ${{C}_{1}}$ и ${{С}_{2}}$ (активности в воде поверхностного и подземного стока), полученные на предыдущем временном шаге, использовались для расчета значений этих функций на текущем шаге. В периоды, когда не происходит генерация поверхностного стока, инфильтрационная вода идет на разбавление радионуклидов в слое смешения, и весь массовый поток поступает в грунтовые воды. При образовании поверхностного стока инфильтрационный массовый поток сокращается, но увеличивается интенсивность выноса радионуклидов с водосборной площади склоновыми водами.

Собственно функции C1 и C2 непрерывны, а наблюдаемые колебания концентраций связаны с разбавлением радионуклидов в слое смешения в период дождей. В свою очередь, потоковые компоненты qi(t) дискретны, что связано с периодичностью дождевых выпадений, и имеют многочисленные разрывы – нулевые значения q1 и q2 в засушливые периоды. В эти периоды величина Cq (34) не может быть определена. Однако при осреднении значений Cq за определенный промежуток времени $T = \sum\nolimits_n^N {\Delta {{t}_{n}}} $ (N – количество точек осреднения), включающий как дождливые дни (P > Ia), так и дни без осадков (или когда P < Ia), может быть рассчитано среднее значение активности:

${{C}_{{aver}}} = \frac{{\sum\limits_n^N {({{q}_{{1n}}}{{C}_{{1n}}} + {{q}_{{2n}}}{{C}_{{2n}}})} }}{{\sum\limits_n^N {({{q}_{{1n}}} + {{q}_{{2n}}})} }},$
на которое мы и будем в дальнейшем ориентироваться, задавшись временным диапазоном T.

7.2. Пример

Представленная модель позволяет проиллюстрировать особенности радиоактивного загрязнения природных (поверхностных и подземных) вод в различных ландшафтных условиях, базируясь на CN-концепции. Пусть дождевые осадки характеризуются синтетическим гидрографом, представленным на рис. 4. Параметры почвенного слоя и подстилающего горизонта следующие: $\theta {{d}_{e}}{{R}_{1}} = 1$ м; ${{\phi }_{2}}S_{2}^{0}{{R}_{2}} = 5$ ${{j}^{{tot}}} = ({{q}_{1}} + {{q}_{2}}){{C}_{q}}$; $E(ET) = 0$. Рассматриваются условия нормальной эксплуатации и аварийный режим выброса радионуклидов в атмосферу, соответственно: $N{\text{'}} = {{10}^{{ - 4}}}$ Бк/м2 (НЭ); ${{N}_{0}} = {{10}^{3}}$ Бк/м2 (ЗПА). Для наглядности будем считать, что выброс представлен двумя радионуклидами – Cs-137 и Cs-134, имеющими периоды полураспада 30.17 и 2.06 года ($\lambda $ = 6.29 × × 10–5 сут–1 и $\lambda $ = 9.22 × 10–4 сут–1).

Рис. 4.

Слой атмосферных осадков (P) и эффективные (r) атмосферные осадки (мм/сут): а – CN = 90, б – CN = 70.

Как видно на рис. 4, номер кривой (CN) контролирует величину эффективных осадков, а также их перераспределение в водных потоках. Чем выше значение CN, тем больше влаги участвует в процессе переноса радионуклидов на поверхности почвы.

В случае НЭ при больших CN радионуклиды не успевают накапливаться на почве, и соответствующие кривые C1, C2 и Caver на графиках C(t) (рис. 5а, б) располагаются под кривыми с более низкими значениями CN. При ЗПА рост CN приводит к более быстрому истощению запасов аварийных радионуклидов в почвенном слое (рис. 6а, б), они активно выносятся поверхностным стоком с водосборной площади, что способствует более быстрому появлению радионуклидов в выходном створе бассейна. Наоборот, пониженные значения CN отражают повышенную роль инфильтрационной компоненты в общем водном балансе, поэтому основная масса радионуклидов поступает в грунтовые воды. Вынос радионуклидов из бассейна происходит в основном с подземными водами (при низких CN) с большой временной задержкой (см. рис. 6а, б). Водоносный горизонт становится источником долговременного загрязнения природных вод.

Рис. 5.

Нормальная эксплуатация. а, б – графики концентраций C1 (сплошные кривые) и C2 (пунктирные кривые): а – Cs-137, б – Cs-134; в, г – средние концентрации Caver; в – Cs-137, г – Cs-134. Трехмесячное осреднение. Числа у кривых – номер кривой стока.

Рис. 6.

Последствия аварии. а, б – графики концентраций C1 (сплошные кривые) и C2 (пунктирные кривые): а – Cs-137, б – Cs-134; в, г – средние концентрации Caver: в – Cs-137, г – Cs-134. Трехмесячное осреднение. Числа у кривых – номер кривой стока (CN).

Графики средних концентраций Caver (формула (34б), рис. 5в, г и 6в, г) построены для трехмесячного (квартального) периода осреднения (N = 90, T = 90 сут) (суточные кривые, как отмечалось, существенно дискретны и имеют многочисленные разрывы в связи с нулевыми значениями функций qi в засушливые периоды). Для условий ЗПА (см. рис. 6в, г) восходящая ветвь кривой Caver образуется вследствие загрязнения поверхностного стока, нисходящая ветвь – вследствие миграции радионуклидов в водоносном горизонте.

8. МОДЕЛЬ СВЯЗНЫХ ПОТОКОВ ДЛЯ “ЦЕПОЧЕК” БАССЕЙНОВ (СУБ-БАССЕЙНОВ)

8.1. Принцип суперпозиции

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

Рис. 7.

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

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

(35)
$C(t) = \frac{{\sum\limits_{i = 1}^k {{{j}_{i}}(t)} }}{{Q(t)}},$
где ${{j}_{i}}(t)$ – поток радионуклидов (сформированный поверхностным стоком) на замыкающем створе i-го суб-бассейна:
(36)
${{j}_{i}}(t) = [{{C}_{1}}(t){{q}_{{1i}}} + {{C}_{2}}(t){{q}_{{2i}}}]{{F}_{i}},$
где ${{F}_{i}}$ – площадь i-го суб-бассейна.

Общий расход, Q, в формуле (35) характеризует поток воды, разбавляющий загрязненный радионуклидами поверхностный сток. В зависимости от условий стока из водосборного бассейна и водного объекта, качество которого подлежит оценке, общий расход может быть представлен двумя способами. Если загрязненный поверхностный сток разбавляется в потоке воды, образовавшемся в том же моделируемом водосборном бассейне (или группе бассейнов), то:

(37)
$Q(t) = \sum\limits_i {({{q}_{{1i}}} + {{q}_{{2i}}}){{F}_{i}}} .$

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

(38)
$Q(t) = {{Q}_{s}}(t) + \sum {({{q}_{{1i}}} + {{q}_{{2i}}}){{F}_{i}}} ,$
где ${{Q}_{s}}(t)$ – расход транзитного потока.

Разрывы потоковых компонент ${{q}_{{1i}}}$ и ${{q}_{{2i}}}$ обусловливают необходимость представления результатов расчетов в виде значений средней активности:

${{C}_{{aver}}} = \frac{{\sum\limits_n^N {\sum\limits_{i = 1}^k {{{j}_{{in}}}(t)} } }}{{\sum\limits_n^N {{{Q}_{n}}} + \sum\limits_n^N {\sum\limits_{i = 1}^k {({{q}_{{1in}}} + {{q}_{{2in}}}){{F}_{i}}} } }},$
где N – количество усредняемых активностей.

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

8.2. Пример

В качестве примера рассмотрим участок бассейна р. Дунай, попадающий в 30-километровую защитную зону проектируемой АЭС (рис. 8).

Рис. 8.

Водосборные суб-бассейны в пределах 30-километровой зоны (полигоны, ограниченные сплошными черными линиями). Пунктирный круговой контур – граница 30- километровой зоны АЭС; жирная линия – внешняя граница модельной области; цифры на карте – номера суб-бассейнов; серый цвет – SW радиоактивный след.

8.2.1. Геоморфологическая и гидрологическая характеристика района. Правобережная водосборная площадь бассейна р. Дунай с несколькими своими притоками ограничена Среднезадунайским горным массивом. Левобережье р. Дунай граничит с бассейном р. Тиса. Водораздел между двумя реками весьма слабо выражен в рельефе. Здесь у притоков Дуная нет хорошо выраженных речных долин. Гидросеть состоит в основном из искусственных водотоков. Средние модули стока с бассейнов меняются от 1.1 до 2.2 л/с/км2.

В 30-километровую зону (см. рис. 8) попадает участок бассейна р. Дунай (1 порядок) с ее главным притоком р. Шио (2 порядок). В пределах этой зоны суб-бассейны, связанные с потоками 3 порядка, могут быть сгруппированы в цепочки в соответствии с конфигурацией и взаимосвязью дренирующих эти бассейны водотоков (рис. 9).

Рис. 9.

Структура путей миграции радионуклидов в 30-километровой зоне в соответствии со структурой связности суб-бассейнов. Числа в прямоугольниках – номера суб-бассейнов. Закрашенная область – SW радиоактивный след.

Траектории (структура связей) путей миграции техногенных радионуклидов в пределах рассматриваемой области будет различной для разных сценариев радиоактивных выпадений. Для НЭ, когда вся 30-километровая зона вокруг АЭС загрязнена, любая из контрольных точек, связанная с каким-либо водотоком, интегрирует массовые потоки, формирующиеся в пределах строго определенного набора суб-бассейнов, имеющих различную плотность загрязнения, контролируемую розой ветров. Для аварийных сценариев (ПА и ЗПА) траектория маршрута радионуклидного транспорта зависит от геометрии и расположения радиоактивного следа. Для заданной контрольной точки в зависимости от расположения радиоактивного следа эта траектория будет разной. Пример выделения цепочки бассейнов для SW радиоактивного следа представлен на рис. 9, где показаны контрольные точки. Отвечающие им гидрографы характеризуют следующие средние расходы потоков (м3/с): р. Дунай – 2265, р. Шио – 8.5, канал Надор – 3.8.

8.2.2. Исходные данные для расчетов. В методе номерных кривых, являющимся элементом бассейновой модели радионуклидного транспорта, используется классификация почв по гидрологическим типам (hydrologic soil groups), в соответствии с которой можно установить минимальную величину инфильтрации. Она достигается на обнаженной почве после продолжительного ее увлажнения. Гидрологические типы почв – A, B, C и D, являются одним из элементов методики определения номера кривой поверхностного стока [16]. Почвенные группы объединяют почвы, имеющие различный потенциал поверхностного стока, в порядке возрастания от A до D.

Каждый из 9 суб-бассейнов в пределах SW радиоактивного следа был разбит на полигоны (площадью Fi), характеризующиеся в зависимости от условий землепользования локальными номерами кривой стока CNi (табл. 2).

Таблица 2.

Характеристики водосборных бассейнов в пределах SW радиоактивного следа

№ суб-бассейна Тип
почвы
Земли населенных пунктов, Fi, км2 CNi Земли промыш-ленности, Fi, км2 CNi Сельхоз зем-ли, Fi, км2 CNi Сады,
Fi, км2
CNi
2 D 1.95 84 0 93 33.23 87 0 82
11 D 0.22 84 0.18 93 24.49 87 0 82
12 D 0.97 84 0.2 93 38.60 87 0.44 82
13 D 0.96 84 0 93 30.06 87 0 82
14 D 0.7 84 0 93 59.80 87 0.11 82
27b D 0.85 84 0 93 52.11 87 0.88 82
30b B 3.53 69 0.94 88 77.86 78 0 65
33 C 1.55 79 0 91 49.69 83 0.7 76
39b C 0 79 0.9 91 32.13 83 0.5 76
№ суб- бассейна Тип
почвы
Лес, кустарник,
Fi, км2
CNi Луга, Fi, км2 CNi Болота, Fi, км2 CNi Водная поверхность, Fi, км2 CNi
2 D 15.52 77 0 78 0 78 2.5 0
11 D 11.22 77 0 78 0 78 0 0
12 D 22.34 77 0 78 0 78 0.9 0
13 D 4.36 77 0 78 0 78 1.1 0
14 D 12.64 77 0 78 0 78 0.9 0
27b D 7.95 77 0 78 0 78 2 0
30b B 12.5 56 1.5 58 0.75 58 0.7 0
33 C 10.5 70 1.52 71 0.3 71 0 0
39b C 7.2 70 0 71 0 71 0 0

Средневзвешенные номера кривых (CN) для 9 суб-бассейнов (модельных водосборных бассейнов) в пределах SW радиоактивного следа оценивались в зависимости от площади определенного вида землепользования и почвенного покрова в процентном отношении к общей площади бассейна (табл. 3).

(39)
${\text{CN}} = \sum\limits_{i = 1}^{\text{9}} {{\text{C}}{{{\text{N}}}_{i}}{{F}_{i}}} {\text{/}}\sum\limits_{i = 1}^{\text{9}} {{{F}_{i}}} .$
Таблица 3.

Средние значения модуля стока (QA) и номера кривых стока (CN) для суб-бассейнов в пределах SW радиоактивного следа

№ суб-бассейна i F, км2 QA, л/с/км2 CN
2 1 53.19 1.7 80
11 2 36.11 1.7 84
12 3 63.45 1.7 82
13 4 36.48 1.5 83
14 5 74.15 1.7 84
27b 6 63.79 1.3 83
30b 7 97.78 1.2 74
33 8 64.26 1.1 80
39b 9 40.73 1.05 81

При отсутствии прямых (экспериментальных) определений сорбционной константы Kd (коэффициент сорбционного распределения) рекомендуется использовать корреляционные связи между Kd , емкостью катионного обмена CEC, содержанием глинистых частиц (CC) и pH почв [15]. Нами были использованы скрининговые (экспертные) значения, приведенные в табл. 4. Они заимствованы из работы [7, табл. 5.6 и 5.13 ] и отвечают минимальным значениям Kd для грунтов с показателями CEC = 3–10 мг-э/100 г и СС = = 4–20%, причем ${{K}_{d}} = {{K}_{{d1}}} = {{K}_{{d2}}}$. Полагалось θ = ϕ2 = 0.2, ${{\rho }_{b}} = 2.05$ г/см3, ${{d}_{e}}$ = 5 см, $S_{2}^{0} = 5$ м.

Таблица 4.

Характеристики аварийного выброса (ТЗПА, сценарий DEC2); средние значения коэффициента распределения Kd

Радио-нуклид T1/2, годы M, ТБк N0, Бк/м2 Kd, см3
Cs-137 30.17 8.4 1.59E+04 70
Cs-134 2.06 11.4 2.15E+04 70
Sr-90 28.79 0.83 1.57E+03 15

M – суммарная активность выброса; N0 – расчетная поверхностная активность осажденных радионуклидов (при $\sum\nolimits_{i = 1}^{\text{9}} {{{F}_{i}}} $ = 529.9 км2).

Прогнозы выполнялись для средних расходов водотоков (см. табл. 2). Потенциальная эвапотранспирация (PET) 1000 мм (среднегодовое значение).

Наконец, анализ данных метеонаблюдений в пределах 30-километровой зоны действующей АЭС выполнен для периода с 1950 по 2018 гг. Достаточно типичным является гидрограф осадков P (мм/сут) за последний 10-летний период (2009–2018 гг.), который использовался нами в качестве базового в модельных расчетах.

На рис. 10 показаны суммы осадков по декадам. На самом деле, при моделировании использовались суточные распределения, пересчитанные по данным 6-часовых замеров. В качестве примера на рис. 11 показан график функции P за 2010 г., который свидетельствует о том, что внутригодовое распределение атмосферных осадков неоднородно. Как видно, достаточно продолжительные засушливые периоды сменяются относительно кратковременными периодами выпадения довольно интенсивных осадков. Для долгосрочных оценок с горизонтом 60 лет этот 10-летний период адаптировался к данному (прогнозному) масштабу времени посредством операции цикличного повторения.

Рис. 10.

Декадные суммы осадков (P), мм (2009–2018 гг.)

Рис. 11.

Суточные распределения осадков за 2010 г.

8.2.3. Результаты расчетов. В соответствии с разд. 2–6, основным гидрологическим механизмом, определяющим загрязнение речных вод, является склоновый сток, формирующийся при выпадении жидких атмосферных осадков. Радионуклиды с поверхности радиоактивно загрязненной почвы поступают в слой воды, стекающей по почве, и переносятся этим потоком к руслу реки. Часть радионуклидов переносится инфильтрационными водами из почвенного слоя в грунтовые воды. При определении удельной активности ${{C}_{{aver}}}$ в качестве расчетной модели может использоваться решение (38а), где $Q$ – суммарный расход (38), причем ${{Q}_{s}}(t)$ – расход транзитного потока Дуная или канала Надор.

Расчеты, выполненные на основе разработанной модели, иллюстрируют значимость различных физических механизмов в радиационном воздействии на водотоки первого (р. Дунай) и второго (канал Надор) порядка. При заданном гидрографе атмосферных выпадений, абсолютные активности ПД в воде во многом определяются степенью разбавления радионуклидов склоновых и подземных вод в речной воде. В Дунае она значительно выше, и поэтому активность загрязненной воды невысока. На абсолютные значения активности влияет также сорбция: благодаря пониженным значениям Kd, несмотря на то, что плотность загрязнения водосборов Sr-90 заметно ниже (почти на десятичный порядок) плотности выпадения Cs-134 и Cs-137 (см. табл. 4), пиковые значения активности в воде всех изотопов в первые годы не слишком сильно различаются. Кроме того, сравнение графиков активностей относительно долгоживущих изотопов Cs-137 и Sr-90 (рис. 12 и 13) показывает, что скорость падения активности Sr-90 в поверхностных водах выше, что объясняется его меньшей сорбируемостью в почвенном слое и породах водоносного горизонта (соответственно, повышенной способностью к вымыванию из почвы и породы). Весьма заметно влияние радиоактивного распада: относительно малый период полураспада Cs-134 (по сравнению с Cs-137 и Sr-90) способствует практически полной очистке территории от этого радионуклида уже примерно через 12 лет (рис. 14).

Рис. 12.

Результаты прогнозов загрязнения поверхностных вод Cs-137 при ЗПА аварии (DEC2) в контрольных (расчетных) точках: а – 1D (р. Дунай), б – 3 (канал Надор) (юго-западный радиоактивный след SW). Тонкие кривые – осреднение 90 сут; жирные кривые – среднегодовые значения.

Рис. 13.

Результаты прогнозов загрязнения поверхностных вод Sr-90 при ЗПА аварии (DEC2) в контрольных (расчетных) точках: а – 1D (р. Дунай), б – 3 (канал Надор) (юго-западный радиоактивный след SW). Тонкие кривые – осреднение 90 сут; жирные кривые – среднегодовые значения.

Рис. 14.

Результаты прогнозов загрязнения поверхностных вод Cs-134 при ЗПА аварии (DEC2) в контрольных (расчетных) точка: а – 1D (р. Дунай), б – 3 (канал Надор) (юго-западный радиоактивный след SW). Тонкие кривые – осреднение 90 сут; жирные кривые – среднегодовые значения.

В целом, полученные результаты показывают, что степень радиоактивного загрязнения р. Дунай (главную водную артерию региона) за счет смыва аварийных (ЗПА, DEC2) радионуклидов с поверхности (с учетом их частичного поступления в подземные воды) относительно невысокая: активность радионуклидов в речной воде не превышает 10-15 мБк/л (среднеквартальные значения, см. рис. 12а–14а). Степень загрязнения радионуклидами канала примерно на порядок выше: объемные активности Cs-134, Cs-137 и Sr-90 достигают нескольких сотен мБк/л (см. рис. 12б–14б).

ЗАКЛЮЧЕНИЕ

Оценка воздействия на окружающую среду при проектировании и строительстве АЭС в соответствии с российским и международным экологическим законодательством, в частности, с требованиями МАГАТЭ, предполагает прогнозирование миграции ПД в водных объектах – поверхностных и подземных водах. Трудности, связанные с реализацией соответствующих требований на практике, обусловлены недостаточной разработанностью научно-методической базы, которая учитывала бы специфику функционирования рассматриваемых инженерных объектов. Представленная работа развивает методические и теоретические подходы к прогнозам миграции ПД в водных объектах, подверженных выбросам ПД в атмосферу на АЭС, как в условиях их нормальной эксплуатации, так и при гипотетических авариях разного уровня (проектных и запроектных).

Математические модели различаются особенностями задания типа граничных условий на почве – поверхности осаждения ПД из атмосферы. Для НЭ зона воздействия и плотность выпадения ПД определяется данными статистического осреднения погодных условий в многолетнем разрезе времени. Модель должна учитывать накопление радионуклидов в почве с последующим их вымыванием атмосферными осадками и талыми водами снежного покрова в течение периода времени, измеряемого многими десятилетиями (для российских АЭС – это 60 лет). Для прогнозирования последствий аварийных выброса (ПА и ЗПА) исходят из предпосылки об импульсном характере площадного источника загрязнения поверхности и вероятностной (наперед не заданной) геометрии и ориентации радиационного следа.

Основной принцип задания плановых границ моделей – бассейновый. Прогнозы, связанные как с НЭ станции, так и с оценкой последствий ПА и ЗПА, строятся дифференцированно применительно к каждому из водосборных бассейнов, попадающих в 30-километровую зону влияния АЭС. Преимущества такого бассейнового подхода обусловлены возможностью привлечения гидрографов замыкающих створов бассейна для определения величин компонентов поверхностного (склоновый сток) и подземного (водоносный горизонт) стока.

Зона потенциального воздействия выбросов ПД с АЭС, измеряемая десятками и сотнями км2, является, как правило, недостаточно изученной для выполнения детальных прогнозов распространения ПД на основе сеточных моделей поверхностного стока и фильтрации с распределенными параметрами. В этих условиях более реалистичные (приближенные, экспертные) оценки могут быть получены на основе моделей с сосредоточенными параметрами. В статье предложены аналитические решения и алгоритмы для численно-аналитических расчетов миграции ПД, описывающие образование слоя воды, насыщенного осажденными ПД, и поступление радиоактивно загрязненных вод в поверхностный водоток и грунтовые воды, которые также движутся к контуру дренирования водосборного бассейна. Расчетные модели интегрированы с методом нумерованных кривых стока [16], что позволяет предсказать цикличность выноса радионуклидов с загрязненных территорий, основываясь на гидрографе дождевых осадков, используя при этом обширную базу данных обобщенных параметров. Эти параметры учитывают природные почвенно-ландшафтные особенности территорий и условия их промышленного, городского и сельскохозяйственного освоения. Номер кривой стока (CN) отражает также перераспределение между двумя компонентами стока – поверхностной (быстрой) и подземной (медленной), предсказывая характер развития процесса радиоактивного загрязнения природных вод в целом.

Список литературы

  1. Шестаков В.М., Поздняков С.П. Геогидрология. М.: ИКЦ “Академкнига”, 2003. 176 с.

  2. Arnold J.G., Allen P.M., Bernhardt G. A comprehensive surface-groundwater flow model // Journal of Hydrology. 1993. V. 142. P. 47–69.

  3. Ahuja L.R. Release of a soluble chemical from a soil to runoff. 1982. V. 25. P. 948–953.

  4. Boughton W.C. A Review of the USDA SCS Curve Number Method // Aust. J. Soil Res. 1989. V. 27. P. 511–523.

  5. Chow V.T., Maidment D.R., Mays L.W. Applied Hydrology. New York: McGraw-Hill, 1988. 572 p.

  6. de Rooij G.H. Aquifer-scale flow equations as generalized linear reservoir models for strip and circular aquifers: Links between the Darcian and the aquifer scale // Water Res. Res. 2013. V. 49. P. 8605–8615.

  7. EPA. Understanding variation in partition coefficient, Kd, values. Washington: Environmental Protection Agency, 1999. V. II.

  8. Erickson T, Stefan H.G. Groundwater recharge from a changing landscape. Project report 490. St Paul: Minnesota Pollution Control Agency, 2007. 112 p.

  9. Gao B., Walter M.T., Steenhuis T.S. et al. Rainfall induced chemical transport from soil to runoff: theory and experiments // J. of Hydrol. 2004. V. 295. P. 291–304.

  10. Geetha K., Mishra S.K., Eldho T.I., Rastogi A.K., Pandey R.P. SCS-CN-based Continuous Simulation Model for Hydrologic Forecasting // Water Resource. Manag. 2008. V. 22 (2). P. 165–190.

  11. Michel C., Vazken A., Perrin C. Soil conservation service curve number method: how to mend a wrong soil moisture accounting procedure // Water Res. Research. 2005. V. 41 (2). P. 1–6.

  12. Nalbandyan A., Ytre-Eide M.A., Thørring H. Potential consequences in Norway after a hy-pothetical accident at Leningrad nuclear power plant. Østerás: Norwegian Radiation Protection Authority, 2012. 30 p.

  13. Plummer A., Woodward D.E. The origin and derivation of Ia/S in the runoff curve number system // Proc. of the International Water Resources Engineering Conference. 1998. P. 1260–1265.

  14. Rumynin V. Overland Flow Dynamics and Solute Transport. Springer International Publishing, 2015. 287 p.

  15. Sheppard S., Lonh J., Sanipelli B. Solid/liquid partition coefficients (Kd) for selected soils and sediments at Forsmark and Laxemar-Simpevarp. Canada: Gustav Sohlenius Geological Survey of Sweden (SGU), 2009.

  16. USDA. Natural Resources Conservation Service. Technical Release 55. Urban Hydrology for Small Watersheds, 1986.

  17. Zhang X.C., Norton L.D., Lei T. et al. Coupling mixinf zone concept with convection-diffusion equation to predict chemical transfer to surface runoff // Transaction of the ASAE. 1999. V. 42 (4). P. 987–994.

Дополнительные материалы отсутствуют.

Инструменты

Геоэкология. Инженерная геология, гидрогеология, геокриология