Н.А. Боголюбов - ОММ Методическое пособие (1133432), страница 11
Текст из файла (страница 11)
Таким11образом, решение uij2 на промежуточном слоево всех точках сетки.23) Пологая k = 0 , решаем методом прогонки при каждом1 ≤ i ≤ Nx − 1 системуuijk+1−τ2k+ 12uij= Λ1uijk+ 2 + Λ2uijk+1 + f (xi, yj, tk+ 12 ) при1граничных условиях:ui0 = μ(xi,0, tk+1),uiNy = μ(xi, ly, tk+1) .157Решение при i = 0 и i = Nx находится из граничных условий.Значит, найдено решение uij1 на слое k = 1 во всех точках сетки.4) Вычислив характеристики полученного решения, увеличиваем номерслоя на единицу k = k + 1 и повторяем пункты 2) и 3) , пока номерслоя не достигнет заданного значения или не будет выполнен другойкритерий окончания счета.Итак, метод прогонки дает выполнение условия: число операцийпропорционально числу узлов сетки, безусловная устойчивость проверяетсяспектральным критерием устойчивости.Данную схему нельзя обобщить на трехмерный случай, так схема будетусловно устойчива и не будет экономичной.
В трехмерном случае мы будемиспользовать локально одномерную разностную схему, которую мырассмотрим чуть позже.Теперь рассмотрим численный метод решения задачи Дирихле для уравненияПуассона.Пусть требуется решить разностную задачу Дирихле:Λ1uij + Λ2uij = − fij, 0 ≤ i ≤ Nx − 1, 0 ≤ j ≤ Ny − 1 при граничных условиях:ui0 = μ(xi,0), 0 ≤ i ≤ Nx,uiNy = μ(xi, ly), 0 ≤ i ≤ Nx,u0j = μ(0, yj), 0 ≤ j ≤ Ny − 1,uNx j = μ(lx, yj), 0 ≤ j ≤ Ny − 1 .158Поставленная задача аппроксимирует дифференциальную задачу:L1u + L2u = − f (x, y),{uΓ= μ(x, y), (x, y) ∈ Γ .
Для вычисления решений многих стационарных задач математическойфизики, описывающих равновесные состояния, которые рассматриваются ,как результат установления развивающегося во времени процесса, расчеткоторого оказывается проще, чем прямой расчет равновесного состояния.Рассмотрим вспомогательную задачу - нестационарную задачу ораспространении тепла:∂u∂tu= L1u + L2u + f (x, y),Γ= μ(x, y), u(x, y,0) = u0(x, y),где f (x, y) имеет прежний смысл , u0(x, y) выбирается произвольно.Так как источники тепла f (x, y) и температура на границе μ(x, y) независят от времени .
Следовательно логично предположить, что и решениеu(x, y, t) с течением времени будет меняться всё медленнее и распределениетемператур в пределе при t → ∞ превратится в равновесноераспределение температуры u(x, y) , описываемое исходной задачей.Задача решается до того момента, пока ее решение не перестанет менятьсяв пределах интересующей нас точности.159Замечание:Если известны границы спектра разностных операторов −Λ1; − Λ2πhyπhxπhx444222δ1 = c1 2 sin, Δ1 = c2 2 cos, δ2 = d1 2 sin,hx2lxhx2lxhy2lyπhy42Δ2 = d2 2 cos, и δ = min(δ1, δ2), Δ = max(Δ1, Δ2) , то в качествеhy2ly2оптимального τ берется τopt ≈.δΔВ модельной задаче τopt ≈h2sin(πh).
Количество иттераций ,необходимое для уменьшения первоначальной нормы погрешности в ξ раз:p≈N 1ln.2π ξДанное замечание является лишь ориентировочным и их следуетиспытать, иными словами, проверить: скорость сходимости при τ большихи меньших τopt , указанного выше.Теперь перейдем к рассмотрению локально одномерной разностной схемы.Итак, трехмерный случай ut = uxx + uyy + uzz .Рассмотрим его на примере решения уравнениятрех пространственных переменных.∂u= Δu + f в случае∂tДля данного уравнения безусловную устойчивость ( дляэкономичности разностной схемы ) обеспечивает локально одномернаяразностная схема.160Если ее решать методом прогонки , то будет выполнено второе условиеэкономичности разностной схемы: число арифметических операций будетпропорционально числу узлов сетки.В трехмерном случае, схема переменных направлений непременимапотому что она будет неявной только по одному направлению, и уже явнойпо двум остальным. Таким образом нарушается условие устойчивостиэкономичных разностных схем!Заменим дифференциальный оператор рассматриваемого уравнения егоразностным аналогом:33∂2 uΔu =→ Λv =Λv .∑ ∂xi2∑ ii=1i=1То есть: оператор 3 пространственных переменных →cумма операторовпроизводные и разностипо одной переменной.Основная идея: переход с s -го на s + 1 cлой в три этапа, на каждом изкоторых учитывается оператор по одной пространственной переменной.Введем два массива промежуточных значений u и u :u−uτu−uτu^ − uτ= Λ1u + f= Λ2 u= Λ3u^каждое из уравнений не аппроксимируетисходное дифференциальное, но хорошорешается методом прогонки.Просуммируем их:u^ − u= Λ1u + Λ2 u + Λ3u^ + f .τ161Полученное уравнение аппроксимирует исходное дифференциальноеO (τ + hx2 + hy2 + hz2) .уравнение с порядком апроксимации _Схема абсолютно устойчива и прогонка дает нам выполнение условия:число арифметических операций пропорционально числу узлов сетки.Следовательно локально - одномерная схема является экономичной !Схема чисто неявная.
Для повышения точности , можноиспользовать схему с весами δ = 0,5 .Пример практического задания по схеме переменных направлений.Условие : используя метод переменных направлений, решить краевую задачу:∂u∂t∂u∂x∂u∂yu= Δu + cos(x) ⋅ cos(y) ⋅ e −t, 0 < x < π, 0 < y < π, t > 0, (1)x=0y=0t=0=∂u∂x=∂u∂yx=πy=π= 0, (2)= 0, (3)= cos(x) ⋅ cos(y) .
(4)Решение:Будем решать задачу численно. Введем разностную сетку в областиG ⊗ [0, T ] :G = {(x, y); 0 ≤ x ≤ π, 0 ≤ y ≤ π},162wh =11xi = ihx, i = 0, . . . , Nx . hx = , yj = jhj, j = 0, . . . , Ny, hy =,NxNy }{wτ =1tk = kτ; k = 0, . . . , Nt, τ =,Nt }{Nx - число узлов вдоль оси x , τ - шаг по времени.Сеточная функция: vi,k j = u(xi, yj, tk) .Введем разностную аппроксимацию оператора Лапласа:Λv = Λx v + Λyv;Λx v =vi−1, j − 2vi, j + vi+1, jhx2,Λyv =vi, j−1 − 2vi, j + vi, j+1hy2.Значение искомой функции известно из начальных условий.Решать поставленную задачу следует методом переменныхнаправлений. Введем обозначения :1^vi,k j ≡ vi, j , vi,k+j 2 ≡ vi, j , vi,k+1j ≡ vi, j .В дальнейшем для сокращения записи неизменяющиеся индексы мыбудем опускать.Начальное условие получаем из (4) с помощью заменыx → xi = ihx .Конкретно в нашей задаче мы имеем граничные условия Неймана.
Ихможно аппроксимировать с помощью односторонней разностнойпроизводной:v1, j − v0, j vNx, j − vNx−1, j== 0, j = 0, . . . , Ny,hxhx163vi,1 − vi,0hy=vi,Ny − vi,Ny−1hy= 0, i = 0, . . . , Nx .Однако, порядок аппроксимации в этом случае будет лишь O(hx) иO(hy) .Для рассматриваемой задачи, явная и неявная схемы имеют одинпорядок точности.
При использовании явной схемы число действий Qyavна k + 1 - ом слое пропорционально числу узлов сетки :Qneyav = O(Nx Ny) .Неявная схема для определения сеточной функции v k+1 разрешаетсистему линейных уравнений , число которых пропорционально числуузлов сетки, то есть число действий возрастает до Qneyav = O (Nx Ny) .()2Кроме того, она абсолютно устойчива. Используемая схемапеременных направлений , как мы уже ранее говорили, сочетает в себедостоинства явной и неявной схемы, то есть Q = O(Nx Ny) и абсолютноустойчива .
Следовательно, она экономична.Запишем разностную аппроксимацию для данного случая:1v k+ 2 − v k11= Λx v k+ 2 + Λyv k + f k+ 2 ,0,5τfk+ 121v k+1 − v k+ 211= Λx v k+ 2 + Λyv k+1 + f k+ 2 ,0,5τ= cos(x) ⋅ cos(y) ⋅ e−tk+ 12.1Таким образом, мы вводим промежуточный дробный слой k + , в2результате переход со слоя k на слой k + 1 будет осуществляться в два164этапа: сначала решается первое разностное уравнение , неявное понаправлению x и явное по направлению y , а затем решается второеразностное направление явное по направлению x и неявное понаправлению y .Запишем явный вид разностных операторов на примере первогополуслоя:vk+ 12k+ 12vi−1k−v=0,5τ−k+ 122vihx2+k+ 12vi+1+kkvj−1− 2vjk + vj+1hy21+ f k+ 2 .Собираем коэффициенты при неизвестных v , и проделываяаналогичное для второго полуслоя , мы приходим к следующей задаче:0,5τ k+ 12vhx2 i−1, j− [1 +k+ 12v0, j=k+ 12v1, jk+ 12Fi, j=0,5=τk+ 12vhx2 ] i, jk+ 12vNx−1, j=+0,5τ k+ 12vhx2 i+1, jk+ 12vNx, j= 0,kkv+vi,j−1i,j+1] + [1 −h y2 [0,5τ k+ 12vh y2 i, j−1− 1+[τk+ 12vh y2 ] i, j+=k+ 12vhx2 [ i−1, j0,5+k+ 12vi+1, j]τkvi,h y2 ] j0,5τ k+ 12vh y2 i, j+1k+1k+1k+1k+1vi,0= vi,1= vi,N=v= 0,−1i,NyyFi,k+1j=−k+ 12Fi, j ,+ [1 −1+ 0,5τf k+ 2 ;= − Fi,k+1j ,τk+ 12vhx2 ] i, j1+ 0,5τf k+ 2 .Рассмотрим первый полуслой, будем решать полученную задачуметодом прогонки :пусть имеется уравнение заданного вида с условиями:Ai yi−1 − Ci yi + Bi yi+1 = − Fi,y0 = χ1y1 + μ1,yN = χ2 yN−1 + μ2 .Ci > Ai + Bi, 0 ≤ χα ≤ 1, α = 1,2, i = 1, .
. . , N .165166167168169170171ГЛАВА 7Иттерационные методы.Схемы бегущего счета.При изучении высокотемпературных процессов необходимо учитыватьзависимость коэффициентов теплоемкости и теплопроводности оттемпературы.Процесс распространения тепла описывается квазилинейным уравнениемтеплопроводности:∂u∂∂uc(u)=k(u)+ f (t), k(u) > 0, c(u) > 0 . ∂t∂t (∂t )Cделаем замену:uv = c(u)du ,∫0и получим:∂u∂∂u=k(u)+ f (t) .()∂t∂t∂tИспользование явных разностных схем не выгодно, например, еслиk(u) - быстроменяющаяся функция ( к примеру, степенная ) .k(u) , c(u) зависят от пространственных переменных .
Тогдаустойчивость требует более меньшего шага по времени:h2.τ≤2maxk(u)172k(u) и c(u) могут быть разрывными.Будем использовать неявную разностную схему:^vm − vm1= 2 am+1(^v)(^vm+1 − ^vm) − am(^v)(^vm − ^vm−1) + fm .]h [τv - значение сеточной функции на s - ом слое ,^v - значение сеточной функции на s + 1- ом слое,vi−1 + viai(y) = k.( 2 )Рассматриваемая схема является нелинейной , следовательно, методпрогонки не применим !Вместо этого, мы воспользуемся заявленными иттерационнымиметодами.Основная идея: при переходе с s - го на s + 1 - ый слой , мы совершаемнесколько повторений расчета , при каждом фиксированном am(y) наоснове результатов предыдущей иттерации.При этом на каждой иттерации , решаемое уравнение мы сводим клинейному, которые в дальнейшем мы уже умеем решать методом прогонки.Нулевое приближение ^v[0] на s + 1 - ом слое является значениемсеточной функции v на s - ом слое.173Первое приближение:111]− ^v[m ]) − am(v)(^v[m ] − ^v[m−1)] + fm ..................................................................................................^v[m1] − vmτ=1]1am+1(v)(^v[m+12h [k - ое приближение:^v[mk] − vm1k]kkk]= 2 am+1(^v[k−1])(^v[m+1− ^v[m ]) − am(^v[k−1])(^v[m ] − ^v[m−1+ fm .)[]hτДалее, полученные линейные уравнения решаем методом прогонки.Условие окончания иттерационного процесса - условие малостиизменения нормы:^v[k+1] − ^v[k] = max0<m<N^v[mk+1] − ^v[mk] < ε .Теперь составим схемы бегущего счета для уравнения переноса.Рассмотрим следующую начально - краевую задачу:∂u∂t+C∂u∂x= f,ux=0= μ(t),ut=0= φ(x) .Так как коэффициент не зависит от u , оно нелинейное ! Данноеуравнение называется уравнением переноса .Если C(u) , то уравнение является квазилинейным.174Введем равномерную сетку по пространственной и временнойкоординатам:w h : {xm = hm, m = 0,1, .
. . }, w τ : {tn = τn, n = 0,1, . . . },w hτ = w h × w τ .Составим 4 разностные аппроксимации рассматриваемого уравненияпереноса:1)vm,n+1 − vm,nτ+Cvm,n − vm−1,nh= fm,n .n+1m−12)vm,n+1 − vm,nτ+Cvm+1,n+1 − vm,n+1hmn= fm,n+1 .n+1m3)vm,n+1 − vm,nτ+Cvm,n+1 − vm−1,n+1hm+1n= fm,n+1 .n+1m−1175mn4)+vm+1,n+1 − vm+1,nτvm,n+1 − vm,nτ+C+Cvm+1,n − vm,nhvm+1,n+1 − vm,n+1h+n+1= fm+1,n + fm,n+1 .mm+1На примере построенной разностной схемы 1) , проверимнеобходимые ( спектральный критерий Неймана ) и достаточные ( критерийКуранта ) . А также рассмотрим геометрический критерий устойчивостирассматриваемой разностной схемы.Итак, поставим задачу:vm,n+1 − vm,nτ+Cvm,n − vm−1,nh= 0,vm,0 = e −iwm .Необходимые условия устойчивости:Будем искать решение поставленной задачи в следующем виде:vm,n = λne −iwm .Подставляем общее решение в разностное уравнение и получаем:λ−1(1 − e+Cτh−iw)= 0.Выразим λ :Cτλ =1−(1 − cosw + isinw) = 0 .h176nλ21 Cτ 2 wCτ=1−sin1−= 0.4 h2(h )>0λ ≤≥1Cτ- необходимое условие устойчивости выполнено!hДостаточные условия устойчивости:vm,n+1 = (1 −vm,0 = φm .Cτvh ) m,n+Cτvh m−1,n+ τ( fm,n),Можарантно оценим левую часть рассматриваемого уравнения,Cτучитывая 1 −≥1 :(h )CτCτvm,n+1 ≤ 1 −v+v+ τ fm,n ≤(h ) m,nh m−1,nCτCτ≤ 1−v +v +τ f .(h ) nh nТак как :f = max fm,n ,m,nvn = max vm,n .mТо:vn+1 ≤ vn + τ f- достаточные условия Куранта.Таким образом:vn ≤ v0 + nτ f ≤ φ + T f ,где T - величина интервала времени , на котором мы ищем решение.Данная оценка является устойчивостью решения по начальнымданным и правой части.
![](https://s.studizba.com/z.php?f=/templates/si/i/noavatar.png&w=100&h=100&t=1"){kind=avatar}
