Моделирование распределения потенциала в МДП-структуре

U^p+1_mn|г = U’_mn|г = (s_mn)

U⁰_mn = ₀(x_m,y_n)

Будем считать, что ₀(x_m,y_n) по уже известному U^p={U^p_mn} для схемы (4) оссуществляется по уже явным формулам.

Вычисление U^p+1 = {U^p+1_mn} по схеме (5) требует решения задачи :

_xxU^p+1_mn + _yyU^p+1_mn - U^p+1_mn = (x_m,y_n) - U^p_mn

  (7)

U^p+1_mn|г = (s_mn)

Вычисление U^p+1 = {U^p+1_mn} по уже известным U^p = {U^p_mn} по схеме (6) осуществляется прогонками в направлении оси OX для вычисления решений {U’_mn} одномерных задач при каждом фиксированом n, а затем прогонками в направлнии оси OY для вычисления решений {U^p+1_mn} одномерных задач при каждом фиксированом m.

Для каждой из двух разностных схем (4) и (6) рассмотрим разность для счёта погрешностеи вычислений:

^p_mn = U^p_mn - U_mn

между сеточной функцией U^p = {U^p_mn} и точным решением U = {U_mn} задачи (1).

Решение {U_mn} задачи (1) удовлетворяет уравнениям:

U^p_mn - U_mn = _xxU_mn - (x_m,y_n)



U_mn|г = (s_mn)

U⁰_mn = U_mn

Вычитая эти равенства из (4) почленно, получим для погрешности ^p_mn следующую разностную задачу:

^p+1_mn - ^p_mn = _xx^p_mn + _yy^p_mn



^p+1_mn|г = 0 (9)

⁰_mn = ₀(x_m,y_n) - U_mn

Сеточная функция ^p_mn при каждом p (p=0,1,...) обращается в ноль на границе Г.

Метод переменных направлений

Рассмотрим двумерное уравнение теплопроводности:

dU = LU + f(x,t) , xG₀₂ , t[0,t₀]

dt

U|г = (x,t) (1)

U(x,0) = U₀(x)

LU = LU = (L_ +L₂)U , где L_U = d²U , =1,2

dx²

Область G_0 =G₀ = {0<= x_ <=l_ , =1,2} -прямоугольник со сторонами l₁ и l₂, Г - граница G₀ = G₀ + Г.

В G₀ построили равномерную по xa сетку _h с шагами h₁ = l₁/N₁ , h₂ = l2/N₂. Пусть _h - граница сеточной области _h, содержащая все узлы на сторонах прямоугольника, кроме его вершин, _h_h + _h.

Оператор L_ заменим разностным оператором _:

_y = yx_x_ , __

В случае одномерного уравнения теплопроводности неявная схема на каждом слое приводит к разностной краевой задаче вида:

A_iy_i-1 - C_iy_i + B_iy_i+1 = -F , i=1,...,N-1

y₀=₁ (2)

y_n=_N

A_i > 0, B_i > 0, C_i > A_i + B_i

которая решается методом прогонки.

Рассмотрим теперь нашу двимерную задачу в прямоугольнике. Сетку _h можно представить как совокупность узлов, расположенных на строках i₂=0,1,2,...,N₂, или как совокупность узлов расположенных на столбцах i₁=1,2,...,N₁. Всего имеется N₁+1 столбцов и N₂+1 строк. Число узлов в каждой строке равно N₁+1, а в каждом столбце N₂+1 - узлов.

Если на каждой строке (или столбце) решать задачу вида (2) методом прогонки при фиксированом i₂(или i₁), то для отыскания решения на всех строках (или столбцах), т.е. во всех узлах сетки, понадобится О(N₁N₂) арифметических действий. Основная идея большинства экономичных методов и состоит в сведении перехода со слоя на слой к последовательному решению одномерных задач вида (2) вдоль строк и вдоль столбцов.

