Задачи изучения дисциплины 5 4 Перечень дисциплин, усвоение которых необходимо для изучения данной дисциплины: 6

СОДЕРЖАНИЕ

1 Введение

1.1 Порядок изучения дисциплины

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

Отчетность: 5 семестр – экзамен;

6 семестр – дифференциальный зачет.

1.2 Цель преподавания дисциплины

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

1.3 Задачи изучения дисциплины

В результате изучения студенты должны:

знать: принципы проведения вычислительного эксперимента, характеристики вычислительных задач, источники погрешностей вычислений, основные методы и алгоритмы решения стандартных вычислительных задач;

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

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

1.4 Перечень дисциплин, усвоение которых необходимо для изучения данной дисциплины:

математика;
дискретная математика;
теория вероятности и математическая статистика;
информатика;
алгоритмические языки и программирование.

2 Программа лекционного курса

2.1 Введение в численные методы

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

Характеристики вычислительных задач. Устойчивые и неустойчивые, корректные и некорректные задачи. Примеры некорректных задач. Требования к вычислительным методам. Устойчивость, корректность, сходимость. Пример неустойчивого алгоритма.

2.2 Погрешности округления в ЭВМ

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

Погрешности округления при выполнении арифметических операций в ЭВМ. Погрешности суммы двух и нескольких чисел. Зависимость погрешности от порядка суммирования. Погрешности произведения двух и нескольких чисел. Алгоритм вычисления произведения чисел. Правила выполнения арифметических операций в ЭВМ. Статистические оценки погрешностей. Примеры организации вычислений.

2.3 Вычисление значений функций

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

2.4 Методы решения систем линейных алгебраических уравнений

Классификация и характеристики методов решения систем линейных алгебраических уравнений (СЛАУ). Прямые и итерационные методы. Методы Крамера, обратной матрицы. Метод Гаусса (схема единственного деления). Метод Гаусса с выбором главного элемента. Вычисление определителя и обратной матрицы методом Гаусса. Метод прогонки.

Погрешности решения СЛАУ. Нормы векторов и матриц. Оценка погрешностей. Число обусловленности. Оценка числа обусловленности.

Итерационные методы решения СЛАУ. Метод простых итераций (метод Якоби). Условия сходимости. Оценка числа итераций. Метод Зейделя.

2.5 Решение нелинейных уравнений с одним неизвестным

Прямые и итерационные методы решения. Число корней нелинейных уравнений. Отделение корней.

Методы уточнения корней. Метод дихотомии. Метод хорд. Метод Ньютона (касательных). Условия и скорость сходимости метода Ньютона. Модифицированный метод Ньютона и метод секущих.

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

2.6 Решение систем нелинейных уравнений

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

2.7 Моделирование случайных величин

Основные характеристики случайных величин. Получение случайных величин на ЭВМ. Генераторы случайных чисел. Метод Монте-Карло. Применение метода Монте-Карло для вычисления определенных интегралов.

2.8 Аппроксимация функций

Понятие приближения функций. Применение аппроксимации функций в САПР. Критерии близости функций. Оптимальная аппроксимация. Классификация задач аппроксимации.

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

Свойства интерполяционных моделей. Погрешность интерполяции. Многочлены Чебышева. Оптимальный выбор узлов интерполяции. Сходимость интерполяции. Теорема Фабера. Локальная интерполяция. Применение глобальной и локальной интерполяции.

Интерполяция с помощью сплайнов. Понятие сплайна. Построение кубического сплайна. Основные соотношения. Применение сплайнов.

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

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

Аппроксимация методом разложения в степенной ряд. Многочлен Тейлора. Погрешность приближения многочленом Тейлора.

2.9 Численное дифференцирование и интегрирование

Прямое вычисление производных. Левая, правая и центральная разностные производные. Ошибки численного дифференцирования. Применение интерполяции.

Численное интегрирование. Формулы прямоугольников, трапеций, Симпсона. Ошибки численного интегрирования. Выбор шага интегрирования.

Учебно-методические материалы по дисциплине

Основная литература

Бабак Л.И. Вычислительные методы. Курс лекций (части 1 и 3).  Томск: ТУСУР, 2002.
Черкашин М.В. Вычислительные методы. Курс лекций (часть 2).  Томск: ТУСУР, 2002.
Демидович Б.П., Марон И.А. Основы вычислительной математики. – М.: Наука, 1966. – 664 с.
Вержбицкий В.М. Численные методы (линейная алгебра и нелинейные уравнения): Учебное пособие для ВУЗов. – М.: Высшая школа, 2000. – 266 с.
Вержбицкий В.М. Численные методы (математический анализ и обыкновенные дифференциальные уравнения): Учебное пособие для ВУЗов. – М.: Высшая школа, 2001. – 382 с.
Турчак Л.И. Основы численных методов. – М., Наука, 1987. – 320с.
Тихонов А.Н., Костомаров Д.П. Вводные лекции по прикладной математике. – М., Наука, 1984.
Самарский А.А., Гулин А.В. Численные методы. – М., Наука, 1989. – 432с.
Боглаев Ю.П. Вычислительная математика и программирование. – М., Высшая школа, 1990. – 544с.

Дополнительная литература

Волков Е.А. Численные методы. – М., Наука, 1987. – 248с.
Мак-Кракен Д., Дорн У. Численные методы и программирование на Фортране. – М., Мир, 1977.
Каханер Д., Моулер К., Нэш С. Численные методы и программное обеспечение. – М., Мир, 1988. –575с.
Копченова Н.В., Марон И.А. Вычислительная математика в примерах и задачах. – М.: Наука, 1972. – 368 с.

Методические указания и прочие учебно-методические материалы

Потемкин В.Г. Система MATLAB: Справочное пособие – М., Диалог –МИФИ, 1997. –352с.
Потемкин В.Г. MATLAB 5 для студентов. – М., Диалог –МИФИ, 1998. –314с.
Черкашин М.В. Система для математических и инженерных расчетов MATLAB: Учебное пособие.  Томск: ТУСУР, 2002.
Ракитин В.И., Первушин В.Е. Практическое руководство по методам вычислений с приложением программ для персональных компьютеров: Учебное пособие. – М.: Высшая школа, 1998, - 383 с.

3 Контрольные работы

3.1 Общие указания

Студентам необходимо выполнить две контрольных работы по курсу.

Вариант контрольной работы выбирается согласно формуле

N_в=(N*k) div 100,

где N_в – искомый номер варианта;

N – общее число возможных вариантов (N=10);

k – две последних цифры пароля (число от 00 до 99);

div – операция целочисленного деления.

Если получается N_в =0, то берется 10 вариант.

Студент выполняет вариант с одним и тем же номером во всех контрольных работах.

Ответы даются письменно или в электронном варианте с использованием тестового редактора MS Word.

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

3.2 Контрольная работа № 1

Контрольная работа № 1 выполняется после изучения разделов 2.1–2.4 и содержит 6 вопросов.

3.2.1 Пример решения типового варианта

Задача 1

Округлить по правилу симметричного округления число а^*=3,1415926 последовательно до тысячных, сотых и десятых. Найти предельные абсолютную и относительную погрешности каждого результата.

Решение:

При симметричном округлении за предельную абсолютную погрешность приближенного числа a принимается половина единицы последнего сохраняемого разряда числа. Таким образом, при симметричном округлении числа а^*=3,1415926 получим следующий результат:

а₁=3,1 – предельная абсолютная погрешность равна (а₁)=0,05;

а₂=3,14 – предельная абсолютная погрешность равна (а₂)=0,005;

а₃=3,142– предельная абсолютная погрешность равна (а₁)=0,0005;

Предельную относительную погрешность числа можно определить по формуле:

Для нашей задачи получим:

;

Задача 2

С каким числом верных знаков в узком смысле следует взять число arcsin(1), чтобы относительная погрешность была не более 0,1%?

Решение:

Вначале вычислим:

а^*=arcsin(1)=/21,57079 …

По значению предельной относительной погрешности

(а^*)=0,1 %=0,001

найдем предельную абсолютную погрешность числа а*:

(a^*) = a^*(а^*) = 1,57079…0,001  0,00158.

(по правилу округления предельной погрешности, последнюю сохраняемую цифру числа увеличиваем на единицу).

Представим число а^* в виде

а^* = 1,57079… = 110⁰ + 510^–1 + 710^–2 + 010 ^–3+ 710^–4+ 910^–5+ … и составим таблицу

Цифра

…

Единица разряда

0,1

0,01

0,001

0,0001

0,00001

…

Половина единицы разряда

0,5

0,05

0,005

0,0005

0,00005

0,000005

…

Сравнивая для каждой цифры половину единицы разряда с абсолютной погрешностью числа (a^*) = 0,00158, получим, что а^* следует записать в виде а=1,57.

Проверка: < 0,1 %.

Ответ: arcsin(1)=1,57 – число верных знаков (в узком смысле) равно 3.

Задача 3

Найти допустимые абсолютные погрешности аргументов и , которые позволяют вычислить функцию f с точностью (f) = 0,2 %:

Решение:

Значение функции при x₁=2,4563 и x₂=0,8473 равно f(x₁, x₂)  2,59833. Предельную абсолютную погрешность (f) можно оценить следующим образом:

(f) = | (f)  f(x₁, x₂) | = 0,0022,59833  0,0051967 = 0,006

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

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

;

Таким образом, величины аргументов x₁=2,4563 и x₂=0,8473 следует брать с предельными абсолютными погрешностями (x₁) = 0,005 и (x₂) = 0,0085.

Ответ: (x₁) = 0,005; (x₂) = 0,0085.

Задача 4

Для матрицы вычислить число обусловленности cond(A).

Решение: число обусловленности можно найти по формуле

cond(A)=||A||||A^–1||,

где ||A|| – норма матрицы А; ||A^–1|| – норма обратной матрицы A^–1.

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

Обратную к матрицу A^–1 можно получить по формуле

где  – определитель матрицы А;

А_ij – алгебраическое дополнение элемента a_ij.

А_ij = (–1)⁽ⁱ⁺^j⁾M_ij,

где M_ij – минор элемента a_ij, т.е. определитель (n–1)-го порядка, получаемый из определителя  вычеркиванием i-й строки и j-го столбца.

Для нашей задачи =14 – 23= –2.

А₁₁=(–1)¹⁺¹4 = 4; А₁₂=(–1)¹⁺²3= –3; А₂₁=(–1)²⁺¹2= –2; А₂₂=(–1)²⁺²1=1.

Отсюда получим обратную матрицу .

Используя разные нормы матриц (-норму и 1-норму), найдем число обусловленности матрицы А.

а) -норма: – максимальная сумма модулей элементов строк. Для нашей задачи:

||A||_= max{ | 1 |+| 2 |, | 3 |+| 4 | } = max{ 3, 7 } = 7;

||A^–1||_= max{ | –2 |+| 1 |, | 1,5 |+| –0,5 | } = max{ 3, 2 }=3.

б) 1-норма: - максимальная сумма элементов столбцов. Для нашей задачи:

||A||₁= max{ | 1 |+| 3 |, | 2 |+| 4 | } = max{ 4, 6 } = 6;

||A^–1||₁= max{ | –2 |+| 1,5 |, | 1 |+| –0,5 | } = max{ 3,5, 1,5 }=3,5.

По полученным нормам рассчитаем число обусловленности матрицы

cond(A) = ||A||_||A^–1||_= 73 = 21;

cond(A) = ||A||₁||A^–1||₁= 63,5 = 21.

Ответ: число обусловленности равно cond(A)=21.

Задача 5

Решить систему методом Гаусса с выбором главного элемента в строке:

Решение:

1 Прямой ход

Первый шаг:

а) выполним выбор главного элемента в 1-м столбце, получим новую систему

перестановкой 1-й и 2-й строк в исходной СЛАУ.

б) для исключения переменной x₁ из второго уравнения умножим первое уравнение на коэффициент и сложим его со вторым уравнением:

2 x₁–x₁+ x₂–x₂+ x₃–x₃= 0 – .

Теперь упростим:

0 x₁+x₂+x₃=0+2;

0 x₁+0,6 x₂+1,4 x₃=2.

в) для исключения x₁ из третьего уравнения умножим первое уравнение на коэффициент и сложим его с 3-м уравнением:

–2 x₁+x₁+ x₂+x₂–2 x₃+x₃=1+;

0 x₁+x₂+x₃=1–2;

0 x₁+1,4 x₂– 2,4x₃= –1.

В результате этих преобразований после первого шага получим:

Видно, что коэффициенты при x₁ во втором и третьем уравнении равны нулю.

