6. Циклические операции. Программирование

Простейший вариант организации циклических вычислений – построение рекурсивного процесса – процесса итераций. Возможность построения итерационного процесса была заложена уже в самых первых версиях пакета Mathcad, работавших в системе DOS. Для организации процесса итераций задается начальное значение (или несколько начальных значений) вычисляемой величины и закон, по которому вычисляется каждое новое значение на основе предыдущего (или нескольких предыдущих). Подчеркнем, что Mathcad позволяет организовать итерационный процесс, как для скалярных, так и для векторных величин. Рассмотрим примеры вычислений с использованием процесса итераций.

Пример 6.1.

Числа Фибоначчи – последовательность чисел, каждое из которых (начиная с третьего) равно сумме двух предыдущих. Итак, зададим количество чисел в последовательности (n+1), первые два числа: f₀ и f₁, индекс i – номер очередного числа и формулу, определяющую очередное число. Результатом расчетов является вектор f, включающий девять чисел Фибоначчи.

Пример 6.2.

Построим собственную формулу для расчета факториала числа n:

В соответствии с данной формулой f(n) последовательно записывается в виде произведения:

f(n) = n f(n-1) = n (n-1) f(n-2) = … = n (n-1) (n-2) … 1

(Напомним, что в условии n=0 должен использоваться знак логического равенства, для написания которого нужно нажать клавиши Ctrl+= .)

П ример 6.3.

Напишем еще один вариант датчика чисел Фибоначчи. С помощью этого датчика можно найти число Фибоначчи с произвольным номером или построить вектор, включающий последовательность чисел.

Пример 6.4.

Решение задачи Коши методом Эйлера. Задача Коши состоит в нахождении решения дифференциального уравнения , удовлетворяющего начальному условию . Требуется найти решение уравнения на отрезке . Введем на оси x сетку: , где Будем искать значения неизвестной функции в узлах сетки. Обозначим: В соответствии с алгоритмом Эйлера значения неизвестной функции в узлах сетки вычисляются по формуле Фактически на каждом шаге приращения аргумента нужно вычислить две величины: и . Поэтому алгоритм вычислений на каждом шаге можно описать векторным соотношением:

(6.1)

Найдем решение уравнения на отрезке [0, 1.5] с шагом h=0.1 и с начальным условием .

Программа расчетов в среде Mathcad показана на (рис. 6.1. Первая строка программы содержит исходные данные. Вторая строка включает вектор начальных данных Y и вектор-функцию F(V), входящую в правую часть соотношения (6.1). Само это соотношение записано в третьей строке.

Р езультатом работы программы является матрица Y, строки которой содержат значения x и y в узлах сетки. Известно, что точным решением данной задачи Коши является экспоненциальная функция . Для сравнения найденных приближенных значений с точными добавим в матрицу Y строчку, включающую значения экспоненты. Итоговая матрица Y содержит в первой строке значения x, во второй – приближенные, а в третьей – точные значения y.

Ф ормулы, описывающие метод Эйлера, записаны в общем виде. Для того чтобы применить эту программу к другой задаче, достаточно изменить определения n, h, f(x,y)и вектор начальных данных Y. Очевидно придется также изменить и точное решение.

Подчеркнем, что в приведенной программе итерационный процесс описан в векторной форме. Изменив описания f, F и Y, можно применить эту программу к системе дифференциальных уравнений. Размерность векторов в этом случае будет равной N+1, где N – число уравнений.

Существенным недостатком организации итерационных процессов в системе Mathcad является невозможность создания вложенных циклов. Для организации вложенных циклов нужно использовать механизм программирования. Операторы, используемые при создании программ, включены в панель инструментов Programming, показанную на рисунке. Клавишные эквиваленты команд приведены в таблице 6.1.

Таблица 6.1. Команды, используемые при создании программ

Название оператора	Команда	Название оператора	Команда
Добавление строки (add line)	]	Локальное присваивание	{
if	}	otherwise	Ctrl+}
for	Ctrl+”	while	Ctrl+]
break	Ctrl+{	continue	Ctrl+[
return	Ctrl+|	on error	Ctrl+’

Программы, написанные на языке Mathcad, имеют своеобразный вид. Вертикальные линии выделяют вложенные группы операторов. Рассмотрим построение программ на примерах.

Пример 6.5.

На рис. 6.2 приведена программа вычисления квадратного корня с помощью процесса Герона. В соответствии с алгоритмом Герона задается начальное приближенное значении корня y₀>0, а затем это исходное значение уточняется в процессе итераций; при этом вычисления проводятся по формуле .

Программа определяет функцию rt (x,), где x – число, из которого извлекается корень,  – требуемая точность. К данной функции можно обращаться как к любой встроенной функции. Обращение к функции возможно справа от прямоугольной области, заключающей программу, и во всех последующих строках. Если аргументы функции заданы неправильно, вычисления прерываются. Если в этом случае щелкнуть мышкой на обращении к функции, появится сообщение об ошибке. Примеры обращения к функции rt(x,) приведены на рисунке. На примере вычисления квадратного корня из четырех видим, что требуемая точность обеспечивается с запасом.

П ример 6.6.

На рисунке 6.3 приведена программа нахождения наибольшего по модулю элемента вектора.

Пример 6.7.

Составим собственную программу перемножения матриц по правилу “строка на столбец”. Программа определяет функцию PM(A,B). Программа и примеры обращения к функции M(A,B) приведены на рис. 6.4. В начале программы проверяется совпадение числа столбцов первой матрицы и числа строк – второй. При наличии совпадения вычисляется произведение. Для этого используются вложенные циклы for. Если совпадения нет, вычисления прерываются.

Обращение к функции M(A,B) доступно в строке, содержащей основной знак присваивания (M(A,B):=), а также во всех последующих строках. В случае, когда размеры матриц не согласованы, можно получить сообщение об ошибке, щелкнув мышкой на обращении к функции.

Пример 6.8.

Адаптивный квадратурный алгоритм, описанный в справочной системе Mathcad.

Алгоритм включает две программы, приведенные на активной вставке – рис. 6.5. Первая программа Simpson(f,a,b,N) вычисляет определенный интеграл по формуле Симпсона. Здесь f –

предварительно определенная подинтегральная функция, [a,b] – отрезок интегрирования, N – количество подотрезков, на которые разбивается отрезок [a,b]. Обозначим: h=(b-a)/N – длина подотрезков, x₀=a, x_i=x₀+h*i. В этих обозначениях формула Симпсона имеет вид:

Точность вычислений растет с увеличением числа N, т.е. с уменьшением длины подотрезков h.

Вторая программа Adapt(f,a,b) использует программу Simpson и служит для вычисления определенных интегралов с заданной точностью (в программе установлено значение допустимой погрешности eps=10^-8). Вначале вычисляется два значения интеграла: при делении исходного отрезка [a,b] на 5 и на 10 подотрезков. Разность между найденными значениями служит оценкой погрешности. Если эта разность превышает значение eps, длина подотрезков уменьшается вдвое, и вычисления повторяются. Весь процесс повторяется до тех пор, пока не будет достигнута требуемая точность.

Н а рисунке приведены также результаты вычислений интеграла . Точное значение этого интеграла равно . Значение интеграла, найденное с помощью программы Adapt, совпадает с числом  вплоть до пятнадцатого знака после запятой. Повторное обращение к программе Simpson при возрастающих значениях N позволяет видеть процесс последовательных приближений к точному значению интеграла. Видно, что для данного интеграла уже при N=3 требуемая точность обеспечивается с запасом.

Пример 6.8.

Напишем еще одну программу вычисления машинного эпсилон: