Программирование в MathCAD

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

И, тем не менее, встроенный язык программирования MathCADпозволяет решать самые различные, в том числе и довольно сложные, задачи и является серьезным подспорьем для расчетов. Поэтому последние версии MathCAD имеют не очень мощный, но весьма удобный собственный язык программирования.

Для вставки программного кода в документы MathCAD имеется специальная панель инструментов

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

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

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

 

Внутри программы можно использовать глобальные переменные документа, но изменить их значение внутри программы никак нельзя. Можно создать в программе другие переменные, доступ к которым может осуществляться только из самой программы. Эти переменные называются локальными переменными. Локальные переменные «не видны» извне. Локальная переменная создается с помощью знака локального присвоения с панели Programming. Для оператора локального присваивания, так же как и для операторов обычного и глобального присваивания, можно изменить внешний вид так, чтобы он выглядел как обычный знак равенства. Для этого достаточно вызвать контекстное меню этого оператора и в нем выбрать команду ViewDefinitionAs/Equal.

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

Вставить строку программного кода в уже созданную программу можно в любой момент с помощью той же самой кнопки AddLine. Для этого следует поместить на нужное место внутри программы линии ввода.

Если вертикальная линия ввода находится в начале текущей строки, то нажатие кнопки приведет к появлению новой строки перед текущей строкой, если вертикальная линия ввода находится в конце строки, то нажатие кнопки приведет к появлению новой строки после текущей. Можно выделить только часть текущей строки, что повлияет на положение новой строки в программе (например, если необходимо для некоторого условия выполнить не одно действие, а несколько). Таким образом, основной принцип создания программных модулей заключается в правильном расположении строк кода. Ориентироваться в их действии довольно легко, т.к. фрагменты кода одного уровня сгруппированы в программе с помощью вертикальных черт.

Проверка условий в программах.

Программы в MathCAD могут быть не только линейными, но и разветвленными. Одним из вариантов ветвления в программах является проверка условия. Условия могут проверять значения как локальных, так и глобальных переменных, а также выражений, содержащих эти переменные.

Для проверки условий в программах MathCAD служит оператор if . В поле ввода справа нужно ввести условие. Для ввода условий служит панель Boolean. На этой панели есть кнопки, предназначенные для проверки условий , а также кнопки, предназначенные для вставки логических операций (логическое отрицание, операция «и», операция «или», операция «исключающее или»). Вторые позволяют создавать сложные условия. В поле ввода слева нужно ввести строку программы, которая должна выполняться, если введенное условие истинно.

Если для условия «истинно» необходимо выполнение нескольких строк программы, надо воспользоваться кнопкой

.

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

Следует помнить, что если в программе введено подряд несколько строк с оператором if

, то выражение слева от оператора otherwise будет выполнено только в том случае, если не выполняются условия, заданные во всех операторах if.

Примеры:

Создание циклов.

Другим очень важным элементом при создании нелинейных программ являются циклы. Циклы позволяют повторять несколько раз выполнение одного и того же программного блока. Для создания циклов в MathCAD предусмотрено два оператора for и while.

Цикл for – цикл со счетчиком.

В таких циклах создается некоторая переменная-счетчик, значение которой изменяется после каждого выполнения тела цикла. Выход из цикла происходит при достижении этой переменной заданного значения. Этот оператор вводится с панели Programming клавишей for . В поле ввода после слова for следует указать имя переменной – счетчика. Это может быть любое имя, которое не использовалось ранее в программе. Внутри цикла можно использовать эту переменную в любых выражениях, нельзя только присваивать ей никакого значения. В поле ввода после знака следует указать диапазон значений переменной-счетчика. Вводить диапазон в данном случае следует так же, как и при создании ранжированной переменной. Вместо диапазона в данном поле ввода можно указатьимя некоторого массива (вектора или матрицы). В этом случае переменная-счетчик будет последовательно принимать значения всех элементов этого массива. Возможность перебора элементов массива не может быть реализована с помощью цикла while, поэтому именно в таких случаях цикл for и является незаменимым. В поле ввода под словом for следует ввести тело цикла.

Пример: Заполнить вектор числами от xнач до xкон с шагом h. Затем определить сумму элементов этого вектора и найти их среднее арифметическое значение.

Цикл while – цикл, который выполняется до тех пор, пока выполняется определенное условие.

В поле ввода справа от слова while следует ввести условие. Это условие строится по тем же правилам, что и в операторе if. Оно будет проверяться после каждого выполнения тела цикла и в тот момент, когда условие перестанет выполняться, повторение тела цикла прекратится. В поле ввода ниже слова while следует ввести тело цикла (для ввода нескольких строк в теле цикла надо воспользоваться кнопкой AddLine).

Пример 1: Вычислить сумму с точностью е.

Пример2: Найти первый элемент, превышающий определенный порог.

Использование операторов break и continue.

Иногда возникает необходимость повлиять на выполнение цикла некоторым образом, например, прервать его выполнение по какому-либо условию или выполнять некоторые итерации не так, как другие. Для этого и служат операторы break и continue.

Оператор break, если он расположен внутри цикла, означает немедленное прекращение выполнения текущей итерации и выход из цикла. Если есть необходимость прекратить выполнение цикла по какому-либо условию, то следует использовать конструкцию следующего вида breakif (условие).

Пример: Выделить из массива все элементы от начала и до первого вхождения в него заданного числа.

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

Пример: Требуется заполнить элементы квадратной матрицы в шахматном порядке.

Оператор return (возврат значения).

Как было сказано ранее, результат выполнения программного модуля помещается, как правило, в последней его строке. Но можно прервать выполнение программы в любой ее точке (например, с помощью условного оператора) и выдать некоторое значение, применив оператор return. В этом случае при выполнении указанного условия значение, введенное в поле ввода после return, возвращается в качестве результата, и никакой другой код больше не выполняется. Вставляется в программу оператор return с помощью одноименной кнопки панели Programming . Пример:

Обработка ошибок.

Система MathCAD предоставляет пользователю некоторый контроль над ошибками, которые могут возникнуть при вычислении выражений или при выполнении программ. Для этой цели служит оператор onerror, который можно вставить с помощью одноименной кнопки панели Programming . В поле ввода справа следует ввести выражение или программу, которые необходимо вычислить (известно, что это выражение может содержать ошибку при определенных значениях входных параметров). В поле ввода слева следует ввести выражение, которое будет выполнено вместо правого выражения, если при выполнении последнего возникнет ошибка. Пример: Если аргументу функции присвоено нулевое значение, то в программе возникает ошибка – деление на нуль. Но за счет оператора onerror сообщение не выводится, а функции в этой точке присваивается значение, указанное слева от оператора onerror – значение машинной бесконечности.

В поле ввода слева может быть введено текстовое выражение, сообщающее об ошибке

Конечно, если программа содержит только простейшие математические операции, то никаких ошибок, кроме деления на нуль или превышения наибольшего допустимого числа в ней возникнуть не может. В таких программах оператор onerror используется редко. Но для более сложных программ, которые содержат функции решения дифференциальных уравнений, аппроксимации или другие сложные функции MathCAD, использование onerror может избавить от многих трудностей.

Иногда может возникнуть ситуация обратная той, которая была описана выше, т.е. необходимо, чтобы при определенных условиях результатом выражения было сообщение об ошибке, хотя в действительности при этом не возникает ни одной стандартной ошибки MathCAD. Для таких случаев в MathCAD предусмотрена встроенная функция error. В качестве аргумента этой функции нужно в кавычках указать текст сообщения об ошибке, который должен быть выведен. Таким образом, если необходимо, чтобы программа возвращала ошибку при определенном условии, то следует использовать конструкцию вида: error («текст ошибки») if (условие).

Пример:

Для того, чтобы иметь возможность нормально вводить текст на русском языке в аргумент функции error (а также во все другие функции со строками), следует изменить шрифт, который используется во встроенном стиле Constant. Для того чтобы этот стиль правильно отображал русские буквы, установите курсор на любом числе или строковом выражении в формульном блоке. При этом в поле на панели инструментов Formatting, отображающем текущий стиль, должно быть написано – Constant. Теперь выберите из раскрывающегося списка шрифтов шрифт, поддерживающий кириллицу.

Примеры программирования.

1. Вычислить функцию sin(x) с точностью е.

2. Даны массивы А(5) и В(5). Получить массив С, в который записаны сначала элементы массива А в порядке возрастания, а затем элементы массива В порядке убывания.

3. По введенным значениям коэффициентов А, В, С определить корни квадратного уравнения

4. Дан массив натуральных чисел В(10). Определить, есть ли в нем 4 последовательных числа (например, 1, 2, 3, 4, и т.п.). Напечатать такие группы чисел.

Вопросы

1. Какая панель служит для вставки программного кода в документ MathCAD? Можно ли операторы программирования набрать с клавиатуры?

2. С какой команды начинается создание программного блока? Как с ее помощью можно создавать разветвленный программный блок?

3. Что такое определение программного блока? Обращение к программному блоку?

4. Что такое глобальные и локальные переменные для программного блока? Что может содержать последняя строка программного блока?

5. Как работает оператор if в программном блоке? Приведите пример.

6. Создание цикла с параметром в программном блоке. Приведите пример.

7. Создание цикла while в программном блоке. Приведите пример.

8. Для чего служат операторы break, continue в программном блоке? Приведите примеры.

9. Как работает оператор return в программном блоке? Приведите пример.

10. Как осуществляется обработка ошибок в программном блоке? Приведите пример.

Лекция №3

Численные методы решения задач.

Обработка экспериментальных данных средствами MathCAD

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

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

Численными методами часто приходится решать следующие математические задачи:

1. решение нелинейных (алгебраических и трансцендентных) уравнений;

2. вычисление определенных интегралов;

3. решение обыкновенных дифференциальных уравнений;

4. решение дифференциальных уравнений в частных производных;

5. решение задач оптимизации;

6. обработка массивов числовых данных.

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

Решение нелинейных уравнений

Обычно нелинейные уравнения делят на трансцендентные и алгебраические. Нелинейные уравнения, содержащие тригонометрические функции или другие специальные функции, например, lg(x) или ex, называются трансцендентными. Методы решения нелинейных уравнений такого типа делятся на аналитические и численные.

Первые позволяют найти решение непосредственно с помощью формул и всегда обеспечивают получение точного решения. Известным примером такого рода является формула корней квадратного уравнения. В численных методах задается процедура решения в виде многократного применения некоторого алгоритма. Задача отыскания корней нелинейного уравнения f(x) = 0 считается решенной, если мы сумеем определить корни с нужной степенью точности.

Для решения нелинейных уравнений известны следующие численные методы: метод половинного деления (метод дихотомии), метод хорд, метод касательных (Ньютона), метод секущих, метод простой итерации. Рассмотрим метод половинного деления.

Графическая интерпретация метода показана на рис.1.

Рис.1 Графическая интерпретация метода половинного деления

В этом методе отыскание корня уравнения f(x) = 0 проходит в два этапа. На первом этапе необходимо отделить корень, т.е.выделить интервал на оси абсцисс, на котором функция f(x) меняет свой знак. Для отделения корня следует провести вычисление функции f(x) в точках, расположенных через равные интервалы по оси x, до тех пор, пока не будут найдены два последовательных значения функции f(xn) и f(xn+1), имеющие противоположные знаки.

На втором этапе производится уточнение корня. Найденный интервал [ xn, xn+1 ], содержащий корень, делится пополам

Затем по разности знаков функции на концах интервала определяем, на каком из полученных двух интервалов находится корень уравнения. Найденный интервал снова делится пополам и т.д. В результате, интервал, на котором находится корень, сужается. Процесс повторяется до тех пор, пока f (x ср) не станет достаточно близким к нулю. Блок-схема алгоритма метода показана на рис.2.

Рис.2 Блок-схема алгоритма метода половинного деления

 
 


Численное интегрирование

К численному интегрированию обращаются, когда нельзя через элементарные функции аналитически записать первообразную интеграла

или когда подобная запись имеет сложный вид.

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

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

1. Методы Ньютона - Котеса. Эти методы требуют, чтобы значения x были заданы с постоянным шагом. Они основаны на полиномиальной аппроксимации подынтегральной функции. Алгоритмы методов просты и легко поддаются программной реализации.

2. Методы Гаусса – Кристоффеля – методы наивысшей алгебраической точности. Эти методы используют неравно отстоящие узлы, расположенные по алгоритму, обеспечивающему минимальную погрешность интегрирования. Требуют большего объема памяти, чем методы первой группы.

3. Методы Монте-Карло. В этих методах узлы выбираются с помощью датчика случайных чисел. Ответ носит вероятностный характер. Методы эффективны при вычислении большой кратности.

 
Рассмотрим методы первой группы. По этому методу интервал разбивается на n равных отрезков, длина каждого из которых . Для вычисления каждой элементарной площади подынтегральную функцию на данном отрезке заменяют с некоторой степенью точности более простой функцией, интеграл от которой можно вычислить, используя только ординаты на концах отрезка. Метод Ньютона -Котеса дает три формулы для приближенного вычисления определенного интеграла

1. При замене f(x) на каждом отрезке прямой, параллельной оси x (рис.3,а), получим формулу прямоугольников

2. При замене f(x) на каждом отрезке прямой, соединяющей ординаты концов отрезка (рис.3,б), получим формулу трапеций

3. При замене f(x) дугой параболы, проведенной через концы трех ординат (рис.3,в), получим формулу Симпсона. При этом число отрезков должно быть четным

Точность вычисления интеграла тем выше, чем больше n и меньше h, и в пределе при n ®¥ h ® 0 указанные формулы дадут точную величину определенного интеграла.

Рис.3 Численное определение интеграла

Блок схема алгоритма численного определения интеграла с заданной степенью точности представлена на рис. 4. На блок-схеме a,b,e,n – исходные данные, е – заданная точность, n – число разбиений, i1 – начальное значение интеграла, его можно задать равным нулю, i- вычисленное значение интеграла.

Для достижения заданной точности число разбиений удваивают до тех пор, пока будет удовлетворено условие |i – i1|<e.< p=""></e.<>

Рис.4 Блок схема алгоритма численного определения интеграла с заданной степенью точности по формуле прямоугольников

Обработка экспериментальных данных средствами MathCAD.

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

Для анализа экспериментальных зависимостей и приведения их к наглядному виду часто используются такие средства, как интерполяция, сглаживание и аппроксимация. В некоторых задачах для анализа зависимости бывает необходимо найти ее Фурье-компоненту. Все перечисленные средства анализа можно реализовать средствами MathCad.

Интерполяция.

Интерполяцией называют заполнение отрезками кривых промежутков между заданными точками по тому или иному закону. Для проведения интерполяции в первую очередь должна быть задана экспериментальная зависимость в виде набора точек на плоскости. Для этого должны быть заданы два одномерных массива (вектора) – vx и vy, содержащие соответственно значения координат x и y каждой точки. При этом важно, чтобы значения в векторе vx были заданы в порядке возрастания.

Система MathCad позволяет проводить линейную интерполяцию и сплайн-интерполяцию наборов экспериментальных точек.

Простейшим вариантом интерполяции является линейная интерполяция. Она заключается в простом соединении точек между собой отрезками прямых. Для реализации такой интерполяции в MathCad существует встроенная функция linterp(vx,vy,x), где vxvy – уже известные векторы, содержащие координаты последовательности точек, x – координата точки, в которой нужно вычислить значение интерполирующей функции. Пример построения линейной интерполяции приведен на рис.5

Рис.5 Линейная интерполяция

На практике линейная интерполяция применяется редко.

