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

СОДЕРЖАНИЕ: ИНСТИТУТ ЭКОНОМИКИ, УПРАВЛЕНИЯ И ПРАВА (г. КАЗАНЬ) Кафедра высшей математики КУРСОВОЙ ПРОЕКТ ПО ВЫСШЕЙ МАТЕМАТИКЕ НА ТЕМУ: ЧИСЛЕННЫЕ МЕТОДЫ Конечно-разностный метод решения краевых задач

ИНСТИТУТ ЭКОНОМИКИ, УПРАВЛЕНИЯ И ПРАВА (г. КАЗАНЬ)

Кафедра высшей математики

КУРСОВОЙ ПРОЕКТ

ПО ВЫСШЕЙ МАТЕМАТИКЕ

НА ТЕМУ:

ЧИСЛЕННЫЕ МЕТОДЫ

Конечно-разностный метод решения краевых задач

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

Выполнил:

студент группы № _____

дневного (заочного) отделения

факультета менеджмента и маркетинга

специальность: управление качеством

_______________________________

номер зачетной книжки __________

контактный телефон: ____________

Руководитель:

доц. Лебедев В.Н.

Казань – 2010 г.

СОДЕРЖАНИЕ

Ведение ……………………………………………………………………………….3

1. Общая теоретическая часть………………………………………………………..4

1.1. Действия с приближёнными величинами……………………………………....4

1.2. Основные численные методы……………………………………………………7

1.2.1. Решение алгебраических и трансцендентных уравнений…………………...7

1.2.2. Интерполяция функций………………………………………………………..10

1.2.3. Метод наименьших квадратов и его применения……………………………11

1.2.4. Численное интегрирование……………………………………………………14

2. Конечно-разностный метод решения краевых задач для обыкновенных

дифференциальных уравнений……………………………………………………...16

3. Расчетная часть…………………………………………………………………….18

3.1. Поиск действительных корней уравнения……………………………………..18

3.2. Приближеннее вычисление интеграла………………………………………….20

3.3. Построение интерполяционных многочленов Лагранжа и Ньютона………...21

3.4. Оценки параметров линейной и квадратичной моделей……………………...22

Результаты и выводы…………………………………………………………………24

Список использованной литературы………………………………………………...24

ВВЕДЕНИЕ

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

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

Вторая крупная проб­лема – освоение космоса – связана с созданием летатель­ных аппаратов и решением для них многих задач аэроди­намики и баллистики (например, расчет движения ракеты и управление ее полетом). Здесь также имеется комплекс сложных задач механики, физики и техники, которые мо­гут быть решены только с использованием численных методов.

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

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

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

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

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

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

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

Па третьем этапе осуществляется программирование вычислительного алгоритма для ЭВМ и на четвертом эта­пе – проведение расчетов па ЭВМ. Деятельность по программированию должна быть тес­но связана с разработкой конкретных численных алго­ритмов.

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

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

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

1. Общая теоретическая часть

1.1. Действия с приближенными величинами.

Прикладные математические методы оперируют числами. Но при изучении сложных природных систем встречаются такие величины, точное значение которых неизвестно. В таких случаях принято говорить, что используются приближенные числа или величины. При работе с ними можно использовать разное число десятичных знаков. Ни к чему записывать результат измерения или вычисления с точностью, превышающей точность исходных данных или точность самого метода. Так, если расстояние между двумя точками 120м измерено с точностью 1 метр, а время 11с – с точностью 1с, то значение скорости 120/11 = 10,9м/с неправильно – следует результат округлить и записать V = 11м/с.

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

Приближенные числа всегда имеют некоторую погрешность, или ошибку. Различают абсолютные и относительные ошибки. Абсолютной ошибкой e приближенной величины называют разность между ее точным (A) и приближенным (X) значениями:

e = A – X

Размерность e совпадает с размерностью исследуемой величины. Относительной ошибкой e называют отношение абсолютной ошибки к точному значению величины. Ее размерность – доли единицы или проценты. Так как точное значение, как правило, неизвестно, вместо него используют приближенное:

e = e/X

Из сказанного следует практический вывод.

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

При измерениях физических величин неизбежны погрешности измерения

Погрешность измерения – оценка отклонения величины измеренного значения величины от её истинного значения. Погрешность измерения является характеристикой (мерой) точности измерения.

Поскольку выяснить с абсолютной точностью истинное значение любой величины невозможно, то невозможно и указать величину отклонения измеренного значения от истинного. Возможно лишь оценить величину этого отклонения, например, при помощи статистических методов. При этом за истинное значение принимается среднестатистическое значение, полученное при статистической обработке результатов серии измерений. Это полученное значение не является точным, а лишь наиболее вероятным. Поэтому в измерениях необходимо указывать, какова их точность. Для этого вместе с полученным результатом указывается погрешность измерений. Например, запись Х = 2,8 ± 0,1 означает, что истинное значение величины Х лежит в интервале от 2,7 до 2,9 с некоторой оговоренной вероятностью.

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

Абсолютная погрешность – X является оценкой абсолютной ошибки измерения. Величина этой погрешности зависит от способа её вычисления, который, в свою очередь, определяется распределением случайной величины X_ИЗМ . При этом неравенство:

X | X_ИСТ X_ИЗМ |,

где X_ИСТ – истинное значение, а X_ИЗМ – измеренное значение, должно выполняться с некоторой вероятностью близкой к 1. Если случайная величина X_ИЗМ распределена по нормальному закону, то обычно, за абсолютную погрешность принимают её среднеквадратичное отклонение. Абсолютная погрешность измеряется в тех же единицах измерения, что и сама величина.

Относительная погрешность – отношение абсолютной погрешности к тому значению, которое принимается за истинное:

Относительная погрешность является безразмерной величиной, либо измеряется в процентах.

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

,

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

По причине возникновения:

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

Методические погрешности – погрешности, обусловленные несовершенством метода, а также упрощениями, положенными в основу методики.

Субъективные / операторные / личные погрешности – погрешности, обусловленные степенью внимательности, сосредоточенности, подготовленности и другими качествами оператора.

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

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

Обобщенной характеристикой средств измерения является класс точности, определяемый предельными значениями допускаемых основной и дополнительной погрешностей, а также другими параметрами, влияющими на точность средств измерения; значение параметров установлено стандартами на отдельные виды средств измерений. Класс точности средств измерений характеризует их точностные свойства, но не является непосредственным показателем точности измерений, выполняемых с помощью этих средств, так как точность зависит также от метода измерений и условий их выполнения. Измерительным приборам, пределы допускаемой основной погрешности которых заданы в виде приведенных основных (относительных) погрешностей, присваивают классы точности, выбираемые из ряда следующих чисел: (1; 1,5; 2,0; 2,5; 3,0; 4,0; 5,0; 6,0)*10ⁿ , где показатель степени n = 1; 0; 1; 2 и т. д.

По характеру проявления:

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

Систематическая погрешность – погрешность, изменяющаяся во времени по определенному закону (частным случаем является постоянная погрешность, не изменяющаяся с течением времени). Систематические погрешности могут быть связаны с ошибками приборов (неправильная шкала, калибровка и т. п.), неучтёнными экспериментатором.

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

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

По способу измерения:

Погрешность прямых измерений – погрешностьизмеряемой величины.

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

При решении математических задач приближёнными (численными) методами также возникают погрешности вычислений.

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

1) математическое описание задачи является неточным, в частности неточно заданы исходные данные описания;

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

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

Погрешности, соответствующие этим причинам, называют:

1) неустранимой погрешностью,

2) погрешностью метода,

3) вычислительной погрешностью.

Часто неустранимую погрешность подразделяют на две части:

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

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

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

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

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

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

4. Относительная погрешность n-ой степени приближенного числа в n раз больше относительной погрешности основания (как для целых, так и для дробных n).

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

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

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

1. При сложении и вычитании приближённых чисел в результате следует сохранять столько десятичных знаков, сколько их в приближённом данном с наименьшим числом десятичных знаков.

2. При умножении и делении в результате следует сохранять столько значащих цифр, сколько их имеет приближённое данное с наименьшим числом значащих цифр.

3. При возведении в квадрат или куб в результате следует сохранять столько значащих цифр, сколько их имеет возводимое в степень приближённое число (последняяцифра квадрата и особенно куба при этом менее надежна, чем последняя цифра основания).

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

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

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

1.2. Основные численные методы

1.2.1. Решение алгебраических и трансцендентных уравнений.

Общая постановка задачи . Найти действительные корни уравнения f(x) = 0, где f(x) алгебраическая или трансцендентная функция.

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

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

1) отделение (локализация) корня;

2) приближённое вычисление корня до заданной точности.

Отделение корня. Отделение действительного корня уравнения f(x) = 0 – это нахождение отрезка [a; b], в котором лежит только один корень данного уравнения. Такой отрезок называется отрезком изоляции (локализации) корня.

Наиболее удобным и наглядным является графический метод отделения корней:

1) строится график функции y = f(x), и определяются абсциссы точек пересечения этого графика с осью OX, которые и являются корнями уравнения f(x) = 0;

2) если f(x) – сложная функция, то её надо представить в виде так, чтобы легко строились графики функций . Так как . Тогда абсциссы точек пересечения этих графиков и будут корнями уравнения

Уточнение корня. Если искомый корень уравнения отделён, т.е. определён отрезок [a; b], на котором существует только один действительный корень уравнения, то далее необходимо найти приближённое значение корня с заданной точностью.

Такая задача называется задачей уточнения корня.

Уточнение корня можно производить различными методами:

1) метод половинного деления (бисекции);

2) метод итераций;

3) метод хорд (секущих);

4) метод касательных (Ньютона);

5) комбинированные методы.

Метод половинного деления (бисекции).

Отрезок изоляции корня можно уменьшить путём деления его пополам.

Такой метод можно применять, если функция f(x) непрерывна на отрезке [a; b] и на его концах принимает значения разных знаков, т.е. выполняется условие^

(1)

Разделим отрезок [a; b], пополам точкой , которая будет приближённым значением корня .

Для уменьшения погрешности приближения корня уточняют отрезок изоляции корня. В этом случае продолжают делить отрезки, содержащие корень, пополам.

Из отрезков [a; c₁ ] и [c₁ ; b] выбирают тот, для которого выполняется неравенство (1).

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

Достоинство метода: простота (достаточно выполнения неравенства (1)).

Недостаток метода: медленная сходимость результата к заданной точности.

Метод хорд (секущих).

Этот метод применяется при решении уравнений вида f(x) = 0, если корень уравнения отделён, т.е. и выполняются условия:

1) (функция f(x) принимает значения разных знаков на концах отрезка.

2) производная сохраняет знак на отрезке [a; b], т.е. функция f(x) либо возрастает, либо убывает на отрезке [a; b].

Первое приближение корня находится по формуле:

.

Для следующего приближения из отрезков [a; х₁ ] и [х₁ ; b]выбирается тот, на концах которого функция f(x) имеет значения разных знаков.

Если , то второе приближение вычисляется по формуле:

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

Метод касательных (Ньютона).

Этот метод применяется, если уравнение f(x) = 0 имеет корень , и выполняются условия:

1) (функция принимает значения разных знаков на концах отрезка [a; b];

2) производные сохраняют знак на отрезке [a; b], т.е. функция f(x) либо возрастает, либо убывает на отрезке [a; b], сохраняя при этом направление выпуклости.

На отрезке [a; b] выбирается такое число х₀ , при котором f(x₀ ) имеет тот же знак, что и , т. е. выполняется условие . Таким образом, выбирается точка с абсциссой x₀ , в которой касательная к кривой на отрезке [a; b], пересекает ось OX. За точку x₀ сначала удобно выбирать один из концов отрезка.

