Дифференциальные уравнения с частными производными и краевыми условиями (краевые задачи)

 

Для уравнений эллиптического типа

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

 

, (42)

 

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

1. Первая краевая задача. На контуре Г, ограничивающем область G (рис. 4), задана непрерывная функция . Требуется найти функцию , удовлетворяющую внутри G уравнению (42) и принимающую на границе заданные значения , т.е. должны быть выполнены условия при ; при .

 

 

Рис. 4. Область исследования

 

2. Вторая краевая задача. На контуре Г, ограничивающем область G (рис. 4), задана непрерывная функция . Требуется найти функцию , удовлетворяющую внутри G уравнению (42), нормальная производная которой на Г принимает заданные значения , т.е. требуется чтобы при ; при .

3. Третья краевая задача. На контуре Г, ограничивающем область G (рис. 4), задана непрерывная функция . Требуется найти функцию , такую, чтобы при ; при , где .

Третья краевая задача может быть рассмотрена как общая. Действительно, при и получаем первую краевую задачу, а при и получаем вторую краевую задача. Если область G ограничена, то соответствующая краевая задача называется внутренней, а в противном случае – внешней. Для уравнения Лапласа первая краевая задача называется задачей Дирихле, вторая – задачей Неймана и третья – смешанной краевой задачей.

 

Уравнения Лапласа в конечных разностях

 

Для получения конечно-разностного уравнения, соответствующего уравнению Лапласа (частный случай уравнения эллиптического типа)

 

, (43)

 

достаточно, выбрав шаг h>0, заменить производные отношениями конечных разностей по формулам:

 

и .

 

Тогда после преобразования будем иметь

 

. (44)

 

Применение метода сеток для решения задачи Дирихле

 

Идея метода сеток (метода конечных разностей) для численного решения краевых задач для двумерных дифференциальных уравнений заключается:

1) в плоскости G, в которой разыскивается решение, строится сеточная область , состоящее из одинаковых ячеек и приближающая данную область G;

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

3) на основании граничных условий устанавливаются значения искомого решения в граничных узлах области .

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

Сеточная область может состоять из квадратных, прямоугольных, треугольных и других клеток. От выбора основного размера клетки h зависит величина погрешности. Обычно задача решается сначала при большом значении h (малом числе клеток), а затем после грубого решения переходят к более мелкой сетке или во всей рассматриваемой области, или в какой-нибудь ее перспективной части.

Покажем применение метода сеток для построения решения задачи Дирихле

 

при и при , (45)

 

где – заданная непрерывная функция.

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

Выбрав шаг h, построим квадратную сетку , с таким расчетом, чтобы узлы (, ) сетки или принадлежали области G, или отстояли от ее границы Г на расстояние меньше, чем h.

Узлы сетки называются соседними, если они удалены друг от друга в направлении оси Оx или Oy на расстояние, равное шагу сетки h. Узел называется внутренним, если он принадлежит области G, а все четыре соседних с ним узла – множеству , в противном случае он называется граничным (например, узлы и сетки ). На рис. 5 внутренние узлы обозначены светлыми кружками, а граничные темными кружками и треугольниками.

 

 

Рис. 5. Узлы сетки

Граничный узел сетки называется узлом первого рода, если он имеет соседний внутренний узел этой сетки (например, узел на рис. 5); в противном случае граничный узел называется узлом второго рода (узел на рис. 5). Внутренние узлы и граничные узлы первого рода сетки называются расчетными точками. Граничные узлы второго рода не входят в вычисление и могут быть изъяты из сетки (на рис. 5 граничные узлы второго рода обозначены темными треугольниками).

Относительно сетки предположим, что множество ее расчетных точек «связное», т. е. любые две расчетные точки можно соединить цепочкой узлов, каждые два смежных элемента которой являются соседними узлами. Кроме того, будем считать многоугольную сеточную область выбранной так, чтобы ее геометрическая граница , возможно ближе примыкала к границе Г области G. Заметим, что узловые точки контура могут лежать как внутри, так и вне области G.

