Моделирование

МИНИСТЕРСТВО ОБРАЗОВАНИЯ И НАУКИ РОССИЙСКОЙ ФЕДЕРАЦИИ _____________ ФЕДЕРАЛЬНОЕ АГЕНТСТВО ПО ОБРАЗОВАНИЮ ____________ НАЦИОНАЛЬНЫЙ ИССЛЕДОВАТЕЛЬСКИЙ УНИВЕРСИТЕТ МОСКОВСКИЙ ЭНЕРГЕТИЧЕСКИЙ ИНСТИТУТ КУРС ЛЕКЦИЙ дисциплины вариативной части обязательных дисциплин Б1.В.ОД.7 «Моделированием» Направление подготовки 09.03.01 Информатика и вычислительная техника Профиль Вычислительные машины, комплексы, системы и сети Авторский коллектив: Профессор кафедры ВМСиС _________________ В.Н.Балашов Доцент кафедры ВМСиС ____________________ А.Г. Гольцов «_15_»_ноября_2016 г. Москва 2016 НИУ «МЭИ» 1 УДК МИНИСТЕРСТВО ОБРАЗОВАНИЯ И НАУКИ РОССИЙСКОЙ ФЕДЕРАЦИИ 621.398 Б202 ___________ ФЕДЕРАЛЬНОЕ АГЕНТСТВО ПО ОБРАЗОВАНИЮ ___________ НАЦИОНАЛЬНЫЙ ИССЛЕДОВАТЕЛЬСКИЙ УНИВЕРСИТЕТ МОСКОВСКИЙ ЭНЕРГЕТИЧЕСКИЙ ИНСТИТУТ В.Н. БАЛАШОВ, А.Г.ГОЛЬЦОВ МОДЕЛИРОВАНИЕ Учебное пособие по курсу "Моделирование" для студентов, обучающихся по направлению "Информатика и вычислительная техника" Москва Издательство МЭИ 2016 2 Курс лекций предназначен для "Моделирование". В курсе студентов, изучающих рассмотрены принципы курс математического моделирования применительно к техническим устройствам, изучаемым по направлению "Информатика и вычислительная техника". В курсе рассмотрены общие принципы создания математических моделей, формулировки целей моделирования и основы проведения вычислительного эксперимента. В рамках раздела курса «Системы массового обслуживания» Приведены основы теории Марковских процессов и вывод системы уравнений Колмогорова. Приведена аналитическая теория одноканальных и многоканальных систем массового обслуживания. В разделе «Имитационное моделирование» изложены подходы численного моделирования на примере имитационного моделирования систем массового обслуживания. В курсе кратко рассмотрены вопросы моделирования по технологии системной динамики и агентного моделирования, а также модели, основанные на сетях Петри. Лекции 1 – 6 написаны Балашовым В.Н., лекции 7 – 13 – Балашовым В.Н. и Гольцовым А.Г., лекции 14 – 18 – Балашовым В.Н. Предназначено для студентов, обучающихся по направлению 230100 "Информатика и вычислительная техника". © Национальный исследовательский университет «МЭИ» 2016 3 Оглавление Лекция 1. МОДЕЛИРОВАНИЕ И МОДЕЛИ............................................ 8 1.1. Методы математического моделирования .......................................... 8 1.2. Вычислительный эксперимент ........................................................... 10 1.3. Имитационное моделирование ........................................................... 13 1.4. Классификация моделей .................................................................... 16 1.4. Цель моделирования ........................................................................... 19 Лекция 2. РАЗРАБОТКА МАТЕМАМИЧЕСКОЙ МОДЕЛИ................ 23 2.1. Этапы разработки математической модели ...................................... 23 2.2. Системный подход к построению моделей....................................... 24 2.2.1.Структура математической модели ................................................. 24 2.1.2. Стадии разработки моделей ............................................................. 27 2.2. Оценка результатов моделирования ................................................. 31 Лекция 3. ПРОИЗВОДЯЩИЕ И ХАРАКТЕРИСТИЧЕСКИЕ ФУНКЦИИ ............................................................................................................. 32 3.1. Производящая функция....................................................................... 32 3.2. Характеристическая функция ............................................................. 33 Лекция 4. ДИСКРЕТНЫЙ МАРКОВСКИЙ ПРОЦЕСС......................... 40 4.1. Матрица переходных вероятностей ................................................... 41 4.2. Конечномерное распределение вероятностей .................................. 43 4.3. Матричный метод нахождения переходных вероятностей ............. 45 4.5. Марковский процесс с регулярной матрицей переходных вероятностей .......................................................................................................... 47 Лекция 5. НЕПРЕРЫВНЫЙ МАРКОВСКИЙ ПРОЦЕСС ..................... 51 4 5.1. Переходные вероятности непрерывного Марковского процесса ... 51 5.2. Уравнения Колмогорова ..................................................................... 56 Лекция 6. НЕПРЕРЫВНЫЙ МАРКОВСКИЙ ПРОЦЕСС ..................... 58 6.1. Пример Марковского процесса с непрерывным временем ............. 58 6.2. Стационарный режим для непрерывного Марковского процесса.. 60 Лекция 7. СИСТЕМЫ МАССОВОГО ОБСЛУЖИВАНИЯ .................. 63 7.1. Основные понятия и определения ..................................................... 63 7.2. Структура системы массового обслуживания .................................. 64 7.3. Дисциплина обслуживания ................................................................. 66 7.4. Характеристики системы массового обслуживания ....................... 67 7.5. Потоки заявок и потоки обслуживания ............................................. 68 7.6. Аналитическая модель системы массового обслуживания ............ 69 7.7. Состояния и интенсивность потока событий.................................. 69 Лекция 8. СИСТЕМЫ МАССОВОГО ОБСЛУЖИВАНИЯ ................ 71 8.1. Одноканальная СМО с ограниченной очередью .............................. 71 8.2. Одноканальная СМО с бесконечной очередью ................................ 73 8.3. Многоканальная СМО с очередью..................................................... 77 8.4. Система без очереди ............................................................................ 79 Лекция 9 СИСТЕМЫ МАССОВОГО ОБСЛУЖИВАНИЯ ................... 82 9.1. Система с ограниченной очередью .................................................... 82 9.2. Система с бесконечной очередью ...................................................... 84 Лекция 10. ИМИТАЦИОННОЕ МОДЕЛИРОВАНИЕ .......................... 89 10.1. Моделирование с постоянным шагом по времени ......................... 90 10.2. Моделирование по особым состояниям .......................................... 91 5 Лекция 11. ИМИТАЦИОННОЕ МОДЕЛИРОВАНИЕ .......................... 98 11.1. Случайные и псевдослучайные числа ............................................. 98 11.2. Период последовательности псевдослучайных чисел. .................. 99 11.3. Оценка качества последовательности случайных чисел ............. 100 11.4. Генерация равномерно распределенных случайных чисел ......... 105 11.4.1. Метод мультипликативного датчика. ......................................... 105 11.4.2. Метод середин квадратов............................................................. 106 11.5. Генерация произвольно распределенных случайных чисел ...... 107 11.5.1. Метод обратной функции ............................................................ 107 11.5.2. Метод исключения (режекции) ................................................. 109 Лекция 12. ИМИТАЦИОННОЕ МОДЕЛИРОВАНИЕ ......................... 112 12.1. Списки текущих и будущих событий ............................................ 112 12.2. Моделирование СМО как Марковского процесса ....................... 113 12.3. Моделирование СМО по событиям .............................................. 117 Лекция 13. ИМИТАЦИОННЫЕ И ФИЗИЧЕСКИЕ МОДЕЛИ ........... 126 13.1. Оценка результатов численного моделирования СМО ............... 126 13.2. Моделирование на основе физических моделей .......................... 130 Лекция 14.ДИНАМИЧЕСКАЯ СИСТЕМА И КОНЕЧНЫЕ АВТОМАТЫ ........................................................................................................ 134 14.1. Модель Динамическая Система ..................................................... 134 14.2. Модель Конечный Автомат ............................................................ 137 Лекция 15. СИСТЕМНАЯ ДИНАМИКА ............................................... 142 15.1. Элементы модели системной динамики ........................................ 142 15.3. Системы обратной связи ................................................................. 145 6 15.4. Нелинейность математической модели ......................................... 147 15.5. Возможность моделирования нетривиальных ситуаций ............. 148 Лекция 16. АГЕНТНО – ОРИЕНТИРОВАННЫЕ МОДЕЛИ.............. 150 16.1. Методы агентного моделирования ................................................ 150 16.2. Разработка агентно ориентированных моделей ........................... 153 Лекция 17. МОДЕЛИ, ОСНОВАННЫЕ НА СЕТЯХ ПЕТРИ .............. 157 17.1. Определение сети Петри ................................................................. 158 17.2. Объекты, образующие сеть Петри ................................................. 159 17.3. Маркировка сети Петри .................................................................. 160 17.4. Пространство состояний сети Петри ............................................. 163 17.5. Моделирование параллельных процессов .................................... 165 17.6. Моделирование процессора с конвейерной обработкой ............. 168 Лекция 18. МОДЕЛИ, ОСНОВАННЫЕ НА СЕТЯХ ПЕТРИ .............. 172 18.1. Кратные функциональные блоки компьютера ............................. 172 18.2.Взаимно исключающие параллельные процессы.......................... 173 18.3. Системы с флагами .......................................................................... 175 18.4. Анализ сетей Петри ......................................................................... 176 18.5. Дерево достижимости сети Петри ................................................. 180 7 Лекция 1. МОДЕЛИРОВАНИЕ И МОДЕЛИ 1.1. Методы математического моделирования Возникновение значительно и расширить совершенствование компьютеров позволяет интеллектуальные возможности человека, раскрепостить его внутренние силы. Впервые за всю историю человечества специалисты, владевшие ранее лишь пером и бумагой, получили в свои руки совершенный инструмент, отвечающий требованиям научно – технического прогресса. Однако необходимо помнить, что компьютер – всего лишь инструмент, возможности которого раскрываются только в сочетании со всеми существующими методами исследований. Это в основном новый научный метод, получивший название математическое моделирование. При математическом моделировании имеют дело не с самим явлением, а с некоторым теоретическим «слепком» с него, с моделью, выражающей в математической форме основные закономерности, которым она подчиняется. Как показывают многочисленные удачные примеры, этот «слепок» может быть выполнен весьма искусно и совершенно. В результате исследователь, проводя вычислительный эксперимент, испытывает как бы саму природу, задавая ей вопросы и получая строгие и относительно полные ответы. Возможность замены исходного объекта его математической моделью и в дальнейшем диалог с ней дает большие преимущества и означает серьезные изменения идеологии научных исследований и технических разработок. Применение компьютерных технологий позволяет существенно повысить производительность труда при проектировании новой техники, повысить качество проекта и технологические возможности оборудования. 8 При использовании математического моделирования необходимо учитывать следующие обстоятельства: Развитие цивилизации все больше идет в условиях ограниченности материальных и людских ресурсов. Поэтому для достижения положительного результата в условиях конкуренции необходимо выбирать наиболее экономичные решения. Как правило, необходимо выбирать решение, позволяющее максимально экономить время. Проблемы, которые необходимо решать, укрупняются, становятся все более сложными для решения. Для их решения необходимо применять системный подход, состоящий в учете не только основных, но и сопутствующих факторов, что приводит к существенному усложнению методов исследований и математических моделей. Мир нелинеен, то есть законы, описывающие поведение неживой и живой природы, поведение экономических и социальных структур, являются нелинейными. Это означает, что результаты развития не однозначно определяются настоящим. Эти результаты нельзя предсказать, опираясь только на предшествующий опыт. Поэтому для получения достоверных результатов необходимо проводить решение задач в максимально точной постановке. Необходимость решения все более сложных проблем ставит задачу изменения всей методологии научной и инженерной работы, позволяющей добиться оперативности, гибкости и универсальности принятия решений. Новая методология должна базироваться на математическом моделировании с использованием современных компьютеров. Разработка дешевых высокопроизводительных компьютеров привело к возникновению новой цивилизации, основанной на информационных технологиях. Условно можно выделить три основные области применения компьютеров: 9 1. Создание и широкое применение информационных систем, включающих системы передачи и обработки информации в реальном времени, создание баз данных, создание глобальных систем электронных. Создание глобальных электронных систем оптовой и розничной торговли, кооперации в промышленности, электронного заказа продукции прямо на предприятии – изготовителе. Информационные системы облегчают доступ к образованию, облегчают быт и доступ к системам развлечений и досуга. 2. Автоматизация различных видов человеческой деятельности, в частности управленческой, научной, инженерной. Эти основанные на информационных системы системах автоматизированные автоматизированного автоматизированного управления методы управления, проектирования, производством, системы создание гибких производственных систем. 3. Математическое разнообразной природы. моделирование Методы сложных математического процессов моделирования, позволяют получать достаточно точные количественные характеристики сложных технических изделий и производственных процессов. Эти три области применения компьютеров на практике неразделимы и их дальнейшее объединение неизбежно. Разработка модели, алгоритма и компьютерной программы является обычным в каждой из этих областей. 1.2. Вычислительный эксперимент Развитие методов математического моделирования привело к понятию вычислительного эксперимента, позволяющего подробно и детально изучать свойства сложных технических и экономических систем, существующих только в форме компьютерной программы. Выделим основные этапы проведения вычислительного эксперимента: 10 • На первом этапе создается математическая модель изучаемого объекта. Для этого выделяются основные свойства объекта с точки зрения решаемой задачи (цели моделирования). Связи, существующие в исследуемом объекте, описываются при помощи математических уравнений, отражающих, как правило, фундаментальные законы природы. Идеализация, упрощение исходной системы является основой любого научного исследования. Однако возможности современных компьютеров позволяют создавать сложные и поэтому точные модели, недоступные для традиционных методов анализа. Традиционные упрощенные аналитические модели, отражающие только отдельные свойства изучаемой системы, оказываются пригодными только на начальных (предварительных) этапах моделирования. Традиционная математика работает в основном с линейными задачами. Линейность в математическом смысле означает, что сумма (суперпозиция) решений линейной задачи также является решением. Линейность в более широком смысле означает, что поведение части объекта схоже с поведением всего объекта, что возможно предсказание поведения объекта по поведению его части. Принцип суперпозиции широко и без достаточных оснований используется в науке и технике. Однако потребности решения практических задач методами вычислительного эксперимента требует глубокой разработки методов решения нелинейных задач. Основная сложность здесь – доказательство существования математической и единственности модели, решения преодоление нелинейных трудностей, уравнений связанных с некорректность постановки задачи и связанной с этим неустойчивостью счета. • На втором этапе вычислительного эксперимента проводится изучение поведения математической модели при различных значениях входящих в нее параметров. Для этого используется аппарат вычислительной 11 математики, позволяющей разработать конечно – разностные аналоги дифференциальных уравнений математической модели, оценить точность аппроксимации и устойчивость счета, построить вычислительный алгоритм решения задачи. • На третьем этапе составляется компьютерная программа, реализующая алгоритм вычислительного эксперимента. Современные реализации языков программирования позволят существенно сократить трудоемкость составления и отладки программы, позволяют формально и без технических трудностей написать программы, реализующие сложные вычислительные алгоритмы, характерные для вычислительного эксперимента. Следует отметить, что явно проявилась тенденция разработки сложных проблемно – ориентированных пакетов программ, позволяющих решать широкий класс задач математического моделирования. Это позволяет в ряде случаев отказаться от разработки собственных программ, предотвратить дублирование в разработках. • На четвертом этапе определяются свойства среды, входящие в математическую модель в форме числовых коэффициентов. Для этого проводится натурный эксперимент или проводится вспомогательный вычислительный эксперимент, базирующийся на экспериментальных данных и соответствующих математических моделях. Важным здесь является калибровка математических моделей, основанная на эксперименте. • На пятом этапе проводится экспериментальных исследований. Дело экспериментах всегда измерить не можно в обработка том, что величины, результатов в реальных входящие в математическую модель. Измеряются другие величины, связанные с основными величинами при помощи своей математической модели. Проблема обработки результатов натурного эксперимента в рамках математического моделирования является очень актуальной. 12 • На шестом этапе проводится вычислительный эксперимент. На этом этапе в подготовленные программы подставляются исходные данные, полученные на 4 и 5 этапах, и проводится счет. Во многом вычислительный эксперимент похож на натурный. В нем также задаются различные варианты исходных данных, многопараметрических а результаты графиков счета или представляются массивов чисел, в форме описывающих поведение объекта. 1.3. Имитационное моделирование Процесс создания новой техники неразрывно связан с проведением расчетов. Некоторые из этих расчетов достаточно просты, требуют только знания формул и проводятся на калькуляторе. Другие расчеты требуют использования сложных компьютерных программ и выполняются на высокопроизводительных компьютерах. При проведении расчета каждому набору исходных данных ставится в соответствие результат. Поэтому, проведя большое количество расчетов с разными наборами исходных данных и оценив полученные результаты, разработчик может получить некоторое представление о работоспособности проектируемого устройства. Такой подход называется методом вариантного счета. При вариантном счете забота о рациональном выборе исходных данных и оценке результатов лежит на разработчике. При большом числе вариантов возрастает трудоемкость и сложность оценки результатов, при малом числе вариантов можно не найти подходящее (оптимальное) решение. Попытки повысить эффективность принятия решений на основе использования компьютеров привели к понятию имитационного моделирования. Существенным в моделировании является то, что выбор исходных данных и анализ результатов вариантного счета поручается компьютеру. В простейшем случае определяется область допустимого 13 изменения исходных данных, а затем расчетные варианты формируются перебором с заданным шагом. Для анализа результатов счета применяется показатель эффективности (оптимальности) системы. Обычно это аналитическое выражение (формула), в которое результаты счета входят с весовыми коэффициентами. Один из наборов исходных данных, соответствующих экстремуму показателя эффективности, принимается за оптимальное решение. Моделирование состоит в наблюдении за поведением модели под действием исходных данных или входных воздействий. Часть входных воздействий оказывает существенное влияние на поведение модели (основные факторы), другие влияют несущественно или вообще не влияют. Задача выделения основных факторов, присутствующих в исходных данных, также входит в процесс моделирования. Часть или все входные воздействия могут носить случайный характер. Тогда и результаты счета будут включать в себя случайные компоненты. В итоге результаты моделирования, по сути, оказываются экспериментальными данными, такими же, как и результаты натурного эксперимента. Метод моделирования, часть входных воздействий в котором носит случайный характер, получил название имитационное моделирование. Методика имитационного моделирования совпадает с методикой натурного эксперимента, в котором проводится статистическая обработка результатов измерений, позволяющая определить средние значения полученных величин и оценить точность измерений. Натурный эксперимент также является статистическим, так как на его результаты всегда оказывают влияние случайные воздействия. Имитационное моделирование особенно эффективно на начальном этапе разработки, когда изучается поведение объекта, заданного алгоритмами работы, под действием случайных факторов. 14 Эффективность имитационного моделирования может быть существенно повышена на основе теории статистического эксперимента. В его основе лежит метод статистических испытаний (метод Монте – Карло), позволяющий при сохранении точности значительно сократить число вариантов счета. В этом методе принимается, что входные воздействия являются случайными величинами, распределенными по какому то подходящему закону. Поэтому и в результате эксперимента получается случайная величина. Моделирование по методу Монте – Карло состоит в проведении серии испытаний с входными воздействиями, выбранными с использованием генераторов случайных чисел. Тогда результаты моделирования представляют собой случайную выборку, которую можно обработать статистическими методами и получить оценки интересующих характеристик моделируемого устройства. Теоретической основой метода статистических испытаний является теория вероятностей, в частности, предельные теоремы теории вероятностей, позволяющие оценить точность результатов в зависимости от числа испытаний (прогона математической модели). Каждый прогон отличается случайными последовательностями чисел, образующих варианты входных воздействий. Отметим основные особенности имитационного моделирования, вытекающие из теории статистического эксперимента. Каждый прогон модели можно рассматривать как одно наблюдение при проведении статистического эксперимента. Увеличение продолжительности прогона приводит к уменьшению отклонения измеряемой величины от точного значения, так как система постепенно переходит в стационарное состояние. Влияние переходных процессов можно уменьшить, увеличивая число прогонов. 15 Слишком большая продолжительность прогона не приводит к повышению точности результатов. Результаты имитационного моделирования необходимо рассматривать как экспериментальные данные и применять к ним соответствующие методы обработки и анализа. Имитационное моделирование требует предварительного построения плана эксперимента, включающего ответы на следующие вопросы: • Какова оптимальная продолжительность прогона? • Являются ли результаты статистически независимыми наблюдениями? • Сколько прогонов нужно провести для получения заданной точности? Рациональное планирование эксперимента, адекватно отвечающее поставленной задаче, позволяет значительно повысить эффективность имитационного моделирования. 1.4. Классификация моделей Модели технических систем разделяются на две группы: • Физические модели, и • Математические модели. Физическая модель – это другая, обычно более простая техническая система, которая по всем или по части параметров функционирует подобно оригиналу. Физическая модель и оригинал могут иметь одинаковую или разную физическую природу. К физическим моделям в частности относят натурные модели, масштабные модели и функциональные модели. Натурные модели включают макеты и опытные образцы, созданные в процессе разработки технической системы. Они имеют практически полное совпадение с оригиналом. Испытание макетов и опытных образцов 16 позволяют получить функционировании наиболее системы. По достоверную результатам информацию испытаний о проводится доработка опытного образца с целью устранения выявленных недостатков. Масштабные модели – это технические системы, выполненные в масштабе по сравнению с оригиналом и имеющие ту же физическую природу. Это, например, модели самолетов, предназначенные для продувки в аэродинамических трубах. Основой моделирования при использовании масштабных моделей является теория подобия, широко применяемая в аэродинамике и теплотехнике. При разработке вычислительных систем масштабное моделирование можно применять, например, для изучения условий теплообмена в конструктивах (шкафах, корпусах) компьютеров. Функциональные модели – это обычно электронные аналоговые функциональные преобразователи, которые после соответствующей настройки и соединения в электрическую схему позволяют моделировать систему, описываемую сложной функциональной зависимостью, дифференциальным или алгебраическим уравнением. Такие системы входят в состав вычислительных машин непрерывного действия, которые широко применялись в недалеком прошлом. Функциональные модели можно использовать при моделировании электронных компонентов вычислительной техники. Математические модели представляют собой формализованное описание технической системы или ее частей средствами математики, в частности системами алгебраических, дифференциальных или интегральных уравнений. При создании математических моделей широко применяются физические законы, теория алгоритмов, математическая теория систем, математическая статистика, теория вероятностей и теория случайных процессов. Понятно, что каждой технической системе может быть сопоставлено несколько различных математических моделей, отражающих существенные 17 свойства этой системы. Изучение каждой модели позволяет найти ответ на определенную группу вопросов. Выбор моделей определяется целями и задачами моделирования. Выделим основные виды математических моделей. Это в первую очередь физические математические модели, в основу которых положено достаточно точное математическое описание технической системы. Примерами таких моделей является описание явлений, подчиняющихся классическим дифференциальным уравнениям математической физики (колебания струн и пластин, распространение радиоволн и волн в упругой среде, теплопроводность, диффузия и т.д.). Обыкновенные дифференциальные уравнения второго порядка достаточно точно описывают колебания и переходные процессы в электрических цепях и поэтому используются в математических моделях схемотехнического моделирования в аналоговой и цифровой электронике. На первоначальных этапах моделирования проводятся значительные упрощения (идеализация свойств элементов технической системы), позволяющие построить математическую модель на основе линейных дифференциальных уравнений. Иногда линейные уравнения удается решить и получить решение в аналитическом виде в виде формулы, позволяющей строить диаграммы и графики. Такая математическая модель называется аналитической. При более точном моделировании технической системы используются математическое описание, основанное на нелинейных уравнениях. Нелинейные уравнения, как правило, не решаются аналитически и поэтому для их решения применяют численные методы, реализуемые в форме программ для компьютера. По результатам счета также строятся наглядные диаграммы и графики, позволяющие получить наглядное описание результатов моделирования. Такая математическая модель называется численной. 18 Следует отметить, что по мере роста производительности компьютеров, разработке и освоению численных методов решения и языков программирования становится применение наиболее численных целесообразным. Уже математических разработаны моделей и широко применяются сложнейшие численные модели аэродинамики самолетов на всех режимах полета, модели работы ядерных реакторов и других сложнейших технических систем. Счет по этим моделям на суперкомпьютерах позволяет исследовать реальные системы с точностью до единиц процента. Имитационная модель строится как описание технической системы в форме алгоритмов и правил ее функционирования под влиянием внешних воздействий. Поэтому имитационные модели могут быть созданы для очень широкого класса систем, в том числе информационных и вычислительных. При создании имитационной модели часто применяют специализированный язык программирования программирования MATLAB). (например, (например, Имитационное пакет GPSS) или расширения моделирование графический системы Simulink проводится язык в форме вычислительного эксперимента, в процессе которого на входы модели подаются последовательности входных воздействий, и проводится анализ поведения модели в ответ на эти воздействия. В качестве входных воздействий, как правило, используются последовательности случайных чисел, а анализ поведения системы проводится методами математической статистики. В результате вычислительный эксперимент имитирует реальное поведение технической системы. 1.4. Цель моделирования Одной из важнейших вспомогательных задач, которую необходимо решить в процессе моделирования, является формулировка цели 19 моделирования. Цель моделирования позволяет выбрать одну из моделей, которые могут быть сопоставлены исследуемой технической системе, определить необходимую точность моделирования и оценить затраты на моделирование. Например, при моделировании архитектуры компьютера в качестве цели можно производительности при принять решении достижение заданного максимальной класса задач. При моделировании электронной схемы функционального узла компьютера в качестве цели можно принять достижение максимального быстродействия или (и) получение минимальной стоимости или (и) получение минимальной потребляемой одновременное мощности. достижение Понятно, этих что целей в большинстве нереально. случаев Поэтому цель моделирования может быть определена в форме подходящего формального критерия эффективности, позволяющего учесть важность противоречивых целей. Правильная формулировка цели моделирования во многих случаях представляет собой очень сложную задачу. Очень часто первоначально поставленная цель оказывается неэффективной или ее не удается формализовать. Поэтому в процессе моделирования часто приходится изменять или уточнять критерий эффективности и выбранную для моделирования математическую модель. Критерий эффективности работы системы строится в форме, адекватной форме выбранной математической модели. Для аналитической модели критерий эффективности строится в форме непрерывной функции или функционала, при этом оптимальное решение можно найти методами вариационного исчисления или решением задачи на условный экстремум (метод множителей эффективности можно Лагранжа). определить Для в численной рамках модели задач критерий линейного или динамического программирования. Для имитационного моделирования эффективность определяется в рамках теории статистических решений. Не 20 следует отбрасывать и методы экспертных оценок, так как очень часто в процессе моделирования хороший специалист быстро находит приемлемое техническое решение. Для определения эффективности технической системы в аналитической форме формулируются частные критерии качества yi, зависящие от ряда параметров {xj}. В результате получаем ряд аналитических выражений y 1  1 ( x1 ,...x j ,...xm );  y i  i ( x1 ,...x j ,...xm );  y n  n ( x1 ,...x j ,...xm ); (частных критериев качества технической системы), значения которых должны находиться в области допустимых значений. yi min  y i  y i max Значения параметров xi также должны находиться в допустимых пределах x j min  x j  x j max Дальнейший ход определения эффективности может проводиться при помощи ряда математических методов параметрической оптимизации. Так, если частные критерии качества yi являются линейными функциями относительно параметров xj, получаем задачу линейного программирования, если нелинейными – задачу динамического программирования. Если один из критериев выбирается главным и параметрическая оптимизация проводится по достижению экстремума только этого критерия, то получаем задачу однокритериальной оптимизации. Если из частных критериев составляется 21 обобщенный критерий эффективности, то приходим к задаче многокритериальной оптимизации. В результате под целью моделирования обычно принимается определение оптимального набора параметров, при котором выбранный критерий эффективности принимает экстремальное значение. Это является задачей оптимизацией технической системы. Важно определить и необходимую точность моделирования. В результате формулировки цели моделирования появляется возможность оценить затраты на моделирование, которые должны окупиться в результате эксплуатации оптимизированной технической системы по сравнению с не оптимизированной. Таким образом, цель моделирования определяется контекстом решаемой задачи. Цель моделирования включает в себя широкий класс понятий. Это, например, простейший перебор нескольких допустимых вариантов технических решений, обеспечивающих максимальную производительность вычислительной системы на заданных тестовых задачах, или сложные задачи многопараметрической оптимизации, обеспечивающие автоматический поиск оптимальных, часто совершенно не очевидных технических решений. 22 Лекция 2. РАЗРАБОТКА МАТЕМАМИЧЕСКОЙ МОДЕЛИ 2.1. Этапы разработки математической модели В процессе разработки математической модели технической системы можно выделить следующие этапы. Это разработка концепции, разработка модели и разработка программы, реализующей математическую модель. Концептуальное (содержательное) описание модели включает в себя словесное описание технической системы и ее составных частей, описание свойств и законов взаимодействия отдельных частей системы, описание входных воздействий взаимодействия и входных технической сигналов системы и системы, описание окружающей среды. Концептуальную модель может разработать специалист, детально знающий объект моделирования. Основная задача, решаемая при разработке концепции моделирования, состоит в выделении основных элементов и явлений, которые должны быть включены в модель и выделении второстепенных с точки зрения конкретной цели моделирования элементов и явлений, которые необходимо отбросить. Основной сложностью в разработке концептуальной модели является интуитивная оценка адекватности технической системы и ее математической модели. В сложных случаях этот вопрос может быть окончательно решен только значительно позже, в процессе моделирования и сравнения результатов моделирования и экспериментальных параметров технической системы. Опыт и интуиция разработчика на данном этапе является основой для достижения успеха. Важным этапом в разработке концептуальной модели является ее детализация. Необходимо правильно выделить отдельные подсистемы и блоки, сделать их описание в форме формул, алгоритмов или законов функционирования, составить описание законов взаимодействия этих блоков 23 и подсистем. Необходимо оценивать сложность (размерность) модели, так как дальнейшая работа с очень сложной детальной моделью может оказаться непосильной для коллектива разработчиков. Целесообразно на первых этапах существенно упрощать модель, а в дальнейшем усложнять ее в тех подсистемах, которые не обеспечивают достижение необходимой точности моделирования. Следующим этапом разработки концептуальной модели является описание входных и выходных сигналов, связывающих математическую модель и окружающую среду. Техническая система связана с окружающей средой (в широком смысле слова) при помощи датчиков, интерфейсов, исполнительных и периферийных устройств, пультов управления и других подобных устройств. Эти устройства обычно не включаются в модель и заменяются подходящими эквивалентами. Входные воздействия на техническую систему моделируются при помощи различных генераторов входных сигналов, выходные величины обычно анализируются в форме, естественно представленной в математической модели. 2.2. Системный подход к построению моделей 2.2.1.Структура математической модели Сложные технические системы имеют иерархическую структуру, включающую подсистемы различных уровней. Система S - множество взаимодействующих подсистем, действующие для достижения заданной цели. Внешняя среда Е для системы S, это множество элементов любой природы, взаимодействующие с системой S. В зависимости от цели могут рассматриваться разные варианты взаимодействия между системой S и внешней средой Е. 24 Системный подход позволяет наиболее точно отобразить процесс взаимодействия системы с окружающей средой. При системном подходе к моделированию систем необходимо, прежде всего, четко определить цель моделирования. Определение цели позволяет выделить подсистемы, и определить, какие из них войдут в математическую модель М. Этот отбор проводится на основе критерия, сформулированного на основе цели моделирования. При системным подходе к моделированию необходимо определить структуру системы, учитывающую связи между элементами с учетом их взаимодействия. Структура системы может исследоваться как в общем виде, с учетом взаимодействия отдельных подсистем, так и в частном виде, когда исследуются отдельные свойства подсистем с учетом возможности достижения цели моделирования. В соответствии с этим применяются различные подходы к исследованию структуры системы, в частности структурный и функциональный подходы. При структурном подходе определяется состав подсистем системы S и связи между ними. Совокупность подсистем и связей между ними позволяет судить о структуре системы. Структура, в зависимости от цели исследования может быть описана на разных иерархических уровнях. Обычно применяется топологическое описание, позволяющее определить составные части системы, хорошо формализуемое на базе теории графов. При функциональном подходе рассматриваются отдельные функции, т. е. алгоритмы поведения системы, и реализуется функциональный подход, оценивающий функции, которые выполняет система, причем под функцией понимается свойство, приводящее к достижению цели. Поскольку функция отображает свойство, а свойство отображает взаимодействие системы S с внешней средой Е, то свойства могут быть выражены в виде либо некоторых характеристик элементов Si(j) и подсистем Si,- системы, либо системы S в целом. 25 Проявление функций системы во времени S(t), т. е. функционирование системы, означает переход системы из одного состояния в другое, т. е. движение в пространстве состояний Z. При эксплуатации системы S весьма важно качество ее функционирования, определяемое показателем эффективности и являющееся значением критерия оценки эффективности. Существуют различные подходы к выбору критериев оценки эффективности. Система S может оцениваться либо совокупностью частных критериев, либо некоторым общим интегральным критерием. Создаваемая модель М с точки зрения системного подхода также является системой, т. е. S' = S'(M), и может рассматриваться по отношению к внешней среде Е. Построение математической модели. Подход к изучению взаимосвязей между отдельными частями модели предусматривает рассмотрение их как отражение связей между отдельными подсистемами объекта. Процесс создания модели М представлен на рис. 2.1. Реальный объект, подлежащий моделированию, разбивается на отдельные подсистемы, т. е. выбираются исходные данные для подсистем и ставятся частные цели, соответствующие работе подсистем объекта моделирования. На основе этого формируется компонента математической модели, соответствующая рассматриваемой подсистеме. Совокупность компонент объединяется в модель М. В результате разработка модели состоит в объединении отдельных компонент в единую модель, причем каждая из компонент решает свою частную задачу. 26 Подсистемы реального объекта Частная цель моделирования математическая модель Частная цель 1 Компонента модели 1 Частная цель k Компонента модели k Набор данных 1 Набор данных 2 Набор данных 3 Набор данных i Набор данных j Рис. 2.1. Процесс создания математической модели. Отметим особенность рассмотренного метода построения модели:  Модель строится по частям, соответствующим частным целям моделирования;  Модель строится путем объединения отдельных компонент и не учитывает явно возможность возникновения системного эффекта. Сложные математические модели создаются путем объединения более простых моделей, при этом явления, связанные с взаимодействием подсистем между собой и с окружающей средой учитываются итерационным методом путем добавления в модель новых подсистем и связей. 2.1.2. Стадии разработки моделей При разработке моделей выделяются две основные стадии проектирования: • макропроектирование; • микропроектирование. 27 На стадии макропроектирования на основе данных о реальной системе S и внешней среде Е строится модель внешней среды, выявляются ресурсы и ограничения для построения модели системы, выбирается модель системы и критерии, позволяющие оценить адекватность модели М реальной системы S. Построив модель системы и модель внешней среды, на основе критерия эффективности функционирования системы в процессе моделирования выбирают оптимальную стратегию управления, что позволяет реализовать возможности модели по воспроизведению отдельных сторон функционирования реальной системы S. Стадия микропроектирования в значительной степени зависит от конкретного типа выбранной модели. В случае имитационной модели необходимо обеспечить создание информационного, математического, технического и программного обеспечений системы моделирования. На этой стадии можно установить основные характеристики созданной модели, оценить время работы с ней и затраты ресурсов для получения заданного качества соответствия модели процессу функционирования системы S. Независимо от типа используемой модели М при ее построении необходимо руководствоваться рядом принципов системного подхода: • пропорционально-последовательное продвижение по этапам и направлениям создания модели; • согласование информационных, ресурсных и других характеристик; • правильное соотношение отдельных уровней иерархии в системе моделирования; • целостность отдельных обособленных стадий построения модели. Модель М должна отвечать заданной цели ее создания, поэтому отдельные части должны компоноваться взаимно, исходя из единой системной задачи. Цель может быть сформулирована качественно, тогда она будет обладать большей содержательностью и длительное время может 28 отображать объективные возможности данной системы моделирования. При количественной формулировке цели возникает целевая функция, которая точно отображает наиболее существенные факторы, влияющие на достижение цели. Построение модели относится к числу системных задач, при решении которых синтезируют решения на базе огромного числа исходных данных, на основе предложений больших коллективов специалистов. Использование системного подхода в этих условиях позволяет не только построить модель реального объекта, но и на базе этой модели выбрать необходимое количество управляющей информации в реальной системе, оценить показатели ее функционирования и тем самым на базе моделирования найти наиболее эффективный вариант построения и выгодный режим функционирования реальной системы S. В основе любого вида моделирования лежит некоторая модель, имеющая соответствие, базирующееся на некотором общем качестве, которое характеризует реальный объект. Объективно реальный объект обладает некоторой формальной структурой, поэтому для любой модели характерно наличие некоторой структуры, соответствующей формальной структуре реального объекта, либо изучаемой стороне этого объекта. Одним из наиболее важных аспектов построения систем моделирования является проблема цели. Любую модель строят в зависимости от цели, которую ставит перед ней исследователь, поэтому одна из основных проблем при моделировании — это проблема целевого назначения. Подобие процесса, протекающего в модели М, реальному процессу является не целью, а условием правильного функционирования модели, и поэтому в качестве цели должна быть поставлена задача изучения какой-либо стороны функционирования объекта. Если цель моделирования ясна, то возникает следующая проблема, а именно проблема построения модели М. Построение модели оказывается 29 возможным, если имеется информация или выдвинуты гипотезы относительно структуры, алгоритмов и параметров исследуемого объекта. Для правильно построенной модели М характерным является то, что она выявляет лишь те закономерности, которые нужны исследователю, и не рассматривает свойства системы S, не существенные для данного исследования. Следует отметить, что оригинал и модель должны быть одновременно сходны по одним признакам и различны по другим, что позволяет выделить наиболее важные изучаемые свойства. В этом смысле модель выступает как некоторый «заместитель» оригинала, обеспечивающий фиксацию и изучение лишь некоторых свойств реального объекта. Последним этапом разработки модели является предварительная (умозрительная) проверке ее адекватности исследуемой технической системе. Эту проверку обычно проводят эксперты. Форма математической модели определяет метод моделирования. Если в основу модели положены физические законы в форме аналитических зависимостей, то применяются численные методы и реализующие их программы. Это позволяет достичь максимальной точности моделирования. Проверка адекватности является важной и очень сложной задачей, однако в некоторых частных случаях проверка достаточно проста и очевидна. Классическим методом является моделирование простой по геометрии системы, для которой известно аналитическое решение. Тогда по точности совпадения аналитического и численного решения можно судить о точности моделирования. Вторым классическим методом проверки адекватности является моделирование (контрольный расчет) существующей технической системы (аналога), подобной той, для которой разработана математическая модель. В этом случае результаты моделирования сравниваются с экспериментальными данными. При проверке адекватности проверяется не только вся модель, но и ее отдельные части и блоки. 30 Проверка адекватности необходима всегда, так как иначе по результатам моделирования могут быть приняты неверные технические решения. В проверку адекватности входит и анализ области допустимого изменения входных сигналов и параметров. При необходимости проводится корректировка модели, позволяющая повысить, а иногда и понизить точность моделирования. 2.2. Оценка результатов моделирования Современные компьютеры позволяют без труда представлять результаты моделирования в графическом виде, удобном для анализа человеком. Обычно строятся трехмерные семейства графиков, позволяющие детально представить особенности функционирования технической системы. Сопоставляя в графическом виде входные воздействия, внутренние переменные и выходные величины модели удается определить реакцию модели на входные воздействия и оценить возможности и достижимые технические характеристики проектируемой системы. Результаты моделирования позволяют определить оптимальную конфигурацию, найти оптимальные технические решения и оценить технические возможности проектируемой технической системы. Вопросы: 1. Сколько моделей может быть у моделируемого объекта? 2. Чем определяется выбор модели для моделируемого объекта? 3. В чем состоит процесс согласования целей при моделировании? 4. Как проверяется адекватность математической модели? 5. Как проводится вычислительный эксперимент? 31 Лекция 3. ПРОИЗВОДЯЩИЕ И ХАРАКТЕРИСТИЧЕСКИЕ ФУНКЦИИ 3.1. Производящая функция Производящие и характеристические функции образуют удобный математический аппарат теории вероятностей, позволяющий достаточно просто находить числовые характеристики случайных величин. Рассмотрим дискретную случайную величину ξ, принимающую целочисленные значения 0, 1, 2, …, k с вероятностями p0, p1, p2, …, pk. Функция переменной z (0 < z ≤ 1) k (z)   pi z  p 0  p1z  p 2 z    p k z i 2 k i0 называется производящей функцией дискретной случайной величины ξ. Производящая функция представляет собой конечный или бесконечный степенной ряд, коэффициенты которого pi являются вероятностями в распределении дискретной случайной величины. Зная производящую функцию, легко получить основные числовые характеристики случайной величины. Дифференцируя производящую функцию n раз по z и приравнивая z = 0, получим вероятность pn состояния Sn случайной величины ξ: pn  n 1 d ( z ) n! dz n z  0 Дифференцируя производящую функцию один раз по z и приравнивая z = 1, получим математическое ожидание случайной величины Mξ. M  d( z)  p1  2p2 dz z 1  3p3    kpk 32 Дифференцируя производящую функцию два раза по z и приравнивая z=1, получим второй нецентральный (начальный) момент случайной величины Mξ2 2 d ( z ) 2 2 z 1  2p 2  3  2p3    k ( k  1) p k  M  ( M) 2 dz Из этой формулы получаем выражение для дисперсии 2 2 2 d ( z ) 2 d( z ) d( z ) D  M  (M)  2 z 1  dz z 1 ( dz z 1 ) dz Пример: Рассмотрим дискретную случайную величину ξ со следующей функцией распределения Значение ξ Вероятность p 1 2 3 0.1 0.3 0.2 0.4 Производящая функция для нее имеет следующий вид 2 3 (z)  0.1  0.3p  0.2p  0.4p Математическое ожидание M  0.3  2  0.2  3  0.4  1.9 Второй начальный момент 2 M  2  0.2  3  2  0.4  2.8 Дисперсия 2 D  2.8  1.9  (1.9)  1.09 3.2. Характеристическая функция Характеристической функцией ψ(t) непрерывной случайной величины ξ называется преобразование Фурье ее плотности распределения p(x).   itx ( t )  e p( x )dx  Характеристической функцией ψ(t) дискретной случайной величины ξ называется ряд Фурье следующего вида 33 n ( t )   e  itx j j1 Дискретная случайная pj величина ξ принимает значения с xj вероятностями pj. Для произвольного закона распределения случайной величины F(x) характеристической функцией называется интеграл Фурье – Стилтьеса следующего вида  e ( t )   itx dF( x )  Свойства характеристической функции совпадают со свойствами преобразования Фурье и поэтому подробно описаны в математической литературе. Приведем только основные свойства преобразования Фурье. Преобразование Фурье базируется на формуле интеграла Фурье, записанного в комплексной форме.  f (x)  1 2   d  f (t )e  i ( t  x ) dt   Функция f(t) во внутреннем интеграле должна быть абсолютно интегрируемой. Внешний интеграл понимается в смысле главного значения. Интеграл Фурье можно разделить на два интеграла. Назовем прямым преобразованием Фурье функции p(x) следующее выражение  ( t )  e  itx p( x )dx  Тогда обратным преобразованием Фурье является следующий интеграл  p( x )  1 2  (t )e dt  itx 34 Преобразование Фурье каждой абсолютно интегрируемой функции p(x) ставит в соответствие другую функцию ψ(t), определенную на всей числовой оси. Обратное преобразование Фурье устанавливает связь между функцией ψ(t) и функцией p(x). Необходимо обратить внимание на сходство формул прямого и обратного преобразований Фурье, отличающихся только знаком в показателе и множителем 1/2π. Однако эти две формулы по существу различны. В первой формуле участвует обычный интеграл, зависящий от параметра t. Вторая формула, по сути, является другой записью интеграла Фурье и содержит утверждение, что стоящий справа интеграл равен исходной функции p(x). p( x )  e Пример. Пусть |x| - четная функция. Преобразование Фурье этой функции равно  ( t )  e   x  itx e  dx     2 e  x  x (cos( tx )  i sin( tx ))dx   cos( tx )dx  e 2 2 t  2 ; Аппарат характеристических функций упрощает работу с суммами независимых случайных величин. Пусть независимые случайные величины ξ1 и ξ2 имеют плотности распределения p1(x) и p2(x). Пусть ξ – новая случайная величина, равная сумме ξ1 и ξ2 Для нее функция плотности распределения вероятности r(x) равна свертке функций распределения p1(x) и p2(x).    1   2 . r ( x )   p1 (s)p 2 (s  x )ds 35 Пусть ψ1(t) и ψ2(t) – преобразования Фурье плотностей распределений p1(x) и p2(x).  1 ( t )   p1(x )e  itx  dx и 2 (t)   p 2 ( x )e  itx dx   Пусть R(t) – преобразование Фурье свертки r(x).   itx R (t)   r ( x )e dx  Тогда по теореме о преобразовании Фурье свертки получим R ( t )  1 ( t ) 2 ( t ) Это свойство преобразования Фурье свертки позволяет существенно упростить вычисления при работе с распределениями сумм случайных величин. Для нахождения функции плотности распределения r(x) суммы нескольких случайных величин ξ1 + ξ 2 + … + ξ n необходимо найти преобразование Фурье ψi(t) для каждой функции плотности распределения pi(x), перемножить преобразования Фурье R(t) = ψ1(t) ψ2(t) … ψn(t), а затем вычислить обратное преобразование Фурье. В результате получим функцию плотности распределения суммы независимых случайных величин. Для одинаково распределенных случайных величин преобразование Фурье R(t) для суммы равно  r(x)  1 2  R (t )e dt ixt  где n R (t)   (t) 36 Правила дифференцирования и интегрирования преобразований Фурье следуют из соответствующих правил дифференцирования и интегрирования интеграла по параметру. Дифференцируя и интегрируя интеграл обратного преобразования Фурье по параметру x, получим  dp ( x ) dx  ix 2  (t )e dt itx   x F( x )   p( x )dx  1 2 ix  (t )e dt itx  Пример. Характеристическая функция для случайной величины ξ с показательным законом распределения. p ( x )  e  x Вычислим преобразование Фурье   ( t )   e      it  e   x itx e dx    e  (   it ) x dx   (   it ) x (d (  it ) x )   (   it ) x     it e Разложим     it характеристическую функцию показательного закона распределения в ряд Тейлора в точке t = 0 по переменной it. (0)  1 (t )  i ; (0)   i ;  (it )2 37 (t )   i2(it )i  2 ; (0)  2 ; (it )4 (it )3 2 В результате получим следующий отрезок ряда Тейлора    it  1 2 i it  2 (it ) 2   2!  1  1 t  1 (it ) 2 2  Коэффициент во втором члене ряда равен математическому ожиданию случайной величины, коэффициент в третьем члене равен дисперсии случайной величины. Mx  1 ; Dx  Эта связь между 1 2 ;  коэффициентами ряда Тейлора для характеристической функции и центральными моментами случайной величины является универсальной. Таким образом, ряд Тейлора для характеристической функции ψ(t) случайной величины ξ позволяет определить центральные моменты для этой случайной величины ( t )  2 (it )  (0)  a  (it )  b  2! (it )  c  3! 3  В этом ряде коэффициенты a, b, c и т.д. являются центральными моментами случайной величины ξ, в частности Mx  a; Dx  b; Техника преобразований Фурье позволяет стандартными способами вычислять обратное преобразование и тем самым по изображению Фурье находить функцию плотности распределения случайной величины. Это можно делать, используя таблицы преобразований Фурье или методы интегралов от функций комплексных переменных (теорию вычетов). Замечание. Характеристические функции могут быть определены аналогичным образом при помощи преобразований Лапласа. Разница между преобразованиями Фурье и Лапласа состоит в замене мнимой величины it 38 комплексной величиной s = σ + it, что при надлежащем выборе значения σ расширяет класс функций, для которых интегральное преобразование существует. 39 Лекция 4. ДИСКРЕТНЫЙ МАРКОВСКИЙ ПРОЦЕСС Рассмотрим процесс случайного блуждания в системе, содержащей конечное число состояний (S0, S1, … , Sn.) = {Si}. Пусть время в системе изменяется шагами и принимает значения (0, 1, … , m ) = {i}. На каждом шаге i система может с вероятностью pi остаться в текущем состоянии Sj или с вероятностью pik перейти в другое состояние Sk. На рис. 4.1 показано дерево процесса случайного блуждания в системе с четырьмя состояниями. На шаге 0 система находится в состоянии S0 . На следующем шаге 1 система может остаться в состоянии Si или перейти с вероятностью pij в одно из двух соседних состояний Si-1 или Si+1. Процесс случайного блуждания реализуется как путь по ветвям дерева возможных переходов. Каждому пути или каждой реализации процесса соответствует функция переходов fi = fi (S0, …Sk). Множество функций переходов f1, f2, … fr называется случайным процессом. S0 S1 S2 S3 S4 1 2 3 4 Рис. 4.1. Дерево возможных переходов Марковского процесса 40 Пусть система на шаге k находится в состоянии Si. Пусть вероятность перехода на шаге k+1 в состояние Sj равна pij .Эта вероятность находится как условная вероятность pij  P{S(k  1)  S j / S(k)  Si } Случайный процесс с дискретным временем и матрицей переходных вероятностей {pij} называется цепью Маркова или Марковским процессом с дискретным временем . Случайный процесс зависит от начального состояния Si(0), вероятность пребывания в котором можно заранее задать. Например, pi = P{S(0) = Si}. 4.1. Матрица переходных вероятностей Для определения Марковского случайного процесса необходимо задать матрицу переходных вероятностей. В качестве примера запишем матрицу для случайного процесса, соответствующего дереву возможных переходов на рис 4.1. 0   p00 p01 0 p 0   10 p11 p12 0 P   0 p 21 p 22 p 23 0   0 0 p32 p33 p34   0 0 p 43 p 44   Случайному процессу можно сопоставить граф состояний, отражающий возможные переходы системы за один шаг по времени. Пусть, например, система находится в состоянии S1 на шаге k. На следующем шаге k+1 система может остаться в состоянии S1 с вероятностью p11, перейти в состояние S0 с вероятностью p10 или перейти в состояние S2 с вероятностью p12. Граф состояний для системы на этом шаге по времени имеет следующий вид (рис. 4.2). 41 p11 S0 S1 p10 S3 p12 Рис.4.2. Граф состояний Марковского процесса Найдем вероятность того, что система с конечным числом состояний {Si} за n шагов перейдет из состояния Sk в состояние Sm. На нулевом шаге система находится в состоянии Si = S(0) и может на следующем шаге перейти в любое состояние Sk, переходная вероятность для которого отлична от нуля. На n-1 шаге система окажется в одном из состояний Sk = S(n-1), каждое из которых является начальным (рис. 3). S0 S1 S2 S4 • • • • • • • • • • • • • • • • • • • • • • • • n-1 S0 S1 p02 n S0 S2 p12 S1 S4 p22 p42 S2 S4 Рис. 4.3. Диаграмма переходов для Марковского процесса с дискретным временем На следующем шаге n система перейдет из одного из этих состояний Sk в состояние Sj. 42 Найдем вероятность этого перехода по формуле полной вероятности p j  P{S(n )  S j}   m {P{S(n)  S j / S(n  1)  Sk }P{S(n  1)  Sk } k 0 Здесь используется факт, что состояния системы {Si} образуют полную систему, так как она обязательно (с вероятностью единица) находится в одном из этих состояний. Из предыдущей формулы следует, что вероятности перехода из состояния Si в состояние Sj за n шагов для Марковского случайного процесса с дискретным временем находятся по следующим рекуррентным формулам p(0)  p j ; p j (n )  n  pk (n  1)  pkj (n  1); k 0 Эти формулы позволяют построить конечномерное распределение вероятностей для Марковского процесса. Другими словами найти вероятности pi нахождения системы в состоянии Si через n шагов. 4.2. Конечномерное распределение вероятностей Построим дерево возможных переходов для Марковского процесса со следующей матрицей переходных вероятностей. 0   1.000  0.500 0.500 0  P 0 0.500 0.500 0   0 0.500 0.500   0 1.000   Эта матрица имеет два поглощающих состояния S0 и S4, для которых вероятности P00 и P44 равны 1. Если система попадает в поглощающее 43 состояние, то уже не может покинуть его. Пусть в начальный момент времени система находится в состоянии S1 (рис. 4.4). S1 0.5 S2 0.5 0.5 S3 S0 0.5 S1 0.5 S0 Рис. 4. 0.5 0.5 0.5 S2 S2 S4 Дерево возможных переходов для матрицы переходных вероятностей P. Начальное состояние S1. На первом шаге система из состояния S1 может перейти в состояния S0 и S2 с вероятностями P10 = 0.5 и P12 = 0.5. Тогда конечномерное распределение вероятностей состояний системы на первом шаге будет иметь следующий вид. 1  (0.5,0,0.5,0,0) На втором шаге система может остаться в поглощающем состоянии S0 или из состояния S2 может перейти в состояния S1 или S3. Вероятности этих переходов находятся путем перемножения соответствующих переходных вероятностей на ветвях дерева. Получаем на втором шаге 2  (0.5,0.25,0,0.25,0) На третьем шаге система из состояния S1 может перейти в состояния S0 или S2, а из состояния S3 – в состояния S2 или S4. Перемножая вероятности на ветвях и складывая вероятности для одинаковых состояний, получим 3  (0.625,0,0.25,0,0.125) Продолжая подсчеты вероятностей по ветвям дерева можно получить конечномерное распределение вероятностей состояний системы на любом шаге по времени. 44 Вычисления конечномерных распределений вероятностей для Марковских процессов с дискретным временем удобно проводить при помощи матричных операций. 4.3. Матричный метод нахождения переходных вероятностей Пусть в начальный момент времени 0 система находится в состоянии Si. Этот факт может быть записан в форме вектора – строки: P (0)  ( p1 ,  p i ,  p n ), в которой одна из вероятностей pi = 1, а остальные вероятности равны нулю. На следующем шаге по времени t1 система будет находиться в состоянии Si(1) с вероятностью Pi(1) P(1)  P(0)  P, где P = {pij} – матрица переходных вероятностей для Марковского процесса, определенная выше. Аналогично, в момент времени n вероятности состояний системы будут равны P(n )  P(0)  P Рассмотрим степени n матрицы переходных вероятностей, соответствующих дереву на рис. 4.4. 0   1.000  0.500 0.250 0.250 0  2  P  0.250 0.500 0.250   0 0.250 0.250 0.500   0 1.000   45  1.000  0.625 0.250 3  P  0.250 0.250 0.250  0.125 0.250  0  0  0.125  0.250  0.625  1.000  Сравним вторые строки этих матриц, соответствующие начальному состоянию системы S1, с конечномерными распределениями вероятностей π2 и π3. Эти строки совпадают. Степени переходной матрицы дают удобный способ вычисления конечномерных распределений по мере развития Марковского процесса во времени. Так, для шага по времени n найдем матрицу Pn путем возведения в степень n исходной матрицы переходных вероятностей P. Первая строка этой матрицы даст конечномерное распределение вероятностей системы на шаге n для начального состояния S0. Вторая строка – для начального состояния S1, третья – для начального состояния S2 и так далее. Это свойство степеней переходной матрицы является утверждением фундаментальной теоремы в теории Марковских процессов. Найдем конечномерное распределение вероятностей в рассматриваемом примере системы для достаточно большого промежутка времени. Так, для n = 50 получим P 50  1.000  0.750   0.500  0.250  0  0  0.250  0.500  0.750  1.000  Результат показывает, что за конечное число шагов система переходит в стационарное состояние, в котором матрица Pn с ростом n перестает изменяться. Переходная матрица P имеет два невозвратных состояния S0 (P00 = 1) и S4 (P44 =1). С течением времени система из любого начального состояния переходит в одно из невозвратных состояний. Вероятность 46 перехода в невозвратное состояние определяются начальным состоянием системы. Так, для начальных состояний S0 и S4 система останется в этих состояниях. Для начального состояния S1 вероятность оказаться в невозвратном состоянии S0 равна 0.75, а в состоянии S4 – 0.25. Для начального состояния S2 вероятности оказаться в невозвратных состояниях S0 и S4 равны 0.5. 4.5. Марковский процесс с регулярной матрицей переходных вероятностей Рассмотрим второй пример Марковского процесса с регулярной матрицей переходных вероятностей. Переходная матрица P называется регулярной, если при некотором n матрица Pn не содержит нулевых элементов. Все состояния Марковского процесса с регулярной переходной матрицей являются достижимыми и возвратными за конечное число шагов. Рассмотрим в качестве примера Марковский процесс со следующей матрицей переходных вероятностей 0   0.500 0.500  0.300 0.400 0.300 0   P 0.300 0.400 0.300 0   0 0.300 0.400 0.300   0 0.500 0.500   После возведения этой матрицы в четвертую степень получим матрицу, не содержащую нулевых элементов  0.295  4  0.238 P  0.131  0.046  0.014  0.396 0.347 0.240 0.130 0.077 0.219 0.240 0.258 0.240 0.219 0.077 0.130 0.240 0.347 0.396 0.014  0.046  0.131  0.238  0.295  47 Таким образом, переходная матрица P является регулярной. Изучим поведение степеней регулярной переходной матрицы Pn при возрастании степени n. В эргодической теореме теории цепей Маркова с конечным числом состояний доказано, что степени переходной матрицы Pn при возрастании n стремятся к вероятностной матрице A, определяющей стационарную конечномерную функцию распределения вероятностей {P*i}, соответствующую состояниям {Si}. Так, для рассматриваемого примера, при n > 50 вероятностная матрица A имеет следующий вид AP 50 Напомним, стационарное  0.143  0.143   0.143  0.143  0.143  что 0.238 0.238 0.238 0.238 0.238 строки конечномерное 0.238 0.238 0.238 0.238 0.238 0.238 0.238 0.238 0.238 0.238 вероятностной распределение 0.143  0.143  0.143  0.143  0.143  матрицы определяют вероятностей состояний системы {Si }для различных начальных состояний. Так, начальному состоянию S0 соответствует первая строка матрицы A, начальному состоянию S1 – вторая строка и так далее. Строки матрицы равны между собой, что означает, что регулярная цепь Маркова с конечным числом состояний из разных начальных состояний сходится к одному и тому же стационарному распределению. В рассматриваемом случае стационарное распределение вероятностей имеет следующий вид Состояния Si S0 S1 S2 S3 S4 Вероятности Pi 0.143 0.238 0.238 0.238 0.143 Степени матрицы переходных вероятностей Марковского процесса с конечным числом состояний являются удобным способом анализа поведения процесса во времени. 48 Найти стационарное распределение вероятностей для регулярного Марковского процесса с конечным числом состояний можно и другим способом. По определению, если конечномерное распределение вероятностей является стационарным, то оно является решением следующей системы уравнений, записанной в матричной форме. p * * * * p1 p 2 p3  p00 p01 0 p p11 p12 0 *  10 p 4   0 p 21 p 22 p 23  0 0 p32 p33  0 0 p 43    p*  0   0*    0   p1  0    p*  2 p34   *  p p 44   *3     p4  Конечномерное стационарное распределение вероятностей позволяет определить, являются ли состояния системы {Si} возвратными или нет. Состояние Si называется возвратным, если Марковский процесс, начавшись в этом состоянии Si, через конечное число шагов вновь вернется в Si. Для возвратного состояния вероятность возвращения в него равна единице. В теории Марковских процессов доказана теорема о возвратных состояниях. Состояние Si является возвратным тогда и только тогда, когда  сумма ряда  pii (n)   n 1 Для регулярной матрицы переходных вероятностей P ответ о сходимости этого ряда очень прост. В самом деле, за конечное число шагов n Марковский процесс достигнет стационарного конечномерного распределения вероятностей. В этом случае переходная вероятность pii(n) = p*i и остается положительной константой при дальнейшем увеличении n. Это противоречит необходимому признаку сходимости ряда с 49 положительными членами и, следовательно, ряд расходится. Следовательно, конечномерный Марковский процесс с регулярной матрицей переходных вероятностей возвращается в каждое состояние через конечное число шагов. 50 Лекция 5. НЕПРЕРЫВНЫЙ МАРКОВСКИЙ ПРОЦЕСС Система с конечным числом состояний {Si} может случайно переходить из одного состояния в другое в любой момент времени t. В этом случае последовательность переходов для процесса случайного блуждания S(t0), S(t1), …S(tn), …определяется функцией S(t) и является случайным процессом с непрерывным временем. Обозначим τ = tn – tn-1 – интервал времени до очередного перехода системы в новое состояние. Заметим, что для Марковского процесса вероятность перехода из состояния S(tn-1) в следующее состояние S(tn) не зависит от предыстории процесса, то есть от последовательности переходов от состояния S(t0) до состояния S(tn-1). Из этого замечания следует, что разности τ = tn – tn-1 (n = 0, 1, … n) образуют последовательность независимых случайных величин. 5.1. Переходные вероятности непрерывного Марковского процесса Переходные вероятности pij(t) для Марковского процесса на отрезке времени от s до t определяются следующим выражением pij (t )  P{S(s  t )  S j / S(s)  Si } Необходимо задать и начальное распределение вероятностей в момент времени t = 0. pi  P{S(0)  Si } Обозначим pj(t) – вероятность перехода системы в состояние Sj за время t. Эта вероятность находится по формуле полной вероятности n p j ( t )   pi pij ( t ) i 0 Пусть система в момент времени s находится в состоянии Si. Через время t система переходит в состояние Sj. Вероятность этого перехода равна 51 pij (s  t )  n  pik (s)pkj (t) k 0 В этой формуле полной вероятности записан переход из состояния Si в состояние Sj через одно из промежуточных состояний Sk, входящее в полную систему {Sn}. В начальный момент времени t = 0 переходные вероятности  1 при i  j pij (0)    0 при i  j определяются так T=0 S1 ••• Si Si+1 • • • Sn-1 T=s S1 ••• Sk Sk+1 • • • T=s+t S1 ••• SJ-1 SJ Sn Sn-1 Sn • • • Sn-1 Sn Рис. 5.1. Диаграмма переходов из состояния Si в состояние Sj в соответствии с формулой полной вероятности Определим вид закона распределения переходных вероятностей, допустимого для Марковского процесса. Пусть процесс в некоторый случайный момент времени τ переходит из состояния Si в состояние Sj. Найдем вероятность φ(t) того, что этот случайный момент времени τ больше текущего времени t. ( t )  P{  t}, t0 Если τ > s , то система будет находиться в том же состоянии Si. 52 (s)  P{  s}, s0 Если τ > s + t, то и для этого случая система будет находиться в состоянии Si. Для этого случая вероятность φ(s + t) находится по формуле условной вероятности (s  t )  (t )  P{  s  t /   s} Переходные вероятности φ(s + t) = φ(t) по свойству Марковского случайного процесса, так как эти вероятности не зависят от времени. По формуле полной вероятности для времен t, s и s + t запишем P{  s  t}  P{  s  t /   s}  P{  s}  ( t )  (s) Из этой формулы следует следующее фундаментальное свойство переходных вероятностей Марковского процесса с непрерывным временем: при любых t и s . (s  t )  ( t )  (s) Прологарифмируем эту формулу ln((s  t ))  ln((s))  ln((t )) Из предыдущей формулы следует, что функция ln(( t )) линейна. Поэтому запишем ln( ( t ))  t t0 Постоянная λ называется плотностью (интенсивностью) перехода из состояния Si в другое состояние. Поэтому переходная вероятность для Марковского процесса распределена по показательному закону с плотностью распределения вероятности p(t). ( t )  P(  t )  p( t )  e  t Математическое ожидание для случайной величины τ равно  Mt   te  t dt  1 53 Это среднее время ожидания перехода из состояния Si в состояние Sj . Если λ = 0, то p(t) = 0 и переход системы в новое состояние невозможен. Функция распределения для показательного закона имеет следующий вид t F( t )    e  t dt  1  e  t ; Тогда обратная функция F-1(t) для показательного закона t ln( 1 F) ;  График обратной функции для показательного закона приведен на рис. 5.2. Генератор случайных чисел с показательным законом распределения, работающий по методу обратной функции, реализует следующий алгоритм: t t(ξ) ξ 1 F Рис. 5.2. График обратной функции для показательного закона распределения вероятности 54 Возьмем в качестве аргумента обратной функции F-1 случайные числа ξ, равномерно распределенные на отрезке [0, 1]. Тогда значения обратной функции t(ξ) являются случайными числами, распределенными по показательному закону. Для Марковского процесса случайное число t(ξ) является временем до ближайшего перехода системы в новое состояние. График на рис. 5.2 показывает, что случайные числа t(ξ) в основном являются малыми и группируются вблизи начала координат. Большие случайные числа t(ξ) или большие интервалы между переходами, возникают относительно редко. Определим вероятность перехода pij из состояния Si в состояние Sj за малый интервал времени Δt. По теореме о среднем эта вероятность равна P{t1  t  t1  t}   t1  t  e  t dt  t t1 Таким образом, вероятность перехода для Марковского процесса из состояния Si в состояние Sj за малое время Δt равна pij (t )  ij t Рассмотрим в качестве примера систему из четырех состояний. Матрица переходных вероятностей для нее имеет следующий вид  p11 p  21  p31 p  41 p12 p 22 p32 p 42 p13 p 23 p33 p 43 p14   0 12 13 14  p 24    21 0  23  24    t p34   31 32 0  34  p 44    41  42  43 0  Переходные вероятности pij зависят от времени t, интенсивности λij для стационарного Марковского процесса от времени не зависят. Для Марковского процесса с непрерывным временем обозначают переходные вероятности с одинаковыми индексами {pii(t)} = {pi(t)}. Это вероятности пребывания системы в каждом состоянии. 55 5.2. Уравнения Колмогорова Рассмотрим стационарный процесс случайного блуждания в системе с непрерывным временем, содержащей конечное число состояний (S0, S1, … , Sn.) = {Si}. В момент времени t + Δt система может с вероятностью pi(t) остаться в текущем состоянии Si или с вероятностью pik(t + Δt) перейти в другое состояние Sk. В этот же момент времени t + Δt система может перейти из любого другого состояния Sj с вероятностью pji(t + Δt) в состояние Si. Таким образом, в момент времени t + Δt вероятность pi нахождения системы в состоянии Si может не измениться, может увеличиться за счет перехода системы из одного из других состояний Sj или уменьшиться за счет ухода системы в одно из состояний Sk. Эти простые соображения позволяют на основе формулы полной вероятности записать дифференциальные уравнения для вероятностей {pi(t)} нахождения системы в одном из состояний {Si}. n n k 1 k i j1 j i pi ( t  t )  pi ( t )   pik (t )p k ( t )  pi ( t ) pij (t ); Вероятность pik(Δt) перехода системы из состояния Si в состояние Sk за малый интервал времени Δt для случайной величины с показательным законом распределения вероятностей равна pij (t )  ij t. Это позволяет переписать предыдущую формулу в следующем виде p i ( t  t )  p i ( t ) t n   p k ( t )ik  pi ( t ) ij ; k 1 k i j1 j i После предельного перехода при Δt → 0 получим систему дифференциальных уравнений Колмогорова. 56 dpi ( t ) dt Система n n k 1 k i j1 j i   p k ( t )ik  pi ( t ) ij ; должна быть дополнена начальным распределением вероятностей, например: p1(0) =1, p2(0) = 0, …, pn(0) = 0. В правой части системы уравнений со знаком ( + ) записана сумма произведений вероятностей pk всех состояний {S}, из которых система может перейти в данное состояние Si, умноженные на интенсивность λik. Со знаком ( – ) записан суммарный поток интенсивностей λij, выходящих из этого состояния, умноженный на вероятность рассматриваемого состояния pi. 57 Лекция 6. НЕПРЕРЫВНЫЙ МАРКОВСКИЙ ПРОЦЕСС 6.1. Пример Марковского процесса с непрерывным временем Рассмотрим систему с семью состояниями со следующей матрицей интенсивностей. 0  0 0 0 0  0  0 0 0   0 2 0  0 0  0 0 3 0  0  0 0 0 3 0   0 0 0 0 3 0  0 0 0 0 0 3  0  0 2 0 0 0 0 0  1.1 0 2 0 0 0  0   0 2.2 0 2 0 0 0    0 0 3.3 0 2 0 0   0 0 0 3.3 0 2     0 0 0 0 3.3 0 0   0 0 0 0 0 3.3 0 0 0 0 0 2 0  Случайный процесс в этой системе можно представить при помощи следующего графа, в котором λ = 2, μ = 1/0.9 = 1.1111. λ λ S0 S1 μ λ S2 2μ λ S3 3μ λ λ S4 3μ S5 3μ S6 3μ Рис. 6.1. Граф состояний системы с непрерывным временем Обозначим {Pi(t)} – вероятности пребывания системы в состояниях {Si} и запишем систему дифференциальных уравнений Колмогорова. Правые части дифференциальных уравнений по графу состояний составляются по следующему правилу: Каждое уравнение системы соответствует одному из состояний системы; 58 В правую часть уравнения входят члены, соответствующие стрелкам переходов. Если система может уйти из состояния Si по стрелке перехода в другое состояние Sj, то знак перед этим членом ( – ). Если система может попасть в состояние Si из Sj, то знак перед этим членом ( + ). dP0 dt dP1 dt dP2 dt dP3 dt dP4 dt dP5 dt dP6 dt  p0  p1;  p0  (  )p1  2p 2 ;  p1  (  2)p 2  3p3 ;  p 2  (  3)p3  3p 4 ;  p3  (  3)p 4  3p5 ;  p 4  (  3)p5  3p6  p5  3p6 ; При решении системы дифференциальных уравнений необходимо задать начальные условия. Например, в момент времени t = 0 система находится в состоянии S0 с вероятностью P0(0) = 1. Тогда P1(0) = 0, P2(0) = 0, P3(0) = 0, P4(0) = 0, P5(0) = 0, P6(0) = 0. Замечание: Начальные вероятности Pi(0) можно задавать произвольно, учитывая, что их сумма должна равняться единице. Система уравнений Колмогорова решается численно, например, средствами пакета Matlab. В результате получим следующие графики (рис. 6.2). Графики пребывания на рисунке системы в 6.2 показывают состояниях {Si} изменение Марковского вероятности процесса с непрерывным временем. Видно, что на начальном участке графиков (t < 4 c) происходит переходный процесс, а затем вероятности перестают изменяться. Это отражает общее свойство Марковских процессов с конечным числом состояний. 59 P P0 P1 P2 P5 P6 P3 P4 t Рис. 6.2. Решения системы уравнений Колмогорова 6.2. Стационарный режим для непрерывного Марковского процесса Доказано, что в системах с конечным числом состояний, в которых из одного состояния в любое другое можно перейти за конечное число шагов, существует стационарное состояние, в котором переходные вероятности становятся постоянными и не зависят от времени. Стационарное состояние достигается при длительной работе системы и не зависит от распределения переходных вероятностей в начальный момент времени. Для исследования стационарного режима необходимо в системе уравнений Колмогорова прировнять производные в левой части к нулю. В результате получим алгебраическую систему уравнений Колмогорова. 60 n n k 1 k i j1 j i  pk (t )ik  pi (t ) ij  0; Система должна быть дополнена условием, что сумма вероятностей для состояний должна равняться единице. n  pi  1; i 1 В результате система уравнений Колмогорова, дополненная этим условием, оказывается переопределенной и поэтому одно любое уравнение из системы может быть удалено. Запишем стационарный вариант системы уравнений Колмогорова для рассмотренного выше примера Марковского процесса с графом, представленным на рис. 6.1.  p0  p1  0; p0  (  )p1  2p 2  0; p1  (  2)p 2  3p3  0; p 2  (  3)p3  3p 4  0; p3  (  3)p 4  3p5  0; p 4  (  3)p5  3p6  0; p5  3p6  0; p0  p1  p 2  p3  p 4  p5  p6  1; В этой алгебраической системе одно любое уравнение, кроме последнего, может быть отброшено. Решение системы проводится методом подстановки. В результате последовательных подстановок из первых шести уравнений получим  )p ; p1  (  3 ) p ; ( p3  1 6 p5  1 54 5 ) p ; ( p2  1 2 2 ) p ; ( 4 1 () p ; p 4  18  6 1 () p ; p 6  162  61 Подставим выражения для вероятностей в последнее условие и получим. 1 (  )  1 (  )  1 (  ) )  1; p0 (1    12 (  )  16 (  )  18  54  162  2 3 4 5 6 Из этой формулы находим вероятность p0 , а затем по предыдущим формулам находим остальные вероятности. Для рассматриваемого примера получим стационарное распределение переходных вероятностей для Марковского процесса. p0 p2 p4 p6  0.1552; p1  0.2794;  0.2515; p3  0.1509;  0.0905; p5  0.0457;  0.0272; Стационарный вариант дифференциальных уравнений Колмогорова широко применяется в качестве математической модели системы массового обслуживания. 62 Лекция 7. СИСТЕМЫ МАССОВОГО ОБСЛУЖИВАНИЯ 7.1. Основные понятия и определения Многие аспекты работы вычислительных машин и сетей моделируются на основе теории систем массового обслуживания (СМО). Приведем примеры простейших СМО: Справочная телефонная служба. На телефон справочной службы поступает запрос (заявка). Телефонистка в течение определенного времени ищет ответ на запрос и отвечает абоненту (выполняет заявку). Следующая заявка может поступить после выполнения предыдущей заявки (телефонистка ждет поступления заявки) или до окончания ее выполнения (заявка отклоняется, телефон занят). Если в справочной службе работает одна телефонистка, то имеет место одноканальная СМО без очереди (СМО с потерями заявок). Если в справочной службе установлен специальный телефон, который не отключает следующего абонента, а удерживает его до окончания обслуживания предыдущего, то получаем одноканальную СМО с очередью (СМО с ожиданием и конечным или бесконечным количеством мест в очереди). Если в справочной службе установлен специальный многоканальный телефон и работают несколько телефонисток, то получаем многоканальную СМО с очередью. Многопроцессорная вычислительная система с общей памятью. Каждый процессор вычислительной системы по мере необходимости обращается к общей памяти (посылает заявку на обслуживание). Если в это время память не взаимодействует с другим процессором, то происходит запись информации в память или считывание информации из памяти (обслуживание заявки). Если память занята - то заявка помещается в очередь на обслуживание. 63 Компьютерная сеть. По линиям связи передается информация (заявки). Один из компьютеров (центр обслуживания) выполняет заявку. Реальные СМО имеют различную природу и обслуживают различные по форме и содержанию заявки. Понятно, что количество заявок, поступающих на вход реальной СМО, и время их обслуживания определяется конкретной ситуацией, которая имеет место в данный момент времени. Сравнение результатов работы реальных СМО с разной структурой (например, результатов работы кассовой службы двух универсамов, расположенных в различных городах) в условиях обслуживания реальных покупателей методически некорректно. Поэтому для анализа работы СМО создают математические модели и определяют их характеристики путем аналитического расчета или имитационного моделирования. 7.2. Структура системы массового обслуживания Рассмотрим в качестве примера структуру многоканальной СМО с очередью: Каналы обслуживания X Очередь Y ……… Рис. 7.1. Структура многоканальной СМО с очередью СМО содержит несколько каналов обслуживания, соединенных по входам параллельно. Между входом СМО X и каналами обслуживания 64 помещен блок "очередь", в который помещаются заявки, если все каналы обслуживания заняты. Поток заявок. Поток заявок X, поступающий на вход реальной СМО, определяется конкретной ситуацией и в принципе неуправляем. При построении математической модели СМО необходимо построить и модель потока заявок. В качестве модели потока заявок принимается случайный процесс (последовательность случайных величин с заданным распределением), определяющий случайные моменты времени поступления заявок. Каналы обслуживания. Заявки поступают на входы свободных каналов обслуживания. Каждая заявка обслуживается в одном канале определенное время, называемое временем обслуживания. СМО может включать один канал обслуживания (одноканальная СМО) или несколько каналов (многоканальная СМО). Поток обслуживания. Время обслуживания каждой заявки определяется потоком обслуживания, который с математической точки зрения аналогичен потоку заявок. Очередь. В структуру СМО может входить блок "очередь". Очередь характеризуется числом мест ожидания для заявок. В модели СМО в очереди запоминаются моменты времени поступления заявки. Очередь может быть конечной, тогда в случае переполнения заявки теряются (СМО с потерями) или бесконечной. Существуют СМО без очереди. Структура СМО обычно характеризуется последовательностью из четырех символов {А, В, n, m}, где А и В - соответственно поток заявок и поток обслуживания, n и m - соответственно число каналов и число мест в очереди. Состояния. Математическая модель СМО отличается от схемы, представленной на рис. 3.1. Модель СМО имеет конечное число дискретных состояний S0, S1, …, Sn+m . 65 Состояние S0 означает, что в системе нет заявок. Состояние S1 означает, что в системе обслуживается одна заявка, очередь пуста. Состояние Sn - обслуживается n заявок, очередь пуста. Состояние S n+1, - обслуживается n заявок и в очереди одна заявка. Состояние Sn+m, - обслуживается n заявок и в очереди m заявок. Система полностью заполнена заявками, поэтому очередная заявка будет отброшена. Состояния S0, S1, …, Sn+m являются дискретной случайной величиной. Каждому состоянию Si соответствует вероятность пребывания в этом состоянии Pi. Замечание. Модели СМО соответствует математическая модель Марковский процесс с непрерывным временем. 7.3. Дисциплина обслуживания В СМО с очередью существуют различные дисциплины обслуживания (порядок обслуживания) заявок, поступивших на вход: Обслуживание в порядке поступления заявок (FIFO, First In First Out), Обслуживание в порядке, обратном порядку поступления заявок (LIFO, Last In First Out). Обслуживание по приоритетам. Каждой заявке в очереди заранее присваивается приоритет (определенное число). Первой обслуживается заявка из очереди, имеющая максимальный приоритет. Приоритет бывает абсолютным или относительным. Если в очередь СМО с абсолютным приоритетом поступает заявка с приоритетом, превышающим приоритет обслуживаемой заявки, то обслуживание прекращается и начинается обслуживание поступившей заявки. Если в очередь СМО с относительным приоритетом поступает заявка с приоритетом, превышающим приоритет обслуживаемой заявки, то обслуживание поступившей заявки начинается после окончания обслуживания текущей заявки. 66 Случайный порядок обслуживания. Выбирается случайным образом одна из заявок в очереди. Обслуживание в первую очередь заявок с минимальным временем обслуживания. Это очень эффективная дисциплина обслуживания, позволяющая минимизировать среднее время пребывания заявок в системе. Однако в большинстве практических случаев время обслуживания заранее не известно. 7.4. Характеристики системы массового обслуживания В процессе моделирования необходимо определить количественные характеристики СМО, позволяющие оценить их работу. Поток заявок и поток обслуживания в модели СМО - последовательности случайных величин. Поэтому для определения характеристик СМО применяют методы теории вероятностей и математической статистики. Обычно представляют интерес следующие характеристики СМО: Pотк - вероятность потери заявки из-за занятости очереди (вероятность отказа), А - пропускная способность - среднее количество обслуженных заявок в единицу времени, kср - среднее количество занятых каналов обслуживания, rср - средняя длина очереди tож - среднее время ожидания заявкой обслуживания в очереди. 67 7.5. Потоки заявок и потоки обслуживания В аналитических моделях СМО применяют простейшие потоки, в которых интервалы между поступлениями заявок подчинены экспоненциальному или Пуассоновскому закону. В экспоненциальном потоке интервалы между заявками подчиняются закону с плотностью распределения p(t) = λe F(t) = 1 - e - λt - λt и функцией распределения . Средняя длина интервала между заявками (математическое ожидание) для экспоненциального потока равно   M( t )   tp ( t )dt   te t dt  1 /  , где λ - интенсивность потока. Дисперсия для экспоненциального потока равна  D( t )   ( t  M( t )) 2 p( t )dt  1 / 2 . Перечислим свойства экспоненциального потока заявок. Суммарный поток заявок, полученный объединением конечного числа экспоненциальных законов с различными интенсивностями λi , является экспоненциальным потоком с интенсивностью N    i ; i 1 Поток заявок, полученный разрежением экспоненциального закона, при котором заявки исключаются из потока с вероятностью p, является экспоненциальным потоком с интенсивностью λp. Важность экспоненциального потока в моделировании СМО определяется не только возможностью на его основе получать аналитические решения, но и хорошим совпадением с практикой. Многие реально наблюдаемые потоки заявок в различных СМО близки к экспоненциальному потоку. 68 7.6. Аналитическая модель системы массового обслуживания Многие важные результаты, связанные с функционированием систем массового обслуживания (СМО), получены при помощи аналитических моделей. Наиболее распространены и изучены Марковские модели, основанные на математическом аппарате Марковских случайных процессов с непрерывным временем. 7.7. Состояния и интенсивность потока событий Будем считать, что СМО имеет конечное число дискретных состояний S1, S2, S3,…, Sn , при этом процесс перехода из одного состояния в другое происходит мгновенно в любой случайный момент времени (время непрерывно). Случайный процесс перехода системы из одного состояния S i в другое состояние Sj называется Марковским, если переходная вероятность pij зависит только от текущего состояния системы и не зависит от прошлых состояний. Для построения Марковской модели вычислительной системы необходимо определить все возможные состояния системы и показать допустимые переходы из одного состояния в другое. Для этого строится граф состояний, в котором состояния обозначаются кружками, а переходы стрелками. Для каждого перехода указана интенсивность потока событий λ ij, (измеряется в 1/сек), определяющая переход системы из состояния Si в состояние Sj.  ij ( t )  lim( pij ( t , t  t )) / t; t  0 69 где pij (t, t + Δt) - вероятность перехода из состояния i в состояние j за время (t, t + Δt). Если случайный процесс является стационарным, то интенсивность потока событий не зависит от времени и определяется соотношением  ij  lim pij (t ) / t; t  0 70 Лекция 8. СИСТЕМЫ МАССОВОГО ОБСЛУЖИВАНИЯ 8.1. Одноканальная СМО с ограниченной очередью Пусть на вход одноканальной СМО (n = 1 - число каналов, m - число мест в очереди) поступает экспоненциальный поток заявок. Время обслуживания заявки в СМО - случайная величина, распределенная по экспоненциальному закону. Интенсивность потока заявок на входе СМО равна λ (1/сек). Среднее время обслуживания заявки обозначим - tобсл. Обратная величина к среднему времени обслуживания определяет интенсивность потока обслуживания μ = 1/ tобсл. Эта величина характеризует среднее число заявок, которое может обслужить система в единицу времени. Определим состояния СМО Si : S0 - система свободна, заявок нет, очереди нет; S1 - система обслуживает одну заявку, очереди нет; S2 - система обслуживает одну заявку, в очереди одна заявка; ………………………………………… S m+1 - система обслуживает одну заявку, в очереди - m заявок. λ S0 λ S1 μ λ S2 μ λ λ S3 μ S4 μ μ Рис. 8.1. Граф состояний одноканальной СМО с очередью Верхние стрелки на этом графе показывают процесс прихода заявок, нижние стрелки - процесс обслуживания заявок. 71 Запишем уравнения Колмогорова для этой системы. Для записи уравнений необходимо воспользоваться следующим правилом: Вероятность того, что система останется в состоянии S1 увеличивается за счет вероятности перехода из состояния S2 (знак + ) и уменьшается за счет ухода в состояние S1 (знак - ).  p 0  p1  0;  p1  p1  p 0  p 2  0;  p 2  p 2  p1  p 3  0;  p 3  p 3  p 2  p 4  0;  p 0  p1  p 2  p 3  p 4    1; После решения первых уравнений этой системы, получим k  p k    p0 ;   Подставим эти вероятности в последнее условие и получим p0  (  )p0  (  ) p0    (  ) 2 m 1 p0  1; Из этой формулы находим p0  p0  1/ 1    (  )    (  ) 2 m 1 ; Коэффициент загрузки и коэффициент простоя СМО. Вероятность p0 соответствует состоянию S0, а это по определению вероятность того, что система свободна и не обслуживает ни одной заявки. Тогда p0 - вероятность простоя СМО, или коэффициент простоя. Поэтому величина pзаг = 1 - p0 вероятность загрузки СМО или коэффициент загрузки. 72 8.2. Одноканальная СМО с бесконечной очередью Пусть очередь в СМО бесконечна (m = ∞). Тогда знаменатель в предыдущей формуле для p0 является геометрической прогрессией. Если (λ/μ) < 1, то ряд сходится и его сумма равна 2   ( )   ( ) S  1    Это позволяет коэффициент простоя. записать выражение m    11 ;  для вероятности p 0. Это p0  1 / S  1   ; Тогда коэффициент загрузки равен p заг  1  p0    1; Вероятность отказа. В СМО с конечной очередью существует вероятность состояния, в котором очередная пришедшая заявка не может быть обслужена. Это последнее состояние Sm+1. Поэтому вероятность отказа pотк равен вероятности пребывания СМО в этом состоянии.  pотк  p m 1  (m  1)  .  m 1 ; В СМО с бесконечной очередью вероятность отказа равна нулю, если p заг    1. Замечание. Условие pзаг < 1 является фундаментальным в теории СМО. Если принять, что λ = μ или pзаг = 1, то вероятность отказа pотк  (m  1)  1  1, превышает единицу, что делает СМО неработоспособной. Пропускная способность СМО. Для СМО с любой структурой пропускная способность А - количество обслуженных заявок, покидающих 73 систему, в единицу времени. Пропускная способность равна количеству входящих заявок в единицу времени, умноженное на вероятность того, что заявка будет обслужена A    (1  p отк ); Среднее число заявок, находящихся в СМО. В системе с конечной очередью в каждый момент времени может быть одно из следующих состояний: S0: в системе 0 заявок; S1: в системе 1 заявка в канале обслуживания и 0 заявок в очереди; S2: в системе 1 заявка в канале обслуживания и 1 заявка в очереди; S3: в системе 1 заявка в канале обслуживания и 2 заявки в очереди; …………………………………………………………………………… S m+1: в системе 1 заявка в канале обслуживания и m заявок в очереди; Номер состояния S k соответствует числу k заявок, находящихся в СМО. Вероятность этого состояния равна pk. 1. Число заявок, находящихся в СМО, является дискретной случайной величиной со следующей функцией распределения вероятностей Значение с.в. 1 2 Вероятность pk p0 p1 p2 … m+1 pn+1 Среднее число заявок в СМО nср равно математическому ожиданию этой случайной величины n ср  0  p0  1  p1  2  p 2    (m  1)  p m 1  m 1  k  pk ; k 1 или 2 n ср  1    2  (  ) 3  3  (  ) m 1  (m  1)  (  )   k  .   ; m 1 k k 1 74 Для системы с бесконечной очередью ( m → ∞) эта формула превращается в убывающую геометрическую прогрессию. В этом случае среднее число заявок в СМО равно n ср  n    /(1   ); Среднее число обслуживаемых заявок. В канале обслуживания одноканальной СМО может обслуживаться одна заявка или ни одной заявки. Этому соответствует следующая дискретная случайная величина Значение с.в. 1 Вероятность pk p0 p1 Среднее число обслуживаемых заявок nоз равно математическому ожиданию этой случайной величины n оз  0  p0  1  p1   ; Средняя длина очереди. Средняя длина очереди rcp равна среднему числу заявок в системе nср за вычетом среднего числа обслуживаемых заявок. rср  n ср  n оз  2 2  (  ) 3  3  (  ) m 1  (m  1)  (  )   k  .   ; m 1 k k2 В одноканальной СМО с бесконечной очередью средняя длина очереди равна 2 rcp  (  ) /(1   ); Средняя длина очереди резко возрастает при стремлении коэффициента загрузки pзаг к единице. pзаг    1; Для pзаг больше единицы работа СМО невозможна. 75 Замечание. Мы рассматриваем только практически значимый стационарный режим работы СМО, в котором переходные вероятности перестают зависеть от времени и система может работать сколь угодно долго. Для исследования переходного режима необходимо численно решить систему дифференциальных уравнений Колмогорова и проанализировать графики изменения переходных вероятностей во времени, как это было показано выше. Среднее время пребывания заявки в СМО. В стационарном режиме работы СМО ( pзаг <1 ) интенсивности входного и выходного потоков равны между собой и равны λ. Пусть в стационарном режиме среднее число заявок в системе в единицу времени равно nср. Тогда согласно теореме Литтла среднее время пребывания заявки в СМО - Tсист равно Tсист  n ср / ; Теорема Литтла доказана для любого распределения потока заявок и при любом распределении времени обслуживания. Для СМО с бесконечной очередью эту формулу можно переписать в следующем виде Аналогично Tсист   (11 /  ) ; можно записать формулу для среднего времени пребывания заявки в очереди Tоч. В общем случае Tоч  rср / ; Для СМО с бесконечной очередью получим Tоч   ;  (1  / ) 2 Эти формулы позволяют определить среднее время обслуживания заявки СМО - tобсл из следующего соотношения 76 Tсист  t обсл  Tоч ; 8.3. Многоканальная СМО с очередью Пусть на вход многоканальной СМО (n – число каналов, m – число мест в очереди) поступает экспоненциальный поток заявок. Длительность обслуживания заявки в каждом канале СМО - случайная величина, также распределенная по экспоненциальному закону. Интенсивность потока заявок одинакова для каждого канала и равна λ (1/сек), среднее время обслуживания заявки - tобсл. Обратная величина к среднему времени обслуживания определяет интенсивность потока обслуживания μ = 1/ tобсл, характеризующее среднее число заявок, которое может обслужить система в единицу времени. При наличии свободного канала приходящая заявка сразу отправляется на обслуживание. Если все каналы заняты, то заявка ставится в очередь. Если очередь заполнена, то заявка получает отказ. При освобождении канала и наличии заявок в очереди в освободившийся канал сразу помещается первая заявка из очереди. Обслуженные заявки покидают систему. Определим состояния Si для многоканальной СМО. S0 - система свободна, заявок нет, очереди нет; S1 - система обслуживает заявку в одном канале, очереди нет; S2 - система обслуживает заявки в двух каналах, очереди нет; ………………………………………… Sn - система обслуживает заявки в n каналах, очереди нет; Sn+1 -система обслуживает заявки в n каналах, одна заявка в очереди; ………………………………………… Sn+m - система обслуживает заявки в n каналах, в очереди находятся m заявок. Длина очереди m может быть неограниченной. 77 λ λ S0 S1 μ λ S2 2μ λ λ S3 3μ S4 4μ 4μ Рис. 8.2. Граф состояний четырехканальной СМО с очередью Граф состояний многоканальной СМО подобен графу СМО (рис. 8.1) и отличается только интенсивностями потока обслуживания (нижние стрелки на графе). При поступлении очередной заявки (если есть свободные каналы) система переходит в соседнее правое состояние. Интенсивность перехода вправо определяется интенсивностью входного потока λ. По-другому обстоит дело с потоком обслуживания. Если система находится в состоянии S1, то работает один канал и интенсивность обслуживания равна μ. В состоянии S2 параллельно работают 2 канала и поэтому интенсивность обслуживания равна 2μ. Рост интенсивности обслуживания прекращается после того, как все каналы оказываются занятыми. Для системы с n каналами максимальное значение интенсивности обслуживания равно nμ. Запишем уравнения Колмогорова для стационарного режима работы 4-х канальной СМО.  p 0  p1  0;  p1  p1  p 0  2p 2  0;  p 2  2p 2  p1  3p3  0;  p3  3p3  p 2  4p 4  0;  p 4  4p 4  p3  4p5  0;  p 0  p1  p 2  p3  p 4    1; 78 В этом примере видно, что в первых 4-х уравнениях, связанных с процессом заполнения каналов обслуживания, прослеживается одна закономерность, а в остальных уравнениях - другая. Выражая последовательно вероятности p1, p2, … , pn+m через p0, получим решение алгебраической системы уравнений Колмогорова k 1   1 pk    p0 ;    k!   1 p1    p 0 ;    1! n   1 pn    p0 ;    n! для k < n,   p n 1      n 1 1 1 n n! k p0 ;   1 p n i    k  n p 0 ;    n n!  pn  m      nm 1 m n n! p0 ; для k > n, Вероятность p0 находится из последнего уравнения после подстановки в него предыдущих формул  p0 1  .   .      2 1  n 1  n 1 1  nm 1     .  .   . m  2!  n!  nn !  n n!   1; На следующем шаге находим вероятности p1, p2, … p n+m. Проведем на основе этих выражений анализ основных вариантов функционирования многоканальной СМО. 8.4. Система без очереди Вероятность отказа. Очередная заявка, поступающая в СМО, не обслуживается, если все каналы обслуживания заняты. Это значит, что система находится в состоянии Sn. Вероятность нахождения системы в этом состоянии pn является вероятностью "отказа" в обслуживании. 79 n p отк   1  pn    p0 ;    n! Вероятность p0 соответствует состоянию S0, означающему, что в СМО нет заявок. Это коэффициент простоя p0. Величина (1 - p0) характеризует факт, что СМО занята обслуживанием заявок. Это коэффициент загрузки pзаг = 1 - p0. A    (1  p отк ); Пропускная способность A. Это количество обслуженных заявок, покидающих систему в единицу времени. Пропускная способность равна количеству входящих заявок в единицу времени, умноженное на вероятность того, что заявка будет обслужена. Среднее число заявок nср, находящихся в СМО. В текущий момент времени в системе может быть 0, 1, 2, …, n заявок, при этом все заявки обслуживаются. Вероятность того, что в системе находится k заявок, равна pk. Среднее число заявок nср равно математическому ожиданию дискретной случайной величины со следующей функцией распределения вероятностей Значение с.в. 1 2 … n Вероятность pk p0 p1 p2 … pn n n ср  0  p 0  1  p1  2  p 2    n  p n   k  p k ; k 1 Среднее число обслуживаемых заявок nоз - математическое ожидание дискретной случайной величины "количество занятых каналов" со следующей функцией распределения вероятности Значение с.в. 1 2 … n Вероятность pk p0 p1 p2 … pn n оз  0  p 0  1  p1    n  p n ; 80 В рассматриваемом случае в СМО очереди нет, nср = nоз. Средняя длина очереди равна нулю. Система работоспособна, если  n  1. 81 Лекция 9 СИСТЕМЫ МАССОВОГО ОБСЛУЖИВАНИЯ 9.1. Система с ограниченной очередью Вероятность отказа. Очередная заявка, поступающая в СМО, не обслуживается, если все каналы обслуживания и все места в очереди заняты. Это значит, что система находится в состоянии Sn+m. Вероятность нахождения системы в этом состоянии pn+m является вероятностью "отказа" в обслуживании. p отк  p n  m      nm 1 m n n! p0 ; Вероятность p0 соответствует состоянию S0, означающему, что в СМО нет заявок. Это коэффициент простоя p0. Величина (1 - p0) характеризует факт, что СМО занята обслуживанием заявок. Это коэффициент загрузки pзаг = 1 - p0. Пропускная способность A. Это количество обслуженных заявок, покидающих систему в единицу времени. Пропускная способность равна количеству входящих заявок в единицу времени, умноженное на вероятность того, что заявка будет обслужена. A    (1  p отк ); Среднее число заявок, находящихся в СМО. В текущий момент времени в системе может быть 0, 1, 2, …, n + m заявок, при этом заявки 1, … n могут обслуживаться, а n + 1, … n + m могут находиться в очереди. Вероятность того, что в системе находится k заявок, равна pk. Среднее число заявок nср равно математическому ожиданию дискретной случайной величины со следующей функцией распределения вероятности 82 Значение с.в. 1 2 … n+m Вероятность pk p0 p1 p2 … pn+m n ср  0  p0  1  p1  2  p2    (n  m)  pn  m  nm  k  pk ; k 1 Среднее число обслуживаемых заявок nоз - математическое ожидание дискретной случайной величины "количество занятых каналов" со следующей функцией распределения вероятности (в состояниях Sn, … Sn+m обслуживается n заявок) Значение с.в. 1 … n … n Вероятность pk p0 p1 … pn … pn+m n оз  0  p 0  1 p1    n (p n  p n 1    p n  m ); Средняя длина очереди rоч. - математическое ожидание дискретной случайной величины со следующей функцией распределения (в первых состояниях S0, … Sn очереди нет) Значение с.в. Вероятность pk 1 2 … m pn+1 pn+2 … pn+m rоч  1 p n 1  2  p n  2    m  p n  m ; Средняя длина очереди резко возрастает при стремлении коэффициента загрузки pзаг к единице. p заг  1  p 0  1; Для pзаг больше единицы работа СМО в стационарном режиме невозможна. 83 9.2. Система с бесконечной очередью Вероятность отказа. В бесконечной очереди могут разместиться все пришедшие заявки. Поэтому вероятность отказа равна нулю. Вероятность p0 соответствует состоянию S0, означающему, что в СМО нет заявок. Это коэффициент простоя p0. Величина p0 находится из последнего уравнения системы уравнений Колмогорова, в которое входит сумма бесконечного числа слагаемых  p0 1  .   .      2 1   1;  n 1  n 1 1  nm 1    .  .   .   m 2!  n!  nn!  n n!  В этом уравнении конечное число первых членов суммы со степенями 1,… n соответствуют n каналам обслуживания СМО. Остальные члены со степенями n + 1, …n + m, … соответствуют очереди и образуют убывающую геометрическую прогрессию со знаменателем членом a 0   / n. q  [  ] n1  1 и первым Подставим в предыдущую формулу выражение для суммы членов геометрической прогрессии и получим формулу для вычисления p0.   p0 1  .   .     2 1  n 1  .  n 1 1  . 1    1;    . 2!  n!  nn!  1[  ] n1     Эта формула позволяет найти вероятность p0. Затем находятся остальные вероятности p1, p2, … pn+m, … по формулам, приведенным выше. В системе с бесконечной очередью бесконечное число состояний Si и соответственно бесконечная последовательность вероятностей { pi }. Последовательность вероятностей { pi } для состояний, соответствующих очереди, монотонно убывает. Поэтому продолжаем вычисление этих вероятностей до тех пор, пока соответствующее значение p N окажется меньше ε, где ε - принятая погрешность. 84 Величина (1 – p0) характеризует факт, что СМО занята обслуживанием заявок. Это коэффициент загрузки pзаг = 1 – p0. Пропускная способность A. Это количество обслуженных заявок, покидающих систему в единицу времени. В системе с бесконечной очередью все поступившие заявки будут обслужены. Поэтому A  ; Среднее число обслуживаемых заявок nоз - математическое ожидание дискретной случайной величины "количество занятых каналов" со следующей функцией распределения вероятности (в состояниях Sn, … Sn+m обслуживается n заявок) Значение с.в. 1 … n n … n … Вероятность pk p0 p1 … pn pn+1 … pn+m … n оз  0  p0  1 p1  n(pn  pn 1  pn  m  ); В этой формуле в скобках записана бесконечно убывающая геометрическая прогрессия с первым членом a1 a1  p 0  nn!  (  ) и знаменателем n  1 1  / n  q   / n. В результате получим формулу для вычисления среднего числа обслуживаемых заявок в СМО.  n оз  p 0 (  )  (  )  (  ) 2 31 2  (  ) 4 1  n 1 1  n 1   ( )  ( ) ( 1 23  ( n  2)!  ( n 1)! 1  / n ; Средняя длина очереди rоч. - математическое ожидание дискретной случайной величины со следующей функцией распределения (в первых состояниях S0, … Sn очереди нет) 85 Значение с.в. Вероятность pk 1 2 … m … pn+1 pn+2 … pn+m … Очередь в рассматриваемом случае бесконечна, однако, средняя длина очереди - конечное число. rоч  1  p n 1  2  p n  2   k  p n  k   Для вычисления rоч преобразуем бесконечный ряд следующим образом rоч  2  p n 1  3  p n  2   (k  1)  p k    p n 1  p n  2   p k   Разобьем этот ряд на два. Сумму положительных членов обозначим S1, а сумму отрицательных - S2 rоч  S1  S2 . Вычислим сумму первого ряда S1. Подставим в ряд формулы для pk, получим S1 ( x )  где n  P0 n 2  n! 3 k (2  x  3  x  4  x  (k  1)  x  ), x   / n. Ряд S1(x) сходится равномерно и поэтому его сумма является непрерывной функцией от x. Этот ряд можно почленно интегрировать и дифференцировать. Вычислим интеграл от суммы ряда  S1(x)dx  n  p0 n  n! 2 3 4 (x  x  x   x k 1  ); В этой функцией в скобках стоят члены геометрической прогрессией с первым членом 2 x  ( / n) 2 и знаменателем x   / n. Запишем сумму этого ряда  S1(x)dx  n  p0 x 2 ; n  n! 1 x 86 Продифференцируем эту формулу по х и получим сумму ряда S1(x). n  p 0 2x  x 2 S1 ( x )  n ; 2  n! (1 x ) Второй ряд S2(x) - геометрическая прогрессия с первым членом и знаменателем x   / n. Сумма ряда S2 равна S2  n  p0 x n  n! 1 x  n(1 / n) ; Объединяя формулы для S1 и S2, получим rоч  n 1  ; n 1  nn !(1  / n ) Суммируя среднее число обслуживаемых заявок и среднее число заявок в очереди, получим среднее число заявок, находящихся в СМО nср. n ср  n оз  rоч . Среднее время пребывания заявки в СМО. В стационарном режиме работы СМО ( pзаг <1 ) интенсивности входного и выходного потоков равны между собой и равны λ. Пусть в стационарном режиме среднее число заявок в системе в единицу времени равно nср. Тогда согласно теореме Литтла среднее время пребывания заявки в СМО - Tсист равно Tсист  n ср / ; Аналогично можно записать формулу для среднего времени пребывания заявки в очереди Tоч. В общем случае Tоч  rср / ; Эти формулы позволяют определить среднее время обслуживания заявки СМО - tобсл из следующего соотношения Tсист  t обсл  Tоч ; 87 В заключение отметим, что простые аналитические решения получены только для простейших систем массового обслуживания с входными потоками заявок экспоненциального или Пуассоновского типов. Вопросы: 1. Как определяется состояние СМО? 2. Сколько состояний у многоканальной СМО без очереди? 3.Сколько состояний у многоканальной СМО с конечной очередью? 3. При каких соотношениях интенсивностей потока заявок и потока обслуживания одноканальная СМО работоспособна? 4. При каких соотношениях интенсивностей потока заявок и потока обслуживания многоканальная СМО работоспособна? 5. Какова пропускная способность СМО с бесконечной очередью? 88 Лекция 10. ИМИТАЦИОННОЕ МОДЕЛИРОВАНИЕ Имитационная модель является частным случаем математической модели и основана на алгоритмическом описании работы технической системы, подверженной имитационного воздействию моделирования случайных является факторов. статистический Основой эксперимент. Описание поведения технической системы во времени является основой имитационной модели. Определим основные понятия имитационного моделирования. Процесс – логически связанный набор действий, выполняемый технической системой. Считается, что работа системы разделена на отдельные шаги и развивается во времени. Событие – мгновенное изменение некоторого параметра технической системы. Заявка – это некоторое сообщение (транзакт), поступающее на вход технической системы извне. Прохождение заявки по имитационной модели соответствует процессу обслуживания заявки в СМО. Реальное и модельное время. В процессе имитационного моделирования необходимо учитывать три представления о времени: Реальное время – время, в котором функционирует реальная техническая система, Модельное (системное) время, связанное с реальным временем масштабным коэффициентом. В модельном времени организуется работа имитационной модели, Машинное время – время работы компьютера, необходимое для реализации процесса моделирования. При организации модельного времени следует обращать внимание на процессы, происходящие одновременно (параллельные процессы) и на процессы с существенно различным масштабом времени. Понятно, что 89 реализация реальных параллельных процессов в модели, реализуемой на однопроцессорном компьютере, невозможна. Поэтому такие процессы моделируются последовательно, результаты моделирования запоминаются, а потом проводится их обработка. При моделировании разномасштабных процессов необходимо шаг приращения времени соразмерять с временем протекания самого быстрого процесса. В этом случае для моделирования медленного процесса потребуется слишком много шагов и следовательно больших затрат машинного времени. 10.1. Моделирование с постоянным шагом по времени Алгоритм моделирования с постоянным шагом представлен на рис. 10.1. Начало Установка начальных состояний Продвижение модельного времени Tm = Tm +Δt Tm > TMax ДА Конец НЕТ Обработка нового состояния НЕТ Событие есть? ДА Обработка событий Рис. 10.1. Алгоритм моделирования с постоянным шагом по времени. 90 Моделирование времени с постоянным шагом применяется, если события появляются регулярно и достаточно равномерно распределены во времени. Интервал изменения времени выбирается исследователем, при этом считается, что событие наступает на границе двух соседних интервалов. Процесс в этом случае, как правило, развивается во времени. После начала алгоритма проводится начальная установка модельного времени и установка начального состояния системы. На следующем шаге модельное время Tn продвигается на величину заданного шага Δt и затем проверяется условие "конец интервала модельного времени А". Если интервал модельного времени не исчерпан, проводится переход в следующее состояние системы по времени. В этом случае очередное событие может наступить, а может не наступить. Если событие наступило, проводится его обработка, а затем проводится переход на шаг алгоритма, продвигающего модельное время. Если очередное событие не наступило, то проводится продвижение модельного времени. Работа прекращается после исчерпания интервала модельного времени. 10.2. Моделирование по особым состояниям В этом случае модельное время изменяется шагами, строго соответствующими интервалам времени до наступления очередного события. События обрабатываются в порядке их поступления, а одновременными считаются только те события, которые одновременны на самом деле. Алгоритмы моделирования по особым состояниям являются более сложными, так как приходится проводить планирование очередности наступления событий (составлять и обрабатывать списки текущих и будущих событий). Алгоритм моделирования по особым состояниям представлен на рис. 10.2. 91 Моделирование по особым состояниям является основным методом моделирования СМО. Это связано с тем, что потоки событий в типичных СМО распределены по законам, близким к экспоненциальному или Пуассоновскому, соответствующим потокам редких событий. Кроме того, в моделях элементов вычислительных систем приходится организовывать обработку параллельных (одновременных) событий, что проще делать при моделировании по особым состояниям. Рассмотрим более подробно моделирование систем массового обслуживания (СМО) по особым состояниям (рис. 10.2). Событием в СМО является продвижение транзакта (заявки на обслуживание) по системе или изменение состояния транзакта (активен, заблокирован). В текущий момент времени транзакт в СМО может находиться в одном из каналов в процессе обработки (находиться в активном состоянии) или в одной из позиций очереди на обработку (в состоянии задержки). Для моделирования СМО по особым состояниям формируются три списка транзактов (заявок на обслуживание): Список текущих событий. В этом списке находятся события, время наступления которых меньше или равно текущему модельному времени. События со временем наступления меньше текущего являются задержанными (заблокированными) в очереди. Список будущих событий. В этом списке находятся события, время наступления которых произойдет позже текущего момента модельного времени. Время наступления будущих событий в процессе моделирования легко рассчитать, так как время обработки каждого события заранее известно. 92 Начало Установка начального времени Прогноз ближайшего состояния T = Tm Tm > TMax ДА Конец НЕТ Продвижение времени Одновременные события НЕТ ДА Обработка одновременных событий Обработка одиночного события Рис. 10.2. Алгоритм моделирования по особым состояниям Список будущих событий. В этом списке находятся события, время наступления которых произойдет позже текущего момента модельного времени. Время наступления будущих событий в процессе моделирования легко рассчитать, так как время обработки каждого события заранее известно. Список прерываний. В этом списке находятся заблокированные события (транзакты), обработка которых может быть возобновлена при выполнении некоторых условий. Прерывания в обслуживании могут 93 наступать для транзактов с приоритетами. Если поступает заявка с высоким приоритетом, то обработка текущей заявки может прерываться и возобновляться в тот момент, когда ее приоритет окажется максимальным. Рассмотрим временную диаграмму работы модели одноканальной СМО с продвижением модельного времени по событиям (рис. 10.3). Предположим, что канал обслуживания СМО свободен. Первая заявка (транзакт) 1 приходит на вход СМО в момент времени t0 и начинает обслуживаться в течение времени t1 обсл . В течение этого времени канал обслуживания занят. Следующая заявка 2 приходит позже, когда канал обслуживания уже свободен. Эта заявка начинает обслуживаться в момент времени t2 и заканчивает обслуживание в момент времени t4. Следующая заявка 3 приходит в момент времени t3, когда канал обслуживания занят, и ждет в очереди до момента времени t4. Эта заявка начинает обслуживаться в момент времени t4 и заканчивает – в момент времени t6. Очередная заявка 4 приходит в момент времени t5, ждет в очереди до времени t6 и обслуживается до момента времени t7 Событиями в системе являются или поступление очередной заявки (транзакта), или изменение состояния заявки (транзакта) (обслуживаемый транзакт активен, транзакт в очереди заблокирован). Моменты наступления событий отмечены на оси модельного времени t. Для этого примера в список текущих событий для моментов времени t0, … tn, отмеченных на оси модельного времени t, помещаем только одно событие, происходящее в этот момент времени (рассматривается одноканальная СМО с очередью, параллельных процессов обслуживания нет). 94 Заявка 4 Заявка 3 Заявка 2 Заявка 1 t0 t1 t 2 t3 t4 t5 t6 t7 t0 – время прихода первой заявки; t1 – время окончания обслуживания первой заявки; t2 – время прихода второй заявки; t3 – время прихода третьей заявки; t4 – время окончания обслуживания второй заявки; t5 – время прихода четвертой заявки; t6 – время окончания обслуживания третьей заявки; t7 – время окончания обслуживания четвертой заявки. t1 – t0 – время обслуживания первой заявки; t2 – t1 – время простоя СМО, заявок нет; t4 – t2 – время обслуживания второй заявки; t4 – t3 – время ожидания в очереди третьей заявки; t6 – t4 – время обслуживания третьей заявки; t6 - t5 – время ожидания в очереди четвертой заявки; t7 – t6 – время обслуживания четвертой заявки. Рис. 10.3. Временная диаграмма работы одноканальной СМО с продвижением модельного времени по событиям В список будущих событий для момента времени ti помещаются события, находящиеся правее его на оси модельного времени в порядке их наступления (слева направо). Понятно, что при моделировании, в отличие от реальной СМО, моменты наступления всех будущих событий могут быть рассчитаны. Глубина списка будущих событий определяется из алгоритма моделирования и обычно включает одно ближайшее событие. Рассмотрим временную диаграмму работы модели трехканальной СМО с продвижением модельного времени по событиям. В этом случае обслуживание транзактов проводится параллельно и возможно появление одновременных событий. 95 Заявка 6 Заявка 5 Заявка 4 Заявка 3 Заявка 2 Заявка 1 Канал 3 Канал 2 Канал 1 Канал 3 Канал 2 Канал 1 t 0 t1 t 2 t 3 t4 t 5 t6 t7 t 8 t9 t10 t11 t0 – время прихода первой заявки в первый канал; t1 – время прихода второй заявки во второй канал; t2 – время окончания обслуживания первой заявки; t3 – время прихода третьей заявки в третий канал; t4 – время прихода четвертой заявки в первый канал; t5 – время прихода пятой заявки в СМО; t6 – время окончания обслуживания второй заявки и начало обслуживания пятой заявки во втором канале; t7 – время прихода шестой заявки в СМО; t8 - время окончания обслуживания третьей заявки в третьем канале и начало обслуживания шестой заявки;. t9 – время окончания обслуживания четвертой заявки; t10 – время окончания обслуживания пятой заявки; t11 – время окончания обслуживания шестой заявки; t6 – t5 – время ожидания пятой заявки в очереди; t8 – t7 – время ожидания шестой заявки в очереди; Рис. 10.4. Временная диаграмма работы трехканальной СМО с продвижением модельного времени по событиям Предположим, что каналы обслуживания СМО свободен. Первая заявка (транзакт) приходит на вход СМО в момент времени t0 и начинает обслуживаться в первом канале до момента времени t2. Первый канал остается свободным до момента времени t4. Вторая заявка попадает в свободный второй канал в момент времени t1 и обслуживается в нем до момента времени t6. 96 Третья заявка попадает в свободный третий канал, а четвертая заявка – в освободившийся первый канал. К моменту поступления пятой заявки t5 все каналы обслуживания заняты. Поэтому пятая заявка ждет в очереди до момента освобождения второго канала t6. Шестая заявка также попадает в очередь. 97 Лекция 11. ИМИТАЦИОННОЕ МОДЕЛИРОВАНИЕ 11.1. Случайные и псевдослучайные числа Основным инструментом имитационного моделирования являются генераторы случайных чисел. Потоки случайных чисел имитируют окружающую среду и создают модель внешних воздействий на объект моделирования. Понятия случайных событий и случайных чисел определяется в Теории вероятностей. Это чисто математические, абстрактные понятия. Предполагается, что проводится некий эксперимент, результат которого (событие) может произойти или не произойти. Случайные числа являются случайными событиями с точки зрения теории вероятности. Так, генератор случайных чисел может выдать очередное число 2.5, а может выдать и другое число из множества допустимых чисел. Теория вероятностей изучает количественные характеристики случайных чисел. Поэтому в определение случайного числа входит закон распределения вероятности появления этого значения числа на некотором интервале. Понятие случайное число и закон распределения вероятности появления этого числа эквивалентны и неразрывны. Например, равномерному закону распределения вероятности на интервале (a, b) соответствует поток случайных чисел, отличных от нуля только на этом интервале и имеющих одинаковую вероятность появления. В математическом понятии «поток случайных чисел» предполагается, что ни одно из чисел не повторяется дважды. Создать аппаратный генератор случайных чисел можно, например, фиксируя время между приходом заряженных частиц или гамма – квантов из космоса при помощи физических счетчиков. Физический генератор случайных чисел устроен просто и 98 содержит счетчик Гейгера и усилитель. Однако закон распределения полученных таким образом случайных чисел не известен, и его можно оценить только в статистическом эксперименте. Задать произвольный закон распределения, необходимый для построения имитационной модели, в физическом генераторе случайных чисел нельзя. Существует большое число алгоритмов, позволяющих формировать последовательности случайных чисел на компьютере. Обычно программным образом формируется последовательность равномерно распределенных случайных чисел, а затем эта последовательность преобразуется в другую последовательность с заданным законом распределения. Строго говоря, получаемые программным путем последовательности чисел являются псевдослучайными, так как при сохранении в алгоритме начальных значений (параметров и "затравок") повторно генерируется та же последовательность чисел. 11.2. Период последовательности псевдослучайных чисел. Программные генераторы случайных чисел реализуются на компьютерах, обладающих ограниченной разрядной сеткой. Это приводит к тому, что в компьютере можно использовать большое, но конечное число различных рациональных чисел. Последовательность случайных чисел является выборкой из конечного числа рациональных чисел и поэтому может содержать ограниченное количество различных чисел. Поэтому любая последовательность случайных чисел, программно реализуемая компьютером, будет повторяться, начиная с некоторого числа, другими словами, будет последовательностью псевдослучайных чисел. Периодом генератора случайных чисел называют количество неповторяющихся чисел в псевдослучайной последовательности {xi}. Для 99 имитационного моделирования важно, что бы псевдослучайная последовательность была бы достаточно длинной. Период генератора случайных чисел определяется экспериментально в ходе вычислительного последовательность эксперимента. псевдослучайных Следует чисел имеет учитывать, что апериодическую (начальную) и периодическую части. Длина апериодической части заранее не известна. Поэтому при определении периода генератора случайных чисел необходимо зафиксировать одно из целых чисел, полученных в результате работы процедуры генератора, и подсчитать количество вызовов процедуры генератора до того момента, когда в результате получится такое же целое число. В процессе эксперимента по определению периода последовательности псевдослучайных чисел необходимо постепенно увеличивать номер числа, для которого подсчитывается число элементов последовательности до повторного появления этого числа. Начиная с некоторого момента, число промежуточных чисел в последовательности до его повторения перестанет изменяться. Это число и принимается за период последовательности. 11.3. Оценка качества последовательности случайных чисел Оценка качества последовательности случайных чисел проводится методами математической статистики. Это стандартная задача проверки статистической гипотезы о виде функции распределения вероятности появления случайного числа. Пусть генератор создает последовательность псевдослучайных чисел с равномерным на интервале [0, 1) законом распределения вероятности (статистическая гипотеза). Определим вероятность этого предположения. Для проверки гипотезы необходимо выполнить следующие действия: 1. разбить интервал [0, 1) на m частей ( m = 20 – 50 ); 100 2. определить число vj случайных чисел ui , попадающих в каждый интервал; 3. определить теоретическое число wj равномерно распределенных случайных чисел, которые должны попасть в эти интервалы; 4. вычислить значение критерия Пирсона χ2; 5. по таблице распределения χ2 найти вероятность правильности статистической гипотезы. Рассмотрим пример решения задачи по проверке статистической гипотезы. Допустим, что на интервале [0, 1) последовательность чисел является случайной последовательностью с равномерной функцией распределения. Разобьем интервал 0 ≤ x < 1 на m = 5 частей (Это демонстрационный пример, в последовательности всего N = 25 чисел). Построим гистограмму попаданий случайных чисел экспериментальной и теоретической последовательности в эти интервалы. v v2 = 7 v1 = 5 v4 = 5 v5 = 5 w1 – w5 = 5 v3 = 3 x 0.2 0.4 0.6 0.8 1 v1 – v5 – экспериментальная гистограмма, w1 – w5 – теоретическая гистограмма. Рис. 11.1. Гистограммы попаданий случайных чисел в интервалы График функции равномерного распределения плотности вероятности на интервале [0, 1) показан на рис. 11.2. Функция распределения отлична от нуля только на интервале [0, 1). 101 p(x) 1 x 0.2 0.4 0.6 0.8 1 Рис. 11.2. Функция плотности вероятности для последовательности равномерно распределенных случайных чисел Функция распределения вероятностей связана с функцией плотности распределения вероятностей следующим соотношением. x x F( x )   p( t )dt   1  dt  x; График функции равномерного распределения вероятности представлен на следующем рис. 11.3. F(x) 1 0.8 0.6 0.4 0.2 x 0.2 0.4 0.6 0.8 1 Рис. 11.3. Функция распределения вероятности для последовательности равномерно распределенных случайных чисел. Вероятность попадания случайной величины в каждый интервал на оси x равна 102 Pi  Fi 1  Fi  0.2; Тогда число случайных чисел последовательности длиной N, которые теоретически должны попасть в каждый интервал, равно wi  N  Pi  25  0.2  5; Гистограмма теоретического распределения числа попаданий равномерно распределенных случайных чисел в интервалы показана на рисунке 1 пунктиром. Понятно, что чем ближе теоретическая и экспериментальная гистограммы, тем ближе распределение вероятности случайных чисел экспериментальной и теоретической последовательности. Это простое соображение положено в основу критерия согласия Пирсона (критерия χ2 ). Количественную оценку близости теоретического и экспериментального распределения вероятности случайных величин по критерию Пирсона находят следующим образом. Вычисляем значение χ2 по следующей формуле m   2 i 1 ( vi  w i ) wi Для нашего примера значение χ2 2 ; = 1.6; Дальше необходимо взять таблицу распределения χ2 и по ней определить вероятность выдвинутой статистической гипотезы. В эту таблицу, кроме значения χ2 входит второй параметр – число степеней свободы n. Число степеней свободы n на единицу меньше числа интервалов m на гистограмме. n  m  1; Для рассматриваемого примера число степеней свободы n = 4. В результате по таблице находим, что вероятность выдвинутой статистической гипотезы равна 0.81 (81%). Как относиться к этому результату? Понятно, что если вероятность правильности статистической 103 гипотезы мала (обычно меньше 70%), то результат нельзя считать хорошим. Аналогично обстоит дело и с очень высокой вероятностью правильности статистической гипотезы (обычно больше 95%). В этом случае, возможно, что последовательность чисел не является случайной. Проверку статистической гипотезы по критерию χ2 удобно проводить средствами программы Microsoft Excel. Для этого необходимо загрузить числа из экспериментальной и теоретической гистограмм в ячейки программы Excel и вызвать функцию ХИ2ТЕСТ. Для нашего примера получим: Рассмотренная схема проверки статистических гипотез работает и для случайных величин Аналогичным образом с произвольным строится распределением экспериментальная и вероятности. теоретическая гистограммы и вычисляется значение χ2 . При построении теоретической гистограммы число случайных величин wi , которые должны попасть в каждый интервал, вычисляется по формуле. wi  N  Pi  25  0.2  5; Вероятность попадания числа в каждый экспериментальный интервал ΔPi определяется по следующей формуле Pi  Fi 1  Fi  0.2; В общем случае эти вероятности разные. При выполнении лабораторной работы число случайных чисел в последовательности N имеет порядок 15 – 20 тыс. Число интервалов m на гистограмме необходимо выбирать порядка 20 – 50. Удобнее всего выбирать 104 m из ряда чисел 21, 31, 41 и т.д., что позволяет получить удобные значения для числа степеней свободы n = 20, 30, 40 и т.д. Эти числа обычно указываются в таблице распределения χ2 . 11.4. Генерация равномерно распределенных случайных чисел 11.4.1. Метод мультипликативного датчика. В этом методе очередное случайное число xi+1 с равномерным законом распределения получается как остаток от деления двух больших целых чисел (деление чисел по модулю M). Результат деления по модулю определяется следующей формулой x  (axi  b)  mod M ; i1 где a, b, M - заданные целые числа,. xi – предыдущее целое случайное число. Слагаемое b используется, чтобы устранить потенциальную опасность получения непрерывной последовательности нулей, начиная с определенного шага. Если M и a — простые числа, то такой опасности не существует. При работе мультипликативного датчика вырабатывается случайная последовательность целых чисел. Последовательность случайных чисел, принадлежащих интервалу [0,1), получается путем нормировки целочисленной случайной последовательности. ui  xi / M ; Понятно, что для начала работы генератора случайных чисел необходимо задать затравку x0 (целое число). 105 11.4.2. Метод середин квадратов В этом методе также генерируется последовательность равномерно распределенных псевдослучайных чисел. Генератор, реализующий метод середин квадратов, хорошо работает на разрядных сетках в 64 или 128 из двоичных чисел. На первом шаге алгоритма задается произвольное число "затравка" x0 половинной разрядности (32 или 64 двоичных разрядов). На каждом шаге алгоритма целое число xi возводится в квадрат, при этом разрядность произведения увеличивается в два раза. Из числа – произведения выделяются средние разряды, которые и принимаются за следующее случайное число разрядностью 32 или 64. Числа, входящие в случайную последовательность целых чисел, делятся (нормируются) на самое большое целое число, которое может появиться в случайной последовательности. Это или N  232 или N  264. 1 2 В результате получаем последовательность псевдослучайных чисел с равномерным законом распределения на интервале [0,1). Рассмотрим алгоритм метода середин квадратов на примере для калькулятора. Зададимся 8-разрядной десятичной сеткой, пусть начальное число x0 = 4824. Будем генерировать случайные числа в формате 0.хххх, то есть с четырьмя десятичными знаками после запятой. Получим следующую последовательность чисел x0=4824; x02=23270976; u0=0.4824; x1=2709; x12=07338681; u1=0.2709; x2=3386; x22=11464996; u2=0.3386; x3=4649; x22 = 21613201; u3=0.4649; ... 106 Нельзя априори утверждать, что полученная последовательность чисел ui обладает необходимыми свойствами. Необходимо проверить статистическую гипотезу о том, что последовательность u0 является последовательностью случайных чисел с равномерным распределением. Недостатком метода середин квадратов является возможность получить на очередном шаге алгоритма нулевого значения ui. Тогда и все последующие числа в последовательности также будут нулями. (Например: 56002034  0020  00000400  0004  00000016  0000 Поэтому при программной реализации метода середин квадратов необходимо проводить "проверку на ноль". Если в последовательности {ui} появляется нулевой элемент, то алгоритм необходимо выполнять сначала, задав новую затравку. 11.5. Генерация произвольно распределенных случайных чисел В качестве основы при построении генераторов случайных чисел с заданным законом распределения берется хорошо работающий генератор случайных чисел с равномерным законом распределения. На практике широко применяются два метода получения последовательности случайных чисел с заданным законом распределения: метод обратной функции и метод режекции. 11.5.1. Метод обратной функции Пусть функция распределения вероятности F(x) задана аналитически. Метод обратной функции основан на взаимно – однозначном соответствии значений функции F(x) точкам случайной переменной, равномерно распределенной на отрезке [0, 1]. 107 Рассмотрим в качестве примера генератор псевдослучайных чисел с экспоненциальным законом распределения. F(x)  1  e  x ; График этой функции имеет следующий вид F(x) 1 0.8 F(t) = 1 - e – λt λ=1 0.6 0.4 0.2 2 4 6 t Рис. 11.4. Функция распределения вероятности для экспоненциального закона Для получения обратной функции аналитически преобразуем формулу закона распределения, положив в качестве независимой переменной F. Получим x  (ln(1  F)) / ; Независимая переменная F имеет смысл вероятности. Она задана на интервале 0 ≤ F < 1. В методе обратной функции в качестве независимой переменной F принимаются случайные числа с равномерным законом распределения на интервале [0, 1). Таким образом, для получения последовательности случайных чисел с экспоненциальным законом распределения необходимо взять генератор случайных чисел с равномерным законом распределения и преобразовать случайные числа в числа, распределенные по заданному закону. 108 График обратной функции представлен на следующем рисунке 10.5. x 5 4 3 x = (ln(1 – F))/λ λ=1 2 1 0.5 1 F Рис. 11.5. График обратной функции для экспоненциального закона распределения вероятности. В качестве аргумента F в обратной функции используется поток случайных чисел с равномерным распределением на интервале [1,0). На второй оси x отображаются случайные числа, распределенные по экспоненциальному закону. 11.5.2. Метод исключения (режекции) Для реализации метода режекции необходимо знать функцию распределения плотности вероятностей случайной величины p(x). Случайная величина должна быть задана на ограниченном отрезке [a, b]. . Пусть функция плотности распределения случайной величины p(x) задана графически (рис. 10.6). Площадь фигуры плотности распределения вероятности по определению равняется единице. Это свойство позволяет найти значение параметра h на рисунке 10.6. Для этого найдем площадь фигуры на рис. 10.6 и приравняем ее единице. S  h(3  1) / 2  h(4  3)  1; Из этой формулы найдем величину h h  1/ 2; 109 Генератор РРСЧ [0, h] y h h/2 x 1 2 3 4 Генератор РРСЧ [1, 4] Рис. 11.6. Функция плотности распределения случайной величины, заданная графически. Для создания генератора случайных чисел, работающих по методу режекции, необходимо взять два независимых генератора равномерно распределенных случайных чисел. Первый генератор РРСЧ создает случайные числа {xi } на отрезке [1, 4]. Второй генератор РРСЧ создает случайные числа { yi } на отрезке [0, h]. На рис. 11.6 представлен график функции плотности распределения вероятности, подготовленной для реализации метода исключения. Далее с этими генераторами проводится вычислительный эксперимент. Генераторы вырабатывают два случайных числа, после чего проверяется, пересекаются ли линии, соответствующие этим числам, внутри фигуры функции плотности распределения или вне ее (рис.11.6). 110 На этом рисунке линии пересекаются внутри фигуры. Поэтому случайное число xi принимается за результат работы генератора случайных чисел с функцией распределения плотности вероятностей, показанной на рис 11.6. Если линии не пересекаются, то очередное число xi отбрасывается. 111 Лекция 12. ИМИТАЦИОННОЕ МОДЕЛИРОВАНИЕ 12.1. Списки текущих и будущих событий Во время имитационного моделирования по событиям время прихода заявок (транзактов) и время обработки транзактов (заявок) определяется генераторами псевдослучайных чисел (ГСЧ), реализованных в программе имитационного моделирования. Если при очередном запуске ГСЧ не изменить начальные установки, то генератор повторит ранее созданную последовательность псевдослучайных чисел. Этим фактом широко пользуются при проведении имитационного моделирования. Значительно проще еще раз запустить ГСЧ, чем сохранять огромные массивы сгенерированных псевдослучайных чисел. Основой алгоритма имитационного моделирования является создание списков текущих и будущих событий. Анализ диаграмм на рис. 10.3 и 10.4 позволяет установить алгоритм. В самом деле, зафиксируем любое событие на оси времени. Это текущее событие. События, расположенные справа от текущего на оси времени – будущие события. Ближайшее будущее событие к текущему событию наступит первым и поэтому это событие необходимо взять в качестве текущего на следующем шаге алгоритма. Если в качестве модели СМО используется Марковский процесс, то в процессе моделирования достаточно вычислять и помнить только одно следующее по времени будущее событие. Это следует из свойств Марковского процесса. Пример списка текущих и будущих событий представлен в следующей таблице. 112 Пример списка текущих и будущих событий СМО. Модельное время ti t0 t1 t2 t3 t4 t5 t6 t7 t8 t9 t10 t11 Список текущих событий Zi Z1 Z2 Z3 Z4 Z5 Z6 Z7 Z8 Z9 Z10 Z11 Z12 Список будущих событий Zi+1 Z2, Z3, Z4 Z3, Z4, Z5 Z4, Z5, Z6 Z5, Z6, Z7 Z6, Z7, Z8 Z7, Z8, Z9 Z8, Z9, Z10 Z9, Z10, Z11 Z10, Z11, Z12 Z11, Z12 Z12 — 12.2. Моделирование СМО как Марковского процесса Вероятность состояний СМО можно найти в процессе численного моделирования Марковского процесса. Рассматриваемый далее алгоритм может применяться для любых Марковских процессов с дискретными состояниями и непрерывным временем. Моделирование возможно, если случайное время перехода системы из одного состояния в другое распределено по экспоненциальному закону. В результате численного моделирования можно получить только вероятности состояний СМО. Остальные характеристики СМО определяются дополнительно по формулам, рассмотренным выше. В общем случае моделирования Марковского процесса размерность квадратной матрицы равна количеству состояний m и она не является трехдиагональной. Графу состояний трехканальной СМО можно сопоставить следующую матрицу интенсивностей переходов Марковского процесса. 113 0  0 0  0  0 Λ   ij  0 2 0  0 0 3 0 Рис. 12.1. Матрица интенсивностей переходов для СМО. Будем для определенности считать, что строки и столбцы матрицы нумеруются с 0. Тогда ij — интенсивность перехода в состояние j при условии, что система находится в состоянии i. Логично задать начальное состояние системы i=0 (нет заявок). Тогда, если i — номер текущего состояния, алгоритм моделирования произвольного Марковского процесса работает следующим образом: 1. Установить текущее состояние i = 0, обнулить счетчик модельного времени T и счетчики времени пребывания в k-том состоянии Tk, k = 0, …m-1, а также счетчик числа обслуженных заявок N. 2. Для каждого состояния j = 0, … m – 1 генерируется случайная величина – время до перехода из текущего состояния i в состояние j. Это время рассчитывается по формуле t ij   1 ln u  ij где u – случайное число с равномерным законом распределения на интервале [0, 1). Это генератор случайных чисел с экспоненциальным законом распределения с интенсивностью ij , построенный по методу обратной функции. 3. Среди полученных времен tij находим минимальное. Другими словами находим индекс j такой, что для него время tij минимально при изменении 0 ≤ j ≤ m-1. Найденное минимальное время определяет действительное время пребывания системы в текущем состоянии i и показывает новое действительное состояние системы j. 114 4. Увеличиваем счетчик модельного времени T и счетчик времени пребывания системы в состоянии i на величину tij . T = T + tij: 5. Ti = Ti + tij; Если j < i, то увеличить число в счетчике обслуженных заявок N на единицу. 6. Если N меньше заданного числа обслуженных заявок, то перейти в следующее состояние j. ( i = j ), а затем перейти к п.2. Если N равно заданному числу заявок – СТОП. Алгоритм определяет новое состояние, время перехода в которое минимально, и определяет, находится ли новое состояние выше или ниже главной диагонали. Для первого случая принимается решение, что пришла очередная заявка, для второго – обслуженная заявка покинула систему. На каждом шаге алгоритма расчет времени до ближайшего перехода производится заново. Это допустимо в рамках модели Марковского процесса, в которой условная вероятность перехода из одного состояния в другое равно безусловной вероятности. Предысторию в Марковском процессе учитывать не нужно. Вероятности состояний рассчитываются как отношение времени, в течение которого система находилась в данном состоянии, к общему времени моделирования, т.е. pi  Ti T Для случаев с конечным количеством состояний m накапливать модельное время необязательно — ведь T равно сумме всех Ti. Очевидно, что для случая СМО возможен переход процесса из текущего состояния i только в два соседних состояния: j = i + 1 и j = i – 1. Это позволяет упростить алгоритм моделирования. Рассмотрим два более простых алгоритма моделирования на примере СМО с бесконечной очередью. 115 При численном моделировании СМО с бесконечной очередью возникает ряд трудностей, связанных с бесконечной размерностью матрицы интенсивностей переходов. Рекомендация 1 Заметим, что матрица интенсивностей переходов для СМО – трехдиагональная (рис. 11.1). В такой матрице все отличные от нуля элементы размещены на двух диагоналях, прилегающих к главной. Таким образом, если j = i + 1, то переход происходит к элементу j, лежащему выше главной диагонали (интенсивность равна ). Если j = i – 1, то переход происходит к элементу j , лежащему ниже главной диагонали. Эти простые условия позволяют описать матрицу интенсивностей переходов условным оператором языка программирования. Необходимо вычислять случайное время перехода из текущего состояния i только в два соседних состояния: j = i + 1 и j = i – 1. Переходы в другие состояния невозможны. Необходимо подсчитывать время пребывания системы только в тех состояниях, которые достигались в действительности. Суммарное время пребывания системы во всех действительных состояниях равно модельному времени. Рекомендация 2. Необходимо ввести искусственные ограничения на длину очереди и решить конечную задачу, подобрав в эксперименте достаточно большую (например, 40 или 100) максимальную длину очереди. В процессе моделирования необходимо проверять, что бы очередь ограниченной длины полностью не заполнялась. Тогда результаты моделирования СМО с конечной очередью будут полностью соответствовать системе с бесконечной 116 очередью. Если последнее состояние в очереди достигается, необходимо увеличить ее длину. Характеристики СМО, полученные в результате моделирования Марковского процесса необходимо сравнить с характеристиками, полученными в результате аналитического расчета. 12.3. Моделирование СМО по событиям Компьютерное моделирование СМО на основе Марковского процесса возможно, если потоки заявок и потоки обслуживания распределены по экспоненциальному закону. В общем случае для произвольных распределений этот метод не применим. Для моделирования СМО с произвольными распределениями потоков заявок и потоков обслуживания применяется метод имитационного моделирования по событиям. Событиями в СМО называются моменты времени, когда очередная заявка приходит в систему и когда обслуженная заявка покидает систему. Другими словами, это моменты изменения состояния системы. Эти моменты образуют список событий, в котором можно выделить текущее событие и списки прошедших и будущих событий. Метод имитационного моделирования включает последовательное формирование списка событий и продвижение процесса моделирования от одного события к другому. Модельное время при имитационном моделировании является зависимой переменной и изменяется скачками, равными интервалам времени между последовательными событиями. Схема алгоритма имитационного моделирования по событиям приведена на рис. 12.2. 117 Принятые обозначения: Входные параметры: N — максимальное число входящих заявок (условие окончания моделирования). NK — количество каналов. LMAX — максимально допустимая длина очереди. Накапливаемая статистика: T — текущее модельное время, изменяющееся скачком между моментами возникновения событий в системе. TA — момент прихода очередной заявки (время, прошедшее с начала моделирования). K — счетчик пришедших заявок. KS — счетчик обслуженных заявок. L — текущая длина очереди. LOS — счетчик отказов (заявок, поучивших отказ в обслуживании). Массивы: OCP[i] — признак занятости i-го канала (0 — канал свободен, 1 — канал занят). TD[i] — ожидаемый момент выхода заявки из i-го канала (время, прошедшее с начала моделирования). TOS[i] — счетчик суммарного времени занятости i-го канала — сколько единиц модельного времени в течение всего процесса моделирования канал был занят. TL[M] — суммарное время пребывания системы в состоянии, когда в ней ровно M заявок. Вспомогательные величины: MIN — ближайший момент выхода обслуженной заявки из канала, считая от текущего модельного времени. 118 S — номер канала, который в текущем состоянии системы освободится первым (в момент времени MIN). DTA — время между приходами заявок (генерируемая в процессе моделирования случайная величина). DTS — время обслуживания заявки в канале (генерируемая в процессе моделирования случайная величина). Для рассматриваемой СМО 1 DTA   ln u1 ,  DTS   1  ln u 2  t обсл. ln u 2 где u1, u2 — равномерно распределенные на интервале [0, 1) случайные величины (рекомендуется использовать для их генерации два независимых датчика равномерно распределенных случайных чисел на интервале [0,1) ). В процессе моделирования в любой момент времени известны моменты наступления событий для всех возможных источников (момент прихода очередной заявки и моменты освобождения каждого из каналов). Если канал свободен, то считается (с целью упрощения алгоритма), что ожидаемое время его освобождения бесконечно велико; в предлагаемой схеме алгоритма в качестве заведомо большого значения взято число 999999. Установка начальных значений (блок 2) состоит в заполнении массива TD большими числами (999999), а всех переменных и других массивов — нулевыми значениями. Это соответствует ситуации, когда все каналы свободны и модельное время равно 0. Блок 3 определяет момент прихода очередной заявки, при этом увеличивается счетчик пришедших заявок K, модельное время не изменяется. Блок 4 — проверка условия окончания моделирования. Блоки 7-11 реализуют определение номера канала, который освобождается первым, и момента его освобождения. При начальном запуске системы все каналы свободны, поэтому номер канала S не важен (вообще-то 119 это будет последний канал), а в качестве момента его освобождения MIN будет обнаружено заведомо большое число. Далее происходит проверка: какое из событий наступит раньше, приход заявки в момент времени TA или освобождение канала в момент времени MIN (блок 12). При начальном запуске ближайшим событием окажется приход заявки (так как все каналы свободны) и произойдет переход на блок 19. Блок 19 осуществляет накопление статистики по пребыванию системы в состоянии с M заявками. Очевидно, что M равно сумме количества занятых каналов (суммы всех элементов массива OCP[i]) и количества заявок в очереди L. Если текущее модельное время T, а переход в новое состояние произойдет в момент времени TA по приходу заявки, то в состоянии М система пробудет (TA+T) единиц модельного времени. Блок 20 изменяет модельное время. Теперь состояние системы можно сформулировать как "очередная заявка только что пришла в систему". Блоки 21-24 организуют цикл поиска свободного канала. Если свободный канал найден, то происходит переход на блок 25 и досрочный выход из цикла поиска. Блоки 25-26 реализуют действия алгоритма моделирования в случае, когда в системе есть свободный канал. В этом случае канал объявляется занятым (блок 25) и сразу же разыгрывается время обслуживания заявки в канале DTS и момент его освобождения TD[i], равный текущему модельному времени T плюс DTS. Сразу же в массиве TOS[i] учитывается и время, которое занимаемый канал потратит на обслуживание заявки. Если в системе нет свободных каналов, то заявка помещается либо в очередь, либо, если очередь уже имеет максимально допустимую длину, получает отказ. Эта логика реализуется блоками 27-29. 120 После обработки прихода очередной заявки в блоках 19-29 происходит переход на блок 3, где происходит определение момента прихода очередной заявки. 121 Рис. 12.2. Схема алгоритма моделирования многоканальной СМО по событиям 122 Рассмотрим теперь вариант, когда текущее состояние системы таково, что приход очередной заявки случится позже, чем освобождение одного из каналов. При этом в результате выполнения поиска первого освобождающегося канала в блоках 7-11 будет определен номер S канала, освобождающегося первым, и момент его освобождения MIN. Переход из блока 12 произойдет на блок 13. Блок 13 идентичен блоку 19 и осуществляет накопление статистики по пребыванию системы в состоянии с M заявками. Блок 14 имитирует окончание обслуживания заявки в канале S. Счетчик обслуженных заявок увеличивается на единицу, а модельное время становится равным времени MIN выхода заявки из канала S. В блоке 15 осуществляется проверка наличия заявок, ожидающих в очереди. Если заявки в очереди есть, то ожидающая заявка сразу направляется в только что освободившийся канал S, генерируется время ее обслуживания DTS и определяется момент окончания обслуживания TD[S]. Сгенерированное время DTS учитывается также в статистике времени занятости канала TOS[S] (блоки 16-17). Если очередь пуста, то канал S объявляется свободным и в качестве момента его освобождения TD[S] устанавливается заведомо большое число. Для систем с бесконечной этот алгоритм обычно тоже применим, при этом в LMAX должно быть установлено достаточно большое число (например 100). Вычисление характеристик СМО. По окончании моделирования основные характеристики СМО могут быть получены следующим образом: 1. Вероятности пребывания в состоянии, когда в системе ровно i заявок, определяется как доля времени, в течение которого в системе было i заявок, относительно общего времени моделирования: 123 pi  TL[i] T 2. Вероятность отказа — доля заявок, получивших отказ, в общем количестве заявок: p отк .  LOS N 3. Пропускная способность — количество обслуживаемых в единицу времени заявок: A KS T 4. Среднее количество занятых каналов можно определить как отношения суммарного времени занятости всех каналов к общему времени моделирования: NK k ср.   TOS[i] i 1 T Среднее время ожидания в очереди можно определить как 5. суммарное время ожидания заявок в очереди, отнесенное к количеству входящих заявок. Если B — это время, проведенное всеми заявками в системе, и C — время, проведенное всеми заявками в каналах, то их разность — время, проведенное всеми заявками в очереди. Тогда B LMAX  NK  m  TL[m] , m0 C NK  TOS[i] , i 1 t ож.  BC N ; 6. Средняя длина очереди может быть вычислена как математическое ожидание количества заявок в очереди в процессе моделирования: NK  LMAX rср.  NK  LMAX  i  NK 1 TL[i]  (i  NK )  T  TL[i]  (i  NK) i  NK 1 T 124 Полученные по этим формулам характеристики СМО необходимо сравнить с характеристиками, полученными в результате аналитического расчета. Программы численного моделирования СМО могут быть написаны на любом алгоритмическом языке высокого уровня. Вопросы: 1. В чем отличие методов математического и имитационного моделирования? 2. В чем отличие методов продвижения модельного времени с постоянным шагом и по особым состояниям? 3. Как строятся списки текущих и будущих событий? 4. Как вычисляется вероятность пребывания СМО в одном из состояний при численном моделировании? 5. Как определяется состояние, в которое перейдет СМО при численном моделировании на основе Марковского процесса? 125 Лекция 13. ИМИТАЦИОННЫЕ И ФИЗИЧЕСКИЕ МОДЕЛИ 13.1. Оценка результатов численного моделирования СМО Имитационное моделирование систем массового обслуживания проводится путем прямого моделирования алгоритма работы СМО по событиям или на основе модели Марковского процесса. Внешняя среда для имитационной модели создается при помощи генераторов случайных чисел, результаты моделирования, включающие основные характеристики СМО, вычисляются методами математической статистики. Для получения статистически обоснованных результатов имитационного моделирования необходимо создать потоки псевдослучайных чисел (поток заявок на обслуживание и поток обслуживания), содержащих более 20 тыс. псевдослучайных псевдослучайных чисел, чисел. соответствующих Поток заявок времени до состоит из поступления следующей заявки. Поток обслуживания состоит из псевдослучайных чисел, соответствующих времени обслуживания заявки в канале обслуживания. Запуск программы имитационного моделирования на компьютере, выполняющей генерацию потоков случайных чисел и обработку результатов моделирования, называется прогоном модели. Перед каждым прогоном изменяются параметры СМО, например интенсивности потока заявок и потока обслуживания, или число каналов в СМО. Каждый прогон модели является вычислительным экспериментом, позволяющим экспериментально оценить выходные характеристики Yiвых СМО. Перечислим обычно решаемые при моделировании задачи: 1. Оценка вероятности pi нахождения системы в состоянии Si, т.е. когда в системе находится ровно i заявок. Фактически в каждом прогоне 126 вычисляется относительная частота pi(k) нахождения системы в состоянии Si , вычисляемая по формуле (k ) pi  (k ) Ti T ; (k ) где Ti(k) – среднее время пребывания системы в состоянии Si в прогоне k, T(k) - время прогона k. 2. Оценка вероятности отказа pотк - доля не обслуженных заявок в общем количестве заявок. В каждом прогоне k вычисляется относительная частота отказа (k ) p отк  N отк N  N отк , где Nотк – количество не обслуженных заявок, N – число обслуженных заявок. 3. Оценка пропускной способности А - количество заявок, обслуживаемых в единицу времени. В каждом прогоне k вычисляется относительная частота A(k) A 4. (k )  (k) N (k) T . Оценка среднего количества занятых каналов kср можно определить как отношение суммарного времени нахождения системы в состояниях S1, S2, … Sm (m – максимальное число каналов обслуживания) к общему времени моделирования. Относительная частота k(k)ср , полученная в k-м прогоне, равна (k ) k ср  m (k )  Ti i 1 (k) . T 127 5. Оценка времени, проведенного всеми заявками в системе. Можно найти время B(k), проведенное заявками в системе в k-м прогоне. (k ) B  m n  (k ) i  Ti , i 1 где m и n - число каналов обслуживания и число мест в очереди, Ti время пребывания системы в состоянии i. 6. Оценка времени, проведенного всеми заявками в каналах обслуживания. Можно найти время С(k), проведенное заявками в каналах обслуживания в k-м прогоне C (k ) m   i  Ti (k ) , i 1 7. Оценка времени, проведенного всеми заявками в очереди. Можно найти время D(k), проведенное заявками в очереди в k-м прогоне mn (k ) (k ) D  (i  m)  Ti . i  m 1 Эти оценки позволяют найти несколько важных характеристик СМО.  8. Оценка среднего времени обслуживания заявки в СМО. Находим среднее время обслуживания в k-м прогоне (k ) t обсл  9. (k) B N Оценка среднего времени обслуживания заявки в каналах обслуживания СМО. Находим среднее время обслуживания в k-м прогоне (k ) (k) t обсл _ к  CN . 10. Оценка среднего времени ожидания заявки в очереди. Находим среднее время ожидания в k-м прогоне 128 (k ) t ож 11.  (k) D N . Оценка средней длины очереди. Находится средняя длина очереди в k-м прогоне как математическое ожидание (k ) L оч n   i  p m i . i 1 Также просто могут быть встроены в программу и найдены другие интересные для конкретных задач моделирования оценки выходных характеристик Y(k)i вых СМО в каждом прогоне. При проведении имитационного моделирования проводится достаточное число прогонов модели с разными реализациями входных потоков заявок на обслуживание. В результате получаются достаточно большие выборки выходных характеристик, которые могут быть подвергнуты статистической обработке. Отметим следующие статистические задачи, связанные с оценкой результатов моделирования: • Статистическая оценка выходных характеристик; • Определение необходимого числа прогонов; • Аппроксимация полученных в результате моделирования методами математической распределений аналитическими формулами. Эти задачи решаются стандартными статистики. На первом этапе строятся гистограммы, затем находятся оценки среднего значения, дисперсии и других моментов, затем высказывается гипотеза о действительной величине оцененного параметра или о виде распределения и находится вероятность оправданности этой гипотезы. 129 13.2. Моделирование на основе физических моделей Математические модели физических систем в большинстве случаев отражают поведение моделируемого объекта во времени. На практике используются математические модели с непрерывным и дискретным временем. Такие математические модели обобщенно называются динамическими системами. Математическая модель динамической системы с непрерывным временем базируется на дифференциальном уравнении или системе дифференциальных уравнений. Дифференциальным уравнением называется уравнение, в котором в качестве неизвестного входит функция одной или нескольких переменных и производные от нее. Если неизвестная функция зависят от многих переменных, то уравнения называются уравнениями в частных производных. Если неизвестная функция зависит от одной переменной, то получаем обыкновенное дифференциальное уравнение. В качестве примера моделей простейших динамических систем с непрерывным временем рассмотрим процесс колебаний математического маятника и процесс изменения энергии электромагнитного поля в колебательном контуре. Математической моделью будет дифференциальное уравнение второго порядка. Схема математического маятника показана на рис. 13.1, схема колебательного контура без потерь – на рис. 13.2. Процесс малых колебаний маятника описывается обыкновенным дифференциальным уравнением mL2  d 2  t  dt 2   mgL  t   0; ,   где m, L - масса и длина подвеса маятника; g - ускорение свободного падения; 130   t  - угол отклонения маятника в момент времени t. θ m Рис.13.1. Математический маятник L q C Рис. 13.2. Колебательный контур без потерь Решение этого уравнения моделирует процесс колебаний идеального математического маятника. Аналогично, процессы в электрическом колебательном контуре описываются обыкновенным дифференциальным уравнением L  d 2q  t  dt 2    q  t  C   0; ,   131 где L, C - индуктивность и емкость конденсатора; q(t ) - заряд конденсатора в момент времени t. Сделаем замену переменных дифференциальных уравнениях физических моделей: h  mL2  Lx ; h  0; h  mgL  1/ C;  (t )  q(t )  z(t ); 2 1 После подстановки новых параметров и переменных в уравнения колебаний математического маятника и колебательного контура получим одинаковое дифференциальное уравнение второго порядка. Это дифференциальное уравнение гармонических колебаний в системах без h  d 2 z  t  dt 2   h  dz  t  dt   h z  t   0; 0 1 2  где h , h и h  безразмерные коэффициенты дифференциального 0 1 2 уравнения, z(t )  независимая переменная. В результате мы получили важный факт, что разные по природе физические объекты могут иметь одинаковые математические модели. Одинаковые математические модели имеют другие физические объекты. Так, например, электродинамика (наука о свойствах и распространении электромагнитных полей) и механика сплошных сред (наука о движении жидкостей) имеют одинаковые математические модели, основанные на дифференциальных уравнениях в частных производных. Усложним однородное дифференциальное уравнение гармонических колебаний. Добавим в правую часть дифференциального уравнения вынуждающую силу дифференциальное x ( t ). В уравнение результате второго получим порядка, неоднородное учитывающее взаимодействие системы и окружающей среды. Уравнение гармонических колебаний системы, взаимодействующей с окружающей средой, имеет следующий вид: h0  d 2 z  t  dt 2   h  dz  t  dt   h z  t   x  t  ;   1 2 132 В математической теории систем правая часть уравнения гармонических колебаний x(t ) является входным воздействием (сигналом) и называется управлением. Переменная z(t ) является выходной переменной. Эту переменную можно принять за внутреннее состояние системы. Две физически разные динамические системы (маятник и колебательный контур) имеют одинаковую линейную математическую модель. Замечание. Приведенные простейшие линейные модели являются приближенными. Модели реальных динамических систем нелинейны. Соответствующие системы дифференциальных уравнений аналитических решений не имеют. 133 Лекция 14.ДИНАМИЧЕСКАЯ СИСТЕМА И КОНЕЧНЫЕ АВТОМАТЫ 14.1. Модель Динамическая Система Любой объект, совершающий какие либо действия в ответ на входные воздействия, может быть назван динамической системой. Будем рассматривать динамическую систему как структуру, в которую вводится нечто (вещество, энергия, информация и т.д.) и из которой в какие-то моменты времени выводится что-то. Математическое определение Динамической Системы включает в себя множества входных воздействий, выходных значений и множество внутренних состояний. Понятие динамическая система включает в себя также вспомогательное непрерывное или дискретное множество моментов времени. Обозначим множество входных воздействий векторной функцией U t  с r входами U t  {u ,u , . . u. ,r } ,множество выходных величин векторной 1 2 функцией Yt с m выходами Y t  { y , y ,..., ym}. 1 2 Множество внутренних состояний динамической системы обозначим векторной функцией X t с n состояниями X t  {x1, x2 ,..., xn}. Введем два отображения, определяющие поведение динамической системы во времени: dX / dt  F t ( X t ,U t , t ) Y t  Gt ( X t ,U t , t). Это Входная и Выходная функции Динамической Системы. Динамическая Система называется стационарной, если входная и выходная функции явно не зависят от времени: dX / dt  F t ( X t ,U t ) Y t  Gt ( X t ,U t ). 134 Схема Динамической Системы показана на рис. 14.1. Рис. 14.1. Динамическая Система Вектор внутреннего состояния динамической системы X называется точкой фазового пространства. Компоненты этого вектора x , x ,..., xn 1 2 называются фазовыми координатами, а функции x (t ), x (t ),..., xn (t )  1 2 фазовыми траекториями. Моделирование динамических систем базируется на решении уравнений для внутренних состояний при заданных начальных условиях и заданных входных воздействиях, входящих в правую часть системы уравнений. Главной технической сложностью при этом является обеспечение высокой точности вычислений при решении нелинейных дифференциальных уравнений. Стандартные численные методы, например, метод Рунге – Кутта или прогноз – коррекция, обычно не приводят к хорошим результатам. Приходится методами вычислительной математики разрабатывать конкретную схему решения дифференциальных уравнений и получать конечно – разностные аналоги исходных уравнений. При выборе конечно – разностных аналогов исходного уравнения необходимо добиваться, что бы для них выполнялись законы сохранения энергии и импульса, а также другие инварианты, присущие исходной задаче. Значительно проще обстоит дело с решением линейных и особенно стационарных систем. Здесь возможны аналитические решения, а принцип 135 суперпозиции позволяет изучать общие свойства линейной динамической системы. Особый подход к решению уравнений Динамической Системы развит в Теории управления. Традиционной задачей Теории управления является задача поиска входных воздействий Ut , обеспечивающих движение Динамической Системы по заданной выходной траектории Y t . Частным случаем задачи Теории Управления являются задачи теории Следящих Систем. Блок схема системы автоматического управления общего вида представлена на (рис. 14.2). На этой схеме: X (t ) - вектор входных воздействий;  (t )  вектор возмущающих воздействий; h' (t )  вектор сигналов ошибки; h'' (t )  вектор управляющих воздействий, Y (t )  вектор выходных переменных. υ h’ X h’’ Управляющая система Y Объект управления Рис 14.2. Структура системы автоматического управления Система автоматического управления под действием входных воздействий X (t ) получает новое значение выходной величины Y (t ). Разность между X (t ) и Y (t ) образует ошибку управления h' (t ), которая должна быть минимальной. 136 При проектировании и эксплуатации систем автоматического управления необходимо выбрать такие параметры системы, которые обеспечили бы требуемую точность управления, а также устойчивость системы в переходном процессе. Если система устойчива, то она представляют практический интерес. Для такой системы изучают поведение системы во времени, максимальное отклонение выходной переменной Y (t ) в переходном процессе, время переходного процесса и т. п. Выводы о свойствах систем автоматического управления различных классов можно сделать в процессе моделирования, например, численно решая систему дифференциальных уравнений при заданных начальных условиях. 14.2. Модель Конечный Автомат Автоматом называется техническое устройство, на вход которое подаются дискретные входные сигналы, а на выходе регистрируются дискретные выходные сигналы. Выходные сигналы зависят от входных сигналов и внутренних состояний автомата. Конечным автоматом называется автомат, у которого множество внутренних состояний и входных сигналов (следовательно, и множество выходных сигналов) являются конечными множествами. Модель Конечного Автомата является частным случаем модели Динамическая Система. Все переменные, Входная и Выходная функции задаются на дискретном множестве моментов времени. Поэтому автомат функционирует по тактам. Другими словами, время t в конечном автомате определяется как дискретное множество тактов, синхронизирующих работу устройства. Поэтому такой автомат называется синхронным. 137 Математическая модель конечного автомата включает алфавит и пару логических функций, определяющих его функционирование: А – входной алфавит {a1, a2, … ai, ….. an}; В – выходной алфавит {b1, b2, … bj, … bm}; Q – внутренний алфавит {q1, q2, … qk, … qp}; Переменные ai и bj называются входными и выходными переменными, переменные qk называются внутренними состояниями автомата. Автомат Мили – синхронный конечный автомат, функционирование которого задается парой функций: t=1,2,3,...  q (t )  ( q (t  1), a (t )) i l  k  b (t )  ( q (t  1)ai (t ))  l  j Здесь qk(t-1), qk(t) – соответственно предшествующее и последующее внутреннее состояния автомата. Заметим, что в автомате Мили внутреннее состояние определяется входной переменной ai(t) на такте t и внутренним состоянием qk(t-1) на предыдущем такте (t – 1). Выходная переменная bj определяется входной переменной на такте t и внутренним состоянием на предыдущем такте (t – 1). Автомат Мура – синхронный конечный автомат, функционирование которого задается парой функций:  qk (t )   ( ql (t  1), ai (t ))  b j (t )  ( ql (t  1)) t= 1,2,3,... В отличие от автомата Мили, в этом автомате выходная переменная bj зависит только от внутреннего состояния qk на предыдущем такте (t – 1). Структурная схема конечного автомата представлена на рис. 14.3. 138 a(t) b(t) Комбинационная часть автомата Φ(…),Ψ(…) q(t) q(t-1) Память Рис. 14.3. Структурная схема конечного автомата Работа конечного автомата обычно описывается при помощи таблиц истинности для соответствующих логических функций. Структура Таблицы переходов показана на следующем рисунке. В этой таблице записаны внутренние состояния q(t) автомата на текущем такте в зависимости от значения входной переменной a(t) и внутреннего состояния q(t-1) на предыдущем такте. Таблица переходов, задающая функцию переходов q (t )   (q (t  1), a (t )) в описании автомата Мили имеет следующий вид: Внутренние состояния Входные переменные Таблица переходов q0(t-1) q1(t-1) q2(t-1) ... qp(t-1) a0(t) q2(t) q0(t) q3(t) … q6(t) а1(t) q3(t) q2(t) q3(t) … q2(t) а2(t) q1(t) q6(t) q2(t) … q6(t) … … … … … an(t) q2(t) q1(t) q5(t) … q4(t) 139 Таблица выходов, задающая функции выходов b(t )   (q (t  1), a (t )) в описании автомата Мили: Внутренние состояния Входные переменные Таблица выходов q0(t-1) q1(t-1) ... qp(t-1) a0(t) b1(t) b3(t) … b5(t) а1(t) b1(t) b5(t) … b2(t) … … … … … an(t) b3(t) b4(t) … b7(t) В этой таблице записаны выходные переменные b(t) автомата на текущем такте в зависимости от значения входной переменной a(t) и внутреннего состояния q(t-1) на предыдущем такте. На основании этих таблиц строится автоматная таблица, показывающая внутреннее состояние автомата q(t) и значение выходной переменной b(t) в зависимости от внутреннего состояния q(t-1) на предыдущем такте и входной переменой a(t). Внутренние состояния Входные переменные Автоматная таблица q0 q1 ... qn a0 q2b1 q0b3 … q0b5 а1 q3b1 q2b5 … q2b3 … … … … … am q2b3 q1b4 … q4b7 Аналогично стоятся таблицы функционирования и для автомата Мура. Существует второй удобный на практике способ описания работы конечных автоматов. Это аппарат графов, вершины которых соответствуют 140 состоянию автомата, а ребра - переходам автомата из одного внутреннего состояния в другое. Модель моделью конечного автомата является удобной математической для описания широкого класса реальных объектов в цифровой технике. В качестве примера таких объектов, в первую очередь следует назвать цифровую схемотехнику, элементы и узлы компьютеров. Все эти объекты имеют дискретные состояния и работают по тактам во времени 141 Лекция 15. СИСТЕМНАЯ ДИНАМИКА 15.1. Элементы модели системной динамики Системная динамика является разделом моделирования и изучает поведение сложных систем, развивающихся во времени. В системной динамике предполагается, что структура системы определяет ее поведение. Другими словами, после построения математической модели, отражающей связи в системе, можно дополнить модель исходными данными и провести моделирование поведения системы при этом варианте исходных данных. Поведение модели системной динамики должно быть подобно поведению реальной системы, при этом в работе модели используются те же причинно – следственные связи, что и в реальности. В системной динамике применяются качественные и количественные модели. Качественная модель строится на основе структурной схемы системы, показывает взаимную зависимость элементов системы и общую структуру системы. В структуру системы в рамках модели системной динамики включаются цепи обратной связи, образующие «петли». Петли обратной связи могут быть двух видов: «Усиливающая» обратная связь и «Ослабляющая» обратная связь. В технических терминах усиливающая обратная связь называется положительной, а ослабляющая – отрицательной. Модель реальной системы содержит достаточно много элементов и большое количество петель обратной связи. Различные петли обратной связи имеют разный «вес», другими словами, по-разному влияют на поведение системы. Для моделирования в рамках системной динамики необходимо дополнить структурную модель конкретными данными и закономерностями, присущими данной системе. Закономерности поведения в моделях системной 142 динамики представляются в форме взаимосвязанных дифференциальных уравнений, по форме близким к уравнениям типовых звеньев в классической теории автоматического регулирования. Наиболее широко в системной динамике применяются модели колебательного звена с обратной связью (дифференциальные уравнения второго порядка), и инерционные звенья, определяющие задержку во времени при протекании моделируемых процессов. В количественной модели системной динамики применяются следующие элементы: • Stock – Резервуар. Резервуар аккумулирует рабочее тело модели (товары, деньги, объемы нефти и т.д.) и характеризуется Уровнем рабочего тела, который может увеличиваться или уменьшаться. Уровень является характеристикой объема рабочего тела, находящегося в резервуаре. • Flow – Поток. Элемент системы, влияющий на состояние резервуара. Уровень рабочего тела в резервуаре под действием потока увеличивается, уменьшается или остается неизменным. • Pace – Темп. Элемент системы, определяющий скорость изменения потока. Темп не накапливается (не аккумулируется) при работе модели. • Connection – Соединение. Элемент системы, устанавливающий связь между резервуарами, Потоками и Переменными. • Delay – Задержка. Временная задержка между процессами, которая имеет существенное влияние на работу системы. Графическое изображение элементов схемы модели системной динамики приведено на рис. 15.1. В правой части рисунка показана информационная цепь обратной связи, по которой информация о Уровне рабочего тела в Резервуаре передается в элемент Темп, который управляет интенсивностью Потока, заполняющего Резервуар. 143 Поток Stoc Уровень Flow Обратная связь Изменение потока Pace Основные элементы Схема управления потоком по цепи обратной связи Рис. 15.1. Основные элементы модели системной динамики Модель системной динамики управляется петлями обратной связи. Обратная связь бывает положительной и отрицательной. Положительная обратная связь увеличивает Поток при увеличении уровня в Резервуаре. Отрицательная обратная связь уменьшает Поток при увеличении уровня в Резервуаре. Поток интегрируется и изменяет уровень в Резервуаре. Уровень в одном Резервуаре может влиять на Уровень в другом Резервуаре только через Темп – скорость изменения потока. Темпы в разных частях модели также не взаимодействуют прямо. Важной особенностью модели системной динамики является возможность создания динамической границы, выделяющей необходимую для моделирования часть структурной схемы системы. Оставшаяся за пределами границы часть структурной схемы принимается за окружающую среду. 15.2. Потоки в модели системной динамики При построении модели системной динамики принято, что движение рабочего тела в Потоке начинается из источника (из «Облака»), лежащего за 144 динамической границей системы. Примеры потоковых диаграмм показаны на следующем рисунке. Граница Ресурсы Поток потребления ресурсов Граница «Облако» источник «Облако» поглотитель Ресурсы Поток поступления ресурсов Граница «Облако» поглотитель Поток потребления ресурсов 15.2. Примеры простейших потоков в моделях системной динамики Модель потока дополняется задержками распространения, задаваемыми при помощи дискретных величин. Возможность указания задержек в потоковых диаграммах и управление ими в процессе моделирования является существенным отличием системной динамики от других имитационных моделей. 15.3. Системы обратной связи При полагается построении на математической причинно-следственные модели связи. разработчик Для этого обычно задается математическая модель, на ее вход подаются исходные данные и после вычислений получают и анализируют результаты. Однако такой подход не учитывает многих особенностей протекания реальных процессов. 145 Рассмотрим производственную систему, на вход которой поступают ресурсы и которая на выходе через время производственного цикла вырабатывает полезный продукт. Допустим, что продукт пользуется успехом на рынке, что позволяет повысить его цену и получить дополнительную прибыль. Часть этой прибыли производитель может вложить в ресурсы и повысить выработку дополнительно полезного произведенного продукта. В результате продажи продукта производитель получит дополнительную прибыль и так далее. Процесс роста производства может продолжаться до момента насыщения рынка. Продажи и прибыль перестают расти, уровень выпуска продукции стабилизируется. Управление работой рассмотренной производственной системой происходит по цепи обратной связи. Сигнал на изменение выпуска продукции определяется разностью между спросом и предложением. Если спрос превышает предложение, то нужно увеличивать выпуск продукции. Это положительная обратная связь. Если наоборот, предложение превышает спрос – то нужно уменьшать выпуск продукции. Это отрицательная обратная связь. Пример потока с обратной связью показан на рис. 15.3. Граница Продукция «Облако» Обратная источник связь Поток производства продукции Рис. 15.3. Поток продукции с обратной связью. Заметим, что в рассмотренном примере цепь обратной связи меняет знак. На этапе роста производства обратная связь положительная, на этапе 146 прекращения роста – отрицательная. Таким образом, цепь обратной связи играет роль управления и стабилизации в модели. 15.4. Нелинейность математической модели Повышение точности компьютерных вычислений при решении задач математического моделирования на физических моделях показало, что окружающий нас мир нелинеен. Линейные модели не обеспечивают достаточной точности согласования результатов моделирования и эксперимента. Показано, что нелинейными являются и достаточно сложные экономические и социологические модели. В моделях системной динамики нелинейность появляется при композиции потоков, охваченных обратными связями. Точность вычислений в физических моделях обеспечивается выполнением законов сохранения и инвариантов, присущих исходной задаче. В экономических моделях роль законов сохранения играют уравнения баланса и условия сохранения товарных потоков. Учет законов сохранения повышает точность вычислений. Однако методы системной динамики часто отказываются от законов сохранения, что позволяет расширить диапазоны изменения параметров моделирования и проанализировать обобщенную ситуацию. Полученные в результате приближенных вычислений законы поведения технической, экономической или социологической модели системной динамики позволяют выделить фундаментальные закономерности, присущие исходной задаче. Эти закономерности позволяют лучше понять функционирование изучаемой системы и уточнить математическую модель. 147 15.5. Возможность моделирования нетривиальных ситуаций Функционирование нелинейной модели системной динамики позволяет получить неожиданные результаты, не предсказуемые при анализе на базе обычных моделей. Возможность достаточной простой коррекции моделей по результатам моделирования делает метод системной динамики очень удобным при изучении реальных технических, экономических и социологических систем. Для проведения решения моделей системной динамики на компьютерах разработан ряд программных пакетов, позволяющих проводить анализ функционирования достаточно сложных практических задач. В настоящее время существуют следующие программные продукты предназначенные для создания и работы с моделями системной динамики: AnyLogic - программа для имитационного моделирования бизнеспроцессов, разработанная российской компанией «Экс Джей Текнолоджис». Использует язык Java для разработки моделей. Помимо методов системной динамики, позволяет использовать дискретно-событийное моделирование и агентно ориентированное е моделирование. Powersim - программа системно-динамического и дискретно- событийного моделирования. Имеет тесную интергацию с SAP. Vensim - программа системно-динамического моделирования. iThink - программа моделирования, основанная на методологии системной динамики. Пакет iThink Analysis, разработанный фирмой High Performance Systems и предназначенный для работы по технологии System Dynamics. Пакет реализует функции, необходимые для построения и анализа моделей системной динамики и графического представления результатов моделирования. В пакете предусмотрена конвертация модели в программу, 148 написанную на языке высокого уровня, что расширяет возможности моделирования и обработки результатов. Замечание 1. Моделирование сложных систем по технологии системной динамики реализует принцип моделирования «сверху вниз». на первом этапе строится обобщенная, достаточно простая модель, при помощи моделирования определяется ее адекватность, а затем модель уточняется и усложняется путем добавления промежуточных аккумуляторов и потоков, связанных цепями обратной связи. Замечание 2. Технология системной динамики универсальна и применяется при построении физических и имитационных моделей. 149 Лекция 16. АГЕНТНО – ОРИЕНТИРОВАННЫЕ МОДЕЛИ 16.1. Методы агентного моделирования Агентно ориентированное моделирование – один из новых подходов к моделирования в технике, экономике и социологии. Несмотря на очень простую идеологию агентно ориентированное моделирование пока находит ограниченное применение. Это связано с трудностями по созданию модели основного элемента процесса моделирования – агента. Системы моделирования, работающие с простыми агентами, не позволяют учитывать все важные свойства системы. Системы со сложными агентами требуют для моделирования гигантских вычислительных ресурсов. Агентное моделирование широко применяется в социологии. В качестве примера посмотрим результаты моделирования поведения толпы, покидающей помещение при возникновении чрезвычайной ситуации (рис. 16.1)[ ]. Рис. 16.1. Моделирование поведения толпы. Толпа в начальный момент равномерно заполняет помещение. Каждый агент в толпе содержит индивидуальную жесткую модель поведения, в соответствии с которой агент стремится двигаться к ближайшему выходу из помещения с учетом противодействия окружающих его других агентов. В 150 результате на некотором дискретном шаге модельного времени образуется расположение агентов, показанное на рис. 16.1. Агентным моделированием поведения групп людей занимаются достаточно давно и поэтому здесь наработаны достаточно точные феноменологические модели поведения агентов. Агент принимает решение о своем поведении на основе локального анализа ситуации в своем ближайшем окружении. Рассмотрим основные свойства агентно ориентированной модели: 1. Модель агента автономна. Агенты не зависимы, внешняя система явно не управляет их поведением. Модель каждого агента содержит макро- и микроуровни поведения. Модели макроуровня являются одинаковыми для всех агентов. Модели микроуровня у агентов могут быть различны, что отражает индивидуальные различия людей в толпе. Взаимодействие агентов определяется поведением на микроуровне, однако результат этого взаимодействия влияет на поведение агентов на макроуровне. 2. Модель агента неоднородна. Модели агентов отличаются друг от друга, что позволяет учесть индивидуальные различия агентов. 3. Модель агента ограничена. Агент может принимать только те решения, варианты которых предусмотрены в макромодели. 4. Модели агентов расположены в ограниченном модельном пространстве. Это может быть плоская площадка с границами или ограниченное трехмерное пространство. Шаг перемещения агента – дискретный. При решении некоторых задач задание геометрических размеров агентов не требуется и тогда они являются точками в модельном пространстве. Главным отличием агентно ориентированных моделей является метод их организации «снизу». Строится модель поведения каждого агента, правила их взаимодействия между собой и границами модельного пространства. Модель содержит большое количество взаимодействующих 151 агентов, при этом модели каждого агента могут быть разными. Эти отличия являются принципиальными, что позволяет выделить объектно ориентированные модели в отдельный класс методов имитационного моделирования. Отметим основные преимущества агентно ориентированного подхода в моделировании. 1. точные Агентно ориентированный подход позволяет создать более модели ряда реальных объектов. Преимущества агентно ориентированного подхода особенно проявляются при моделировании сложных, многоэлементных систем с различными свойствами составляющих их элементов. Сложность модели, по сути, ограничивается только возможностями компьютера. Модели, определяющие поведение агента являются локальными и поэтому достаточно простыми, часто имитирующими поведение человека при заданных обстоятельствах. Такие модели успешно применяются в социологии при моделировании коллективного поведения групп людей в различных обстоятельствах. Однако существуют примеры применения агентно ориентированных моделей к производственным и транспортным системам. 2. Практика агентно ориентированного моделирования показала, очень часто полученные в результате моделирования результаты являются неожиданными и непредсказуемыми. В этом проявляются закономерности в свойствах больших систем, когда в них проявляются закономерности, отсутствующие в их элементах или подсистемах. 3. Агентно ориентированное моделирование позволяет создавать модели сложных систем, законы функционирования которых не известны. Для построения такой модели достаточно знать правила работы или поведения отдельных элементов системы. Свойства системы целиком можно изучить по результатам моделирования. 152 4. Агентно ориентированное моделирование является очередным шагом развития методов математического моделирования, позволяющим создавать модель, содержащую множество других независимых алгоритмов и моделей. Результаты моделирования не являются простым объединением свойств отдельных подсистем модели, но позволяют получить новые свойства, присущие сложной моделируемой системе. 16.2. Разработка агентно ориентированных моделей Будем считать, что создаваемая агентно ориентированная модель является инвариантной относительно предметной области. Для этого необходимо создать абстрактную сеть, поместить в нее программных агентов с общими свойствами, уточнить топологию сети и поместить эту абстрактную модель в среду моделирования. Среда моделирования должна иметь ресурсы для простого создания модели агента, а также для объединения агентов в сеть, проведение моделирования и наглядного графического вывода результатов. Модель агента должна содержать текущее внутреннее состояние, правила поведения, механизм принятия решения, механизм коррекции элементов модели в процессе имитационного моделирования. Для создания модели необходимо: • определить границы применения модели; • найти преимущества агентного подхода по сравнению с другими подходами к моделированию; • построить математические модели; • построить алгоритмы поведения и взаимодействия агентов. При построении модели агента необходимо реализовать онтологический подход, состоящий в необходимости получения и хранения текущей информации о предметной области, о других агентах и окружающей 153 среде. Эта задача является достаточно сложной, так как целесообразно создавать универсальную модель, инвариантную к предметной области. Иначе затраты на моделирование могут оказаться чрезмерными. Программная модель агента должна создаваться в параметризованном виде и иметь возможность простой настройки при переходе в новую предметную область. Блок схема Агента приведена на рис. 16.2. Аппаратный блок Входы Контроллер Память Программный блок Механизм принятия решений Выходы Математическая модель Диспетчер Интерфейсы Диспетчер задач Рис. 16.2. Блок схема Агента На рисунке представлена Условно Графическая реализация Агента. На практике аппаратные функции контроллера, памяти и интерфейсов реализуются программно, в результате чего Агент является программным модулем, входящим в среду моделирования. Часть входов Агента подключаются к выходам других агентов, что обеспечивает их взаимодействие. Оставшаяся часть входов выполняет роль датчиков (рецепторов), обеспечивающих связь Агента с окружающей средой. Блок Диспетчер задач через свой вход/выход обеспечивает связь Агента с центральным блоком управления процессом моделирования. Блок схема системы агентного моделирования показана на рис. 16.3. Основу системы составляет среда моделирования, блок управления и набор модулей Агентов. 154 Содержание модели Среда моделирования Модули Агентов Система управления Агентным моделированием Система обработки результатов моделирования Рис. 16.3. Система агентно ориентированного моделирования Система агентно ориентированного моделирования содержит среду моделирования, представляющую собой программный продукт, позволяющий определить границы области моделирования, установить в области моделирования необходимое количество Агентов, установить связи между агентами и провести их настройку. Например, модель среды в агентной модели динамики населения страны, построенной в программе AnyLogic, включает наличие и расположение жилья и рабочих мест, транспортную инфраструктуру, наличие и расположение объектов образования и здравоохранения и т.д. Модель агента оформлена в форме карты состояний, в окнах которой записываются параметры модели. Карты состояний позволяют графически определить возможные состояния агента, переходы между ними, события, вызывающие эти переходы, временные задержки и действия, совершаемые агентом на протяжении своей жизни. Такие конструкции, как вложенные состояния, позволяют задавать режимы функционирования агента. Агент может иметь несколько параллельно активных и взаимодействующих карт состояний, 155 каждая из которых отвечает за какой-либо аспект его жизни: например, работу, образование и семейное положение. Система управления агентным моделированием выполняет дискретное продвижение модельного времени и фиксирует состояние имитационной модели на очередном шаге моделирования. Система обработки результатов осуществляет построение графических изображений и обеспечивает диалоговый режим взаимодействия с моделью. Отметим особенности процесса имитационного моделирования в моделях системной динамики и агентно ориентированного моделирования. Модель системной динамики строится «сверху вниз», при этом основу модели составляет ее структура и система обратных связей. Заявки, проходящие в модели, являются пассивными объектами. Агентно ориентированная модель строится «снизу вверх». Агенты, находящиеся на нижнем уровне модели, определяют ее функционирование. В результате моделирования получаем глобальные характеристики и свойства реальной системы. Оба подхода в моделировании можно реализовать в программе AnyLogic. Это относительно новый инструмент, который использует компания XJ Technologies. Существенным преимуществом про AiiyLogic является возможность быстрого создания профессиональных агентных моделей в той же самой графической среде, широко встраивая в модель программные блоки, написанные на основных языках программирования. Существует также возможность использовать разные подходы для разных частей модели. 156 Лекция 17. МОДЕЛИ, ОСНОВАННЫЕ НА СЕТЯХ ПЕТРИ В современных вычислительных системах особое внимание уделяется организации параллельных процессов. Примером параллельного процесса в работе любого компьютера являются вычисления, проводимые процессором, и одновременный фоновый вывод информации на принтер. Параллельно работающие системы могут взаимодействовать между собой различным образом, а могут и вообще не взаимодействовать. Параллельные процессы независимыми. могут Подчиненный быть ведущими, параллельный подчиненными процесс и управляется (синхронизируется) ведущим процессом. Синхронные процессы могут обмениваться информацией. Независимые процессы являются асинхронными и могут в процессе работы вообще не обмениваться информацией. Отметим особенности вычислительных многопроцессорных системах протекания параллельных и Понятно, системах сетях. можно процессов что организовать только в в параллельный вычислительный процесс. Вычисления в каждом процессоре могут быть синхронизированы или не синхронизированы между собой. Однако работа процессоров всегда согласуется при обращении к общим ресурсам вычислительной системы или сети. Использование общих ресурсов вычислительной системы может организовываться по схеме взаимного исключения, когда процесс, первым захвативший общий ресурс, запрещает доступ к нему других параллельных процессов. Эта простейшая схема может быть усовершенствована путем присвоения процессам различных абсолютных или относительных приоритетов. Синхронизация параллельных процессов может проводиться на аппаратном уровне путем использования общего тактового генератора или на программном (системном) уровне путем обмена сигналами по 157 установленному протоколу. Синхронизация параллельных процессов в вычислительной системе организуется на этапе разработки и включает в себя специальные аппаратные средства и программы операционной системы. Моделирование параллельных процессов проще всего проводить в рамках специализированных языков моделирования (например, GPSS) и на алгоритмических языках реального времени (например язык АДА). Однако перед написанием программы моделирования необходимо описать алгоритм реализации параллельных процессов, что сделать традиционными методами очень даже не просто. Одним из распространенных и эффективных способов описания параллельных процессов являются сети Петри. Сети Петри могут быть определены в графической форме в виде специфичного графа или в формально аналитической форме, позволяющей легко перейти к реализации в виде компьютерной программы. 17.1. Определение сети Петри Определение: Сеть Петри С включает четыре объекта {P, T, I, O}:  конечное множество позиций P = {p1, p2, … pn};  конечное множество переходов T = {t1, t2, … tm};  входную функцию I: T ← P, отображение из позиций в переходы;  выходную функцию O: T → P, отображение из переходов в позиции. Позиция pi является входной позицией для перехода tj в том случае, если p i  I(t j ) ; позиция pi является выходной позицией для перехода tj в том случае, если pi  O(t.j ) Пример 1. Примером аналитического описания структуры сети Петри является следующий текст: 158 17.2. Объекты, образующие сеть Петри C = {P, T, I, O}; Множество позиций P = {p1, p2, p3, p4, p5} Множество переходов T = {t1, t2, t3, t4, t5} Входная функция Выходная функция I(t1) = { p2, p4 }; O(t1) = { p1 }; I(t2) = { p4 }; O(t2) = { p2 }; I(t3) = { p1 }; O(t3) = { p3 }; I(t4) = { p3, p5 }; O(t4) = { p4 }; I(t5) = { p3 }; O(t5) = { p5 }; Определим расширенную входную функцию I: P ← T, отображение из переходов в позиции; и расширенную выходную функцию O: P → T, отображение из позиций в переходы. Эти функции по форме являются обратными к ранее введенным входной и выходной функциям. Для сети Петри, определенной выше эти функции имеют следующий вид: Расширенная входная Расширенная выходная функция функция I(p1) = { t1 }; O(p1) = { t3 }; I(p2) = { t2 }; O(p2) = { t1 }; I(p3) = { t3 }; O(p3) = { t4, t5 }; I(p4) = { t4 }; O(p4) = { t1, t2 }; I(p5) = { t5 }; O(p5) = { t4 }; 159 Для наглядного описания сети Петри применяется граф, обладающий двумя типами узлов: Кружок является позицией, черта – переходом. Стрелки соединяют позиции и переходы, при этом стрелка, направленная от позиции к переходу определяет входную позицию для этого перехода, а стрелка, направленная от перехода к позиции – выходную позицию для этого перехода. В качестве примера рассмотрим сеть Петри с графом, представленным на рис. 17.1 а. p1 t3 p3 t5 p5 t1 t4 p2 p4 t2 Рис. 17.1а. Граф сети Петри. Начальная маркировка сети p1 t3 p3 t5 p5 t1 t4 p2 p4 t2 Рис. 17.1б. Граф сети Петри. Маркировка после первого такта Сеть Петри, представленная в форме графа на рис. 17.1, соответствует аналитическое описание, приведенное выше в примере 1. 17.3. Маркировка сети Петри Позициям в сети Петри могут быть присвоены метки или фишки, служащие для определения процесса работы сети. Маркировка сети Петри на 160 рис. 6.1 представлена в виде точек (меток, фишек), находящихся в некоторых позициях сети. 2. В аналитической форме маркировка μ сети Петри определяется как отображение множества позиций P в множество неотрицательных целых чисел N μ: P → N. Маркировка может быть задана и в форме вектора μ = { μ1, μ2, … μn} (μi = 0, 1, 2, …), где индекс i при компоненте вектора μi совпадает с номером соответствующей позиции pi. Вообще говоря в каждой позиции может быть ноль, одна или несколько фишек. Приведем маркировку сети для примера 1, показанную на рис. 7.1а. Маркировка сети μ = { μ1, μ2, μ3, μ4, μ5} μ1 = 1; μ2 = 0; μ3 = 0; μ4 = 1; μ5 = 1; В результате аналитическое определение маркированной сети Петри включает пять объектов {P, T, I, O, μ}. Вообще говоря число фишек в каждой позиции сети Петри не ограничено. Поэтому каждая сеть Петри может иметь бесконечное число маркировок. Сеть может быть "оживлена", при этом фишки будут переходить из одной позиции в другую, изменяя маркировку сети. Фишки могут перемещаться по сети только через "активные" (разрешенные) переходы. Переход активен, если в каждой для него входной позиции имеется хотя бы одна метка. 161 В примере на рис. 17.1а разрешены (активны) переходы t2, t3 и t4, во входных позициях которых p4, p1 и p5 находятся фишки. В результате срабатывания перехода число фишек в его входных позициях уменьшается на единицу, а число фишек в его выходных позициях увеличивается на единицу. Новая маркировка сети, возникшая после срабатывания активных переходов на рис. 17.1а, показана на рис. 17.1б. В результате фишки оказываются в позициях p2, p3 и p5, при этом активным становится переход t4. На следующем такте работы сети произойдет срабатывание перехода t4 , после чего изменится маркировка сети (рис. 17.1в). p1 t3 p3 t5 p5 t1 t4 p2 p4 t2 Рис. 17.1в. Граф сети Петри. Маркировка после второго такта p1 t3 p3 t5 p5 t1 t4 p2 p4 t2 Рис. 17.1г. Граф сети Петри. Маркировка после третьего такта 162 В результате активными становятся переходы t1 и t2 , так как в их входных позициях p2 и p4 имеется по одной фишке. После срабатывания этих переходов маркировка сети изменится еще раз (рис. 17.1г). 17.4. Пространство состояний сети Петри Пространство состояний сети Петри определяется ее маркировкой. Каждое срабатывание активных переходов изменяет маркировку сети и поэтому изменяет ее состояние. Для рассмотренного примера сети Петри маркировке на рис. 1а соответствует состояние S1, маркировке на рис. 1б – состояние S2, маркировке на рис. 1в – состояние S3, маркировке на рис 1г – состояние S4 и так далее. Число различных состояний сети Петри, содержащей n позиций и не более N фишек в каждой позиции не превышает число Nn. При очередном срабатывании активных переходов возникает новая маркировка сети Петри. Поэтому в процессе работы возникают две последовательности:  последовательность маркировок μ0, μ1, μ2, μ3, μ4, …  и последовательность активных переходов {ti0}, {ti1}, {ti2}, {ti3}, {ti4}, … Эти две последовательности связаны между собой при помощи функции следующего состояния δ(μk, tik) = μk+1, определяющей маркировку сети на следующем такте работы. Если хотя бы один переход tik в маркировке μk активен, то в результате его срабатывания возникнет новая маркировка μk+1, непосредственно достижимая из маркировки μk. Понятно, что в результате последовательности переходов из одной маркировки в другую непосредственно достижимую маркировку (из состояния Si в состояние Si+1 ) возникает множество достижимости R(S, μ0) 163 для сети Петри, определяемое начальной маркировкой μ0 . Множество достижимости совпадает с наименьшим множеством всех маркировок (всех состояний), достижимых из одного начального состояния сети Петри. Приведем пример сети Петри, моделирующую работу простейшего компьютера, содержащего устройство ввода, устройство вывода и процессор (рис. 17.2). Поступление очередного задания на вход компьютера моделируется при помощи фишки, поступающей в позицию p1 через переход t1. Если процессор свободен (в позиции p2 находится фишка), то переход t2 окажется активным. После его срабатывания фишки в позициях p1 и p2 исчезнут, и появится фишка в позиции p3 , что означает обработку задания в процессоре. В результате оказывается активным переход t3. После срабатывания этого перехода фишки окажутся в позициях p2 и p4. Из позиции p4 фишка через активный переход t4 покинет сеть Петри. t1 t1 – задание помещается во входную очередь, p1 – задание ждет, p2 – процессор свободен, t2 – начало выполнения задания в процессоре, p3 – задание выполняется, t3 – завершение выполнения задания, p4 – задание ждет вывода, t4 – задание выводится. p1 t2 p2 p3 t3 p4 t4 Рис. 17.2. Сеть Петри - модель простого компьютера Этот пример показывает, как в сети Петри происходят события, т.е. действия, связанные с переходом сети из одного состояния в другое. Для 164 того, что бы событие произошло необходимо, что бы выполнилось некоторое условие (фишки находятся во всех входных позициях одного из переходов). Эти условия называют предусловиями для события. Так, в предыдущем примере для события, состоящего в срабатывании перехода t2 необходимо выполнить условие, состоящее в наличии фишек в позициях p1 и p2. 17.5. Моделирование параллельных процессов В сети Петри два невзаимодействующих события могут происходить независимо друг от друга. Это позволяет моделировать два или несколько асинхронных независимых процессов. Два связанных между собой события должны происходить одновременно, что также легко моделируется в сети Петри. Считается, что событие происходит мгновенно, за нулевое время, поэтому вероятность одновременного возникновения двух событий принимается равной нулю. Такое событие называется примитивным. Примитивные события не могут быть одновременными. Непримитивными событиями называются события, длительность которых конечна и больше нуля. Такие события могут пересекаться во времени и соответствуют реальным событиям в реальных системах. Ясно, что сеть Петри с примитивными событиями не совсем точно отражает реальные системы. Однако на практике при построении модели на основе сети Петри это обстоятельство легко обходится. Так непримитивное событие заменяется двумя примитивными: "начало непримитивного события", "конец непримитивного события" и условия: "непримитивное событие происходит". начало условие конец Рис. 17.3. Моделирование непримитивного события 165 Вторым широко применяемым способом моделирования непримитивного события является обозначение соответствующего перехода прямоугольником. непримитивный переход входные выходные позиции позиции Рис.17.4. Моделирование непримитивного события В качестве примера покажем сеть Петри для моделирования компьютера, в работы котором простого время работы процессора учитывается как непримитивное t1 p1 событие (рис. 17.5). На этой схеме при p2 поступлении фишки в позицию p1 происходит срабатывание непримитивного перехода, в результате чего через время Δt фишка окажется в позиции p3. Обратим внимание на важную p3 особенность сетей Петри. Они в принципе моделируют независимые параллельные t2 процессы как асинхронные. Рис. 17.5. Сеть Петри для моделирования работы простого компьютера на основе непримитивного события 166 Синхронизация событий возникает автоматически в точках объединения параллельных процессов, в тот момент, когда во всех входных позициях общего перехода возникают фишки. Эта ситуация показана на рис. 17.6. Применение прямоугольника для моделирования отдельных переходов приводит к значительному упрощению сети Петри, что важно при моделировании сложных систем. Таким образом может быть свернут целый функциональный блок, что эквивалентно применению подпрограмм в программировании. Два параллельных процесса объединяются в переходе t5, что и приводит к синхронизации. t1 p1 t2 p2 t5 t3 Рис. p3 t4 p4 17.6. Автоматическая синхронизация двух процессов при их объединении на переходе t5 параллельных Достаточно просто моделируется и ситуация "захвата ресурсов", когда один из параллельных процессов блокирует работу другого при возникновении некоторого условия (рис. 17.7). На этом рисунке фишка, проходящая по верхнему параллельному процессу, первой достигает позиции p2 и через разрешенный переход t5 перемещается дальше. В результате срабатывания перехода t5 фишка из позиции p5 удаляется и поэтому нижний параллельный процесс остановится в позиции p4. 167 t1 p1 t2 p2 t5 p5 t6 t3 p3 t4 p4 Рис. 17.7. Блокировка одного параллельного процесса другим Дальнейшее продвижение фишек через переходы t5 и t6 может произойти только после попадания фишки в "разрешающую" позицию p5 . Заметим, что при наличии фишек в позициях p2 и p4 приход фишки в "разрешающую" позицию p5 синхронизирует работу переходов t5 и t6. 17.6. Моделирование процессора с конвейерной обработкой Большинство современных процессоров имеют конвейерную организацию вычислительного процесса. Конвейер состоит из ряда операций, которые могут выполняться одновременно. Каждая операция конвейера получает исходные данные от предыдущей операции и передает результат последующей. 3. Наибольший выигрыш времени вычислений в процессоре с конвейером достигается при выполнении однородных операций, например, при сложении ряда чисел. Каждая операция сложения чисел с плавающей запятой может быть представлена как следующая последовательность элементарных операций: 1. сравнить экспоненциальные части чисел и найти разность между ними, 2. если необходимо, переместить точку в мантиссе меньшего числа, 168 3. сложить мантиссы, 4. нормализовать сумму, 5. вычислить новое значение экспоненты, 6. проверить это значение на переполнение, 7. сформировать формат суммы. Каждая операция может быть выполнена специальным аппаратным блоком. После выполнения первой операции результат передается во второй блок, при этом первый блок оказывается свободным. Поэтому в нем одновременно с работой второго блока по сложению первых двух чисел может быть выполнена первая операция по сложению третьего и четвертого числа. Если последовательность слагаемых достаточно велика, то в конвейере может одновременно складываться семь чисел, при этом результаты сложения выдаются в темпе один такт вместо семи тактов при обычной организации вычислений. Синхронизация работы блоков конвейера выполняется различным образом. При синхронной организации информация из блока в блок передается в темпе работы самого медленного блока. При более сложной асинхронной организации информация после обработки в блоке n поступает в его выходной регистр и ждет там до тех пор, пока не освободится входной регистр следующего блока n + 1. В асинхронном режиме скорость вычислений в среднем выше, хотя организация вычислений более сложна. 169 Блок n – 1 конвейера Выходной регистр блока n – 1 Входной регистр блока n Блок конвейера n Выходной регистр блока n Входной регистр блока n + 1 Рис. 17.8. Организация конвейерной обработки чисел в процессоре Таким образом при асинхронной конвейерной обработке для каждого вычислительного блока n возникают следующие состояния сети Петри:  входной регистр блока n заполнен,  входной регистр блока n свободен,  выходной регистр блока n заполнен,  выходной регистр блока n пуст,  вычислительный блок n занят,  вычислительный блок n свободен. Построим модель асинхронной конвейерной обработки в форме сети Петри. Устройство управления процессом конвейерной обработки в процессоре (рис. 17.9) находится в исходном состоянии и готово к работе. При поступлении второй фишки в переход t1 он срабатывает и устанавливает фишку в позицию p5. После обработки информации в этом блоке фишка переходит в позицию p9, обозначая окончание работы блока n – 1 конвейера. 170 В результате оказывается разрешенным переход t6 и после его срабатывания фишка оказывается в позиции p6 , обозначая передачу информации во входной регистр следующего блока n. После этого блок n начинает обработку информации, а блок n – 1 освобождается и оказывается готовым к работе. Выходной регистр n – 1 пуст Входной регистр n p2 заполнен Выходной p3 регистр n пуст Входной p4 регистр n + 1 заполнен t1 p1 t2 p5 t5 блок n - 1 занят t6 p6 пересылка из n - 1 в n t3 p7 t7 p11 блок n занят t4 p8 t8 пересылка из n в n + 1 Выходной p9 регистр n – 1 заполнен Входной p10 регистр n пуст Выходной регистр n заполнен p12 Входной регистр n + 1 пуст Рис.17.9. Модель асинхронной конвейерной обработки информации в процессоре 171 Лекция 18. МОДЕЛИ, ОСНОВАННЫЕ НА СЕТЯХ ПЕТРИ 18.1. Кратные функциональные блоки компьютера Для повышения производительности компьютера в состав процессора часто вводятся дополнительные регистры и другие функциональные блоки (блоки сложения и умножения чисел, логические блоки), позволяющие выполнять несколько инструкций программы одновременно. Организация параллельных вычислений на основе кратных функциональных блоков должна быть такой, что бы результат работы по линейной программе и с использованием кратных блоков был бы одинаков. Пусть в программе требуется последовательно выполнить две операции a и b. Если для выполнения операции b не требуются данные, полученные при выполнении операции a , то операцию b можно выполнить одновременно или раньше операции a. Управление программами такого типа основано на использовании таблиц резервирования. В такую таблицу заносятся блоки и регистры, необходимые для выполнения данной операции. Если эти блоки свободны, то такую операцию можно выполнять. Процессор резервирует эти блоки и выдает команду на выполнение этой операции. Если хотя бы один необходимый блок занят, то такая операция не выполняется и процессор ждет, когда возникнут условия, необходимые для выполнения этой операции. Такая ситуация может быть промоделирована сетью Петри (рис. 18.1). Поставим в соответствие каждому функциональному блоку и регистру свою позицию. Если регистр свободен, то пометим его фишкой. В кратные регистры поместим соответствующее число фишек. 172 Команда процессора Блок K свободен следующая команда процессора Выдача команды Регистр 1 свободен Блок K, работающий с регистрами 1, 2, 3. Регистр 2 свободен Регистр 3 свободен Команда выполнена Рис. 18.1. Часть сети Петри, моделирующей устройство управления процессора с кратными функциональными блоками Перед тем, как выдать очередную команду процессор проверяет, что необходимые для ее выполнения регистры свободны (в них находятся фишки). Если это так, то команда начинает выполняться и подается команда к проверке готовности системы к выполнению следующей команды. Приведенная схема очень проста и показывает только идеологию параллелизма, основанного на введении кратных функциональных блоков. Понятно, что некоторые функциональные блоки могут использовать общие регистры, а содержимое регистров может управлять порядком выполнения машинных команд. 18.2.Взаимно исключающие параллельные процессы Очень часто два параллельных вычислительных процесса используют общий аппаратный ресурс, например, общую часть оперативной памяти. Для правильной работы вычислительной системы необходимо обеспечить 173 поочередный доступ к общему ресурсу и блокировать доступ к занятому ресурсу. Эта ситуация изображена на следующем рисунке. Процесс 1 Процесс 2 p1 p6 t1 t4 p2 p7 t2 p3 t5 p4 p8 t3 t6 p5 p9 Рис. 18.2. Сеть Петри для взаимно исключающих процессов На этой схеме позиции p3 и p8 соответствуют критическим сегментам программы, в которых они обращаются к общему аппаратному ресурсу. Этот ресурс захватывает тот процесс, в котором раньше появится фишка в позиции p2 или p7. Пусть появилась фишка в позиции p2. Тогда переход t2 оказывается активным и фишки из позиции p2 и "разрешающей" позиции p4 поступают в позицию p3, соответствующую общему аппаратному ресурсу для двух процессов. В результате переход t3 оказывается активным и после его срабатывания фишки поступают в позицию p5 и "разрешающую" позицию p4. Все время, пока первый процесс находится в позиции p3, критическая позиция p8 оказывается заблокированной для второго процесса. 174 Замечание. В рассматриваемых сетях Петри рассматриваются только примитивные, мгновенно происходящие события, и соответствующие им примитивные переходы. Поэтому вероятность возникновения одновременных событий принимается равной нулю. Замечание. Сети Петри моделируют только логику работы программы или технического устройства. 18.3. Системы с флагами Одним из наиболее распространенных механизмов синхронизации является использование флагов (семафоров), указывающих на возможность использования в программе соответствующего блока. Флаг – это элемент данных, принимающих два значения: "0" – блок свободен, "1" – блок занят. Такие системы очень просто моделируются сетями Петри. p1 Входной "разрешающий" переход Функциональный блок Выходной "разрешающий" переход t1 p2 p3 t2 p4 Рис. 18.3. Моделирование системы с флагом 175 При поступлении фишки в позицию p1 переход t1 оказывается разрешенным. После срабатывания этого перехода фишка попадает в позицию p2 , при этом переход t1 запрещен (нет фишки в позиции p3 ). Переход t2 разрешен и после его срабатывания фишки попадают в позиции p3 и p4 . Таким образом в позиции p2 может находиться только одна фишка (только одна заявка находится в функциональном блоке, моделируемым позицией p2 ). 18.4. Анализ сетей Петри Построение модели вычислительной системы в форме сети Петри является только предварительным этапом моделирования. Для изучения свойств моделируемой системы необходимо провести анализ сети. На практике для этого необходимо как правило написать компьютерную программу и задавая различные исходные данные изучить работу моделируемой системы. Широкие возможности для моделирования открывает система Simulink пакета Matlab. Определим основные свойства сетей Петри, которые являются основными задачами моделирования. Безопасность. Позиция сети Петри называется безопасной, если в ней может быть не более одной фишки. Сеть Петри называется безопасной, если все ее позиции безопасны. Понятно, что безопасная сеть Петри позволяет упростить техническую реализацию моделируемой системы. Такая система реализуется в обычной бинарной логике на базе обычной элементной базы (триггеров, логических микросхем). Безопасные позиции соответствуют логическому условию (есть фишка – условие истинно, нет фишки – условие ложно). Если в позиции может быть несколько фишек, то такая логическая интерпретация невозможна. Таким образом, только в безопасной сети возможно простое моделирование событий и условий. 176 Сеть Петри, являющаяся безопасной в начальный момент времени, может быть сделана принудительно безопасной (кратных дуг в сети быть не должно). Рассмотрим способ обеспечения безопасности на примере. Структура сети Петри не гарантирует безопасность, так как при поступлении в начальную позицию нескольких фишек они могут продвигаться по автоматически разрешенным переходам или накапливаться в некоторых позициях. p1 t1 p2 t2 p3 t3 p4 Рис. 18.4. Пример небезопасной сети Петри Переход t2 на этом рисунке размножает фишки, которые накапливаются в позиции p2. Безопасная сеть может быть создана из небезопасной различными способами. Например, можно ввести дополнительные условия, состоящие в принудительном запрещении попадания в любую позицию более одной фишки. Для принудительного обеспечения безопасности, например, позиции p3 достаточно ввести дополнительную позицию p5, реализующую условие: в позиции p3 нет фишки (рис. 18.5). 177 p5 p1 t1 p2 t2 p3 t3 p4 Рис. 18.5. Сеть Петри, которая, начиная с позиции p3, становится безопасной Понятие безопасности сети Петри естественным образом расширяется на понятие ограниченности, состоящее в том, что в позиции сети не может быть больше чем m фишек. Безопасная позиция может быть аппаратно реализована триггером, а ограниченная – счетчиком. Неограниченную сеть Петри реализовать нельзя. Сохраняемость. Для некоторых вариантов сетевых моделей важным является понятие сохраняемость. Сеть Петри называется сохраняющей, если в процессе работы сохраняется число фишек, перемещающихся в сети. Это важно, когда фишкам соответствуют какие либо ресурсы, например, блоки памяти или устройства ввода и вывода информации. 178 p1 p4 t1 t3 p2 p3 p5 t2 t4 Рис. 18.6. Пример сохраняющей сети Петри В процессе функционирования модели эти блоки не могут возникать или уничтожаться. Условие сохраняемости является очень жестким. В качестве примера рассмотрим сеть, моделирующую взаимно исключающие процессы (рис. 18) (сравни с рис. 14). Если в позицию p1 фишка попадет раньше, чем в позицию p4 , то сработает переход t1 и в позицию p2 попадут две фишки (фишки из позиций p1 и p3 будут удалены). Затем через переход t2 фишки вернутся в позиции p1 и p3 . Ситуация, когда в позиции p1 и p4 фишки попадут одновременно невозможна, так как переходы примитивны. Не сохраняющая сеть Петри показана на рис. 18.7. p1 p4 t1 t3 p2 p3 t2 p5 t4 Рис. 18.7. Не сохраняющая сеть Петри 179 Активность. В процессе функционирования сети Петри могут возникать ситуации, когда некоторые фишки попадают в тупики, в результате чего некоторые переходы сети перестают быть активными. Возможно ситуация, когда неактивны все переходы сети, в результате чего движение фишек в ней прекращается. Понятно, что это ситуация зависит от начального расположения меток в сети (от начальной маркировки). Достижимость. С понятием активность близко связано понятие достижимость. Под достижимостью маркировки μj из начальной (из другой) маркировки μk понимается существование последовательности маркировок, приводящих сеть из состояния с маркировкой μk в состояние с маркировкой μ j. 18.5. Дерево достижимости сети Петри Основным методом анализа сетей Петри является построение дерева достижимости, позволяющий найти множество маркировок, достижимых из заданной начальной маркировки сети Петри. Построение дерева достижимости легко выполняется путем счета по компьютерной программе, написанной для данного варианта сети. Рассмотрим процесс построения дерева достижимости на конкретном примере. p2 t1 p1 t3 t2 p3 Рис. 18.8. Сеть Петри, для которой строится дерево достижимости 180 Начальная маркировка сети μ0 = (1, 0, 0) – фишка находится в позиции p1. В этот момент разрешены два перехода: t1 и t2. Дерево достижимости должно охватывать все возможные варианты эволюции сети. В результате срабатывания перехода t1 будет достигнута маркировка μ1 = (1, 1, 0), а в результате срабатывания перехода t2 - маркировка μ2 = (0, 1, 1). Фрагмент дерева достижимости для первого шага имеет следующий вид: (1, 0, 0) (1, 1, 0) (0, 1, 1) Рис. 18.9. Дерево достижимости для первого шага работы сети Из маркировки μ1 = (1, 1, 0) можно запустить переход t1 и получить маркировку μ3 = (1, 2, 0).и запустить переход t2 и получить переход μ4 = (0, 2, 1). Из маркировки μ2 = (0, 1, 1) можно запустить переход и получить маркировку μ5 = (0, 0, 1). В результате получим дерево достижимости для второго шага работы сети (рис. 18.10). (1, 0, 0) t1 t2 (1, 1, 0) t1 (1, 2, 0 (0, 1, 1) t2 (0, 2, 1) t3 (0, 0, 1) Рис. 18.10. Дерево достижимости для второго шага работы сети 181 Для третьего шага работы сети получим следующее дерево достижимости (1, 0, 0) t1 t2 (1, 1, 0) t1 t2 (1, 2, 0 t1 (1, 3, 0) (0, 1, 1) (0, 2, 1) t3 (0, 0, 1) t2 t3 (0, 3, 1) (0, 1, 1) Рис. 18.11. Дерево достижимости для третьего шага работы сети Продолжая таким образом, заметим, что в результате для конечной сети Петри строится бесконечное дерево достижимости. Приведение такого дерева к конечному виду проводится следующим способом. Введение новых маркировок ограничивается на каждом шаге граничными вершинами. Это:  пассивные маркировки, из которых никакие переходы невозможны (такой является маркировка μ5 = (0, 0, 1), полученная в результате срабатывания переходов t2 t3, и  ранее встречавшиеся маркировки, из которых повторяются уже построенные фрагменты дерева достижимости. Такой является маркировка μ2 = (0, 1, 1), полученная в результате срабатывания последовательности переходов t1 t2 t3 . Для этого дерева только маркировка μ5 = (0, 0, 1) является пассивной. Для всех остальных маркировок некоторые переходы сети оказываются разрешенными, что приводит к постепенному накоплению фишек в позиции 182 p2. Это означает, что сеть имеет бесконечное множество маркировок, что отмечено на конечном дереве достижимости для этой сети. (1, 0, 0) t1 t2 (1, ∞, 0) t1 t2 (1, ∞, 0 t1 (1, ∞, 0) (0, 1, 1) (0, ∞, 1) t2 (0, ∞, 1) t3 (0, 0, 1) t3 (0, ∞, 1) Рис. 18.12. Конечное дерево достижимости для работы сети Если последовательность маркировок ограничена (символ ∞ отсутствует в дереве достижимости), то число состояний сети Петри конечно. Дерево достижимости является графом состояний для этой сети. В этом практически важном случае находить множество всех достижимых маркировок (множество состояний сети) можно простым перебором. Понятно, что только ограниченная сеть Петри является безопасной, в которой число фишек в позициях не может превышать некоторое конечное число. Такая сеть аппаратно осуществима. Пример ограниченной сети Петри приведен на следующем рисунке. 183 t6 t1 p1 p2 t2 t4 t3 p4 t5 p3 Рис. 18.13. Ограниченная сеть Петри Дерево достижимости для этой сети имеет следующий вид: t1 (1, 0, 0, 0) t2 t3 (0, 1, 0, 0) (0, 0, 1, 0) (0, 1, 1, 0) t4 t4 (0, 0, 1, 0) (0, 0, 2,0) t5 (0, 0, 0, 1) t6 (1, 0, 0, 0) Рис. 18.14. Дерево достижимости для сети на рис. 18.17 Граф достижимости показывает, что сеть Петри ограничена и безопасна. Дерево достижимости позволяет легко проверить, является ли сеть сохраняющей или нет. Для этого достаточно сложить количество фишек 184 для каждой маркировки (для каждой вершины дерева достижимости) и сравнить суммы. Если они одинаковы, то сеть сохраняющая. Сеть на рис. 18.13 не является сохраняющей. Анализ дерева достижимости для ограниченной сети Петри позволяет определить все его свойства, перечисленные выше. Если сеть неограниченна (в разметке дерева достижимости присутствует символ ∞ ), то анализ сети затруднен. Наличие символа ∞ означает потерю части информации о свойствах сети. Основным недостатком сетей Петри являются трудности моделирования времени работы (времени задержки) переходов. Поэтому классические сети Петри позволяют эффективно моделировать логику работы программ и технических устройств, но имитационное моделирование временных характеристик в них затруднено. Поэтому разработаны различные варианты расширения возможностей сетей Петри, приспособленных к решению конкретных задач. 185 Литература 1. Худсон Д. Статистика для физиков. М: Из-во МИР, 1970 г. - 295 с. 2. Бородин А.Н. Элементарный курс теории вероятностей и математической статистики. СПб: Из-во Лань, 1999 г. - 224 с. 3. Розанов Ю.А. Случайные процессы. М.: Издательство ПАУКА, 1971г. 286 с. 4. Розанов Ю.А. Лекции по теории вероятностей. М.: Издательство НАУКА, 1968 г. - 120 с. 5. Кемени Дж., Снелл Дж. Конечные цепи Маркова. М. Издательство НАУКА, 1970г. - 271 с. 6. Гнеденко Б.В., Коваленко И.Н. Введение в теорию массового обслуживания. М.: Издательство НАУКА, 1966 г. - 206 с. 7. Матвеев В.Ф., Ушаков В.Г. Системы массового обслуживания. Издательство московского университета, 1984 г. - 239 с. 8. Гультяев А.К. Имитационное моделирование в среде Windows. С-Пб Издательство КОРОНА принт, 2001 г. - 400 с. 9. Стренг Г., Фикс Дж. Теория метода конечных элементов. М.: Мир, 1977 г. - 349 с. 10. Бухарева Р.Г. Вероятностные автоматы. Казань, 1970 г. -235 с. 11. Котов В.Е. Сети Петри. — М.: Издательство Наука, 1984. - 160 с. 12.Питерсон Дж., Теория сетей Петри и моделирование систем. Пер. с англ.М.: Мир, 1984. - 264 с., ил. 13. А.Л. Бекларян. Фронт выхода в модели поведения толпы при чрезвычайных ситуациях. Вестник ТГУ, т.20, вып.4,2015 – стр. 28 -32. 186