ЖУРНАЛ ВЫЧИСЛИТЕЛЬНОЙ МАТЕМАТИКИ И МАТЕМАТИЧЕСКОЙ ФИЗИКИ, 2015, том 55, № 6, с. 928-932
УДК 519.653
О НЕУСТОЙЧИВОСТИ ПОРЯДКА СИММЕТРИЧНЫХ ФОРМУЛ ЧИСЛЕННОГО ДИФФЕРЕНЦИРОВАНИЯ И ИНТЕГРИРОВАНИЯ
© 2015 г. В. М. Вержбицкий
(426069 Ижевск, ул. Студенческая, 7, ИжГТУ) Поступила в редакцию 29.07.2014 г.
Рассматривается вопрос об устойчивости порядка симметричных формул аппроксимации производных, применяемых в конечноразностных методах решения дифференциальных уравнений. Вводится определение устойчивости порядка формулы численного дифференцирования в процессе смещения точки, в которой она используется. Описываются условия и приводятся некоторые результаты вычислительных экспериментов по выявлению характера поведения порядка простейших симметричных формул аппроксимации первых и вторых производных в указанном процессе. На примерах показывается неустойчивость максимального порядка этих формул. Аналогично исследуется семейство квадратурных формул прямоугольников и демонстрируется неустойчивость второго порядка квадратурной формулы средней точки. Библ. 2. Табл. 7.
Ключевые слова: численное дифференцирование, аппроксимация производной, квадратурная формула прямоугольников, порядок остаточного члена, устойчивость порядка формулы.
БО1: 10.7868/80044466915060149
При построении конечноразностных методов решения граничных задач для дифференциальных уравнений на системах равноотстоящих узлов с шагом Н важная роль отводится симметричным формулам аппроксимации производных. Типичным примером может служить широко известная схема Кранка—Николсон — частный случай семейства разностных схем с весами а е [0,1], который выделяется заданием значения ст = 0.5, обеспечивающего второй порядок аппроксимации уравнения теплопроводности по временной переменной (см. [1], [2]). Именно при этом значении параметра а производная по времени аппроксимируется симметричным разностным отношением с остаточным членом, имеющим порядок квадрата шага.
Рассмотрим сначала вопрос об устойчивости остаточного члена формулы симметричной аппроксимации производной.
Определение 1. Остаточный член г0(к) приближенной формулы /(к)(х0) ~ у0Хк) называется устойчивым в точке х0, если
ф,5):= /(%0 + 5) - у0к)(Л) т. (1)
Непосредственно из (1) видно, что для устойчивости остаточного члена приближенной формулы / (к)(х0) ~ у^Хк) достаточно непрерывности в окрестности точки х0 самой к-й производной. Следовательно, для представляющего особый интерес остаточного члена формулы симметричной аппроксимации первой производной, имеющего вид
т : = /(х) - /(Х0 + к- /(х0 - к) = -^0)
2к 6
(в предположении, что /(х) е С (х0 - к, х0 + к)), устойчивость имеет место.
Теперь зададимся вопросом, что можно сказать о порядке остаточного члена, если формула симметричной аппроксимации в точке х0 применяется для приближенного вычисления производной в близкой к х0 точке х0 + 5. Будет ли порядок устойчив?
Для ответа на этот вопрос откажемся от привычного представления о порядке аппроксимации как о целой величине и дадим следующее
Определение 2. Пусть к -я производная функции /(х) при малых к > 0 аппроксимируется в точке х0 значением у0(к)(К), причем имеет место равенство /^к)(х0) = у0\к) + С0(к)кр0, где р0 > 1 — не зависящая от к постоянная. Порядок аппроксимации р0 назовем устойчивым в точке х0, если остаточный член в приближенном равенстве /^(х) ~ у^)(к) при любом х е (х0 - к, х0 + к) имеет вид С(х, к)кр(х), где |С(х, к)| < да и р(х) ^ р0 при к ^ 0.
Определение 3. Порядок аппроксимации р0 назовем устойчивым на отрезке [а, Ь], если он устойчив в любой точке х е (а, Ь).
Не имея подходящего аппарата для проведения аналитического исследования поставленной задачи, представим некоторые результаты вычислительных экспериментов.
Для выявления поведения порядка симметричной формулы аппроксимации первой производной
у < (к) = /х + к) - /(х0 - к) 2к
в процессе смещения точки применения этой формулы из позиции х0 в позицию х0 +8 эксперименты организованы следующим образом.
У заданной функции/(х) в заданной точке х0 вычисляем значения у'0(к1), у'0(к2),..., у'0(к„) на некоторой системе шагов к1, к2, ..., кп, уменьшающихся по правилу к = ук{-1 с у е (0,1). Согласно теории при к ^ 0 в случае непрерывности /'"(х) при х е (х0 - кх, х0 + к1) величины г0(к¡) := / '(х0) - у'(>(к1) должны стремиться к нулю по закону С(к(. Полагая, что имеет место функциональная связь вида г0(к) = Ск^, на основе множества вычисленных значений г0(кг) методом наименьших квадратов находим значения параметров С0 и р0. Убедившись, что значение р0 = 2 достаточно четко определяется при варьировании значений п е N и у е (0,1), приступаем к следующему этапу.
Закладываем процесс смещения точки из позиции х 0 в позиции х0 + 5 ], в которых вычисляются значения производной, т.е. для которых считается /'(х0 + 5¡) ~ у0(к). Значения смещения 5]
связываем с величиной шага аппроксимации: 5, := , где I е N — варьируемый параметр, а
I
]придаются значения 0,1,..., I. Вычисляя значения г](к) := /'(х0 + §у-) - у'0(к1), получаем статистику поведения остаточного члена при использовании одних и тех же значений, находимых по формуле симметричной аппроксимации, для приближения производной в смещенных (в пределах шага к) точках. К полученным данным применяем метод наименьших квадратов для нахождения при каждом ] параметров С } и р] в предполагаемой функциональной зависимости вида
г (к) = Скр.
Хорошо известно, насколько существенную роль в образовании полной погрешности результата численного дифференцирования играют погрешности, содержащиеся в значениях функции. Например, для рассматриваемого случая симметричной аппроксимации первой
8 М 2
производной в точке х0 полная погрешность оценивается величиной g(к)::--+--к , где
к 6
М : = тах |/"'(х)|, а е — уровень абсолютных погрешностей используемых значений функ-
х^[ж0-к, х01+к]
ции /(х) (см. [2, с. 527]). Чтобы устранить влияние этих погрешностей на результаты численных экспериментов, для значений функции применялись числа типа шр!7 из библиотеки GMP с точностью до 1000 знаков. Эксперименты проводились со степенными функциями с целыми и с дробными показателями.
Ниже приводятся некоторые численные результаты.
В табл. 1 и 2 для двух функций показывается, насколько точно определяется значение порядка р0 при различных сочетаниях параметров численного эксперимента.
Таблица 1
й0 = 0.1 у = 0.5 и = 7 и = 100 р = 2.0000203565 р = 2.0000001261
/(X) = X5, у = 0.1 и = 7 и = 100 р = 2.0000052045 р = 2.0000000289
Х0 = 3 й0 = 0.001 у = 0.5 и = 7 и = 100 р = 2.0000000020 р = 2.0000000000
у = 0.1 и = 7 и = 100 р = 2.0000000005 р = 2.0000000000
Таблица 2
й0 = 0.1 у = 0.5 и = 7 и = 100 р = 2.0000890964 р = 2.0000005519
/X) = 4х, у = 0.1 и = 7 и = 100 р = 2.0000227804 р = 2.0000001267
Х0 = 3 й0 = 0.001 у = 0.5 и = 7 и = 100 р = 2.0000000089 р = 2.0000000001
у = 0.1 и = 7 и = 100 р = 2.0000000023 р = 2.0000000000
Таблица 3
У /X) = X5, Х0 = 3 /(X) = X5, X,) = -3 /(X) = 4х, Xo = 3
РУ РУ ру
0 0 2.0000000000 2.0000000000 2.0000000000
1 0.0001 0.9999995373 1.0000004621 1.0000002311
2 0.0002 0.9999997902 1.0000002097 1.0000001049
3 0.0003 0.9999998840 1.0000001160 1.0000000580
4 0.0004 0.9999999380 1.0000000620 1.0000000310
5 0.0005 0.9999999762 1.0000000238 1.0000000119
6 0.0006 1.0000000064 0.9999999937 0.9999999968
7 0.0007 1.0000000320 0.9999999680 0.9999999840
8 0.0008 1.0000000548 0.9999999452 0.9999999726
9 0.0009 1.0000000758 0.9999999243 0.9999999621
10 0.001 1.0000000954 0.9999999047 0.9999999523
Табл. 3 показывает динамику изменения порядка при смещении точки, к которой применяется формула симметричной аппроксимации производной при следующих значениях параметров: Н0 = 0.001, у = 0.1, и = 100, I = 10.
Такая динамика сохраняется и при более мелком шаге смещения 5 ]. Например, при I = 100 и тех же значениях Н0 = 0.001, у = 0.1, и = 100 в точке х0 + 81 = 3.00001 для функции /(х) = х5 получаем р = 0.9999951979, для функции /(х) = л/Х получаем р = 1.0000023738. Налицо резкий, скачкообразный сброс значения показателя порядка симметричной формулы с максимального на минимальный, говорящий о его неустойчивости.
Таблица 4
У $ /(х) = х5, х0 = 3 /х) = х5, Х0 = -3 /х) = 4х, х0 = 3
Р Р Р
0 0 2.0000000000 2.0000000000 2.0000000000
1 0.0001 0.9999846964 1.0000001493 1.0000001867
2 0.0002 0.9999939069 1.0000000604 1.0000000755
3 0.0003 0.9999975628 1.0000000244 1.0000000305
4 0.0004 0.9999998565 1.0000000016 1.0000000020
5 0.0005 1.0000016090 0.9999999841 0.9999999801
6 0.0006 1.0000030915 0.9999999693 0.9999999616
7 0.0007 1.0000044198 0.9999999560 0.9999999450
8 0.0008 1.0000056515 0.9999999436 0.9999999295
9 0.0009 1.0000068186 0.9999999319 0.9999999148
10 0.001 1.0000079402 0.9999999205 0.9999999007
Таблица 5
/(х) ] 0 1 2 3 10
$ 0 0.0000001 0.0000002 0.0000003 0.0000001
х5 г(к) -9.0е- 011 5.4е-005 1.1е—004 1.62е-004 5.4е-004
4х г(к) -4.0е- 015 -4.8е-009 -9.6е-009 -1.4е-008 -4.8е-008
Проведенные в тех же описанных выше условиях эксперименты с симметричной формулой аппроксимации второй производной второго порядка
у'Х) = /(*0 - к) - 2ДХо) + Д*0 + к) - к^у к2 12
приводят к совершенно аналогичным результатам, частично отраженным в табл. 4.
Заметим, что в каждом из многочисленных экспериментов с применением симметричных формул аппроксимации первой и второй производных наблюдалось монотонное изменение значений Су, подсчитываемых наряду со значениями ру в процессе увеличения величины смещения 5у из точки х0. Это можно расценивать как аргумент в пользу выполнимости фигурирующего в
определении 2 предположения о сохранении вида зависимости С(х, к)крх) для выражения остаточного члена при использовании симметричных формул в смещенных точках.
Скачкообразное снижение порядка симметричной формулы аппроксимации производных можно связать с соответствующим поведением ошибки аппроксимации при смещении точки, в которой она используется. Типичный пример такого поведения величины г (К) := у (х0 + 5) - у '(к) в точке х0 = 3 при к = 0.000001 показан в табл. 5.
При просмот
Для дальнейшего прочтения статьи необходимо приобрести полный текст. Статьи высылаются в формате PDF на указанную при оплате почту. Время доставки составляет менее 10 минут. Стоимость одной статьи — 150 рублей.