Лабораторная работа 1. Определение погрешностей

ОГЛАВЛЕНИЕ

1. ПОГРЕШНОСТИ ВЫЧИСЛЕНИЙ.. 4

Лабораторная работа 1. Определение погрешностей при вычислении функций методом разложения их в степенной ряд. 5

2. ЧИСЛЕННЫЕ МЕТОДЫ РЕШЕНИЯ СИСТЕМ ЛИНЕЙНЫХ АЛГЕБРАИЧЕСКИХ УРАВНЕНИЙ.. 12

Лабораторная работа 2. Решение систем линейных алгебраических уравнений методом Гаусса. 25

Лабораторная работа 3. Решение систем линейных алгебраических уравнений методом итераций. 36

3. РЕШЕНИЕ НЕЛИНЕЙНЫХ УРАВНЕНИЙ.. 37

Лабораторная работа 4. Решение не линейных уравнений. 41

4. ИНТЕРПОЛИРОВАНИЕ ФУНКЦИЙ.. 42

Лабораторная работа 5. Методы интерполирования функций. 58

5. ЧИСЛЕННОЕ ИНТЕГРИРОВАНИЕ. 61

Лабораторная работа 6. Методы численного интегрирования. 71

6. ЧИСЛЕННОЕ РЕШЕНИЕ ОБЫКНОВЕННЫХ ДИФФЕРЕНЦИАЛЬНЫХ УРАВНЕНИЙ.. 77

Библиографический список. 86

ПОГРЕШНОСТИ ВЫЧИСЛЕНИЙ

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

Абсолютная погрешность ε– это разность между истинным и приближенным значением величины:

.

Относительная погрешность – это отношение величины абсолютной погрешности к приближению:

.

Для величин, близких по значению к единице, абсолютная и относительная погрешности имеют один порядок. Для очень больших и очень маленьких величин эти погрешности существенно отличаются. Например, если точное значение некоторой величины равно 0,000 06, а приближенное значение 0,000 05, то , а (20 %). Если точное значение величины равно 100 500, а приближенное 100 000, то , а .

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

1. Сложение:

; .

2. Умножение:

; .

3. Деление:

; .

Для уменьшения погрешности в практической деятельности рекомендуется соблюдать следующие правила:

1) по необходимости избегать вычитания двух почти равных чисел;

2) использовать формулы типа

a(b-c) = ab -ac ; (b-c)/a = b/a - c/a.

Причины возникновения погрешностей:

1) неточное представление исходной информации;

2) округление;

3) ограничение бесконечного математического процесса конечным числом операций.

Задача лабораторной работы – определение погрешностей третьего вида.

Лабораторная работа 1. Определение погрешностей

При вычислении функций методом разложения их в степенной ряд

Работа выполняется с использованием палитры программирования, системы автоматизации математических вычислений Mathcad (табл. 1.1).

Задание.Вычислить значение функции по приближенной формуле с заданной точностью:

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

2. Вычислить значение функции с помощью программы на Mathcadе.

Таблица 1.1

Вариант	Функция f (x)	Сумма (приближенное значение) Q (x)	Значение аргумента x	Точ-ность ( )
			0,5
			0,8
			0,5
			0,5
			3,5
			1,5
			0,8

Окончание табл. 1.1

			0,5
			0,5
			0,5
			2,5
			3,5
			3,5
			2,5
			0,6
			0,5
			0,5
			0,5
			0,8
			0,7
			0,8

Требования к оформлению отчета:

1) название и цель работы;

2) математическую формулировку задачи (представление заданной функции степенным рядом, вывод рекуррентного соотношения и алгоритм решения);

3) тексты программы (в системе Mathcad);

4) результаты, полученные в процессе выполнения работы;

5) выводы.

Вспомогательные материалы

Пример вычисления функции f(x) = Sin(x)с заданной точностью

eps = 0,000 001

X = 0,5

.

1. Алгоритм решения задачи.

1.1 . В сумму, предварительно обнуленную, заносится первый член ряда.

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

1.3 . С помощью рекуррентной формулы вычисляется следующий член ряда и добавляется к сумме.

1.4 . Переход на 1.2.

Рис.2.1 а

Решение:Приведем систему к виду (2.10)

.

Запишем последовательность итераций

. (2.19)

Для приведенной матрицы достаточное условие ходимости выполняется по m-норме:

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

Число итераций для достижения заданной точности определяем из неравенства (2.13) , которое запишем так:

, действительно:

.

; , так как

то ; .

Вычислим теперь три последовательных приближения по формулам (2.19) и оценим погрешность каждого результата, используя неравенство (2.18) в виде:

. .

Первое приближение:

.

.

Следовательно, дает значение корня ξ с погрешностью, не превышающей величины .

Далее последовательно находим:

;

.

Третья итерация:

;

.

Заданная точность достигается за пять шагов. Точное решение .

Ниже приведена блок-схема алгоритма решения системы линейных алгебраических уравнений методом итераций.

Рис. 2.2 а

		Рис 2.2 б

Блок-схема алгоритма решения системы линейных уравнений алгебраических уравнений приведена на рис. 2.2 а, рис. 2.2 б.

Вспомогательные материалы

Пример разработки программы для решения системы линейных алгебраических уравнений в системе Mathcad (рис. 2.3 – рис. 2.8).

Объявление переменных.

Определим переменные и зададим матрицу системы линейных уравнений:

Создадим новый лист в MathCAD и для начала напишем оператор ORIGIN:=1.Это нужно для того, чтобы MathCAD начинал нумерацию в матрицах с единицы, а не с нуля.

Теперь определим матрицу А.

Обратите внимание, что в MathCAD строчные и прописные буквы различны, следите за регистром переменных.

Наберите A:=нажмите меню Insert è Matrixзадайте размер-ность матрицы 3 на 3 и нажмите OK. Заполните матрицу так, как показано на рис. 2.3. Тем же способом задайте матрицу-столбец B.

Объединим все в одну расши-ренную матрицу matr с помощью оператора augment.

Присвоим переменной nзначе-ние 3. Она будет определять коли-чество корней системы.

Чтобы посмотреть значение любой переменной, достаточно написать ее имя и поставить знак равенства. Выведем значение перемен-ной matr.

Переменные определены, теперь перейдем к написанию функции проверки матрицы на линейно-зависимые строки.

Begin

d1:=0;

da:=0;

lin:=False;

For i:=0 to n-2 do

Forj:=i+1 to n-1 do

Begin

For k:=0 to n-1 do

If matr[j,k]<>0 then d1:=d1+matr[i,k]/matr[j,k];

d1:=d1/n;

For k:=0 to n-1 do

If matr[j,k]<>0 then If d1=matr[i,k]/matr[j,k] then da:=da+1;

If da=n then lin:=True;

d1:=0;

da:=0;

end;

Iflin then Proverka:=True Else Proverka:=False;

End.

Проверка на линейно-зависимые строки в Mathcadе. Откройте панель, изображенную на рис. 2.4. С помощью меню view è Math Palette, щелкните на третьей строке во втором столбике иконки и перед вами откроется еще одна панель инструментов (см. рис. 2.5), отвечающая за встроенный язык программирования. Так же в процессе решения нам понадобятся еще некоторые панели инструментов (см. рис. 2.6, 2.7). Найдите их на панели Math Palette (см. рис. 2.4).

Оператор Add Line обозначает в MathCAD программный модуль: вертикальную черту, справа от которой последовательно записываются операторы. Оператор ç означает присваивание: F ç 2+B.

Otherwise – это оператор аналогичный Else в Pascal. С остальными операторами вы знакомы.

Напишем имя функции, оператор присваивания и вставим программный модуль, несколько раз нажав кнопку Add Line на панели программирования (см. рис. 2.8).

Определим локальные переменные:

 

Еще одной особенностью языка MathCAD является то, что он сам определяет тип переменных:

d1 –хранит результат деления строк матрицы;

da –хранит количество совпавших элементов деления;

