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

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

В численном анализе , Экстраполяция Ричардсона является ускорение последовательности метод , используемая для улучшения скорости сходимости в виде последовательности оценок некоторого значения . По сути, учитывая значение для нескольких значений , мы можем оценить , экстраполируя оценки на . Он назван в честь Льюиса Фрая Ричардсона , который представил эту технику в начале 20-го века, хотя идея была уже известна Христиану Гюйгенсу в его вычислении π . По словам Биркгофа и Роты , «его полезность для практических вычислений трудно переоценить».

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

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

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

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

где и - два различных размера шага.

потом

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

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

Общая формула

Позвольте быть приближением (точное значение), которое зависит от положительного размера шага h с формулой ошибки вида

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

Кроме того, k 0 - это поведение размера шага ведущего порядка ошибки усечения, как

Искомое точное значение может быть выражено как

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

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

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

который можно решить, чтобы дать

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

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

Общее рекуррентное соотношение, начинающееся с, может быть определено для приближений формулой

где удовлетворяет

.

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

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

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

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

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

Следующий псевдокод в стиле MATLAB демонстрирует Экстраполяция Ричардсона , чтобы помочь решить ОДУ , с методом трапеций . В этом примере мы уменьшаем вдвое размер шага на каждой итерации, и поэтому в приведенном выше обсуждении это будет . Ошибка трапециевидного метода может быть выражена в терминах нечетных степеней, так что ошибку на нескольких шагах можно выразить в четных степенях; это приводит нас к возведению во вторую степень и к степеням псевдокода. Мы хотим найти значение , которое имеет точное решение, поскольку точное решение ОДУ равно . Этот псевдокод предполагает, что существует вызываемая функция, которая пытается вычислить , выполняя трапециевидный метод для функции с начальной точкой и размером шага . Trapezoidal(f, tStart, tEnd, h, y0)y(tEnd)fy0tStarth

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

tStart = 0          % Starting time
tEnd = 5            % Ending time
f = -y^2            % The derivative of y, so y' = f(t, y(t)) = -y^2
                    % The solution to this ODE is y = 1/(1 + t)
y0 = 1              % The initial position (i.e. y0 = y(tStart) = y(0) = 1)
tolerance = 10^-11  % 10 digit accuracy is desired

maxRows = 20                % Don't allow the iteration to continue indefinitely
initialH = tStart - tEnd    % Pick an initial step size
haveWeFoundSolution = false % Were we able to find the solution to within the desired tolerance? not yet.

h = initialH

% Create a 2D matrix of size maxRows by maxRows to hold the Richardson extrapolates
% Note that this will be a lower triangular matrix and that at most two rows are actually
% needed at any time in the computation.
A = zeroMatrix(maxRows, maxRows)

%Compute the top left element of the matrix. The first row of this (lower triangular) matrix has now been filled.
A(1, 1) = Trapezoidal(f, tStart, tEnd, h, y0)

% Each row of the matrix requires one call to Trapezoidal
% This loops starts by filling the second row of the matrix, since the first row was computed above
for i = 1 : maxRows - 1 % Starting at i = 1, iterate at most maxRows - 1 times
    h = h/2             % Halve the previous value of h since this is the start of a new row.
    
    % Starting filling row i+1 from the left by calling the Trapezoidal function with this new smaller step size
    A(i + 1, 1) = Trapezoidal(f, tStart, tEnd, h, y0)
    
    for j = 1 : i     % Go across this current (i+1)-th row until the diagonal is reached
        % To compute A(i + 1, j + 1), which is the next Richardson extrapolate, 
        % use the most recently computed value (i.e. A(i + 1, j)) and the value from the
        % row above it (i.e. A(i, j)).
     
        A(i + 1, j + 1) = ((4^j).*A(i + 1, j) - A(i, j))/(4^j - 1);
    end
    
    % After leaving the above inner loop, the diagonal element of row i + 1 has been computed
    % This diagonal element is the latest Richardson extrapolate to be computed.
    % The difference between this extrapolate and the last extrapolate of row i is a good
    % indication of the error.
    if (absoluteValue(A(i + 1, i + 1) - A(i, i)) < tolerance)  % If the result is within tolerance
        print("y = ", A(i + 1, i + 1))                         % Display the result of the Richardson extrapolation
        haveWeFoundSolution = true
        break                                                  % Done, so leave the loop
    end
end

if (haveWeFoundSolution == false)   % If we were not able to find a solution to within the desired tolerance
    print("Warning: Not able to find solution to within the desired tolerance of ", tolerance);
    print("The last computed extrapolate was ", A(maxRows, maxRows))
end

Смотрите также

Рекомендации

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

Внешние ссылки