Из всех видов интерполяции наиболее часто используется интерполяция, где экспериментальные точки попарно соединяются отрезками полиномов. Чаще всего для этого выбирают полиномы третьей степени (поэтому такая кривая и называется кубическим сплайном). Для того чтобы найти коэффициенты этих полиномов, очевидно, недостаточно того условия, что кривая должна проходить через экспериментальные точки. Поэтому на сплайн накладываются дополнительные условия сшивки – первая и вторая производные слева и справа от каждой экспериментальной точки должны быть равны между собой. Но и после этого количество условий остается на два меньше, чем количество неизвестных коэффициентов. Дополнительные два условия должны быть наложены в начальной и конечной экспериментальных точках, поскольку в них нет условий сшивки. Эти условия можно выбрать по-разному. В MaqthCad существуют три различных функции для построения кубических сплайнов с различными дополнительными условиями.

· lspline(vx,vy) – в начальной и конечной точках накладывается условие линейности, т.е. вторая производная от функции равна нулю. Первая буква в названии функции – l, означает linear (линейный).

· pspline(vx,vy) – на первом и последнем интервале кривая является параболой, т.е. полиномиальный коэффициент при x3 равен нулю. Буква p означает parabolic (параболический).

· cspline(vx,vy) – полиномиальные коэффициенты при x3 на первых двух интервалах равны между собой точно так же, как на последних двух интервалах. Буква c означает cubic (кубический).

Результатом каждой из перечисленных функций является вектор, содержащий значения вторых производных от интерполяционной кривой во всех точках, заданных в массиве vx. Для того чтобы исходя из этого вектора построить кривую, нужно воспользоваться встроенной функцией interp(v,vx,vy,x), где vx и vy – массивы экспериментальных точек, v – массив, полученный как результат одной из трех функций, перечисленных выше, x- координата, в которой нужно вычислить значение интерполяционной кривой. Пример интерполяции кубическим сплайном приведен на рис.6.

Рис.6 Интерполяция кубическим сплайном

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

Кубический сплайн является эффективным средством построения интерполяционной кривой в подавляющем большинстве случаев. Но иногда использование кубического сплайна может привести к нежелательным результатам. Чаще всего это происходит в тех случаях, когда данные очень неравномерно распределены вдоль оси x. В таких случаях на кривой кубического сплайна могут появляться острые экстремумы в промежутках между экспериментальными точками. В некоторых подобных случаях получить лучщую интерполяционную кривую помогает использование другого вида интерполяции – В-сплайна. Основное отличие В-сплайна от всех описанных выше методов – сшивка отрезков кривых происходит не в экспериментальных точках, а между ними, в специально заданных точках. В MathCad для реализации интерполяции В-сплайном служит функция bspline(vx.vy,u,n), где vx,vy – векторы, содержащие координаты экспериментальных точек, u – вектор, содержащий координаты точек сшивки, n – порядок полинома. Результатом функции bspline является вектор, который далее следует использовать как аргумент функции interp. В-сплайн в MathCad можно построить из отрезков прямых, парабол или кубических парабол, т.е. допустимые значения параметра n – 1,2 или 3. Количество точек сшивки не является произвольной величиной и должно быть всегда на n-1 меньше, чем количество экспериментальных точек. Также на координаты точек сшивки накладываются следующее условие: первая точка сшивки должна быть не правее первой экспериментальной точки, а последняя не левее последней экспериментальной точки. Остальные точки сшивки могут произвольным образом располагаться внутри отрезка. Пример использования В-сплайна приведен на рис.7.

Экстраполяция.

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

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

Рис.7 Сравнение эффективности кубического сплайна и В-сплайна

Также в MathCad есть одна встроенная функция для проведения экстраполяции – predict(v,m,n), где v – вектор значений функции на том отрезке, где она известна (значения аргумента в данной функции не задаются, и считается, что точки распределены равномерно), m – количество элементов вектора v, на основании которых проводится экстраполяция (естественно, выбираются точки, ближайшие к правой границе), n – количество точек в рассчитываемом векторе. Пример экстраполяции от осциллирующей функции приведен на рис.8

k,101+i

Рис.8 Экстраполяция с помощью функции predict

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

Регрессия

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

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

Линейная регрессия

