Отправляет email-рассылки с помощью сервиса Sendsay

Всё о работе в Интернет

  Все выпуски  

Занятие 57. Вычисление экспоненты посредством разложения её в степенной ряд.


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

Объявление. С наступающим Новым годом вас, уважаемые подписчики!

ЭКСПОНЕНТА

Экспонентой называют функцию y = ex, в которой основанием степени является иррациональное число = 2,7182818. Обозначение e для этого числа было введено членом Петербургской Академии наук математиком Леонардом Эйлером (1707–1783), швейцарцем по происхождению.

Значения этой функции можно вычислять с помощью её разложения в степенной ряд

ex = 1 + x / 1! + x2 / 2! + x3 / 3! +  + xi / i! +  .

Теоретически этот ряд можно использовать для вычисления значений экспоненты при любых значениях x.

Легко обнаружить, что первый член этого ряда можно получить при = 0, так как, по определению, x0 = 1 и 0! = 1. Именно это позволяет использовать в качестве общего члена ряда формулу Ai = xi / i!, = 0, 1, 2, .  

Проанализируем значения членов ряда. Так как факториал растёт очень быстро, то члены ряда столь же быстро убывают по модулю, даже не смотря на то, что в числителе находится также достаточно быстро возрастающая степенная функция (например, при x = 10 пятый член ряда ещё велик и приближённо равен 105 / 120 » 103, но уже сотый член практически исчезающе мал: 10100 / 100! » 10100 / 10158 » 10–58).

Быстрое убывание значений членов ряда приводит к тому, что достаточно большая сумма, накопленная за счёт “ближних” членов последовательности, мало изменяется при добавлении к ней “дальних” членов. И хотя теоретически сумма всё же возрастает, то практически ситуация зависит от того, каким числом разрядов ограничена запись числа. Предположим, например, что запись числа ограничена пятью разрядами. Попробуем в таких условиях первоначальное значение некоторой суммы S = 12345 увеличить на величину слагаемого A = 0,1. Вполне очевидно, что “теоретическое” значение новой суммы S = 12345,1, ограниченное пятью разрядами, не изменяется и остаётся прежним S = 12345.

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

Задача C.2.4. “Экспонента”.  Последовательность задана формулой общего члена Ai = xi / i!, = 0, 1, 2, , где x – заданное число. Построить подпрограмму вычисления суммы членов этой последовательности с наибольшей возможной точностью. Определить также, сколько слагаемых для этого потребовалось.

   Procedure Exponent (x:Extended; Var S:Extended; Var N:LongInt );
      Var A,Sp: Extended;
   Begin
      S := 0; N := 0; A := 1;
      Repeat
         Sp := S; S := S + A;
         Inc(N); A := A * x / N;
      Until Sp = S;
   End;

В подпрограмме C.2.4 (процедура Exponent) вычисление суммы членов заданной числовой последовательности организовано на основе циклического процесса, критерием окончания которого является равенство значений новой и предыдущей частичной сумм. При этом наиболее удобным оказалось использовать оператор цикла с постусловием (начальное значение частичной суммы S вычисляется перед циклом и запоминается как Sp уже в цикле, следующее значение частичной суммы S вычисляется тоже в цикле, и только после этого выполняется их сравнение).

Для вычисления членов последовательности применяется искусственное рекуррентное соотношение Ai = Ai–1  x / i при > 0, где Ai = 1 при i = 0. Это соотношение несложно получить, взяв отношение Ai / Ai–1 соседних членов последовательности.

Основные особенности подпрограммы C.2.4 таковы:

1) подпрограмма вычисляет значение экспоненты ex; теоретически никаких ограничений на значение переменной x не накладывается; однако, таковые могут проявиться в процессе вычислений, поскольку значения чисел, представляемых в компьютере, ограничены, а функция ex быстро растёт;

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

3) перед началом цикла устанавливаются начальные значения для частичной суммы S:=0 и для количества N:=0 слагаемых, а также в соответствии с искусственным рекуррентным соотношением вычисляется начальный член последовательности A:=1;

4) в цикле с помощью переменной Sp в первую очередь запоминается предыдущее значение частичной суммы Sp:=S, а затем, поскольку член последовательности A уже известен, вычисляется новое значение частичной суммы S:=S+A;

5) далее в цикле увеличивается на единицу счётчик количества слагаемых (он же порядковый номер члена последовательности) Inc(N) и вычисляется очередной член последовательности в соответствии с искусственным рекуррентным соотношением A:=A*x/N;

6) условием завершения цикла является равенство предыдущей и новой частичных сумм Sp=S.

Результаты тестирования подпрограммы C.2.4 таковы.

При x = 50 применение стандартной функции Паскаля даёт результат Exp(50)=5.18470552858707E+21. Подпрограмма даёт точно такое же значение при количестве слагаемых N=128. Если же вместо типа Extended использовать тип Real, то, очевидно, за счёт уменьшения количества слагаемых (N=109) результат подпрограммы получается несколько хуже: e50 = 5.18470552863559E+21.

Попытки вычислить значение экспоненты при очень больших значениях аргумента x³11400 не удаются вследствие возникновения переполнения. Значения результата получаются настолько большим, что превышают возможности типа Extended. Указанное граничное значение установлено приближённо, и при необходимости его следует уточнить.

Для x = –10 получаем e–10 = 4.53999297623457E–5 при N=64. Стандартное значение Exp(10)=4.53999297624849E5. Появилось незначительное отличие.

Для x = –50 получаем e50 = 5.66922128275782 при N=173. В то же время стандартное значение Exp(50)=1.92874984796392E22. Отличие полученного результата от стандартного очень большое.

Для x = 30 получаем e30 = –5.73715361616713E8 при N=130. Этот результат абсурден вообще, так как значения функции ex не могут быть отрицательными в принципе.

Вывод. Применять процедуру Exponent следует только при неотрицательных значениях аргумента x.

В чём же причины неудовлетворительных результатов работы подпрограммы C.2.4 при вычислении значений функции ex для x < 0?

Внимательное рассмотрение членов суммируемой последовательности приводит к следующему выводу: если x >= 0, то суммируемые члены последовательности только положительны, если же x < 0, то среди суммируемых членов последовательности есть не только положительные, но и отрицательные. Следовательно, причиной потери точности есть вычитание близких по величине чисел.

Выход из этой ситуации довольно прост. Если x >= 0, то обращаться к подпрограмме C.2.4 можно непосредственно. Если же x < 0, то обращаться к ней нужно со значением аргумента, равным x, а в качестве результата брать обратную величину (на основании свойства степени a–b = 1 / ab).

Всё вышесказанное реализовано в функции Universal_Exponent. Её тестирование показало абсолютное совпадение с результатами использования стандартной функции Exp(x).

   Function Universal_Exponent ( x: Extended ): Extended;
      Var y: Extended; N: LongInt;
   Begin
      If x >= 0 Then Exponent (x,y,N)
                Else Begin Exponent(-x,y,N); y:=1/y End;
      Universal_Exponent := y
   End;

На следующих занятиях, уважаемые подписчики, мы рассмотрим ряд других циклических алгоритмов.

Уважаемые подписчики! При необходимости задать вопрос, проконсультироваться, уточнить или обсудить что-либо обращайтесь через Гостевую книгу моего персонального сайта http://a-morgun.narod.ru. При этом настоятельно рекомендую пользоваться браузером Internet Explorer.

С уважением, Александр.

 


В избранное