Значение искомой функции и=и(х,у) в точках обозначим через . Следуя общей схеме, для каждой внутренней точки сетки заменяем дифференциальное уравнение (45) конечно-разностным уравнением

 

, (46)

 

где – расчетные точки.

В граничных узлах первого рода полагаем ,

где B – ближайшая к точка границы Г.

Система (46) является неоднородной системой, причем число неизвестных (т.е. число внутренних узлов сетки) равно числу уравнений. Решив систему (46), получим приближенные значения искомой функции и=и(х,у) в узлах сеточной области . Тем самым будет найдено приближенной решение задачи Дирихле для области .

Если число узлов сетки велико, а также с учетом криволинейности области G для решения системы (46) прибегают к итерационным методам, с одновременным исправлением граничных значений. Так согласно процессу усреднения Либмана, выбрав начальное приближение , последовательные приближения для внутренних узлов сетки определяем по формуле

 

(k=1,2,…). (47)

 

 

Для уравнений параболического типа

 

В качестве примера уравнения параболического типа рассмотрим уравнение теплопроводности для однородного стержня 0 ≤ x ≤ l

 

, (48)

 

где u=u(х,t) – температура и t – время. В дальнейшем для простоты будем полагать а= 1 (к такому случаю всегда можно прийти путем введения нового времени τ=a²t).

Итак, рассмотрим уравнение

 

. (49)

 

Пусть в начальный момент времени t=0 задано распределение температуры u(x,0)=f(x) и законы изменения температуры в зависимости от времени (тепловые режимы) на концах стержня х=0 и x=l: u(0,t)=φ(t), u(l,t)=ψ(t).

Требуется найти распределение температуры u=u(х,t) вдоль стержня в любой момент времени t. Решим эту смешанную задачу методом сеток. Для этого рассмотрим пространственно-временную систему координат (х,t) (рис. 6).

 

Рис. 6. Прямоугольная сетка

 

В полуполосе t ≥0, 0 ≤ x ≤ l построим прямоугольную сетку x=ih (i=0, 1,…, n), t=jk (j=0, 1, 2,…), где h=l/n (n – целое) – шаг вдоль оси Ох и k²=σh² – шаг вдоль оси Ot (σ – постоянная), вообще говоря, различны.

Величина σ будет выбрана ниже. Введя обозначения xi =ih, t_j =jk, u_ij=u(x_i, t_j) и заменяя уравнение (49) конечно-разностным уравнением, получим:

. (50)

 

Отсюда

 

u_i,j+1= σ u_i-1,j+(1-2σ)u_ij+ σ u_i+1,j. (51)

 

Из рассмотрения формулы (51) ясно, что зная значения функции u(х,t) в точках j -го слоя t=jk, с помощью этой формулы можно вычислить значения u(х,t) в точках следующего (j+1) -го слоя t=(j+1)k (риc. 7). При вычислении пользуются четырьмя соседними узлами – явная схема вида (рис. 7).

 

 

Рис. 7. Явная схема с использованием одного предыдущего слоя j

 

Таким образом, исходя из начального слоя t=0, значения u(х,t) для которого определяются из начального условия u(x_i,0)=f(x_i), (i=0, 1,…, n), и используя значения функции u(х,t) в крайних узлах (0, t_j), (1, t_j) (j=0, 1, 2,…), определяемые граничными условиями u(0,t_j)= φ(t_j), u(l,t_j)= ψ(t_j), по формуле (51) последовательно вычисляем: u(x_i,t₁), u(x_i,t₂), u(x_i,t₃),… (i=0, 1,…, n), т. е. находим значения искомой функции u(х,t) во всех узлах полуполосы.

Остается разумно выбрать величину σ. При этом будем исходить из требования, чтобы ошибка при замене дифференциального уравнения (49) конечно-разностным уравнением (50) была наименьшей. Введем обозначения:

 

