На правах рукописи

Сагдеева Юлия Альбертовна

МЕТОД ЧИСЛЕННОГО ОПРЕДЕЛЕНИЯ

ОСРЕДНЕННЫХ ХАРАКТЕРИСТИК КОМПОЗИТОВ

НА ОСНОВЕ ВЕЙВЛЕТ-ПРЕОБРАЗОВАНИЯ И

МЕТОДА КОНЕЧНЫХ ЭЛЕМЕНТОВ

05.13.18 — Математическое моделирование, численные методы

и комплексы программ

АВТОРЕФЕРАТ

Казань – 2007

Работа выполнена в Институте прикладной механики Уральского отделения РАН (г.Ижевск).

Научный руководитель кандидат физико-математических наук, старший научный сотрудник Копысов Сергей Петрович

Официальные оппоненты: доктор физико-математических наук, доцент Даутов Рафаил Замилович, кандидат физико-математических наук Мартыненко Сергей Иванович

Ведущая организация Институт математического моделирования РАН, г. Москва.

Защита диссертации состоится «25» октября 2007г. в на заседании диссертационного совета Д 212.081.21 в Казанском государственном университете по адресу: 420008, г.Казань, ул. Кремлевская, 18, корп.2.

С диссертационной работой можно ознакомиться в научной библиотеке им. Н.И. Лобачевского Казанского государственного университета.

Автореферат разослан «20» сентября 2007г.

Ученый секретарь диссертационного совета Д 212.081.21, к.ф.-м.н., доцент Задворнов О.А.

Общая характеристика работы

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

Существующие аналитические оценки свойств композитов (например, оценки Хашина-Штрикмана и оценки Фойгта-Рейсса для упругих констант, тепловых и фильтрационных свойств), как правило, дают достаточно широкий диапазон возможных значений свойств материала и могут использоваться только для грубой оценки.

В настоящее время разработаны численные методы получения эффективных характеристик материалов с периодической структурой в задачах линейной упругости, теплопроводности, диффузии и др. — это методика асимптотического осреднения, описанная Н.С. Бахваловым, Г.П. Панасенко, Б.Е. Победрей, Э. Санчес-Паленсией и А. Бенсусаном.

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

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

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

Одно из главных преимуществ вейвлет-преобразования заключается в возможности получать представление величин на интересующем уровне масштаба. С помощью вейвлет-преобразования получают осредненное представление функции (грубый масштаб – «низкое разрешение») и выделяют ее локальные компоненты (мелкий масштаб – «высокое разрешение»). Данное свойство преобразования позволяет ввести многомасштабный анализ исследуемой функции или анализ с переменным разрешением. Свойства вейвлетов позволяют предположить, что вейвлет-преобразование будет полезным и при осреднении решений уравнений в частных производных.

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

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

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

Научная новизна работы состоит в следующем:

— разработана методика численного осреднения линейных эллиптических краевых задач второго порядка с быстро осциллирующими коэффициентами для вычисления эффективных характеристик с помощью одномерного и двумерного вейвлет-преобразования Хаара и метода конечных элементов в областях периодической и непериодической структуры;

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

— численным моделированием получены осредненные значения модуля Юнга, коэффициента Пуассона, коэффициентов теплопроводности и проницаемости в одномерном и двумерном случаях для композитов со случайной структурой, материалов с включениями разной формы (квадратные, круглые и ромбические включения), разной объемной долей (объемная доля включений составляла от 7% до 50%), различного взаимного расположения включений (симметричное и несимметричное), когда свойства составляющих компонент композита различаются на порядки.

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

