Введение

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

Моделирование необходимо для изучения сложных природных процессов, оно заключается в выделении основных особенностей и закономерностей, а также последующей идеализации изучаемого явления. Такое упрощенное представление об исследуемом объекте позволяет получить его описание с дальнейшим прогнозом поведения – например, деформации материала и разрушения. Адекватность моделей должна подкрепляться экспериментальными данными, а важным средством эффективного нахождения численного решения задач являются современные вычислительные методы и их реализация в программном обеспечении. Идеализация – сложный этап процесса изучения объекта, особенно материала сложной структуры – например, слоистой [1–3] или ячеистой [4]. Лед и композиционные материалы на его основе применяют в качестве конструкционного материала для построения инфраструктуры в арктических зонах и в регионах с холодным климатом: дорог, взлетно-посадочных полос, площадок разгрузки и погрузки и др. [5–7]. При использовании льда желательно, чтобы его поведение было прогнозируемым. Хотя большое внимание уделяется его экспериментальному исследованию [8–11], для дальнейшего развития работ необходимо проводить моделирование льда, дающее возможность объяснять и предсказывать поведение конструкций в различных условиях эксплуатации.

Получение численных моделей льда – нетривиальная задача с учетом его варьируемых прочностных свойств, хрупкости, сильной чувствительности свойств льда к заморозке и условиям эксплуатации [11]. Лед – особый материал и на него нельзя напрямую распространить методы моделирования для остальных классов материалов. В композиционных материалах на основе льда изучение свойств усложняется вследствие взаимодействия разных по структуре фаз. Применение традиционных методов часто затруднено и требует предварительных исследований.

Количество работ по численному моделированию механики льда незначительно, но число их растет в последние годы, что подтверждает интерес к данной теме. В данном обзоре рассмотрены основные работы по моделированию льда, поэтому в настоящее время критический анализ и сопоставление разных подходов, представленных в научно-технической литературе, актуально.

Метод конечных элементов

Метод конечных элементов (МКЭ) широко распространен для решения задач прикладной механики конструкций. Лед можно считать упругим только при очень малых деформациях, поэтому в этом случае допустимо использовать простейший вариант – модель линейно-упругого тела [12].

Моделирование более сложного процесса нелинейного деформирования и хрупкого разрушения требует использования расширенных моделей [13, 14], к которым можно отнести нелинейные модели МКЭ изотропных материалов Джонсона–Холмквиста, Леметра и др.

В работе [15] рассмотрены и проанализированы с помощью МКЭ различные модели материалов и проведено моделирование в программном комплексе LS-Dyna. Экспериментальные данные взяты из открытых научно-технических литературных источников и собственных данных автора по пресноводному поликристаллическому льду. Кроме применения различных критериев разрушения, автором для получения сходящегося решения применяется метод удаления элементов из модели при достижении предельного шага интегрирования.

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

В результате моделирования динамики контактного взаимодействия твердого тела с моделью льда получено, что наиболее согласующиеся с экспериментальными данными результаты демонстрируют модели упругопластического материала с кинематическим упрочнением и вязкопластический материал. В работе [15] для моделирования взята упругопластическая модель материала Прандтля–Рейса с изотропным упрочнением и критерием текучести Мизеса:

                                                              (1)

где  – компоненты приращения пластических деформаций; λs_ij – компоненты девиатора напряжения; σ_i – интенсивность напряжений; σ_т – предел текучести.

В зоне ползучести эта модель может применяться без критерия разрушения. Для зоны разрушения более применима гидродинамическая модель поведения разрушенных элементов.

Рассмотренные модели являются наиболее универсальными, но требуют разделения их применения по исследуемым процессам деформации льда.

