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

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

Если одна из переменных имеет физический смысл времени t , то сетку обычно строят так, чтобы среди ее линий (или гиперплоскостей) были линии t = t m . Совокупность узлов сетки, лежащих на такой линии или гиперплоскости, называют слоем.

На каждом слое выделяют направления, вдоль которых меняется только одна пространственная координата. Например, для переменных x , y , t есть направления x (t = const , y = const ) и направление y (t = const , х = const ).

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

Определение. Узлы, в которых разностная схема записана на шаблоне, называются регулярными, а остальные – нерегулярными.

Нерегулярными являются обычно граничные узлы, а иногда также лежащие вблизи границы узлы (такие, что взятый около этого узла шаблон выходит за границу области).

Составление разностной схемы начинается с выбора шаблона. Шаблон не всегда определяет разностную схему однозначно, но существенно влияет на ее свойства; например, далее мы увидим, что на шаблоне рис. 26.2b нельзя составить хорошей разностной схемы для задачи теплопроводности (26.1). Для каждого типа уравнений и краевых задач требуется свой шаблон.

Явные и неявные разностные схемы

Обсудим вопрос о фактическом вычислении разностного решения. Большая часть физических проблем приводит к уравнениям, содержащим время в качестве одной из переменных. Для таких уравнений ставится обычно смешанная краевая задача, типичным случаем которой является задача теплопроводности (26.1).

К подобным задачам применяют послойный алгоритм вычислений. Рассмотрим его на примере схем (26.2) и (26.4).

В схеме (26.4) на исходном слоеm = 0 решение известно в силу начального условия. Положим m = 0 в уравнениях (26.4). Тогда при каждом значении индекса n уравнение содержит одно неизвестное ; отсюда можно определитьпри
Значенияиопределяются по краевым условиям (26.3). Таким образом, значения на первом слое вычислены. По ним аналогичным образом вычисляется решение на втором слое и т.д.

Схема (26.4) в каждом уравнении содержит только одно значение функции на следующем слое; это значение нетрудно явно выразить через известные значения функции на исходном слое, поэтому такие схемы называются явными.

Схема (26.2) содержит в каждом уравнении несколько неизвестных значений функции на новом слое; подобные схемы называются неявными. Для фактического вычисления решения перепишем схему (26.2) с учетом краевого условия (26.3) в следующей форме

(26.5)

На каждом слое схема (26.5) представляет собой систему линейных уравнений для определения величин
; правые части этих уравнений известны, поскольку содержат значения решения с предыдущего слоя. Матрица линейной системы трехдиагональная, и решение можно вычислить алгебраической прогонкой.

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

В этих обозначениях явная и неявная разностные схемы принимают соответственно следующий вид

Невязка. Рассмотрим операторное дифференциальное уравнение общего вида (не обязательно линейное)

Au = f , или Au – f = 0.

Заменяя оператор А разностным оператором A h , правую часть f – некоторой сеточной функцией , а точное решениеu – разностным решением y , запишем разностную схему

или
. (26.6)

Если подставить точное решение u в соотношение (26.6), то решение, вообще говоря, не будет удовлетворять этому соотношению
. Величину

называют невязкой.

Невязку обычно оценивают при помощи разложения в ряд Тейлора. Например, найдем невязку явной разностной схемы (26.4) для уравнения теплопроводности (26.1а). Запишем это уравнение в каноническом виде

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

Разложим решение по формуле Тейлора около узла (x n , t m ), предполагая существование непрерывных четвертых производных по х и вторых по t

(26.7)

где

Подставляя эти разложения в выражение невязки и пренебрегая, в силу непрерывности производных, отличием величин
от (x n , t m ) найдем

(26.8)

Таким образом, невязка (26.8) стремится к нулю при
и
Близость разностной схемы к исходной задаче определяется по величине невязки. Если невязка стремится к нулю приh и стремящихся к нулю, то говорят что такая разностная схема аппроксимирует дифференциальную задачу. Аппроксимация имеетр -й порядок, если
.

Выражение (26.8) дает невязку только в регулярных узлах сетки. Сравнивая (26.3) и (26.1б), легко найдем невязку в нерегулярных узлах

Замечание 1. Решение задачи теплопроводности с постоянным коэффициентом (26.1) в области непрерывно дифференцируемо бесконечное число раз. Однако учет пятых и более производных в разложении в ряд Тейлора (26.7) прибавит к невязке (26.8) только члены более высокого порядка малости поиh , т.е. по существу, не изменит вида невязки.

Замечание 2. Пусть по каким-либо причинам решение исходной задачи дифференцируемо небольшое число раз; например, в задачах с переменным коэффициентом теплопроводности, гладким, но не имеющим второй производной, решение имеет лишь третьи непрерывные производные. Тогда в разложении в ряд Тейлора (26.7) последними будут члены
не точно компенсирующие друг друга. Это приведет к появлению в невязке (26.8) члена типа
т.е. невязка будет иметь меньший порядок малости, чем для четырежды непрерывно дифференцируемых решений.

Замечание 3. Преобразовав выражение невязки с учетом того, что входящая в него функция u (x ,t ) есть точное решение исходного уравнения и для нее выполняются соотношения

Подставляя это выражение в (26.8), получим

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

Схема Куранта - Изаксона - Риса (КИР), которую иногда также связывают с именем С.К. Годунова, получается при , . Ее порядок аппроксимации . Схема КИР условно устойчива, т.е. при выполнении условия Куранта . Приведем разностные уравнения для схемы Куранта - Изаксона - Риса во внутренних точках расчетной области:

Эти схемы, имеющие также название схемы с разностями против потока (в англоязычной литературе - upwind) могут быть записаны в виде

Их преимущество состоит в более точном учете области зависимости решения. Если ввести обозначения

то обе схемы можно записать в следующих формах:

(потоковая форма разностного уравнения);

(здесь явно выделен член со второй разностью, придающий устойчивость схеме);

(уравнение в конечных приращениях).

Рассмотрим также метод неопределенных коэффициентов для построения разностной схемы правый уголок первого порядка точности для уравнения переноса

Схему можно представить в виде

Схема Куранта - Изаксона - Риса тесно связана с численными методами характеристик . Дадим краткое описание идеи таких методов.

Две последние полученные схемы (при разных знаках скорости переноса) можно интерпретировать следующим образом. Построим характеристику , проходящую через узел (t n + 1 , x m ), значение в котором необходимо определить, и пересекающую слой t n в точке . Для определенности считаем, что скорость переноса c положительна.

Проведя линейную интерполяцию между узлами x m - 1 и x m на нижнем слое по времени, получим