Линейная регрессия является наиболее простой, но, тем не менее, используется чаще любого другого вида регрессии. Она заключается в нахождении таких значений параметров a и b, чтобы прямая y = a+bx наилучшим образом аппроксимировала заданной набор точек. Для проведения линейной регрессии по методу наименьших квадратов в MathCad существует функция line(vx,vy). Результатом функции line будет вектор, содержащий значения параметров a и b для построения регрессионной прямой. Пример линейной регрессии представлен на рис.9

Рис.9 Линейная регрессия с помощью функции line

На рис.9 для анализа эффективности линейной регрессии построена псевдоэкспериментальная последовательность точек. Для этого к точным значениям линейной функции прибавлены случайные числа, сгенерированные с помощью функции rnorm.Затем проведена линейная регрессия полученного набора точек. Ее результат можно сравнить с исходной прямой. Как видно из рис., параметры регрессионной прямой заметно отличаются от исходных, но на графике обе прямые проходят достаточно близко.

Для линейной регрессии в MathCad реализован также метод медиан с помощью функции medfit(vx,vy). Результатом этой функции является вектор, аналогичный результату line. Нельзя утверждать, что один из двух методов регрессии более точен. Метод наименьших квадратов является наиболее универсальным, поэтому функция line считается в MathCad основной функцией для проведения линейной регрессии.

Полиномиальная регрессия.

Кроме аппроксимации неизвестной функции с помощью прямой, широкое применение находит и аппроксимация с помощью полиномов различной степени. Для этой цели в MathCad существует функция regress(vx,vy,n). Последний аргумент данной функции задает степень полинома. Можно использовать полином любой степени, но не большей, чем число точек в выборке минус один. При n = 1 получится линейная регрессия. На практике наибольшее применение находит полиномиальная регрессия от второй до пятой степени.

Результатом функции regress является вектор, содержащий коэффициенты аппроксимирующего полинома. Эти коэффициенты располагаются в векторе, начиная с четвертого элемента в порядке возрастания степеней. Первые три элемента данного вектора являются служебными и используются для того, чтобы результат функции regress можно было использовать как первый аргумент функции interp по аналогии со сплайн-интерполяцией. Пример полиномиальной регрессии представлен на рис.10.

Рис.10 Полиномиальная регрессия

На рис.10 для иллюстрации полиномиальной регрессии построена псевдоэкспериментальная последовательность точек. В качестве теоретической функции был использован полином третьей степени с коэффициентами 0,1,-2,1. Как видно из примера, коэффициенты, рассчитанные функцией regress, значительно отличаются от коэффициентов исходного полинома. Тем не менее в области экспериментальных точек обе кривые достаточно близки, но за пределами этой области резко расходятся.

Аппроксимация набора точек различными элементарными функциями

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

· expfit(vx,vy,vg) – регрессия экспоненциальной функцией y = a*eb*x+c.

· sinfit(vx,vy,vg) – регрессия синусоидальной функцией y = a*sin(x+b)+c.

· pwrfit(vx,vy,vg) - регрессия степенной функцией e = a*xb +c.

Перечисленные функции используют трехпараметрическую аппроксимирующую функцию, нелинейную по параметрам. При вычислении оптимальных значений трех параметров регрессионной функции по методу наименьших квадратов возникает необходимость в решении сложной системы из трех нелинейных уравнений. Такая система часто может иметь несколько решений. Поэтому в функциях MathCad, которые проводят регрессию трехпараметрическими зависимостями, введен дополнительный аргумент vg. Данный аргумент – это трехкомпонентный вектор, содержащий приблизительные значения параметров a,b и c, входящих в аппроксимирующую функцию. Неправильный выбор элементов вектора vg может привести к неудовлетворительному результату регрессии. На рис.11 приведен пример проведения экспоненциальной регрессии с помощью функции expfit, регрессия проведена для двух различных значений вектора vg.

Рис.11 Экспоненциальная регрессия

Аппроксимация набора точек функциями, заданными пользователем.