s1 – для хранения сообщения и результата проверки (в Pascal эта перемен-ная lin).

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

Запишем это так:

 

i будет обозначать делимое, а j делитель.

Важно!Чтобы вставить последовательность 1..n-1, нажмите на панели инструментов.

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

Рис 3.

.

 

 

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

Обратите внимание! Что в MathCAD оператор if записывается в отличие от записи в Pascal следующим образом: сначала действие, которое произойдет, если выполнится условие.

Следующим шагом найдем среднее от частного, разделив на число слагаемых n и присвоив это значение d1.

 

 

Теперь проверим, сколько слагаемых равны d1. Переменную da будем

увеличивать на единицу, если условие верно.

Так же, как и в предыдущем случае не стоит допускать деления на ноль, поэтому здесь используется вложенный оператор if.

После выполнения цикла, проверим, равно ли da = n. Если так, то эти строки линейно зависимые.

Пример:

Допустим, матрица выглядит так: .

d1 будет равно = 1/2+2/4+3/6= 1,5

d1 =1,5 / n = 1,5 / 3 = 0,5

d1 = 1/2 = 2 /4 = 3/6 = 0,5 , то da= 3 и равно n, значит, строки линейно-зависимые.

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

Возвратим результат функции. Если значение строковой переменной s1 так и не изменилось, то возвратим матрицу matr, в противном случае выведем сообщение о наличии линейной зависимости.

3. Нахождение максимального элемента в столбце.

Для реализации модифицированного метода Гаусса нам понадобится функция нахождения максимального элемента в столбце и возвращение номера строки этого элемента.

Алгоритм этой функции достаточно прост и не требует подробных объяснений.

4. Перестановка строк в матрице.

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

Определим переменные: i1будет использоваться в цикле, а temp для временного хранения элементов матрицы. Далее цикл и стандартный алгоритм замены переменных.

Возвратим матрицу matr.

5. Прямой ход.

Все дополнительные функции определены, теперь можно перейти непосредственно к поиску решения. Выполним прямой ход в методе Гаусса:

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

Первым делом найдем строку с максимальным элементом через ранее объявленную функцию maxcol и присвоим это значение переменной max. Чтобы застраховать себя от частного случая, проверим, не является ли максимальный элемент равным нулю. Если это так, то система не имеет единственного решения,и выйдем из цикла, используя команду break. В Pascal сразу произведем замену Zamena (i1,MaxCol(i1)) и реализуем алгоритм проверки в теле функции maxcol:

В MathCAD:

Если все в порядке, то поменяем строки местами.

Получим разрешающее уравнение, разделив каждый элемент на первый в строке i1:

Присвоим единицу коэффициенту при x_{i1 .}

Умножим разрешающее уравнение i1 на первый коэффициент уравнения j и вычтем его из уравнения j.

Присвоим ноль коэффициенту при х_{j,i1 .}

После выполнения всех итераций матрица примет треугольный вид, теперь несложно найти корни системы.

Не забудьте возвратить преобразованную матрицу, если система не имеет единственного решения, соответствующее сообщение.

6. Обратный ход

Теперь не составит труда найти корни уравнения, выполнив обратный ход. Напишем еще одну функцию oh:

Выполним обратный цикл, задав диапазон значений от n до2, где n >= 2.

Тем самым пройдем по строкам в обратном порядке.

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

Задача решена, осталось только вывести результаты.

7. Вывод ответа

Создадим матрицу-столбец a из последнего столбца матрицы matr, затем возвратим ее.

Соберем все воедино в функции otvet:

1. Выполним:

1) проверку на линейную зависимость;

2) прямой ход, если значение matr удовлетворяет условию;

3) обратный ход с аналогичной проверкой;

2. Присвоим matr полученный результат;

3. Вернем значение функции (им может быть или матрица-столбец с корнями системы или одно из двух сообщений).

В заключении выведем ответ, набрав otvet =.Ответ может быть трех типов:

Пример программы в MathCAD

Контрольные вопросы:

1. Какой вид имеет матрица коэффициентов после окончания прямого хода в методе Гаусса ?

