Бураго Н.Г. Вычислительная механика (1185926), страница 41

Просмтор этого файла доступен только зарегистрированным пользователям. Но у нас супер быстрая регистрация: достаточно только электронной почты!

Текст из файла (страница 41)

Расчет подвижных границ разделаГраничне условия учитываются в дифференциальнойформулировке включением в исходные уравнения интеральныхчленов с дельта-функциями. [C.S.Peskin, 1977].. При численнойреализацииучаствующаявуравненияхдельта-функцияаппроксимируется обычной функцией-шапочкой, при этом границераздела приписывается конечная толщина., что позволяетобеспечить устойчивочть и гладкость решения в окрестностиграницы.20.7. Граничные дискретные маркерыВ методе граничных дискретных маркеров на границераздела сред вводится подвижная сетка поверхностных ячеек.

Расчетдвижения сред как и в подходе непрерывных маркеров проводитсяна эйлеровой сетке, покрывающей область движения сред.Впервые сетка поверхностных дискретных маркеров была введена вработе [Нох, 1964]. Подробный анализ возможностей этого подходавыполнен в работе [Unverdi, Tryggvason, 1992].Отслеживание изменений топологии фаз при разделении илислиянии фрагментов среды в двумерном случае еще удаетсяреализовать, хотя уже в этом случае логика алгоритмов генерацииповерхностных сеток для границ раздела становится чрезмерносложной. В трехмерном случае реализация переменной топологиинаталкивается на непреодолимые трудности реализации.20.8.

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

Принципиальныхтрудностей в определении этих параметров геометрии границраздела нет, но такое определение связано с анализом соседствамаркеров, что потребует очень большого объема вычислений.233Глава 20. Расчет подвижных границ разделаСовместное использование дискретных маркеров и маркерфункций позволяет улучшить описание границ раздела и учестьсложные граничные условия..20.9.

Метод гладких частицГидродинамический метод гладких частиц (SPH - Smoothedparticle hydrodynamics) предложен в работе [Gingold, Monaghan,1977] является лагранжевым методом частиц используемым дляописания гидродинамических явлений в областях решения,имеющих переменную связностьи подвижные границы.Подробное обсуждение метода можно найти в обзоре Монагана(Monaghan, 1992).

Чтобы дать представление о методе гладкихчастиц ниже приведен современный вариант расчетной схемыметода. Выкладки по выводу формул опущены, их можно найти вцитированных источниках.Рассматривается течение идеального газаdρ= −ρ∇ ⋅ vdtdv1= − ∇pdtρdUp= − ∇⋅vρdtp = ( γ − 1)ρU(1)(2)(3)(4)где использованы традиционные обозначения.В методе гладких частиц газ моделируется системойлагранжевых частиц, каждая из которых имеет свою массу mi ,скорость v i и внутреннюю энергию U i . Эти дискретные параметрыf i связаны со значениями соответствующих непрерывных искомыхфункций f следующими аппроксимационными формулами.f i ≈< f > (x i ) = ∫ f (x)W (x − xi , h)dx =∑ ∆fj∈ωii, j(5)d 1  −|x|2 / h2W (x, h) =  eh π ∆f i , j ≈ m jf ( xi )W ( xi − x j , h)ρ( xi )234(6)Глава 20.

Расчет подвижных границ раздела(7)где ωi - номера соседей частицы i, h – радиус сглаживания, d -число пространственных координат, W - гауссова нормированнаясферически симметричная функция ядра (возможны другиеопределения). обладающая следующими свойствамиW ( x, h ) = W ( − x, h )(8)∫ W (x, h)dx = 1(9)Таким образом, плотностьопределяется соотношениемρi =∑ m W (xj∈ωijiвокрестности− x j , h)точкиxi(10)Сверткой уравнений движения и уравнения энергии с ядромW исходные уравнения приводятся к дискретным уравнениямлагранжевых траекторий частиц и к уравнениям эволюциивнутренней энергии частицd 2 xi= − ∑ m j p*[Vi ,2j (hi )∇W (xi − x j , 2hi ) +2dtj∈ωi+Vi ,2j (h j )∇W (xi − x j , 2h j )](11)dU idx= − ∑ m j ([ pv]* − p* i )[Vi ,2j (hi )∇W (xi − x j , 2hi ) +dtdtj∈ωi+Vi ,2j (h j )∇W (xi − x j , 2h j )](12)где полагается, что радиус сглаживания h не очень сильно меняется**для соседних частиц, pij и [ pv ]ij - значения в середине ребра (i,j):si*, j , s – локальная координата на ребре (i,j).Поскольку слагаемые в правых частях дискретныхуравнений антисимметричны по номерам частиц, то присуммировании по всем частицам их вклады взаимно уничтожаются.Поэтому при отсутствии физических источников дискретныеуравнения обеспечивают сохранение массы, импульса и полнойэнергии для рассматриваемой системы частиц235Глава 20.

