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

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

объекта, его параметры и внутренние и внешние связи. Для построения математической модели необходимо следующее 1 тщательно проанализировать реальный объект, процесс или систему 2 выделить наиболее существенные черты и свойства 3 определить переменные параметры, значения которых влияют на основные черты и свойства объекта, процесса или системы 4 описать зависимость основных свойств объекта, процесса или системы от значений переменных с помощью математических соотношений 5

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

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

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

S Обратная задача кинематики S Необходимость отыскания корней нелинейных уравнений встречается при расчетах систем управления и регулирования, собственных колебаний машин и механизмов, задачах кинематического анализа и синтеза плоских и пространственных механизмов и т.д. Дано fx1 Необходимо решить это уравнение. Если функция fx имеет вид многочлена , то уравнение 1 называется алгебраическим и имеет m корней. Пример . Если fx включает в себя тригонометрические или экспоненциальные

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

сужение интервала изоляции до величины, равной . Процесс отделения корней проводят, как правило, исходя из физического смысла прикладной задачи графически, с помощью таблиц, при помощи специальной программы. Считают, что интервал a b при решении прикладных задач задан. Для алгебраических и трансцендентных уравнений пригодны одни и те же методы уточнения корней 1 метод половинного деления 2 метод итераций 3 метод Ньютона касательных 4 метод хорд 5 метод секущих и др.

Метод половинного деления Дано уравнение fx0. Необходимо решить его с заданной точностью , если известен интервал изоляции a b. Словесный алгоритм. делим интервал a b пополам. координата середины отрезка в качестве нового интервала изоляции принимаем ту половину интервала a b, на концах которой функция имеет разные знаки. Для этого a вычисляем значение функции в точках a и t b проверяем если faft 0, то корень находится в этой половине интервала. Тогда b t. Если условие не выполняется, то корень находится на

интервале t b и a t. В обоих случаях получаем новый интервал a b, длина которого в 2 раза меньше предыдущего интервала. итерационный процесс повторяем, начиная с первого пункта до тех пор, пока длина a b не станет равной или меньше , то есть пока не выполнится условие . Блок-схема алгоритма. Метод итераций Для применения данного метода исходное уравнение fx0 должно быть преобразовано к виду x x. Пример В качестве начального приближения x0 принимаем любую точку интервала

a b. Далее итерационный процесс поиска строится по схеме x1 x0 x1 x0 x2 x1 x3 x2 В результате схождение к корню реализуется рекуррентной формулой xn1 xn Итерационный процесс поиска прекращается, как только выполняется условие . Для использования этого метода то есть чтобы последовательность x1, x2, сходилась к корню необходимо выполнение условия . Блок-схема алгоритма. Метод Ньютона метод касательных

Метод основан на замене функции fx на каждом шаге итерационного процесса поиска касательной, пересечение которой с осью абсцисс дает приближение корня. В качестве точки начального приближения примем точку b. Проведем в этой точке касательную на первом шаге поиска при x x0. Получим первое приближение корня x1 x0 h, где . В результате итерационный процесс схождения к корню реализуется рекуррентной формулой 1. Итерационный процесс продолжается, пока не выполнится условие 2.

Подставив 1 в 2 получаем . Метод обеспечивает быструю сходимость, если . Первую касательную рекомендуется проводить в той точке интервала a b, где знаки функции и ее кривизны совпадают. Блок-схема алгоритма. Модифицированный метод Ньютона. Этот метод заключается в том, что вместо вычисления производной на каждом шаге поиска находится ее приближенное значение с помощью конечных разностей .

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

В первом случае итерационный процесс схождения к корню реализуется формулой . Во втором случае . Процесс поиска продолжается, пока не выполнится условие xn1 xn . Метод обеспечивает быструю сходимость, если fx0 fx0 0, поэтому хорды фиксируются в том конце интервала a b, где знаки функции и ее кривизны совпадают. Блок-схема алгоритма Моделирование линейных систем со многими степенями подвижности.

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

В общем виде система линейных уравнений имеет вид где aij коэффициенты системы, bi свободные члены, xj неизвестные, i номер строки, j номер столбца, n порядок. В матричной форме СЛУ имеет вид или Численные методы решения СЛУ можно разбить на две группы 1 точные прямые методы 2 приближенные итерационные методы. Методы первой группы позволяют получить решение за конечное число арифметических действий.

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

треугольному виду. На втором этапе обратный ход из системы треугольного вида последовательно в обратном порядке находятся корни. Для простоты рассмотрим систему 4-го порядка. 1 Прямой ход 1-й шаг, k 1. Находим из 1-го уравнения системы 1 x1. ведущий элемент первой строки. Если 0, то Обозначив , получим , 2 где и вычисляются по рекуррентным формулам Подставим выражение 2 для x1 во 2-е, 3-е и 4-е уравнения системы 1.