Далее перенесем вдоль характеристики значение u n (x") без изменения на верхний слой t n + 1 , т.е. положим . Последнее значение естественно считать приближенным решением однородного уравнения переноса. В таком случае

или, переходя от числа Куранта снова к сеточным параметрам,

т.е. другим способом пришли к уже известной схеме "левый уголок", устойчивой при . При точка пересечения характеристики , выходящей из узла (t n + 1 , x m , с n - м слоем по времени расположена левее узла (t n , x m - 1 ). Таким образом, для отыскания решения используется уже не интерполяция , а экстраполяция, которая оказывается неустойчивой.

Неустойчивость схемы "правый уголок" при c > 0 также очевидна. Для доказательства этого можно использовать либо спектральный признак, либо условие Куранта, Фридрихса и Леви. Аналогичные рассуждения можно провести и для случая c < 0 и схемы "правый уголок".

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

Схема Лакса - Вендроффа возникает при . Порядок аппроксимации схемы Лакса - Вендроффа есть . Схема устойчива при выполнении условия Куранта .

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

При выводе выражения для главного члена погрешности аппроксимации использовано следствие исходного дифференциального уравнения переноса

Которое получается путем дифференцирования исходного уравнения (3.3) сначала по времени t , затем по координате x и вычитанием одно из другого получившихся соотношений.

Далее, заменяя вторую производную во втором слагаемом в правой части с точностью до O(h 2) , получим новую разностную схему, аппроксимирующую исходное дифференциальное уравнение с точностью . Сеточные уравнения для схемы Лакса - Вендроффа во внутренних узлах расчетных сеток есть

Неявная шеститочечная схема возникает при q = 0 ; при ее порядок аппроксимации , при .

Пример 1. Разностная схема для уравнения Пуассона, относящегося к эллиптическому типу.

Рассмотрим построение разностной схемы для первой краевой задачи для уравнения А и = f(x,y) в области, представляющей собой прямоугольник со сторонами, параллельными осям координат. Пусть с этим прямоугольником связана равномерная сетка с шагами h x и h y .

Краевую задачу

можно записать в операторной форме:

Отметим, что в эту запись включены и граничные условия.

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

которые аппроксимируют исходное дифференциальное уравнение со вторым порядком 0(h 2 + h 2) точности и действуют во всех внутренних точках области.

Разностные аналоги краевых условий будут иметь вид

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

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

где в L/, включены как разностное уравнение, так и разностное краевое условие:

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

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

Относительно системы уравнений отметим два обстоятельства.

• 1. Система имеет очень высокую размерность (М - 1) х (N - 1), а традиционные методы точного решения (например, метод Гаусса) требуют для решения число алгебраических операций, пропорциональное третьей степени размерности системы.
• 2. В матрице системы много нулевых элементов (неплотная матрица). Это обстоятельства позволяет разработать экономичные методы приближенного решения.

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

Рис. 2.8.

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

Рассмотрим следующую задачу:

Отмстим, что в случае параболического уравнения имеем открытую область. При построении разностной схемы возникают несколько вариантов связи между разностными производными по пространству и по времени.

Проинтегрируем уравнение в пределах одного временного шага:

В зависимости от того, какую квадратурную формулу мы примем для вычисления интеграла в правой части, получим разные разностные схемы (рис. 2.9).

Связывая разностную производную по времени с пространственной производной, определенной на п -м временном слое, получим

явную ’разностную схему

Это эквивалентно приближенному вычислению интеграла, стоящего в правой части (2.12), но способу левых прямоугольников.

Рис. 2.9. Сетка и шаблоны для уравнения теплопроводности: а - область и сетка; б - шаблон явной схемы; в - шаблон неявной схемы; г - шаблон семейства шеститочечных схем; д - шаблон схемы

«чехарда»

В приведенной формуле заключен и метод решения сеточных уравнений:

Значение сеточной функции на следующем временном слое

определяется через известные значения гф на предыдущем. Перемещаясь последовательно слоями от начального условия и(х , 0) = у(х), можно найти решение во всей расчетной области. Разностный шаблон для этой схемы приведен на рис. 2.9, б .

Оценивая интеграл через значение подынтегральной функции па слое п + 1, используем разностный шаблон вида рис. 2.9, б, а разностный аналог дифференциального уравнения примет вид

Для того чтобы найти значения сеточной функции на следующем временном слое, при использовании этой разностной схемы необходимо решить совместно столько уравнений вида (2.14), сколько внутренних узлов расположено на п - 1-1 -м временном слое. С учетом краевых условий = / п+1 , Мд Г +1 = m n+1 система позволяет построить решение на следующем временном слое при известных значениях сеточной функции на предыдущем. Передвигаясь от начальных значений слоями, на каждом из которых необходимо решать систему уравнений, можно построить приближенное решение во всей области.

Рассмотренная разностная схема представляет собой пример неявной разностной схемы, ее называют схемой с опережением или чисто неявной схемой.

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

При а = 0 имеем явную схему, при а = I - неявную с опережением, при а > 0 - неявную. При а - 0,5 получаем широко известную в вычислительной практике симметричную схему Кранка Николсона.

Приведенные схемы, разумеется, не исчерпывают всего многообразия разностных схем, основанных на разностной аппроксимации дифференциальных операторов. Вот пример явной разностной схемы, основанной на центрировании временной производной, схемы, использующей сеточную функцию на трех временных слоях:

Разностный шаблон захватывает три временных слоя. Схема имеет второй порядок аппроксимации как по времени, так и по пространственной переменной и является явной. Эта схема обладает рядом существенных недостатков, от большей части которых можно избавиться, заменив и ” в аппроксимации пространственной производной средним значением по двум временным слоям:

Полученная таким образом явная трехслойная схема

называется схемой Дюфортпа - Франкела , а отсутствие в ней значения сеточной функции в центральном узле объясняет название «чехарда», которое иногда применяется для схем такого рода.

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

Используя шаблон для каждого внутреннего узла области решения апроксимируется уравнение теплопроводности

Отсюда найдем:

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

Затем с помощью соотношений

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

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

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

Поэтому рассматриваемая разностная схема, использующая этот шаблон, называется явной разностной схемой . Точность её имеет порядок .

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

Явная разностная схема является условно устойчивой . Если условие не выполняется, то небольшие погрешности вычислений, например, связанные с округлением данных компьютера приводит к резкому изменению решения. Решение становится неприемлемым для использования. Это условие накладывает весьма жесткие ограничения на шаг по времени, что может оказаться неприемлемым из-за значительного увеличения времени счета решения этой задачи.

Рассмотрим разностную схему, использующую другой шаблон

Метод 36

Неявная разностная схема для уравнения теплопроводности.

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

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

Схема решения задачи следующая:

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

Отличие от явной схемы - значения на очередном временном уровне вычисляются не непосредственно с помощью готовой формулы, а находится путем решения системы уравнений, т.е. значения неизвестных находятся неявно путем решения СЛАУ. Поэтому разностная схема называется неявной. В отличие от явной неявная является абсолютно устойчивой.

Тема №9

Задачи оптимизации.

Эти задачи являются одними из важнейших задач прикладной математики. Под оптимизацией понимают выбор наилучшего варианта из всех возможных решений данной задачи. Для этого необходимо сформулировать решаемую задачу как математическую, придав количественный смысл понятиям лучше или хуже. Обычно в процессе решения необходимо найти оптимизируемые значения параметров. Эти параметры называют проектными. А число проектных параметров определяет размерность задачи.

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

- целевая функция.

Наиболее просты случаи, когда целевая функция зависит от одного параметра и задаётся явной формулой. Целевых функций может быть несколько.

Например, при проектировании самолёта требуется одновременно обеспечить максимальную надежность, минимальные вес и стоимость и т.д. В таких случаях вводится система приоритетов . Каждой целевой функции ставится в соответствие некоторый целевой множитель в результате получается обобщенная целевая функция(функция компромиссов).

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

Теория и методы решения задач оптимизации при наличии ограничений составляют предмет исследований одного из разделов прикладной математики – математического программирования.

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

Требуется найти значения на при которых целевая функция имеет максимальное значение. Если целевая функция задана аналитически и может быть найдено выражение для её производных, то оптимальное решение будет достигаться либо на концах отрезка, либо в точках в которых производная обращается в ноль. Это критические точки и . Необходимо найти значения целевой функции во всех критических точках и выбрать максимальное.

В общем случае для нахождения решения применяют различные методы поиска. В результате происходит сужение отрезка содержащего оптимальное решение.

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

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

То есть, делаем пока выполняется (*)

При этом производится вычислений целевой функции.

Требуется найти такое значение, чтобы выражение (*) было получено при наименьшем

числе вычислений .

Метод 37

Метод половинного деления.

Рассмотрим метод поиска при . Он называется методом половинного деления, так как на каждом шаге отрезок, содержащий оптимальное решение уменьшается в два раза.

Эффективность поиска можно повысить путём специального выбора точек, в которых вычисляется целевая функция на определённом шаге сужения.

Метод 38

Метод золотого сечения.

Одним из эффективных способов является метод золотого сечения. Золотым сечением отрезка называется точка для которой выполняется условие

Таких точек две: =0,382 +0,618

0,618 +0,382 .

Отрезок делится точками и а после находится точка, целевая функция в которой максимальна. В результате чего находится изменённый отрезок длинною 0,618( - ) .

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

Метод 39

Метод покоординатного подъёма (спуска).

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

Математика и математический анализ

Решение разностной схемы называется приближенным решением дифференциальной задачи. Характеристика неявной разностной схемы Рассмотрим одномерное дифференциальное уравнение параболического типа с начальным и граничными условиями: 4.7 записана на n 1ом шаге по времени для удобства последующего изложения метода и алгоритма решения неявной разностной схемы 4. В разделе Порядок аппроксимации разностной схемы было отмечено что разностная схема 4.

8 вопрос: Разностные схемы: явная и неявная схемы:

Разностная схема — это конечная система алгебраических уравнений, поставленная в соответствие какой-либо дифференциальной задаче, содержащей дифференциальное уравнение и дополнительные условия (например краевые условия и/или начальное распределение ). Таким образом, разностные схемы применяются для сведения дифференциальной задачи, имеющей континуальный характер, к конечной системе уравнений, численное решение которых принципиально возможно на вычислительных машинах. Алгебраические уравнения, поставленные в соответствие дифференциальному уравнению получаются применением разностного метода , что отличает теорию разностных схем от других численных методов решения дифференциальных задач (например проекционных методов, таких как метод Галёркина ).

Решение разностной схемы называется приближенным решением дифференциальной задачи.

Характеристика неявной разностной схемы

Рассмотрим одномерное дифференциальное уравнение параболического типа с :

Запишем для уравнения (4.5) неявную разностную схему :

Запишем :

Аппроксимация граничных условий (4.7) записана на (n метода и алгоритма решения неявной разностной схемы (4.6).
В разделе " " было отмечено, что разностная схема (4.6) имеет такой же порядок аппроксимации , как и соответствующая ей явная разностная схема (4.2) , а именно:

В разделе " Доказательство абсолютной устойчивости неявной разностной схемы " было доказано, что неявная разностная схема (4.6) абсолютно устойчива, т.е. вне зависимости от выбора интервала деления на разностной сетке (или, иначе говоря, выбора расчётного шага по независимым переменным) погрешность решения неявной разностной схемы в процессе вычислений возрастать не будет. Отметим, что это, безусловно, является достоинством неявной разностной схемы (4.6) по сравнению с явной разностной схемой (4.2) , которая устойчива только при выполнения условия (3.12) . В то же время явная разностная схема имеет достаточно простой метод решения , а метод решения неявной разностной схемы (4.6), называемый методом прогонки , более сложен. Прежде чем перейти к изложению метода прогонки , необходимо вывести ряд соотношений , используемых этим методом.

Характеристика явной разностной схемы.

Рассмотрим одномерное дифференциальное уравнение параболического типа с начальным и граничными условиями :

Запишем для уравнения (4.1) явную разностную схему :

Запишем аппроксимацию начального и граничных условий :

Аппроксимация граничных условий (4.3) записана на (n + 1)-ом шаге по времени для удобства последующего изложения метода и алгоритма решения явной разностной схемы (4.2).
В разделе " Порядок аппроксимации разностной схемы " было доказано, что разностная схема (4.2) имеет порядок аппроксимации :

В разделе " Доказательство условной устойчивости явной разностной схемы " было получено условие устойчивости данной разностной схемы, накладывающее ограничение на выбор интервала деления при создании разностной сетки (или, иначе говоря, ограничение на выбор расчётного шага по одной из независимых переменных):

Отметим, что это, безусловно, является недостатком явной разностной схемы (4.2). В то же время она имеет достаточно простой метод решения .