Расчет подвижных границ раздела2d  1  dxi ∑i mi dvi / dt = 0 и ∑i mi dt  2  dt  + U i  = 0Масса системы частиц сохраняется так как масса каждой частицы иих число неизменны.Подобно схеме Годунова для вычисления потоков в среднихточках ребер между соседними частицами i и j величины p* и[ pv ]* в эволюционных уравнениях для скорости и внутреннейэнергии определяются из решения задачи Римана о распадепроизвольного разрыва.Подобно схемам коррекции потоков для устойчивого расчетаразрывов вводятся монотонизирующие ограничители градиентовфизических величин.

Численные эксперименты показывают, чтомонотонизирующий ограничитель требуется только для градиентовскорости. Таким образом, ограничение по монотонности имеет вид ∂v   ∂v  ∂v   ∂v   =   = 0 при  ∂s   ∂s  < 0 i   j ∂s i  ∂s  jгде s – координата вдоль ребра (i,j). Численная схема должна иметьпервый порядок по пространству на скачке. Соответствующийпереключатель имеет вид ∂f   ∂f   =   = 0 при cshock ei , j ⋅ ( v j − v i ) > min(cs ,i , cs. j ) ∂s i  ∂s  jгде f = ρ, p, v и cshock является счетной постоянной, отвечающейчислу частиц в окрестности скачка, ei , j - единичный вектор вдольребра (i,j), cs ,i и cs , j - скорости звука в точках xi и x j .

Обычнополагается cshock = 3 .Радиус сглаживания определяется формулой1/ dm hi = η  *i  ρi где функция ρ* более гладкая, чем ρ при csmooth > 1 :ρ*i =∑ m W (xj∈ωijihi* = hi csmooth− x j , hi* ) ,236Глава 20. Расчет подвижных границ разделаЧисленные эксперименты показывают, что η ≈ 1 с csmooth = 2дают хорошие результаты. Эффективное число соседей околокаждой частицы зависит от отношения радиуса сглаживания ксреднему расстоянию до соседней частицы. Например, взаимнымивкладамичастицiиjможнопренебрегать,если| xi − x j |> 3min(hi , h j ) , так как входящий в функцию ядра множительexp(−32 ) ≈ 1.23410 − 4 мал.

Таким образом число соседей равно236ηCsmooth при d=1, 28η2Csmoothпри d=2 и 128η3Csmoothпри d=3.237Глава 21. Метод граничных элементовГлава 21. Метод граничных элементов21.1. Граничные интегральные уравненияРассмотрим метод граничных интегральных уравнений напримере смешанной краевой задачи для уравнения Лапласа∆ϕ = 0 внутри Ωс граничными условиямиϕ = f* или ∂ϕ / ∂n = g* на ∂ΩУравнение Лапласа имеет сингулярное (фундаментальное)решение 1/(4π r ( p, q )) . Здесь r ( p, q ) является расстоянием междупроизвольными точками p и q в области решения. Это решениеотвечает единичному возмущению (точечному источникуединичной интенсивности), приложенному в точке q :∆ϕ = δ ( r ( p, q ))гдеδ (r ( p, q ))соотношением-дельта-функцияДирака,определяемая∫ Φ( p)δ (r ( p, q))d Ω( p) = Φ(q)ΩУбедиться в справедливости этого решения для точек с r ( p, q ) ≠ 0можно непосредственной его подстановкой в уравнение Лапласа,которое в полярной системе координат для центральносимметричного и зависящего только от радиальной координатыr =| x p − xq | решения принимает вид∂ 2u * 2 ∂u *+− δ (r ) = 0∂r 2 r ∂rПри r ( p, q ) = 0 решение имеет особенность и трактуется в смыслеобобщенного решения в соответствии с определением Дирака.Умножая уравнение Лапласа на произвольную функцию w иинтегрируя его по частям, получаем238Глава 21.