Применение изотропного упругопластического материала с изотропным упрочнением приводится в работах [16, 17] при моделировании ударного взаимодействия упругого цилиндра со льдом. Расчеты проводятся в программном комплексе LS-Dyna, а также приводятся результаты для авторской программы для двухмерной постановки задачи. Указывается и основной недостаток проведения подобного рода моделирования – точность определения свойств льда из экспериментальных данных. Так потребовалась корректировка методики обработки экспериментальных данных из-за завышенных значений кривой зависимости предела текучести при сжатии от скорости деформации. Получена удовлетворительная корреляция экспериментальных и расчетных данных в широком диапазоне скоростей удара (118–550 м/с). В настоящее время в работах [16, 17] указывается недостаток данных по деформированию разрушаемого льда, что необходимо для определения и расширения границ применимости физико-механических моделей динамического поведения материала.

В работах [18–20] подробно рассмотрено моделирование взрывного воздействия на лед. Коллективом авторов разработано собственное программное обеспечение для расчета такого воздействия на лед в двухмерной постановке задачи и используется модель упругопластического течения Прандтля–Рейса с критерием текучести Мизеса. Получена хорошая сходимость результатов расчетов с большим количеством собственных экспериментальных данных, накопленных авторами, а также взятых ими из научно-технических источников. Специфика работ авторов – взрывное воздействие на лед для расчистки русла рек.

Данные работы посвящены скоростному нагружению льда. Для изучения поведения ледяного покрова, ледяных переправ, льдин, припойного льда необходимо исследование механики квазистатического нагружения менее интенсивного характера. В отечественных работах [21–23] представлен опыт моделирования ледяной балки при четырехточечном изгибе. Лед описывается упругопластическим материалом с нелинейным критерием текучести Виллама–Варнке, разработанного для моделирования бетона. Выбор параметров для определения поверхности текучести по данному критерию представлен в работе [24]. Исследовано армирование, выполненное со стальным каркасом. Нагружение проведено только до начала разрушения, т. е. до активного растрескивания ледяной матрицы и наступления условной «пластичности» нагружаемой балки. Проведено моделирование и продемонстрирована хорошая корреляция с результатами эксперимента. На основе полученных моделей делается попытка вариативного исследования других материалов для армирования. Прочие типы нагружения не принимались во внимание. Данное исследование отличается тем, что рассматривается армированный лед, применение которого позволяет существенно повысить прочностные и деформационные свойства материала.

В зарубежной работе [25] также исследуется разрушение ледяной балки при четырехточечном изгибе. Проводится оптимизация модели материала с раздельным поведением при растяжении и сжатии, с критерием текучести по Мизесу. В качестве критерия разрушения выбран критерий для разрушения композиционных материалов Цая–Ву. При нагружении четырехточечным изгибом получена хорошая корреляция с экспериментальными данными, но применение модели потребовало проведения оптимизационного поиска параметров модели материала.

В работах [26, 27] представлено необычное применение метода когезии, который заключается в моделировании между конечными элементами промежуточных виртуальных интерфейсных элементов. При достижении предельных условий интерфейсные элементы разъединяются. Основная особенность метода заключается именно в заранее определенной зоне с интерфейсными элементами, т. е. возможное направление роста трещин определяется заранее. Традиционная область применения метода – механика композиционных материалов, граница сред, моделирование клеевых соединений, учет адгезии арматуры в железобетоне. Однако авторы адаптировали данный метод для хрупкого разрушения льда, для чего интерфейсные элементы были добавлены между всеми конечными элементами модели. Представлены результаты трехмерного моделирования для следующих типов нагружения: растяжение, сдвиг, растяжение при раскалывании (brazilian test), падение с высоты на упругую и жесткую поверхность. В целом получена удовлетворительная корреляция экспериментальных и расчетных данных, особенно в процессах с преобладанием одной трещины при растрескивании объемного тела. Но существует проблема контакта интерфейсных элементов и внешних тел – например, при моделировании падения образца на жесткую поверхность получены на порядок большие контактные силы. С увеличением количества элементов увеличилось количество степеней свободы, соответственно, значительно возросли требования к вычислительным ресурсам. Авторами работ [26, 27] получена зависимость результата от параметров разбиения объемного тела. Параметры контакта должны быть подобраны каждый раз в зависимости от сетки модели. Данные исследования заявляются первой попыткой в применении данного метода в механике хрупкого разрушения. Постановка задачи напоминает переход уже к бессеточным методам.

Применение МКЭ в моделировании хрупкого разрушения льда обобщено в работе [28]. В отличие от ранее указанных работ исследуется и пресноводный, и морской лед – менее прочный и более пластичный. Исследование проводилось на разномасштабных моделях, при этом демонстрируется неуниверсальность существующих моделей материалов.

Автор работы [28] использует разные типы моделей материалов МКЭ для хрупкого характера разрушения, в основном разработанных для таких материалов, как бетон и гранит (применен программный комплекс LS-Dyna). На основе анализа ряда работ приводится выбор подходящих для хрупкого разрушения моделей материалов. Проведено моделирование с применением следующих моделей материалов [29]:

– разрушаемый пеноматериал (MAT_063, MAT_CRUSHABLE_FOAM в LS-Dyna);

– модель материала с повреждаемостью Леметра (MAT_153, MAT_DAMAGE_3);

– модель материала с хрупким разрушением Говинджи–Кея–Симо (MAT_096, MAT_BRITTLE_DAMAGE);

– модель материала Джонсона–Холмквиста (MAT_111, AT_JOHNSON_ HOLMQUIST_CONCRETE).

Применение моделей пеноматериалов применительно ко льду рассматривалось также в работах [30–32]. Данная модель изотропная, одномерная. При растяжении происходит разрушение. Модуль упругости постоянен, напряжения при нагрузке для каждого шага пересчитываются по линейному закону, разгрузка пеноматериала происходит по эластичному закону.

Расчет напряжений после шага нагружения (n) с промежуточной оценкой (n+1/2) производится согласно формуле

                                                  (2)

где σ_ij – тензор напряжений; E – модуль упругости; ^·ε_ij– тензор скоростей деформации; Δt – расчетный шаг времени.

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

Модель материала с хрупким разрушением Говинджи–Кея–Симо и модель бетона Джонсона–Холмквиста применялись ранее только для бетона и гранита, для льда они применены и адаптированы в работе [28].

Модель с хрупким разрушением Говинджи–Кея–Симо предполагает, что разрушение (накопление повреждений) в материале начинается, когда главные напряжения превышают некоторое пороговое значение (подробное описание дано в работе [29]). Повреждаемость учитывается путем уменьшения модуля упругости материала по экспоненциальному закону. Преимуществом данной модели является минимальный набор параметров, описывающих материал: массовые, жесткостные, вязкостные параметры, а также предельные напряжения.

Модель Джонсона–Холмквиста разработана для моделирования бетона, подверженного большим деформациям, высоким скоростям деформации и высокому давлению. Уравнения определяются в нормированных величинах.

Текущая прочность материала на расчетном шаге времени (n+1) с учетом повреждаемости равна (здесь и далее – знак * обозначает нормированный параметр)

                                            (3)

где  – прочность неповрежденного материала;  – прочность для полностью разрушенного материала;  – скорость пластической деформации;  Ρ^*– гидростатическое давление; T, B, C, m – константы материала.

Модель аккумулирует повреждения как от эквивалентной пластической деформации, так и от пластической объемной деформации:

                                                       (4)

где D₁, D₂ – константы, описывающие характер разрушения конкретного материала.

Более подробное описание этой модели можно найти в работе [33].

Для полного определения большого количества констант для конкретного материала требуются целый ряд предварительных испытаний: на неограниченное сжатие, на вытягивание, на трехосное сжатие, а также серии испытаний на сжатие при различных скоростях деформации, что затрудняет использование модели материала.

Численные эксперименты по исследованию льда совместно с водой на основе модели повреждаемости Леметра можно найти в работе [34]. Совместные усилия авторов работ [25–27] и [34] представлены в международном проекте Евросоюза SAFEWIN для анализа судоходного движения судов во льдах Балтики. Как видно в публикации по проекту [35], авторами была выбрана и применялась именно модель Леметра с накоплением повреждаемости для анализа давления ледяного покрова на корпус судна. Проведенные испытания в специальном бассейне демонстрируют хорошую корреляцию моделирования с экспериментальными данными при данном типе разрушения.