В результате получим 4 вычисляются по рекуррентным формулам 5 В результате выполнения всех шагов прямого хода система 1 приводится к системе треугольного вида, которая получается объединением всех выражений для x1, x2, xn. Наша система принимает вид 6 где Вычисления по всем шагам прямого хода при построении алгоритма организуем в цикле по . Последний, n-й шаг прямого хода выводим из цикла, так как на этом шаге не нужны вычисления

aij и bi, а реализуется только одно выражение . На каждом шаге прямого хода все преобразования коэффициентов и свободных членов проводятся по рекуррентным формулам 3 и 5. Обратный ход. В процессе обратного хода из системы 6 неизвестные находятся в обратном порядке цикл по . Используем рекуррентную формулу bi bi xjaij xi bi. Поиск ненулевого ведущего элемента. Рассмотренный простейший метод

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

Если в строке kn этого столбца такой элемент найден, то строки k и kn меняются местами. Для этого используем рабочую ячейку R. Блок-схема алгоритма. Блок 2 ввод данных С помощью двух вложенных циклов с управляющими переменными i и j организуем ввод коэффициентов aij и свободных членов bi исходной системы. Для того, чтобы в дальнейшем можно было выполнить в блоке 9 проверку результата, в алгоритме предусматривается

сохранение значений aij и bi исходной системы с помощью переприсвоений cij aij di bi. Блок 3. Организуем в цикле по k, где k 1, n 1, вычисления по всем шагам прямого хода. Последний n-й шаг прямого хода выводим из цикла. Блок 4. Поиск ненулевого ведущего элемента. Поиск ненулевого ведущего элемента ведется в следующем порядке а на каждом шаге прямого хода, в k-й строке, элемент ak, k сравнивается с нулем. б если в k-й строке имеется

нулевой ведущий элемент ak, k, то в k-м столбце, начиная со следующей строки, в цикле по k1 k 1, n осуществляется поиск ненулевого ведущего элемента. в если в строке kn такой ненулевой элемент akn, k найден, то строки k и kn поэлементно, в цикле по j 1, n, меняем местами. Для перестановки элементов строк используется рабочая ячейка R. г если ненулевой ведущий элемент не найден, то коду ошибок

KO присваивается единица и расчет прекращается. Блок 5 шаг прямого хода. Проводим прямой ход исключения неизвестных путем преобразования коэффициентов и свободных членов системы по полученным ранее рекуррентным формулам k 1, n 1 i k 1, n j k 1, n Блок 6. В этом блоке выполняется последний шаг прямого хода. Блок 7. В процессе обратного хода метода Гаусса из системы треугольного вида последовательно в обратном

порядке в цикле i n 1, 1 находим неизвестные системы по рекуррентной формуле. Блок 8. Корни системы можно вывести в виде массива x, либо поэлементно в цикле. Блок 9 . Этот блок в алгоритме рекомендуется использовать в процессе отладки метода в этом случае можно объединить блоки 8 и 9. В дальнейшем при использовании метода Гаусса, при решении других прикладных задач блок 9 можно опустить.

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

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

и Ньютона. Дана СНУ n-го порядка Необходимо решить эту систему, т.е. найти с точностью . В качестве начального приближения принимаем вектор начальных решений . Выбор начального приближения вектора начальных решений во многом определяет эффективность итерационных методов. Если начальное приближение лежит за пределами области сходимости системы, то решение получить нельзя. Метод простых итераций. Для применения этого метода исходная система 1 должна быть преобразована

к виду 2 Далее, выбрав вектор начальных решений, итерационный процесс строим по схеме , k шаг. Итерационный процесс прекращается как только выполняется условие Условие сходимости итерационного процесса , т.е. сумма частных производных каждого преобразованного уравнения по всем переменным по модулю должна быть меньше 1. Блок-схема алгоритма Метод Ньютона Метод Ньютона наиболее распространенный метод для решения систем

нелинейных уравнений. Он обеспечивает более быструю сходимость, чем другие методы. В основе метода Ньютона лежит представление всех n уравнений исходной системы 1 в виде рядов Тейлора R1, , Rn члены более высоких порядков. Если приращение hi таковы, что все функции fi принимают значения, близкие к корню, то будем считать, что все левые части этих уравнений обращаются в нули. Отбросив все Ri, получим систему линейных уравнений, в которой все приращения hi являются неизвестными.

3 Запишем эту систему в матричной форме 4 Матрица частных производных , i 1, n j 1, n называется матрицей Якоби или якобианом. В методе Ньютона система 4 решается на каждом шаге итерационного процесса поиска. Найденные значения приращения hi на каждом шаге используются как поправки к предыдущему решению. Итерационный процесс прекращается как только выполнится условие i 1, n 5. Если хотя бы по одному приращению условие не будет выполнено, то поиск продолжается, пока не выполнится