Второй шаг:

а) выполним выбор главного элемента во втором столбце (первую строку не трогаем), получим новую систему

перестановкой 2-й и 3-й строк системы.

в) теперь с помощью второго уравнения можно исключить переменную x₂ из третьего уравнения системы. Умножим второе уравнение на коэффициент и сложим его с 3-м уравнением:

0,6 x₂+x₂+1,4 x₃+x₃=2+;

0 x₂+2,4286 x₃=2,4286.

В результате получаем следующую систему с треугольной матрицей коэффициентов:

2 Обратный ход

Вторая часть метода Гаусса, называемая обратным ходом, заключается в последовательной подстановке найденных переменных x₃ и x₂ во второе и первое уравнения СЛАУ и их решения.

Первый шаг:

а) находим значение переменной x₃ из последнего уравнения

x₃=2,4286/2,4286=1;

б) подставляем полученное значение x₃ во второе уравнение системы:

0 x₁+1,4 x₂–2,41=–1,

решая которое, получим значение x₂=(–1+2,4)/1,4 = 1.

Второй шаг:

Подставляем найденные значения x₃ и x₂ в первое уравнение системы:

5 x₁+11–11= –5;

x₁= (–5)/5= –1;

x₁= –1.

В результате получаем решение заданной системы линейных уравнений x=[–1, 1, 1]^t.

Проверка:

Ответ: вектор неизвестных равен x=[–1, 1, 1]^t.

Задача 6

Решить систему итерационным методом Гаусса-Зейделя, сделать 3 шага итерационного процесса. (При необходимости для выполнения условий сходимости сначала преобразовать систему).

Решение:

Приведем исходную систему уравнений к виду, необходимому для выполнения итерационного процесса. Для этого выразим неизвестные x₁, x₂ и x₃ из первого, второго и третьего уравнений СЛАУ. Получим:

В матричном виде приведенная система имеет вид x = x + ,

где ; .

Проверим условие сходимости итераций. Исходная система не отвечает условиям диагонального преобладания

; . (*)

Следовательно, достаточное условие сходимости не выполняется.

Исследуем также сходимость метода, используя нормы матрицы для приведенной системы:

а) -норма: =max{(0++); (2+0+4); (0+2+0)}=

= max{; 6; 2} = 6.

б) 1-норма: =max{(0+2+0); (+0+3); (+4+0)}=

=max{2; ; } = 4,33.

Так как обе нормы больше 1, то достаточное условие сходимости итерационного процесса || || < 1 не выполняется.

Преобразуем исходную систему таким образом, чтобы получить матрицу А с преобладающими диагональными элементами. Для этого переставим местами 2 и 3 строки:

Проверим условия сходимости (*):

для строк матрицы А: | 3 | > | 1 | +| 1 |; | –2 | > | 0 | +| 1 |; | 4 | > | 2 | +| –1 |;

для столбцов матрицы А: | 3 | > | 0 | +| 2 |; | –2 | = | 1 | +| –1 |; | 4 | > | 1 | +| 1 |.

Так как одна из групп неравенств (*) (для строк) выполняется со строгим знаком неравенства, система удовлетворяет достаточному условию сходимости.

Приведем полученную СЛАУ к итерационному виду. Получим следующую систему

В матричном виде приведенная система имеет вид x = x + ,

где ; .

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

а) -норма: =max{(0++); (0+0+); (++0)}=

= max{; ; } = 0,75.

б) 1-норма: =max{(0+0+); (+0+); (++0)}=

=max{; ; } = 0,83.

Как и следовало ожидать, в данном случае условие сходимости выполняется и по -норме и по 1-норме.

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

1-е приближение:

=max{2–1,5875; 1,5–1,375; –0,25+0,7}=

= 3max{0,4125; 0,125; 0,45}=1,35.

2-е приближение:

=3max{1,773–1,5875; 1,15–1,375; 0,8511+0,7}=

= 3max{0,1898; 0,2250; 0,1511}=0,675.

3-е приближение:

=3max{1,9014–1,773; 1,0744–1,15; –0,9321+0,8511}=

= 3max{0,1241; 0,0756; 0,081}=0,3723.

Итак, полученное решение х = [1,9014, 1,0744, –0,9321]^Т. Ошибка решения ₃0,3723.

Проверка: подставим найденный вектор х в исходную СЛАУ:

Полученный вектор свободных членов близок к заданному вектору .

Замечание: Точное решение данной системы уравнений равно х^*= [2, 1, –1]^Т.

Ответ: решение исходной системы уравнений равно

х = [1,9014, 1,0744, –0,9321]^Т.

3.2.2 Варианты контрольной работы № 1

Вариант №1

Представить число 8912,342 в виде десятичной дроби (разложить по степеням числа 10).

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

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

Для матрицы А вычислить обратную матрицу (использовать алгебраические дополнения элементов матрицы). Убедиться, что .

Решить СЛАУ методом простой итерации, сделать 3 шага итерационного процесса (при необходимости для выполнения условий сходимости сначала преобразовать систему).

Решить СЛАУ методом Гаусса с выбором главного элемента в столбце:

Вариант №2

Определить абсолютную погрешность (а) приближенного числа а по его относительной погрешности (а):

Округлить сомнительные цифры числа а, оставив в нем верные знаки в узком смысле:

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

Для матрицы А вычислить нормы , , :

Решить СЛАУ методом обратной матрицы. Найти число обусловленности системы уравнений. Оценить возможное отклонение от полученного решения, если допустить относительную ошибку в левой части (использовать любую из норм).

Решить СЛАУ методом Зейделя, сделать 3 шага итерационного процесса (при необходимости для выполнения условий сходимости сначала преобразовать систему).

Вариант №3

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

Округлить сомнительные цифры числа а, оставив в нем верные знаки в узком смысле, если

и предельная относительная погрешность

Определить относительную погрешность вычисления объема цилиндра V, если погрешность измерения радиуса основания R равна 1% и погрешность измерения высоты H равна 1%.

Найти произведение , если:

а) ;

б) .

Решить систему уравнений: а) обычным методом Гаусса (схема единственного деления); б) методом Гаусса с выбором главного элемента. Все вычисления производить с точностью до 4-х значащих цифр. Вычислить число обусловленности системы уравнений (использовать любую из норм), сделать выводы.

Показать, что для СЛАУ процесс решения по методу Зейделя сходится. Определить необходимое число итераций для нахождения корней системы с точностью до .

Вариант №4

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

Определить, какое из измерений выполнено точнее – 80 км с ошибкой 20 м или 8 см с ошибкой 2 мм (сравнить относительные погрешности измерений).

При измерении радиуса R круга с точностью до 0,5 см получилось число 12 см. Найти предельные абсолютную и относительную погрешности при вычислении площади круга.

Найти произведение , если:

Вычислить определитель методом Гаусса:

Вариант №5

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

С каким числом верных знаков в узком смысле следует взять числа arctg(6) и ln(30), чтобы относительная погрешность была не более 0,1%?

Каждое ребро куба, измеренное с точностью до 0,02 см, оказалось равным 8 см. Найти предельные абсолютную и относительную погрешности при вычислении объема куба.

Решить СЛАУ по методу обратной матрицы:

Вычислить обратную матрицу методом Гаусса:

Вариант №6

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

С каким числом верных знаков в узком смысле следует взять числа arctg(9) и ln(45), чтобы относительная погрешность была не более 0,1% ?

Найти предельные абсолютную и относительную ошибки при вычислении напряжения U по закону Ома (), если ток I равен и сопротивление R равно 10 Ом с точностью 10%.

Решить систему уравнений методом прогонки:

Вариант №7

Определить абсолютную погрешность (а) приближенного числа а по его относительной погрешности (а):

Определить, какое из измерений выполнено точнее – 100 км с ошибкой 30 м или 9 см с ошибкой 4 мм (сравнить относительные погрешности измерений).

Реактивное сопротивление емкости (в Омах) задается формулой

где f – частота в герцах,

С – емкость в фарадах.

Указать границы возможных значений для Гц и

4. Для матрицы А вычислить нормы , , :

Вариант №8

С какой точностью следует определить радиус основания R и высоту H цилиндрической банки, чтобы ее вместимость можно было вычислить с точностью 1% ?

Найти произведение , если:

;

Вариант №9

Представить число 6834,148 в виде десятичной дроби (разложить по степеням числа 10).

Найти допустимые абсолютные погрешности аргументов и , которые позволяют вычислить значения функции f с 4 верными знаками:

Найти произведение , если:

Решить СЛАУ методом Гаусса с выбором главного элемента в столбце:

Вариант №10

Округлить сомнительные цифры числа а, оставив в нем верные знаки в узком смысле:

С каким числом верных знаков в узком смысле следует взять значения аргумента х, чтобы получить значение функции f с точностью 0,1%:

Решить СЛАУ по методу обратной матрицы:

Вычислить определитель методом Гаусса:

3.3 Контрольная работа № 2

Контрольная работа № 2 выполняется после изучения разделов 2.5–2.9 и содержит 5 вопросов.

3.3.1 Пример решения типового варианта

Задача 1

Вычислить методом секущих корень уравнения в интервале [–3, 0]. Сделать три итерации.

f(x) = x²+2x–1=0.

Решение:

Выполним отделение корней на заданном интервале. Для этого рассчитаем значение заданной функции в интервале от –3 до 0 с шагом 0,5. Результат показан в таблице ниже.

–3,0

–2,5

–2,0

–1,5

–1,0

–0,5

f(x)=

x²+2x–1

2,0

0,25

–1,0

–1,75

–2,0

–1,75

–1,0

Из таблицы видно, что искомый корень лежит в интервале [–2,5, –2,0]. Дальнейшее уточнение корня выполним методом секущих. Итерационная формула метода секущих имеет следующий вид

В качестве начального приближения возьмем левую границу и середину найденного интервала, т.е. х⁽⁰⁾= –2,5, х⁽¹⁾= (–2,5+(–2,0))/2 = –2,25.

1-я итерация:

2-я итерация:

3-я итерация:

Проверка: f(x⁽⁴⁾) = (–2,4142)² + 2(–2,4142) – 1  –3,8410^–5.

Ответ: найденный корень равен x–2,4142.

Задача 2

Функция y = f(x) задана таблично. Найти коэффициенты интерполирующего полинома непосредственно из условий интерполяции. Рассчитать значение первой производной в точке 2,5.

–5

–8

–12

Решение:

Интерполирующий полином имеет следующий вид: P_n(x) = C₀+C₁x + C₂x². Для нахождения коэффициентов C₀, C₁, C₂ составим следующую систему уравнений, используя условие интерполяции P_n(x_i)= y_i, :

C₀+C₁x₀ + C₂x₀²= y₀

C₀+C₁x₁ + C₂x₁²= y₁

C₀+C₁x₂ + C₂x₂²= y₂

или в матричном виде

 AC=B,

где – вектор искомых коэффициентов полинома.

Подставляем значения x_i и y_i из таблицы, получим

, т.е. .

Решение данной системы можно выполнить методом обратной матрицы. Вектор искомых коэффициентов найдем по формуле C=A^–1B:

Таким образом, интерполирующий полином будет иметь следующий вид

P_n(x) = –3 – 1,5x – 0,5x².

Проверка: P_n(x₀) = P_n(1) = –3 – 1,51 – 0,51²= –5;

P_n(x₁) = P_n(2) = –3 – 1,52 – 0,52²= –8;

P_n(x₂) = P_n(3) = –3 – 1,53 – 0,53²= –12.

Видно, что полученный полином проходит через заданные точки (x_i, y_i).

Получим значение первой производной: P'_n(x) = – 1,5 – x.

Значение производной в точке x=2,5 равно P'_n(2,5) = – 1,5 – 2,5 = –4.

Ответ: значение производной равно (2,5)  P'_n(2,5) = –4.

Задача 3

Для функции y = f(x), заданной таблицей, с помощью метода наименьших квадратов (МНК) найти коэффициенты аппроксимирующей функции в виде P_n(x) = C₀+C₁x

Решение:

В методе наименьших квадратов минимизируется сумма квадратов разностей между значениями функции = C₀+C₁x и измеренными значениями y_i в точках x_i, :

В нашем случае число точек равно 3, т.е. n=2. Тогда

Дифференцируя функцию ошибки по переменным С₀, С₁, получим систему уравнений следующего вида:

В матричном виде эту систему можно записать как

или AC=B.

Для решения данной системы используем правило Крамера:

, ,

где  – определитель матрицы коэффициентов A=;

₀ – определитель матрицы , полученной из исходной матрицы коэффициентов A заменой первого столбца вектором свободных членов B.