Модель Леметра строится на концепции зависимости повреждаемости D от величины пластической деформации. Влияние поврежденности среды на прочность выражается уравнением

                                                           (5)

Изменение меры поврежденности определяется соотношением

                                                            (6)

где  – эквивалентная скорость пластической деформации; S и t – константы, характеризующие материал; Y – плотность энергии деформации.

Полное описание модели приведено в работе [20].

В работе [28] рассмотрены две расчетные схемы – ударное нагружение бойком и контакт объемной пластины с твердым телом – аналог контакта «льдина–судно» или «льдина–опора моста» [36]. В отличие от работ [30–32] исследование в работе [28] показало, что модель пеноматериала с разрушением дает неудовлетворительный результат с завышенными на порядок напряжениями прежде всего для ударного нагружения, что можно объяснить более подробной моделью в работах [30–32], в которой использовалось разбивание ледового объекта на отдельные блоки с корректировкой характеристик для каждого блока.

Результаты для моделей, учитывающих повреждаемость, имеют желаемый физичный характер разрушения: кривая зависимости силы от перемещений для всех моделей материалов имеет пилообразный вид, адекватно описывает скачкообразно меняющийся характер разрушения льда. Порядок результирующих величин (реакций и перемещений) близок к экспериментальным данным. Однако следует отметить неустойчивость решения к изменению входных параметров моделей (варьирование приводит к несходящемуся решению) – определение параметров требует существенных сопутствующих материальных затрат на валидацию моделей, особенно для модели Джонсона–Холмквиста. Хрупкое разрушение с помощью приведенных моделей получить проще, чем пластичное деформирование льда. Моделирование потребовало значительных вычислительных ресурсов – проблемы со сходимостью решения в зоне разрушения приводили к вынужденному сгущению сетки. Отмечается недостаток экспериментальных данных с учетом вариативности свойств льда. В работе [28] указывается на возможное значительное влияние анизотропии морского льда на характер его деформирования, что пока подробно не исследовано. Отмечаются непростые условия применения МКЭ, приспособленного для континуальных сред при решении задач хрупкого разрушения. Отдельно предполагается, что предпочтительнее использовать дискретные бессеточные методы, такие как, например, гидродинамика сглаженных частиц.

Гидродинамика сглаженных частиц

Гидродинамика сглаженных частиц (Smoothed Particle Hydrodynamics – SPH) представляет собой бессеточный метод, в котором среда заменяется набором частиц, т. е. параметры среды распределяются между частицами [37]. Первоначально разработанный для астрофизики метод применяется во многих областях и прежде всего в вычислительной гидродинамике. Активно ведутся работы по адаптации метода для механики твердых тел [38–40].

На рис, 1, а приведена схема дискретизации области Ωсо смещением частицы i на вектор u. Для i-й частицы выделена область распространения функции W(r, h), которая представлена на рис. 1, б. Положение в пространстве меняется с r₀до r_d (радиус-векторы i-й частицы), расстояние h называют сглаживающей длиной.

Рис. 1. Гидродинамика сглаженных частиц:

а – схематичное представление дискретизации области Ω и ее деформации; б – схематичное представление нормированной сглаживающей функции ядра W(r, h)

Уравнение плотности тела в формулировке гидродинамики сглаженных частиц записывается через сглаживание параметров соседних частиц:

                                                     (7)

Уравнение состояния упругого тела, с учетом уравнения неразрывности среды, имеет вид [30]:

                                                      (8)

где  ν^a– а-я компонента скорости;  σ^ab– тензор напряжений;  g^a– объемные силы.

Уравнение движения для i-й частицы примет вид:

                                      (9)