Метод граничных элементов0 = ∫ ∇ 2ϕ wd Ω = ∫ ∇ ⋅ (∇ϕ w)d Ω − ∫ ∇ϕ ⋅∇wd Ω =ΩΩΩ= ∫ n ⋅ ∇ϕ wdS − ∫ ∇ ⋅ (ϕ∇w) d Ω + ∫ ϕ∇ 2 wd ΩΩSΩили∫ n ⋅∇ϕ wds − ∫ ϕ n ⋅ ∇wd Ω + ∫ ϕ∇ wd Ω = 02∂Ω∂ΩΩПодставляя сюда вместо w сингулярное решение получаемисходное интегральное тождество метода граничных интегральныхуравнений:ϕ ( p) =14π∂ 11∫  g (Q) r ( p, Q) − f (Q) ∂n  r ( p, Q)  ds(Q)∂Ωгде p ∈ Ω , Q ∈ ∂Ω , g (Q ) = n ⋅ ∇ϕ Q , f (Q ) = ϕ Q . Устремлениемp → P ∈ ∂Ω это тождество сводится к интегральному граничномууравнениюf ( P) =12π1∂  1    ds (Q )Q  r ( P, Q )  ∫  g (Q) r ( P, Q) − f (Q) ∂n∂Ωгде f ( P ) = lim ϕ ( p ) при p → P .

Это уравнение устанавливаетсвязь между f и g , а также обеспечивает соответствие обеих этихфункций одной и той же гармонической функции ϕ . В каждойточке границы задана одна из функций f или g , а другаяподлежит определению. После решения этого уравнения решениеϕ в области Ω определяется с помошью исходного интегральноготождества.21.2. Численная реализацияГраница представляется набором N граничных элементов(отрезков). Значения функций f и g ищутся в классе кусочнопостоянных функций, принимающих постоянные значения накаждом из отрезков.

Интегральное уравнение записывается длякаждого из отрезков. На каждом отрезке искомым является значение239Глава 21. Метод граничных элементоводной из функций f или g , в то время как другое заданограничными условиями.Основанная на этих предположениях запись интегральногоуравнения в точке Pi имеет вид1δ ij +2π 1∂  1  ds (Q )  f j = ∂n(Q)  r ( Pi , Q ) 2π∂Ω j∫1i∫ ( r ( P , Q) ds(Q)  g∂Ω jjгде ∂Ω j - граничный элемент, i, j = 1,..., N . Краткая форма записисистемы уравнений имеет видAij f j = Bij g jЗдесь матрицы Aij и Bij размера N × Nсодержат интегралы,указанные в записи интегрального уравнения выше.

На каждомграничном элементе одно из значений f i и gi задано в соответствиис граничными условиями задачи, а другое определяется выписаннойсистемой алгебраических уравнений. Искомая функция ϕ ( p ) влюбой внутренней точке p определяется затем при помощиподстаноски полученных граничных значений f i и gi в исходноеинтегральное тождествоϕ ( p) =14π∂ 11∫  g (Q) r ( p, Q) − f (Q) ∂n  r ( p, Q)  ds(Q)∂ΩПроизводные от решения определяются непосредственнымдифференцированием интегрального тождества (дифференцируютсяподинтегральные выражения, содержашие r ( p, Q ) ).Описанный способ решения носит название методаграничных интегральных уравнений (метода ГИУ) или методаграничных элементов (МГЭ).

Для большинства задач он хорошоработает при малом числе N даже для самых простыхаппроксимаций функций f и g . Повышение точности достигаетсяза счет увеличения числа элементов N и за счет примененияаппроксимаций более высокого порядка точности (линейных,квадратичных, кубических и так далее). Интегралы в выраженияхдля коэффициентов алгебраических уравнений метода определяютсячисленно с использованием квадратурных формул. Матрицы МГЭявляются полностью заполненными с диагональным преобладанием.С ростом числа элементов N их обусловленность ухудшается, что240Глава 21. Метод граничных элементовможет потребовать предобусловливания (умножения системыуравнений на приближенную обратную матрицу системыуравнений) для обеспечения устойчивости процесса решения.Ключевой особенностью метода ГИУ является тообстоятельство, что для решения задачи надо аппроксимироватьтолько граничную поверхность.

Характеристики

Список файлов книги