Г. Е. КРАСНИКОВ, О. В. НАГОРНОВ
Национальный исследовательский ядерный университет «МИФИ»
МОДЕЛИРОВАНИЕ ТЕЧЕНИЯ ЛЕДНИКОВ
МЕТОДОМ КОНЕЧНЫХ ЭЛЕМЕНТОВ
Разработана численная двухмерная модель расчета скоростей течения льда в горном леднике, основанная на методе конечных элементов. Модель реализована в системе Comsol Multiphysics 3.2 (FEMLAB). Проведено тестирование модели на точность решения задач подобного класса. Модель используется для расчета поля скоростей в ледниках Григорьева и Академии Наук при различных сценариях изменения климата.
Изучение динамики льда в горных ледниках необходимо для составления прогнозов эволюции ледников в условиях возможного потепления климата. Работы по исследованию течения льда в горных ледниках ведутся уже более 20 лет. Pattyn, Hindmarsh, MacAeyal разработали различные модели течения льда, основанные на конечно-разностном методе [1–3]. Были созданы модели течения ледника Григорьева (Тянь-Шань) [4, 5]. В нашей работе описывается модель расчета скоростей течения в различных горных ледниках, основанная на методе конечных элементов. Эта модель может быть основой для составления детальных прогнозов о дальнейших изменениях ледниковых систем.
Течение горного ледника может быть описано с помощью уравнений движения сплошной среды:
u – скорость течения льда; σ – тензор напряжений; σ’ – девиатор напряжений; 
– тензор скоростей деформаций; η – вязкость; n – показатель Глена (обычно n = 3) [6]. Величина A зависит от температуры льда. Область ледника ограничена поверхностями hb снизу и hs сверху.
Граничные условия выбираются из следующих соображений. На ложе ледника выполняется условие прилипания – скорости течения здесь равны нулю. На свободной границе давление равняется атмосферному:
где u и v – горизонтальная и вертикальная компоненты скорости соответственно; p – давление. После преобразований система уравнений сводится к системе трех дифференциальных уравнений с тремя неизвестными u, v, p:
Поскольку ледник представляет собой замкнутую область со сложным профилем, то удобно решать систему уравнений методом конечных элементов. Система конечно-элементных расчетов Comsol Multiphysics (FEMLAB), основанная на этом методе, довольно удобна для решения краевых задач со сложными областями.
Для тестирования системы FEMLAB выбраны данные международного эксперимента «Benchmark experiments for higher-order and full Stokes ice sheet models» [7]. В данном эксперименте производился расчет скоростей течения в тестовом леднике с заданным профилем для различных численных моделей. Расчеты были произведены более, чем для 20 моделей, разработанных разными авторами.
Сравним скорости течения, полученные в системе FEMLAB, со скоростями течения, полученными в других моделях согласно пункту B эксперимента (двумерный ледник). Профиль ледника обладает следующими особенностями. Свободная граница ледника – это прямая, наклоненная под некоторым углом α к горизонту, а основание ледника параллельно свободной поверхности и имеет синусоидальный вид с некоторым периодом. Область ледника задается следующими выражениями:
Расчеты проводятся для значений L = 5, 10, 20, 40, 80, 160 км; H = 1000 м, A = 500 м. Для моделирования бесконечного ледника была задана область с размером в несколько раз превышающим характерный размер профиля ледника (период функции zb(x, y)). Интересующая нас область будет находиться в центре такого ледника. Будем рассматривать профили скоростей на протяжении одного периода функции, задающей форму основания. В эксперименте сравниваются максимальные скорости на поверхности ледника.
Сравнение результатов расчетов модели Femlab и результатов расчетов моделей Pattyn и других представлено на рис. 2. Относительная разница результатов, полученных в Femlab, и результатов, полученных другими моделями, довольно низкая, не превышает 5 %. При этом она понижается с ростом периода функции, задающей профиль основания L. Если при L = 5 км расхождение в полученных скоростях у нашей модели с остальными достигает 5 %, то уже при L = 10 км максимальное отклонение составляет 4 %. При дальнейших значениях L расхождение результатов нашей модели с результатами большинства остальных экспериментов незначительно (1–2 %). Исходя из этого, пакет Femlab можно считать пригодным для моделирования течения льда в ледниках. Однако стоит учитывать, что при увеличении неровности основания ледника точность может уменьшаться.
С помощью модели были рассчитаны скорости течения льда в ледниках Академии Наук и Григорьева. Ледник Академии Наук располагается на острове Комсомолец архипелага Северная Земля в Северном Ледовитом океане (координаты 80–81º с. ш., 92–96º в. д.) и является крупнейшим ледником в российской Арктике. По этой причине прогнозы эволюции ледника при возможном глобальном потеплении представляют значительный интерес. Ледник занимает преимущественно южную и центральную части острова Комсомолец, составляя около 60 % его территории. Общая площадь ледника на 2002 г. составляла 5575 км2. Максимальная возвышенность ледника около 800 м, максимальная толщина льда 819 м [8].
На рис. 3 представлены результаты модели расчета скоростей в системе FEMLAB. Наибольшие значения скоростей течения льда находятся в местах наибольшего угла наклона поверхности. При этом высокие значения скоростей сохраняются вдоль всей глубины ледника, сходя на нуль только у основания. Максимальная скорость льда наблюдается на склоне справа от вершины и составляет ~ 3 м/год. Вертикальная скорость течения на порядок меньше горизонтальной, ее максимальное значение составляет 0,6 м/год.
1) |
2) |
3) |
4) |
5) |
6) |
Рис. 2. Скорости течения льда при значениях L: 1) 5 км,км,км,км,км,км, где по результатам расчетов в Femlab (а); по результатам моделей Pattyn и др. [6] (б)
а) |
б) |
Рис. 3. Скорости течения льда в леднике Академии Наук: а) горизонтальная; б) вертикальная
Ледник Григорьева располагается на южном склоне системы Терскей Ала-Тау, располагающейся в горном массиве Тянь-Шань. Характерные размеры этого ледника: длина – 3,6–3,7 км; ширина – 3 км. Расположен ледник на высоте 4100–4700 м. Толщина ледникового покрова меняется на его протяжении в среднем плавно возрастая от нижнего края и достигая максимального значения в области 
км (
110 м).
а) |
б) |
Рис. 4. Скорости течения льда в леднике Григорьева: а) горизонтальная; б) вертикальная
Горизонтальные скорости находятся в диапазоне от 0 до 5,15 м/год (рис. 4). Наибольшие значения скоростей наблюдаются на втором уклоне. На первом уклоне значения достигают 3,6–3,7 м/год вблизи поверхности и сохраняются на этом уровне до глубин ~ 100 м. Значения вертикальных скоростей лежат в диапазоне от 0 до 1,97 м/год.
Работа выполнена при поддержке ФЦП «Научные и научно-исследовательские кадры инновационной России» на 2009–2013 годы (ГК П943) и АВЦП «Развитие научного потенциала высшей школы», проект 2.1.1/6827.
СПИСОК ЛИТЕРАТУРЫ
1. Hindmarsh, R. C. A. // J. Geophys. Res. 2004. V. 109. F01012, doi:10.1029/2003JF000065.
2. MacAyeal, D. R. // J. Geophys. Res. 1989. V. 94(B4). С. 4071.
3. Pattyn, F. // Ann. Glaciol. 2000. V. 31. С. 211.
4. Konovalov Y. V., Nagornov O. V. // The Cryosphere Discuss. 2009. V. 3. С. 77.
5. , // Материалы гляциологических исследований. 2005. Т. 98. С. 183.
6. Б. Физика ледников: перевод с английского. М.: Мир, 1984.
7. Pattyn F. et al // The Cryosphere Discuss. 2008. V.2. P. 111.
8. Атлас снежно-ледовых ресурсов мира / Под ред. . М.: РАН, Институт географии, 1997.
