Н.В.Абасов, В.В.Ветрова
О технологии построения планетарных предикторов для задач долгосрочного прогнозирования
Несмотря на то, что отношение к долгосрочному прогнозу природно-климатических факторов (например, приземные температуры воздуха, приток воды в водохранилища ГЗС) остается неоднозначным из-за их невысокой точности и надежности (они носят вероятностный характер), интерес к этой проблеме значительно возрос за последнее время, связанный, прежде всего, с глобальным изменением климата. Начавшееся в конце 70-х годов XX века интенсивное глобальное потепление в Байкальском регионе вызвало изменение ранее наблюдавшихся устойчивых связей. Так, например, достаточно высокая корреляционная связь между полезным притоком в озеро Байкал и боковым притоком в Братское водохранилище практически исчезла, а средняя температура за последние 20 лет увеличилась по сравнению с аналогичным предшествующим периодом на 1,2°С (при среднеквадратическом отклонении за 130 лет 1,0°С), Что же является главной причиной изменения глобального климата?
В настоящее время нет однозначного ответа на этот вопрос. Долгое время доминировала точка зрения о преобладающей роли антропогенного влияния (подписание Киотского протокола). На основе анализа тесной корреляционной связи между средней мощностью солнечной активности за 11-летний солнечный цикл и средней многолетней температурой за этот же период делается вывод о завершении глобального потепления в Байкальском регионе и начале постепенного похолодания [1].
Предложенный ранее системный подход к долгосрочному прогнозированию природно-климатических процессов, влияющих на функционирование и развитие энергетических систем (с практической реализацией в виде системы ГИПСАР [2]), позволяет предсказывать в вероятностной форме гидрологическую и метеорологическую обстановку на несколько лет. В основе этого подхода лежат модели выявления внутренних закономерностей исследуемых процессов за 150-летний период наблюдений. Но, с учетом глобального изменения климата (его неустойчивой фазы), многие найденные ранее закономерности либо перестали срабатывать, либо значительно ослабли, что неизбежно приводит к увеличению ошибок долгосрочных прогнозов и снижению их надежности.
Для повышения точности прогноза, А.П.Резниковым [3] оценивается влияние положения планет-гигантов, сходящихся в одном секторе пространства, а также лунно-приливных сил на геофизические и метеорологические явления Земли. Одним из перспективных способов исследований в долгосрочном прогнозировании является поиск связей различных астрономических показателей с природно-климатическими процессами. В работе Н. Н. Завалишина [4] рассмотрена динамика конфигураций планет, в которой с помощью модели с постоянными угловыми коэффициентами, выявлены циклические компоненты (сидерические и синодические), представленные как интервалы времени, через которые выбранная группа планет занимает одинаковое положение относительно либо звезд, либо Земли. В другой его работе [5] проведена проверка гипотезы влияния характеристик взаимного положения "Земля-Солнце" на режимы среднемесячных температур и осадков на территории Евразии, а также на циркуляцию атмосферы. Для проведения подобных исследований необходимо получить астрономические параметры для различных временных периодов. Выше перечисленные факторы послужили причиной для создания модели движения объектов в Солнечной системе - "5о1аг15".
Уравнения движения планет
Если исключить влияние различных внешних сил на всю Солнечную систему (т.к. они малы по сравнению с ее внутренним гравитационным взаимодействием), то движение каждого объекта можно описать законом всемирного тяготения через систему уравнений (1) в некоторой декартовой системе координат (с началом в центре масс):
и г/г,. • г. (7) . _ , _
Г, (,) = О. £ 1Л[ ЛЧ к, -r.it)- г, (г) , / = 1,
Н Кц
(1)
Ш У- [t )
где в - гравитационная постоянная, ' - масса ¡-ой планеты, '4 ' - ее радиус-вектор, п - количество рассматриваемых объектов (для влияния на земные процессы нас будут интересовать: Солнце, Луна - спутник Земли и все известные планеты Солнечной системы). В скалярной форме систему уравнений (1) можно записать в виде:
>1 Ки ./=1 кп Н Ко
)*> ]*> ^2)
а ее порядок для 11 объектов, наиболее сильно оказывающих влияние на земные процессы, составляет 3*11=33. Учитывая нелинейность дифференциальных уравнений и сложность их интегрирования, в астрономических расчетах часто используют разложение на основное движение по эллиптической орбите тела с координатами (х.у.г) и возмущающее воздействие от других тел (3)
/2 /1 \ * дф х ~~к (1 + т )-т + —
V V дх
.2, v у ЭФ у = -к (1 + /И) — + — г ду
i2 /1 \ z дф -к (1 + т)~г + —
V V &
т.
1 XX j + у у. + ZZ
R,
\ ч
(3)
где Ф- возмущающая функция [6].
Для более точных оценок гравитационного влияния планет требуется учет релятивистского запаздывания сигнала, распространяющегося со скоростью света, как принять в общей теории относительности (например, время прохождения гравитационного сигнала от планеты Плутон до Земли составляет более 4-х часов, а от Солнца до Земли - около 8-ми минут), Учет релятивистского запаздывания существенно усложняет не простую для непосредственного интегрирования систему уравнений (2). Тогда система уравнений (2) будет иметь вид:
уу1 . у(Л — у ({ ~ т \ _
т = у1< ;-2-А, = (4)
(Л Ки
где с - скорость света в вакууме, а гу. - определяет временную задержку гравитационного сигнала от /-ой массы к /-ой.
Для большинства статистических оценок связей не требуется очень высокая точность вычисления траекторий планет, поэтому можно воспользоваться линеаризованными моделями (5), выведенными из уравнений (3) [7-8], в которых отсутствуют производные, а орбитальные параметры достаточно точно описываются в виде явно задаваемых функций от времени. Такой подход позволяет осуществить быстрый доступ к различным параметрам небесных тел, и, в то же время, получаемые значения обладают точностью, достаточной для проведения исследований. Так, для ХХ-ХХ! века, погрешность вычислений траекторий планет (кроме Плутона и Луны, для которых точность значительно ниже, особенно на период более 100 лет) по уравнениям (5) составляет менее 0,2% от наблюдаемых характеристик, зафиксированных в астрономических альманахах. Погрешность нарастает до 2 угловых минут при вычислении данных на 1000 лет. Ниже приведены уравнения движения для 9 планет и Луны в эклиптической системе координат, в центре которой находится Солнце.
Земля: п(/) = 0; /(/)- 0;
eo(t) = 282.9404 + 4.70935 * 10° * t\ a(t) = 1.0;
e(t) = 0.016709-1.151*10""' *?; M(t) = 356.0470 + 0.9856002585 */;
Луна:
n(i) = 125.1228-0.0529538083 */; /(/)= 5.1454;
co(t) = 318.0634 + 0.1643573223 */; ait) = 60.2666; e(t) = 0.054900;
M{t) « 115.3654+ 13.0649929509 */; Меркурий:
Cl(t) = 48.3313+ 3.24587 *10~5 */; i(t)= 7.0047 + 5.00£ -8 */; co(t) = 29.1241 + 1.01444*10°*/; a(t) = 0.387098; e(t) - 0.205635 + 5.59 * 10"10 */; M(t) = 168.6562 + 4.0923344368 *t
Юпитер:
а(г) = 100.4542+2.76854 * 10"5 */; /(/) = 1.3030- 1.557* 1 о"7 */; û>(t)= 273.8777+1.64505 *10"5 */; a(t) = 5.20256; <?(/) = 0.048498+4.469 * 10"9 */; M (t) = 19.8950 + 0.0830853001 * /;
Сатурн:
a(t) = 113.6634+2.38980 *10"5*/; i(t) = 2.4886- 1.081 *10"7*/; ®(t)= 339.3939 +2.97661 *10"5 *t\ a(t) = 9.55475; e(t) = 0.055546 -9.499 * 10~9 *t\ M (t) = 316.967 + 0.0334442282 * /;
Уран:
n(f) = 74.0005+1.3978 *10"5 *t\ i(t) = 0.7733+ 1.9* 10"8*/; û)(t)= 96.6612 + 3.05650*10° *i; a(t) = 19.18171-1.55* 10"5 * t; e(t) = 0.047318 +7.450* 10"9*/; M(t)= 142.5905+ 0.0117258060 * /;
Венера:
П(/) - 76.6799 + 2.46590 * 10-5 */; /(*) = 3.3946 + 2.75 * Ю"8 * Г, со(() = 54.8910 + 1.38374*10"^ а (г) = 0.723330; е(/) = 0.006773 -1.302* 10"9 *(; М(() = 48.0052 + 1.6021302244*/;
Марс:
п(/) = 49.5574 + 2.11081* 10"5 */;
/(/) - 1.8497- 1.78*10"8*/; (о 286.5016 -и 2.92961 *10'5 *(;
а(!) = 1.523688;
е(/) = 0.093405 + 2.516* 10'9 */; М 18.6021 + 0.5240207766 */;
Нептун:
п(/) = 131.7806+3.01730 * Ю""* -5 *1\ /(/) = 1.7700-2.550 *10"7*/; со (1) = 272.8461 - 6.02700 * 10"6 *Г; а(/) = 30.05826 +3.313 * 10"8 * е(с) = 0.008606+2.150*10"9*/; М(х)= 260.2471+ 0.0059951470 *
Плутон:
5(/) = 5= 50.03 + 0.033459652 * С, Я (г) = р = 238.95 + 0.003968789 * /; ср(() = 238.9508 + 0.00400703 * (
- 19.799 * яп(Р) + 19.848 * соз(Р)
+ 0.897 * яп(2 *Р) - 4.956 * со$(2* Р) + 0.610 * 51п(3*Р) + 1.211 * со$(Ъ*Р)
- 0.341 * вт(4 *Р) - 0.190 * соз(4 * Р) + 0.128 * гт(5 * Р) - 0.034 * сов(5 * Р)
- 0.038 * ип(6 * Р) + 0.031 * со8(6*?) + 0.020 * эт(Б-Р) - 0.010 * соэ(5-Р);
^0)= -3.9082
- 5.453 * $т(Р) - 14.975 * сов(Р)
+ 3.527 * ып(2 * Р) + 1.673 * соз(2 * Р)
- 1.051 * Бт(3*Р) + 0.328 * со.ч(3 * Р) + 0.179 * 5т(4 * Р) - 0.292 * С08(4*Р) + 0.019 * эт(5 * Р) + 0.100 * соз(5 * Р)
- 0.031 * зт(6 *Р) - 0.026 * соз(6 * Р)
+ 0.011 * 008(5" - Р);
г (/) = 40.72
+ 6.68 * зт(Р) + 6.90 * соя(Р) - 1.18 * ът(2*Р) - 0.03 * со8(2 * Р) + 0.15 * ип(3*Р) - 0.14 * 008(3 * Р);
Рис. 1. Эклиптическая система коорлинат: а - орбитальные параметры; б - плоскость эклиптики
Здесь П - долгота восходящего узла; /' - наклонение орбиты; со - аргумент перицентра; а - большая полуось; е - эксцентриситет; М - средняя аномалия; -широта, ^/-долгота; 5 и Р - вспомогательные параметры (рис.1); г - радиус вектор; / - время. Для определения точного соответствия астрономического времени его календарным показателям (год, месяц, день, час, минута, секунда, доля секунды) используется формула пересчета (6),
3 6 7 * - 7
7 +
(т + 9 )
1 2
/4 + 275
т ~9~
+ d -730530
(6)
Здесь используются обозначения: у -год; m - месяц; d - день, а также операции выделения целой части из вещественного числа, определяемые квадратными скобками. Шкала времени начинается с 1 января 2000 года в 0:0 часов по всемирному времени (TDT), Часы, минуты и секунды включаются соответствующими долями в показатель дня d.
Программная реализация модели Солнечной системы "SolarlS"
В программную модель включены 9 основных планет, Солнце, а также спутник Земли - Луна. За основу принята декартова гелиоцентрическая эклиптическая система координат, где главной плоскостью является плоскость эклиптики, и Солнце всегда находится в начале координат, а за единичный отрезок расстояния принята одна астрономическая единица. Ось абсцисс в данной системе направлена на точку весеннего равноденствия, ось аппликат перпендикулярна к плоскости эклиптики, а ось ординат образует с ними правую систему (рис, 16). Для вычисления координат планет относительно плоскости эклиптики, используя элементы орбиты, применяются соответствующие геометрические преобразования (7).
X = г * [cos (Q) * cos (и + со) - sin (Q) * sin (о + со) * cos (/)] у = г * [sin (п) * cos (ü + Ú)) + cos (Q) * sin (ü + со) * cos (/)] (7)
z = г * sin (í;+ <£>)* sin (/)
V
arctg(yc! хс), хс = a • (cos E - e), yc - a- (vi - e2 - sin E)
где О, о, со, а, Е, е - орбитальные параметры (рис.1).
Основу модели составляет вычислительное ядро, включающее в себя вышеприведенные уравнения (3). Для получения различных параметров движения небесных тел (или их производных) необходимо задать период времени, за который будут вычисляться данные; период осреднения (день, месяц, квартал, год и др.) для вычислений интегрированных показателей; систему координат; точку начала отсчета (Земля, Солнце, центр масс), Система 5о1аг1$ позволяет выводить во внешние файлы динамические параметры с целью их использования в других программных комплексах, а также строить графики временнбй изменчивости различных производных индексов, Она включает в себя также блок трехмерной визуализации таких процессов, как движение планет относительно Солнца или Земли, Солнца - вокруг центра масс и плоскости вращения лунной орбиты (18,6-летний цикл). Кроме того, в систему включен блок вычисления времени начала и длительности солнечных и лунных затмений.
Следует отметить, что компьютерная реализация модели (2) позволит проигрывать варианты влияния других космических тел (астероидов, комет), находящихся в пределах Солнечной системы, на глобальные климатические изменения на Земле в прошлом и будущем.
Что касается уравнений (4), учитывающих временную задержку распространения гравитационного взаимодействия, то ее компьютерное решение значительно сложнее, Учет запаздывания от Плутона (а возможно и более дальних, например, Седна и других неизвестных пока планет) требует хранения траекторий движения этих планет, по крайней мере, за этот период запаздывания. При достаточно малом шаге интегрирования потребуется очень большая память для хранения их последних фрагментов. При практическом интегрировании можно воспользоваться заданием лишь небольшим набором опорных точек траектории и интерполяцией их, например, кубическим сплайном.
О поиске статистических связей
Реализованная модель Солнечной системы позволяет оперативно формировать различные индексы, с которыми можно статистически сопоставлять природно-климатические процессы на Земле, Например, можно вычислять траекторию движения Солнца относительно центра масс Солнечной системы по формулам:
м (о м (О м (0 (о
x(t)
где х(1), у(!), г(1)- координаты центра масс; М - суммарная масса всех объектов; т - временной интервал. На основе уравнений (8) траектория Солнца имеет достаточно сложную кривую с удалением от центра масс на расстояния более 1,5 млн. км. Для наглядности на рис.2 представлена траектория движения центра масс (в плоскости эклиптики) за период с 1900 по 2011 годы относительно неподвижного Солнца, Несмотря на сложный вид траектории, прямая корреляционная связь за указанный период между среднегодовым боковым притоком воды в Братское водохранилище и среднегодовым отклонением Солнца от центра масс составляет примерно (0.4 ± 0.1), что говорит о не случайности этой связи.
1500
1000 -
500 - 1
-500 -
-1000
-150С
-1500
000
1500 Тыс. км
Рис. 2. Траектория лвижения центра масс относительно Солнца
Согласно уравнениям (5), движение Луны осуществляется не только по примерно эллиптической орбите вокруг Земли, но и сама ее плоскость постоянно вращается, делая полный оборот вокруг Земли за 18,6 лет. Кроме того, система позволяет формировать показатели, отражающие взаимодействия отдельных планет, например, гравитационное воздействие по оси I, вызываемое наклонами плоскостей вращения планет относительно плоскости эклиптики, косинус угла их расположения относительно Земли и многие другие.
Основным методом при оценке связей временных процессов является взаимо-корреляционная функция, учитывающая запаздывание по времени от синтезированного процесса у(1), который можно определить в виде:
}>({) = F(R^(tl...,Rfl(t)),
где Л, (?) - определяют траектории /-ой планеты или ее производных параметров, например, показателей магнитного
поля. Тогда взаимокорреляционную функцию можно определить в виде (для заданного набора дискретных периодов времени):
л-1 , _, _
р{х(1),у(1-т))= £ ^ _
ЩхМ-*«) дЯ'-г-О-о'М) (9!
V /=0 /=0
В формуле (9) x(t),y(t -г) -исследуемые процессы; т- сдвиг во времени процесса у относительно х;
x(t), y(t)- средние значения величин х и у на рассматриваемых отрезках времени. В качестве x(t) можно рассматривать некоторый гидрологический или метеорологический процесс, а y(t-r) - процесс, который, по предположению, оказывает на него влияние (положение центра масс, изменение лунной орбиты); t, т - дискретные значения времени с некоторым заданным периодом (например, год),
Вместо формулы (9) можно воспользоваться некоторой мерой близости, например среднеквадратичным отклонением для нормализованных временных рядов, приведенных к одному масштабу.
Одной из задач является синтез некоторого производного индекса от параметров движения планет с наиболее тесной корреляцией с исследуемыми процессами, который определяется формулой (10):
F(R](t),...,Rn(t),P]\...,P;\ (10)
где (Р;,..., Р;) = arg mm(M(X(t), F(R] (/),..., Rn (0, Pl > • • ■ ,РЯ )) : P, e Pf, i = Ü).
Здесь ju(X(t),F(...)) определяет заданную меру близости между исследуемым процессом на Земле и синтезируемым индексом F(...) через оптимизацию выбранного вектора параметров Pi е Р;, / = 1,.у (Р; - определяет допустимое множество варьирования параметра Pt). Задача (10) не является простой, особенно в случае нелинейной функции F и большой размерности вектора параметров Р,
Для решения задачи (10) можно воспользоваться распределенными и параллельными вычислениями, описанными в работе [9]. Для этого нужно синтезировать распределенную систему на некотором специализированном кластере, где основными вычислителями являются ядро модели SolarlS, вычисляющее параметры планет и их преобразования к выбранной системе координат, а также процедуры вычисления функций F, которые могут иметь не только разный набор параметров, но и, возможно, саму структуру. Некоторые опытные численные эксперименты показали перспективность такого направления, когда синтезируемый индекс уже на линейной функции F от параметров
PI е Q. ,/ -l,s (s ^15) позволяет получать корреляционные связи более 0,7 не только на выбранном временном
интервале (около 100 периодов и более), но и, самое главное, на различных верификационных выборках. К сожалению, подавляющее большинство методов моделирования сложных природно-климатических процессов, какими являются приземные температуры воздуха или гидрологические (полезный среднегодовой или квартальный приток в озеро Байкал; боковой в Братское водохранилище и др.) дают неприемлемую ошибку именно на некоторых верификационных выборках, когда принципиально изменяются глобальные климатические процессы. В данном случае решение задачи (10) можно рассматривать в качестве специализированного метода математического моделирования. Вопрос же реальной физической интерпретации найденных результатов остается открытым, Устойчиво-найденным связям (хорошей корреляции), проявляющимся не случайно и на различных независимых подвыборках, безусловно, можно найти достаточно обоснованную физическую интерпретацию.
Рассмотренный выше подход к выявлению предикторов на основе учета динамики различных планетарных показателей и реализованного ядра моделирования линеаризованных моделей движения планет Солнечной системы позволяет выдвигать различные гипотезы в виде некоторых задаваемых функций от координат планет и их производных показателей с определением оптимального вектора параметров, входящих в состав этих функций, Система SolarlS позволяет наглядно визуализировать произвольно заданные конфигурации в выбранном ракурсе с отображением необходимых производных параметров, а также точные временные промежутки лунных и солнечных затмений. Кроме того, система позволяет формировать в выходной файл временные ряды основных или производных параметров для их дальнейшей обработки в других программных комплексах.
В данной работе описан подход с компьютерной реализацией компонентов, основанный на гипотезе наличия существенных связей между земными и планетарными взаимодействиями, Модели движения планет, рассмотренные в работе, могут быть существенно изменены в зависимости от результатов многочисленных компьютерных экспериментов.
Библиографический список
1, Башкирцев B.C., Машнич Г.П. Переменность Солнца и климата Земли II Солнечно-земная физика,-2004,-Вып.6,-С. 135-137.
2, Абасов Н.В., Бережных Т.В., Резников А.П. Долгосрочный прогноз природообусловленных факторов энергетики в информационно-прогностической системе ГИПСАР II Известия РАН, Энергетика,-2000,-№6,-С, 22-30.
3, Резников А.П. Предсказание естественных процессов обучающейся системой (физические, информационные, методологические аспекты). - Новосибирск: Наука, 1982.-287 с.
4, Завалишин H.H. Циклические компоненты в динамике планетных конфигураций II Сибирское отделение ВАСХНИЛ 1983.-60 с.
5, Виноградова Г, М., Завалишин Н. Н. Один из факторов формирования аномалий погоды II Сибирское отделение ВАСХНИД 1988. - 167 с.
6. Абалакин В,К. Основы эфемеридной астрономии, - М,; Наука, 1979, - 448 с,
7, Schlyter Раи!, "How to compute planetary positions", http://docs.avramishin.com/other/astronomy/ ppcomp.html
8. Jean Meeus, Astronomical algorithms - Virginia: Willman Bell Inc., 1991, - 220pp,
9, Абасов H.В, Подход к организации пакетно-ориентированных распределенных вычислений II Тр. X Байкальской Всерос. конф. «Информационные и математические технологии в науке, технике и образовании». - Иркутск: ИСЭМ СО РАН, 2005. - 4,1. - С.109.