условие 5 по всем приращениям или пока не станет ясно, что процесс расходится и получить решение не удастся. Если значения корней hi значительно отличаются друг от друга, то следует пользоваться нормированными приращениями i 1, n 6 В методе Ньютона на каждом шаге поиска необходимо формировать матрицу Якоби, при этом частные производные вычисляются аналитически, либо все частные производные можно заменить конечно-разностными значениями либо . Величина области сходимости метода

Ньютона обратно пропорциональна степени сложности и числу уравнений. С увеличением числа неизвестных область сходимости уменьшается и в случае больших систем сходимость обеспечивается если начальное приближение X0 очень близко к исходному решению. Блок-схема алгоритма Обработка опытных данных с помощью ЭВМ Любому специалисту в своей практической и теоретической деятельности приходится изучать зависимости

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

получена зависимость, но не в виде функции, а виде таблицы или говорят в результате эксперимента получена табличная функция. Вслед за этим возникают различные задачи по обработке опытных данных. Интерполирование функции Дана табличная функция yi fxi, i 0, n, т.е. дана таблица, в которой для некоторых дискретных значений аргумента xi, расположенных в порядке возрастания заданы соответствующие значения функции. ixy0x0y01x1y12x2y2nxnynТочки с координатами xi, yi называются узловыми точками или узлами.

Количество узловых точек N равно n 1. На графике табличная функция представляет собой совокупность узлов. Задача интерполирования функции интерполяции состоит в том, чтобы найти значение функции, заданной таблично, в любой промежуточной точке с аргументом xk, причем xk не является узлом табличной функции и лежит в промежутке xk x0 xn. Задача экстраполирования функции экстраполяции состоит в том, чтобы найти значение функции, заданной таблично, в точках xl, лежащих вне отрезка x0 xn.

Эта задача называется задачей прогноза. Обе эти задачи решаются с помощью выражения некоторой функции Px, которая приближала бы заданную табличную функцию т.е. в узловых точках xi принимала бы значения yi. Для определенности задачи искомую функцию Px будем выбирать из класса алгебраических многочленов. Степень многочлена n зависит от числа узловых точек и равна n N 1. Общий вид многочлена . Назовем этот многочлен интерполяционным.

Интерполирование с помощью алгебраических многочленов называется параболическим интерполированием. Построение интерполяционного многочлена Для этого построения необходимо найти его коэффициенты ai. Для нахождения коэффициентов необходимо построить систему линейных уравнений, которая может быть получена на основании того, что многочлен проходит через все узловые точки. Pnxi yi. В результате имеем систему Порядок системы равен n x0, x1, , xn, y1, y2, ,yn известны и заданы

в табличной функции. Неизвестными системы являются коэффициенты a0, a1, , an. Далее полученный многочлен используется для решения задачи интерполирования. Интерполяционный многочлен может быть построен и другими методами, например, при помощи специальных интерполяционных формул Лагранжа, Ньютона, Бесселя и др. Построение интерполяционного многочлена по формуле

Лагранжа Интерполяционный многочлен по формуле Лагранжа имеет вид Многочлен LN является интерполяционным многочленом т.е. в узловых точках он принимает значения таблицы. Свернем формулу Лагранжа , где 1 Можно избавиться в формуле от нулевых индексов, т.к. не все языки программирования допускают построение цикла с 0 , где . D значение аргумента в промежуточной точке, для которой необходимо решить задачу интерполирования,

N количество узловых точек. Алгоритм метода Лагранжа не предусматривает получение многочлена в явном виде, а сразу находит значение в промежуточных точках. Блок-схема алгоритма Построение интерполяционного многочлена по методу Ньютона x1, x2, xn узлы. y1, y2, yn значения функции или разделенные разности нулевого порядка. разделенные разности 1-го порядка разделенные разности 2-го порядка.

Разделенная разность k-го порядка для k 1-й точки вычисляется через разделенную разность k 1-го порядка . Лемма. Пусть x1, x2, xn произвольные попарно несовпадающие узлы, в которых известны значения функции y1, y2, yn. Тогда алгебраический многочлен n-1-го порядка, который записывается в виде является интерполяционным многочленом. Например x1, x2, x3, x4, x5, y1, y2, y3, y4, y5. xixixi, xi1xi, xi1, xi2xi, xi1, xi2, xi3xi, xi1, xi2, xi3, xi4x1y1x2y2yx1, x2x3y3yx2, x3yx1, x2, x3x4y4yx3, x4yx2, x3, x4yx1, x2, x3, x4x5y5yx4,

