Разностные схемы для решения многомерного уравнения теплопроводности

Мы рассмотрим основные методы решения многомерных задач на примере построения разностных схем для двумерного уравнения теплопроводности.

Чтобы применить метод конечных разностей в этой ситуации, поступим совершнно так же, как и в одномерном случае. С этой целью построим множество равноотстоящих точек х_п (п = 0,…, N) на отрезке 0<�х<1_хсх₀= 0,x_N= l_x,х"+ — х_п = /г, а также множество равноотстоящих точек у_П1 {т = 0, …, М) на отрезке 0 < у < 1_у с у₀ = 0, у_м = l_y, у_т+ - у_т = h. После этого область, где ищется решение, покрывается прямоугольной пространственной сеткой: через каждую точку х_п проводятся прямые параллельные оси г/, а через каждую точку у_т — прямые параллельные оси х (рис. 11.12).

Рис. 11.12. Двумерная пространственная сетка

В дополнение введем множество точек по времени ^ (k = 0, 1,…) с tk+1 - tk = т. Тогда типичная точка разностной сетки имеет координаты вида (х_п, у_т, fy). После этого входящие в уравнение (11.25) частные производные заменяются соответствующими конечно-разностными аппроксимациями.

Для более компактной записи разностных схем примем следующие обозначения для разностных операторов, аппроксимирующих вторые производные по х и у:

Тогда двумерный аналог схемы (11.21) записывается как.

Шаблон, схематически представляющий эту схему, показан на рис. 11.13.

Если у > 0, то схема (11.26) неявная; для каждого значения k необходимо решать разностное уравнение вида

относительно неизвестных значений в момент времени ?*₊₁. Эффективное решение такого уравнения возможно только для довольно узкого класса задач (параграф 11.3). С другой стороны, явная схема (11.43) имеет довольно сильное огра-.

Рис. 11.13. Шаблон для схемы (11.43).

ничение на шаг, но времени (для квадратной сетки р < ¼), что также может потребовать большого объема вычислений. Поэтому мы не будем останавливаться на подробном анализе этой схемы, а перейдем к рассмотрению более экономичных разностных схем построенных на основе метода расщепления.

Существуют различные способы построения разностных схем, основанных на уравнениях (11.40) и (11.41). Мы рассмотрим лишь несколько примеров. Если использовать явную схему (11.20) для решения задач (11.40) и (11.41), то получим явную схему расщепления

+ разностная форма граничных условий,.

которая представляется в виде шаблона, показанного на рис. 11.14.

При записи этой схемы мы ввели обозначения w^k_n>m = v_{n m}

k+ — k+1 т г и w_{n m} = и_{п т}. Параметры n_s> п_е, m_s, т_е зависят от типа граничных условий:

Рис. 11.14. Шаблон для схемы (11.44)

Вычисления по схеме (11.44) производятся следующим образом. На первом шаге, при каждом значении т, решаются одномерные задачи по х, что дает промежуточное решение v = {v" _т}. Затем, решая, при каждом значении п, одномерные задачи по у, мы получаем сеточную функцию i/⁺¹, которая и считается приближенным решением в момент времени 4₊₁.

Для анализа аппроксимации запишем сначала схему (11.44) в виде одного уравнения. Для этого выразим v_nm из первого уравнения и подставим это выражение во второе уравнение. В результате получим.

Подставляя в это выражение точное решение и учитывая, что операторы Д_хи_{п т} и и_{п т} аппроксимируют вторые производные по х и у (11.11), приходим к дифференциальному приближению в точке (.г," у_т, Г/₍).

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

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

Подставим это решение в (11.45) и получим (функция //…) при этом не учитывается, так как она не влияет на устойчивость):

Решение этого уравнения относительно X дает выражение для собственного значения оператора перехода на один шаг по времени:

Поэтому.

и для выполнения критерия фон Неймана необходимо |1 — 4р| < 1, откуда следует условие устойчивости [3 < ½. Это говорит о том, что мы можем выбрать вдвое больший шаг по времени по сравнению с явной схемой (11.43).

Если мы используем для решения задач (11.40) и (11.41) неявную схем}' (11.21) с у = 1, то получим неявную схему расщепления.

которая представляется в виде шаблона, показанного на рис. 11.15.

Рис. 11.15. Шаблон для схемы (11.44).

Для последовательного вычисления решения по схеме (11.47) на первом шаге, при каждом значении т, решаются системы уравнений вида (11.28), где.

Затем, при каждом значении п, решаются системы уравнений вида (11.28), где.

В результате мы получаем требуемое решение в момент времени 4+н Для нахождения порядка аппроксимации выразим v_nm из второго уравнения и подставим результат в первое уравнение. Это дает следующее разностное уравнение:

Анализ этого уравнения показывает, что схема расщепления (11.47) аппроксимирует дифференциальную задачу (11.42) с точностью до членов 0(т + А²). Подстановка решения (11.46) в уравнение (11.48) дает следующее выражение для собственных значений оператора перехода:

Тогда |Ца_1? сх.2) | < 1 для всех значений сц, а₂> и значит схема (11.47) устойчива при любых значениях р > 0.

Рассмотренный нами способ расщепления не является единственно возможным. Например, уравнение (11.38) можно представить в следующей форме:

На основе такого расщепления можно записать следующую схему:

которая имеет шаблон, показанный на рис. 11.16.

Эта схема получила название метод переменных направлений (ADI, Alternating Direction Implicit). Так как эта схема является неявной, то процедура расчета по этой схеме аналогична процедуре расчета по неявной схеме расщепления (11.49). Более подробная запись, анализ аппроксимации и устойчивости схемы (11.49) оставляется читателю в качестве упражнения.

В заключение обсудим расчет граничных условий. Разностная форма граничных условий первого рода практиче;

Рис. 11.16. Шаблон для схемы (11.49).

ски одинакова для всех рассмотренных нами схем. Например, условие и ( О, у, t) = J (у> t) в случае явной схемы (11.26) записывается как.

В случае схем расщепления (11.44) и (11.48) граничное условие имеет вид.

и так как на втором шаге расчета по этим схемам это значение не изменяется, то можно положить Uq '"] = В случае граничных условий второго рода их разностная форма зависит от вида основной разностной схемы. В качестве примера рассмотрим граничное условие.

Для явной схемы (11.43) применение метода фиктивных областей дает следующую разностную форму этого граничного условия:

Представление граничных условий для схем расщепления (11.44) и (11.47) непосредственно следует из выражения (11.26). Поэтому для схемы (11.44) можно записать.

Аналогично, для схемы (11.47) получим.

Пример 11.3 (нагрев пластины при воздействии на нее теплового потока) Рассмотрим квадратную пластину размером / и толщиной d. Пусть при 2 = 0 она теплоизолирована, а при z = d на нее воздействует распределенный тепловой поток вида.

В начальный момент времени t = 0 пластина имеет температуру щ. Края пластины теплоизолированы. Тогда распределение температуры в пластине подчиняется уравнению (11.25) с f/x, у у t) = g (x_y у, t) и условиям.

Здесь под и (х, у, t) понимается средняя по толщине пластины температура.

На рис. 11.17 показана зависимость температуры в центре пластины от времени, при к=1,/=1,/₀=1,с=10иг/₀=0. Точное решение этой задачи имеет вид.

где.

Рис. 11.17. Решение уравнения (11.42) с условиями (11.50):

• —точное решение;——приближенное решение
• (h = 0.1,т = 0.005)

Приближенное решение вычисляется с помощью явной схемы расщепления (11.44).