дга яитегрироваиия неочевиден. Тогда имеет смысл использовать 2)лее сложные программы, в которых автоматически осуществляется

[бор шага, а иногда и квадратурной формулы. Такие программы азывают адаптивными квадратурными программами [32]. В них используется одна или несколько квадратурных формул и определяются величины Подынтервалов \Xj, Xj+tJ так, чтобы вычисленный результат удовлетворял предписанной точности. В разных частях HHTepBaiT интегрирования могут использоваться сетки разной густо-.jbi. Пользователь подобной программы указывает интервал [а, Ь], составляет подпрограмму-функцию, вычисляющую / (х), и выбирает допустимую погрешность е. Адаптивная программа пытается вычислить величину у так, чтобы выполнялось условие

lN-{[ix)dx

< е.

Программа может сделать вывод, что предписанная точность недостижима, получить наилучший возможный для нее результат и выдать оценку реально достигнутой погрешности. Пример такой Программы подробно рассмотрен в книге [32].

В заключение отметим, что при вычислении иесобствеииых интегралов в задачах о полупространстве возникает вопрос о выборе пределов интегрирования. Разумеется, приходится задавать конечный интервал [а, Ь], но предварительно или в самой программе проверять, чтобы отбрасываемая часть интеграла ие превышала допустимую погрешность.

§ 1Л. ПРОГРАММНАЯ РЕАЛИЗАЦИЯ

ТОЧНОГО АНАЛИТИЧЕСКОГО РЕШЕНИЯ ОДНОМЕРНЫХ ЗАДАЧ

Рассмотрим пример программы для расчета одномерного нестационарного температурного поля пластины по точному решению (2,13). Исходными данным» являются, во-первых, параметры, входящие в постановку задачи (2.1)-(2.3): толщина /, теплопроводность X, температуропроводность а, коэффициент теплоотдачи а, начальный перегрев во-вторых, массивы координат {лг;}/, и моментов времени {т;}-.для которых требуется вычислить температуру, п. наконец, число членов ряда N, учитываемых в решении (2.13), и Исло интервалов разбиения отрезка [О, /] при вычислении интегра-• (2.10) для изображения W.

Распределение объемной плотности мощности (х) задается с Помощью подпрограммы-функции с именем QV (X). Для примера

а программе задана следующая функция:

*7Л-)=*7«оехр1 -(Зх )].

Вообще говоря, число членов ряда Л можно выбирать в npoin-ме автоматически, оценивая величины слагаемых. В этом случае о? жется, что для малых моментов времени (прн числе Fo< 0,05) д достижения погрешности порядка единиц процентов понадобит? несколько десятков членов ряда, а при Fo > 0,5 будет дост;поч}, нескольких слагаемых. В данной учебной программе для упрощени. логики число Л фиксировано в исходных данных н должно иодб* раться путем пробных расчетов, что весьма часто н делают па црак тике.

Для определения собственных чисел р,, используется под11п)грам-ма KORNI из §2.2, а для вычисления интегралов W,, - Huirmn. грамма SIMPS нз §2.3.

Так как суммирование членов ряда требуется выполнять для каж-дого Xi и Xj, то чтобы сократить затраты машинного времени, целесообразно предварительно вычислить все комплексы, зависящие от собственных чисел р„, но не зависящие от xi и т ,-, и записать н.\ в соответствующие массивы длиной /V. При числе N порядка нескольких десятков это не потребует ощутимых затрат памяти, но в дальнейшем сократит машинное время прн вычнсленнн и (x;,Tj) по формуле (2.13). В случае двойных и тем более тройных рядов данная рекоменлацир уже становится проблематичной.

В рассматриваемой программе до суммирования ряда вычисляют следующие массивы:

А (N) - собственные числа р„; В1 (N) и В2 (N) - комплексы в решениях (л-, т) н (х, т):

0 I sin \in

\ln Win