Первое приближение корня определяется по формуле: .

Второе приближение корня определяется по формуле: .

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

Достоинства метода: простота, быстрота сходимости.

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

Комбинированный метод хорд и касательных.

Если выполняются условия:

1) ,

2) и сохраняют знак на отрезке [a; b],

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

Схема решения уравнения методом хорд и касательных

1. Вычислить значения функции и .

2. Проверить выполнение условия . Если условие не выполняется, то неправильно выбран отрезок [a; b].

3. Найти производные .

4. Проверить постоянство знака производных на отрезке [a; b]. Если нет постоянства знака, то неверно выбран отрезок [a; b].

5. Для метода касательных выбирается за х₀ тот из концов отрезка [a; b], в котором выполняется условие , т.е. и одного знака.

6. Приближения корней находятся:

а) по методу касательных: ,

б) по методу хорд: .

7. Вычисляется первое приближение корня: .

8. Проверяется выполнение условия: , где - заданная точность.

Если условие не выполняется, то нужно продолжить применение метода по схеме 1 – 8.

В этом случае отрезок изоляции корня сужается и имеет вид . Приближённые значения корня находятся по формулам:

и .

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

В простейших случаях используют метод простых итераций, вычисляя последовательно значения функции у = f(x), изменяя значения х, чтобы у 0 и его модификацию – метод вилки, изменяя величину х так, чтобыесли f(x₁ ) 0, то f(x₂ ) 0и сужая интервал [х₁ ; х₂ ].

1.2.2. Интерполяция функций.

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

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

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

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

Линейная интерполяция – интерполяция алгебраическим двучленом Р₁ (x) = ax + b функции f(x), заданной в двух точках x₀ и x₁ отрезка [a, b].

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

Интерполяционная формула Ньютона применяется, если узлы интерполяции равноотстоящие и упорядочены по величине, так что x_{i + 1} x_i = h = const, то есть x_i = x₀ + ih. Тогда интерполяционный многочлен можно записать в форме Ньютона.

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

Короткая форма интерполяционной формулы Ньютона

В случае равноудаленных центров интерполяции, находящихся на единичном расстоянии друг от друга, справедлива формула:

где – обобщенные на область действительных чисел биномиальные коэффициенты.

Прямая интерполяционная формула Ньютона

где

а выражения вида ^k y_i – конечные разности.

Обратная интерполяционная формула Ньютона

где .

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

Для n + 1 пар чисел , где все x_i различны, существует единственный многочлен L(x) степени не более n, для которого L(x_i ) = y_i . В простейшем случае (n = 1) – это линейный многочлен, график которого – прямая, проходящая через две заданные точки.

Лагранж предложил способ вычисления таких многочленов:

где базисные полиномы определяются по формуле:

l_j (x) обладают следующими свойствами:

а) являются многочленами степени n; б) l_j (x_j ) = 1; в) l_j (x_i ) = 0 при .

Отсюда следует, что L(x), как линейная комбинация l_j (x), может иметь степень не больше n, и L(x_j ) = y_j .

1.2.3. Метод наименьших квадратов и его применения.

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

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

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

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

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

Но использовать в качестве критерия близости сумму отклонений не имеет смысла, т. к. при сложении разности будут компенсировать друг друга. Поэтому, учитывая также и то, что величина погрешности в экспериментальных точках может быть разной, необходимо минимизировать среднее значение суммы абсолютных погрешностей в заданных точках. Если приближаемая функция y = f(x) задана таблицей своих значений: y_j = f(x_j ), j = 1, 2, ..., n, и имеется некоторая приближающая функция Ф(х), определенная для всех значений x_j , то данный критерий запишется следующим образом:

Это условие было предложено Эджвортом. В современной литературе этот способ аппроксимации носит название равномерное приближение. Однако приближение функций по этому способу в широкое употребление не вошло.

Вместо среднего значения модуля отклонения используется среднее квадратическое отклонение эмпирической и теоретической величины в соответствии с выражением:

Если же приближаемая функция y = f (x) задана аналитически, т. е. она считается известной в любой точке x отрезка [a; b], то близость между y и приближающей функцией Ф(x) понимается в интегральном смысле:

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

Метод наименьших квадратов был предложен в начале XIX столетия К. Гауссом (1794–95) и независимо от него А. Лежандром (1805–06). Первоначально этот метод использовался для обработки результатов астрономических и геодезических наблюдений. Строгое математическое обоснование и установление границ содержательной применимости метода наименьших квадратов были даны А. А. Марковым и А. Н. Колмогоровым. Сейчас этот метод представляет собой один из важнейших разделов математической статистики и широко используется для статистических выводов в различных областях науки и техники.

Сущность обоснования метода наименьших квадратов (по Гауссу) заключается в допущении, что «убыток» от замены точного (неизвестного) значения физической величины её приближённым значением X , вычисленным по результатам наблюдений, пропорционален квадрату ошибки:

(X – )² , где – оцениваемая величина.

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

Рассмотрим частный случай зависимости наблюдаемой величины от искомых параметров :

.

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

Обозначим:

.

Теперь для определения искомых x, y, z будем иметь систему n уравнений:

.

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

.

Построим систему нормальных уравнений:

,

,

.

Выполнив дифференцирование, получим

,

,

.

Имеем систему, состоящую из трех линейных уравнений с тремя неизвестными, которую легко тем или иным способом решить.

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

, и т.д.

В этих обозначениях нормальные уравнения примут вид

,

,

.

Основное свойство нормальных уравнений - симметричность матрицы системы:

Решить эту систему можно используя, например, формулы Крамера, или один из матричных способов.

Покажем численный пример применения МНК для аппроксимации результатов эксперимента линейным уравнением первой степени.

В результате эксперимента получены семь значений искомой функции Y при семи значениях аргумента X. Используя метод наименьших квадратов, найти функциональную зависимость между Х и У в виде линейной функции у = ах + b. Построить график этой функции. Отметить экспериментальные значения.

Х	1	2	3	4	5	6	7
У	6,5	7,0	5,1	5,8	4,5	4,9	3,0

Решение : Построим корреляционное поле по данным семи наблюдений.

Получили точки определяющие линию регрессии У на Х. Характер изменения функции У позволяет предположить, что зависимость между У и Х близка к линейной вида у = ах + b.

Коэффициенты a и b уравнения найдём, выполнив необходимые вычисления:

Расчеты сведем в таблицу. В правом столбце этой таблицы записаны суммы по её строкам.

Номер измерения	1	2	3	4	5	6	7	Сумма
Х	1	2	3	4	5	6	7	28
У	6,5	7,0	5,1	5,8	4,5	4,9	3,0	36,8
х2	1	4	9	16	25	36	49	140
хy	6,5	14	15,3	23,2	22,5	29,4	21	131,9
	6,90	6,35	5,80	5,26	4,71	4,16	3,62	36,80
	1,54	3,04	0,02	0,29	0,57	0,13	5,09	10,70
	0,16	0,42	0,50	0,29	0,04	0,54	0,38	2,34

Нормальные уравнения для этого случая будут иметь вид:

Из этих уравнений имеем:

Т.о. линейная зависимость У от Х имеет вид: .

По этому уравнению вычислены теоретические значения величины .

Прямая линия зависимости У от Х построена вместе с ломаной линией по двум точкам:

При х = 1. у = 7,443 – 0,5464·1 = 6,897.

При х = 7. у = 7,443 – 0,5464·7 = 3,618.

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

Это означает, что зависимость переменной у от переменной х определяется полученным уравнением регрессии на 78,1%, а вариация величины у на оставшиеся 21,9% объясняется другими факторами.

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

Численное интегрирование – вычисление значения определённого интеграла (как правило, приближённое), основанное на том, что величина интеграла численно равна площади криволинейной трапеции, ограниченной осью абсцисс, графиком интегрируемой функции и отрезками прямых x = a и x = b, где a и b – пределы интегрирования.

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

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

