Диссертация (Моделирование пространственных течений в газовых трактах с использованием адаптивных сеток), страница 6

Если его эффективный размер равен 0 – требуемая триангуляцияпостроена.4. Если координаты центра описанной окружности p0 принадлежаттреугольнику, делим треугольник на три, добавлением этой точки.5. В противном случае определяем ребро самое длинное ребро. Еслиребро не принадлежит границе, делим его пополам, тем самымпреобразуя два образующих его треугольника в четыре. Если ребропринадлежит границе, делим ребро пополам, тем самым преобразуяодин приграничный треугольник в два.6.

После всех этих операций вычисляем эффективный размер новыхтреугольников и добавляем их в отсортированный список всехтреугольников триангуляции и удаляем старые треугольники из списка.Все рёбра, принадлежащие новым треугольникам триангуляции,проверяем на возможность переворота. Переходим к пункту 1.35Число треугольников < NБерём первый треугольник из отсортированного массиваЭффективный размер треугольника > 0Центр описанной окружности p0 принадлежит треугольникуНаибольшее ребро l принадлежит границеДелимтреугольник натри,добавлениемточки p0Удаляем ребро l и двасоседнихтреугольника.Разбиваем получившийсятреугольникна4добавлениемточкицентра ребра lРазбиваем ребро l иприлежащий треугольникна дваПроверяем все рёбра, принадлежащие новым треугольникам, на возможностьпереворотаТриангуляция завершенаРисунок 2.1.10 - Блок-схема алгоритма триангуляции с ограничениями2.7.

Алгоритмы преобразования триангуляции для уменьшенияразрешающей способности сеткиДанный алгоритм позволяет увеличить минимальный размер рёбер втриангуляции.Вкачествеисходныхданныхзадаётсяминимальнаяэффективная длина ребра.Все рёбра, принадлежащие триангуляции, сортируются по величинеэффективной длины. Отсортированные рёбра хранятся в виде красно-чёрногодерева, позволяющего добавлять и удалять элементы со сложностью log(n).Рёбра, образующие границу области, сохраняются отдельно. Алгоритмпреобразования делится на две части: сначала преобразуются внутренниерёбра, после – граничные.36Алгоритм преобразования внутренних рёбер:1.

Берётся первое внутреннее ребро из отсортированного списка рёбер (снаименьшей эффективной длиной).2. Если эффективная длина больше минимально допустимой – переходимк преобразованию граничных рёбер.3. Проверкавозможностиудаленияребраизтриангуляциибезпересечений рёбер и вырожденных треугольников.4. Если возможно, удаляем ребро, изменённые рёбра проверяются навозможность переворота и добавляются в отсортированный массив.Возврат к пункту 1.5. Если удаление ребра невозможно, устанавливаем размер этого ребрабольше минимально допустимого;6. Для каждого изменённого ребра проверяется возможность переворота,заново вычисляется эффективная длина и определяется положение вотсортированном списке.Алгоритм преобразования приграничных рёбер:1.

Берётся первое ребро из списка граничных рёбер.2. Если эффективная длина больше минимально допустимой, алгоритмуспешно завершён.3. Для выбранного ребра и каждого из соседних граничных рёбервысчитывается расстояние между прямой, проходящей через крайниеточки и центральной точкой. Если расстояние меньше допустимойточности, а ячейки, образованные этими рёбрами, имеют общее ребро,производится слияние ячеек. В противном случае ребру присваиваетсяэффективный размер больший минимально допустимого.4. Для изменённых рёбер проверяется возможность переворота и занововычисляется эффективная длина.372.8.Интерполяция величин на неструктурированной сеткеИнтерполяция на неструктурированной сетке позволяет получитьзначение величины в узлах сетки по значениям в центрах конечных объёмов.Интерполяция производится локальными методами.

Рассмотрим набор точекp i в пространстве D. Пусть f i - значения функции в этих точках. В этомслучае значения функции в точке p может быть определено каксредневзвешенное значение: ff q  iii, где i i1p  pi2.(2.10)iДля более точной интерполяции используется модифицированныйалгоритм Шепарда [45].

В этом случае интерполяция осуществляется поформулам:f x    f x  ,ii(2.11)i R  p  pii   x R x x  xi2.9.2 , Rx  max p  pi .(2.12)Аппроксимация оператора градиента на неструктурированнойсеткеПри решении задачи Навье-Стокса, уравнений диффузии, припостроении схем высокого порядка точности и адаптации расчётной сетки кособенностямтечениятребуетсяопределитьзначенияфизических величин. На неструктурированной сеткеградиентовданная задачаосложняется неравномерностью ячеек и сложностью расширения шаблона.Интересный способ вычисления градиента функции на произвольнойсетке представлен в работе [43], в которой предлагается задавать величины в38виде комплексных чисел, действительной частью которых являетсяфизическая величина, а комплексной – малое возмущение. Все операциипроизводятся в пространстве комплексных чисел.

В этом случае, как следуетиз формулы Тейлора, значение производной вычисляется со вторымпорядком точности по формуле:f ' x  Im f x  hi  O h2 .hi (2.13)Данный метод имеет высокую точность, но требует повышенногорасхода машинных ресурсов в связи с отсутствием аппаратной поддержкиопераций над комплексными числами.Вторым часто используемым методом для вычисления градиентасеточной функции является метод Грина-Гаусса [23,44,47,48], основанный наиспользовании теоремы о дивергенции: dV   dS .(2.14)SВ этом случае значение оператор градиента  в центре ячееквычисляется по формуле: 1VfSf ,(2.15)fгде V – объём контрольного объёма, S – площадь грани.Согласно работам [23, 47], данный метод при использовании нанеструктурированных сетках может иметь низкую точность и даётсущественно немонотонное решение в случае большого отношения размеровэлементов.

Однако, данный метод показывает хорошие результаты в задачахадаптации расчётной сетки к особенностям течения, в которых используетсяне сам градиент, а его модуль.39Значения градиентов, полученных методом Грина-Гаусса, может бытьуточнено при помощи метода Фринка [44]. Для этого значения величининтерполируются в узловые точки сетки:Niri 1Ni1i 1 ri,(2.16)где ri – расстояние до центра ячейки i.Далее для расчёта градиентов применяется метод Грина-Гаусса, вкотором интегрирование потока через границы ячейки производится припомощи метода трапеций 2-го порядка точности: 1 f1 f 2 2V f  S f ,(2.17)где  f 1 ,  f 2 - значения в узловых точках начала и конца ребра.Такая модификация значительно увеличивает точность определенияградиентовзасчётрасширенногошаблона,нотребуетбольшихвычислительных ресурсов.Такжеградиентынанеструктурированнойсеткемогутбытьрассчитаны методом наименьших квадратов, основанном на построенииградиента в точке центра масс ячейки P, который минимизирует суммаквадратов разностей между значениями в ячейках, принадлежащих шаблону,и значениями, полученными экстраполированием из точки P:  r    P   P r  rP  .(2.18)Линейная система уравнений для получения градиента физическойвеличины  xyzTметодом МНК в трёхмерном случае имеет вид:404  4 4 4222xxyxz k2 x k  kkkkkkk x  k kyzk 1k 1k 1 k41444 222yyz k2 y k  k ,   k x k y k kkkkkxyzk 1k 1k 1 k 41444 222xzxzz k2 z k  kkkkkk x  k k kyzk 1k 1k 1k 1(2.19)гдеωk – весовой коэффициент, значения которого вычисляются согласноследующим формулам:1) ωk=1;(2.20)12)  k  r f2;(2.21) f S fr f  rL3)  k  2 ,  f  r  r  r  r , J P  якобиан .J P r ffRfL(2.22)Сравнение точности метода МНК с различными коэффициентамиприведено в работах [41], [44], [46], [47], [48].

В данной работе в качествевесового коэффициента используется формула (2.20) в случае равномернойсетки и (2.21) для сетки с адаптацией.Сравнение методов Гаусса и МНК для неструктурированной сеткиприведено в работах [23], [43], [47], [48]. На основании этих работ можносделать вывод, что оригинальный метод Гаусса применим только для«хороших» сеток, соседние ячейки которых не имеют большой разницы вразмерах.

Методы Фринка и МНК показывают схожие результаты. В даннойработе метод Гаусса используется для расчёта градиентов при проведенииадаптации сетки, а МНК – для вычисления градиентов при решенииуравнений Навье-Стокса.Значениявеличинвцентревнутреннихинтерполяцией величин из смежных ячеек:41гранейопределяется f   L 1      R  .(2.23)Для невзвешенной интерполяции интерполяционный фактор ω = ½.r f  rLДля линейной интерполяции     .r f  rR  r f  rLДляполучениявторогопорядка(2.24)интерполяциииспользуетсяразложение функции в ряд Тейлора в окрестности центра ячейки (9).Интерполяция второго порядка точности достигается при использованииформулы: f   L 1      R    L 1      R    r f  rL 1     rR(2.25)Необходимо отметить, что выбор метода реконструкции градиентазависит как от решаемой задачи, так и от используемой сетки.

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

Вто же время значения, получаемое при помощи этого метода, в ряде случаевоказывается более точными в связи с использованием расширенногошаблона.2.10. Примеры построения расчётных сетокНарисунке2.12приведенатриангуляциярасчётнойобласти,представляющей из себя береговую линию озера, а так же гистограммыраспределения параметров B и минимального угла. Данный пример являетсяклассическимтестомдляопределениятриангуляции с заданным минимальным углом.42возможностипостроенияРисунок 2.12 – Триангуляция береговой линии озера.На рисунке 2.13 приведён пример адаптации расчётной сетки кособенностямтечениявтрактевоздухозаборногоустройства.Рассчитывалось стационарное течение, возникающее при сверхзвуковомобтекании ВЗУ. Невязка рассчитывалась как разница втекающих ивытекающихпотоковчерезграницырасчётнойобласти.Сеткаперестраивалась каждый раз, при достижении величиной невязки значенияменьше1 10 5 ,коэффициентадаптациивнутрирасчётнойобластирассчитывался как модуль градиента плотности, а максимальная величинакоэффициентаадаптацииувеличиваласьпропорциональнономеруадаптации.Рисунок 2.13– Расчётная сетка, адаптированная к параметрам течения.43Полученные результаты сравниваются с результатами, полученными врезультате работы общераспространённых общераспространённых программTriangle [30], netgen [51] и gmsh [52].