Выражение (9) дополнено искусственной вязкостью П_ij для стабилизации решения (δ^ab – символ Кронекера).

Если в сеточных методах применительно к моделированию твердых тел существенен вопрос – как разделить элементы сетки при разрушении, то в бессеточных методах – как удержать частицы вместе.

Метод позволяет моделировать большие деформации и разрушение среды путем разделения скопления частиц. Это объясняет возросший интерес к данному методу в области моделирования хрупкого льда. В работе [41] представлено решение задачи разрушения образцов льда в двухмерной постановке задачи при четырехточечном изгибе и одноосном сжатии. Для пластической деформации используется критерий Друкера–Прагера, для валидации модели упругого тела применяется сравнительный расчет упругих колебаний консольно закрепленной балки.

Качественно характер разрушения точно совпадает с экспериментальными данными. Количественно же потребовалась корректировка модели, авторами работы [41] вводится интерполяция полей деформации по первой производной во времени для устранения неустойчивых колебательных процессов в частицах – данная особенность метода может вносить нефизичные эффекты в деформацию твердых тел.

В работе [42] метод получил развитие: проведено двухмерное моделирование разрушения ледяной балки при трехточечном изгибе, прогрессирующее разрушение балки об препятствие и разрушение равномерного ледяного покрова при прохождении судна с разными скоростями (но без моделирования воды). Эксперимент для сопоставления проводили в бассейне на масштабной модели. Получено хорошее соотношение экспериментальных и расчетных данных. В вычислительном эксперименте получена более консервативная оценка нагрузок, действующая на корпус модели корабля. Метод гидродинамики сглаженных частиц с точки зрения авторов работы [42] показывает удовлетворительную точность и устойчивость решения к изменению входных параметров. Предполагается дальнейшее развитие моделей, дополнение их водной средой и расчет системы «лед–судно–вода».

Метод гидродинамики сглаженных частиц больше подходит для моделирования водной среды. В ряде работ [43–46] лед рассматривается как упругий объект, движущийся на поверхности воды. Данный подход больше подходит для крупномасштабного исследования, анализа переноса больших масс льда, оценки давления жидкости на ледяной покров, без рассмотрения локальной деформации и разрушения.

Метод дискретных элементов

Другим бессеточным методом, применяемым для анализа хрупкого разрушения льда, является метод дискретных элементов (МДЭ). В нем физическое состояние системы складывается из физических состояний большого числа дискретных элементов. Набор элементов образует исследуемую среду, но каждый из них является независимым и взаимодействует с соседними элементами. Общее (макроскопическое) поведение системы является результатом взаимодействия отдельных частиц. В 1956 г. был разработан метод для моделирования частиц в молекулярной динамике [47], в 1971 г. метод был адаптирован для исследования динамики горных пород и назван методом дискретных элементов [48]. Теоретическое обоснование метода дано в работе [49] и показано, что его можно считать обобщением МКЭ. В настоящее время этот метод активно применяется для моделирования сыпучих сред [50].

Метод дискретных элементов, как бессеточный метод, хорошо адаптирован под описание быстро протекающих процессов, связанных с переносом вещества. Такие процессы легче описать в виде, когда каждая частица имеет отдельное физическое состояние, а эволюция всей среды обусловлена локальными физическими взаимодействиями между частицами. Применение метода дискретных элементов для моделирования непрерывных сред потребовало выработки специальных подходов моделирования, нетипичных для исходного метода. Метод дискретных элементов был дополнен моделями связей в виде балок, соединяющих частицы. В зарубежной научно-технической литературе выделяют модели для объемных тел, подчеркивая их обособленность от исходного метода, – модель из связанных частиц (рис. 2). Связь между частицами представляет собой виртуальную балку, которая ограничивает тангенциальное t и нормальное n смещение частиц относительно друг друга (рис. 3).

                                                        (10)