x5yx3, x4, x5yx2, x3, x4, x5yx1, x2, x3, x4, x5Интерполяционный многочлен Ньютона используется для неравных промежутков. Блок-схема алгоритма Аппроксимация опытных данных Проведен эксперимент, в результате которого получена табличная функция. ixy0x0y01x1y12x2y2nxnynxi, yi координаты узловых точек. На графике функция представляется совокупностью узловых точек, их количество равно n 1.

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

Pnx a0xna1xn-1 an. В этом случае аппроксимирующая кривая проходит через все узловые точки. Однако этот метод имеет недостатки 1 точность аппроксимации гарантируется лишь в небольшом интервале x0 xn при количестве точек не более 7 8. 2 значение функции в узлах должно быть задано с большой точностью. Однако известно, что как бы точно не проводился эксперимент, в результатах неизбежны погрешности. Разброс результатов неизбежен в результате эксперимента.

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

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

МНК. На практике чаще всего применяется другой метод аппроксимации опытных данных МНК. Сущность метода состоит в том, что опытные данные аппроксимируют кривой Fx, которая необязательно должна проходить через все узлы, а должна сгладить все случайные помехи табличной функции. При этом аппроксимирующую кривую Fx стремятся провести так, чтобы все ее отклонения от табличной функции уклонения были минимальными. 0 Fx0 - y0 1 Fx1 - y1 i

Fxi - yi Избавимся от знака уклонения возведем его в квадрат. i2 Fxi - yi2 МНК. Его суть. Для табличных данных, полученных в результате эксперимента, отыскать аналитическую зависимость кривой Fx, сумма квадратов уклонений которой по всем узлам была бы минимальной. Интерполяционный многочлен Для определенности задачи будем выбирать функцию Fx из класса алгебраических многочленов. Степень многочлена m.

Аппроксимирующий многочлен не проходит через все узлы, поэтому степень многочлена не зависит от количества узлов. Как правило, в задаче степень аппроксимирующего многочлена задается, при этом всегда m n. Если m1, то задача называется линейная регрессия. Если m2, то задача называется квадратичная аппроксимация. Если m3, то задача называется кубическая аппроксимация.

Уточним задачу аппроксимации. Для табличных данных, полученных в результате эксперимента fxi yi, i 0, n, построить аппроксимирующий многочлен , причем m n, для которого выполняется условие 1 Изменим вид многочлена Перепишем условие 1 Необходимым условием существования минимума функции S является равенство нулю ее частных производных Дифференцируя S по каждой переменной, получаем систему линейных уравнений 2

Порядок системы равен m 1. k 0, 2m j 0, m aj j 0, m коэффициенты аппроксимирующего многочлена неизвестные. Решая полученную систему, определяем значения коэффициентов аппроксимирующего многочлена и в результате сможем построить аппроксимирующий многочлен степени m. Для решения полученной системы 2 используем метод Гаусса. Перепишем систему 2 3 , k 1, m1, j 1, m1 , k 1, m1 aj коэффициенты многочлена вида 4

Задачу решаем в следующем порядке 1. Задаемся степенью m. 2. Строим систему 3. Т.к. матрица коэффициентов симметрична относительно главной диагонали, то достаточно определить только наддиагональные элементы. 3. Решаем систему 3 методом Гаусса и находим неизвестные коэффициенты аппроксимирующего многочлена. 4. Строим аппроксимирующий многочлен и определяем его значение в каждом узле.

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

P am1 P aj xiP, j m, 1, -1 Блок-схема алгоритма Блок 3 Блок 4 Блок 6 Примечание. Вывод результата блок 7 может быть выполнен в блоке 8 или отдельно. В последнем случае для значений многочлена уклонений необходимо формировать массивы. Блок 8. Решение задач интерполяции или прогноза для заданных точек выполняется по схеме Горнера в случае необходимости. Выбор оптимального варианта аппроксимирующего многочлена

Выбор оптимального варианта аппроксимирующего многочлена, т.е. выбор его оптимальной степени m, зависит от того, что мы положим в основу выбора критерия оптимальности. Можно положить в эту основу минимальные значения D Dmin. В этом случае поиск оптимальной степени m будет вестись в пределах mнач mкон. Схема алгоритма Окончательно ответить на вопрос, что должно быть положено в основу выбора критерия

оптимальности многочлена, может только практика. В нашем случае это решение задач интерполяции или прогноза и сравнение результатов расчетов с результатами экспериментов. Численные методы интегрирования Дана функция y fx. Необходимо найти интеграл на отрезке a b. Если подынтегральная функция fx задана аналитически, если эта функция непрерывна на a b, если известна первообразная