, L_h[u]=1/h²[(u_i+1,j –2u_ij+u_i-1,j)-1/ σ(u_i,j+1 –u_ij)],

 

где L_h[u] – конечно-разностный оператор, соответствующий дифференциальному оператору L[u].

Разность R_h[u]=L_h[u]-L[и], называемая ошибкой аппроксимации, есть погрешность, которая получается при замене оператора L[и] оператором L_h[u]. Для L_h[u] можно записать

 

. (52)

 

Тогда выберем число σ так, чтобы первая скобка формулы (52) обратилась в нуль, т. е. положим σ/2=1/12 и, следовательно, σ=1/6. При этом значении σ будем иметь

 

. (53)

 

При выполнении равенства R_h[u]=L_h[u] при таком выборе σ для погрешности R_h[u] получаем оценку R_h[u]=О(h⁴), тогда как при другом выборе числа σ имеем R_h[u]=О(h²). В этом смысле значение σ =1/6 является для расчетной схемы 1 наилучшим.

Соответствующая расчетная формула (51) при таком выборе σ окончательно принимает вид

 

u_i,j+1=1/6(u_i-1,j+4u_ij+u_i+1,j). (54)

 

Отметим, что оценка ошибки аппроксимации R_h[u] в общем случае для граничных узлов (x_i, t_j) не годится.

 

Для уравнений гиперболического типа

 

Рассмотрим простейшее уравнение гиперболического типа:

 

. (55)

 

Будем искать решение уравнения (55) при заданных начальных и краевых условиях:

 

u(х,0)=f(х), u_t(х,0)=F(х) (0 ≤ х ≤ l) (56)

 

u(0,t)= φ(t), u(l,t)= ψ(t) (0 ≤ t < ∞). (57)

 

Как и в случае параболического уравнения, решая задачу методом сеток, покроем полуполосу (0 ≤ х ≤ l), (0 ≤ t <∞) прямоугольной сеткой x_i=ih (i=0, 1,…, n), t_j=jk (j=0, 1, 2,…), где ∆x_i=x_i+1-x_i=h=l/n (n – целое) и ∆t_j=t_j+1-t_j=k. На сетке x_i, t_j приближенно заменим дифференциальное уравнение (55) соответствующим конечно-разностным уравнением.

Пользуясь симметричными формулами для производных, будем иметь

 

, (58)

 

При k=h/а уравнение (58) упрощается и принимает вид

 

u_i,j+1=u_i+1,j+u_i-1,j-u_i,j-1 . (59)

 

Из уравнения (59) видно, что для получения значений u(х,t) в (j+1) -м слое используются значения u(х,t) в двух предыдущих слоях: j -м и (j-1)- м (рис. 8).

 

 

Рис. 8. Схема с использованием двух предыдущих слоев j-1 и j

 

Для начала вычисления по формуле (59) также необходимо знать значения u(х,t) на двух слоях, в то время как начальные условия (56) задают нам значения u(х,t) лишь на нулевом слое j=0. Однако, используя начальные условия, можно определить значения u(х,t) на фиктивном слое с номером j=-1. Для этого заменим производную во втором начальном условии конечно-разностным отношением.

 

, или u_i,-1=u_io –kF_i, (60)

 

где F_i=F(X_i). Отсюда u_i,-1=u_io –kF_i. Теперь, зная значения u(х,t) на слое j=-1, определяемые с помощью формулы (60), можно начать вычисления. Краевые условия (57) используются для получения значений u_oj и u_nj.

Вместо определения значений u(х,t) на слое j=-1 можно вычислить значения u(х,t) на слое j=1. Это достигается, например, с помощью формулы Тейлора:

 

. (61)

 

Учитывая вид уравнения (55), а также из начальных условий (56) предполагая, что , получив , , , перепишем формулу (61) в другом виде, а именно:

 

. (62)

 

Очевидно формулу (62) целесообразно применять в том случае, когда функция f(х) задана аналитическим выражением.