где

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

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

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

x_i = a + ih; где h = (b a) / n – шаг сетки; n – число узлов сетки, а индекс узлов . Слагаемое – погрешность метода, которая может быть найдена разными способами. Для нечетных погрешность может быть найдена интегрированием погрешности интерполяционного полинома подынтегральной функции.

Частными случаями формул Котеса являются: формулы прямоугольников (n = 0), формулы трапеций (n = 1), формула Симпсона (n = 2), формула Ньютона (n = 3) и т. д.

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

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

Формула б) выражает площадь ступенчатой фигуры, состоящей из «выходящих» прямоугольников, также называемая формулой правых прямоугольников.

Чем меньше длина отрезков, на которые делится отрезок , тем точнее значение искомого интеграла, вычисляемое по этим формулам.

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

Учитывая большую точность последней формулы при том же объеме и характере вычислений её называют формулой прямоугольников

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

Площадь трапеции на каждом отрезке:

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

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

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

.

Если разбить интервал интегрирования на 2N равных частей, то имеем

Описанные выше методы используют фиксированные точки отрезка (концы и середину) и имеют низкий порядок точности (1 – методы правых и левых прямоугольников, 2 – методы средних прямоугольников и трапеций, 3 – метод парабол (Симпсона)). Если можно выбирать точки, в которых вычисляется значения функции f(x), то можно при том же количестве вычислений подынтегральной функции получить методы более высокого порядка точности. Этот способ реализован в методе Гаусса .

Так для двух (как в методе трапеций) вычислений значений подынтегральной функции, можно получить метод уже не 1-го, а 3-го порядка точности:

.

В общем случае, используя n точек, можно получить метод с порядком точности 2n 1. Значения узлов метода Гаусса по n точкам являются корнями полинома Лежандра степени n.

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

Метод Гаусса-Конрада позволяет оценить точность значения интеграла, вычисленного методом Гаусса.

В численном интегрировании имеются труды Чебышева, интегрирование при бесконечных пределах рассмотрены в методе Самокиша, существуют методы Монте-Карло , применяемые в многомерных случаях, методы Рунге-Кутта .

2. Конечно-разностный метод решения краевых задач.

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

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

с граничными условиями, заданными на концах отрезка [a; b]:

Следует найти такое решение у(х) на этом отрезке, которое принимает на концах отрезка значения у₀ , у₁ . Если функция линейна по аргументам , то задача поиска этой функции – линейная краевая задача, в противном случае – нелинейная..

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

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

где – такие числа, что

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

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

- метод стрельбы (пристрелки);

- конечно-разностный метод.

Используя конечно-разностный метод, рассмотрим двухточечную краевую задачу для линейного дифференциального уравнения второго порядка на отрезке [а; b].

Введем разностную сетку на отрезке [а; b]:

Решение задачи будем искать в виде сеточной функции:

предполагая, что решение существует и единственно.

Введем разностную аппроксимацию производных следующим образом:

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

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

Для этой системы уравнений при достаточно малых шагах сетки h и q(x_k ) 0выполнены условия преобладания диагональных элементов:

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

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

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

Пример .Решить краевую задачу:

с шагом 0,2.

Здесь р(х) = х; q(x) = 1; f(x) = 0; N = 5; x₀ = 0; x₁ = 0,2; x₂ = 0,4; x₃ = 0,6; x₄ = 0,8; x₅ = 1;

Во всех внутренних узлах отрезка [0; 1]после замены производных их разностными аналогами получим:

На левой границе y₀ = 1, на правой границе аппроксимируем производную односторонней разностью 1-го порядка:

С помощью группировки слагаемых, приведения подобных членов и подстановки значений x_k , а также с учётом у₀ = 1,получим систему линейных алгебраических уравнений:

.

В результате решения системы методом Крамера в Excel, получим:

Решением краевой задачи является табличная функция:

k	0	1	2	3	4	5
xk	0	0,2	0,4	0,6	0,8	1,0
yk	1,0	0,772	0,583	0,431	0,313	0,223

