Содержание к диссертации
Введение
ГЛАВА 1. Обзор методов решения обыкновенных дифференциальных уравнений и вариационных задач
1.1. Аналитические методы решения 10
1.2. Численные методы решения 19
1.3. Алгоритмы генетического программирования 25
Выводы 31
ГЛАВА 2. Алгоритм генетического программирования и методы его применения для решения задачи Коши в символьном виде для обыкновенных дифференциальных уравнений
2.1. Постановка задачи и предлагаемые подходы к решению 33
2.2. Алгоритм генетического программирования для решения задачи 37
2.3. Полученные результаты 43
Выводы 70
ГЛАВА 3. Алгоритм генетического программирования и методы его применения для решения в символьном виде задач вариационного исчисления
3.1. Постановка задачи и предлагаемые подходы к решению 73
3.2. Алгоритм генетического программирования для решения краевых задач и вариационных задач
3.3. Результаты численных экспериментов 81
Выводы 98
ГЛАВА 4. Программная реализация алгоритмов генетического программирования для решения задач Коши, краевых задач, вариационных задач
4.1. Программная система для решения задачи Коши 101
4.2. Программная система для решения вариационной задачи 109
4.3. Реализация программ для работы алгоритма ГП на распределенных и кластерных вычислительных системах
4.4. Реализация алгоритма генетического программирования для 121
решения в символьном виде задачи Коши для систем обыкновенных
дифференциальных уравнений
Выводы
Заключение
Список литературы
Список публжаций автора
- Численные методы решения
- Полученные результаты
- Алгоритм генетического программирования для решения краевых задач и вариационных задач
- Реализация программ для работы алгоритма ГП на распределенных и кластерных вычислительных системах
Численные методы решения
Обыкновенные дифференциальным уравнением (ОДУ) [2,21,49] /7-го порядка называется соотношение вида F(x,y ,у ,...,у(п)) = 0 (1.1.1) между независимым переменным х, его функцией у и производными О) у ,у ,-,г Функция у = р(х) называется решением этого дифференциального уравнения, если после замены у на (р{х), у на (р\х), ..., у на qrn\x) оно (уравнение (1.1.1)) обращается в тождество [27]. Замечание. Будем считать, что рассматриваемые величины принимают только действительные значения, а рассматриваемые функции однозначны. Таким образом, в ОДУ неизвестная функция зависит только от одного аргумента, в отличие от уравнений с частными производными.
Основной задачей при решении ОДУ является нахождение функций, являющихся решениями заданного ОДУ. Конечная последовательность элементарных действий над известными функциями и интегрирований этих функций (с целью получить решение) называется сведением к квадратурам. Процесс нахождения решения называется интегрированием дифференциального уравнения [30].
Рассмотрим некоторые распространенные ОДУ и методы их решения. Одними из наиболее часто встречающихся уравнений являются уравнения с разделяющимися переменными. ю
Дифференциальными уравнениями с разделяющимися переменными называются уравнения вида [44] = fl(x)f2(y). (1.1.2) ах Для таких уравнений доказана теорема существования и единственности решения, гласящая о том, что через каждую точку (xQ,yo) прямоугольника Q a x b, c y d проходит одно и только одно решение этого уравнения, при выполнений следующих условий: функции f\{x) и /2(.у) непрерывны, причем /2(.у) нигде (в Q) не обращается в ноль.
Для получения решения уравнения, учитывая обозначение у = р{х), нужно преобразовать уравнение (1.1.2) к виду /2 ( ?( )) и проинтегрировать обе части последнего равенства по х. Пусть F2(y) - какая-нибудь первообразная от , и F\(x) - какая fiiy) нибудь первообразная от /\(х). Учитывая монотонность, так как f2{y) нигде не обращается в ноль, решение уравнения (1.1.2) можно записать в виде: ф) = F l ([F2 (у0 ) + F{ (х) - Fl ( о )], где F - обратная функция для функции F2 (у).
Замечание. Если не выполняется условие, о том что /2(у) нигде не обращается в ноль, единственность решения может нарушаться и такой случай требует дополнительного исследования. її Нужно также отметить, что функция у = р(х) может быть не представима элементарными функциями. В таком случае решение уравнения представляется в виде, содержащем неопределенный интеграл, а на практике используются значения, полученные численным решением. Однородными дифференциальными уравнениями называют уравнения вида [42] :г=/()- о-1-3) ах у х Если функциями) определена при а и Ь, то функция /(—) будет определена У X в углах, образованных такими точками (х,у), для которых а — Ь, эту У область обозначим через G.
Для однородных уравнений справедлива теорема существования и единственности решения: через каждую точку(xQ,yo) из G проходит одно и только одно решение уравнения (1.1.3), если выполняются следующие требования: функция f(u) при а и Ь и f(u) u всюду непрерывна на интервале (а, Ь). Для решения необходимо сделать замену у=их, тогда уравнение (1.1.3) перепишется в виде хи л-и- f(u), откуда получается уравнение с разделяющимися переменными ay _ f(u) - и dx х которое интегрируется описанным выше способом. При этом получается следующее решение: 1пЫ = Ф(-) + С, У где Ф(и) - некоторая первообразная от /(и) - и Как и в случае с однородными дифференциальными уравнениями, решение многих других уравнений сводится к решению простейших уравнений с разделяющимися переменными. (1.1.4) Линейные дифференциальные уравнения 1-го порядка - уравнения вида: dy dx = а(х)у + b(x). Пусть функции а(х) и Ь(х) непрерывны в интервале а х Ь. Тогда через каждую точку (х$,уо) полосы, определенной неравенствами а х Ъ и - х _у оо, проходит одно решение. Для получения решения нужно рассмотреть сначала линейное однородное уравнение: dy ( л — = а(х)у, ах которое является уравнением с разделяющимися переменными. Последнее уравнение имеет единственное решение (для точки (xQ,yo)): у(х) = у0 ехр (х \a{t)dt Затем решается первоначальное уравнение (1.1.4), например, методом вариации постоянной. Результатом решения (с учетом решения однородного уравнения) является следующая функция [1]: у(х) = у0 ехр \a(t)dt JC (X \ + J&(s)exp \a(t)dt ds, XQ \S J которая является единственным решением уравнения (1.1.4). із К линейным уравнениям сводится уравнение Бернулли: & = а(х)у + Ь(х)уп, (л 1)[30]. ах Уравнение Риккати —— = а(х)у + Ь{х)уп + с(х) - также можно свести к ах уравнению Бернулли, но лишь в том случае, если известно какое-либо частное решение уравнения (в общем случае уравнение Риккати не сводится к квадратурам) [8]. Линейными уравнениями с постоянными коэффициентами называются уравнения вида: а0УМ+а1у(п-Г +... + апу = Ях). (1.1.5)
Решается такое уравнение в два этапа: сначала решается однородное уравнение (/(л;)=0), затем - неоднородное уравнение. Для решения однородного уравнения составляется характеристическое уравнение, корни которого определяют вид общего решения. Решение ищется в виде: у = С\у\ + С2У2 + + СпУп- Если коэффициенты уравнения (1.1.5) вещественные, то решение можно написать в вещественной форме, даже если корни характеристического уравнения комплексные. Решение однородного уравнения (для каждой пары комплексных сопряженных корней Я = а ± fii) в таком случае можно записать следующим образом:
Полученные результаты
Рассмотрим следующую задачу Коши для обыкновенных дифференциальных уравнений в общем виде. Пусть дано уравнение в виде: (х,у ,у ,...,у(п))=0 (2.1.1) и начальные условия у(х0) = о ,у ( о) = о ,у\ о) = о,..., У(п 1) ( о) = У о 1 »(2Л-2) где у, у , у", ... , у(п) - неизвестная функция со своими производными, до п-го порядка включительно, х0 - начальная точка интервала, на котором определена функция, Уо Уо УО --- УО " заДанные вещественные числа, необходимые для определения единственного решения задачи Коши.
Требуется найти функцию, удовлетворяющую уравнению (2.1.1) и начальным условиям (2.1.2). Для поставленной задачи доказана следующая теорема:
Теорема существования и единственности решения для уравнения п-го порядка. Существует единственное решение дифференциального уравнения п-го порядка у = f(x, у, у ,..., у п ), удовлетворяющее условиям У( о)= У0 У (х0) = у 0,у"(х0) = у о,..., yi"- [\x0) = у%-\ если в окрестности начальных значений (xQ,yQ,yQ,...,yQ ) функция f является функцией всех своих аргументов и удовлетворяет условию Липшица по всем аргументам, начиная со второго. Условия Липшица по аргументам можно записать следующим образом: Данная теорема дает основания считать возможным получить единственное решение (при удовлетворении условий теоремы).
Замечание. Здесь и дальше будем полагать, что требования теоремы выполнены для искомых функций.
Для решения поставленной задачи предлагается два подхода. Возможны некоторые другие варианты решения задачи [65,66, 68, 85], но они не рассматриваются в рамках работы в силу того, что необходимо решить задачу в символьном виде (найти точное решение) при том, что структура результата (количество узлов, набор функций) априорно неизвестны.
Прямой метод решения
Для решения поставленной задачи предлагается метод решения, основанный на использовании алгоритма генетического программирования (алгоритм описан ниже). При этом в качестве входных данных используется уравнение (2.1.1) и начальные условия (2.1.2). Метод предполагает прямую подстановку решений в заданное уравнение, а также проверку начальных условий.
Точным решением задачи (при использовании прямого метода решения) будем считать функцию, удовлетворяющую уравнению и начальным условиям в некоторых точках заданного интервала. Количество и метод выбора таких точек зависит от длины интервала и характера задачи, если задача представляет реальный процесс, выбираемые точки должны качественно его описывать. Необходимые для подстановки в ОДУ производные в точках могут быть получены из символьной функции, которая подставляется в уравнение.
В частном случае, если ОДУ является уравнением нулевого порядка, прямой метод решения задачи Коши можно рассматривать как задачу символьной регрессии. Для существования и единственности решения поставленной задачи старшая производная должна быть выражена явно У -/(х,У у ,---,У ) предполагается, что это возможно, но для решения задачи прямым методом уравнение не обязательно приводить к специальному виду.
Гибридный метод решения Задачу отыскания решения задачи Коши в символьном виде можно решать, используя комбинацию численного метода и алгоритма генетического программирования. Для этого предварительно ОДУ решается численно, например, методом Рунге-Кутта с четвертым порядком точности. Такой подход является одним из вариантов получения символьных решений ОДУ [96-98]. Решение численным методом предполагает предварительное видоизменение уравнения в пригодную для решения форму, вследствие чего некоторые частные решения должны быть рассмотрены отдельно. Если уравнение (1) является уравнением первого порядка, преобразуем его, если это возможно, к виду пригодному для итерационной процедуры численного метода:
В случае наличия производных более высокого порядка, такое уравнение сводится к системе аналогичных дифференциальных уравнений 1 го порядка. Если в общем виде приведение уравнения к приемлемому для решения виду невозможно, применяют частные варианты решения таких уравнений. где (хп,уп) - аргумент и значение сеточной функции (решения задачи),/ функция из правой части уравнения (2.1.3), т - шаг сетки.
Результатом численного метода является числовая таблица, представляющая собой значения искомой функции в заданных точках. Эта таблица рассматривается как входные данные для алгоритма ГП, который в этом случае решает задачу символьной регрессии в обычном смысле, подбирая символьное выражение, наиболее точно описывающее числовые данные.
Такой подход позволяет ускорить процесс работы алгоритма ГП - нет необходимости вычислять производные в точках выборки. При этом теоретическая обоснованность (слабое место всех стохастических процедур) метода повышается - доказаны существование и единственность решения для результата, полученного методом Рунге-Кутта (а алгоритм ГП представляет собой процедуру поиска интерполирующей функции). Однако возникают недостатки, свойственные численным методам: не все уравнения можно решить таким образом (для некоторых уравнений нужно применять специальные схемы решений), в процессе приведения часть решений теряется — их приходится рассматривать отдельно, в некоторых задачах не выполняется условие непрерывной зависимости от входных данных. Кроме того, возникает недостаток, связанный с комбинацией численного метода с алгоритмом ГП - наложение вычислительных погрешностей (особенно для степенных и экспоненциальных функций).
Алгоритм генетического программирования для решения краевых задач и вариационных задач
На третьем этапе производится адаптация каждого индивида из популяции. Индивид, лучше адаптирующийся в условиях поставленной задачи, имеет более высокую пригодность, а, следовательно, более высокую вероятность быть отобранным для порождения потомков. С точки зрения теории ОДУ адаптация индивида означает изменение его частей так, чтобы решение имело меньшую погрешность, то есть означает оптимизацию решения - изменение вещественных коэффициентов и/или функционального набора бинарного дерева. Если зафиксировать все вершины дерева кроме какого-то одного узла, то его (дерево) можно рассматривать как функцию одной переменной.
Теория оптимизации содержит много методов поиска экстремума [3, 34, 37]. Методы отличаются главным образом количеством вычислений целевой функции и порядком старшей производной, которая необходима для вычислений (и связанным с ним видом функции, для которого метод наиболее эффективен). В силу того, что априорно ничего не известно о функции (решении ОДУ), предпочтение целесообразно отдавать методам нулевого порядка. Предлагается воспользоваться методом Хука-Дживса, отличающимся простотой и эффективностью. Бинарное дерево (решение ОДУ) может в общем случае содержать много вещественных коэффициентов и узлов-операций, это делает целесообразным при оптимизации проведение лишь небольшого числа итераций. Для ускорения поиска решения можно варьировать и функциональный набор, содержащийся в дереве.
Большое количество коэффициентов в дереве (и большой функциональный набор решения) может значительно затруднять оптимизацию решения задачи. В таком случае можно воспользоваться генетическими алгоритмами [13,25] или вероятностными генетическими алгоритмами [11,38-41]. Коэффициенты, а возможно, и функциональный набор из дерева, представляющего решение задачи, может быть закодирован в бинарную строку. Так, используя генетические алгоритмы, можно оптимизировать все коэффициенты решения задачи Коши для ОДУ. Хотя для оптимизации с помощью генетических алгоритмов необходимо гораздо больше ресурсов. Поэтому такой подход целесообразно применять лишь в тех случаях, когда традиционные методы оптимизации не обеспечивают поиск оптимальных коэффициентов в дереве решения задачи.
На четвертом этапе применяются генетические операторы: селекция, рекомбинация (скрещивание), мутация. Учитывая сложность задачи, генетические операторы должны предотвращать преждевременную сходимость и захват локальных минимумов. Поэтому необходима установка средних настроек. Селекция не должна осуществлять слишком жесткий отбор, иначе велик риск захвата определенного решения задачи (определенно структуры дерева). Коэффициент мутации должен обеспечивать изменение одного гена (узла дерева) в среднем. Слишком большой коэффициент скрещивания будет подвергать популяцию слишком частым изменениям, что не позволит оптимизировать коэффициенты при структуре близкой к искомой, слишком маленький - не позволит изменять структуру дерева. Поэтому настройки целесообразно выбирать средние, если нет специфичных требований к решению задачи.
Далее проверяется условие остановки алгоритма: выполнено заданное количество вычислений (прошло заданное количество времени) или достигнута заданная точность аппроксимации. Если условие выполняется, задача считается решенной, иначе - происходит переход ко второму этапу.
Отметим, что параметры, касающиеся алгоритма генетического программирования: коэффициенты скрещивания, мутации, максимальная глубина дерева, размер популяции и др. можно не изменять от задачи к задаче - эти параметры настраиваются, исходя из общей концепции (высокая надежность, быстрая сходимость и т.д.). Оставшиеся параметры: границы интервала, объем выборки, коэффициент штрафа за несоответствие начальным условиям, порядок аппроксимации производных и др., должны задаваться, исходя из конкретной задачи.
Полученные результаты Рассмотрим тестовую задачу, у"—2-sin(x), у(0)=0.0, у (0)=2.0. В качестве функционального множества будем использовать только элементарные арифметические операции (+,-, /), что существенно усложнит алгоритму ГП поиск решения. Результатом решения такой задачи является выражение: у(х)=1.999611801-х-0.3330926297-х3+0.01662624753-х5-0.0003886344121-х7+0.4403-10-5-х9.
Легко видеть, что это выражение является представлением точного решения (y=2sin(x)) с помощью формулы Тейлора с остаточным членом достаточно высокого порядка малости.
Были рассмотрены задачи Коши с различными типами уравнений: уравнения с разделяющимися переменными, однородные, линейные уравнения первого порядка (приводимые к ним) и другие [47]. Задачи были решены прямым методом и гибридным. Задачи запускались многократно, данные усреднялись. Ниже представлены некоторые из решенных задач,
Реализация программ для работы алгоритма ГП на распределенных и кластерных вычислительных системах
Алгоритм генетического программирования для решения поставленной задачи должен оперировать бинарными деревьями (в простейшем варианте), представляющими собой потенциальные решения поставленной задачи. Каждое бинарное дерево состоит из элементов функционального множества (+, -, ,/, sin, cos, exp, log и т.д.) и элементов терминального множества (х или компоненты вектора X, а также вещественные коэффициенты). В случае решения специфических задач эти множества могут быть дополнены.
В общем виде схему работы алгоритма можно представить следующим образом (схема в целом аналогична схеме из главы 2, хотя в силу различия задач есть некоторые отличия):
На первом шаге начальная популяция инициализируется случайным образом: до заданной начальной глубины случайно выбираются элементы функционального множества и записываются в бинарное дерево, там выбираются унарные операторы, дальнейшее строительство дерева продолжается по одной ветви; на заданной глубине в дерево записываются случайно выбранные элементы терминального множества. Количество индивидов в популяции также является параметром работы алгоритма и должно в некоторой степени зависеть от поставленной задачи, так как слишком большая популяция требует много ресурсов для обработки, слишком маленькая не дает точный и надежный результат.
На втором шаге происходит оценка каждого индивида в условиях поставленной задачи, при этом учитываются различные факторы: во-первых, соответствие функции граничным условиям, во-вторых, нарушение тождества при подстановке проверяемого решения задачи в уравнение (2), в-третьих, сложность полученной функции.
Соответствие функции граничным условиям определяет значительно более узкий круг функций, подходящих в качестве решения задачи. При этом возможны три варианта: 1 - через граничные точки не проходит ни одно решение, 2 - проходит единственное решение, 3 - проходит несколько решений. В первом случае вариационная задача не имеет решения при поставленных краевых условиях. Второй случай является наиболее желанным, когда решение существует и единственно. В третьем случае нужно найти несколько решений, этот случай будет отмечен ниже. Самым простым способом учитывать граничные условия является вычисление отклонения от них с некоторым коэффициентом. В более сложных случаях может быть введена нелинейная зависимость функции пригодности от разницы между функцией и точками на границе.
Примечание. В случае, когда решение не единственное и существует конечное число решений, можно последовательно находить решения поставленной задачи, попутно штрафуя (например, добавив при вычислении пригодности дополнительное слагаемое) те решения, которые сходятся к уже найденным (или известным) решениям.
Наличие тождества в уравнении (3.1.2) определяет функцию являющуюся возможным решением задачи (3.1.1). Если функционал (на определенной функции) достигает минимума, то уравнение Эйлера превращается в тождество. Проверив все решения уравнения (3.1.2) - если решение не единственное, найдем (локальный) минимум (максимум) функционала (глобального экстремума функционал может и не достигать).
Нужно учитывать и сложность полученной функции, она должна быть наименьшей - при идентичности графиков кривых, предпочтение отдается функции, которая имеет более короткую, лаконичную запись. Поэтому функцию пригодности можно представить следующим образом:
где Fit(Pk) - значение функции пригодности k-го индивида, Е(Рк) - ошибка аппроксимации, вычисляемая по всем точкам выборки (N - объем выборки), К1 - коэффициент штрафа за сложность дерева, D(Pk) - число вершин дерева Рь К2 - коэффициент штрафа за граничные условия, Y0,YN - заданные граничные условия, y(Xj) - значение функции в точке Xj, Pk( ,) - отклонение от нуля (например, среднеквадратическое) функции F из уравнения (2), получаемое при подстановке решения Рк в ОДУ в выбранных точках ( ,) -точки из интервала могут быть выбраны равномерно, случайно или каким-либо специальным образом. Здесь при вычислении ошибки соответствия решения ОДУ требуется вычисление производных в точках выборки. Возможны два способа вычисления: численный (с определенным порядком малости) и аналитический - построить дерево производной функции и вычислить значения этой функции в точках выборки. Так как дифференцирование (в отличие от интегрирования) проводится всегда по строгим правилам, аналитически вычислить производную не представляет особого труда. Нетрудно заметить, что получение численной оценки производной требует меньше машинной памяти и времени, однако аналитически вычисленная производная может быть точнее (при не слишком большом количестве операторов в функции), поэтому метод вычисления производных должен выбираться с учетом решаемой задачи.
Замечание. При решении задачи (3.2.1) прямым методом несколько изменится способ вычислением функции пригодности, поскольку пригодность оценивается по значению функционала. Если решается задача минимизации функционала, то функцию пригодности можно представить следующим образом: Fit(Pk) = l + E(Pk) + K\-D(Pk) + K2- І у(х,)-, 7=0 где E(Pk) - значение функционала на заданном интервале (остальные обозначения соответствуют формулам (3.2.1), (3.2.2)).
В случае, когда через граничные точки проходит не единственное решение, алгоритм может сойдется к какому-то одному решению. Чтобы найти остальные решения, в функцию пригодности (а именно, в знаменатель выражения (3.2.1)) может быть добавлен дополнительный компонент, например, такой:
