CC BY
28
УДК 622:510.67:523.98 А.Н. Якубович
МОДЕЛИРОВАНИЕ СРЕДНЕГОДОВОЙ ДИНАМИКИ ПРИХОДА СОЛНЕЧНОЙ РАДИАЦИИ НА ТЕРРИТОРИИ, НАРУШЕННЫЕ ПРИ ПРОИЗВОДСТВЕ ГОРНЫХ РАБОТ
Рассмотрен метод получения численных оценок среднегодовой динамики прихода солнечной радиации на территории региона с использованием модели, учитывающей динамику годового движения Солнца по эклиптике, соответствующую рассматриваемому участку территории, состав и строение идеальной атмосферы, состояние реальной атмосферы, наиболее типичное для рассматриваемого участка территории, особенности ландшафта территории — ее альбедо и ориентация в пространстве участков земной поверхности.
Ключевые слова: солнечная радиация, моделирование годового движения Солнца, самовосстановление почвенного покрова, ГИС-технологии.
Герритория природно-технических комплексов Крайнего Северо-Востока России отличается повышенной уязвимостью к техногенным воздействиям и низкими темпами самовосстановления после нарушения почвенного покрова при добыче полезных ископаемых. При этом только на территории Магаданской области расположено более 2000 разведанных месторождений россыпного золота, около 50% которых в настоящее время активно разрабатывается. Эксплуатация россыпных месторождений связана с полным разрушением почвенного покрова, и, соответственно, уничтожением растительного слоя почв. Период самовосстановления растительного слоя по окончании горных работ в условиях Крайнего Севера может составлять более 100 лет.
Одной из главных природно-климатических особенностей территории Крайнего Северо-Востока России является явно выраженный дефицит приходящей солнечной радиации. Именно ограниченность солнечной радиации является основным фактором, обуславливающим уязвимость формирующихся на почвенном покрове экосистем и значительную продолжительность периода их самовосстановления. В связи с этим представляется обоснованным использовать среднегодовую динамику прихода солнечной радиации в качестве основного показателя, характеризующего самовосстано-вительный потенциал территории при наличии нарушений почвенного покрова техногенного характера. Для получения численных
оценок среднегодовой динамики прихода солнечной радиации на территорию предлагается использовать модель, учитывающую следующие основные факторы:
^ динамику годового движения Солнца по эклиптике, соответствующую рассматриваемому участку территории;
> состав и строение идеальной атмосферы;
^ состояние реальной атмосферы, наиболее типичное для рассматриваемого участка территории;
> особенности ландшафта территории - ее альбедо и ориентация в пространстве участков земной поверхности.
Моделирование годового движения Солнца по эклиптике предполагается выполнять в первой экваториальной системе координат, когда координатами места светила являются 5 (солнечное склонение) и t (часовой угол). Началом счета экваториальных координат является точка весеннего равноденствия (20-21 марта), когда 5 = 0°; продолжительность тропического года (промежутка между двумя весенними равноденствиями) составляет 365 суток 5 часов 48 минут 46 секунд, и в настоящей работе обозначается как Туеаг = 31 556 926 с [1].
Необходимо отметить некоторую неравномерность движе-ния Земли по орбите, и, соответственно, неравномерность перемещения Солнца по эклиптике. В результате промежутки времени между последовательными характерными точками тропического года (равноденствиями и солнцестояниями) неодинаковы, что необходимо учитывать при построении модели движения Солнца по эклиптике.
Принимая за основу периодический характер изменения солнечного склонения во времени (с периодом, равным Туеаг), учитывая фактические отклонения от строго периодического закона изменения 5 в пределах одного периода, принимая за начальный момент времени 0 часов 00 минут 00 секунд 1 января и предполагая наступление момента весеннего равноденствия в 0 часов 00 минут 00 секунд 21 марта (через 6 825 600 с после начального момента времени), получаем зависимость, моделирующую изменение склонения во времени:
5 = S • ko = £ • sin
2^
(AT - 6 825 600) • kv
T
year
(1)
где АТ - количество секунд с 0 часов 00 минут 00 секунд 1 января; kv - поправочный коэффициент, учитывающий неравномерность движения Земли по орбите.
В качестве опорных точек для определения kv принимаются характерные точки тропического года (равноденствия и солнцестояния). Зависимость kv от времени принята в виде полинома 5 степени; окончательно выражение для поправочного коэффициента:
kv = 1 + 0.001(-7.07 + 64.4662 -155.9004
( Л 2
2 •AT
+117.7263
f V 2 AT
T
у year J
-35.0415
Tyear year
( V 2 AT
Tyear у year J
+
Tyear у year J
+3.5476
f V 2 AT
Tyear у year J
)
Вторая координата, определяющая положение Солнца на небесной сфере в экваториальной системе координат - часовой угол t - изменяется равномерно и пропорционально времени. Возможные значения часового угла: 0° < t < 360° ; за 1 час значение t изменяется на 15°. Время, измеряемое часовым углом Солнца (солнечное время), отсчитывается от момента наивысшего положения светила (верхней кульминации), когда t = 0° . Поскольку начало гражданских суток отнесено к полуночи, гражданское время на 12 часов (180° часового угла) больше солнечного времени. С учетом вышеизложенного, для моделирования изменения часового угла Солнца во времени предлагается зависимость: mod(АT, 86 400)
t = ж + -
-•2ж
(2)
86 400
где mod - остаток от деления AT на количество секунд в сутках (86 400 с).
Таким образом, при моделировании годового движения Солнца по эклиптике в качестве единственного независимого аргумента используется AT - количество секунд с 0 часов 00 минут 00 се-
кунд 1 января. При любом значении АТ координаты 5 и t, характеризующие положение Солнца в первой экваториальной системе координат, определяются с помощью (1) и (2). Данные зависимости учитывают основные закономерности движения Земли по орбите вокруг Солнца, в том числе и неравномерность скорости движения в пределах тропического года, и сочетают простоту численной реализации с достаточной для целей настоящего исследования точностью.
Предлагаемая модель поглощения и рассеяния солнечной радиации в атмосфере основывается на учете наиболее значи-мых факторов, влияющих на перенос солнечной радиации: поглощение кислородом, водяным паром и атмосферным озоном, а также моде-кулярное рассеяние и рассеяние на атмосферных аэрозолях [2].
Основные полосы поглощения молекулярного кислорода приходятся на ультрафиолетовую область спектра с длинами волн X е [0,13 мкм; 0,24 мкм\. Поскольку солнечная радиация с длинами волн X < 0,3 мкм полностью поглощается атмос-ферным озоном, влияние данных полос поглощения кислорода в настоящей работе не учитывается.
В красной части видимой области спектра, около X = 0,76 мкм , расположена узкая, но сильно дифференцированная по интенсивности поглощения полоса, в отдельных линиях которой наблюдается полное поглощение солнечного излучения. Моделирование пропускания атмосферы данной части спектра осуществлялось путем аппроксимации полосы поглощения кубическими сплайнами. При этом полоса поглощения разбивалась на отдельные участки, ограниченные соседними экстремумами функции пропускания солнечной радиации р(Х) . В пределах каждого участка, при п известных точках, использовалось п — 1 сплайнов вида
pi (Хсм ) = k1,i ^ k2,г' Хсм ^ k3,г Хсм ^ k4, г Хсм
Где г - порядковый номер сплайна в пределах очередного участка; Хсм = 7600 — X - смещенная длина волны.
Всего в пределах полосы пропускания было выделено 40 участков с количеством известных (опорных) точек от 4 до 14 [2]. Определение значений коэффициентов сплайнов производилось, от-
дельно для каждого участка, решением систем линейных уравнений, содержавших от 12 до 52 уравнений.
Основные полосы поглощения водяного пара приходятся на длины волн (центр полосы) 0,72, 0,84, 0,94, 1,14, 1,38, 1,87, 2,70,
3,20 мкм. Зависимость поглощения радиации р от содержания водяного пара в столбе воздуха Qп весьма близка к логарифмической; соответственно, при моделировании величин поглощения в качестве аргумента использовалось Q0 = 1п^п ) .
Моделирование зависимости атмосферного поглощения от содержания водяного пара осуществлялось отдельно для каждой из 8
полос поглощения. Аппроксимация значений р^0) производилась с помощью кубических сплайнов; в пределах каждой полосы поглощения использовалось 7 сплайнов.
Поглощение атмосферным озоном учитывалось для 2 основных участков спектра: полос Хюггинса и Шаппюи. Определение значений параметров для каждой зависимости моделируемых значений десятичного коэффициента поглощения (X)
производилось из условия ^ (аг — X ))2 ^ min , где аг -
опорные точки в пределах данной полосы поглощения [2].
Влияние полосы Хюггинса учитывалось в диапазоне X е [0,3 мкм; 0,34 мкм\. С целью повышения точности аппроксимации полоса поглощения разбивалась на 2 участка: X е [0,3 мкм; 0,32 мкм\ и X е [0,32 мкм; 0,34 мкм\; на каждом участке использовалась своя интерполяционная кривая. При количестве параметров в моделируемой зависимости fj (X) менее 5 не
удалось обеспечить достаточную точность аппроксимации (со средней погрешностью не более 2%). Лучшая из полученных зависимостей fj (X), со средней ошибкой аппроксимации 1,37%, и максимальной ошибкой 3,13%, имеет вид:
- при 0,3 мкм < X < 0,32 мкм :
а^) = —7,5304 + 19,7845е 49,9855—1^з89^ — 6,8402^)
- при 0,32 мкм < X < 0,34 мкм :
а(Л) = 7,0413 + 20,8732е15,9137-58,8520Л + 6,83961п(Л) (3)
Полосы поглощения Шаппюи в диапазоне Л є [0,45 мкм; 0,75 мкм] описываются сложной двухмодальной зависимостью а (Л). Вследствие сложного вида зависимости а (Л) не удалось с удовлетворительной точностью аппроксимировать ее единым полиномом на всем диапазоне значений Л. С целью повышения точности аппроксимации было использовано раздельное моделирование значений а(Л) на 4 участках:
- при 0,45 мкм < Л < 0,5767 мкм:
fI (Л) = -0,0257 + 0,0689е р(Л), р(Л) = 183,05Л —1001,9Л2 + 1508ДЛ3 - 2,4774Л4 - 958,9Л5
- при 0,5767 мкм < Л < 0,5883 мкм :
/я (Л) = —878,23 + 4085,1Л — 5502,6Л2 + 51,2Л3 + 2230,95Л4 + 0,0009Л5
- при 0,5883 мкм < Л < 0,6 мкм :
fIII (Л) = -1289,8 + 9935,5Л-23982,8Л2 + 3717Л3 + 19905,8Л4 -20322,7Л5
- при 0,6 мкм < Л < 0,75 мкм :
fV (Л) = -0,0015 + 0,0163е Р(Л), р(Л) = -249,69 + 832,37Л- 305,58Л2 - 1269,4Л3 + 1047,8Л4
Средняя ошибка аппроксимации на всем интервале значений Л составила 3,6%.
Оптическая толщина атмосферы для монохроматического потока при молекулярном рассеянии [3]:
ч ,» 8п\п2 -1)2 6 + 3d ,
тк (Л, г0) =1 ----------------------dz
^ 3Ш4 6 - 7d
где п - показатель преломления воздуха, d = 0,035 - фактор деполяризации молекул воздуха, N - количество молекул воздуха в единице объема.
Учитываем, что (п 2 -1)2 = (п -1)2 (п +1)2 , и с погрешностью
не более 0,01% полагаем (п + 1)2 ~ 4 . Принимаем линейную зависимость показателя преломления от удельного количества молекул и строим данную зависимость по двум точкам:
(4)
> при нормальных условиях: N = 2,68-1025 м 3,
п = 1,000292;
> в вакууме: N = 0 , п = 1.
Тогда получаем п = 1 + 2 92-10 -4----N-------, и соответственно
2,68 -1025
(п2 -1)2 = 4
( N
2,92-10-4---------
= 4,7485-10-58 N 2.
Для идеального газа N = ——Р;
RT
принимая NA = 6,0225-1023 моль_1 (число Авогадро) и
R = 8,3144 —Д------- (универсальная газовая постоянная) оконча-
моль - К
тельно получим:
„ Л 3,0169-10-33 г р(х)^ ^
тЕ (Я, х0) =-------- ------------------------------------------I -dz (5)
«V , о; я ^ т(х) ' '
Поскольку первообразной подинтегрального выражения (5) не существует, необходимо прибегать к численному интегрированию. При численной реализации описываемой модели использовался метод Ромберга, модифицированный экстраполяционной формулой Ричардсона [4].
Атмосферный аэрозоль в большинстве случаев является смесью частиц пыли, льда и мельчайших капелек воды. Основным фактором, определяющим его наличие в атмосфере, в настоящей работе предлагается считать влажность воздуха, повышенные значения которой с большой вероятностью сопровождаются дымкой и туманами.
Влажность воздуха при моделировании аэрозольного рассеяния рассматривается как случайная величина т, распределенная
по нормальному закону с параметрами распределения т (среднее значение) и о т (дисперсия). В качестве условной границы, разделяющей область образования воздушной дымки и область образования тумана, принято значение ттах = 95% . Нижней границей участка возможного образования воздушной дымки считается 444
2
ттп = т - 3от . Данный участок разбивается на отдельные интервалы последовательными значениями т0 = тт;п , т1 = т0 +1, ... тп = ттах . Для каждого интервала ] т1 -1; т1 ] вероятность попадания в него реализации т определяется как рг (тг-1 < т < тг) = F(тг) - F(тг-1)
В качестве дискретной числовой характеристики реализации т используются медианы возможных значений на каждом из участков трасч:
F (тг) - F (трасч) = F (трасч) - F (тг-1) = 0,5р(тг-1 < т < т1) Удельный показатель рассеяния ^ (т, х) при длине волны
Я0 = 0,514 мкм в пределах каждого атмосферного слоя аппроксимируется зависимостью
^ (Щ х) = ах,0 + ах,1т + ах,2т2
где значения коэффициентов ах }- определяются по 3 известным
точкам: ^ (т = 0, х) = 0, а ^ (т = 85%, х) и ^ (т = 90%, х')
принимаются по осредненным данным многолетних климатических наблюдений [5]. Оптическая толщина атмосферы для каждого показателя влажности воздуха определяется как
?гЛ = 1 kЯ>(mг, х] )А/]
]
где ] - порядковый номер атмосферного слоя, А/] - длина пути
солнечного луча в атмосферном слое, км.
Средневзвешенное значение рассеяния монохроматического потока солнечной радиации с длиной волны Я0 = 0,514 мкм при влажности воздуха тг :
1д = 10(1 - е р,
Средневзвешенное значение рассеяния монохроматического потока во всех возможных направлениях в результате образования воздушной дымки определяется как сумма величин рассеяния при всех возможных значениях т :
і=1 і=1 і=1
С учетом симметрии индикатрисы рассеяния относительно горизонтальной плоскости, вероятностное значение оптической толщины атмосферы, обусловленной возможным рассеянием на воздушной дымке монохроматического потока с А0 = 0,514 мкм :
Зависимости (1) и (2) могут служить основой моделирования годового движения Солнца по эклиптике для любого участка Земной поверхности. С помощью (3) - (6) и с учетом результатов моделирования поглощения для молекулярного кислорода и водяного пара осуществляется построение оператора пропускания солнечной радиации атмосферой Земли. В результате становится возможным количественно оценить годовую динамику прихода солнечной радиации на территорию как одного из важнейших факторов, характеризующих ее самовосстановительный потенциал.
В качестве основы для практической реализации модели оценки среднегодовой динамики прихода солнечной радиации предлагается использовать ГИС-технологии, включающие в себя средства фиксации и обработки территориально-привязанной информации, в том числе и визуализацию результатов моделирования. Практическое значение предлагаемой модели заключается в возможности контроля и поддержки принятия научно-обоснованных управленческих решений в сфере планирования природопользования с учетом требований к сохранению экологического равновесия при производстве горных работ.
1. Кранц П. Сферическая тригонометрия / пер. с нем. А. Цейтлина - 2-е изд. - М.: Изд-во ЛКИ, 2007. - 93 с.
2. Хриган А.Х. Физика атмосферы (в 2-х т.). - Л.: Гидрометеоиздат, 1978.
3. Матвеев Л.Т. Физика атмосферы. - 3-е изд. - СПб:Гидрометеоиздат, 2000. - 777 с.
4. Якубович А.Н. Численные методы. - Магадан: Кордис, 2000. - 63 с.
(
(6)
СПИСОК ЛИТЕРА ТУРЫ
5. Научно-прикладной справочник по климату СССР. Серия 3, части 1-6, выпуск 33. - Д.: Гидрометеоиздат, 1990. - 584 с. ШШ
A.N. Yakubovich
THE MODELING OF THE ANNUAL AVERAGE DYNAMIC OF THE SOLAR RADIATION ARRIVAL ON THE TERRITORY DESTURBED AFTER THE MINING WORKS
Is the method of receipt of numeral estimations of average annual dynamics of arrival of sun radiation considered on territory of region with the use of model, taking into account next basic factors: dynamics of annual motion of a Sun on an ecliptic, proper the examined area of territory; composition and structure of ideal atmosphere; state of the real atmosphere, most typical for the examined area of territory; features of landscape of territory are its reflectance and orientation in space.
Key words: solar radiation, modeling of the annual Sun movement, selfreconstruction of the soil layer, geographic information system technologies.
— Коротко об авторе ----------------------------------------------
Якубович А.Н. - кандидат технических наук, директор политехнического института Северо-Восточного государственного университета, [email protected]