При переборе моментов времени {т ,};(, н пространственных точек {.V;}., во внешнем цикле (операторы 50-65) выбирают моменты времени. При этом для каждого т ,- величины ехр (-р;;Ро;) записывают в рабочий массив С длиной N и благодари этому в цикле по Xf (операторы 57-61) повторно не вычисляют. Таким образом, нрнвы-числеинн температур d{Xi, Xj) путем суммирования ряда рассчитывают выражения:

{Xt, 2 {lnC„(Tj)+&2«[l-C„(T;)l}COS(p„Xi ).

Описание входных н выходных параметров программы дано в t ментариях к тексту.

ГЛАВА

КОНЕЧНО-РАЗНОСТНЫЕ МЕТОДЫ РЕШЕНИЯ ЗАДАЧ ТЕПЛОПРОВОДНОСТИ

Применение вычислительной техники и численных методов значительно расширяет классы исследуемых полевых задач теплообмена, позволяя гюлучать приближенные решения многомерных, нелинейных, нестационарных задач, для которых использование точных и приближенных аналитических методов ие представляется возможным. При выборе математических моделей, описывающих процессы теплообмена в реальных объектах, границы их допустимой сложности в настоящее время часто определяются не столько возможностями численных методов ]1 ресурсами ЭВМ, сколько недостатком достоверной входной информации для этих моделей.

При опрелеленнн различных пространственно-временных полей необходимо находить решения краевых задач для ди(})ференцнальцых уравнений в члстных производных в заданных областях изменения Пространственных переменных и временных интервалах. Отличительной особенностью применения численных методов является дискретизация пространствошой н времешюй областей на первом же этапе решения задачи. При дискретизации выбираются узловые точки в пространственной и временной областях. На втором этапе составляется система алгебраических уравнении относительно значений искомых функций в этих узловых точках. На третьем - прово-Д1ГГСЯ решение системы и находятся значения исследуемых величин в узловых точках. Отметим, что днскрепгчацня области часто делается и прн расчете па основе аналитических решений, однако в этих случаях она проводится на заключительных этапах, реализуемых уже после получения аналитического решения.

Существуют два основных численных метода решения уравнений в частных прои;водных; метод конечных разностей и метод конечных элементов. Онн отличаются сп-кобамн получения системы уравнении для значений искомых функций в узловых точках. Метод конечных разностей базируется непосредственно на дифференциальном Уравнении и граничных условиях, а метод конечных элементов - На эквивалентной варнацнонной постановке задачи.

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

них какие-то

§ 3.1. ОСНОВНЫЕ ПОНЯТИЯ ТЕОРИИ РАЗНОСТНЫХ СХЕМ

Теории численных методов решения уравнений в частных про. водных представляет собой весьма обширный и достаточно сложны" ра.чдел математики, называемый теорией разностных схем, с котц! рым можно познакомиться, например, по книгам 14, 14. 24, 2б] В данном учебном пособии основное внимание уделяется практичес КИМ вопросам построения и реализации на ЭВМ различных числе. ных методик, а не их теоретическому исследованию и обоснованию Как правило, будем ограничиваться лишь объяснением основных понятий, которые понадобятся в дальнейшем, причем некоторые вопросы рассмотрим не вполне строго с позиции математики.

Разностная схема и разностное решение. Основные нонятия теории разностных схем разберем на примере одномерного нестационарного уравнения теплопроводности для пластины с внутренним источником теплоты

дТ дх

-f 0<.v</, 0<т <т,

(3.1)

На границах пластины заданы граничные условия третьего рода дТ , „-

х = 0. I

(3.2)

а начальное условие имеет вид

Т{х, t)UoT,{x). (3.3)

Искомой в задаче (3.1)-(3.3) является функция Т (х, х),заданная в непрерывной области

й- {О Х<1}Х{0 < X <Тпах}.

При использовании чистенных методов ставится более скромная задача. В пространственной области выбирается некоторое конечное число значений координаты л:,, Xg, хм (узлы пространственной сетки), для временной переме]1Ной также выбирается конечное число значений То, т, %j (узлов временной сеткн). Цель - определение значений температуры Т в узлах пространственной сетки в моменты времени Tj-:

(3.4)

т. е. находятся значения искомой функции в дискретной области л. Лт (рис. 3.1):

Л.Д-с = {Хи Xs}X{Xq, ....

Для упрощения будем считать пространственное и временное разбиения равномерными с шагами h по координате х и Дт по времени:

Х=(Л- 1)/1, hl/{N~\), П=\, N- Tj =/&т, Лт =

= mas/-/, / о, 1, J.

асно что для определения TjJ необходимо иметь для ния и эти уравнения следует получить из основной задачи

Возьмем произвольный момент времени xj (только не начальный) оизвольный узел х,, (только не лежащий на границе, т. е. и ] уу) и обратимся к уравнению (3.1). Поставим задачу выразить оизводные дТ/дт и дТ/дх в точке (л,,, Xj) через значения функции

......

в этой точке и в некоторых со-„нмх узлах сетки Йл,дт. "Используя определение производной, запишем

г 7/ -

= iim ---

X, Дт->0 Дт

дТ дх

Т{.

- + биДт), (3.5)

где б« (Дт) ~ некоторая величина, Рнс. 3.1

стремящаяся к нулю при Дт-vO.

Конечную разность (Т - Т~) называют «разностью назад» или

левой разностью.

Найдем выражение для 6; (Дт), представив в (3.5) 7"/- с помощью разложения в ряд Тейлора в точке (x„, tj):

l дТ V

Тогда

(3.6)

При достаточно малом Дт выполняется неравенство 6i (Дт) I < Лг Дт, const.

Вспомним, что в математике условия типа (3.6) записывают в символическом виде 6 О (Дт).

Понятно, что аналогичным образом можно построить аппроксимацию для временной производной с по.мощью «разности вперед» (•и правой разности):

!L +.0(Дт).

(3.7)