3. Расчетная часть

3.1 . Найти действительные корни уравнения методами простых итераций и касательных (Ньютона) с точностью до 0,00001.

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

Построим в осях ХОУ графики функций:

:

Линии графиков пересекаются в единственной точке с абсциссой х₀ , лежащей в интервале [0,5; 0,6], т.е.

а = 0,5; b = 0,6.

Значение функции на концах интервала:

Т.к. знаки различны, то уравнение имеет единственный корень в интервале [0,5; 0,6].

3.1.1 . Уточнение корня методом простых итераций.

Приведём исходное уравнение к виду:

(x) = x + С·f(x).

Т.к. первая производная заданной функции в этом интервале положительна и численно первая производная на этом участке близка к 1,5, то константу С выбираем из интервала:

Примем С = – 1.

Т.о. итерационная функция приобретает вид:

(x) = x – f(x).

Делаем первую итерацию:

Делаем вторую итерацию:

Делаем третью итерацию:

Делаем четвёртую итерацию:

Делаем пятую итерацию:

Делаем шестую итерацию:

Делаем седьмую итерацию:

Делаем восьмую итерацию:

Делаем девятую итерацию:

Продолжая далее, получаем:

На 19-ой итерации изменение шестого знака после запятой, позволяет утверждать, что пятый знак – после запятой – 5. Т.о. значение корня с заданной точностью:

х₀ = 0,57615

3.1.2 . Уточнение корня методом касательных (метод Ньютона):

Т.к. уравнение то же, то интервал, содержащий искомый корень, оставляем тот же [0,5; 0,6], т.е. а = 0,5; b = 0,6.

Находим первую и вторую производную функции :

Очевидно необходимые условия выполняются, т.к.:

, т.е. сохраняют знак на отрезке .

Выполняем первое приближение (х₀ = 0,5):

Выполняем второе приближение (х₁ = 0,571429):

Выполняем третье приближение (х₂ = 0,576128:

Выполняем четвёртое приближение (х₃ = 0,576146):

В пределах заданной точности f(x₂ ) оказался равен нулю, т.е. требуемая точность достигнута за 4 шага. Значение корня с заданной точностью:

3.2 . Вычислить приближенное значение интеграла , используя формулы:

а) трапеций (n = 10); б) Симпсона (n = 10); в) Гаусса (n = 5).

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

Разобьём интервал (–1; 9) на n = 10 отрезков (h =1) и вычислим значения подынтегрального выражения для начала и конца каждого отрезка.

№	0	1	2	3	4	5	6	7	8	9	10
х	0	1	2	3	4	5	6	7	8	9	10
	2,4495	2,6458	3,7417	5,7446	8,3666	11,4455	14,8997	18,6815	22,7596	27,1109	31,7175

Тогда по формуле трапеций, имеем:

Используя формулу Симпсона (формулу параболических трапеций) в виде:

получим:

Применяя к исходному интегралу квадратурную формулу Гаусса, имеем:

где

Для n = 5, коэффициенты t_i , представляющие нули полинома Лежандра и коэффициента А_i (эти значения табулированы в справочных таблицах) составляют:

i	1	2	3	4	5
ti	–0,9061	–0,5385	0	0,5385	0,9061
A1	0,2369	0,4786	0,5689	0,4786	0,2369
хi	0,4695	2,3075	5	7,6925	9,5305
	2,4705	4,2763	11,4455	21,4756	29,5239

Тогда:

3.3. Построить интерполяционные многочлены Лагранжа и Ньютона по следующим табличным данным:

	2,9	4,4	6,3	9,7
	2,84	4,53	6,04	5,50

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

Решение : Интерполяционный полином Лагранжа для четырёх узлов интерполяции записывается в виде:

Подставим численные значения из заданной таблицы:

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

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

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

Тогда интерполяционный полином Ньютона L_n (x) приобретает следующую форму:

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

По заданным узлам интерполяции х_i значения полинома по этому уравнению составляют:

х	2,9	4,4	6,3	9,7
Ln (x)	2,840133	4,530614	6,041651	5,504897
f(x)	2,84	4,53	6,04	5,50

Расчётные значения практически совпадают с заданными значениями f(x).

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

3.4 . Найти оценки параметров линейной и квадратичной моделей функциональной зависимости величин у и х по результатам наблюдений , приведенным в таблице:

	0,4	2,4	3,4	4,4	5,4
	2,14	2,14	2,24	2,34	2,34

Построить чертеж: на плоскости нанести экспериментальные точки , построить графики полученных эмпирических функций .

Решение : Коэффициенты a₀ и а₁ линейной модели найдём, выполнив необходимые вычисления. Расчеты сведем в таблицу:

Номер наблюдения	1	2	3	4	5	Сумма
х	0,4	2,4	3,4	4,4	5,4	16
у	2,14	2,14	2,24	2,34	2,34	11,2
х2	0,16	5,76	11,56	19,36	29,16	66
хy	0,856	5,136	7,616	10,296	12,636	36,54
	2,108	2,202	2,249	2,297	2,344	11,200
	0,0011	0,0039	0,0001	0,0019	0,0000	0,0069

Тогда:

Т.о. линейная зависимость у = а₀ + а₁ х имеет вид: у = 2,08865 + 0,0473х.

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

Коэффициенты а₀ , а₁ , а₂ квадратичной зависимости найдём, также выполнив необходимые расчёты в таблице:

Номер наблюдения	1	2	3	4	5	S
х	0,4	2,4	3,4	4,4	5,4	16
у	2,14	2,14	2,24	2,34	2,34	11,2
х2	0,16	5,76	11,56	19,36	29,16	66
х3	0,064	13,824	39,304	85,184	157,464	295,84
х4	0,0256	33,1776	133,634	374,81	850,306	1391,95
у·х	0,856	5,136	7,616	10,296	12,636	36,54
у·х2	0,3424	12,3264	25,8944	45,3024	68,2344	152,1
	2,128	2,182	2,230	2,292	2,368	11,200
	0,0001	0,0018	0,0001	0,0023	0,0008	0,0051

Составим систему уравнений:

Решение этой системы методом Крамера даёт:

Т.о. квадратичная зависимость у = а₀ + а₁ х + а₂ х² имеет вид:

у = 2,12433 + 0,00729·х + 0,006996·х² .

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

Эмпирическая ломаная, а также линии линейной и квадратичной модели построены на рисунке.

Результаты и выводы.

1. Т.о. интерполяционный полином Лагранжа и Ньютона, построенный по 4 заданным узлам интерполяции имеет вид:

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

Полученное уравнение позволяет найти приближённые значения функции в любых промежуточных точках от х₁ = 2,9 до х₄ = 9,7.

2. Применение метода минимальных квадратов (МНК) к аппроксимации пяти экспериментальных точек линейной зависимостью вида у = а₀ + а₁ х, т.е. прямой линией и квадратичной зависимостью вида , т.е. параболой дало следующие выражения:

– линейная зависимость реализована уравнением: у = 2,0887 + 0,0473х

– квадратичная зависимость реализована уравнением: у = 2,1243 + 0,0073·х + 0,007·х² .

Судя по остаточной сумме квадратов отклонений, квадратичная зависимость несколько лучше аппроксимирует экспериментальные данные, т.к. для неё остаточная сумма квадратов отклонений меньше, чем для линейной функции.

Список использованной литературы

1. Самарский А.А. Гулин А.В. Численные методы. М. МГУ. 1989 год.

2. Н. С. Бахвалов; Н.П. Жидков; Г.М. Кобельков. Численные методы. М 2003 год;

3. В.А. Буслов, С.Л.Яковлев. Численные методы иисследование функций. СПГУ. Курс лекций. СПБ 2001 г

4. Г.А. Зуева. Метод наименьших квадратов и его применение. Электронное учебное пособие. Иваново, 2009

Скачать архив с текстом документа