Апробация работы. Основные положения диссертационной работы обсуждались на следующих научных конференциях: Всероссийской молодежной школе-конференции «Численные методы решения линейных и нелинейных краевых задач», 27 июня – 1 июля 2003, Казань; Конференции молодых ученых ФТИ УрО РАН «КоМУ-2004», 8–10 декабря 2004, Ижевск; Международной конференции по избранным вопросам современной математики, 4–8 апреля 2005, Калининград; 14 Зимней школе по механике сплошных сред, 28 февраля – 3 марта 2005, Пермь; VI International Congress on Math. Modeling, 20–26 сентября 2004, Нижний Новгород; 14-й Всероссийской школе-конференции молодых ученых, 4 октября – 7 октября 2005, Пермь; Научной конференции «Теория управления и математическое моделирование», 3–8 июля 2006, Ижевск; III научно-практической конференции «Проблемы механики и материаловедения», 14–15 июня 2006, Ижевск; Всероссийской молодежной школе-конференции «Численные методы решения задач математической физики», 26 июня – 1 июля 2006, Казань; IX Всероссийском съезде по теоретической и прикладной механике, 22–28 августа 2006, Нижний Новгород; 15 Зимней школе по механике сплошных сред, 26 февраля – 3 марта 2007, Пермь.

Публикации. По теме диссертации опубликовано 19 печатных работ. Из них в рецензируемых журналах 4 работы, в журналах, рекомендованных ВАК, 1 работа.

Благодарности. Диссертационная работа выполнена при финансовой поддержке РФФИ (проекты №№07-01-96069, 06-07-89015, 02-07и УрО РАН.

Структура и объем работы. Работа состоит из введения, четырех глав, 15 параграфов (нумерация параграфов по главам), заключения, списка литературы и содержит 56 рисунков, 13 таблиц. Объем работы 124 страницы. Библиографический список включает 91 наименование.

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

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

В п. 1.1 рассмотрены преимущества и недостатки экспериментального определения свойств композиционных материалов. В п. 1.2 описаны уравнения с быстро осциллирующими коэффициентами, приведены постановки задач статической теории упругости, стационарной теории теплопроводности и теории фильтрации, а также введено понятие эффективных оценок. В п. 1.3 описаны известные верхние и нижние аналитические оценки для упругих эффективных модулей, коэффициента теплопроводности и проницаемости, использующие информацию только об объемных долях включений и материала. Это оценки ФойгтаРейсса (называемые также смесевыми), получаемые на основе закона сохранения энергии, и более точные оценки Хашина-Штрикмана, получаемые из вариационного принципа Хашина-Штрикмана. В п. 1.4 рассмотрены численные методы вычисления оценок свойств, среди которых наиболее проработанным является асимптотический метод осреднения, предложенный Н.С. Бахваловым, Г.П. Панасенко1 и Б.Е Победрей. Рассмотрены методы давлений и ренормализации получения осредненных фильтрационных свойств, а также многосеточный метод осреднения дифференциальных уравнений.

Вторая глава посвящена описанию многомасштабного анализа на основе вейвлет-преобразования Хаара и метода конечных элементов. В п. 2.1 вводится понятие семейств масштабирующих функций j,k (x) = 2j/2 (2j x k), k, j Z, и вейвлетов j,k (x) = 2j/2 (2j x k), k, j Z. В качестве вейвлет-базиса был выбран базис Хаара, поскольку он обладает свойствами ортогональности и симметрии:

1 Бахвалов Н. С., Панасенко Г. П. Осреднение процессов в периодических средах.

— М.: Наука. — 1984.

В п. 2.2 описана многомасштабная декомпозиция функции. Известно, что пространства, порождаемые вейвлетами и масштабирующими функциями, Vj = span{j,k (x)} и Wj = span{j,k (x)}, j, k Z удовлетворяют соотношению Vj = Vj1 Wj1, где — прямая сумма.

Поэтому функция fj (x) Vj представляется в виде:

Пусть произвольный вектор cj имеет размер 2j. Этот вектор представим как последовательность коэффициентов разложения некоторой функции fj (x) Vj :

Коэффициенты разложения в пространствах Vj, Vj1 и Wj1 связаны между собой: cj1 = Pj cj, dj1 = Qj cj. Вектор cj1 является проекцией вектора cj на пространство Vj1, т. е. его огрубленным представлением. Вектор dj1 соответствует уточняющим коэффициентам.

Матрица Wj = является матрицей вейвлет-преобразования.

Вейвлет-преобразование Хаара ортогонально: Wj1 = WjT.

В п. 2.3 рассмотрено применение вейвлет-преобразования для осреднения решения эллиптических дифференциальных уравнений в комбинации с методом конечных разностей (МКР) или методом конечных элементов (МКЭ) на примере одномерного уравнения dx k(x) du = f, x [0, 1]. В работе вейвлет-преобразование впервые использовалось в комбинации с МКЭ. Методом конечных элементов исходное дифференциальное уравнение Lu = f заменялось аппроксимирующим его на некоторой сетке дискретным уравнением Lh uh = fh, решение которого сводится к решению системы линейных алгебраических уравнений где матрица A имеет размер 2 2, симметрична и положительно определена. Затем к системе (1) применялось вейвлет-преобразование Хаара Wn :

Система (2) имеет блочный вид. Вектор неизвестных разбивается на две составляющих компоненты — вектор осредненных неизвестных и вектор «деталей» (уточняющих компонент):

Вводились следующие обозначения, которые используются далее при описании вычислительных особенностей:

Система (3) принимает вид Исключая u1 из первого уравнения, получим Разрешив (4), получаем искомое осредненное решение u2. Применяя вейвлет-преобразование несколько раз, на каждом шаге имеем решение с разной степенью осреднения.

В п. 2.4 вводится расширение вейвлет-представления с одномерного случая на двумерный в пространстве L2 (R2 ). Базисные функции образуются за счет комбинаций тензорных произведений базисных функций для одномерного случая. Разложение системы уравнений проводится аналогично одномерному случаю. Преобразованный вектор неизвестных состоит из четырех компонент — три компоненты содержат информацию о быстро осциллирующих составляющих решения по трем возможным направлениям, а одна компонента описывает осредненное поведение решения.

В п. 2.5 предложена и обоснована численная методика получения осредненных характеристик композитов. Методика основана на совместном применении МКЭ, вейвлет-осреднения системы, описанного в пп. 2.3 и 2.4, и решении краевых задач при определенных граничных условиях:

а) Модуль Юнга (одномерный случай). Рассматривалась модельная задача: одноосное напряженное состояние упругого стержня где p — давление, l — половина длины стержня.

Для однородной и изотропной среды с модулем Юнга E(x) = E = const, решение задачи имеет вид Предполагалось, что композитный стержень реагирует на внешнюю нагрузку как однородный, поэтому из (6) следует, что E = u p.

Эффективный модуль Юнга E определялся на основе численного решения методом конечных элементов и вейвлет-осреднения задачи (5) по соотношению Ei — модуль упругости в узле сетки xi, uw = u(xi ) — значение перемеi щения в узле сетки xi, осредненное с помощью вейвлет-преобразования, Nn — число узлов сетки.

В соотношении (7) используются координаты узлов новой сетки осредненного поля перемещений. Координаты этих узлов для стержня определяются согласно схеме на рис. 1. Вектор неизвестных системы (1), сформированной в МКЭ, дает значения неизвестной функции в узлах сетки. В случае вейвлет-преобразования Хаара на каждом шаге число неизвестных в системе сокращается в два раза, а каждый столбец (строка) матрицы осредненной системы уравнений получается c помощью преобразования двух соседних столбцов (строк) матрицы исходной системы. Две соседние переменные в векторе неизвестных исходной системы преобразуются в одну неизвестную осредненной системы.

Переменные в векторе неизвестных упорядочены по номеру узла, то есть на каждом шаге осредняются значения, соответствующие соседним узлам сетки. Полученные осредненные значения уже нельзя сопоставить какому-либо узлу исходной сетки, т. к. они сочетают в себе значения для двух узлов. Новый узел помещается посередине между двумя осредняемыми узлами. Недостатком данного подхода является то, что осредненное решение определено лишь на части исходной области, поскольку узлы новой сетки смещаются к центру. Шаг новой сетки h1 будет отличаться от старого h0 в два раза: h1 = 2h0. Этот недостаток преодолевается с помощью расширения физической области на фиктивную область. Например, решение продолжается за пределы новой сетки нечетным образом относительно своих значений на границе.

б) Модуль Юнга и коэффициент Пуассона (двумерный случай). Эффективные значения двух констант — модуля Юнга E и коэффициента Рис. 1. Схема изменения сетки при осреднении (одномерный случай) Пуассона — в плоскости пластины вычислялись по следующей схеме.

Рассматривалась краевая задача с граничными условиями ij (x)nj (x)|x=±l = p, ij (x)nj (x)|y=0,y=l = 0.

Эта задача в силу симметрии эквивалентна задаче на четверти области с краевыми условиями Методом конечных элементов и вейвлет-осреднением вычисляется осредненное поле перемещений в узлах конечно-элементной сетки.

Обозначим через u1 и u2 составляющие x и y осредненного пеij ij ремещения u в узле (xi, yj ). В условиях одноосного растяжения вдоль оси x модуль Юнга и коэффициент Пуассона в направлении оси x в узле (xi, yj ) определяются по формулам Используя (10), находим значения констант в каждом узле и затем вычисляем эффективное значение Координаты новой сетки в двумерном случае определяются аналогично одномерному случаю. За один шаг осредняются значения сразу в четырех узлах, причем число неизвестных в системе сокращается в четыре раза.

Таким образом, методика вейвлет-осреднения для получения упругих эффективных свойств имеет вид:

Шаг 1. Рассматривается задача (8), (9) в перемещениях для одноосного напряженного состояния композитной пластины или стержня;

Шаг 2. Методом конечных элементов формируется матрица системы линейных алгебраических уравнений (СЛАУ) (1);

Шаг 3. К полученной СЛАУ применяется вейвлет-преобразование Хаара (2); в результате имеем систему для осредненного вектора (4);

Шаг 4. Существует два варианта дальнейших вычислений: а) система (4) разрешается относительно вектора неизвестных и в качестве результата имеем вектор осредненных перемещений, б) к системе (4) вновь применяется вейвлет-преобразование; таким образом, могут быть получены осредненные значения перемещений на нескольких масштабах;

Шаг 5. Осредненные модуль Юнга пластины (или стержня) и коэффициент Пуассона пластины определяются из (11).

в) Коэффициент теплопроводности (проницаемости). Эффективный коэффициент теплопроводности определялся на основе вейвлет-осреднения и решения методом конечных элементов двумерной стационарной краевой задачи теплопроводности Обозначим через Tij — осредненное с помощью вейвлет-преобразования поле температуры в узле сетки (xi, yj ). Эффективный коэффициент теплопроводности рассчитывался как среднее значение по всем узлам Для задачи фильтрации при условиях выполнения закона Дарси, коэффициент проницаемости K определяется по такой же схеме, как отношение произведения начальной скорости фильтрации на координату узла к значению давления в узле.

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