Наряду с основными значениями искомой сеточной функции y(x,t), т.е. с y = yⁿ и y` = yⁿ⁺¹ вводится промежуточное значение y = y^n+Ѕ , которое можно формально рассматривать как значение при t = t_n+Ѕ = _n+Ѕ . Переход от слоя n на слой n+1 совершается в два этапа с шагами 0.5t .

y^n+Ѕ - yⁿ = ₁y^n+Ѕ + ₂yⁿ + ⁿ (3)

0.5t

yⁿ⁺¹ - y^n+Ѕ = ₁y^n+Ѕ + ₂yⁿ⁺¹ +ⁿ (4)

0.5t

Эти уравнения пишутся во всех внутренних узлах x = x_i сетки _h и для всех t=t_h > 0.

Первая схема неявная по направлению х₁ и явная по х₂, вторая схема явная по х₁ и неявная по х₂. К уравнениям (3),(4) надо добавить начальные условия:

y(x,0) = U₀(x) , x_h (5)

и разностно краевые условия, например, в виде:

yⁿ⁺¹ = ⁿ⁺¹ при i₁=0, i₂=N₂ (6)

y^n+Ѕ =  при i₁=0, i₂=N₁ (7)

где  = 1 (ⁿ⁺¹ + ⁿ) -  L₂(ⁿ⁺¹ - ⁿ) (8)

2 4

Т.о. , разностная краевая задача (3)-(8) соответствует задаче (1). Остановимся на методе решения этой задачи. Пререпишем (3) и (4) в виде:

2 y - ₁ y = F , F = 2 y + ₂ y + 

• 9)

2y` - <sub>2</sub> y` = F’ , F = 2 y + ₁ y + 

• 

Введём обозначения:

x_i = (i₁h₁ , i₂h₂)

F = F_i1,i2

y = y_i1,i2

при этом, если в уравнении один из индексов фиксирован, то его не пишем. Тогда (9) можно записать в виде (2), т.е.:

1 y_i1-1 - 2 1 + 1 y_i1 + 1 y_i1+1 = - F_i1

h²₁ h²₁  h²₁

i1 = 1,...,N₁-1 (10)

y = при i1 = 0,N₁

1 y`<sub>i2-1</sub> - 2 1 + 1 y`_i2 + 1 y`_i2+1 = - F_i2

h²₂ h²₂  h²₂

i2 = 1,...,N₂-1 (11)

y` = ` при i2 = 0,N₂

Пусть задано у=уⁿ. Тогда вычисляем F, затем методом прогонки вдоль строк i₂=1,...,N₂-1 решаем задачу (10) и определим y’ во всех узлах сетки _h, после чего вычисляем F и решаем задачу (11) вдоль столбцов i₁=1,...,N₁-1, определяя y`=yⁿ⁺¹. При переходе от слоя n+1 к слою n+2 процедура повторяется, т.е. происходит всё время чередование направлений.

Построение разностных схем

Для каждой области МДП - структуры построим консервативную разностную схему, учитывая при этом заданные условия.

Разобьём данную МДП - структуру на несколько областей следующим образом:

L M N

y

K₀

K₁

x

I : j_k0,y = U_n

 . ^k+Ѕ_i-1,y + 1 +  +  . ^k+Ѕ_ij -  . ^k+Ѕ_i+1y = _ij

2h^*_ih_i 2h^*_ih_i+1 2h^*_i2h_i 2h^*_ih_i+1

_k1,y = U_n

где _ij = ^k_ij +  (_y^k_ij + f^k_ij )

2

_y = 1 ^k_ij+1 - ^k_ij - ^k_ij - ^k_ij-1

r^*_j r_j+1 r_j

II: _ij=U₃

 . ^k+Ѕ_i-1,j + 1 +  +  . ^k+Ѕ _ij -  ^k+Ѕ_i+1,j =

2h^*_ih_i 2h^*_ih_i+1 2h^*_ih_i 2h^*_ih_i+1

^k_ij +  _y^k_ij

2 , 0 < i < k0-1 L< j <M

_ok . ^k+Ѕ _i-1,j + - _nn - _ok . ^k+Ѕ _ij + _n . ^k+Ѕ _i+1,j = ^*_ij , i=k0

h^*_i-1 h^*h_i h^*h_i-1