Алгоритм Ромберга
После сделанных вводных замечаний приведем основные идеи итерационного алгоритма Ромберга, который применяется в системе Mathcad для выполнения операции численного интегрирования.
Сначала строится несколько интерполирующих полиномов, которые заменяют на интервале интегрирования подынтегральную функцию f (х). В качестве первой итерации полиномы вычисляются по 1, 2 и 4 интервалам. Например, как уже отмечалось выше, первый полином, построенный по 1 интервалу, – это просто прямая линия, проведенная через две граничные точки интервала интегрирования, второй – квадратичная парабола и т. д.
Интеграл от каждого полинома с известными коэффициентами легко вычисляется аналитически. Таким образом, определяется последовательность интегралов от интерполирующих полиномов: I1,I2,I4,… Например, по правилу трапеций: I1=(b-a) (f (a) + f (b)) / 2 и т. д.
Из-за интерполяции по разному числу точек вычисленные интегралы I1, I2,… несколько отличаются друг от друга. Причем чем больше точек используется для интерполяции, тем интеграл от интерполяционного полинома ближе к искомому интегралу I, стремясь к нему в пределе бесконечного числа точек. Поэтому определенным образом осуществляется экстраполяция последовательности i1,i2,i4,… до нулевой ширины элементарного интервала. Результат этой экстраполяции j принимается за приближение к вычисляемому интегралу.
Осуществляется переход к новой итерации с помощью еще более частого разбиения интервала интегрирования, добавления нового члена последовательности интерполирующих полиномов и вычисления нового (N -го) приближения Ромберга JN.
Чем больше количество точек интерполяции, тем ближе очередное приближение Ромберга к вычисляемому интегралу и, соответственно, тем меньше оно отличается от приближения предыдущей итерации. Как только разница между двумя последними итерациями | JN -JN-1 1 становится меньше погрешности TOL или меньше TOL|JN |, итерации прерываются, и JN появляется на экране в качестве результата интегрирования.