Fx, то На практике, при решении прикладных задач численные методы интегрирования применяются в следующих случаях 1 Функция задана таблично yi fxi, i 0, 1, n. 2 Функция задана аналитически, но ее первообразная Fx не выражается через элементарные функции. 3 Функция задана аналитически, имеет первообразную, но ее определение слишком сложно. В численных методах интегрирования не используется нахождение первообразной.

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

функция fx аппроксимируется интерполяционным многочленом 2 на практике чтобы не иметь дело с полиномами высоких степеней весь участок a b делится на части 3 строят интерполяционный многочлен, но не для всего отрезка a b, а для отдельных его частей 4 вычисляются приближенно площади частичных трапеций 5 в результате приближенное значение интеграла равно сумме площадей частичных трапеций. Нахождение приближенного значения интеграла называется квадратурой, а формулы для приближенного вычисления

интеграла называются квадратурными формулами или квадратурными суммами. Величина называется остаточным членом или погрешностью квадратурной формулы. Из всех численных методов интегрирования наибольшее распространение получили методы 1 прямоугольников 2 трапеций 3 Симпсона. Если в каждой части деления отрезка a b подынтегральная функция аппроксимируется многочленом нулевой степени прямой, параллельной оси абсцисс, то квадратурная формула называется формулой

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

Симпсона. Метод прямоугольников Рассмотрим два случая 1 функция fx задана аналитически 2 функция fx задана таблично. Разделим отрезок a b на n равных частей где n количество частей деления участка. N n 1 На каждом участке через узловые точки проводим многочлен нулевой степени. В результате вся кривая подынтегральной функции fx на отрезке a b заменяется ломаной линией, отрезки которой в каждой части деления параллельны x, а вся площадь под кривой fx заменяется суммой площадей

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

вид Точность этого метода имеет порядок h, т.е. чем меньше h, тем точнее результат. Блок-схема алгоритма I вариант. Подынтегральная функция задана аналитически, алгоритм с автоматическим выбором шага. C 0 метод левых прямоугольников C 1 метод правых прямоугольников S1 предыдущее значение интеграла S текущее значение интеграла Внутри цикла S сумма игреков II вариант. Функция задана таблично в неравноотстоящих точках.

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

Точность этого метода имеет порядок h2. Блок-схема алгоритма. 1 Функция задана аналитически алгоритм с автоматическим выбором шага. 2 Функция задана таблично в неравноотстоящих точках. Квадратурная формула имеет вид . Метод Симпсона В этом методе вся подынтегральная функция fx на отрезке a b заменяется кусочно-непрерывной линией, состоящей из отрезков квадратичных парабол.

В результате значение интеграла равно сумме всех площадей под квадратичными параболами. Т.к. для построения квадратичной параболы необходимо учесть три узла, то каждая часть деления в методе Симпсона равна 2h. Тогда количество частей деления N2 N 2. N всегда четное число. Определим S1 под квадратичной параболой на отрезке x0, x2 1 Неизвестные коэффициенты определим из условия прохождения параболы через 3 узла.

Составим систему линейных уравнений 2 Решая систему 2, найдем коэффициенты a0, a1, a2 и подставим их в уравнение 1. В результате площадь под квадратичной параболой на отрезке x0, x2 равна . На участке xi 1, xi 1 площадь равна , где номер части деления длиной 2h k 1, N2. Суммируя площади под квадратичными параболами на всех участках деления длиной 2h, получим квадратичную формулу по методу Симпсона . Точность метода имеет порядок h3 h4.

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

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

поэтому большое распространение получили численные методы решения дифференциальных уравнений. Численные методы решения дифференциальных уравнений Общий вид дифференциального уравнения Fx, y, y 0, где y yx неизвестная функция от x. Нормальная форма дифференциального уравнения y fx, y. Дифференциальное уравнение вида y fx, y, в котором неизвестная функция y зависит от одного аргумента

x называется обыкновенным дифференциальным уравнением, если от нескольких уравнением в частных производных. Общим решением дифференциального уравнения y fx, y является семейство функций y yx, C. При решении прикладных задач обычно ищут частное решение. Выделение частного решения из семейства общих осуществляется с помощью начальных условий. Нахождение частного решения дифференциального уравнения, удовлетворяющего условию yx0 y0, называется

задачей Коши. Постановка задачи Коши. Дано y fx, y yx0 y0 a b, h Решить дифференциальное уравнение. В численных методах задача Коши ставится следующим образом найти табличную функцию yi fxi, которая удовлетворяет заданным начальным условиям на участке a b с шагом h. ixy01x1y12x2y2nxnynНа графике решение задачи Коши численными методами представляется в виде совокупности узловых точек.

Метод Рунге-Кутта Наиболее эффективными и часто применяемыми методами решения задачи Коши являются методы Рунге-Кутта. Они основаны на аппроксимации искомой функции y в пределах каждого шага h многочленом, который получен при помощи разложения функции y в ряд Тейлора. Разложим функцию yx в окрестности шага h каждой i-й точки 1 Усекая ряд Тейлора в различных точках Рунге и Кутт обобщили, усовершенствовали и развили методы решения

дифференциальных уравнений. Метод Рунге-Кутта первого порядка метод Эйлера Если шаг h мал, то члены ряда, содержащие h2, h3, h4 и т.д. являются малыми более высоких порядков. Отбросим в 1 члены ряда, содержащие h2, h3 и т.д тогда получим формула Эйлера. Точность метода Эйлера имеет порядок h2. Блок-схема алгоритма Геометрический смысл метода Эйлера yi1 yi hfxi, yi y fxi, yi tg i tg i тангенс угла наклона касательной

к искомой функции y в точке xi yi. yi1 yi htg i. tg 0 fx0, y0 tg 1 fx1, y1 В методе Эйлера на графике искомая функция yx заменяется ломаной линией, каждый отрезок которой на шаге h линейно аппроксимирует искомую функцию. Поэтому метод Эйлера получил название метода ломаных. Основной недостаток метода Эйлера систематическое накопление ошибок. Модифицированный метод

Эйлера метод Рунге-Кутта второго порядка Отбросим в 1 члены ряда, содержащие h3, h4 и т.д. Получим ряд 2 Чтобы сохранить член ряда, содержащий h2, нужно определить член ряда, содержащий вторую производную. Будем аппроксимировать вторую производную разделенной разностью 2-го порядка. 3 Подставляя 3 в 2, получим модифицированная уточненная формула Эйлера. Для определения предварительного значения yi1 воспользуемся формулой

Эйлера, тогда модифицированная формула Эйлера примет вид , где . Точность метода имеет порядок h3. Обозначим Блок-схема алгоритма Геометрический смысл модифицированного метода Эйлера Т.к. fxi, yi yxi tg i fxi1, yi1 yxi1 tg i1, то модифицированная формула Эйлера примет вид . P1 ошибка по методу Эйлера P2 ошибка по уточненному методу

В модифицированном методе Эйлера модифицированная функция yx на каждом шаге аппроксимируется не одной, а двумя прямыми, каждая из которых аппроксимирует функцию yx на полушаге. Метод Рунге-Кутта 4-го порядка В методе Эйлера используется член ряда Тейлора 1, содержащий h. Ошибка метода имеет порядок h2 . В модифицированном методе Эйлера используется член ряда 1, содержащий h2.

Ошибка этого метода составляет h3. Математики Рунге и Кутт обобщили и развили эти методы. Согласно их классификации, метод Эйлера получил название метода Рунге-Кутта 1-го порядка. Рунге и Кутт предложили методы для решения задачи Коши, в которых более высокая точность решения дифференциальных уравнений может быть достигнута, если

в методах сохранять члены ряда Тейлора, содержащие h3, h4 и т.д. Наибольшее распространение получил метод, в котором искомая функция yx аппроксимируется рядом Тейлора 1, в котором сохранен член ряда, содержащий h4. Метод получил название метода Рунге-Кутта 4-го порядка, но в литературе он известен больше как просто метод Рунге-Кутта. Точность этого метода имеет порядок h5

Для получения значения функции yi1 по методу Рунге-Кутта необходимо определить y и y. Определяем эти производные используя разделенные разности. Окончательно для определения yi1 выполняется следующая последовательность вычислительных операций Блок-схема алгоритма Решение дифференциальных уравнений высоких порядков Методы Рунге-Кутта можно использовать не только для решения дифференциальных уравнений 1-го порядка

y fx, y, но и для решения дифференциальных уравнений более высоких порядков ym fx, y, y, y, , ym 1 1 Любое дифференциальное уравнение m-го порядка можно свести к системе, состоящей из m дифференциальных уравнений 1-го порядка. Для этого заменим y1 y y2 y y1 y3 y y2 ym ym ym 1. Мы получили систему дифференциальных уравнений 1-го порядка 2 Получили систему 2 дифференциальных уравнений 1-го порядка.

Порядок системы равен m. Неизвестными системы являются m функций y, y1, y2, ym. В прикладных задачах чаще всего решаются дифференциальные уравнения второго порядка. Дифференциальное уравнение второго порядка y fx, y, y сведем к системе, состоящей из двух уравнений первого порядка Пример 3 4 Решением полученной системы являются функции y и y1. Задача Коши для системы дифференциальных уравнений второго порядка должна содержать два начальных условия.

В прикладных задачах начальные условия имеют четкий физический смысл. Сформулируем задачу Коши для системы дифференциальных уравнений 2-го порядка. Необходимо решить систему дифференциальных уравнений на участке a b с шагом h. В численных методах задача Коши для системы дифференциальных уравнений 2-го порядка ставится следующим образом найти табличные функции yxi и y1xi, i 1, n. ixyy0x0y0y10x0 a1x1y1y11xn b2x2y2y12n b x0 hnxnyny1nНа

графике табличные функции y и y1 представлены совокупностью узловых точек. Для решения системы дифференциальных уравнений используются те же методы, что и для решения одного дифференциального уравнения 1-го порядка, но необходимо соблюдать условие на каждом шаге интегрирования, т.е. в точках x1, x2, xn все уравнения системы нужно решать параллельно. Рассмотрим систему 4. Необходимо решить ее с условиями 1 1,5 h 0,01 x0 1 y0 0,4401 y0 0,3251.

Блок-схема алгоритма Решение дифференциального уравнения m-го порядка методом Рунге-Кутта Любое дифференциальное уравнение m-го порядка 1 приводится к системе 2. Численным решением системы 2 являются m табличных функций. Каждую из m табличных функций определяют на участке a b с шагом h, она имеет n узловых точек. Таким образом, численным решением дифференциального уравнения m-го порядка является матрица порядка

n m 5 m порядок дифференциального уравнения n b a h i 1, n j 1, m Каждая i-я строка матрицы 5 это совокупность решений m табличных функций на одном i-м шаге интегрирования. Каждый столбец матрицы 5 это совокупность решений одной j-й табличной функции по всем n шагам интегрирования. При построении алгоритма задачи будем, как и ранее, расчет вести в цикле по i 1, n. На каждом i-м шаге интегрирования будем определять решения дифференциальных уравнений и тут же их печатать.

Тогда нет необходимости формировать матрицу 5, а можно ограничиться формированием массива решений на одном шаге интегрирования. Y y1, y2, ym 6 Массив 6 соответствует одной строке матрицы 5. При этом в массиве 6 используется соотношение y1 y y2 y y3 y и т.д. С учетом этих обозначений система 4 запишется как Тогда подпрограмма PRAV будет иметь вид F массив длиной m, куда на каждом шаге интегрирования нужно

заносить значения правых частей уравнений системы. PRAVm, x, Y, F Пример Программу решения системы дифференциальных уравнений m-го порядка методом Рунге-Кутта реализуем в виде трех программных модулей 1 Основная программа, в которой реализуется вычислительный процесс по всем шагам интегрирования 2 Подпрограмма RGK, в которой реализуется метод Рунге-

Кутта 4-го порядка для каждого i-го шага. И на этом шаге вычисляется массив Y длиной m 3 Подпрограмма PRAV, обращение к которой осуществляется из подпрограммы RGK для вычисления массива F. Блок-схема алгоритма Подпрограмма RGK В подпрограмме RGK на каждом шаге интегрирования вычисляется массив Y длиной m. При этом формулы для вычисления массива

Y методом Рунге-Кутта имеют тот же вид, что и при вычислении решений дифференциальных уравнений 1-го порядка, но с учетом поправки на массив. T рабочая матрица порядка 4 m, куда заносятся значения T1, T2, T3, T4 для всех элементов массива Y. Y1 рабочий массив длиной m. Организуется для сохранения решений дифференциального уравнения предыдущего шага интегрирования. Блок-схема алгоритма Оптимизация в задачах проектирования

Проектирование технических объектов, систем и технологических процессов ТОСТП можно представить себе как выбор оптимальных вариантов этих ТОСТП. Поиск наилучших вариантов требует создания и исследования большого числа дорогостоящих опытных образцов. Применение ЭВМ в проектировании сокращает сроки, стоимость проектирования, улучшает качество изделия и значительно сокращает сроки его внедрения.

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

автоматического проектирования САПР. Математическая модель в задачах оптимального проектирования. Для использования ЗВМ задача проектирования должна быть переведена на формальный математический язык, т.е. для реальных ТОСТП должна быть построена математическая модель. Построение математической модели включает в себя следующие задачи 1. Анализ проектируемых ТОСТП и выделения вектора переменных , где xi - набор конструктивных параметров,

значения которых определяют вариант конструкции ТОСТП. Меняя значения переменных можно менять свойства и качество ТОСТП. Пример. Необходимо спроектировать гидроцилиндр, который должен удовлетворять каким-то техническим требованиям. Значение n определяет мерность задачи. 2. Определяется критерий качества, позволяющий указать лучший из вариантов проектируемых

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

целевой функции.Допустим, если речь идет о цене 3. Формирование ограничивающих условий, исходя из требований, предъявленных к ТОСТП и формализация этих ограничений в виде системы равенств и неравенств . Необходимо определить цилиндр, обладающий заданным тяговым усилием и минимальными габаритами. В этом случае тяговое усилие-функция ограничена, а функция, определяющая габариты целевая функция.

Вектор может быть представлен точкой в n-мерном евклидовом пространстве En. Тогда ограничения hi 0 и gj 0, накладываемые на возможные изменения вектора , зададут область G в n-мерном пространстве En. Задача оптимального проектирования может быть сформулирована следующим образом найти вектор G En, доставляющий минимум или максимум целевой функции, т.е. arg min F при ограничениях в виде равенств hi 0 i 1, n или неравенств gj 0 j m 1, p.

Пример. Спроектировать контейнер объемом 1 м3. Желательно, чтобы при изготовлении контейнера было затрачено как можно меньше материала, чтобы контейнер удобно было брать погрузчиком его ширина должна быть не менее 1,5 м. Построим модель. 1 x1, x2, x3 2 Для того, чтобы использовать меньше металла, должна быть минимальна площадь боковой поверхности. S 2x1x2 x2x3 x1x3 3 x1x2x3 1 x1 1,5 Ограничение-равенство позволяет уменьшить размерность задачи.

Целевая функция здесь нелинейная, поэтому сформулированная задача является задачей нелинейного программирования, задача является двумерной, однокритериальной. Геометрическая интерпретация задач нелинейного программирования Рассмотрим интерпретацию для двумерной задачи. F x12 x2 g1 x12 x22 9 0 g2 x1 x2 1 0 g1 x12 x22 9 g2 x1 x2 1 Первое ограничение отсекает плоскость круга. Второе ограничение в плоскости круга отбрасывает одну полуплоскость.

Построим линии уровней множества точек с координатами x1, x2, при которых целевая функция имеет постоянное значение целевой функции. x12 x2 C x2 C x12 Линии уровней квадратичные параболы, симметричные относительно оси x2. C 0 x2 x12 C 1 x2 1 x12 C 3 x2 3 x12 0 3 Методы нелинейного программирования Построение математической модели первый этап в задаче оптимального проектирования. Построив математическую модель, необходимо выбрать или разработать достаточно эффективный метод поиска

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

программирования делятся на два основных класса 1 задачи линейного программирования 2 задачи нелинейного программирования. Если целевая функция и функция ограничения линейны, то соответствующая задача поиска экстремума является задачей линейного программирования. В противном случае нелинейного. Математические модели в задачах проектирования реальных ТОСТП должны отражать реальные физические процессы, которые в них протекают, а они, как правило, нелинейны.

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

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

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

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

Необходимо найти вектор , доставляющий минимум или максимум заданной функции на отрезке a, b с заданной точностью . Словесный алгоритм. 1 Делим отрезок a, b пополам x a b2 2 Вычисляем значение функции в окрестности точки X. F1 Fx F2 Fx 3 Сравниваем F1 и F2 и в результате отбрасываем одну из половин отрезка a, b . В случае, если ищем минимум, то если F1 F2 , то отбрасываем x, b, b x, если

F1 F2, то отбрасываем a, x, a x. 4 деление отрезка a, b продолжаем, пока длина отрезка a, b не станет меньше . Блок-схема алгоритма Метод золотого сечения. Этот метод также является методом прямого поиска. В методе золотого сечения на каждом шаге итерационного процесса поиска отрезок a, b делится на две неравных части по правилу золотого сечения отношение длины всего отрезка к длине большей его части равно отношению

длины большей части к длине меньшей части отрезка. При этом золотое сечение отрезка a, b производит две симметрично расположенные точки x1 и x2. Если золотое сечение отрезка a, b производит точка x1, то 1, если золотое сечение отрезка a, b производит точка x2, то 2. При этом точка x1 производит золотое сечение отрезка a, x2 3, а точка x2 производит золотое сечение отрезка x1, b 4. Из выражения 1 имеем

Найдем координаты точек x1 и x2. Из выражения 2 имеем Из выражения 3 находим координаты точки x1 имеем Словесный алгоритм метода поиск минимума. 1 Делим отрезок a, b точками x1 и x2 по правилу золотого сечения 2 F1Fx1 F2Fx2 3 Сравниваем F1 и F2 и отбрасываем либо интервал a, x1, либо интервал x2, b. Если F1 F2, то отбрасываем x2, b и делаем переприсвоение b x2, иначе отбрасываем отрезок a, x1 и делаем

переприсвоение a x1. В результате интервал a b уменьшился в раз. На оставшемся интервале уже есть одна точка, производящая золотое сечение, поэтому на каждом шаге итерационного процесса поиска, начиная со второго шага, достаточно вычисления лишь одной точки, производящей золотое сечение. 4 Поиск продолжается до тех пор, пока b a . Метод гарантирует нахождение минимума максимума в самых неблагоприятных условиях, но обладает медленной

сходимостью. Блок-схема алгоритма