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

Размещено на http://www.allbest.ru/

Введение

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

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

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

уже не устойчиво ни для какого положительного запаздывания [5]. Другим примером проблемы неучтённого запаздывания может служить модель системы автоматического регулирования («идеальный предсказатель»), где значение входного сигнала в будущий момент времени полностью определяется значением выходного сигнала в настоящий момент времени , что противоречит как здравому смыслу, так и принципу причинности [5].

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

1) уравнение Маккея-Гласса, описывающее концентрацию белых кровяных клеток в организме человека:

,

где и - параметры;

2) уравнение кинетики ферментов

,

где и - параметры;

3) общая модель роста капитала Солоу

где ;

4) логистическое уравнение с запаздыванием (уравнение Хатчинсона или уравнение Райта)

,

где и - параметры.

Другие интересные примеры использования ОДУ с запаздывающими аргументами можно найти, например, в книге [5].

Модели, как правило, содержат набор параметров, которые их характеризуют. Эти параметры - неизвестны, определяются в каждом случае отдельно по некоторому массиву наблюдаемых значений (значений, полученных в ходе эксперимента). Для уравнения Маккея-Гласса (1.3), например, параметрами выступают переменные и , а наблюдаемыми значениями - величины концентрации лейкоцитов в моменты времени .

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

оптимизационной задачи

ввиду отсутствия в явном виде.

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

Цели дипломной работы:

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

2. Реализовать алгоритм в виде библиотеки на языке программирования MATLAB, а также программы с графическим интерфейсом пользователя

3. Апробировать полученное решение на некоторых реальных примерах.

Постановка задачи

Рассмотрим систему обыкновенных дифференциальных уравнений (ОДУ) с запаздывающим аргументом:

,

где - независимая переменная, обычно - время,- вектор параметров размерности : , - неизвестная вектор-функция независимого аргумента и параметра размерности , - запаздывания, - известная вектор-функция размерности .

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

,

где - начальная точка времени, - конечная.

Отметим, что некоторые элементы вектора-параметра могут входить только в функцию , некоторые только в .

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

,

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

Задача оценки параметров системы ОДУ с запаздывающим аргументом заключается в нахождении такого значения параметра , что решение системы (2.1)-(2.2) некоторым образом приближает данные .

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



Задача оценки параметров (2.1)-(2.2) при использовании МНК может быть записана следующим образом:

при условии, что

Отметим, что вид функции решения в общем случае неизвестен (см. главу «Введение») и точные значения заменяются их приближениями , полученные в результате численного решения ОДУ (2.6).

Задача оптимизации с ограничениями

В настоящей работе был разработан эффективный и устойчивый алгоритм для оценки неизвестных параметров ОДУ. Метод наименьших квадратов (МНК) использовался для определения функции цели.

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

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

Введём замену переменных (внесём в вектор неизвестных):

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

Запишем лагранжиан, связанный с задачей (2.8)-(2.9):

,

Гессиан лагранжиана имеет вид:

,

,

,

- нулевая квадратная матрица размерности .

Положим, что якобиан ограничений равен , т.е. что

Отметим, далее, что ограничения входят только в выражения якобиана (2.18) и гессиана ограничений (2.17), а, следовательно, только эти матрицы зависят от схемы численного интегрирования ОДУ.

Вычислим их для некоторых численных методов.

Метод Эйлера

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

Тогда ограничения для задачи (2.8)-(2.9) будут равны



Якобиан (2.20) примет вид:

,

где для :

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

и заметим, что в функцию входит .

Гессиан в методе Эйлера примет следующий вид:

где для :

Неявный метод Эйлера

Неявный метод Эйлера характеризуется следующей формулой:

Тогда ограничения задачи (2.8)-(2.9) будут иметь вид:



Якобиан (2.33) будет иметь такую же структуру что и якобиан (2.20), но элементы матрицы примут другие значения, а именно:

Значения гессиана будут следующими:

1. Описание алгоритма

Поиск решения