₁ – определитель матрицы , полученной из исходной матрицы коэффициентов A заменой второго столбца вектором свободных членов B.

Следовательно, решением данной системы уравнений будет

;

Замечание: При преобразовании формул следует правильно выполнять действия с величинами, входящими под знак суммы. В частности, необходимо иметь в виду, что , и т.д. В случае сомнений рекомендуется проверить преобразования, раскрывая знаки сумм при небольшом числе слагаемых. Например, раскрывая первую из приведенных формул при числе слагаемых n=2, легко убедиться, что .

Подставляя в полученные выражения значения x_i и y_i из таблицы, получим:

;

Таким образом, аппроксимирующая функция будет иметь следующий вид: (x) = –0,3333 + 2x.

Проверка: для проверки построим график функции (x) = –0,3333 + 2x и нанесем на него точки (x_i, y_i).

Ответ: аппроксимирующая функция имеет вид (x) = –0,3333 + 2x.

Задача 4

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

0,15

0,20

0,25

0,0034

0,0080

0,0156

Решение:

Аппроксимация первой производной по формуле левых разностных производных (ЛРП) имеет следующий вид:

Отсюда получим

Аппроксимация первой производной по формуле правых разностных производных (ПРП) имеет следующий вид:

Отсюда получим

Аппроксимация первой производной по формуле центральных разностных производных (ЦРП) имеет следующий вид:

Отсюда получим

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

а производные первого порядка y'_i₊₁и y'_i – c помощью левых разностей:

и окончательно получим

Для нашей задачи вторая производная равна

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

Ответ: левая разностная производная: y'(0,2) = 0,092;

правая разностная производная: y'(0,2) = 0,152;

центральная разностная производная: y'(0,2) = 0,122;

вторая производная: y''(0,2) = 1,2.

Задача 5

Вычислить интеграл по формулам прямоугольников, трапеций, Симпсона (взять 5 значений подынтегральной функции, т.е. число подинтервалов n=4). Рассчитать оценки погрешности вычислений интеграла, сравнить с действительными погрешностями.

Решение:

Выполним разбиение интервала интегрирования на 4 части и рассчитаем значение функции в узлах. Результаты представим в виде таблицы.

x_i

0,25

0,5

0,75

1,0

у_i =f(x_i)

2,0

2,2188

3,0

4,5312

7,0

а) Метод средних прямоугольников. Формула вычисления интеграла функции f(x) по методу средних прямоугольников имеет следующий вид

где .

Вычислим значения функции в середине выбранных интервалов и результаты сведем в таблицу

x_ср_i

0,125

0,375

0,625

0,875

у_ср_i =f(x_ср_i)

2,0508

2,5273

3,6602

5,6367

Значение искомого интеграла по формуле прямоугольников будет равно:

I_пр=0,25(2,0508+2,5273+3,6602+5,6367)=3,4688.

Погрешность метода прямоугольников можно оценить по формуле

где .

Для нашего случая имеем: h=0,25; (x)=6x²+ 6x; (x)=12x + 6; . Следовательно

б) Метод трапеций. Формула вычисления интеграла функции f(x) по методу трапеций имеет следующий вид

= 0,25((2,0+7,0)/2+2,2188+3,0+4,5312)=3,5625.

Погрешность метода трапеций можно оценить по формуле

0,0938.

в) Метод Симпсона. Формула вычисления интеграла функции f(x) по методу Симпсона имеет следующий вид (n – должно быть четным числом)

где ₁= f(x₁) + f(x₃) +…+ f(x_n_–1); ₂= f(x₂) + f(x₄) +…+ f(x_n_–2)

Для нашей задачи формула Симпсона примет следующий вид

 =

= 0,25/3(2,0+7,0+4(2,2188+4,5312)+23,0)=3,5.

Погрешность метода Симпсона можно оценить по формуле

где . Четвертая производная заданной функции f(x)=2x³+ 3x²+ 2 равна нулю, т.е. M₄=0.

Проверка: Точное значение интеграла равно

Следовательно, истинные значения ошибок будут равны:

а) метод прямоугольников =| I – I_пр | = | 3,5 – 3,4688 | = 0,0312;

б) метод трапеций = | I – I_тр | = | 3,5 – 3,5625 | = 0,0625;

в) метод Симпсона = | I – I_c | = | 3,5 – 3,5 | = 0.

Таким образом, оценки ошибок близки к действительным ошибкам вычисления интеграла.

3.3.2 Варианты контрольной работы № 2

Вариант №1

Вычислить методом дихотомии (деления отрезка пополам) корень уравнения в интервале . Сделать три итерации.

Функция интерполируется на интервале [2; 4] полиномом

1) Найти оптимальные (чебышевские) узлы интерполяции.

2) Оценить ошибку интерполяции при таких узлах.

3) Выполнить интерполяцию и сравнить оценку с действительной ошибкой в точке .

Экспериментально получена зависимость дальности приема радиосигнала D от интенсивности осадков I:

I, мм/час

D, км

100

С использованием МНК найти приближенную зависимость D(I) в виде квадратного трехчлена.

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

1,00

1,10

1,20

2,7183

3,0042

3,3201

Вычислить интеграл по формулам прямоугольников, трапеций, Симпсона (взять 5 значений подынтегральной функции, т.е. ). Рассчитать оценки погрешности вычислений интеграла, сравнить с действительными погрешностями.

Вариант №2

Вычислить методом хорд корень уравнения в интервале . Сделать три итерации.

Методом интерполяции найти первую производную функции в точке . Функция задана числовыми значениями. Указание: вначале найти интерполирующий полином.

Экспериментально получена зависимость дальности приема радиосигнала D от интенсивности осадков I:

I, мм/час

D, км

100

С использованием МНК найти приближенную зависимость D(I) в виде квадратного трехчлена.

Дан полином

Вычислить при (шаг ):

1) по формулам левой, правой и центральной разностных производных.

2) с помощью интерполяции по 5 узлам.

3) по формуле центральной разностной производной.

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

Определить, с какой точностью можно вычислить , используя 9 значений подынтегральной функции ():

1) по формуле прямоугольников;

2) по формуле трапеций;

3) по формуле Симпсона.

Вариант №3

Вычислить методом Ньютона корень уравнения в интервале . Сделать три итерации.

Функция интерполируется полиномом в точках

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

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

Указание: вначале выполните преобразование аппроксимирующей функции к линейному виду.

1,00

1,10

1,20

0,3679

0,3329

0,3012

Сколько узлов следует взять для вычисления интеграла с точностью при использовании формулы Симпсона?

Вариант №4

Вычислить модифицированным методом Ньютона корень уравнения в интервале . Сделать три итерации.

Функция интерполируется полиномом в точках

1) Найти коэффициенты полинома.

2) Оценить ошибку интерполяции в точке

3) Сравнить оценку с действительной ошибкой.

Экспериментально получена зависимость дальности приема D сигнала от мощности передатчика Рп. С использованием МНК найти приближенную зависимость D(Рп) в виде квадратного трехчлена.

Рп, кВт

D, км

120

190

230

260

0,00

0,12

0,24

0,0000

0,34045

0,76507

Определить, с какой точностью можно вычислить , используя 9 значений подынтегральной функции ():

1) по формуле прямоугольников;

2) по формуле трапеций;

3) по формуле Симпсона.

Вариант №5

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

Дана таблично заданная функция . Пользуясь формулой Лагранжа, найти значение функции в точке

2,5

3,5

Рп, кВт

D, км

100

150

170

200

Дан полином

Вычислить при (шаг ):

1) по формулам левой, правой и центральной разностных производных.

2) с помощью интерполяции по 5 узлам.

3) по формуле центральной разностной производной.

Определить, с какой точностью можно вычислить , используя 9 значений подынтегральной функции ():

1) по формуле прямоугольников;

2) по формуле трапеций;

3) по формуле Симпсона

Вариант №6

Уточнить корень уравнения методом Ньютона и модифицированным методом Ньютона. Сделать три итерации, начальная точка – ближайшее к корню целое число.

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

3,60

3,65

3,70

36,598

38,475

40,447

Указание: вначале выполните преобразование аппроксимирующей функции к линейному виду.

1,00

1,10

1,20

0,8415

0,8912

0,9320

Сколько узлов следует взять для вычисления интеграла с точностью при использовании формулы Симпсона?

Вариант №7

Вычислить методом дихотомии (деления отрезка пополам) корень уравнения в интервале . Сделать три итерации.

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

Рп, кВт

D, км

120

160

200

1,00

1,10

1,20

0,5403

0,4536

0,3624

Вариант №8

Вычислить методом хорд корень уравнения в интервале . Сделать три итерации.

Функция интерполируется на интервале [3; 5] полиномом

Найти оптимальные (чебышевские) узлы интерполяции.

Оценить ошибку интерполяции при таких узлах.

Выполнить интерполяцию и сравнить оценку с действительной ошибкой в точке .

Экспериментально получена зависимость дальности приема радиосигнала D от интенсивности осадков I:

I, мм/час

D, км

150

С использованием МНК найти приближенную зависимость D(I) в виде квадратного трехчлена.

Дан полином

Вычислить при (шаг ):

1) по формулам левой, правой и центральной разностных производных.

2) с помощью интерполяции по 5 узлам.

3) по формуле центральной разностной производной.

Определить, с какой точностью можно вычислить , используя 9 значений подынтегральной функции ():

1) по формуле прямоугольников;

2) по формуле трапеций;

3) по формуле Симпсона.

Вариант №9

Вычислить методом Ньютона корень уравнения в интервале . Сделать три итерации.

Указание: вначале выполните преобразование аппроксимирующей функции к линейному виду.

1,00

1,10

1,20

1,1752

1,3356

1,5095

Сколько узлов следует взять для вычисления интеграла с точностью при использовании формулы Симпсона?

Вариант №10

Вычислить модифицированным методом Ньютона корень уравнения в интервале . Сделать три итерации.

Функция интерполируется полиномом в точках

Экспериментально получена зависимость дальности приема радиосигнала D от интенсивности осадков I:

I, мм/час

D, км

С использованием МНК найти приближенную зависимость D(I) в виде квадратного трехчлена.

Дан полином

Вычислить при (шаг ):

1) по формулам левой, правой и центральной разностных производных.

2) с помощью интерполяции по 5 узлам.

3) по формуле центральной разностной производной.

Сколько узлов следует взять для вычисления интеграла с точностью при использовании формулы Симпсона?

4 Лабораторные работы

4.1 Общие указания

Целью лабораторных работ является закрепление на практике теоретических знаний, полученных в ходе изучения дисциплины «Вычислительные методы».

Предполагается, что студенты уже теоретически освоили соответствующие методы вычислительной математики и имеют опыт программирования на алгоритмических языках. Весь курс работ поставлен в среде системы для математических и инженерных расчетов MATLAB версии не ниже 4.2. Ознакомиться с системой MATLAB можно в [14-16], а также в сети INTERNET на сайтах www.exponenta.ru\matlab и www.matlab.ru.

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

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

В ходе выполнения лабораторных работ для написания программ на языке MATLAB можно использовать любой текстовый редактор, работающий с ASCII-кодами (например, «Блокнот», текстовый редактор файлового менеджера FAR или NORTON COMMANDER). В состав системы MATLAB 5.x и выше входит специальный редактор/отладчик m-файлов MEdit.

В ходе выполнения каждой работы оформляется отчет, который должен содержать:

вариант и формулировку задания (по каждому пункту);

краткие теоретические сведения по изучаемой проблеме;

описание разрабатываемого или исследуемого алгоритма и программы;

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

выводы по каждому пункту работы.

Отчет оформляется согласно требованиям настоящего пособия в печатном или электронном виде с применением редактора MS Word и присылается на проверку в ТМЦ ДО ТУСУР. На титульном листе отчета студент должен указать свои ФИО, идентификатор, пароль и номер варианта.

Вариант задания к лабораторным работам определяется по формуле

N_в=(N*k) div 100,

где N_в – искомый номер варианта;

N – общее число возможных вариантов;

k – две последних цифры пароля (число от 00 до 99);

div – операция целочисленного деления.

Если получается N_в =0, то берется вариант номер N.

Таблицы вариантов заданий для всех лабораторных работ приведены в приложении А в конце данного пособия.

4.2 Лабораторная работа № 1. Знакомство с пакетом для математических и инженерных расчетов MATLAB

4.2.1 Вводные замечания

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

Своим названием (MATrix LABoratory) система MATLAB обязана ориентации на матричные и векторные операции [14-16]. Система содержит также средства работы с графикой, методы численного анализа и оптимизации и многое другое.

