Maple 9.5/10 в математике, физике и образовании

5.10.3. Паде-аппроксимация

Теперь опробуем рациональную аппроксимацию Паде (Pade) функции f(x) степени (4,4). Приближения, по этому разложению, будут аппроксимировать функцию более точно, и потому ошибки округления в вычислениях станут более заметными. Поэтому зададим вычисления с двумя дополнительными знаками точности:

> Digits := 12:
> s := map(evalf, taylor(f(x), x=2, 9)):
> PadeApprox := pade(s, x=2, [4,4]);

PadeApprox := (0.341034792604 + 0.0327799035348x - 0.00612783638188(x-2)² + 0.00452991113636(x-2)³ - 0.000431506338862(x-2)⁴)/( 0.068484906786 + 0.465757546607x+ 0.159149610837(x-2)² + 0.0266813683828(x-2)³ + 0.00346967791444(x-2)⁴)

> PadeApprox := unapply(PadeApprox, x):

Кривая ошибки для интервала [0,4] строится командой

> plot(f - PadeApprox, 0..4,color=black);

и имеет вид, показанный на рис. 5.26.

Рис. 5.26. Кривая погрешности при Паде-аппроксимации степени (4.4)

Как и при аппроксимации рядом Тейлора, ошибка здесь мала вблизи точки разложения и велика вдали от нее. Мы снова видим из графика, что для указанной функции, самая большая ошибка — в левой оконечной точке. Однако, максимальная ошибка в Паде-аппроксимации уже на порядок меньше, чем при аппроксимации полиномом Тейлора:

> maxPadeError := abs(limit(f(x), x=0) - PadeApprox(0));

maxPadeError:=0.000353777322

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

5.10.4. Аппроксимация полиномами Чебышева

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

Разложим функцию f(x) на [0,4] в ряд Чебышева с точностью 1*10^{- 8}. Это означает, что все члены с коэффициентами меньше, чем эта величина, будут опущены. Такая точность обеспечивается полиномом 13 степени:

> evalf(limit(f(x), х=0));

.500000000000

> fproc := proc(x) if x=0 then 0.5 else evalf(f(x)) fi end:
> ChebApprox := chebyshev(fproc, x=0..4, 1E-8);

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

Для последующих вычислений, полезно заметить, что мы можем использовать процедуру для нахождения численных значений f(x), которая будет намного эффективнее, чем прямое определение, которое требует численного интегрирования для каждого значения х. А именно, определим процедуру численной оценки, основанную на разложении в ряд Чебышева степени 13, так как максимальная ошибка при такой аппроксимации меньше, чем 10^-8, и обеспечивает для нашей цели достаточную точность. Мы определим полином Чебышева Т(х) из пакета orthopoly, и затем для эффективной оценки преобразуем его в форму Горнера:

> F := hornerform(eval(subs(T=orthopoly[T],
 ChebApprox)));

F = 0.499999998610 + (0.192405358503 + (-0.163971754264 + (-0.0083861432817 + (0.0277082269676 + (-0.00593172541573 + (-0.00132728874257 + (0.000910057654178 + (-0.000180351181100 + (0.57685696534 10^-5 + ( 0.448885653549 10^-5 + (-0.990274556116 10^{- 6} + (0.925433855729 10^-7 - 0.347161977631 10^-8x)x)x)x)x)x)x)x)x)x)x)x) x

> F := unapply(F, x):

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

5.10.5. Аппроксимация Чебышева-Паде

Теперь рассмотрим еще более точную рациональную аппроксимацию Чебышева-Паде. Это такая рациональная функция r[m, n](х) с числителем степени m и знаменателем степени n такой же, как и для разложения в ряд Чебышева. Функция r[m, n](х) согласуется с разложения в ряд ряда Чебышева f(x) членом степени m+n. Мы вычислим аппроксимацию Чебышева-Паде степени (4, 4), подобную обычной Паде- аппроксимации, успешно выполненной ранее:

> ChebPadeApprox := chebpade(F, 0..4, [4,4]);

Построим кривую ошибок:

> with(orthopoly, Т):
> plot(F - ChebPadeApprox, 0..4,color=black);

Она представлена на рис. 5.27.

Рис. 5.27. Кривая ошибки при Паде-Чебышева рациональной аппроксимации

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

> maxChebPadeError := abs(F(0) - ChebPadeApprox(0));

maxChebPadeError := 0.1236749 10^-5

Мы достигли впечатляющего успеха и остается сделать еще один шаг в направлении повышения точности аппроксимации.

5.10.6. Минимаксная аппроксимация