где F_n, F_t – проекции результирующей силы в направлениях n и t; M_n – нормальный момент; M_t – сдвиговый момент;A=πr²_b  – площадь сечения; r_b – радиус жесткой связи;  – момент инерции сечения; δt – шаг по времени; S_n и S_t – жесткости в нормальном и поперечном направлении; v_n, v_t – скорости в направлениях n и t; ω_n, ω_t – угловые скорости вокруг соответствующих осей n и t.

Рис. 2. Модель связанных частиц тела:

а – дискретные элементы; б – визуализация балочных связей между элементами

Рис. 3. Схема двух частиц со связью в виде виртуальной балки (F_b – результирующая сила)

Разрыв связи между частицами происходит при условиях

                                                     (11)

После разрыва связи частицы взаимодействуют как независимые твердые тела.

Эта модель была представлена в работе [51], разработана она для моделирования бетонных конструкций и скальных пород. На ее основе в настоящее время делаются попытки построить модели для льда. Одной из первых работ по оценке применения МДЭ для моделирования льда можно указать статью [52], в которой проведено моделирование одноосного сжатия и трехточечного изгиба простейших образцов морского льда прямоугольной формы в трехмерной постановке задачи. Полученные откалиброванные мезопараметры материала распространены на модели массивных льдин, воздействующих на опоры шельфовых установок, проведена оценка их воздействия [53]. Материал рассматривается как изотропный. Зависимость температуры, солености и анизотропии не исследуется, но предполагается в будущих исследованиях.

Метод дискретных элементов подходит для моделирования льдин с учетом того, что он адаптирован к независимому движению большого числа элементов. Причем геометрическая форма данных элементов не обязательно должна представлять собой сферу, что является простейшим случаем. В целом ряде работ [54–57] рассматривается подобное моделирование ледяного покрова, его раскалывание и воздействие на судно или опоры. Возможно также использование нескольких методов в единой модели: МДЭ – для хрупкого льда, МКЭ – для конструкций, на которые воздействует лед.

Как видно из данных приведенных работ, использование МДЭ перспективно в моделировании механики льда, но пока разработка этого направления только начинается.

Обсуждение и заключения

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

Наблюдается расширение данного направления: применение новых, ранее не адаптированных методов; корректировка уже созданных моделей. Однако модели пока ограничены, в них не анализируется анизотропия льда, изменяемость свойств во времени и пространстве, соленость и температурные условия. Многие модели разномасштабны и подходят только для определенных постановок задач.

На основе проанализированных работ можно сделать некоторые обобщающие выводы.

Стандартизированная методика моделирования льда и композиционных материалов на его основе в настоящее время не разработана.

При моделировании льда зачастую адаптации подвергаются модели и методики моделирования, разработанные для бетона и горных пород. По характеру деформирования и разрушения эти материалы ближе всего оказываются ко льду.

Добиться с помощью МКЭ качественно верного решения возможно с применением моделей, основанных на повреждаемости тела, но даже при качественно верном решении получается количественное отклонение в оценке напряжений. Зачастую определение полного набора параметров таких моделей затруднительно – например, рассмотренная модель Джонсона–Холмквиста. Более простые модели с учетом повреждаемости демонстрируют более стабильный результат. Модель с повреждаемостью Леметра в нескольких исследованиях показала адекватный и удовлетворительный результат.

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

Специфика методов гидродинамики сглаженных частиц и дискретных элементов с движением отдельных частиц может привести к неустойчивому поведению частиц. Методы нуждаются в стабилизации и настройке для моделирования сплошных тел.

Метод дискретных элементов выглядит перспективным для применения в моделировании деформации и разрушения льда. Решаются задачи и в двухмерной, и в трехмерной постановке задачи. Имеется возможность построения мультифизичных моделей сложных систем «вода–лед конструкция» (объединение вычислительной гидродинамики, МДЭ и МКЭ). При этом в методе имеется возможность применения неоднородных и несферических частиц. В дальнейшем также могут быть усовершенствованы модели связей частиц для более точного моделирования сложных процессов разрушения.

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

Исследование выполнено при финансовой поддержке Российского научного фонда (грант №18-13-00392).