Система MATLAB может работать как в режиме непосредственных вычислений (диалоговый режим), так и в режиме операционной среды, т.е. программу на языке MATLAB можно запускать на выполнение из командного окна среды. Первое знакомство с системой можно получить, запустив в командном окне MATLAB файл demo.m, который показывает основные возможности системы. Для этого, получив приглашение к вводу (символ «), необходимо набрать имя m-файла demo и нажать Enter.

Файлы, содержащие команды на языке MATLAB, имеют расширение *.m. Они бывают двух типов: script-файлы и функции. Script-файл представляет собой самостоятельную программу, написанную на языке системы MATLAB. Функции (function) представляют собой встроенные в систему или написанные пользователем подпрограммы, которые вызываются обращением к ним по имени с указанием входных аргументов и возвращают вектор, матрицу, скалярное или символьное значение. Информацию по каждой из внутренних функций системы MATLAB можно получить с помощью команды help <имя функции>. Команда help без параметра выдает информацию о всех функциях системы.

4.2.2 Задание на лабораторную работу

I часть – Изучение команд и операторов MATLAB

Ознакомиться с основными командами системы MATLAB:

Справочные команды - DEMO, WHATSNEW,VERSION, HELP, HELPWIN,

Управление командным окном – СLC, CLF, CLEAR, CD, DIR HOME, PAUSE, DISP, PATH, QUIT, FLOATS.

Арифметические и логические операторы – +, -, *, /, \, ^, ’, <, >, <=, >=, = =, ~ =, &, |, ~, ANS, PI, I, J.

Операторы ввода, вывода данных и операторы циклов – INPUT, PAUSE, ERROR, FPRINTF, SPRINTF, FOR…END, WHILE…END, IF…ELSE…END, BREAK, MENU.

Основные команды для работы с графикой – FIGURE, CLOSE, HOLD, SUBPLOT, PLOT, LOGLOG, POLAR, SEMILOGX,SEMILOGY, MESHGRID, PLOT3, MESH,SURF, CONTOUR.

Оформление графиков – GRID, TEXT, TITLE, XLABEL, YLABEL, GTEXT, CLABEL, AXIS.

Замечания:

a) В системе MATLAB реализовано два типа арифметических операций:

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

поэлементные операции над массивами.

Чтобы различать эти операции, перед знаком поэлементной операции ставится точка (для операторов *, /, \, ^).

Пример:

1) X=A./B – поэлементное деление матрицы А на матрицу В. Размеры матриц Aи B должны быть одинаковыми. Результатом является массив с элементами .

2) X=А/В – деление матрицы А на матрицу В по правилам линейной алгебры (это эквивалентно X=АВ^–1, где В^–1 – обратная матрица).

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

в) Логические операторы &, |, ~ соответствуют операторам булевой алгебры AND, OR, NOT.

II часть – Работа с системой MATLAB в режиме непосредственных вычислений

Запустить файл demo.m и ознакомится с основными возможностями системы MATLAB. Проделать в режиме непосредственных вычислений следующие операции:

1) Работа с матрицами и векторами.

Ввод векторов, матриц, строковых переменных (из командной строки и с помощью команды INPUT).

Вывод произвольного элемента матрицы (вектора) или нескольких произвольных строк или столбцов. Формирование матрицы меньшего размера из исходной матрицы.

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

Операции поэлементного умножения, деления, возведения в степень векторов и матриц.

Решение систем линейных уравнений вида AX=B с помощью стандартного решателя системы (использование оператора \).

2) Работа с функциями.

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

Создать файл-функцию (m-файл типа function).

Построить график y = f(x). Добавить подписи к осям и заголовок графика.

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

Построить 3D-график произведения заданной функции, зависящей от переменной x, на функцию sin(y)/y (следует пользоваться функциями MESHGRID и MESH).

III часть - Программирование в среде пакета MATLAB

1) Составить программу на языке MATLAB (script-файл), выполняющую матричные операции и работу с заданной функцией (см. п.п. 1 и 2 раздела II). Программа должна выполнять следующие действия:

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

2) Написать программу построения фигуры Лиссажу. Отношение частот k1/k2 указано в таблице вариантов заданий (см. приложение А). Фигура Лиссажу определяется комплексной функцией вида y(x)=cos(k1*x)+j*sin(k2*x), где – мнимая единица. При этом независимая переменная x изменяется в пределах от 0 до 2.

3) Написать программу для построения трехмерной поверхности функции z=f(x, y) следующего вида z = f(x).*sin(y)./y, где f(x) берется из табл. 1 (см. приложение А). Задать цветовую шкалу для значений z, используя команды COLORMAP и COLORBAR.

4) Написать программу для проверки генератора случайных чисел. На одном графическом экране построить гистограммы нормального и равномерного законов распределений (использовать операторы RAND, RANDN, SUBPLOT и HIST).

Отчет по данной работе представляет собой программу (один или несколько m-файлов) на языке MATLAB. Листинг программы в письменном виде или саму программу в виде m-файла следует прислать для проверки в ТМЦ ДО ТУСУР. Программа должна содержать все необходимые комментарии, объясняющие работу отдельных фрагментов и основных операторов. Также в заголовке программы в комментариях студент должен указать свои ФИО, идентификатор, пароль и номер варианта.

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

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

4.3.1 Вводные замечания

Система алгебраических линейных уравнений (СЛАУ) имеет вид:

(1)

Также СЛАУ можно записать в матричной форме

Ax=b,

где x=[x₁,x₂,…,x_n]^T– вектор-столбец неизвестных; b = [b₁,b₂,…,b_n]^T – вектор-столбец коэффициентов правой части системы (вектор свободных членов); А – матрица nxn коэффициентов левой части СЛАУ. Значок ^Т означает транспонирование матрицы (вектора), т.е. замену строки соответствующим столбцом.

В общем случае система (1) будет иметь единственное решение, если выполняются следующие условия:

1) число неизвестных x_i () равно числу уравнений, т.е. матрица коэффициентов А должна быть квадратной;

2) определитель матрицы коэффициентов A не равен нулю, т.е. detA0.

Если detA=0, то в системе (1) есть линейно зависимые уравнения и СЛАУ в общем случае не имеет единственного решения. Система уравнений, у которой определитель матрицы коэффициентов равен нулю (detA=0), называется вырожденной.

Для решения СЛАУ используются два класса методов: прямые и итерационные. Прямые методы являются универсальными и используются для решения СЛАУ сравнительно невысокого порядка (n  100…200). Прямые методы – это методы, которые теоретически позволяют получить точное решение системы уравнений за конечное число действий. Они используют аналитические выражения, которые дают точный результат, если не делать округления чисел в процессе расчетов. К прямым методам относятся, в частности, методы Гаусса обратной матрицы, Крамера и др.

На практике при выполнении вычислений на ЭВМ неизбежны ошибки округления. Поэтому решение СЛАУ на ЭВМ с использованием прямых методов осуществляется с определенной погрешностью, которая возрастает с увеличением порядка системы уравнений.

В данной лабораторной работе студентам предлагается реализовать наиболее распространенный метод решения СЛАУ – метод Гаусса. Этот метод при реализации в виде программы требует гораздо меньшего числа вычислений, чем методы обратной матрицы и Крамера. Метод Гаусса применим только к СЛАУ с невырожденной матрицей коэффициентов (detA0).

4.3.2 Основные положения метода Гаусса

Метод Гаусса основан на приведении матрицы к треугольному виду. Это достигается последовательным исключением неизвестных из уравнений системы. Сначала при использовании первого уравнения в качестве ведущего исключается x₁ из всех последующих уравнений СЛАУ. Затем с помощью второго уравнения исключается x₂ из третьего и всех последующих уравнений и т.д. Этот процесс называется прямым ходом. Он продолжается до тех пор, пока в левой части последнего (n-го) уравнения не останется лишь один член с неизвестным x_n, т.е. пока матрица коэффициентов А системы не будет приведена к треугольному виду.

На произвольном k-м шаге прямого хода формулы для пересчета коэффициентов имеют следующий вид:

; ; .

Это основные формулы для прямого хода метода Гаусса.

Обратный ход метода Гаусса состоит в последовательном вычислении искомых неизвестных: решая последнее уравнение приведенной системы, находим неизвестное x_n. Далее, используя это значение, из предыдущего (n–1)-го уравнения вычисляем x_n_–₁, и т.д. Последним найдем x₁ из первого уравнения.

Общая формула расчета переменных x_k () имеет следующий вид

В процессе исключения неизвестных приходится выполнять деление неизвестных на коэффициенты a₁₁, a₂₂, … и т.д. Поэтому они должны быть отличными от нуля, в противном случае необходимо соответственным образом переставить уравнения СЛАУ. Перестановка уравнений должна быть предусмотрена в вычислительном алгоритме при его реализации на ЭВМ.

Недостаток стандартного метода Гаусса (метода единственного деления). При малой величине главного элемента каждое уравнение системы во время прямого хода необходимо умножать на большое число , что приводит к росту ошибки округления при расчетах. Для уменьшения этой ошибки в качестве главного элемента следует выбирать наибольший элемент в строке, столбце или во всей матрице A. Это делается с помощью соответствующей перестановки столбцов и/или строк матрицы коэффициентов A и вектора свободных членов b.

4.3.3 Стандартные функции системы MATLAB для работы с СЛАУ

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

norm – вычисление нормы векторов и матриц

Синтаксис: n=norm(A,p)

Функция n=norm(A,p) вычисляет p-норму матрицы А. В качестве параметра p можно использовать значения 1, 2 или inf. Матрица А имеет размер nxm.

n=norm(A,1) – вычисляет 1-норму матрицы А: - максимальная сумма модулей элементов столбцов;

n=norm(A,inf) – вычисляет max-норму (-норма) матрицы А: - максимальная сумма модулей элементов строк;

n=norm(A,2) – вычисляет 2-норму (евклидову норму) матрицы А: .

cond – вычисление числа обусловленности матрицы

Синтаксис: k=cond(A)

Функция k=cond(A) возвращает число обусловленности матрицы А. Число обусловленности k – есть мера относительной погрешности при операциях обращения матрицы. Оно характеризует близость матрицы к вырождению. Матрицы с большим по величине числом обусловленности (близкие к вырожденной) дают большие ошибки при решении СЛАУ или обращении матрицы.

rcond – оценка числа обусловленности матрицы

Синтаксис: k=rcond(A)

Функция k=rcond(A) возвращает величину, обратную значению числа обусловленности матрицы А относительно 1-нормы. Если матрица А хорошо обусловлена, то значение k близко к единице, если матрица А плохо обусловлена, то значение k близко к нулю.

det – определитель матрицы

Синтаксис: d=det(A)

Функция d=det(A) вычисляет определитель квадратной матрицы.

inv – вычисление обратной матрицы

Синтаксис: Y=inv(A)

Функция Y=inv(A) вычисляет матрицу, обратную квадратной матрице А. В случае, когда матрица А плохо масштабирована или близка к вырожденной выдается соответствующее предупреждение.

На практике вычисление явной обратной матрицы требуется не часто. Если необходимо искать решение СЛАУ, то в среде MATLAB рекомендуется использовать решатели систем, т.е. операторы вида x=A\b или x=A/b, а не операцию x=inv(A)*b. Расчеты с помощью решателя систем выполняются быстрее и с большей точностью.

/ \ – решатели систем линейных уравнений

Синтаксис: X=B\A

X=B/A

Функция X=B\A находит решение системы уравнений вида AX=B, где А – прямоугольная матрица размера m x n и В – матрица размера n x k.

Функция X=B/A находит решение системы уравнений вида XA=B, где А – прямоугольная матрица размера n x m и В – матрица размера m x k.

lu – треугольное LU разложение матрицы

Синтаксис: [L,U,P]=lu(A)

Функция [L,U,P]=lu(A) находит разложение для квадратной матрицы А в виде трех составляющих – нижней треугольной матрицы L, верхней треугольной матрицы U и матрицы перестановок P, так что PA=LU.

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

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

4.3.4 Задание на лабораторную работу

В ходе данной лабораторной работы студенты должны выполнить следующие пункты.

Реализовать метод исключения Гаусса без выбора главного элемента в виде процедуры на языке системы MATLAB (блок-схема алгоритма приведена ниже).

Протестировать реализованную процедуру решения СЛАУ на следующих тестовых задачах:

а) СЛАУ с хорошо обусловленной матрицей коэффициентов размером 50х50. Для формирования матрицы коэффициентов А и вектора свободных членов b можно использовать генераторы случайных чисел rand и randn (например, для матрицы A=rand(N) и вектора b=randn(N,1)).

