CC BY
128. Новацкий В. Теория упругости. М.: Мир, 1975.
9. Атоян А.А., Саркисян С.О. Изучение свободных колебаний микрополярных упругих тонких пластин // Докл. НАН Армении. 2004. 104, № 4. 287-294.
Поступила в редакцию 19.12.2007
УДК 539.3
ПРИМЕНЕНИЕ МНОГОСЕТОЧНОГО МЕТОДА ДЛЯ РЕШЕНИЯ ЗАДАЧИ О ШИНЕ
К. А. Лопухин1, C.B. Шешенин2
В работе исследуется эффективность одного варианта геометрического многосеточного метода для решения задачи о качении шины и производится подбор его оптимальных компонент. Построенный метод сходится эффективно для достаточно больших с практической точки зрения систем. Особенностью является измельчение сетки только по окружному направлению.
Ключевые слова: многосеточный метод, моделирование шины.
This paper studies one geometric multigrid method in a rolling tire problem, finding components that make it reasonably effective for solving systems that are large enough for real simulation. Its distinct feature is grid refinement in circumferential direction only.
Key words: multigrid method, tire simulation.
При численном моделировании качения шины приходится решать контактные задачи на каждом шаге дискретизации по времени [1]. Для этого используется контактный итерационный алгоритм [2-4], на каждой итерации которого встает необходимость решения системы линейных алгебраических уравнений после дискретизации задачи методом конечных элементов. Размер системы, которую следует решать, составляет от 200 до 500 тысяч уравнений, причем число узлов в окружном направлении выбирается от 60 до 200. Такой размер системы позволяет использовать для решения системы многосеточные методы (ММ) [5]. Контактный алгоритм и внутренний итерационный процесс составляют вместе двухступенчатый итерационный метод, основное достоинство которого — небольшое число внутренних итераций.
Целью настоящей работы является исследование эффективности одного варианта ММ. Особенность данного варианта ММ состоит в измельчении сетки только в окружном направлении. Дело в том, что двумерная сетка в меридиональном сечении является очень сложной, что исключает применение геометрических ММ, когда измельчение идет по всем направлениям. Именно сильной неоднородностью и анизотропией материала обусловлена проблематичность применения геометрического ММ.
При использовании многосеточных методов система уравнений решается на нескольких сетках. Эффективность метода зависит от ряда факторов: функции перевода векторов с сетки на сетку; методов сглаживания и последующей коррекции ошибки; типа применяемых сеток и типа цикла. В то время как для многих модельных задач эти составляющие ММ хорошо исследованы, для нестандартных задач их правильный выбор требует проведения численных экспериментов.
Запишем решаемую линейную систему как
Lu = f. (1)
Простейшим вариантом многосеточного метода является V-цикл, составляющий основу для более сложных и эффективных ММ. Необходимо задать по крайней мере две сетки, грубую и мелкую, причем
1 Лопухин Константин Александрович — студ. каф. механики композитов мех.-мат. ф-та МГУ, e-mail: kostia.lopuhin@gmail.com.
2 Шешенин Сергей Владимирович — доктор физ.-мат. наук, проф. каф. механики композитов мех.-мат. ф-та МГУ, e-mail: shesheni@mech.math.msu.su.
мелкая должна являться измельчением грубой. За один шаг выполняется первоначальное сглаживание решения, потом проводят коррекцию ошибки на грубой сетке и завершающее сглаживание. Целью сглаживания является уменьшение высокочастотных компонент ошибки. Для сглаживания мы применяли метод Гаусса-Зейделя. Низкочастотные компоненты ошибки могут быть эффективно уменьшены на более грубой сетке. Сетки характеризуются шагами h и H соответственно. Будем обозначать сглаживание ш-й итерации um через й™ = SVl (uim,Lh,fh), а весь цикл — через ь^+1 = V(п^1,Lh,fh,vi,v2) = M^. Тогда многосеточный метод, применяемый для решения системы (1), запишется в виде следующего алгоритма [5].
1. Сглаживание U"m = SVl(п™, Lh, fh).
2. Коррекция на грубой сетке
um ^ rm = fh — Lhum ^ rm = Ih rm ^ LHwm = rm ^ wm = IHwm ^ um = um+wm.
3. Сглаживание um+1 = SV2(um,Lh,fh).
Невязка на ш-й итерации обозначена через rm, а поправка приближенного решения um — через wm. Оператор IH осуществляет проектирование с мелкой сетки на грубую, а оператор iH, наоборот, с грубой сетки на мелкую. Через Vi и V2 обозначено число итераций сглаживания в начале и в конце цикла. Верхний индекс у векторов соответствует номеру итерации, нижний — сетке, на которой задан вектор.
При коррекции на грубой сетке нет никакой необходимости решать уравнение точно. Без существенного влияния на скорость сходимости можно заменить точное решение уравнения на приемлемую аппроксимацию. Естественный способ получения этого приближения — опять использовать V-цикл, но уже на более грубых сетках. Эта идея может применяться ко все более грубым сеткам, что приводит к многосеточному методу. На самой грубой сетке уравнение для поправки решается прямым методом, во всех остальных случаях используется несколько итераций V-цикла. Количество этих итераций и схема переходов между сетками определяют тип многосеточного метода.
Огромное влияние на эффективность многосеточных методов оказывает способ измельчения сетки. Как правило, измельчение сетки выполняется по всем направлениям. Однако, как было сказано выше, этот способ неприменим в данной задаче, поэтому измельчение происходит только в окружном направлении. Таким образом, специфика рассматриваемого ММ состоит в том, что в меридиональном сечении двумерные сетки одинаковы, а в окружном направлении различаются лишь количеством узлов. Разбиение в одном направлении часто применяется для решения анизотропных задач, однако данный случай отличается тем, что разбиение идет не по прямолинейной, а по криволинейной оси, а также тем, что здесь при разбиении не уменьшается, а увеличивается анизотропия (неоднородность) задачи, поскольку отношение характерного размера сетки в окружном направлении и в меридиональном сечении возрастает при укрупнении сетки. В результате численных тестов было выявлено, что стандартные алгоритмы сглаживания, применяемые в геометрических ММ, становятся крайне неэффективными. Повысить эффективность удалось путем замены метода Гаусса-Зейделя на его блочный аналог.
Один блок переменных соответствует узлам одного меридионального сечения. Другими словами, обращение каждого блока соответствует решению двумерной задачи теории упругости в сечении. Для реализации этого метода требуется выделить из матриц, хранимых в разреженном формате (причем в силу симметрии матрицы хранится только верхнетреугольная часть), отдельные блоки и провести с ними различные матричные операции. Блочный метод можно записать следующим образом: ЛцХi = bi — BiX—i — CiXi+i, i = 1,...,ш, где ш — число блоков, Лц — блок матрицы, относящийся к блоку переменных Xi; Bi и Ci — подматрицы, расположенные выше и ниже блока Aii; Xi-i и Xi+i — соответственно подвекторы, расположенные выше и ниже Xi. Программа, реализующая многосеточный метод, написана на языке Fortran с использованием специальной высокоэффективной библиотеки для векторных и матричных операций. Матрица решаемой системы линейных уравнений (матрица жесткости) получена с помощью конечно-элементной аппроксимации задачи о деформировании шины при наезде на препятствие. Кроме того, для тестов и подбора параметров метода применялись матрицы, возникающие при решении еще двух задач — о симметричном раздувании внутренним давлением шины и однородного цилиндрического тела. Эти задачи проще вследствие осевой симметрии. Кроме того, задача о раздувании однородного цилиндра приводит к существенно лучше обусловленной системе с меньшим числом уравнений. Далее в статье анализируется поведение V-цикла, так как при некоторых естественных допущениях [5] из его хорошей сходимости следует хорошая сходимость и более эффективных в вычислительном смысле многосеточных методов. Описание результатов вычислительных экспериментов, которым, собственно, и посвящена статья, приводится ниже.
Рассмотрим сначала результаты применения метода для задачи о раздувании однородного цилиндра. V-цикл исследовался на сетках, состоящих из 32 и 16 сечений в окружном направлении. При этом
представляло интерес не только изучение собственно сходимости, но и сравнение блочного и обычного метода Гаусса-Зейделя. Для анализа сходимости использовалась норма невязки г = Ьй — / , где й обозначает приближенное решение. Нормой служила обычная евклидова норма (корень из суммы квадратов).
Оценить сходимость можно исходя из графиков, на которых изображено отношение норм к = ц||Гг^ц (см.
рис. 1). Видно, что даже в простой задаче блочный метод Гаусса-Зейделя существенно повышает скорость сходимости, увеличивая отношение к с 1,45 до 2,66. Однако настоящая проверка эффективности может быть получена только на матрице жесткости для шины, где применение блочного метода является просто необходимым для обеспечения сходимости.
Рис. 1. Сравнение обычного и блочного методов Гаусса-Зейделя для цилиндра: 1 — блочный метод; 2 — обычный метод
В задаче о шине в качестве начального приближения выступал вектор решения на грубой сетке, переведенный оператором IН на мелкую сетку. Его невязка была лишь в полтора раза меньше, чем у нулевого вектора начального приближения, но это позволяло избавиться от нескольких итераций (при нулевом приближении на первой итерации норма вектора невязки увеличивалась в 10 раз, и только со второй итерации начиналась сходимость).
Попытки использовать в задаче о шине уже отлаженный на цилиндре метод приводили к тому, что вычисления очень быстро расходились. Оказалось, что сглаживание методом Гаусса-Зейделя практически не уменьшает невязку: даже очень далекие от точного решения векторы почти не менялись. Предположение о том, что благодаря блочному методу Гаусса-Зейделя можно улучшить сглаживание, полностью подтвердилось в дальнейших экспериментах (см. рис. 2 и 3).
Применение блочного метода Гаусса-Зейделя позволяет добиться сходимости V-цикла для матрицы, соответствующей задаче о шине. Матрица жесткости в этой задаче получена при дискретизации реальной и весьма сложной внутренней структуры шины. Следует заметить, что на скорость и на само наличие сходимости влияет также ряд параметров, например количество итераций сглаживания до и после коррекции: и При (и\= (1,2) сходимость по-прежнему отсутствовала, но уже при проведении четырех итераций сглаживания (и\= (2, 2) метод начинал медленно сходиться. При (^1,^2) = (2, 3) и (VI ,^2) = (2, 4) скорость сходимости возросла и для 5 и 6 итераций сглаживания составила к = 1,62 и 1,79 соответственно (вместо 1,17 для 3 итераций сглаживания). Такая сходимость уже позволяет получить решение за небольшое число итераций, причем видно, что дальнейшее увеличение числа итераций сглаживания почти не влияет на скорость сходимости, увеличивая только время решения.
При использовании сеток с 128 и 64 сечениями скорость сходимости при (и\, = (2, 4) почти не уменьшается, что хорошо согласуется с тем, что для многосеточного метода скорость сходимости не должна зависеть от шага сетки (и от размера системы) [5].
Больший интерес представляет исследование сходимости ММ для несимметричной нагрузки, возникающей при моделировании наезда шины на препятствие. После первых вычислительных тестов выяснилось, что для обеспечения сходимости необходимо существенно увеличить число шагов сглаживания. В результате при 120 узлах сетки в окружном направлении хорошая сходимость достигается при об-
щем числе итераций сглаживания, равном 50 ,и2) = (5, 45)), что увеличивает время решения в 5 раз. По-видимому, это связано с тем, что большие искажения сетки, вносимые препятствием, снижают эффективность блочного метода Гаусса-Зейделя из-за большей зависимости между перемещениями в соседних сечениях. Однако сходимость имеет место, как следует из графиков, представленных на рис. 3.
Рис. 2. Сходимость У-цикла для шины: 1 — на 32 сечениях; 2 — на 120 сечениях при несимметричной нагрузке; 3 — на 32 сечениях при 6 итерациях сглаживания; 4 — на 128 сечениях
Рис. 3. Скорость сходимости У-цикла для шины: 1 — на 32 сечениях; 2 — на 120 сечениях при несимметричной нагрузке; 3 — на 32 сечениях при 6 итерациях сглаживания; 4 — на 128 сечениях
Основные результаты проведенного исследования заключаются в построении многосеточного метода, позволяющего решать очень плохо обусловленные линейные системы, возникающие при численном моделировании наезда шины на препятствие. Отличительными особенностями апробированного варианта ММ являются измельчение сетки только вдоль криволинейной окружной оси, а также применение для сглаживания блочного метода Гаусса-Зейделя. Такой вариант, возникший в результате проведения многочисленных численных тестов, оказался достаточно эффективным в рассматриваемой задаче о шине вследствие того, что переменные в сечении связаны между собой сильнее, чем переменные в соседних сечениях.
Работа выполнена при поддержке гранта РФФИ-ГФЕН, проект № 07-01—92111-ГФЕН.
СПИСОК ЛИТЕРАТУРЫ
1. Шешенин С.В. Трехмерное моделирование шины // Изв. РАН. Механ. твердого тела. 2007. № 3. 12-21.
2. Кравчук А. С. К теории контактных задач с учетом трения на поверхности соприкосновения // Прикл. матем. и механ. 1980. 44, вып. 1. 122-129.
3. Simo J. C, Laursen T. A. An augmented Lagrangian treatment of contact problems involving friction // Comput. and Struct. 1992. 42, N 1. 97-116.
4. Rothert H., Idelberger H., Jacobi W., Laging G. On the contact problem of tires, including friction // Tire Sci. and Technol. 1985. 13, N 2. 111-123.
5. Trottenberg U., Oosterlee C. W., Shuller A. Multigrid. London: Academic Press, 2001.
Поступила в редакцию 27.02.2008
