Экстраполяция Ричардсона

В численном анализе , Экстраполяция Ричардсона является ускорение последовательности метод, используемый для улучшения скорости сходимости в виде последовательности оценок некоторого значения . По сути, учитывая значение для нескольких значений , мы можем оценить , экстраполируя оценки на . Он назван в честь Льюиса Фрая Ричардсона , который представил эту технику в начале 20 века. ^{[1] [2]} По словам Биркгофа и Роты , «его полезность для практических вычислений трудно переоценить». ^[3]${\ displaystyle A ^ {\ ast} = \ lim _ {h \ to 0} A (h)}$${\ displaystyle A (h)}$${\displaystyle h}$${\displaystyle A^{\ast }}$${\displaystyle h=0}$

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

Пример экстраполяции Ричардсона [ править ]

Предположим, что мы хотим приблизиться , и у нас есть метод, который зависит от небольшого параметра таким образом, что${\displaystyle A^{*}}$${\displaystyle A(h)}$${\displaystyle h}$

${\displaystyle A(h)=A^{\ast }+Ch^{n}+O(h^{n+1}).}$

Определим новую функцию

${\displaystyle R(h,t):={\frac {t^{n}A(h/t)-A(h)}{t^{n}-1}}}$где и - два различных размера шага.${\displaystyle h}$${\displaystyle {\frac {h}{t}}}$

потом

${\displaystyle R(h,t)={\frac {t^{n}(A^{*}+C\left({\frac {h}{t}}\right)^{n}+O(h^{n+1}))-(A^{*}+Ch^{n}+O(h^{n+1}))}{t^{n}-1}}=A^{*}+O(h^{n+1}).}$

${\displaystyle R(h,t)}$называется Ричардсон экстраполяцией из А ( ч ), и имеет оценку погрешности более высокого порядка по сравнению с .${\displaystyle O(h^{n+1})}$${\displaystyle A(h)}$