б) СЛАУ с плохо обусловленной матрицей коэффициентов размером 50х50. Для получения плохо обусловленной матрицы можно использовать функции системы MATLAB A=hilb(N) или A=invhilb(N), которые генерируют плохо обусловленную матрицу Гильберта. Матрицы с большим числом обусловленности cond(A)   (близкие к вырожденной) при использовании функции A=hilb(N) или A=invhilb(N) получаются при N  15.

Для каждой сформированной матрицы коэффициентов А вычислить определитель det(A), число обусловленности cond(A) и оценку числа обусловленности rcond(A).

Сравнить полученное решение (вектор x) с результатами расчетов с помощью стандартного решателя системы MATLAB x=A\b.

Оформить отчет по лабораторной работе, который должен включать следующие обязательные разделы:

а) описание алгоритма метода исключения Гаусса (включая блок-схему);

б) текст программы решения СЛАУ методом Гаусса с необходимым числом комментариев;

j=j+1

ввод n, [a_ij] и [b_i]

начало

k=1

i=k+1

a_kk=0

перестановка уравнений

c=a_ik/a_kk

a_ik=0; j=k+1

j<n

нет

да

b_i=b_i–c*b_k

a_ij=a_ij–c*a_kj

да

нет

i<n

i=i+1

k=k+1

k<n–1

да

нет

да

нет

x_n=b_n/a_nn

i=n–1

j=i+1; S=0

S=S+a_ij*x_j

j=j+1

j<n

да

нет

x_i=(b_i–S)/a_ii

i=i–1

i>1

да

нет

вывод

[x_i]

конец

прямой ход

обратный ход

Блок-схема алгоритма исключения Гаусса (без выбора главного

элемента) для решения СЛАУ

в) результаты тестирования процедуры;

для каждой тестовой задачи необходимо привести вектор решений х, полученный своей программой и стандартным решателем MATLAB, значение определителя det(A) и величину числа обусловленности cond(А) матрицы коэффициентов и оценку числа обусловленности rcond(A);

г) выводы о проделанной работе.

Отчет по данной лабораторной работе выполняется с помощью текстового редактора MSWord или письменно. Отчет и написанную программу (m-файл) следует прислать для проверки в ТМЦ ДО ТУСУР. Программа должна содержать все необходимые комментарии, объясняющие работу отдельных фрагментов и основных операторов. Также в заголовке программы в комментариях студент должен указать свои ФИО, идентификатор и пароль.

4.4 Лабораторная работа № 3. Линейная полиномиальная интерполяция функций

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

4.4.1 Теоретические сведения

Интерполяция – это нахождение функции, которая в некоторых заранее указанных точках совпадает с заданной функцией.

Пусть на отрезке [a, b] задана n+1 точка x₀, x₁, … , x_n (узлы интерполяции). Пусть в этих точках известны значения интерполируемой функции f(x):

y_i = f(x_i), .

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

(x_i) = y_i, . (4.1)

Геометрически это означает, что нужно найти кривую y=(x) определенного вида, проходящую через заданную систему точек (x_i, y_i), (см. рис. 4.1). Функция (x) называется интерполирующей функцией.

Рассмотрим задачу линейной интерполяции. При этом интерполирующая функция имеет следующий вид:

, (4.2)

где ₀(x), ₁(x), … , _m(x) – базисные функции, вид которых задан; C_k – искомые коэффициенты.

Рисунок 4.1 – Интерполяция исходной функции f(x) интерполирующей

функцией (x)

y_i

0 x₀=a x₁ x₂ x_i x_i₊₁ x_n=b x

(x)

f(x)

Используя условие (4.1) и выражение (4.2), получаем систему уравнений

(4.3)

Единственное решение системы (4.3) существует при двух условиях:

число точек (x_i, y_i), равно числу коэффициентов С_k, , т.е. m = n.

система уравнений (4.3) должна быть невырожденной, т.е. определитель системы  0.

Таким образом, если выполняются вышеуказанные условия, то через точки (x_i, y_i) проходит единственная функция вида .

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

В случае линейной полиномиальной интерполяции базисные функции имеют следующий вид: ₀(x)= x⁰=1, ₁(x)= x¹=x, ₂(x)= x², … , _m(x)= x^m. Интерполирующая функция при этом имеет вид полинома степени m: (x)= P_m(x)=C₀ + C₁x + C₂x ² + … +C_mx ^m. При m = n система (4.3) примет вид

(4.4)

В матричной форме систему (4.4) можно переписать как А*C=B, где

– матрица Ван дер Монда; ; .

Решением системы (4.4) будет вектор коэффициентов полинома С. Так как определитель матрицы Ван дер Монда всегда отличен от нуля (при x_ix_j), то решение системы (4.4) – единственное. Для решения системы (4.4) необходимо найти обратную матрицу A. В этом случае решением (4.4) будет C=A^–1*B.

Вывод: Таким образом, через заданные на интервале [a, b] точки (x_i, y_i), всегда можно провести единственный интерполяционный полином степени n P_n(x)=C₀ + C₁x + + C₂x² + … +C_nxⁿ, коэффициенты которого находятся в результате решения системы (4.4).

Выражение (4.1) определяет поведение функции (x) только в узлах интерполяции (x_i, y_i), . Между узлами (x) может вести себя произвольным образом, сколь угодно далеко, в принципе, отклоняясь от зависимости f(x). Определить погрешность приближения можно, используя выражение для абсолютной ошибки =| f(x) – (x) |.

Ошибка полиномиальной интерполяции. Лучший способ проверить качество интерполяции – вычислить значения интерполирующей функции в большом числе точек и построить график. Однако в некоторых ситуациях качество интерполянта можно проанализировать. Предположим, что величина y_i представляет собой точные значения известной функции f(x) в точках x_i. Пусть P_n(x) – единственный полином n-й степени, интерполирующий функцию по этим точкам (x_i, y_i), . Предположим, что во всех точках х[a, b] функция f(x) имеет (n+1)-ю непрерывную производную. Тогда можно показать, что абсолютная ошибка интерполяции (x) = | f(x)– P_n(x) | в некоторой точке определяется выражением

, (4.5)

где - максимальное значение (n+1)-й производной функции f(x) на интервале [a, b]; .

Рисунок 4.2 – Интерполяция функции Рунге

полиномом степени n

Теперь посмотрим, что получится, если интерполировать известную функцию f(x) все в большем и большем числе точек на фиксированном интервале. Выражение для погрешности (4.5) состоит из трех разных частей: факториал и произведение разностей _h(x) с увеличением n уменьшают ошибку, но порядок производной при этом растет. Для многих функций величина M_n₊₁ увеличивается быстрее, чем (n+1)!. В результате полиномиальные интерполянты далеко не всегда сходятся к обычной непрерывной функции. Практический эффект выражается в том, что интерполирующий полином высокой степени может вести себя очень плохо в точках, отличных от узлов интерполяции (x_i, y_i),

. Поэтому на практике обычно используют интерполирующие полиномы степени не выше 5-6.

Примером может служить функция Рунге вида f(x)=1/(1+25x²), график которой представлен на рис. 4.2. С увеличением порядка интерполирующего полинома P_n(x) при равномерном распределении узлов интерполяции на интервале [–1, 1] происходит ухудшение качества приближения на краях интервала. Это объясняется тем, что производные функции Рунге, которые фигурируют в выражении для погрешности интерполяции (4.5), быстро растут с увеличением числа n.

Формула (4.5) показывает, что точность приближения зависит не только от числа узлов интерполяции (т.е. порядка интерполирующего полинома), но и от их расположения на интервале [a, b]. В простейшем случае выбирается равномерное расположение точек x_i, на интервале [a, b] с шагом x = (b–a)/n. Однако, как показывает практика, равномерное расположение не является оптимальным с точки зрения лучшего приближения (x) к зависимости f(x). Более оптимальным для полиномиальной интерполяции является расположение узлов на интервале [a, b] по формуле

, . (4.6)

Выражение (4.6) определяет так называемое оптимальное распределение узлов интерполяции на интервале [a, b].

4.4.3 Содержание лабораторной работы

При выполнении лабораторной работы используется программа-макет polinom.m, написанная в среде пакета для математических и инженерных расчетов MATLAB. Программа предназначена для решения задачи линейной полиномиальной интерполяции и позволяет вывести графики исходной f(x), приближающей (x)=P_n(x) функций и ошибки интерполяции (x), коэффициенты полинома P_n(x), значение максимальной ошибки на интервале аппроксимации [a, b].

Студентам предлагается выполнить следующие задания в ходе данной лабораторной работы:

1) Изучить предлагаемую программу-макет polinom.m (текст программы см. в Приложении Б). Исходные интерполируемые функции (см. табл. 2 в приложении А) записать в виде m-файлов со структурой function.

2) Исследовать точность приближения с помощью полинома P_n(x) следующих трех функций (см. табл. 2 в приложении А):

а) функция 1 из таблицы 2;

б) полином P_m(x) из таблицы (выполнить расчет для случаев n<m, n=m, n>m, где n - порядок интерполирующего полинома);

в) функция 2 из таблицы 2.

Исследования выполнить при различном порядке n интерполирующего полинома P_n(x) и способе распределения узлов на интервале. Например, n = 3, 4, 5, 7, 10 – не менее 5-6 значений n для каждой из функций.

3) По результатам испытаний для каждой исходной функции заполнить таблицу или построить график зависимости максимальной ошибки _max от порядка интерполирующего полинома n (числа узлов). Кроме этого, в случае 2б для полинома P_m(x) заполнить таблицу 4.1 и сравнить коэффициенты исходного P_m(x) и интерполирующего P_n(x) полиномов для случаев n<m, n=m, n>m).

Таблица 4.1 – Исследование интерполяции функции P_m(x)

Порядок полинома P_n(x)

Равномерное распределение узлов

Оптимальное распределение узлов

коэффициенты полинома P_n(x)

_max

коэффициенты полинома P_n(x)

_max

n<m

C₀ =

C₁ =

…

C_n =

n=m

C₀ =

C₁ =

…

C_n =

n>m

C₀ =

C₁ =

…

C_n =

4) Оформить отчет по лабораторной работе. Отчет должен содержать следующие обязательные разделы:

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

Вариант задания – номер варианта, приближаемые функции и интервалы интерполяции (из табл. 2 в приложении А).

Графики исходной функции f(x), интерполирующего полинома P_n(x) и ошибки (x) на интервале интерполяции – для всех заданных приближаемых функций f(x) для одного из значений n.

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

Заполненная таблица 4.1.

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

4.5 Лабораторная работа № 4. Применение метода наименьших квадратов для сглаживания и выравнивания экспериментальных данных

Выполнение проектных операций и процедур в системах автоматизированного проектирования (САПР) основано на оперировании математическими моделями проектируемых объектов или технологических процессов.

Для решения задач, связанных с получением и хранением моделей, в состав современных САПР входят специальные информационные подсистемы. Важным разделом математического обеспечения таких подсистем является метод наименьших квадратов, который используется, в частности, для обработки результатов экспериментального определения параметров объектов и процессов при получении их математических моделей.

Цель работы: изучение методики сглаживания и выравнивания экспериментальных данных и получения математических моделей объектов и процессов по методу наименьших квадратов.

4.5.1 Теоретические сведения

Пусть данные некоторого эксперимента представлены в виде таблицы значений переменных x и y:

Таблица 4.2 – Исходные данные для аппроксимации

x_i

x₁

x₂

…

x_n

y_i

y₁

…

y_n

0 x₀=a x₁ x₂ x_i x_i₊₁ x_n=b x

Рисунок 4.3 – Аппроксимация табличных данных (x_i, y_i) функцией (x)

y_i

(x)= С₀+С₁x

_i

где x – независимая, а y – зависимая переменные.

Можно поставить задачу об отыскании аналитической зависимости между x и y, т.е. некоторой формулы y = (x), явным образом выражающей y как функцию от x. Естественно требовать, чтобы график искомой функции y = (x) изменялся плавно и не слишком уклонялся от экспериментальных точек (x_i, y_i). Поиск такой функциональной зависимости называют «сглаживанием» или «выравниванием» экспериментальных данных (см. рис. 4.3):

В частности, рассматриваемая задача встречается при построении математических моделей объектов и процессов, характеристики которых получены путем экспериментальных измерений. В этом случае требуется найти математическую модель y = (x), которая достаточно хорошо будет приближать экспериментальные характеристики (x_i, y_i) реального объекта или процесса. Функцию y = (x) иногда называют эмпирической зависимостью (моделью).