Задача (2.8)-(2.9) - типичная задача нелинейной оптимизации. Тем не менее, общие методы решения оптимизационных задач с ограничениями не могут работать достаточно эффективно с такого рода проблемой, потому что они не учитывают специальную структуру функций . Также, желательно использовать метод, который способен учесть и свойства целевой функции МНК (2.8), а, именно, специальный вид её матрицы Гессе.

Важно заметить, что задача (2.8)-(2.9) обычно имеет достаточно большую размерность, но только с несколькими степенями свободы, определяемых ОДУ с запаздывающими аргументами (2.1), нежели размерами входных данных. Для того чтобы сделать процедуру исключения переменных эффективной, вместо использования исключения по Гауссу в настоящей работе был рассмотрен метод, основанный на редукции переменных.

Также были рассмотрены другие методы решения системы, получаемой на каждом шаге SQP алгоритма, такие как:

1. Метод, основанный на блочной факторизации

2. Прямые и итерационные методы

Ниже (см. Рисунок 1) представлено дерево поиска решения - красными цветом отмечены отброшенные по тем или иным причинам варианты решения, зелёным - методы и способы решения, которые были приняты и применены.



Рисунок 1. Дерево поиск решения

Основываясь на результатах работы [1] для решения задачи (2.8)-(2.9) будет использовано последовательное квадратичное программирование (SQP), на каждом шаге которого функции цели и ограничений заменяются на их квадратичные приближения, и решается следующая подзадача: найти направление , такое, что является решением для квадратичной задачи:

В качестве матрицы (гессиан лагранжиана) может быть использован как полный гессиан:



так и неполный:

,

, а .

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

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

На каждом шаге работы SQP алгоритма решается система уравнений

,

где - гессиан лагранжиана, рассчитанный в точке при множителях Лагранжа (), - якобиан ограничений в точке, - множители Лагранжа на -м шаге алгоритма.

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

Для ОДУ без запаздывающих аргументов структура матрицы (3.5) будет иметь следующий вид (см. Рисунок 2):

Рисунок 2. Структура матрицы подзадачи SQP для ОДУ без запаздывающих аргументов (n=200)

Для ОДУ с запаздывающим аргументом в блоках матрицы (3.5) появятся дополнительные по диагонали, отвечающие запаздываниям (формулы расчёта элементов матриц приведены ниже в главах «Метод Эйлера» и «Неявный метод Эйлера») и структура матрицы примет следующий вид (см. Рисунок 3):



Рисунок 3. Структура матрицы подзадачи SQP для ОДУ с запаздывающими аргументами (n=200, три запаздывания)

Основными свойствами матрицы системы (3.5) являются:

1. Отсутствие положительной определённости

2. Симметричность

3. Разреженность

Далее (глава «Ускорение шага SQP») будут рассмотрены различные методы для решения такой системы линейных уравнений.

Ускорение алгоритма: редукция переменных

В качестве первого способа ускорения работы алгоритма был рассмотрен метод редукции переменных. Основная идея метода - уменьшение числа переменных и нахождение решения на пространстве меньшей размерности. Рассмотрим, например, метод редукции, приведенный в [2]: пусть имеется система разностных уравнений вида:

,

где - искомые векторы размерности , - заданные векторы и - заданная квадратная матрица порядка .

Запишем уравнение (3.6) в точках и :

Складывая (3.8) и (3.9) получим

,

откуда, учитывая, что:

,

придём к уравнению

,

связывающему значения искомого вектора в узлах одинаковой чётности. В частности, если - чётные, то проведено исключение нечётных узлов. Далее этот процесс исключения можно продолжить аналогичным образом. При этом предполагая, что число узлов является степенью двойки (), можно свести систему (3.6) к системе с двумя неизвестными. Случай для произвольного числа подробно рассмотрен в [2]. При проведение исключения переменных на -м этапе исключения получаем систему

,

где матрицы и векторы находятся из рекуррентных соотношений

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

Формулы редукции для ОДУ без запаздывающих аргументов

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

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

Обозначим матрицы преобразования

Запишем уравнения для и для системы (3.17)

,

откуда, выражая и , получаем

Аналогично, для получаем рекуррентную формулу

Введём замену переменных

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

,

Таким образом, получим окончательные формулы для матриц преобразования

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

Случай ОДУ с запаздывающими аргументами

В отличие от системы (3.17) система для ОДУ с запаздывающими аргументами будет иметь более сложную структуру, а именно:

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

Ускорение шага SQP

Так как применить метод редукции переменных не удалось, попытаемся на каждом шаге SQP алгоритма решать систему уравнений максимально быстро. Как уже упоминалось выше, матрица коэффициентов системы (3.5) не является положительно определённой, и достаточно быстрое разложение Холецкого не применимо в данном случае. Попробуем рассмотреть ряд других методов решения системы уравнений и выбрать из них лучший:

1. Асимметричное блочное разложение

2. Прямые и итерационные методы

Сначала, однако, зададимся вопросом о выборе типа гессиана лагранжиана для задачи (3.5).

Выбор типа гессиана

Рассмотрим вопрос выбора гессиана для подзадачи квадратичного программирования (3.1)-(3.2). Полная матрица Гессе (или матрица Гессе Ньютона) будет иметь вид

,

,

.

Неполный же гессиан (или гессиан Гаусса-Ньютона) получается отбрасыванием второго слагаемого в формуле полного гессиана (3.31):

.

При использовании неполного гессиана, говорят, что задача решается методом Гаусса-Ньютона.

Для наглядности структуры матриц Гессе и матриц коэффициентов систем (3.5) приведены на рисунках:



В теории, локальная скорость сходимости метода Ньютона квадратичная. С другой стороны скорость сходимости алгоритма Гаусса-Ньютона зависит от нелинейности задачи. Для линейной модели метод Гаусса-Ньютона имеет высокую скорость сходимости. Детальное сравнение методов можно найти, например, в [6].

Понятно, что от методов, использующих полный гессиан, ожидают более высокую скорость сходимости, чем от методов использующих неполный, из-за того, что первые используют более точное квадратичное приближение. Однако было показано (например, в [7]), что матрица Гаусса-Ньютона предпочтительнее матрицы Ньютона в приложениях, если данные сильно зашумлены. Дополнительным плюсом использования неполного гессиана является быстрота его построения, его структура и большая разреженность (см. Рисунок 2 и Рисунок 3), что приводит к более простому виду системы (3.5) (см. Рисунок 4 и Рисунок 5), а, следовательно, к более быстрому способу её решения.

Основываясь, на результатах работы [7] было решено замерить времена работы двух методов для нескольких тестовых задач и выбрать наиболее быстрый из них.

Рисунок 8. Время работы SQP для полного и неполного гессиана (тестовая задача 4)

Рисунок 9. Время работы SQP для полного и неполного гессиана (тестовая задача 5)

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

Асимметричное блочное разложение

В качестве первого способа решения системы (3.5) рассмотрим метод асимметричного блочного разложения описанный в [3].

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

и пусть матрица системы - положительно определена.

Тогда существует единственное разложение матрицы (разложение Холецкого) вида:

,

где и - множители Холецкого соответственно для матриц и , а .

Запишем два алгоритма нахождения матрицы в виде таблиц:

Таблица 1. Блочное разложение для положительно определённой матрицы

Номер шага	Симметричная схема	Асимметричная схема
1
2
3	как простое произведение матриц	Вместо перемножения сначала вычисляют , затем получают
4
5

В [3] показано, что если матрица имеет специальную разреженную структуру, то число операций для вычисления - меньше, а, следовательно, скорость разложения выше.

Запишем аналогичные формулы для случая, когда матрица не является положительно определённой:

Формула (3.37) - формула LDL разложения - используется вместо разложения Холецкого в случаях не положительно определённости матрицы .

Таблица 2. Блочное разложение для не положительно определённой матрицы

Номер шага	Симметричная схема	Асимметричная схема
1
2
3	как простое произведение матриц	Вместо перемножения сначала вычисляют , затем получают
4
5

Численные опыты показали, что данный подход не приводит к увеличению скорости решения системы (3.5). Такой результат обусловлен следующим фактором: начальная матрица - разреженная и разложение блока также приводит к разреженным матрицам и . Однако, вычисление произведения , как в симметричной схеме, так и в ассиметричной, приводит к заполненной матрице и её разложение занимает много времени.

Выбор метода решения

Последним вопросом, рассмотренным в данной работе, стал вопрос о выборе метода решения системы уравнений (3.5).

Как было показано выше в качестве левого верхнего блока матрицы (3.5) - - был выбран неполный гессиан лагранжиана задачи. Структура матрицы для такого выбора была представлена выше (см. Рисунок 4).

Выше уже был рассмотрен один специальный способов решения системы (3.5), а именно - ассиметричное блочное разложение, который, однако, не привёл к уменьшению времени решения. Так как других специальных методов для решения систем вида (3.5) неизвестно, было принято решения исследовать стандартные методы решения систем с разряженными матрицами и выбрать наиболее быстрые из них. Среди рассмотренных методов были:

1. Прямые методы:

a. LDL разложение (LDL)

b. LU разложение (LU)

2. Итерационные методы:

a. Метод бисопряжённых градиентов (BICG)

b. Квадратичный метод сопряжённых градиентов (CGS)

c. Симметричный LQ метод (SYMMLQ)

d. Обобщённый метод минимальный невязок (GMRES)

e. Метод квази-минимальных невязок (QMR)

f. Симметричный LQ метод (SYMMLQ)

Подробное описание вышеперечисленных методов можно найти, например, в [8], [9], [11].

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

Рисунок 10. Время работы SQP (тестовая задача 4)

Рисунок 11. Время работы SQP (тестовая задача 5)

Основываясь на численных результатах времени работы SQP алгоритма для различных методов решения системы шага SQP, было отобрано два наиболее быстрых метода -LDL разложение и метод бисопряжённых градиентов. Для этих двух методов далее проводилось дополнительное исследование по подбору их параметров для обеспечения минимального времени работы.

LDL разложение - это расширение разложения Холецкого для матриц не являющихся положительно определёнными. Формула разложения имеет следующий вид:

,

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

Как и в случае с методом Холецкого разложение может привести к увеличению числа ненулевых элементов в матрице относительно матрицы , и чтобы избежать этого применяют различные способы переупорядочивания элементов оригинальной матрицы, такие как, например, метод минимальной степени (MD), метод приближённой минимальной степени (AMD) и метод множественной минимальной степени (MMD). Подробное описание LDL метода и методов минимальной степени можно найти в [9] и [10] соответственно.

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

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



Рисунок 12. Время решения системы в зависимости от порогового значения для тестовой задачи 5

Результатом данного исследования стал выбор параметра .

Метод бисопряжённых градиентов

Вторым рассмотренным методом стал метод бисопряжённых градиентов. Этот метод, в отличие от метода сопряжённых градиентов, подходит для не положительно определённых систем линейных уравнений. В его основе лежит идея построения двух взаимно ортогональных последовательностей. Подробное описание алгоритма можно найти в [11].

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

Рисунок 15. Пороговое значение для построения предобуславливателя

Результатом данного исследования стал выбор порога 0.15.

2. Реализация алгоритма

дифференциальный уравнение программирование алгоритм

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

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

Рисунок 16. Скриншот части окна работающей программы

Рисунок 17. Скриншот части окна работающей программы

При написании программы были использованы некоторые встроенные в MATLAB и внешние библиотеки, среди, которых:

1. Optimization Toolbox для непосредственного решения задач оптимизации.

2. Symbolic Math Toolbox для символьного вычисления градиента и гессиана, введённой пользователем функции.

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

4. SqpLab для реализации SQP алгоритма и его модификации.

Общее число строк, содержащееся в программе - более 5000 без учёта внешних библиотек, среди которых около 3000 содержат ядро вычисление - реализация алгоритма решения, 500 - для работы с интерфейсом пользователя, а оставшиеся 1500 - тестовые примеры и unit-тесты.

3. Численные примеры

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

Тестовые примеры

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

В примеры для тестирования входят, как ОДУ с одним или несколькими неизвестными параметрами, так и ОДУ с запаздывающим аргументом.

Таблица 3. Решения для тестовых примеров

№	ОДУ	Начальные условия	Аналитическое решение ОДУ	- точное решение	- численное решение
1
2
3
4		Из аналитического решения

Графики решения для тестовых задач:

Рисунок 18. Решение тестовой задачи № 1

Рисунок 19. Решение тестовой задачи № 2



Рисунок 20. Решение тестовой задачи № 3

Рисунок 21. Решение тестовой задачи № 4

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

Таблица 4. Решения для тестовых примеров, не имеющих аналитического решения

№	ОДУ	Начальные условия	- точное решение	- численное решение
5		Значения, полученные решением ОДУ с точным значением параметра
6

Графики решения для тестовых задач, не имеющих аналитического решения:

Рисунок 22. Решение тестовой задачи № 5



Рисунок 23. Решение тестовой задачи № 6

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

Если в качестве примера взять точную сетку , то можно увидеть, как ведёт себя ошибка метода численного интегрирования (см. Рисунок 25). Известно, что для метода Эйлера ошибка дискретизация пропорциональна шагу , и, уменьшая шаг, можно добиться пропорционального уменьшения ошибки.



Рисунок 24. Разность решения и начальных данных для тестовой задачи 4. Начальные данные с ошибками

Рисунок 25. Разность решения и начальных данных для тестовой задачи 4. Начальные данные без ошибок

Задача демографической динамики

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

Данная модель характеризуется следующим дифференциальным уравнением:

,

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

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

Таблица 5. Оценка параметров демографической модели

Начальные параметры метода	Известные параметры модели	Оцениваемые параметры модели	Сумма квадратов отклонений
, , неполный гессиан, ограничения определяются неявным методом Эйлера, метод шага SQP: LDL,

Рисунок 26. Решение задачи для демографической модели (Pc=4)

Рисунок 27. Решение задачи для демографической модели (Pc=5.2)



Рисунок 28. Решение задачи для демографической модели (Pc=6)

Рисунок 29. Решение задачи для демографической модели (Pc=7)

Рисунок 30. Решение задачи для демографической модели (Pc=9)

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

Рисунок 31. Ошибка решения для задачи демографической модели

4. Безопасность жизнедеятельности

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

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

Цель и содержание БЖД:

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

2. Ослабление действия этих факторов до безопасных пределов или исключение их если это возможно.

3. Ликвидация последствий катастроф и стихийных бедствий.

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

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

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

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

· определение оптимальных условий труда инженера - программиста;

· расчет освещенности;

· расчет уровня шума.

Характеристика условий труда программиста

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

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

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

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

Требования к производственным помещениям

Окраска и коэффициенты отражения

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

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

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

· окна ориентированы на юг: - стены зеленовато-голубого или светло-голубого цвета; пол - зеленый;

· окна ориентированы на север: - стены светло-оранжевого или оранжево-желтого цвета; пол - красновато-оранжевый;

· окна ориентированы на восток: - стены желто-зеленого цвета; пол - зеленый или красновато-оранжевый;

· окна ориентированы на запад: - стены желто-зеленого или голубовато-зеленого цвета; пол зеленый или красновато-оранжевый.

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

· для потолка: 60…70%

· для стен: 40…50%

· для пола: около 30%

· для других поверхностей и рабочей мебели: 30…40%

Освещение

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

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

Существует три вида освещения - естественное, искусственное и совмещенное (естественное и искусственное вместе).

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

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

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

Согласно СанПиН 2.2.2/2.4.1340-03 в помещения вычислительных центров необходимо применять систему комбинированного освещения.

При выполнении работ категории высокой зрительной точности (наименьший размер объекта различения 0,3…0,5мм) величина коэффициента естественного освещения (КЕО) должна быть не ниже 1,5%, а при зрительной работе средней точности (наименьший размер объекта различения 0,5…1,0 мм) КЕО должен быть не ниже 1,0%.

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

Требования к освещенности в помещениях, где установлены компьютеры, следующие: при выполнении зрительных работ высокой точности общая освещенность должна составлять 300лк, а комбинированная - 750лк; аналогичные требования при выполнении работ средней точности - 200 и 300лк соответственно.

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

Параметры микроклимата

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

Вычислительная техника является источником существенных тепловыделений, что может привести к повышению температуры и снижению относительной влажности в помещении. В помещениях, где установлены компьютеры, должны соблюдаться определенные параметры микроклимата. В СанПиН 2.2.4.548-96 установлены величины параметров микроклимата, создающие комфортные условия. Эти нормы устанавливаются в зависимости от времени года, характера трудового процесса и характера производственного помещения (Таблица 6).

Объем помещений, в которых размещены работники вычислительных центров, не должен быть меньше 19,5м3/человека с учетом максимального числа одновременно работающих в смену. Нормы подачи свежего воздуха в помещения, где расположены компьютеры, приведены в табл. 7.

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

Период	Параметр микроклимата	Величина
Холодный	Температура воздуха в помещении	22…24?
	Относительная влажность	40…60%
	Скорость движения воздуха	до 0,1м/с
Теплый	Температура воздуха в помещении	23…25?
	Относительная влажность	40…60%
	Скорость движения воздуха	0,1…0,2м/с

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

Объём помещения на человека	Объемный расход подаваемого в помещение свежего воздуха, м3 /на одного человека в час
до 20м3	Не менее 30
20…40м3	Не менее 20
более 40м3	Естественная вентиляция

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

Шум и вибрация

Шум ухудшает условия труда, оказывая вредное действие на организм человека. Работающие в условиях длительного шумового воздействия испытывают раздражительность, головные боли, головокружение, снижение памяти, повышенную утомляемость, понижение аппетита, боли в ушах и т. д. Такие нарушения в работе ряда органов и систем организма человека могут вызвать негативные изменения в эмоциональном состоянии человека вплоть до стрессовых. Под воздействием шума снижается концентрация внимания, нарушаются физиологические функции, появляется усталость в связи с повышенными энергетическими затратами и нервно-психическим напряжением, ухудшается речевая коммутация. Все это снижает работоспособность человека и его производительность, качество и безопасность труда. Длительное воздействие интенсивного шума [выше 80 дБ(А)] на слух человека приводит к его частичной потере.

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

Таблица 8. Предельные уровни звука, дБ, на рабочих местах

Категория напряженности труда	Категория тяжести труда
	I. Легкая	II. Средняя	III. Тяжелая	IV. Очень тяжелая
I. Мало напряженный	80	80	75	75
II. Умеренно напряженный	70	70	65	65
III. Напряженный	60	60	-	-
IV. Очень напряженный	50	50	-	-

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

Электромагнитное и ионизирующее излучения

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

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

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

Таблица 9. Допустимые значения параметров неионизирующих электромагнитных излучений (в соответствии с СанПиН 2.2.2/2.4.1340-03)

Наименование параметра	Допустимые значения
Напряженность электрической составляющей электромагнитного поля на расстоянии 50см от поверхности видеомонитора	10В/м
Напряженность магнитной составляющей электромагнитного поля на расстоянии 50см от поверхности видеомонитора	0,3А/м
Напряженность электростатического поля не должна превышать:	для взрослых пользователей	20кВ/м
	для детей дошкольных учреждений и учащихся средних специальных и высших учебных заведений	15кВ/м

Для снижения воздействия этих видов излучения рекомендуется применять мониторы с пониженным уровнем излучения (MPR-II, TCO-92, TCO-99), устанавливать защитные экраны, а также соблюдать регламентированные режимы труда и отдыха.

Эргономические требования к рабочему месту

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

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

Эргономическими аспектами проектирования видеотерминальных рабочих мест, в частности, являются:

· высота рабочей поверхности

· размеры пространства для ног

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

· характеристики рабочего кресла

· требования к поверхности рабочего стола

· регулируемость элементов рабочего места

Главными элементами рабочего места программиста являются стол и кресло. Основным рабочим положением является положение сидя.

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

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

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

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

Рисунок 32. Зоны досягаемости рук в горизонтальной плоскости

Зоны, представленные на рис. 32:

а) зона максимальной досягаемости

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

в) зона легкой досягаемости ладони;

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

д) оптимальное пространство для тонкой ручной работы

Оптимальное размещение предметов труда и документации в зонах досягаемости:

· дисплей размещается в зоне а (в центре);

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

· клавиатура - в зоне г/д;

· компьютерная мышь - в зоне в справа;

· сканер - в зоне а/б (слева);

· принтер - в зоне а (справа);

· документация:

o необходимая при работе - в зоне легкой досягаемости ладони - в

o в выдвижных ящиках стола - литература, неиспользуемая постоянно.

Рисунок 33. Размещение основных и периферийных составляющих ПК

На рис. 33 показан пример размещения основных и периферийных составляющих ПК на рабочем столе программиста:

1) сканер

2) монитор

3) принтер

4) поверхность рабочего стола

5) клавиатура

6) компьютерная мышь.

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

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

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

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

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

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

Большое значение придается характеристикам рабочего кресла. Так, рекомендуемая высота сиденья над уровнем пола находится в пределах 420-550мм. Поверхность сиденья мягкая, передний край закругленный, а угол наклона спинки - регулируемый.

Необходимо предусматривать при проектировании возможность различного размещения документов: сбоку от видеотерминала, между монитором и клавиатурой и т.п. Кроме того, в случаях, когда видеотерминал имеет низкое качество изображения, например заметны мелькания, расстояние от глаз до экрана делают больше (около 700мм), чем расстояние от глаза до документа (300-450мм). Вообще при высоком качестве изображения на видеотерминале расстояние от глаз пользователя до экрана, документа и клавиатуры может быть равным.

Положение экрана определяется:

1) расстоянием считывания (0,6…0,7м);

2) углом считывания, направлением взгляда на 20 ниже горизонтали к центру экрана, причем экран перпендикулярен этому направлению.

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

1) по высоте +3 см;

2) по наклону от -10 до +20 относительно вертикали;

3) в левом и правом направлениях.

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

1) голова не должна быть наклонена более чем на 20,

2) плечи должны быть расслаблены,

3) локти - под углом 80…100,

4) предплечья и кисти рук - в горизонтальном положении.

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

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

Существенное значение для производительной и качественной работы на компьютере имеют размеры знаков, плотность их размещения, контраст и соотношение яркостей символов и фона экрана. Если расстояние от глаз оператора до экрана дисплея составляет 60…80 см, то высота знака должна быть не менее 3мм, оптимальное соотношение ширины и высоты знака составляет 3:4, а расстояние между знаками - 15…20% их высоты. Соотношение яркости фона экрана и символов - от 1:2 до 1:15.

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

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

Режим труда

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

В табл. 10 представлены сведения о регламентированных перерывах, которые необходимо делать при работе на компьютере, в зависимости от продолжительности рабочей смены, видов и категорий трудовой деятельности с ВДТ (видеодисплейный терминал) и ПЭВМ (в соответствии с СанПиН 2.2.2/2.4.1340-03 «Гигиенические требования к персональным электронно-вычислительным машинам и организации работы»).

Таблица 10. Время регламентированных перерывов при работе на компьютере

Категория работы с ВДТ или ПЭВМ	Уровень нагрузки за рабочую смену при видах работы с ВДТ	Суммарное время регламентированных перерывов, мин
	Группа А, количество знаков	Группа Б, количество знаков	Группа В, часов	При 8-часовой смене	При 12-часовой смене
I	до 20000	до 15000	до 2,0	30	70
II	до 40000	до 30000	до 4,0	50	90
III	до 60000	до 40000	до 6,0	70	120

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