2. Назовите операции, которые выполняются при обратном ходе?

3. Какие преимущества имеет модифицированный метод Гаусса по сравнению с обычным методом?

РЕШЕНИЕ НЕЛИНЕЙНЫХ УРАВНЕНИЙ

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

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

1. Отделение корней, то есть установление промежутков, в которых содержится один корень уравнения.

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

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

Возможны и другие способы отделения корней, например, графические (рис. 3.1).

После отделения корней для вычисления корня можно применить ниже следующие методы.

Метод половинного деления.

Пусть дано уравнение

f(x) = 0, (3.1)

причем функция f(x) непрерывна на отрезке [a,b] и f(a)f(b) < 0.

Для вычисления корня уравнения (3.1), принадлежащего отрезку [a,b] , найдем середину этого отрезка x₁=(a+b)/2. Если f(x₁) 0 , то для продолжения вычислений выберем ту из частей банного отрезка [a,x₁] или [x_1,b], на концах которой функция f(x) имеет противоположные знаки. Концы нового отрезка обозначим через a₁ и b_1.

Новый суженный промежуток [a₁, b₁] снова делим пополам и продолжаем вычисления по разработанной схеме и т. д. В результате получаем либо точный корень уравнения (3.1) на каком - то этапе , либо последователь-ность вложенных отрезков [a, b] , [a₁, b₁ ] , . . . , [a _n, b_n] , . . таких , что

f(a _n)f(b_n)<0 (n =1,2,...); (3.2)

b_n-a _n=(1/2ⁿ) (b-a).

Число ξ – общий предел последовательности { a _n } и { b _n } – является корнем уравнения f(x) = 0.

Оценку погрешности решения на n-м шаге вычислений можно получить из соотношения (3.3)

0 < ξ - a _n ( 1/ 2 ⁿ) ( b- a ) = b _n - a _n .

Здесь a _n ξ c точностью ε не превышающей ( 1/ 2 ⁿ) ( b- a ).

Рис. 3.1.

Наличие единственного корня уравнения на интервале [a,b]

Вспомогательные материалы

Для разработки программ на Mathcadе можно использовать приемы, описанные в лабораторных работах 1 и 2.

Существует бесчисленное множество линейных функций, имеющих корень на интервале [a,b]. Поясним это на примере.

Пусть a=2,b=7. Пусть корень уравнения равен 5. Тогда функция

F(x) = k1+k2*x, имеющая корень, равный пяти, может иметь следующий вид:

F(x)=10–2*x(один из коэффициентов задается произвольно, другой находится из уравнения F(x) = 0).

Контрольные вопросы:

1. Зачем нужна процедура отделения корней уравнения?

2. Что называется корнем уравнения?

3. Какова точность метода половинного деления?

4. Каким образом исходное уравнение преобразуется к виду, удобному для итераций?

5. Чему равна оценка погрешности k -го приближения?

ИНТЕРПОЛИРОВАНИЕ ФУНКЦИЙ

Постановка задачи интерполирования. На отрезке заданы n+1 точки , которые называются узлами интерполяции, и значение некоторой функции в этих точках

. (4.1)

Требуется построить интерполирующую функцию F(x), принадлежащую известному классу и принимающую в узлах интерполяции те же значения, что и :

(4.2)

Рис. 4.1. Интерполирование функции y = f(x)

7. Сплайн – интерполяция (spline – рейка, планка).

Механические сплайны – гибкие деревянные рейки, закрепленные на концах. В узлах (точках) интерполяции подвешивают грузила. Сплайн принимает форму, минимизирующую его потенциальную энергию. Если сплайн представить функцией S(x) , то S и непрерывны на [ ].

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

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

Пусть отрезок [a, b] разбит на n частей точками .

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

Сплайн 1-й степени – кусочно-линейная функция (непрерывная). Производная терпит разрыв в точках излома.

Задача интерполяции функции на отрезке [a, b] кубическим сплайном (сплайном 3-й степени) состоит в нахождении функции S(x), равной многочлену третьей степени на каждом отрезке ,