Если для получения модели y = (x) использовать глобальную интерполяцию функции в узлах x_i, с помощью полинома, то при большом числе узлов степень интерполяционного полинома будет высока, это приведет к усложнению вычислений и «колебаниям» y = (x) между узлами. Между тем, на практике даже при большом количестве точек обычно нет необходимости в использовании полиномов высокой степени.

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

Задачу сглаживания экспериментальных данных можно решать, используя метод наименьших квадратов (МНК). Согласно МНК, указывается вид аппроксимирующей функции y = (x) (математической модели):

y = (x, С₀, С₁, С₂, … , С_m) =, (4.7)

где С₀, С₁, С₂, … , С_m – неизвестные числовые параметры (коэффициенты модели).

Например, часто используется функция в виде обобщенного полинома (линейной комбинации заданных базисных функций y = _j(x), ):

. (4.8)

Если выбрать систему базисных функций _j(x) = x ^j, , то функция будет представлять собой обычный алгебраический полином:

Аналогично, функция (4.8) может описывать тригонометрические, экспоненциальные полиномы и т.д.

При использовании аппроксимирующей функции в виде (4.8) задача нахождения коэффициентов C_j, называется линейной задачей дискретной среднеквадратической аппроксимации или линейной задачей МНК.

Наилучшими значениями параметров С₀, С₁, С₂, … , С_m считают те, для которых сумма квадратов уклонений функции в точках x_i от экспериментальных значений y_i () является минимальной, т.е. функция

(4.9)

достигает минимума. Таким образом, в методе наименьших квадратов минимизируется сумма квадратов разностей между значениями функции и измеренными значениями y_i в точках x_i:

, (4.10)

где _i – разность в i-ой точке значений функции и y_i. Величина  называется среднеквадратической ошибкой. Выбор функции ошибки в виде (4.10) удобен с практической точки зрения, так как делает функцию , построенную по экспериментальным данным, нечувствительной к случайным ошибкам измерений.

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

 = f (С₀, С₁, … , С_m).

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

. (4.11)

Дифференцируя функцию ошибки (4.9) по переменным С₀, С₁, … , С_m, согласно (4.11), получим систему уравнений следующего вида

(4.12)

Заметим, что в общем случае система (4.12) является нелинейной.

Далее будем полагать, что функция является линейной относительно параметров С₀, С₁, … , С_m, т.е. имеет вид (4.8).

Из (4.8) следует, что . С учетом этого система (4.12) примет следующий вид (общий множитель 2 в уравнениях опустим):

(4.13)

Полученная система уравнений является системой (m+1) линейных алгебраических уравнений с (m+1) неизвестным С₀, С₁, … , С_m и называется нормальной системой. Ее решение определяет искомые коэффициенты функции y=.

При решении систем линейных уравнений на ЭВМ обычно используют стандартные подпрограммы. В этом случае систему уравнений необходимо представить в стандартной матричной форме:

А*С=В, (4.14)

где А – квадратная матрица размером (m+1)(m+1), составленная из коэффициентов системы линейных уравнений; С – вектор-столбец, содержащий (m+1) неизвестных переменных С₀, С₁, … , С_m; В – вектор-столбец свободных членов системы линейных уравнений.

Представим (4.13) в матричном виде (4.14). Для этого сгруппируем слагаемые в уравнениях (4.13) таким образом, чтобы каждое уравнение с номером j () приняло следующий вид:

В результате получим систему уравнений, представляющую собой развернутую запись системы (4.14):

(4.15)

Здесь - элементы матрицы А; первый индекс (j) коэффициента a_jk равен номеру уравнения, а второй (k) – номеру соответствующей неизвестной переменной; – элементы вектора–столбца свободных членов B; .

Система линейных уравнений (4.15) при невысоком порядке может быть решена вручную по правилу Крамера или с помощью метода подстановки. Решение же систем линейных уравнений на ЭВМ чаще всего проводится с использованием метода исключения неизвестных (методы Гаусса или LU-разложения).

4.5.2 Выравнивание экспериментальных данных на основе МНК

Изучим подробнее случай, когда аппроксимирующая зависимость (4.7) имеет два неизвестных коэффициента: y = (x, С₀, С₁). Используя выражения (4.12) и опуская несложные выкладки, получим систему двух уравнений с двумя неизвестными С₀, С₁:

(4.16)

Рассмотрим, в частности аппроксимацию экспериментальных данных с помощью линейной функции (см. рис. 4.3)

y = (x, С₀, С₁) = С₀+ С₁x. (4.17)

В этом случае, очевидно, , и система (4.16) примет вид:

 (4.18)

Систему (4.18) можно записать в матричном виде (4.14), где

Теперь коэффициенты С₀, С₁ легко могут быть найдены из решения (4.18) по методу Крамера.

Среднеквадратичная аппроксимация с помощью линейной функции (прямой линии) (x) = С₀+ С₁x целесообразна, очевидно, тогда, когда экспериментальные данные (x_i, y_i) приближенно описывают линейную зависимость y от x (об этом можно судить непосредственно по расположению точек (x_i, y_i) на координатной плоскости).

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

X=₁(x, y); Y=₂(x, y) (4.19)

так, чтобы преобразованные экспериментальные данные

X_i=₁(x_i, y_i); Y_i=₂(x_i, y_i) (4.20)

в новой системе координат X, Y давали точки (Х_i, Y_i), менее отклоняющиеся от линейной зависимости. Для аппроксимирующей прямой в новой системе координат

коэффициенты и можно определить из уравнений (4.18), где вместо (x_i, y_i) подставляют соответствующие значения (Х_i, Y_i). Нахождение подходящих зависимостей (4.20) называют выравниванием экспериментальных данных.

После нахождения функции находят функцию y = (x, С₀, С₁), выполнив обратную замену переменных. Функция y = (x) определена неявно уравнением

Явное выражение для (x) получают их этого уравнения, решая его относительно y (заметим, что это возможно не всегда).

Пример. Экспериментальная зависимость y от x представлена таблицей:

x_i

y_i

7,1

27,8

62,1

110

161

Установить вид эмпирической формулы, описывающей эту зависимость, используя аппроксимирующую функцию y = (x, С₀, С₁) с двумя параметрами С₀ и С₁. Определить наилучшие значения коэффициентов С₀ и С₁.

Решение. Легко убедиться в том, что экспериментальные точки (x_i, y_i) не располагаются вблизи прямой. Произведем замену переменных X = ln x, Y = ln y и составим таблицу экспериментальных данных в новых переменных (X_i, Y_i)

X_i =ln x_i

0,000

0,693

1,099

1,386

1,609

Y_i =ln y_i

1,960

3,325

4,129

4,700

5,081

Теперь точки (X_i, Y_i) лежат приблизительно на прямой (см. рис. 4.4). Значения коэффициентов и эмпирической зависимости Y=+X (в новых переменных) находятся из системы уравнений (4.18)



Рисунок 4.4 – Исходные точки и график аппроксимирующей

функции в логарифмическом масштабе

Решив эту систему, получим

= 1,97,

= 1,95.

Найдем теперь связь между исходными переменными x и y в виде функции y = (x, С₀, С₁), выполнив обратную замену переменных. Неявное уравнение, выражающее связь между x и y, имеет вид

ln y = +ln x = 1,97 + 1,95 ln x.

Отсюда легко получить явную зависимость между x и y в виде степенной функции, выполнив антилогарифмирование:

Обозначим . В результате получим эмпирическую формулу в виде

y = (x, С₀, С₁) = = 7,16 x^1,95, (4.21)

где , .

Сравнение экспериментальных данных с результатами вычислений по формуле (4.21) в соответствующих точках представлены в виде таблицы

x_i

y_i

7,1

27,8

62,1

110

161

y = 7,16 x^1,95

7,16

27,703

61,081

107,04

165,39

Как видно, эмпирическая формула (4.21) хорошо приближает эти данные.

Заметим, что параметры и в зависимости (4.21) в принципе можно найти из решения нелинейных уравнений (4.16). Однако применение способа выравнивания существенно упрощает вычисления этих параметров. В данном случае , С₁=.

Часто используемые аппроксимирующие зависимости y = (x, С₀, С₁) с двумя параметрами С₀ и С₁ приведены в таблице 4.3. Здесь же указано преобразование переменных для приведения зависимости к линейному виду Y=+X, а также связь коэффициентов С₀ и С₁ с коэффициентами и .

Таблица 4.3 – Вид эмпирических формул для приведения данных к линейной

зависимости

№

Преобразование переменных

Эмпирическая формула

(x, С₀, С₁)

X = x, Y = xy

, С₀=, С₁=

X = x, Y =

, С₁=, С₀=

X = x, Y =

, С₁=, С₀=

X = x, Y = ln y

, С₀=, С₁=

X = ln x, Y = y

y = С₁ ln x + С₀, С₀=, С₁=

X = ln x, Y = ln y

, С₀=, С₁=

Естественно, наряду с шестью предложенными формулами преобразования к переменным X, Y следует проверить возможность применения линейной зависимости непосредственно к исходным данным (x_i, y_i). Условием выбора наилучшей эмпирической формулы является наименьшее уклонение исходных или преобразованных экспериментальных данных от прямой. Уклонение данных от прямой в каждом варианте выравнивания будем определять величиной

. (4.22)

Для наилучшей эмпирической формулы величина d является наименьшей.

4.5.3 Задание на лабораторную работу

Варианты заданий к лабораторной работе представлены в табл. 3 приложения А. Студентам необходимо выполнить следующие действия:

1) нанести на координатную сетку x, y экспериментальные точки (x_i, y_i);

2) выбрать одну из шести предложенных формул преобразования к переменным X, Y так, чтобы преобразованные экспериментальные данные (X_i, Y_i) наименее уклонялись от прямой линии;

3) методом наименьших квадратов найти значения коэффициентов и в уравнении прямой Y=+X;

4) найти явный вид эмпирической формулы y = (x, С₀, С₁) и построить график эмпирической функции на плоскости x, y, оценить максимальную ошибку приближения;

5) оформить отчет по лабораторной работе.

Указания: в ходе выполнения лабораторной работы необходимо написать программу на языке системы MATLAB, реализующую вычисление коэффициентов С₀ и С₁ функции y = (x, С₀, С₁) на основе МНК. Выбор формулы преобразования можно сделать как автоматически (по величине уклонения d (4.22) от прямой), так и вручную – т.е. пользователь сам выбирает тип преобразования на основе оценки расположения точек (x_i, y_i) на плоскости x, y или точек (X_i, Y_i) на плоскости X, Y.

4.5.4 Содержание отчета

Отчет по данной лабораторной работе должен содержать следующие обязательные разделы:

1) цель работы;

2) вариант задания и постановку задачи;

3) математические выкладки при получении системы уравнений, результирующие формулы для вычисления коэффициентов С₀ и С₁;

4) текст программы вычисления коэффициентов С₀, С₁ по заданным экспериментальным точкам (x_i, y_i) с достаточным для понимания количеством комментариев;

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

величина уклонения d для выбранного способа выравнивания;
значения коэффициентов С₀, С₁ для эмпирической функции y = (x, С₀, С₁);
график функции y = (x, С₀, С₁) и точек (x_i, y_i) на плоскости x, y;
график функции Y=+X и точек (X_i, Y_i) на плоскости X, Y;

6) выводы о проделанной работе.

Отчет по лабораторной работе выполняется письменно или в электронном виде с применением текстового редактора MS Word. Программа и сам отчет должны быть присланы в ТМЦ ДО ТУСУР для проверки преподавателем.

5 Курсовое проектирование

Целью курсового проектирования является: углубленное изучение методов и алгоритмов решения различных вычислительных задач; получение навыков в разработке программ для решения вычислительных задач.

5.1 Вводные указания

Выбор языка программирования для реализации программы остается за разработчиком (Pascal, C, C++, Java и др.). Готовая программа должна быть самостоятельным приложением (т.е. должна быть представлена в виде самостоятельно исполняемого exe-файла). Приветствуется реализация в среде ОС Windows с использованием сред визуального программирования Delphi, Visual C++, C Builder, J Builder и пр.

По результатам работы оформляется пояснительная записка (ПЗ) в соответствие с требованиями изложенными ниже. Разработанная программа, исходные тексты и сама ПЗ присылается для проверки в ТМЦ ДО ТУСУР.

Таблица – Варианты заданий на курсовую работу

№ варианта

№ задачи

Вариант задания определяется по формуле

N_в=(N*k) div 100,

где N_в – искомый номер варианта;

N – общее число возможных вариантов (N=17);

k – две последних цифры пароля (число от 00 до 99);

div – операция целочисленного деления с отбрасыванием остатка.

Если получается N_в = 0, то берется последний вариант N.

Задачи 18, 19 и 20 предлагаются студентам для самостоятельного выбора без учета № варианта.

5.2 Темы курсовых проектов

I Линейная алгебра

1) Решение СЛАУ методом Гаусса с выбором главного элемента в строке или столбце.

2) Решение СЛАУ методом простых итераций.

3) Решение СЛАУ итерационным методом Гаусса-Зейделя.

4) Обращение матрицы на основе метода Гаусса с выбором главного элемента в строке или столбце.

Основные требования к проекту:

программа должна быть оформлена в виде самостоятельной процедуры (подпрограммы);
входные данные (размер системы уравнений, матрицы коэффициентов А и В) задаются с клавиатуры и из внешнего текстового файла;
результат решения (значения искомых переменных) выводится на экран и во внешний текстовый файл;
для тестирования вычислительной процедуры написать программу-оболочку, с помощью которой можно быстро изменять входные данные и просматривать результаты;
исследовать решение СЛАУ размером 10х10 и 50х50, при этом для получения матриц А и В следует использовать генератор случайных величин;
проверить работоспособность программы на плохо обусловленной (близкой к вырождению) задаче. Для этого сформировать входную матрицу коэффициентов А размером не менее 30х30 с одним или двумя маленькими диагональными элементами (порядка 10^–6-10^–8);
сравнить результаты, полученные с помощью разработанной программы, с результатами расчетов в системах для математических расчетов MATLAB или MathCAD.

II Приближение функций и аппроксимация данных

5) Аппроксимация экспериментальных данных на основе метода наименьших квадратов – в качестве аппроксимирующей функции использовать полином степени m, т.е. y = (, х) = C₀ + C₁ x + C₂ x² + … + C_m x^m

6) Построение графика функции y=f(x) с использованием кубических сплайнов

Основные требования к проекту:

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

III Численное интегрирование функций

7) Программа расчета определенных интегралов методом прямоугольников (правых или левых и центральных).

8) Программа расчета определенных интегралов методом трапеций.

9) Программа расчета определенных интегралов методом Симпсона (парабол).

Основные требования к проекту:

исходными данными для программы является вид подынтегральной функции, требуемая точность  и пределы интегрирования а и b; результат – значение интеграла с указанной точностью;

процедура интегрирования должна быть оформлена в виде самостоятельной подпрограммы;

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

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

выполнить тестирование подпрограммы на 3-4 функциях, для тестирования можно использовать стандартную таблицу интегралов или пакет для математических расчетов MATLAB или MathCAD;

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

IV Нахождение корней (нулей) нелинейных уравнений

10) Нахождение нулей функции y = f(x) методом дихотомии.

11) Нахождение нулей функции y = f(x) методом хорд.

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

13) Нахождение нулей функции y = f(x) модифицированным методом Ньютона.

14) Нахождение нулей функции y = f(x) методом секущей.

Основные требования к проекту:

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

выходными данными программы будут найденные на указанном интервале все корни функции f(x)=0;

процедура нахождения нулей должна быть оформлена в виде самостоятельной подпрограммы;

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

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

для тестирования использовать 3-4 функции: тригонометрические функции, алгебраический полином и др. Полученные результаты сравнить с решением этой же задачи с использованием пакетов MATLAB или MathCAD;

исследовать зависимость требуемого числа итераций от точности определения корня.

V Решение системы нелинейных уравнений

15) Решение системы нелинейных уравнений (СНУ) методом Ньютона – для расчета элементов матрицы Якоби на каждом шаге использовать аппроксимацию частных производных конечно-разностными соотношениями.

16) Решение системы нелинейных уравнений методом простых итераций.

17) Решение системы нелинейных уравнений многомерным методом секущих.

Основные требования к проекту:

исходными данными к программе являются вид СНУ, начальное приближение и требуемая точность нахождения решения;

выходными данными – значения переменных, в которых СНУ имеет решение, и соответствующие значения функций ;

вычислительная процедура решения СНУ должна быть оформлена в виде самостоятельной подпрограммы;

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

решить 23 системы уравнений (взять из учебника или составить самостоятельно), исследовать зависимость числа итераций от точности нахождения решения.

VI Дополнительные задачи повышенной трудности

18) Нахождение нулей алгебраического полинома – необходимо найти все корни (в том числе и комплексные).

19) Построение графика функции двух переменных z=f(x, y) в виде поверхности на основе триангуляции или двухмерных сплайнов.

20) Аппроксимация периодической функции тригонометрическим рядом Фурье и построение спектра амплитуд и фаз гармоник.

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

Пояснительная записка (ПЗ) к курсовому проекту (работе) должна включать в себя следующие разделы:

1. Титульный лист – пример оформления титульного листа представлен в приложении В.

2. Аннотация (Реферат) – содержит краткое описание (реферат) выполненной работы (2-3 предложения). Перечисляются ключевые слова, указывается количество страниц и приложений. Реферат размещают на отдельной странице. Заголовком служит слово «Реферат», написанное прописными буквами по центру страницы.

3. Задание на проектирование – формулируется по выбранному варианту. Пример оформления задания на проектирование представлен в приложении Б.

4. Содержание ПЗ – нумерованный по страницам список разделов ПЗ. Нумерация страниц ПЗ – сквозная: титульный лист имеет первый номер. На листе с содержанием обычно рисуется рамка для текстовых документов (см. приложение В).

Замечание: по последним требованиям СТП ТУСУР рамки может и не быть, если курсовой проект не связан с разработкой какого-то конкретного устройства или блока.

5. Введение – содержит общую информацию по проекту: краткую характеристику решаемой задачи, назначение программного продукта, необходимость его разработки и пр.

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

7. Постановка задачи. В данном разделе подробно формулируется решаемая задача.

8. Анализ задания – включает в себя анализ технического задания, представленного в разделе «задание на проектирование». Здесь рассматриваются особенности решаемой задачи, вырабатываются требования к программному продукту и интерфейсу пользователя, определяются все функции программы, выбирается среда реализации и пр.

9. Основной раздел ПЗ – содержит описание хода работы над программным продуктом и результаты тестирования готовой программы. Может состоять из нескольких подразделов. Желательно отразить следующие моменты:

1) описание алгоритмов программы и используемых математических моделей (расчетные формулы, блок-схемы используемых алгоритмов, общее описание алгоритмов);

2) описание реализации программы (структурная схема программы, назначение отдельных блоков, структура и типы данных, основные процедуры, их назначение и взаимосвязь, входные и выходные данные программы и отдельных процедур, смысл основных (глобальных) переменных программы, построение файловой системы, используемые средства среды программирования и т.д.). Для объектно-ориентированного программирования приводится описание используемых в программе стандартных и разрабатываемых классов и методов;

3) описание программы для пользователя (описание интерфейса, запуск и порядок работы с программой, состав файлов, входные и выходные данные, экранные формы и пр.) – данный раздел может быть оформлен в виде инструкции для пользователя;

4) тестирование программы (описание тестовых задач, обсуждение результатов тестирования).

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

11. Список литературы – список источников, используемых при работе над проектом. Может содержать не только литературные источники, но и ссылки на различные ресурсы в сети ИНТЕРНЕТ. Пример оформления списка литературы представлен в Приложении В.

12. Приложение – обычно в приложение выносят следующее:

блок-схемы алгоритмов, структурную схему и листинг программы с комментариями, результаты тестирования (если они обширные) и т. д.

Если листинг программы занимает достаточно большой объем (более 5-6 страниц), то в текст ПЗ следует включать не весь листинг целиком, а только тексты процедур, выполняющих основные действия и расчеты.

6.1 Указания к оформлению ПЗ

ПЗ пишется в редакторе MS Word шрифтом Times New Roman, размером 12 пунктов, на листе формата A4. Нумерация страниц должна быть сквозной, первой страницей является титульный лист (номер страницы на титульном листе не ставится). Номер страницы проставляется вверху справа. Заголовки разделов пишутся прописными буквами по середине текста. Заголовки подразделов пишутся с абзаца строчными буквами, кроме первой прописной. В заголовке не допускаются переносы слов. Точку в конце заголовка не ставят. Если заголовок состоит из двух предложений, то их разделяют точкой.

Необходимо обратить внимание на следующее:

1) наличие нумерации и подписей к рисункам, нумерация и заголовки таблиц;

2) наличие заголовков приложений и их сквозная нумерация;

3) выполнение блок-схем алгоритмов и структурных схем программ в соответствии с требованиями ГОСТа (или СТП ТУСУР);

4) обязательное детальное комментирование текста программ;

5) оформление списка литературы в соответствии с требованиями ГОСТа;

6) наличие в тексте ссылок на используемую литературу (в тексте ПЗ в квадратных скобках пишется номер источника).

Оценка за курсовой проект выставляется с учетом качества выполнения программного продукта и пояснительной записки.

Приложение А

Варианты заданий к лабораторным работам

Таблица 1 – Варианты заданий на лабораторную работу № 1

№ варианта

Функция f(x)

Интервал аргумента

x для f(x)

Отношение частот

k1/k2

[-3, 3]

3/5

[2, 5]

2/5

sin(x)cos(x)

[-2, 4]

3/4

[-1, 2]

2/3

[0, 5]

3/5

[0, 100]

5/6

[-2, 3]

2/3

[1, 3]

3/7

[0, 2]

5/8

[1, 7]

3/8

Приложение А (продолжение)

Таблица 2 – Варианты заданий на лабораторную работу № 3

№ варианта

Функция 1

Интервал

Полином P_m(x)

Интервал

Функция 2

Интервал

sin(x)

3x⁵+12x–7

[0, 2]

| x–5 |

[0, 10]

cos(x)

4x⁷–5x⁵+8x³+12

[–2, 3]

| x+3 |

[–8, 5]

cos(x)+sin(x)

7x⁷–9x⁴–7x²+12

[–4, 1]

| 2–x |

[–2, 4]

cos(x)–sin(x)

12x⁶–3x⁵–7x³+11

[–1, 2]

| 2x–2 |

[–5, 8]

sin(x)–cos(x)

3x⁹–12x³–18

[–3, 2]

| 2x+1 |

[–5, 5]

sin(x)cos(x)

5x⁴–4x²+9

[–4, 7]

| 4–2x |

[0, 8]

1+sin(x)

–2x⁴+7x²–8x

[–1, 3]

| x |+5

[–9, 9]

1+cos(x)

–8x⁷+12x⁴–8x²+19

[–3, 6]

| 3x | – 4

[–7, 8]

1–sin(x)

9x⁴–7x²–21

[–4, 5]

| 5–7x |

[–4, 2]

1–cos(x)

–2x⁸+11x⁵–4x²–8

[–2, 3]

| 2x+3 |

[–6, 4]

Приложение А (продолжение)

Таблица 3– Варианты заданий на лабораторную работу № 4

№ варианта

1,1

1,4

1.6

1,7

1,9

1,05

1,55

1,7

1,75

1,8

0,4

0,55

0,13

0,09

0,07

7,5

6,2

5,5

3,5

8,2

5,9

4,9

3,2

7,2

5,9

4,9

3,2

7,1

6,1

4,9

3,1

0,55

0,7

0,77

0,82

0,85

1,1

1,55

1,9

2,3

2,6

1,1

1,55

1,9

2,25

2,5

5,1

4,4

3,2

2,7

2,55

5,1

3,4

3,2

2,7

2,55

1,9

5,5

3,5

3,67

3,75

3,8

0,25

0,09

0,07

0,05

0,04

0,25

0,111

0,071

0,053

0,042

0,20

0,28

0,33

0,36

0,38

4,8

5,76

6,912

8,294

9,95

3,08

4,3

5,16

5,83

0,33

0,5

0,6

0,67

0,71

1,5

1,75

1,83

1,87

1,9

0,2

0,11

0,077

0,059

0,4

0,33

0,31

0,29

2,25

3,37

5,06

7,59

11,4

2,69

3,1

3,39

3,61

Приложение Б

Листинг программы-макета Polinom.m (для системы MATLAB4.x)

% ----------------------------------------------------------------

% Программа линейной полиномиальной интерполяции функций

% Кафедра КСУП ФВС ТУСУР

% Автор: Черкашин М.В.

%-------------------------------------------------------------------------

clc; % --- очистка главного окна

home; % --- перемещение курсора в начало

disp('---------------------------------------------------------------- ');

disp(' ');

disp(' Линейная полиномиальная интерполяция функций ');

ddisp('---------------------------------------------------------------- ');

disp(' ');

disp(' Жми любую клавишу ...')

pause

%-------------------------------------------------------------------------

key=1;

while key~=0, % Общий цикл программы

clc; % Очистить экран

close all; % Закрыть все графические окна

clear all; % Удаление переменных из памяти

% --- входные данные

NN=100; % -- количество точек для построения графиков

% --- запрос имени исходной фукнции

clc;

disp(' Введите имя М-файла с интерполируемой функцией')

disp(' Замечание: Функция должна быть оформлена в виде function ');

disp(' т.е. должен существовать М-файл с Вашей функцией ');

disp(' подробности см. в MATLAB HELP ');

% --- цикл запроса имени функции

k=1;

while k~=0, txtFUN=input(' Введите имя файла ->','s');

if length(txtFUN)==0,

k=1;

else

if (exist(txtFUN)==2 | exist(txtFUN)==5),

k=0;

else

disp([' Ошибка: ',txtFUN,' - такой функции или файла нет... ' ]);

k=1;

end

% --- задание интервала интерполяции

disp(' Задайте интервал интерполяции [a b] ');

predel=input(' Введите [a b] ->');

key1=1;

while key1~=0, % цикл по узлам интерполяции N

% -- очистка переменных для нового цикла расчетов

i=[];j=[];ii=[];x=[];x1=[];x2=[];A1=[];A2=[];n=[];k=[];

C1=[];C2=[];P1=[];P2=[];Error1=[];Error2=[];

f=[];f0=[];f1=[];f2=[];

%-------------------------------------------------------------------------

disp(' Задайте порядок интерполирующего полинома Pn(x) ')

disp(' (Число узлов интерполяции будет N+1): ')

N=input(' Введите порядок полинома N ->');

delta=(predel(2)-predel(1))/(NN-1); % -- Шаг построения графика функции

delta1=(predel(2)-predel(1))/N; % -- Шаг интерполяции

%-------------------------------------------------------------------------

j=1:NN; % - количество точек для построения графика функции

i=1:N+1;

ii=0:N;

x(j)=predel(1)+delta*(j-1); % -- Расчет х для построения графика функции

f=feval(txtFUN,x); % -- расчет значений интерполируемой функции f(x) для

% --построения графика функции

% -- Узлы интерполяции (x-координата)

x1(i)=predel(1)+delta1*(i-1); % -- равномерное распределение

x2(ii+1)=(predel(1)+predel(2))/2+((predel(2)-predel(1))/2)*… cos(((2*ii+1)*pi)/(2*N+2)); % -- оптимальное распределение

% --- расчет функции (имя функции находится в переменной txtFUN)

f1=feval(txtFUN,x1); % -- для узлов с равномерным распределением

f2=feval(txtFUN,x2); % -- для узлов с оптимальным распределением

%-------------------------------------------------------------------------

%------- установки графического окна для отображения графика функции

[a,fig]=figflag('Программа линейной интерполяции: Графики',1);

if a==0

figN1=figure;

end

figure(figN1);

set(figN1,'Name','Программа линейной интерполяции: Графики',...

'NextPlot','add',...

'NumberTitle','off');

% --- Построение графика функции и узлов интерполяции

plot(x,f,'-g',x1,f1,'ow',x2,f2,'+w');

title(' Вид интерполируемой функции f(x) в интервале интерполяции');

xlabel(' o - равномерное х - оптимальное распределение узлов')

disp(' ');

disp(' Жми любую клавишу ...')

pause

%-------------------------------------------------------------------------

% --- Формирование матрицы системы уравнений для расчета коэф-тов

% --- интерполирующего полинома

for i=1:N+1,

for n=0:N,

A1(i,n+1)=x1(i)^n; % -- для равномерного распределения узлов

end

%------------------------------------------------------

for i=1:N+1,

for n=0:N,

A2(i,n+1)=x2(i)^n; % -- для оптимального распределения узлов

end

%-------------------------------------------------------------------------

% --- Решение системы уравнений и нахождение коэф-тов полинома

disp(' ');

disp(' Интерполяция функции f(x) полиномом вида : ');

disp(' P=C(1)*X^0+C(2)*X^1+C(3)*X^2+...+C(N+1)^N, ');

disp(' где N - порядок интерполирующего полинома ');

disp(' ');

%-------------------------------------------------------

C1=inv(A1)*f1'; % -- решение системы уравнений

disp(' Коэф-ты интерполирующего полинома P1(x) ')

disp(' (при равномерном распределении узлов) ')

for k=1:N+1,

fprintf(' C1(%g) = %g \n',k-1,C1(k));

end

%-------------------------------------------------------

C2=A2\f2'; % -- решение системы уравнений

disp(' Коэф-ты интерполирующего полинома P2(x) ');

disp(' (при оптимальном распределении узлов) ');

for k=1:N+1,

fprintf(' C2(%g) = %g \n',k-1,C2(k));

end

%-------------------------------------------------------

disp(' ');

disp(' Жми любую клавишу ...')

pause

% --- расчет интерполирующего полинома для равномерного распределения

for j=1:NN,

P1(j)=0;

for n=0:N,

P1(j)=P1(j)+C1(n+1)*x(j)^n;

end

% --- расчет интерполирующего полинома для оптимального распределения

for j=1:NN,

P2(j)=0;

for n=0:N,

P2(j)=P2(j)+C2(n+1)*x(j)^n;

end

%-------------------------------------------------------------------------

% -- вызов графического окна

[a,fig]=figflag('Программа линейной интерполяции: Графики',1);

if a==0

figN=figure;

else

figN=fig;

end

figure(figN);

set(gcf,'NextPlot','add');

%--- построение интерполирующих полиномов в виде графика

plot(x,f,'-w',x1,f1,'og',x,P1,'-.g',x2,f2,'ob',x,P2,':b');

title(' Интерполяция функции f(x) полиномом P(x)')

xlabel(' белый - f(x) зеленый - P1(x) синий - P2(x)')

pause

%-------------------------------------------------------------------------

% -- Расчет ошибки интерполяции

Error1=abs(f-P1); % -- ошибка для равномерного распределения

Error2=abs(f-P2); % -- ошибка для оптимального распределения

% -- Нахождение максимального значения ошибки

MaxErr1=max(Error1);

MaxErr2=max(Error2);

f0=zeros(1,N+1);

% --- Построение графика ошибки

[a,fig]=figflag('Программа линейной интерполяции: Ошибка',1);

if a==0

figN2=figure;

else

figN2=fig;

end

figure(figN2);

set(figN2,'Name','Программа линейной интерполяции: Ошибка',...

'NextPlot','add',...

'NumberTitle','off');

plot(x,Error1,'-g',x,Error2,'-b',x1,f0,'og',x2,f0,'ob')

title(' Ошибка интерполяции ')

xlabel(' зеленый - P1(x) синий - P2(x)')

pause

disp('Максимальная ошибка интерполяции функции f(x) :')

fprintf(' MaxErr1= %g \n MaxErr2= %g \n',MaxErr1,MaxErr2)

%-------------------------------------------------------------------------

% --- Блок запроса на повторный запуск процедуры расчета коэф-тов

k=1; % -- Начало блока запроса

while k~=0,

t=input(' Изменить число узлов интерполяции (y/n) ? ->','s');

if isempty(t)==1, t='y'; end;

txt=sscanf(t(1),'%s');

if txt=='y' | txt=='n'

k=0;

else

k=1;

end

end % -- Конец блока запроса

if txt=='y'

key1=1;

elseif txt=='n'

key1=0;

end

end % -- Конец цикла по узлам N

%-------------------------------------------------------------------------

% --- Блок запроса на повторный запуск программы

k=1; % -- Начало блока запроса

while k~=0,

t=input(' Ввести новую функцию (y/n) ? ->','s');

if isempty(t)==1, t='y'; end;

txt=sscanf(t(1),'%s');

if txt=='y' | txt=='n'

k=0;

else

k=1;

end

end % -- Конец блока запроса

if txt=='y'

key=0;

elseif txt=='n'

key=1;

end

end % -- Конец общего цикла программы whille key~=0

clc;

close all;

disp(' Вот и все ... ');

% ----------- END PROGRAM

Приложение В

Пример оформления титульного листа ПЗ

Министерство образования Российской Федерации

ТОМСКИЙ ГОСУДАРСТВЕННЫЙ УНИВЕРСИТЕТ СИСТЕМ УПРАВЛЕНИЯ И РАДИОЭЛЕКТРОНИКИ (ТУСУР)

ТМЦ ДО ТУСУР

РЕШЕНИЕ СИСТЕМ ЛИНЕЙНЫХ АЛГЕБРАИЧЕСКИХ

УРАВНЕНИЙ МЕТОДОМ ГАУССА

Пояснительная записка к курсовой работе по дисциплине

"Вычислительные методы"

Выполнил:

студентка гр.587-1

И.В.Никитина

16.05.2001 г.

Домашний адрес

Имя и пароль

Вариант

ТОМСК – 2002

Приложение В (продолжение)

Министерство образования РФ

Томский государственный университет систем управления и радиоэлектроники (ТУСУР)

ТМЦ ДО ТУСУР

Утверждаю

Директор ТМЦ ДО

Уваров А.Ф. _______________

"_____"_____________ 2002 г.

ЗАДАНИЕ

на курсовую работу по дисциплине "Вычислительные методы"

Студентке гр. 587-1 Никитиной Ирине Владимировне

Тема работы: РЕШЕНИЕ СИСТЕМ ЛИНЕЙНЫХ АЛГЕБРАИЧЕСКИХ

УРАВНЕНИЙ МЕТОДОМ ГАУССА

Исходные данные к работе:

1. Разработать программу решения систем линейных алгебраических уравнений (СЛАУ) методом Гаусса

2. Программа решения СЛАУ должна быть оформлена в виде самостоятельной вычислительной процедуры

3. При прямом ходе необходимо осуществлять выбор главного элемента в матрице коэффициентов

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

5. Результат работы процедуры – вектор искомых переменных Х. Вывод осуществляется на экран монитора и во внешний текстовый файл.

4. Выбор среды реализации предоставляется разработчику

5. При тестировании рассмотреть решение СЛАУ с матрицей коэффициентов, близкой к вырожденной, размером не менее 20х20.

Литература:

М.В. Черкашин. Вычислительные методы. Курс лекций (часть 2). Томск, ТУСУР, 2002.
В.М. Вержбицкий. Численные методы (линейная алгебра и нелинейные уравнения): учебное пособие для ВУЗов. – М.: Высшая школа, 2000. –266 с.
Л.И. Турчак. Основы численных методов. – М., Наука, 1987. – 320с.

Задание принял к исполнению _________ Никитина И.В.

1 февраля 2002 г.

Пример оформления задания на курсовой проект (работу)

Приложение В (продолжение)

Пример оформления списка литературы

Список использованных источников

Черкашин М.В. Вычислительные методы. Курс лекций (часть 2).  Томск: ТУСУР, 2003.

Вержбицкий В.М. Численные методы (линейная алгебра и нелинейные уравнения): Учебное пособие для ВУЗов. – М.: Высшая школа, 2000. –266 с.

Турчак Л.И. Основы численных методов. – М., Наука, 1987. – 320 с.

Тихонов А.Н., Костомаров Д.П. Вводные лекции по прикладной математике. – М., Наука, 1984.

Самарский А.А., Гулин А.В. Численные методы. – М., Наука, 1989. – 432 с.

Боглаев Ю.П. Вычислительная математика и программирование. – М., Высшая школа, 1990. – 544с.

Ракитин В.И., Первушин В.Е. Практическое руководство по методам вычислений с приложением программ для персональных компьютеров: Учебн. пособие. – М.: Высшая школа, 1998. - 383 с.

www.exponеnta.ru\matlab - сайт по математическим методам и программам.

www.ishodniki.ru\pascal – сайт, посвященный программированию на языке PASCAL.

Иванова Л.М. Методы решения СЛАУ // «Информатика и образование».  2000.  № 4. – С. 12-23.

Листать вверх Листать вниз

Получить код

Скачать материал (0.64 Мб)

Нравится материал? Поддержи автора!

Ещё документы из категории разное:

Содружества независимых государств

Пушкин а. с. - Природа в романе а. с. пушкина евгений онегин

Нормативы и методические указания по разработке учебных планов

Безопасности всех уровней

Отглагольные имена существительные в современной английской прессе и их перевод

О великой отчественной войне - Суровые будни войны в повести в. быкова его батальон

Международный рынок ценных бумаг