В MathCad существуют средства для проведения регрессии самого общего вида. Это означает, что можно использовать любые функции в качестве аппроксимирующих и находить оптимальные значения любых их параметров, как линейных, так и нелинейных. В том случае, если регрессионная функция является линейной по всем параметрам, т.е. представляет линейную комбинацию жестко заданных функций, провести регрессию можно с помощью встроенной функции linfit(vx,vy,F). Аргумент F – это векторная функция, из элементов которой должна быть построена линейная комбинация, наилучшим образом аппроксимирующая заданную последовательность точек. Результатом работы функции linfit является вектор линейных коэффициентов. Каждый элемент этого вектора – коэффициент при функции, стоящей на соответствующем месте в векторе F. Таким образом, для того чтобы получить регрессионную функцию, достаточно скалярно перемножить эти два вектора. Пример такой аппроксимации представлен на рис.12.

Рис.12 Использование функции linfit

Сглаживание

Система MathCad применяется для обработки различного рода сигналов. При этом очень важной процедурой является очистка сигнала от шумов. Один из вариантов решения данной задачи – это использование алгоритмов сглаживания (smoothig). Существует множество алгоритмов сглаживания данных, причем не всегда можно заранее сказать, какой из них будет наиболее эффективен для той или иной задачи. В MathCad реализовано три алгоритма сглаживания и соответственно существует три функции для их выполнения.

· medsmooth(vy,n) – реализует алгоритм «бегущих» медиан. Параметр n определяет ширину «окна» сглаживания, т.е. количество точек, которые используются при вычислении сглаженного значения в каждой точке. Параметр n должен быть целым нечетным числом. Данный алгоритм лучше всего подходит для сглаживания наборов точек, в которых лишь некоторые точки резко выбиваются из общей гладкой последовательности (см.рис.13).

· ksmooth(vx,vy,b) – алгоритм Гауссового ядра. В данном алгоритме сглаженное значение в каждой точке вычисляется как весовое среднее от всего набора данных с ядром в виде функции Гаусса. Параметр b - это параметр ширины функции ядра. Данный метод наилучшим образом подходит для фильтрации зашумленного сигнала (см.рис.14).

· supsmooth(vx,vy) – в данном алгоритме значение в каждой точке заменятся на значение регрессионной прямой, построенной с использованием некоторого количества близлежащих точек, причем данное количество выбирается по-разному для каждой точки с помощью адаптивного алгоритма.

Рис.13 Сглаживание сигнала с узкими нерегулярностями

Рис.14 Сглаживание зашумленного сигнала

Дискретное преобразование Фурье

Еще одна операция, широко применяемая при обработке сигналов, - вычисление Фурье-спектра. Для проведения данной операции с сигналами, заданными в виде массивов данных, система MathCad содержит численный алгоритм, называемый быстрым преобразованием Фурье (FastFourierTransform FFT). Для реализации данного алгоритма в MathCad существует несколько различных функций.

· В том случае, если набор данных v состоит из 2m элементов, а также все числа в наборе данных действительны, то для проведения Фурье-преобразования следует пользоваться функциями fft(v) или FFT(v). Обе эти функции выполняют Фурье-преобразование, а различие между ними заключается лишь в нормировке результата. Результатом функций fft и FFT будет массив из 2m -1+1 комплексных чисел.

· Если какое-либо из условий, перечисленных выше, не выполняется, т.е. массив содержит комплексные числа или не может быть расширен до размерности 2m, то следует пользоваться функциями cfft(v) и CFFT(v). Результатом функций cfft и CFFT будет массив комплексных чисел той же размерности, что и исходный.

Результатом каждой из перечисленных функций будет массив комплексных чисел. Исходя из этого массива, можно построить амплитудно-частотную (АЧХ) или фазово-частотную (ФЧХ) характеристику сигнала. Для построения АЧХ следует вычислить абсолютное значение каждого элемента в массиве. Пример такого построения приведен на рис.15.

Рис.15 Использование функции fft для построения АЧХ - сигнала

На рис 16 проведена фильтрация сигнала, построенного на рис.15. Конечно, такой метод фильтрации можно использовать, только если уровень шума ощутимо ниже уровня сигнала.

Рис.16 Амплитудная фильтрация сигнала

 


Понравилась статья? Добавь ее в закладку (CTRL+D) и не забудь поделиться с друзьями:  



double arrow
Сейчас читают про: