Н.А. Боголюбов - ОММ Методическое пособие (1133432), страница 10

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

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

Точка x =не2совпадает с точкой разбиения сетки xm ⇒ значение коэффициента k(xm)определено.Обозначим точку, лежащую левее от точки x =Тогда xn+1 - точка, лежащая правее от точки x =N1как xn , n =.221.2Подберем разностную схему с аппроксимацией нашего разностногоуравнения _O (h 2) :kmvm+1 − 2vm + vm−1h2+km+1 − km−1 vm+1 − vm−12h2h= 0.Перепишем разностное уравнение в следующем виде:1111km−1 + km − km+1 vm−1 + (−2km)vm + − km−1 + km + km+1 vm+1 = 0 . ( * )(4)( 4)44Будем искать численное решение разностного уравнения в видекусочно линейной функции:v(xm) =1 − a xm, 0 < xm < 12 ,b(1 − xm), 12 < xm < 1 .144По физическому смыслу решение должно быть непрерывным:limh→0vn = 1 −a1= limh→0vn+1 = b 1 −⇒a+b =2.(22)Теперь рассмотрим уравнение ( * ) для m = n, m = n + 1 .

Будемиспользовать следующие равенства:xn−1 = xn + ah,⇒{xn+2 = xn+1 + bh .13v+ah−2v+v = 0,()n2 n2 n+157v−6v+v − bh) = 0 .n+12 n2 ( n+1Приведем подобные слагаемые:3vn+1 − 3vn + ah = 0,21⇒a=b .{5vn+1 − 5vn + 7bh = 0 .5a + b = 2 ⇒ считаем коэффициенты ⇒ a =v= limh→0vn+1 =x=0,5521., b=13135511−0,5=≠= u(0,5) .()1326 4Cледовательно, решение u(x) рассматриваемой дифференциальной задачиотличается от численного решения:v(x) = limh→0v(xm) .Причиной расхождения численного и аналитического решений являетсяиспользование неконсервативной разностной схемы, которая не учитывает законысохранения!145В данной задаче нужно было учесть закон сохранения теплового потока:q = − k(x)du= const .dxОдним из методов построения консервативных разностных схемявляется интегро - интерполяционный метод или метод баланса , к изложениюкоторого мы сейчас перейдем.Рассмотрим заявленный метод для решения уравнения :∂u∂∂u=k(x)+ f (x, t) .∂t∂x (∂x )k , f - быстроменяющиеся, разрывные функции.Как и прежде, будем искать значение сеточной функции vm на s - омслое в узлах сетки xm , а на s + 1 - ом слое в узлах сетки xm через ^vm .Введем промежуточные точки:xm+0,5 = h + xm−0,5 = h(m + 0,5) .Будем рассматривать в этих точках аналоги потока тепла Wm+0,5 .Запишем уравнение баланса тепла на отрезке [xm−0,5; xm+0,5 ] .∂uПоток тепла будет равен W = − k.∂xПроинтегрируем выражение∂uW=−на отрезке [xm−1, xm] :∂xk146xmum−1 − um =Wdx .∫ kxm−1Тогда мы можем представить исходное рассматриваемое уравнение вследующем виде:∂W∂u= f (x, t) −.∂x∂tПроинтегрируем записанное выражение и получим:xm+0,5Wm+0,5 − Wm−0,5 =∫ (f−xm−0,5∂udx .)∂t∂uи W - непрерывны и мало изменяются на отрезке [xm−1, xm] .∂tДля нахождения приближенного значения интегралов,можно вынести из под интегралов в средней точке участкаинтегрирования.Получаем следующие аппроксимации:xm+1^vm+1 − ^vm ≈ Wm+0,5dx,∫ kxm∂uWm+0,5 − Wm−0,5 ≈ − h∂t147xm+0,5+x=xm∫xm−0,5dx .∂uи W∂tТогда используя последние соотношения, мы можем положить:Wm+0,5 = am+0,5^vm+1 − ^vmh;am+0,5 =hxm+1dx∫ k;xmWm+0,5 − Wm−0,5h^vm − vm1=−+ Fm , где Fm =4τxm+0,5∫fd x .xm−0,5Исключая из полученных соотношений W , мы получим разностноеуравнение в следующем виде:^vm+1 − ^vm^vm − vm 1^vm − ^vm−1=am+0,5− am−0,5+ Fm .τh(hh)am+0,5 и Fm определяются функциями k(x) и f (x, t) .Когда функции k(x) и f (x, t) непрерывны ⇒ вычисление интеграловможно заменить разностной аппроксимацией:xm+11dx1≈,h ∫ k(x) km−0,5xmxm+11dx111≈+.h ∫ k(x) 2 ( km+1 km )xmТак как мы учли закон сохранения тепла , полученная схема являетсяконсервативной .

Для каждой ячейки сетки, схема неявная : абсолютноустойчива и имеет порядок аппроксимации _O (τ + h 2) .148Теперь перейдем от консервативных разностных схем к экономичнымразностным схемам.Определение: экономичной разностной схемой называется схема, котораяявляется:1) безусловно устойчивой ;2) при расчетах: число операций пропорционально числу узлов сетки( основное свойство ) .Экономичные разностные схемы сочетают в себе преимущества явных и неявныхсхем.Основная идея использования экономичных разностных схемзаключается в сведении многомерных задач к цепочке одномерных.Явная схема для уравнения теплопроводности является условноустойчивой и не является экономичной.Рассмотрим уравнение ut = uxx + uyy .