В п. 3.2 рассмотрено пять задач осреднения коэффициентов и решений одномерного дифференциального уравнения вида u(x) Рис. 2. Решения: а) вейвлет-решение, конечно-элементное решение, решение по асимптотическому методу; б) увеличенный фрагмент интервала В Задаче 1 рассматривалось вейвлет-осреднение для однородного материала k(x) = k = const. В этом случае осредненный коэффициент k w естественным образом совпадает с начальным модулем материала k. В Задаче 2 получено осредненное решение уравнения (13) с быстро меняющимся периодическим коэффициентом k(x) ( k(x) принимало значение 0 или 1), u0 = 0, u1 = 0. На рис. 2 показаны графики решения дифференциального уравнения, соответствующие конечноэлементному решению на сетке из Nn = 256 узлов, решению, осредненному с помощью j = 3 шагов вейвлет-преобразования Хаара (число узлов сетки Nn = 32 ), и решению, полученному с помощью асимптотического метода осреднения. Обе модели осреднения хорошо аппроксимируют осциллирующее решение и практически совпадают. В Задачах 3–4 получены осредненные решения уравнения (13) для непериодического коэффициента k ; сделан вывод о том, что в этих задачах вейвлетосреднение описывает поведение решения точнее, чем асимптотическое осреднение. В Задаче 5 исследовалось изменение осредненного модуля Юнга k(x) = E(x) для стержня в зависимости от числа включений nv ( nv = 1,..., 10 ) и шага осреднения j ( j = 1,..., 12 ) при одинаковой объемной доле включений c(v) = 0.3. Согласно графикам на рис.

3 наблюдается незначительный рост значений модуля E на последних шагах осреднения. Также можно заметить, что с ростом числа включений уменьшается диапазон вариаций значений модуля, так как чем меньше характерный размер неоднородности (включения) в сравнении с размером исследуемой области, тем меньше в среднем оказывает влияние эта неоднородность. В примере с одним включением этот модуль изменяется в диапазоне E [1.392, 1.492], тогда как для примера с десятью включениями модуль изменяется в диапазоне E [1.479, 1.49], т.е. практически не зависит от шага осреднения.

В п. 3.3 рассматривались результаты двумерного вейвлет-осреднения для задач теории теплопроводности, фильтрации и упругости. В Задачах 6–8 проводится осреднение коэффициента проницаемости K(x) для уравнения (12) для расчетных областей трех типов: области с квадратным включением (доля включения c(v) = 0.25 ), области шахматной структуры и области, в которой проницаемость распределена по закону d = z ln(a) (рис. 4–6).

лиз для разных отношений проm) ницаемости материала K циенте при сравнении с асим- 1. вышает 1%. Полученный коэффициент удовлетворяет оценРис. 3. Модуль Юнга E в зависимости кам Фойгта-Рейсса и Хашинаот шага осреднения j Штрикмана.

В случае шахматной структуры (задача 7) для материалов с проницаемостью K1 и K2 известна теоретическая эффективная проницаемость среды K g = K1 K2 (среднее геометрическое двух величин). Получены следующие результаты для вейвлет-осреднения на сетке 128 128 узлов, когда различие в модулях материала составляет один и два порядка: при K1 = 100 и K2 = 1 осредненный с помощью вейвлет-преобразования коэффициент K w = 9.5 ( K g = 10, относительная ошибка = 5% ), при K1 = 10 и K2 = 1 осредненный с помощью вейвлет-преобразования модуль K w = 3.08 ( K g = 3.16, относительная ошибка = 2.5% ). Для сравнения при K1 = 100 и K2 = 1 оценки Хашина-Штрикмана задают широкий интервал 2. 34.21, а оценки Фойгта-Рейсса (арифметическое и гармоническое среднее) 1.98 50.5 ; во втором случае при K1 = 10 и K2 = Рис. 4. Квадратное Рис. 5. Шахматная Рис. 6. Случайное оценки Хашина-Штрикмана задают интервал 2.38 K 4.19, а оценки Фойгта-Рейсса 1.82 K 5.5. Коэффициенты, полученные с помощью вейвлет-осреднения, удовлетворяют этим аналитическим оценкам.

В Задаче 8 для распределения коэффициентов материала d = z ln(a), где z — случайно распределенная величина на интервале (0, 1), известно, что значение эффективного коэффициента проницаемости определяется как геометрическое среднее значений данного распределения. Результаты моделирования для разных значений параметра a = 2, 5, 10 приведены в табл. 1. Во втором-четвертом столбцах приведены статистические данные о начальном распределении проницаемости — максимальное значение Kmax, минимальное значение Kmin, среднеквадратическое отклонение. В пятом-восьмом столбцах показаны соответственно арифметическое среднее K a, гармоническое среднее K h, геометрическое среднее K g и значение K w, полученное с помощью вейвлет-осреднения. Результаты моделирования близки к геометричеg ским средним — относительная погрешность |K K g | мала и равна для разных a соответственно 3.5%, 3.6%, 3.7%. Арифметическое и гармоническое средние дают настолько широкий интервал оценок, что их использование не имеет смысла. Оценки Хашина-Штрикмана в данном случае не применимы, так как они могут использоваться только для двухфазного материала.

В Задаче 9 сравнивается распределение температур на разных шагах осреднения для пластины с девятью включениями (рис. 7) (коэффициент теплопроводности матрицы K (m) = 1, включения K (v) = 0.1 ). На рис. 8 показано конечно-элементное решение на сетке 64 64 узла. На Таблица 1. Сравнение аппроксимаций эффективной проницаемости для случайного распределения коэффициента проницаемости композита рис. 9 и 10 показано распределение температуры после двух и четырех шагов осреднения. После четырех шагов, когда сетка становится достаточно грубой, решение сглаживается и описывает поведение температуры в среднем.

В Задаче 10 исследовалось изменение коэффициента теплопроводности в зависимости от изменения формы включения в пластинке (квадратная, круглая и ромбическая осредненные коэффициенты теплопроводноРис. 7. Девять включений сти: а) для включения квадратной формы K as = 1.548 (асимптотическое осреднение), K w = 1.522 (вейвлетосреднение) и K b = 1.598 (метод многосеточного осреднения)2 ; б) для включения круглой формы K as = 1.516, K w = 1.514 и K b = 1.563 ; в) для включения в форме ромба K as = 1.573, K w = 1.537 и K b = 1.608.

Отличие коэффициентов K as и K w мало и составляет от 1 до 2%.

Следует отметить, что теоретические оценки будут давать одинаковый результат в независимости от формы включения. Оценки ХашинаШтрикмана для трансверсально изотропного материала имеют следующие значения: нижняя граница Klhs = 1.514, верхняя — Ku = 2.394.

Полученные с помощью вейвлетов оценки близки к нижней границе и удовлетворяют диапазону [Klhs, Ku ]. Нижняя и верхняя оценки Фойгтаhs Рейсса (гармоническое и арифметическое среднее) равны соответственно Klf r = 1.29 и Klf r = 3.25. Эти оценки также выполнены.

Задача 11 содержит сравнение экспериментальных и численных оценок осредненного коэффициента теплопроводности одиннадцатипустотной области. Вычисленный по методу вейвлет-осреднения эффективный коэффициент теплопроводности равен K = 0.62. В эксперименте K = 0.64. Расхождение расчетных и опытных данных составило 3.6%.

В Задаче 12 вычислялись осредненные значения модуля Юнга E 2 Moulton J. D., Dendy J. E., Hyman J. M. The black box multigrid numerical homogenization algorithm // Journal of Computational Physics. – 1998. – Vol. 142. — P. 80–108.

Конечно-элементное Вейвлет-решение Вейвлет-решение Таблица 2. Зависимость упругих констант от объемной доли включения и коэффициента Пуассона пластинки с одним включением для разной объемной доли включения ( E (m) = 1, (m) = 0.3, E (v) = 0.5, (v) = 0.3 ). Результаты вейвлет-осреднения после j = 6 шагов (область стягивается в точку) представлены в табл. 2. В той же таблице показаны результаты асимптотического осреднения. Различия для значений модуля Юнга при объемной доле c(v) = 0.1 имеют порядок 5%, тогда как при c(v) = 0.5 отклонение близко к 15%. Такое поведение при большой объемной доле объясняется тем, что в этом случае напряженное состояние в пластине уже больше нельзя считать близким к одноосному.