Очень часто гораздо проще получить заданную точность, используя R (h), а не A (h ') с гораздо меньшим h' . Где A (h ') может вызвать проблемы из-за ограниченной точности ( ошибки округления ) и / или из-за увеличения количества необходимых вычислений (см. Примеры ниже).

Общая формула [ править ]

Позвольте быть приближением (точное значение), которое зависит от положительного размера шага h с формулой ошибки вида ${\displaystyle A(h)}$${\displaystyle A^{*}}$

${\displaystyle A^{*}=A(h)+a_{0}h^{k_{0}}+a_{1}h^{k_{1}}+a_{2}h^{k_{2}}+\cdots }$

где a _i - неизвестные константы, а k _i - известные константы, такие что h ^{k _i} > h ^{k _{i + 1}} .

k ₀ - поведение размера шага ведущего порядка ошибки усечения, поскольку${\displaystyle A^{*}=A(h)+O(h^{k_{0}})}$

Искомое точное значение можно определить как

${\displaystyle A^{*}=A(h)+a_{0}h^{k_{0}}+a_{1}h^{k_{1}}+a_{2}h^{k_{2}}+\cdots }$

который можно упростить с помощью нотации Big O до

${\displaystyle A^{*}=A(h)+a_{0}h^{k_{0}}+O(h^{k_{1}}).\,\!}$

Используя размеры шага h и некоторую константу t , две формулы для A :${\displaystyle h/t}$

${\displaystyle A^{*}=A(h)+a_{0}h^{k_{0}}+O(h^{k_{1}})\,\!}$

${\displaystyle A^{*}=A\!\left({\frac {h}{t}}\right)+a_{0}\left({\frac {h}{t}}\right)^{k_{0}}+O(h^{k_{1}}).}$

Умножение второго уравнения на t ^{k ₀} и вычитание первого уравнения дает

${\displaystyle (t^{k_{0}}-1)A^{*}=t^{k_{0}}A\left({\frac {h}{t}}\right)-A(h)+(t^{k_{0}}-1)O(h^{k_{1}})}$

который можно решить, чтобы дать${\displaystyle A^{*}}$

${\displaystyle A^{*}={\frac {t^{k_{0}}A\left({\frac {h}{t}}\right)-A(h)}{t^{k_{0}}-1}}+O(h^{k_{1}}).}$

Таким образом, ошибка усечения была уменьшена до . Это в отличие от , где ошибка усечения находится за тот же размер шага${\displaystyle R(h,t)={\frac {t^{k_{0}}A\left({\frac {h}{t}}\right)-A(h)}{t^{k_{0}}-1}}}$${\displaystyle O(h^{k_{1}})}$${\displaystyle A(h)}$${\displaystyle O(h^{k_{0}})}$${\displaystyle h}$

Благодаря этому процессу мы достигли лучшего приближения к A , вычтя наибольший член ошибки, который был O ( h ^{k ₀} ). Этот процесс можно повторить, чтобы удалить больше ошибок, чтобы получить еще лучшие приближения.

Общее рекуррентное соотношение, начинающееся с, может быть определено для приближений формулой${\displaystyle A_{0}=A(h)}$

${\displaystyle A_{i+1}(h)={\frac {t^{k_{i}}A_{i}\left({\frac {h}{t}}\right)-A_{i}(h)}{t^{k_{i}}-1}}}$

где удовлетворяет${\displaystyle k_{i+1}}$

${\displaystyle A^{*}=A_{i+1}(h)+O(h^{k_{i+1}})}$.

Экстраполяцию Ричардсона можно рассматривать как преобразование линейной последовательности .

Кроме того, общая формула может использоваться для оценки k ₀ (поведение размера шага ведущего порядка ошибки усечения ), когда ни его значение, ни A ^* (точное значение) не известны априори . Такой метод может быть полезен для количественной оценки неизвестной скорости сходимости . Учитывая приближения A из трех различных размеров шага h , h / t и h / s , точное соотношение

${\displaystyle A^{*}={\frac {t^{k_{0}}A\left({\frac {h}{t}}\right)-A(h)}{t^{k_{0}}-1}}+O(h^{k_{1}})={\frac {s^{k_{0}}A\left({\frac {h}{s}}\right)-A(h)}{s^{k_{0}}-1}}+O(h^{k_{1}})}$

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

${\displaystyle A\left({\frac {h}{t}}\right)+{\frac {A\left({\frac {h}{t}}\right)-A(h)}{t^{k_{0}}-1}}\approx A\left({\frac {h}{s}}\right)+{\frac {A\left({\frac {h}{s}}\right)-A(h)}{s^{k_{0}}-1}}}$

которое можно решить численно, чтобы оценить k ₀ для некоторых произвольных выборов h , s и t .

Пример кода псевдокода для экстраполяции Ричардсона [ править ]

Следующий псевдокод в стиле MATLAB демонстрирует Экстраполяция Ричардсона , чтобы помочь решить ОДУ , с методом трапеций . В этом примере мы уменьшаем вдвое размер шага на каждой итерации, и поэтому в приведенном выше обсуждении это будет . Ошибка трапециевидного метода может быть выражена в терминах нечетных степеней, так что ошибку на нескольких шагах можно выразить в четных степенях; это приводит нас к возведению во вторую степень и к степеням псевдокода. Мы хотим найти значение , которое имеет точное решение, поскольку точное решение ОДУ равно . Этот псевдокод предполагает, что существует вызываемая функция, которая пытается вычислить${\displaystyle y'(t)=-y^{2}}$${\displaystyle y(0)=1}$${\displaystyle h}$${\displaystyle t=2}$${\displaystyle t}$${\displaystyle 4=2^{2}=t^{2}}$${\displaystyle y(5)}$${\displaystyle {\frac {1}{5+1}}={\frac {1}{6}}=0.1666...}$${\displaystyle y(t)={\frac {1}{1+t}}}$Trapezoidal(f, tStart, tEnd, h, y0)y(tEnd)путем выполнения метода трапеции для функции fс начальной точкой y0и tStartразмером шага h.

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

tStart = 0 % Время запуска   tEnd = 5 % Время окончания   f = - y ^ 2 % Производная y, поэтому y '= f (t, y (t)) = -y ^ 2    % Решение этого ОДУ: y = 1 / (1 + t)y0 = 1 % Исходное положение (т.е. y0 = y (tStart) = y (0) = 1)   Допуск = 10 ^ - 11 % Требуется 10-значная точность   maxRows = 20 % Не позволять итерации продолжаться бесконечно   initialH = Tstart - Tend % Pick первоначального размера шага     haveWeFoundSolution = false % Удалось ли найти решение в пределах желаемого допуска? еще нет.   h = initialH  % Создайте 2D-матрицу размером maxRows на maxRows для хранения экстраполяций Ричардсона% Обратите внимание, что это будет нижняя треугольная матрица и что на самом деле не более двух строк% требуется в любой момент вычислений.A = zeroMatrix ( maxRows , maxRows )   % Вычислить верхний левый элемент матрицы( 1 , 1 ) = Трапециевидный ( е , Tstart , как правило , ч , у0 )       % Каждая строка матрицы требует одного вызова Trapezoidal% Этот цикл начинается с заполнения второй строки матрицы, так как первая строка была вычислена вышедля i = 1 : maxRows - 1 % Начиная с i = 1, повторять не более maxRows - 1 раз         h = h / 2 % Половина предыдущего значения h, так как это начало новой строки     % Вызовите трапециевидную функцию с этим новым меньшим размером шага ( Я + 1 , 1 ) = Трапециевидный ( е , Tstart , как правило , ч , у0 )           для j = 1 : i % Пересеките строку, пока не достигнете диагонали       % Используйте только что вычисленное значение (например, A (i + 1, j)) и элемент из % строка над ней (т.е. A (i, j)) для вычисления следующей экстраполяции Ричардсона  A ( i + 1 , j + 1 ) = (( 4 ^ j ) . * A ( i + 1 , j ) - A ( i , j )) / ( 4 ^ j - 1 );                конец  % После выхода из указанного выше внутреннего цикла был вычислен диагональный элемент строки i + 1. % Этот диагональный элемент представляет собой последнюю вычисляемую экстраполяцию Ричардсона. % Разница между этой экстраполяцией и последней экстраполяцией строки i является хорошей % индикации ошибки. if ( absoluteValue ( A ( i + 1 , i + 1 ) - A ( i , i )) < толерантность ) % Если результат находится в пределах допуска             print ( "y (5) =" , A ( i + 1 , i + 1 )) % Показать результат экстраполяции Ричардсона        haveWeFoundSolution = true   break % Done, поэтому оставьте цикл  конецконецif ( haveWeFoundSolution == false ) % Если мы не смогли найти решение в пределах желаемого допуска     print ( «Предупреждение: невозможно найти решение в пределах желаемого допуска» , толерантность );  print ( "Последняя вычисленная экстраполяция была" , A ( maxRows , maxRows ))  конец

См. Также [ править ]

• Дельта-квадрат процесс Эйткена
• Такебе Кенко
• Итерация Ричардсона

Ссылки [ править ]

• Методы экстраполяции. Теория и практика К. Брезинского и М. Редиво Загля, Северная Голландия, 1991.
• Иван Димов, Захари Златев, Иштван Фараго, Агнес Хаваси: Экстраполяция Ричардсона: практические аспекты и приложения , Walter de Gruyter GmbH & Co KG, ISBN 9783110533002 (2017). 

Внешние ссылки [ править ]

• Основные методы численной экстраполяции с приложениями , mit.edu
• Ричардсон-Экстраполяция
• Экстраполяция Ричардсона на сайте Роберта Исраэля (Университет Британской Колумбии)
• Код Matlab для обобщенной экстраполяции Ричардсона.
• Код Джулии для общей экстраполяции Ричардсона.