Будем решать его численно.Разностная схема будет безусловно устойчива. Чтобы число узлов былопропорционально числу действий , мы будем использовать схему переменныхнаправлений, которая была впервые предложена в 1955 году ПисменомРэкфордом.149а) Постановка задачи:Рассмотрим двумерное уравнение теплопроводности:∂u∂t= Lu + f (x, y, t), (x, y) ∈ G, t ∈ (0; T ],Lu =∂∂upx,y()∂x (∂x )+∂∂uqx,y,()∂y (∂y )p(x, y) , q(x, y) - достаточно гладкие функции в рассматриваемойобласти.

Такие , что 0 < c1 ≤ p(x, y) ≤ c2, 0 < d1 ≤ q(x, y) ≤ d2,c1 , c2 , d1 , d2 - постоянные.Обозначим A = max(c2, d2) . Пусть G = {0 ≤ x ≤ lx, 0 ≤ y ≤ ly} прямоугольник со сторонами lx и ly , Γ - его граница, G = G + Γ .Необходимо найти решения двумерного уравнениятеплопроводности, удовлетворяющего следующим условиям:u= μ(x, y, t), (x, y) ∈ Γ,{u(x, y,0) = u0(x, y) .ΓЗададим Nx и Ny - число разбиений по x и y . В области G построимlylxравномерную сетку whx hy с шагами hx =и hy =:NxNywhx hy ={(}xi, yj), xi = ihx, jhj, 0 ≤ i ≤ Nx, 0 ≤ j ≤ Ny Введем обозначения:∂∂u∂∂uL1u =p(x, y), Lu=q(x, y).∂x (∂x ) 2∂y (∂y )Имеем уравнение :Lu = L1u + L2u .150.Λ2 :Теперь заменим операторы L1 и L2 их разностными аналогами Λ1 иΛ1uij = pi+ 12 jΛ2uij = qij+ 12ui+1j − uijhx2uij+1 − uijhy2− pi− 12 j− qij− 12uij − ui−1jhx2uij − uij−1hy2,.Тогда:Λuij = Λ1uij + Λ2uij, 1 ≤ i ≤ Nx − 1, 1 ≤ j ≤ Ny − 1 .Пусть M - количество разбиений по t , длина частичного разбиенияTи точки разбиения tk = kτ .τ=MОбозначим uijk ≈ u(xi, yj, tk) .Теперь составим разностную схему для рассматриваемой задачи:б) Двухслойная схема с весами :Будем аппроксимировать нашу задачу следующей разностной схемой:uijk+1 − uijkτ= Λ(δuijk+1 + (1 − δ)uijk) + f (xi, yj, tk),1 ≤ i ≤ Nx − 1, 1 ≤ j ≤ Ny − 1, 0 ≤ k ≤ M − 1 .151ui0k+1 = μ(xi,0, tk+1), 0 ≤ i ≤ Nx,k+1uiN= μ(xi, ly, tk+1), 0 ≤ i ≤ Nx,yu0jk+1 = μ(0, yj, tk+1), 1 ≤ j ≤ Ny − 1,uNk+1= μ(lx, yj, tk+1), 1 ≤ j ≤ Ny − 1 .xjПри достаточной гладкости функций u(x, y, t), μ(x, y, t), p(x, y), q(x, y)полученная разностная схема аппроксимиует исходную дифференциальнуюс порядком аппроксимации O( h2+ τ) .Решение при k = 0 находится из условия:uij0 = u0(xi, yj), 0 ≤ i ≤ Nx; 0 ≤ j ≤ Ny .Рассмотрим два варианта значений параметра δ :1) при δ = 0 мы получаем явную разностную схему и решение вовнутренних точках вычисляется по формуле :uijk+1 = uijk + τ(Λ1uijk + Λ2uijk) + τf (xi, yj, tk),i = 1 ≤ j ≤ Nx − 1, j = 1 ≤ j ≤ Ny − 1, k = 0 ≤ k ≤ M − 1 .Ah 2Схема условно устойчива при τ ≤и общее число операций при4переходе со слоя на слой пропорционально числу узлов сетки h , то есть1O 2 , следовательно, схема не является экономичной.(h )1522) При δ = 1 , мы получаем неявную разностную схему.

Онаустойчива при любых пространственных и временных шагах h и τ .Для определения uijk+1 мы получаем на каждом шаге линейнуюсистему:uijk+1 − τ(Λ1uijk+1 + Λ2uijk+1) = uijk + τf (xi, yj, tk+1),i = 1 ≤ j ≤ Nx − 1, j = 1 ≤ j ≤ Ny − 1, k = 0 ≤ k ≤ M − 1 .Матрица данной системы пятидиагональная и решать систему можно методомматричной прогонки или методом исключения Гаусса , который при учетеспециального вида матрицы требует O(Nx2 Ny2) действий , то есть схема являетсяэкономичной.Теперь перейдем к заявленной схеме переменных направлений. Так какона является экономичной, она сочетает в себе лучшие качества явнойсхемы ( экономичность) и неявной ( устойчивость) .Помимо основных значений сеточной функции на слоях uijk и uijk+1 ,1вводится промежуточный дробный слой uijk+ 2 , которое можноτрассматривать , как значение при t = tk+ 12 = tk + .2Тогда решение задачи сводится к решению двух систем стрехдиагональными матрицами следующего вида:k+ 12uij− uijkτ2= Λ1uijk+ 2 + Λ2uijk f (xi, yj, tk+ 12 ),11 ≤ i ≤ Nx − 1, 1 ≤ j ≤ Ny − 1, 0 ≤ k ≤ M − 1,k+ 12uijk+1 − uijτ2= Λ1uijk+ 2 + Λ2uijk+1 f (xi, yj, tk+ 12 ),11 ≤ i ≤ Nx − 1, 1 ≤ j ≤ Ny − 1, 0 ≤ k ≤ M − 1 .153В граничных узлах решение должно принимать заданные условия1τ и t = (k + 1)τ , поэтому каждому изисходной задачи при t = k +()2наших разностных уравнений добавятся дополнительные условия:1uijk+ 2 − uijkτ2=k+ 12Λ1uijk+ 12ui0k+ 12uiNyk+ 12u0y+ Λ2uijk f (xi, yj, tk+ 12 ) := μ(xi,0, tk+ 12 ), 0 ≤ i ≤ Nx,= μ(xi, ly, tk+ 12 ), 0 ≤ i ≤ Nx,= μ(0, yj, tk+ 12 ), 1 ≤ j ≤ Ny − 1,uNk+x y2 = μ(lx, yj, tk+ 12 ), 1 ≤ j ≤ Ny − 1 .1схема неявная по направлению x и явная по направлению y .uijk+1−τ2k+ 12uij=k+ 12Λ1uij+ Λ2uijk+1 f (xi, yj, tk+ 12 ) :ui0k+1 = μ(xi,0, tk+1), 0 ≤ i ≤ Nx,k+1uiN= μ(xi, ly, tk+1), 0 ≤ i ≤ Nx,yk+1u0y= μ(0, yj, tk+1), 1 ≤ j ≤ Ny − 1,uNk+1= μ(lx, yj, tk+1), 1 ≤ j ≤ Ny − 1 .xyсхема явная по направлению x и неявная по направлению y..154Поэтому , разбив многомерную задачу ( для которой изначально методпрогонки не применим ) на две одномерные задачи, мы можем поотдельности решить их методом прогонки .При достаточной гладкости функций u(x, y, t) , p(x, y) , q(x, y)разностная схема аппроксимирует исходную задачу с порядкомO( h2+ τ 2) .Перепишем нашу разностную схему в следующем виде:k+ 12Aijui−1j12G k+= − uijk −ij−k+ 12Bijuij+k+ 12Cijui+1j=k+ 12G ij ,τΛ2uijk + f (xi, yj, tk+ 12 ) , 1 ≤ i ≤ Nx − 1, 1 ≤ j ≤ Ny − 1, 0 ≤ k ≤ M − 1 .)2(В граничных узлах решение находится из разностной схемы :k+ 12uij− uijkτ2=k+ 12Λ1uij+ Λ2uijk + f (xi, yj, tk+ 12 ) .Прямая и обратная прогонки осуществляются вдоль строк:y1u k ↦ u k+ 2lylx155xТеперь возьмем второе разностное уравнение:uijk+1−τ2k+ 12uij= Λ1uijk+ 2 + Λ2uijk+1 + f (xi, yj, tk+ 12 ) ,1и перепишем его в следующем виде:k+1k+1Aijuij−1− Bijuijk+1 + Cijuij+1= Gijuijk+1,где1Gijuijk+1 = − uijk+ 2 −1τΛ1uijk+ 2 + f (xi, yi, tk+ 12 ) , 1 ≤ i ≤ Nx − 1; 1 ≤ j ≤ Ny − 1; 0 ≤ k ≤ M − 1 .)2(В граничных точках решение находится аналогичным образом,прогонки осуществляются вдоль столбцов:y1lyu k+ 2 ↦ u k+1lxxВ простейшей задаче (p(x, y) ≡ 1, q(x, y) ≡ 1, h x = hy = h) ,Aij = Aij =τ,22hBij = Bij =ττ+1,C=C=.ijij22h2h156Рассмотрим алгоритм схемы переменных направлений:1) из начального условия , мы получаем решение при k = 0 во всехточках сетки:uij0 = u0(xi, yj), 0 ≤ i ≤ Nx, 0 ≤ j ≤ Ny .2) Пусть k = 0 , решаем методом прогонки при каждом 1 ≤ j ≤ Ny − 11cистемуuijk+ 2 − uijkτ2=k+ 12Λ1uij+ Λ2uijk + f (xi, yj, tk+ 12 ) при граничныхусловиях:u0j = μ(0, yj, tk+ 12 ),uNx j = μ(lx, yj, tk+ 12 ) .Решение при j = 0 и j = Ny находится из граничных условий.

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

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