В Задаче 13 исследовалось изменение модуля Юнга E и коэффициента Пуассона для пластинки с nv = 12 включениями, c(v) = 0.27. В этом случае относительная погрешность для E по сравнению с асимпK11 K Рис. 11. Портрет Рис. 12. Матрица Рис. 13. Структура тотическим методом составляет 1.5%. Коэффициент Пуассона практически не меняется.

В Задаче 14 проводилось сравнение аналитических оценок ХашинаШтрикмана и Фойгта-Рейсса для упругих констант (модуль Юнга, коэффициент Пуассона, модуль сдвига, объемный модуль) для пластинки с четырьмя круглыми включениями при разных соотношениях начальных модулей Юнга. Полученные значения попадают внутрь интервала, задаваемого оценками Фойгта-Рейсса и близки к нижним или верхними оценками Хашина-Штрикмана.

В четвертой главе описаны вычислительные особенности, присущие вейвлет-преобразованию. Основными операциями при осреднении являются преобразование (2) матрицы A и вектора правых частей, обращение матрицы K11, матричные и матрично-векторные операции умножения и сложения. В п. 4.1 исследовалась структура матриц A, K11, K12, K22, S, участвующих в вычислении, на примере двумерной задачи теории упругости. Матрица системы A сильно разрежена, а ненулевые элементы компактно размещаются в трех «лентах-полосах».

После применения вейвлет-преобразования структура матрицы меняется — она становится блочной, причем каждый малый блок имеет портрет, близкий к портрету исходной матрицы и несет информацию или об уточнениях (локальные свойства), или об осреднении исходной матрицы. На рис. 11 представлена структура преобразованной матрицы.

Пунктирной линией показано разделение матрицы W AW T на блоки K11, K12, K12, K22. Каждая из подматриц имеет блочную структуру ( K11 – девять блоков, K12 и K21 – три блока).

При обращении матрицы K11 возникает заполнение (рис. 12), для хранения которого требуется дополнительная память. Система уравнений для дополнения Шура также становится более плотной (рис. 13).

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

Распределение таково, что более 70% всех элементов этих матриц по модулю меньше величины 106.

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

Анализ ошибки вычислений при задании порогового значения для всех матриц, участвующих в алгоритме, на примере задачи фильтрации (Задача 6) дает следующие результаты. Вычислялась относительная погрешность получаемого значения осредненного коэффициента проницаемости и полей решений. В качестве «точного» выбирались поля решений и коэффициент проницаемости при пороговом значении t = 1016.

Пороговое значение в расчетах выбиралось равным t = Kmax /10s, где s = 5, 7, 9, 11, 13 и Kmax = max{K (m), K (v) }. Расчеты были проведены для разных значений параметра = 0.01, 0.1, 10, 100, задающего свойства материала (см. условия задачи 6) и значениях Kmax = 1, 1, 10, 100.

Определялись значения относительной погрешности полей давлеNn Nn pi — поле давления в i -ом узле, подсчитанное для порогового значения t = Kmax /10s, p0 — поле давления в i -ом узле, подсчитанное для контрольного порогового значения t = 1016. Для всех выбор порогового значения t c s 7 приводит к погрешности, не превышающей 3%, причем в большинстве случаев погрешность составляет всего лишь несколько десятых процента. При этом от 50 до 90% всех ненулевых элементов имеют модуль в диапазоне от Kmax /107 до 1016.

Таким образом, значительное сокращение занимаемой памяти для хранения элементов не приводит к большой погрешности результатов и существенного накопления ошибки c ростом шага осреднения не происТаблица 3. Относительная ошибка при вычислении коэффициента проницаемости в зависимости от и порогового значения t Kmax Kmax /105 Kmax /107 Kmax /109 Kmax /1011 Kmax /