Модели и методы анализа проектных решений

№ 5267 Л.А. Гладков Н.В. Гладкова Модели и методы анализа проектных решений МИНИСТЕРСТВО ОБРАЗОВАНИЯ И НАУКИ РОССИЙСКОЙ ФЕДЕРАЦИИ Федеральное государственное автономное образовательное учреждение высшего профессионального образования «Южный федеральный университет» Л.А. Гладков, Н.В. Гладкова Модели и методы анализа проектных решений КОНСПЕКТ ЛЕКЦИЙ Таганрог 2014 УДК: 519.8(075) Рецензенты: зав. кафедрой компьютерных образовательных технологий СанктПетербургского государственного университета информационных технологий, механики и оптики, доктор технических наук, профессор Лисицина Л.С.; кандидат технических наук, доцент кафедры САиТ ЮФУ Лебедев В.Б. Гладков Л.А., Гладкова Н.В. Модели и методы анализа проектных решений: Конспект лекций. – Таганрог: Изд-во ЮФУ, 2014. – 105 с. Рассмотрены общие сведения о процессе проектирования. Приведена классификация математических моделей, приведены основные требования к математическим моделям. Дано описание блочно-иерархического подхода к проектированию. Описаны основные математические модели элементов и систем на микро-, макро- и метауровнях. Подробно описаны основные методы получения и анализа математических моделей объектов проектирования. Конспект предназначен для студентов вузов, обучающихся по направлениям «Информатика и вычислительная техника» и «Информационные системы». Учебное пособие может быть полезным для специалистов, занятых разработкой интеллектуальных САПР, поддержки и принятия решений, новых информационных технологий в науке, технике, образовании, бизнесе и экономике. Ил. 21. Табл. 4. Библиогр.: 9 назв.  Л.А. Гладков, 2014  Н.В. Гладкова, 2014  ЮФУ, 2014 3 ВВЕДЕНИЕ Проектирование сложных технических объектов в машиностроении, авиационной и автомобильной промышленности, электронновычислительных и радиоэлектронных устройств и других отраслях является одной из наиболее значимых и актуальных задач, стоящих перед человечеством. Корректное и в то же время эффективное решение данной задачи позволяет получить огромный экономический эффект, позволяет развивать и внедрять новые технологии и материалы. В процессе проектирования инженерам постоянно приходится решать задачу моделирования проектируемого объекта, производить выбор из множества вариантов, анализировать получаемую информацию. От эффективности решения этих задач напрямую зависит результативность проектирования. Учебная дисциплина «Модели и методы принятия решений» посвящена рассмотрению наиболее известных и популярных методов построения и анализа различных видов моделей реальных объектов проектирования из различных областей техники. Конспект лекций является составной частью учебного курса «Модели и методы принятия решений» наряду с учебно-методическими пособиями по проведению лабораторных и практических занятий по данному учебному курсу. Он предназначен для студентов очной и заочной форм обучения, проходящих подготовку по направлениям «Информатика и вычислительная техника» и «Информационные системы». При подготовке конспекта лекций был использован классический учебник Корячко В.П., Курейчик В.M., Норенков И.П. «Теоретические основы САПР: Учебник для вузов». – M.: Энергоатомиздат, 1987. Также были использованы другие учебники и учебные пособия, изданные ранее. Список рекомендуемой для изучения литературы находится в конце данного конспекта. Составители данного конспекта лекций выражают благодарность всем авторам, упомянутым в библиографическом списке, студентам и сотрудникам кафедры систем автоматизированного проектирования Южного федерального университета, помогавшим в подготовке к изданию данного конспекта и будут признательны всем высказавшим свои замечания и дополнения. 4 ГЛАВА 1. МАТЕМАТИЧЕСКИЕ МОДЕЛИ ПРОЕКТИРУЕМЫХ ОБЪЕКТОВ 1.1. ОБЩИЕ СВЕДЕНИЯ О ПРОЕКТИРОВАНИИ При построении новых объектов техники по заданному описанию несуществующего объекта выполняется его материализация в работоспособную надежную конструкцию. Проектирование – это процесс создания описания, необходимого для построения в заданных условиях еще не существующего объекта, на основе первичного описания этого объекта. Процесс создания описания нового объекта может выполняться различными способами. Выделим три основных. Если весь процесс проектирования осуществляет человек, то проектирование называют неавтоматизированным. В настоящее время неавтоматизированное проектирование таких сложных объектов, как электронно-вычислительная аппаратура (ЭВА), практически не применяют. Наибольшее распространение получило проектирование, при котором происходит взаимодействие человека и ЭВМ. Такое проектирование называют автоматизированным. Автоматизированное проектирование, как правило, осуществляется в режиме диалога человека с ЭВМ на основе применения специальных языков общения с ЭВМ. Проектирование, при котором все преобразования описаний объекта и алгоритма его функционирования осуществляются без участия человека, называют автоматическим. Рассмотрим ряд понятий, которые используются, например, при разработке ЭВА. Первичное описание ЭВА, представленное в заданной форме, называется заданием на проектирование. В задании на проектирование ЭВА должны быть сведения о назначении ЭВА, ее параметрах, способах функционирования, конструктивной реализации, изготовления и т. п. Проектным решением называется промежуточное или конечное описание объекта проектирования, необходимое и достаточное для рассмотрения и определения дальнейшего направления или окончания проектирования. Проектное решение или их совокупность, удовлетворяющие заданным требованиям, необходимые для создания объекта проектирования, будут являться результатом проектирования. В заданные требования должны 5 быть обязательно включены требования к форме представляемого проектного решения. Документ, выполненный по заданной форме, в котором представлено какое-либо проектное решение, полученное при проектировании, называется проектным. Совокупность проектных документов в соответствии с установленным перечнем, в котором представлен результат проектирования, называется проектом. Под проектной процедурой понимают формализованную совокупность действий, выполнение которых оканчивается проектным решением. Например, проектными процедурами являются оптимизация, контроль, поиск решения, корректировка, компоновка, проверка правильности трассировки и т.п. Действие или формализованная совокупность действий, составляющих часть проектной процедуры, алгоритм которых остается неизменным для ряда проектных процедур, называется проектной операцией. Примерами проектных операций являются составление таблиц с данными вычисления, вычерчивание топологии, ввод и вывод данных и т.п. Соответственно проектная процедура, алгоритм которой остается неизменным для различных объектов проектирования или различных стадий проектирования одного и того же объекта, называется унифицированной проектной процедурой. 1.2. КЛАССИФИКАЦИЯ МАТЕМАТИЧЕСКИХ МОДЕЛЕЙ Математическая модель (ММ) технического объекта есть совокупность математических объектов (чисел, переменных, матриц, множеств и т.п.) и отношений между ними, которая адекватно отображает свойства технического объекта, интересующие инженера, разрабатывающего этот объект. Выполнение проектных операций и процедур в САПР основано на оперировании ММ. С их помощью прогнозируются характеристики и оцениваются возможности предложенных вариантов схем и конструкций, проверяется их соответствие предъявляемым требованиям, проводится оптимизация параметров, разрабатывается техническая документация и т.п. В САПР для каждого иерархического уровня сформулированы основные положения математического моделирования, выбран и развит 6 соответствующий математический аппарат, получены типовые ММ элементов проектируемых объектов, формализованы методы получения и анализа математических моделей систем. Сложность задач проектирования и противоречивость требований высокой точности, полноты и малой трудоемкости анализа обусловливают целесообразность компромиссного удовлетворения этих требований с помощью соответствующего выбора моделей. Это обстоятельство приводит к расширению множества используемых моделей и развитию алгоритмов адаптивного моделирования. Функциональные и структурные модели. В проектных процедурах, связанных с функциональным аспектом проектирования, как правило, используются ММ, отражающие закономерности процессов функционирования объектов. Такие модели называют функциональными. Типичная функциональная модель представляет собой систему уравнений, описывающих либо электрические, тепловые, механические процессы, либо процессы преобразования информации. В то же время в процедурах, относящихся к конструкторскому аспекту проектирования, преобладает использование математических моделей, отражающих только структурные свойства объекта, например его геометрическую форму, размеры, взаимное расположение элементов в пространстве. Такие модели называют структурными. Структурные модели чаще всего представляются в виде графов, матриц инцидентности и смежности, списков и т.п. Как правило, функциональные модели более сложные, поскольку в них отражаются также сведения о структуре объектов. Однако при решении многих задач конструирования использование сложных функциональных моделей не оправдано, так как нужные результаты могут быть получены на основе более простых структурных моделей. Функциональные модели применяют преимущественно на завершающих этапах верификации описаний объектов, предварительно синтезированных с помощью структурных моделей. Иерархия математических моделей в САПР. Блочно-иерархический подход к проектированию радиоэлектронной аппаратуры (РЭА) включает в качестве своей основы иерархию математических моделей. Деление моделей по иерархическим уровням (уровням абстрагирования) происходит по степени детализации описываемых свойств и процессов, протекающих в объекте. При этом на каждом иерархическом уровне используют свои 7 понятия «система» и «элементы». Так, система k-гo уровня рассматривается как элемент на соседнем более высоком (k – 1) - м уровне абстрагирования. Представим структуру некоторого объекта в виде множества элементов (рис. 1) и связей между ними. Выделим в соответствии с блочноиерархическим подходом в структуре объекта некоторые подмножества элементов и назовем их блоками (на рисунке показаны штриховыми линиями). Пусть состояние каждой связи характеризуется одной фазовой переменной vi, zj или uk. Здесь vi относится к внутренним связям между элементами данного блока, zj и uk относятся к выходам и входам блока соответственно. Теперь поясним важные для функциональных моделей понятия полной модели и макромодели. U1 A 1 V1 V2 2 U3 3 V4 U2 4 V5 7 V3 5 6 Z1 V7 V6 V8 8 9 Z2 Z3 Рис. 1. Представление структуры объекта Полная модель блока есть модель, составленная из моделей элементов с учетом межэлементных связей, т.е. модель, описывающая как состояние выходов, так и состояние каждого из элементов блока. Моделями элементов блока А являются уравнения, связывающие входные и выходные переменные: f1 ( v1 , u1 )  0;   f 2 ( v1 , v 2 )  0;  f3 ( v2 , u3 , v4 )  0; .   f9 ( v9 , v8 , z 2 )  0.  Полная модель блока есть система уравнений (1) 8 F(V, U) = 0, Z = (V, U), (2) где V, Z и U – векторы внутренних, выходных и входных фазовых переменных блока. При большом количестве элементов размерность вектора V и порядок системы уравнений (2) становятся чрезмерно большими и требуются упрощения. При переходе к более высокому иерархическому уровню упрощения основаны на исключении из модели вектора внутренних переменных V. Полученная модель представляет собой систему уравнений (Z, U) = 0 (3) существенно меньшей размерности, чем полная модель (2), и называется макромоделью. Следовательно, макромодель уже не описывает процессы внутри блока, а характеризует только процессы взаимодействия данного блока с другими в составе системы блоков. Модели (2) и (3) относятся друг к другу как полная модель и макромодель на n-м уровне иерархии. На более высоком (п – 1)-м уровне блок A рассматривается как элемент и макромодель (3) становится моделью элемента А. Следовательно, модели (1) и (3) относятся друг к другу как модели элементов соседних иерархических уровней. Из моделей типа (3) может быть составлена полная модель системы на (п – 1)-м уровне. Микро-, макро- и метауровни. В зависимости от сложности объекта при его проектировании используют большее или меньшее число уровней абстракции. Объединение уровней, родственных по характеру используемого математического аппарата, приводит к образованию трех укрупненных уровней – микро-, макро- и метауровня – в иерархии функциональных моделей для большинства проектируемых сложных объектов. На микроуровне используют математические модели, описывающие физическое состояние и процессы в сплошных средах. Для моделирования применяют аппарат уравнений математической физики. Примерами таких уравнений служат дифференциальные уравнения в частных производных – уравнения электродинамики, теплопроводности, упругости, газовой динамики. Эти уравнения описывают поля электрического потенциала и температуры в полупроводниковых кристаллах интегральных схем, напряженно-деформированное состояние деталей механических 9 конструкций и т.п. К типичным фазовым переменным на микроуровне относятся электрические потенциалы, давления, температуры, концентрации частиц, плотности токов, механические напряжения и деформации. Независимыми переменными являются время и пространственные координаты. В качестве операторов F и  в уравнениях (2) фигурируют дифференциальные и интегральные операторы. Уравнения (2), дополненные краевыми условиями, составляют ММ объектов на микроуровне. Анализ таких моделей сводится к решению краевых задач математической физики. На макроуровне производится дискретизация пространства с выделением в качестве элементов отдельных деталей, дискретных электрорадиоэлементов, участков полупроводниковых кристаллов. При этом из числа независимых переменных исключают пространственные координаты. Функциональные модели на макроуровне представляют собой системы алгебраических или обыкновенных дифференциальных уравнений, для их получения и решения используют соответствующие численные методы. В качестве фазовых переменных фигурируют электрические напряжения, токи, силы, скорости, температуры, расходы и т.д. Они характеризуют проявления внешних свойств элементов при их взаимодействии между собой и внешней средой в электронных схемах или механических конструкциях. На метауровне с помощью дальнейшего абстрагирования от характера физических процессов удается получить приемлемое по сложности описание информационных процессов, протекающих в проектируемых объектах. На метауровне для моделирования аналоговой РЭА широко применяют аппарат анализа систем автоматического управления, а для моделирования цифровой РЭА – математическую логику, теорию конечных автоматов, теорию массового обслуживания. Математические модели на метауровне – системы обыкновенных дифференциальных уравнений, системы логических уравнений, имитационные модели систем массового обслуживания. Формы представления моделей. Для представления моделей используют следующие основные формы: Инвариантная форма – запись соотношений модели с помощью традиционного математического языка безотносительно к методу решения уравнений модели. 10 Алгоритмическая форма – запись соотношений модели и выбранного численного метода решения в форме алгоритма. Аналитическая форма – запись модели в виде результата аналитического решения исходных уравнений модели; обычно модели в аналитической форме представляют собой явные выражения выходных параметров как функций внутренних и внешних параметров. Схемная форма, называемая также графической формой, – представление модели на некотором графическом языке, например на языке графов, эквивалентных схем, диаграмм и т.п. Графические формы удобны для восприятия человеком. Использование таких форм возможно при наличии правил однозначного истолкования элементов чертежей и их перевода на язык инвариантных или алгоритмических форм. Модели в алгоритмической и аналитической формах называют соответственно алгоритмическими и аналитическими. Среди алгоритмических моделей важный класс составляют имитационные модели, предназначенные для имитации физических или информационных процессов в объекте при задании различных зависимостей входных воздействий от времени. Собственно имитацию названных процессов называют имитационным моделированием. Результат имитационного моделирования – зависимости фазовых переменных в избранных элементах системы от времени. Примерами имитационных моделей являются модели электронных схем в виде систем обыкновенных дифференциальных уравнений или модели систем массового обслуживания, предназначенные для имитации процессов прохождения заявок через систему. 1.3. ТРЕБОВАНИЯ К МАТЕМАТИЧЕСКИМ МОДЕЛЯМ Основными требованиями, предъявляемыми к математическим моделям, являются требования адекватности, универсальности и экономичности. Адекватность. Модель считается адекватной, если отражает заданные свойства объекта с приемлемой точностью. Точность определяется как степень совпадения значений выходных параметров модели и объекта. Пусть j – относительная погрешность модели по j-му выходному параметру:  j  ( ~y j  y j ) / y j , 11 где ~y j – j-й выходной параметр, рассчитанный с помощью модели; y j – тот же выходной параметр, имеющий место в моделируемом объекте. Погрешность модели m по совокупности учитываемых выходных параметров оценивается одной из норм вектора m = ( 1,  2, …, m), например  M  max |  j | или  M  j 1, , m n   2j . j 1 Точность модели различна в разных условиях функционирования объекта. Эти условия характеризуются внешними параметрами. Если задаться предельной допустимой погрешностью пред, то можно в пространстве внешних параметров выделить область, в которой выполняется условие m < пред. Эту область называют областью адекватности (ОА) модели. Возможно введение индивидуальных предельных значений пред для каждого выходного параметра и определение ОА как области, в которой одновременно выполняются все m условий вида |j |≥ пред. Пример ОА (заштрихована) в двумерном пространстве дан на рис. 2. Здесь qk – k-й внешний параметр. q ОА Рис. 2. Пример области адекватности Определение областей адекватности для конкретных моделей – сложная процедура, требующая больших вычислительных затрат. Эти затраты и трудности представления ОА быстро растут с увеличением размерности пространства внешних параметров. Определение ОА – более трудная задача, чем, например, задача параметрической оптимизации, поэтому для моделей вновь проектируемых объектов ОА не рассчитывают. Однако для моделей унифицированных элементов расчет областей адекватности становится оправданным в связи с однократностью определения ОА и многократностью их использования при проектировании 12 различных систем. Знание ОА позволяет правильно выбирать модели элементов из числа имеющихся и тем самым повышать достоверность результатов машинных расчетов. В общем случае ОА может иметь произвольную форму, сведения о которой выражаются громоздко, и неудобна в использовании, поэтому на практике вместо истинных ОА применяют те или иные их аппроксимации. Наиболее просто представляются и используются сведения об областях, имеющих формулу гиперпараллелепипеда, который задается р двусторонними неравенствами: qk  qk  qk , k  1, p, где р – размерность пространства внешних параметров. В библиотеку моделей элементов наряду с алгоритмом, реализующим модель, и номинальными значениями параметров должны включаться граничные значения внешних параметров q'k и q''k, задающие область адекватности. На рис. 2 дано графическое представление области адекватности и аппроксимирующего ее гиперпараллелепипеда. Такое представление удобно для двумерных случаев. Возможно использование и других аппроксимаций ОА, например областей с линеаризованными границами в виде участков гиперплоскостей, областей в форме гиперсфер и т.п. Универсальность. При определении ОА необходимо выбрать совокупность внешних параметров и совокупность выходных параметров yj, отражающих учитываемые в модели свойства. Типичными внешними параметрами при этом являются параметры нагрузки и внешних воздействий (электрических, механических, тепловых, радиационных и т.п.). Увеличение числа учитываемых внешних факторов расширяет применимость модели, но существенно удорожает работу по определению ОА. Выбор совокупности выходных параметров также неоднозначен, однако для большинства объектов число и перечень учитываемых свойств и соответствующих им выходных параметров сравнительно невелики, достаточно стабильны и составляют типовой набор выходных параметров. Например, для макромоделей логических элементов БИС такими выходными параметрами являются уровни выходного напряжения в состояниях логических «0» и «1», запасы помехоустойчивости, задержка 13 распространения сигнала, рассеиваемая мощность. Если адекватность характеризуется положением и размерами ОА, то универсальность модели определяется числам и составом учитываемых в модели внешних и выходных параметров. Экономичность. Экономичность модели характеризуется затратами вычислительных ресурсов для ее реализации, а именно затратами машинного времени ТM и памяти ПМ. Общие затраты ТM и ПМ на выполнение в САПР какой-либо проектной процедуры зависят как от особенностей выбранных моделей, так и от методов решения. В большинстве случаев при реализации численного метода происходят многократные обращения к модели элемента, входящего в состав моделируемого объекта. Тогда удобно экономичность модели элемента характеризовать затратами машинного времени, получающимися при обращении к модели, а число обращений к модели должно учитываться при оценке экономичности метода решения. Экономичность модели по затратам памяти оценивается объемом оперативной памяти, необходимой для реализации модели. Требования широких областей адекватности, высокой степени универсальности, с одной стороны, и высокой экономичности, с другой, являются противоречивыми. Наилучшее компромиссное удовлетворение этих требований оказывается неодинаковым в различных применениях. Это обстоятельство обусловливает использование в САПР многих моделей для объектов одного и того же типа – различного рода макромоделей, многоуровневых, смешанных моделей и т.п. 1.4. МЕТОДЫ ПОЛУЧЕНИЯ МОДЕЛЕЙ ЭЛЕМЕНТОВ Получение моделей элементов (моделирование элементов) в общем случае – процедура неформализованная. Основные решения, касающиеся выбора вида математических соотношений, характера используемых переменных и параметров, принимает проектировщик. В то же время такие операции, как расчет численных значений параметров модели, определение областей адекватности и другие, алгоритмизированы и решаются на ЭВМ. Поэтому моделирование элементов обычно выполняется специалистами конкретных технических областей с помощью традиционных средств экспериментальных исследований и средств САПР. Методы получения функциональных моделей элементов делят на 14 теоретические и экспериментальные. Теоретические методы основаны на изучении физических закономерностей протекающих в объекте процессов, определении соответствующего этим закономерностям математического описания, обосновании и принятии упрощающих предположений, выполнении необходимых выкладок и приведении результата к принятой форме представления модели. Экспериментальные методы основаны на использовании внешних проявлений свойств объекта, фиксируемых во время эксплуатации однотипных объектов или при проведении целенаправленных экспериментов. Несмотря на эвристический характер многих операций моделирования, имеется ряд положений и приемов, общих для получения моделей различных объектов. Достаточно общий характер имеют методика макромоделирования, математические методы планирования экспериментов, а также алгоритмы формализуемых операций расчета численных значений параметров и определения областей адекватности. Методика макромоделирования. Применение методики состоит из следующих этапов: Определение тех свойств объекта, которые должны отражаться моделью (устанавливаются требования к степени универсальности будущей модели). Сбор априорной информации о свойствах моделируемого объекта. Примерами собираемых сведений могут служить справочные данные, математические модели и результаты эксплуатации существующих аналогичных объектов и т.п. Получение общего вида уравнений модели (структуры модели). Этот этап в случае теоретических методов включает выполнение всех присущих этим методам операций, перечисленных выше. Часто проектировщику модели удобнее оперировать не уравнениями, а эквивалентными схемами, с помощью которых инженеру проще устанавливать физический смысл различных элементов математической модели. Определение численных значений параметров модели. Возможны следующие приемы выполнения этого этапа: а) использование специфических расчетных соотношений с учетом собранных на этапе 2 сведений; б) решение экстремальной задачи, в которой в качестве целевой функции выбирается степень совпадения известных значений выходных параметров объекта с результатами использования модели, а управляемыми 15 параметрами являются параметры модели; в) проведение экспериментов и обработка полученных результатов. Оценка точности полученной модели и определение области ее адекватности. При неудовлетворительных точностных оценках выполняют итерационное приближение к желаемому результату повторением этапов 3 – 5. Представление полученной модели в форме, принятой в используемой библиотеке моделей. Методы планирования экспериментов. Для целей моделирования используют пассивные и активные эксперименты. В пассивных экспериментах нет возможности выбирать условия опыта по своему усмотрению и устанавливать значения факторов на желаемом уровне. В активных экспериментах опыты проводятся по заранее разработанному плану, выражающему количество опытов и значения факторов в каждом опыте. Выбор вида зависимости выходного параметра макромодели у (в общем случае рассматривается вектор выходных параметров Y) от внешних параметров qk, объединенных в вектор факторов Q, осуществляется проектировщиком. Чаще всего в методах планирования эксперимента используются модели линейные y = AQ (4) или квадратичные у = АВ (5) где А – вектор-строка коэффициентов (параметров) модели; В – вектор, включающий факторы qk, те или иные произведения из двух, трех или более факторов и возможно также квадраты факторов qk2; k  1, p ; р – число факторов. Число опытов N, как правило, должно превышать число определяемых параметров вектора А. Параметры рассчитывают по методу наименьших квадратов, т.е. из условия минимизации суммы квадратов отклонений значений ~yl , определенных по уравнению модели (4), и измеренных значений yl . N 2 min  ( ~yl ( A)  yl ) , A где l – номер опыта. l 1 16 В зависимости от способов планирования преимущества активных экспериментов перед пассивными могут выражаться в получении оптимального положения области адекватности, в ее увеличенном объеме, в упрощении оценок точности и т.п. Регрессионный анализ. Связь между у и Q может быть не функциональной, а статистической, что особенно характерно при пассивных экспериментах. Для получения моделей в такой ситуации часто применяют регрессионный анализ. Модель ищется в форме уравнения регрессии (4), в котором роль коэффициентов ak в векторе А выполняют коэффициенты относительной регрессии. Рассмотрим алгоритм вычисления коэффициентов ak. По результатам пассивных экспериментов получаются оценки математических ожиданий Му, Мk, среднеквадратичных отклонений σу, σk соответственно для выходного у и внешних qk параметров, а также коэффициенты корреляции rk между у и qk, образующие вектор R, и коэффициенты корреляции dkj между факторами qk и qj, образующие матрицу D. Далее решается система линейных алгебраических уравнений D = R (6) и полученный вектор  = (1, 2, …, p) используется при расчете относительных коэффициентов регрессии по формуле ak = ky/k. Если факторы qk некоррелированы, то D – единичная матрица и можно обойтись без решения системы (6), так как k = rk. Диалоговое моделирование. Наличие в методике макромоделирования эвристических и формальных операций обусловливает целесообразность разработки моделей элементов в диалоговом режиме работы с ЭВМ. Язык взаимодействия человека с ЭВМ должен позволять оперативный ввод исходной информации о структуре модели, об известных характеристиках и параметрах объекта, о плане экспериментов. Диалоговое моделирование должно иметь программное обеспечение, в котором реализованы алгоритмы статистической обработки результатов экспериментов, расчета выходных параметров эталонных моделей и создаваемых макромоделей, в том числе расчета параметров по методам планирования экспериментов и регрессионного анализа, алгоритмы методов поиска экстремума, расчета областей адекватности и др. Пользователь, разрабатывающий модель, может 17 менять уравнения модели, задавать их в аналитической, схемной или табличной форме, обращаться к нужным подпрограммам и тем самым оценивать результаты предпринимаемых действий, приближаясь к получению модели с требуемыми свойствами. 1.5. МАТЕМАТИЧЕСКИЕ МОДЕЛИ ОБЪЕКТОВ ПРОЕКТИРОВАНИЯ, ИСПОЛЬЗУЕМЫЕ НА МИКРОУРОВНЕ Математические модели деталей и процессов на микроуровне отражают физические процессы, протекающие в сплошных средах и непрерывном времени. Независимыми переменными в этих моделях являются пространственные координаты и время. В качестве зависимых переменных выступают фазовые переменные, такие как потенциалы, напряженности полей, концентрации частиц, деформации и т.п. Взаимосвязи переменных выражаются с помощью уравнений математической физики – интегральных, интегродифференциальных или дифференциальных уравнений в частных производных. Эти уравнения составляют основу ММ на микроуровне, Для получения законченной математической модели, используемой в задачах проектирования, необходимо дополнительно выполнить ряд процедур:  выбрать краевые условия. Краевые условия представляют собой сведения о значениях фазовых переменных и (или) их производных на границах рассматриваемых пространственных и временных областей;  дискретизировать задачу. Дискретизация подразумевает разделение рассматриваемых пространственных и временных областей на конечное число элементарных участков с представлением фазовых переменных конечным числом значений в избранных узловых точках, принадлежащих элементарным участкам;  алгебраизировать задачу – аппроксимировать дифференциальные и интегральные уравнения алгебраическими. Используют два основных подхода к дискретизации и алгебраизации краевых задач, составляющие сущность методов конечных разностей (МКР) и конечных элементов (МКЭ). С помощью любого из этих методов формируется окончательная модель, исследуемая при выполнении различных процедур анализа проектируемого объекта. Пользователь САПР средствами входного языка задает исходную информацию о конфигурации проектируемого объекта, о способе дискретизации – разделения среды на элементы, о физических свойствах 18 участков среды. Формирование модели объекта, т.е. разделение среды на элементы, выбор математических моделей элементов из заранее составленных библиотек, объединение моделей элементов в общую систему уравнений, так же как и решение получающихся уравнений, осуществляется автоматически на ЭВМ. Основные уравнения математической физики, используемые в моделях проектируемых объектов. Процессы, протекающие в техническом объекте при его функционировании, по своей физической природе могут быть разделены на электрические, тепловые, магнитные, оптические, механические, гидравлические и т.п. Каждому типу процессов в математической модели соответствует своя подсистема, основанная на определенных уравнениях математической физики. Например, электрические процессы в современных полупроводниковых приборах с достаточной точностью удается описать с помощью уравнений непрерывности и Пуассона. Уравнения непрерывности выражают скорости изменения концентраций свободных носителей заряда и записываются отдельно для дырок и электронов: p 1   div J p  g p ; t q (7) n 1  div J n  g n . t q (8) где р и п – концентрации дырок и электронов соответственно; q – заряд электрона; gp и gn – скорости процесса генерации-рекомбинации соответственно дырок и электронов; Jp = q(–p p grad  + Dp grad p); (9) Jn = q(–n n grad  + Dn grad n). (10) Jp, Jn – плотности дырочного и электронного токов;  р,  n – подвижности; Dp, Dn – коэффициенты диффузии дырок и электронов; φ – электрический потенциал. Уравнения (7) – (8) показывают, что причинами изменения концентрации носителей могут быть неодинаковость числа носителей, втекающих (и вытекающих) в элементарный объем полупроводника (тогда div J  0), и нарушение равновесия между процессами генерации и рекомбинации носителей. Уравнения (9) и (10), называемые уравнениями плотности тока, характеризуют причины протекания электрического тока в полупроводнике: электрический дрейф под воздействием электрического 19 поля (grad   0) и диффузию носителей при наличии градиента концентрации. Уравнение Пуассона характеризует зависимость изменений в пространстве напряженности электрического поля Е = – grad φ от распределения плотности электрических зарядов ρ: div E =  / (0), где  – относительная диэлектрическая проницаемость среды; ε0 – диэлектрическая постоянная. В качестве краевых условий в моделях полупроводниковых приборов используют зависимости потенциалов на контактах от времени. Значения концентраций носителей на границе между внешним выводом и полупроводником принимают равными равновесным концентрациям p0 и n0. Для границ раздела полупроводника и окисла задаются скоростью поверхностной рекомбинации gs, что определяет величины нормальных к поверхности раздела составляющих плотностей тока Jp и Jn, и т.д. Результат решения уравнений непрерывности и Пуассона при известных краевых условиях – это поля потенциала и концентраций подвижных носителей в различных областях полупроводниковой структуры. Знание этих полей позволяет оценить электрические параметры прибора. В основе моделей диффузионных процессов, используемых, в частности, для описания технологических операций диффузии примесей при изготовлении интегральных схем и полупроводниковых приборов, лежит уравнение диффузии N  div ( D grad N ) , t где N – концентрация примеси; D – коэффициент диффузии. Краевые условия представлены зависимостью распределения примеси N в объеме полупроводника в начальный момент времени и зависимостью поверхностной концентрации от времени. Основные положения метода конечных разностей. Рассмотрим уравнение Lv(Х) = f(X), (12) где L – дифференциальный оператор; v(Х) – зависимая фазовая переменная; X = (x1, x2, x3) – вектор независимых (пространственных) координат; f(X) – заданная функция. Пусть X  R (R – ограниченная область изменения независимых переменных), заданы граничные условия. 20 Дискретизация задачи заключается в покрытии R сеткой и замене множества R конечным множеством точек Хk, являющихся узлами сетки. Сетка может быть прямоугольной, косоугольной, с постоянными или переменными межузловыми расстояниями вдоль координатных осей (величинами шагов). Наиболее часто используют прямоугольную сетку с постоянными величинами шагов. На рис. 3 представлен фрагмент такой сетки для двумерной задачи с величинами шагов h1 и h2 вдоль координатных осей X1 и X2. x1 h2 h1 x2 Рис. 3. Дискретизация пространства с помощью сетки Алгебраизация задачи заключается в замене дифференциального оператора Lv разностным. Это означает, что непрерывная переменная v(X) заменяется конечным множеством значений vk = v(Xk) в узлах сетки, а производные dv / dX аппроксимируются конечноразностными выражениями. Применяется несколько способов выражения производных через значения vk. Вид разностных операторов удобно представлять графически в форме шаблонов. На рис. 4, (а – г) даны примеры шаблонов для одномерных, а на рис. 4, (д, ж) – для двумерных стационарных задач. d h dx xi 1 xi+1 d h dx 1 -1 xi-1 а 2h xi -1 xi-1 d dx б в 1 2 d h2 dx 2 -2 -1 xi-1 1 xi+1 xi h2 2 1 2h2 4 x1,i-1 x1,i x2;j+1 x2;j 1 1 1 1 x1,i-1 x1,i x1,i+1 x2;j x2;j-1 д 1 2 x2;j+1 4 г 1 xi 1 xi+1 x1,i+1 x2;j-1 -1 2 4h2 xy 1 x2;j+1 x2;j -1 1 x1,i-1 x1,i x1,i+1 x2;j-1 е ж Рис. 4. Примеры шаблонов для разностных операторов 21 Шаблон представляет собой часть сетки, включающую множество узлов Хk, значения переменных в которых используются при аппроксимации производных в заданном узле X*. Узлы Хk на рис. 4 показаны темными кружками, а узел X* обведен дополнительной окружностью. В левой части рисунков указан аппроксимируемый дифференциальный оператор, а рядом с узлами сетки записаны значения коэффициентов, с которыми соответствующие величины vk входят в конечноразностное аппроксимирующее выражение. При этом принято h1 = h2 = h. Так, для рис. 4, а) получаем h v  vi1  vi , x x* для рис. 4, г) h2  2v  vi1  2vi  vi1 , x 2 x* (13) для рис. 4, д) h 2 2 v ( x1i ,x2 j )  vij1  vij1  vij 1  vij 1  4vij (14) и т.д. В последнем выражении использовано обозначение vij = v(x1i, x2i). Подстановка выбранных аппроксимаций производных в исходное уравнение (12) преобразует его в систему разностных уравнений F(V) = 0, (15) где V – вектор, элементами которого являются значения фазовой переменной во внутренних узлах сетки. Аналогичная подстановка в исходное нестационарное уравнение (16)  / t  L  f(X, t) преобразует его в систему обыкновенных дифференциальных уравнений (ОДУ) dV / dt = F(V), (17) которая вместе с заданными начальными условиями представляет собой задачу Коши. Следовательно, дискретизация и алгебраизация уравнений в МКР сводит задачу анализа моделей на микроуровне к численному решению систем конечных (15) или обыкновенных дифференциальных (17) уравнений. Следует отметить, что точность аппроксимации растет с уменьшением величин шагов, однако при этом увеличивается порядок систем уравнений (15) или (17). Так, если окажется, что для достижения 22 приемлемой точности рассматриваемую область R нужно делить вдоль каждой из координатных осей на 102 участков, то порядки систем уравнений (15) или (17) в одно-, дву- и трехмерных задачах составляют соответственно около 102, 104 и 106. Очевидно, что решение двумерных и особенно трехмерных задач требует значительных вычислительных ресурсов и тщательного отбора соответствующего математического обеспечения. Методы решения таких уравнений, применяемые в САПР, рассматриваются в следующей главе. Основные положения метода конечных элементов. Рассмотрим применение МКЭ к решению задачи Lv(X) = f(X). Дискретизация исследуемой пространственной области R в МКЭ осуществляется ее разделением на непересекающиеся подобласти – конечные элементы (КЭ). В одномерных задачах КЭ представляют собой отрезки линий, в двумерных имеют форму треугольников, прямоугольников, в трехмерных – тетраэдров, параллелепипедов и т. п. В пределах каждого КЭ выбирают конечное число узловых точек Xk. Непрерывную фазовую переменную v, фигурирующую в модели (12), заменяют конечным числом значений vk этой фазовой переменной в точках Xk. Возможности использования КЭ различной формы, размеров и пространственной ориентации обусловливают легкость дискретизации граничных условий при произвольной форме области R. Это обстоятельство – одно из основных преимуществ МКЭ перед МКР, объясняющее широкое применение конечноэлементных представлений при моделировании процессов в деталях сложной конфигурации. Второе отличие МКЭ от МКР заключается в способе алгебраизации дифференциальных уравнений Lv(X) = f(X). Если в МКР аппроксимируются производные dv/dX, то в МКЭ аппроксимируется решение v(X) некоторой функцией u(Х) с неопределенными коэффициентами. Решение исходной задачи получается путем вычисления этих коэффициентов. В свою очередь задача вычисления коэффициентов формулируется как задача минимизации функционала, характеризующего качество аппроксимации решения v(X) функцией u(Х), а эта задача сводится к решению системы алгебраических уравнений. Очевидно, что чем сложнее применяемые аппроксимирующие функции и чем шире класс этих функций, тем сложнее задача формализации метода и 23 его реализации в САПР. Особенностью МКЭ является выбор аппроксимирующих функций для каждого КЭ в отдельности. Малые размеры КЭ позволяют использовать простые аппроксимирующие функции, причем одного и того же типа для всех КЭ определенной формы. Обычно в качестве u(X) для отдельного КЭ применяют полиномы степени не выше третьей, например, в одномерном случае r u ( x)   ai x i . i 0 В МКЭ u(Х) представляют в форме r u ( x)   qii ( X ) , (18) i 0 где коэффициенты qi имеют вполне определенный физический смысл – это значения аппроксимирующей функции в узловых точках; φi(Х) – функции, называемые координатными (базисными, пробными или функциями формы); r – число узловых точек в конечном элементе. В общем случае аппроксимации вектора V(X) в m–мерном пространстве выражение (18) принимает вид U(X) = NQ, (19) где U – вектор размера m  l; Q – вектор размера (mr)  l; N – интерполяционная матрица размера т  (тr), ее элементами являются координатные функции. Важной составной частью МКЭ является выбор функционала, характеризующего качество используемой аппроксимации. Примером такого функционала может служить следующий: I   ( LU ( X )  f ( X )) 2 dR , (20) R где X  R, R – исследуемая область. Очевидно, что интеграл I, относящийся ко всей области R, получается суммированием интегралов, подобных (20), но относящихся к отдельным КЭ. Иногда математические модели объектов на микроуровне уже в своем исходном виде могут быть представлены в вариационной формулировке, т.е. в виде задачи минимизации функционала. Типичным примером таких моделей служат модели, описывающие статические напряженнодеформированные состояния деталей. В этих моделях в качестве минимизируемого функционала используется выражение полной потенциальной энергии 24 П 1 t  DdR  A, 2 R (21) а в качестве искомой экстремали функционала фигурирует вектор перемещений W(X), XR. Деформации εij связаны с перемещениями i, соотношениями, 1      ij   i  j , 2  x j xi  (22) что можно выразить в матричной форме записью   2 x 1   0   11       22   0  33  1       12  2  x  13   2      23   x3   0  2  x2  x1  x3  0   0     2  x3   0     x1     x2  1   , (2.29)  2 3  (23) или более лаконично  = SW, (24) где S – фигурирующая в (23) матрица–оператор дифференцирования. Подставляя (24) в (21) и заменяя вектор W его аппроксимацией U = NQ в соответствии с (19), получаем П  0,5 Q t B t DBQdR A, R где B = SN. Обозначим K   B t DBdR (25) R и назовем матрицей жесткости. Тогда П = 0,5QtKQ В соответствии с принципом Лагранжа дифференцируем П по вектору Q и приравниваем результат нулю, получаем систему алгебраических уравнений (26) KQ = P, где P = dA/dQ – вектор правых частей, называемый вектором нагрузок. Матрица жесткости K всей исследуемой детали составляется из матриц жесткости Kij отдельных КЭ. Матрицы Kij несут информацию о конфигурации и упругих свойствах материала конечных элементов и подсчитываются по формуле (25), в которой при этом под R понимается 25 подобласть, относящаяся к рассматриваемому КЭ. Система уравнений типа (26) получается при применении МКЭ к решению и других стационарных уравнений. При этом требуется минимизировать некоторый функционал, например (20), что также приводит к решению системы алгебраических уравнений I/Q = 0. В случае нестационарных уравнений основные положения МКЭ – деление на КЭ, подбор аппроксимирующей функции и минимизируемого функционала – по-прежнему применяют по отношению к пространственной области. Тогда получаем систему обыкновенных дифференциальных уравнений (ОДУ) dQ/dt = F(Q), решением которой являются зависимости фазовых переменных Q(t) в узловых точках от времени. 1.6. МАТЕМАТИЧЕСКИЕ МОДЕЛИ ОБЪЕКТОВ ПРОЕКТИРОВАНИЯ, ИСПОЛЬЗУЕМЫЕ НА МАКРОУРОВНЕ Компонентные и топологические уравнения. Для одного и того же объекта (детали) на микро- и макроуровнях используют разные математические модели. На микроуровне ММ должна отражать внутренние по отношению к объекту процессы, протекающие в сплошных средах. На макроуровне ММ того же объекта служит для отражения только тех его свойств, которые характеризуют взаимодействие этого объекта с другими элементами в составе исследуемой системы. Математические модели элементов на макроуровне получают одним из способов, рассмотренных в разд. 1.5. Математические модели систем (ММС) формируют из математических моделей элементов (ММЭ) с помощью методов, излагаемых ниже. Уравнения, входящие в ММЭ, называют компонентными. Наряду с компонентными уравнениями в ММС обязательно входят уравнения, отражающие способ связи элементов между собой в составе системы и называемые топологическими. Топологические уравнения могут выражать законы сохранения, условия неразрывности, равновесия и т.п. В используемых в САПР методах формирования ММС принято моделируемую систему представлять в виде совокупности физически однородных подсистем. Каждая подсистема описывает процессы 26 определенной физической природы, например механические, электрические, тепловые, гидравлические. Как правило, для описания состояния одной подсистемы достаточно применять фазовые переменные двух типов – потенциала и потока. В первых столбцах табл. 1 конкретизированы типы фазовых переменных применительно к ряду встречающихся подсистем. Таблица 1 Параметры простых Фазовая переменная элементов Типа Типа Подсистема потенциала C L R потока I U Механическая Скорость Сила Механическая вращательная Угловая скорость Вращательный момент Электрическая Электрическое напряжение Электрический ток Тепловая Температура Тепловой поток Гидравлическая и пневматическая Давление Расход Механическое Масса Гибкость сопротивление ВращательВращательМомент ное ная инерции сопротивгибкость ление ЭлектриЭлектриЭлектричеческая ческое ская индуктив- сопротивле емкость ность ние Тепловое Теплоем–––– сопротивле кость ние Гидравли- ГидравлиГидравлическая ческое ческая индуктив- сопротивемкость ность ление Особенностью топологических уравнений является то, что каждое из них связывает однотипные фазовые переменные, относящиеся к разным элементам системы. Примером могут служить уравнения законов Кирхгофа, записываемые относительно либо токов, либо напряжений ветвей. Для компонентных уравнений характерно то, что они связывают разнотипные фазовые переменные, относящиеся к одному элементу. Так, уравнение закона Ома связывает ток и напряжение резистора. Формы представления моделей. Элементы подсистем могут быть простыми и сложными. Элемент называют простым, если соответствующая ему ММЭ может 27 быть представлена в виде одного линейного уравнения, связывающего переменную типа потенциала U и переменную типа потока, характеризующие состояние данного элемента. В физически однородных подсистемах различают три типа простых элементов. Это элементы емкостного, индуктивного и резистивного типов. Соответствующие им ММЭ имеют вид CdU /dt = I, LdI / dt = U, U = RI, (27) где С, L, R – параметры элементов, физический смысл которых поясняется в табл. 1. Элементы подсистем в зависимости от числа однотипных фазовых переменных, входящих в ММЭ, делят на двухполюсники и многополюсники. Двухполюсник характеризуется парой переменных типа U и I, определяется так же, как простой элемент, если снять условие линейности уравнения. Многополюсник можно представить как совокупность взаимосвязанных двухполюсников. Для представления математических моделей на макроуровне применяют несколько форм. Инвариантная форма – представление модели в виде системы уравнений, записанной на общепринятом математическом языке, безотносительно к методу численного решения. Применительно к системам обыкновенных дифференциальных уравнений различают две инвариантные формы – нормальную и общую, определяемые тем, в каком виде – явном или неявном относительно вектора производных – представлена система. Ряд форм модели получается при преобразовании ее уравнений на основе формул и требований выбранного численного метода решения. Так, численное решение дифференциальных уравнений, как в частных производных, так и обыкновенных, требует их предварительного преобразования – дискретизации и алгебраизации. Дискретизация заключается в замене непрерывных независимых переменных (времени и пространственных координат) дискретным множеством их значений. Алгебраизованная форма – результат представления дифференциальных уравнений в полученных после дискретизации точках в алгебраизованном виде с помощью формул численного интегрирования. Ряд численных методов решения основан на линеаризации исходных уравнений. Линеаризованная форма модели – представление ее уравнений в линейном виде. Алгебраизация и линеаризации могут осуществляться по 28 отношению ко всем или только избранным переменным, уравнениям или их частям, что увеличивает разнообразие возможных форм представления моделей. Формы представления моделей определяются также используемыми языковыми средствами. Наряду с традиционным математическим языком применяют алгоритмические языки, а также те или иные графические изображения, облегчающие пользователю восприятие модели и приводящие к представлению модели в той или иной схемной форме, например представление моделей в виде эквивалентных схем, графов, к таким формам относится также представление разностных уравнений с помощью шаблонов (см. рис. 4). Рассмотрим особенности представления моделей в виде эквивалентных схем. В разных областях техники применяют специфические системы обозначений элементов на эквивалентных схемах. Будем использовать в дальнейшем единую систему обозначений для элементов всех подсистем, обычно применяемую при изображении электрических эквивалентных схем. При этом элементы представляют собой двухполюсники, которые могут быть пяти различных видов (рис. 5, а). A1 LA1 CA1 O PA1 LD1 B LD 2 LA2 CB U1 U2 A2 CA2 PB PA2 O а б в Рис. 5. Условные обозначения двухполюсных элементов Получение эквивалентных схем – обычная для инженеровсхемотехников операция, выполняемая при анализе функционирования радиоэлектронных устройств. Переход от принципиальной электрической схемы к эквивалентной заключается в замене обозначений электронных приборов обозначениями двухполюсников (рис. 5, а) и добавлении ветвей, отображающих учитываемые паразитные параметры. Не вызывает затруднений и составление на основе электрогидравлических и электротепловых аналогий эквивалентных схем, отражающих гидравлические, пневматические и тепловые процессы в проектируемых 29 устройствах. Составление эквивалентных схем для механических систем начинается с выбора системы координат, начало О которой должно быть связано с инерциальной системой отсчета. Далее формируются n эквивалентных схем, где п – число степеней свободы. В общем случае возможны три эквивалентные схемы, соответствующие поступательным движениям вдоль координатных осей, и три эквивалентные схемы, соответствующие вращательным движениям вокруг осей, параллельных координатным осям. Рассмотрим правила составления эквивалентных схем на примере одной из эквивалентных схем для поступательного движения: 1) для каждого тела At с учитываемой массой Ci в эквивалентной схеме выделяется узел i и между узлом i и узлом О включается двухполюсник массы Сi; 2) трение между контактируемыми телами Аp и Аq отражается двухполюсником механического сопротивления, включаемым между узлами р и q; 3) пружина, соединяющая тела Ар и Аq, а также другие упругие взаимодействия контактируемых тел Ар и Aq, отражаются двухполюсником гибкости (жесткости), включаемым между узлами p и q. В качестве примера на рис. 5, в приведена эквивалентная схема, моделирующая вертикальные скорости и усилия, возникающие в элементах движущегося транспортного устройства, условно изображенного на рис. 5, б в виде платформы В и колес А1 и А2. Здесь учитываются массы платформы СB и колес СА, жесткости колес LA и рессор LD, a также веса РB, РA1, РA2 платформы и колес. Внешние воздействия отражены источниками скорости U. Часто на эквивалентных схемах рядом с обозначением нелинейного элемента указан его тип или записано его компонентное уравнение. Для отражения взаимосвязей подсистем различной физической природы, из которых состоит моделируемая техническая система, в эквивалентные схемы подсистем вводят специальные преобразовательные элементы. Различают три вида связей подсистем. Трансформаторная и гираторная связи выражают соотношения между фазовыми переменными двух подсистем, этим типам связей соответствуют преобразовательные элементы, представляемые парами источников тока или напряжения. Третий вид связи выражает влияние фазовых переменных одной подсистемы на параметры элементов другой и задается в виде зависимостей C, L или R от фазовых переменных. Варианты эквивалентных схем трансформаторной (а) 30 и гираторной (б, в) связей даны на рис. 6, где запись вида А(В) означает, что фазовая переменная А является функцией фазовой переменной B. I1 I2 I1 U1 U 2  I 2  I1  I1 U 2  U2 I 2 U1  U1 а U1  I 2  U2 б U 2  I1  в Рис. 6. Эквивалентные схемы трансформаторной (а) и гираторных (б, в) связей Так, в гидравлическом приводе связь механической и гидравлической подсистем является гираторной и соответствует рис. 6, б, если для источника объемного расхода в гидравлической подсистеме использовать выражение g = SV, а для источника силы в механической подсистеме – выражение F = SP, где V – скорость перемещения поршня; S – площадь поршня; Р – давление жидкости в цилиндре. Примеры математических моделей элементов электронных схем. Для конденсаторов, катушек индуктивности и резисторов чаще всего применяют простые модели (27). Примерами сложных элементов являются транзисторы, диоды, трансформаторы. На рис. 7, а представлена эквивалентная схема биполярного транзистора, используемая во многих программах анализа электронных схем. Этой схеме соответствуют следующие компонентные уравнения: dUCэ / dt = (Iэ – Iг – Iэ.д) / Cэ; dUCк / dt = (Iк – Iг – Iк.д) / Cк, где Сэ, Ск – емкости и Iг, Iэ.д, Iк.д – токи переходов, являющиеся заданными функциями напряжений UСэ и UСк на емкостях, а для rб и r’к используются уравнения закона Ома. δ I К.Д rK' CK IK rδ СЭ Cc.n Id CЗ.U Э а C C IГ I Э.Д CЗ.C K Подложка Cu.n U б Рис. 7. Эквивалентные схемы биполярного (а) и МДП (б) транзисторов На рис. 7, б представлена эквивалентная схема МДП-транзистора, в которой Id является функцией потенциалов на электродах прибора, а 31 емкости между затвором и стоком Сз.с, затвором и истоком Сз.и, истоком и подложкой Си.п, стоком и подложкой Сс.п считаются либо постоянными, либо зависящими от потенциалов электродов. Примеры математических моделей элементов систем неэлектрической природы. Простыми элементами механических поступательных систем являются элементы массы и гибкости (жесткости). Математическая модель массы выражает закон Ньютона: F = mdU / dt, где F – сила; т – масса; U – скорость. Математическая модель упругого стержня получается из закона Гука:  = El / l, (28) где  и l – напряжение и изменение длины стержня в продольном направлении; Е – модуль упругости; l – длина стержня. Так как  = F / S, где S – площадь поперечного сечения стержня, то после дифференцирования (28) по времени имеем dF/dt  U/L M , где LM = l / (SE) –гибкость стержня (напомним, что обратная величина гибкости – жесткость). Сочетания элементов массы, гибкости и механического трения образуют разнообразные сложные элементы многих механических систем, например манипуляторов. В качестве примера на рис. 8 приведена эквивалентная схема плоского сложного элемента «шарнирная связь двух твердых тел», где С1, С2 – массы, а С3, С4 – моменты инерции соединенных тел. Математическая модель представляет собой систему уравнений, отражающих геометрические соотношения, действующие в системе шарнирно связанных тел: Ux1 – a11 – Ux2 + a22 = 0; Uy1 – b11 – Uy2 + b22 = 0, (29) где Ux и Uy – проекции вектора скорости центра масс на оси координат х и у (индексы 1 и 2 относят соответствующую величину к первому или второму телу);  – угловая скорость относительно центра масс тела; ai и bi – коэффициенты, зависящие от угла поворота i-гo тела в выбранной системе координат. Реакции в шарнире Iх и Iy определяются в результате решения системы уравнений, состоящей из компонентных уравнений (29), компонентных уравнений других элементов и топологических уравнений. 32 Реакции Iх и Iy фигурируют в топологических уравнениях, составляемых для внешних узлов эквивалентной схемы. В случае связи тел с помощью упругой тяги в эквивалентной схеме появляется элемент гибкости, при учете трения – элемент механического сопротивления. Iy Ix C1 a1 I x C3 C2 b2 I y b1 I y C2 C1 a2 I x C4 Рис. 8. Эквивалентная схема «шарнирная связь двух твердых тел» 1.7. ПОЛУЧЕНИЕ МАТЕМАТИЧЕСКИХ МОДЕЛЕЙ СИСТЕМ НА МАКРОУРОВНЕ Математическая модель системы (ММС) на макроуровне представляет собой систему обыкновенных дифференциальных уравнений (ОДУ):  dV  F ,V , t   0 ,  dt  (30) где V – вектор зависимых переменных, называемых также базисными или определяющими координатами; t – время. При получении системы (30) исходными являются компонентные и топологические уравнения. Поскольку выбор, как формы исходных топологических уравнений, так и формы итоговой модели неоднозначен, для получения ММС возможно применение ряда методов. В настоящее время используются три основные формы представления ММС (30) на макроуровне: нормальная форма dV / dt = (V, t), ее применение удобно при ориентации на последующее численное решение явными методами; линеаризованная форма D dV  HV  B dt (31) где D и Н – матрицы; В – вектор. Форма (31) используется в программах, открытых по отношению к методам интегрирования; 33 алгебраизованная и линеаризованная форма ЯnVn = Bn, (32) где Я – матрица; В – вектор. Их вид определяется выбранными методами алгебраизации и линеаризации. Обозначения Яn, Vn, Вn используются для значений переменных Я, V, В, соответствующих моменту времени tn. В подавляющем большинстве современных программ анализа применяют форму (32). Для получения ММС в такой форме применяют методы узловых потенциалов (МУП) и табличные методы. В этих методах для алгебраизации реализуют одну из неявных разностных формул численного интегрирования dV dt   nVn   n , (4.41) (33) n где (dV/dt)n и Vn – значения векторов dV/dt и V в момент времени tn; п – номер шага интегрирования; ηn – коэффициент, зависящий от выбранного метода интегрирования и значения шага hn; n – вектор, зависящий также от значений векторов V на р предыдущих шагах интегрирования, р – порядок метода. Линеаризацию выполняют с помощью разложения нелинейных элементов вектора F(dV/dt, V, t) в ряд Тейлора с сохранением в разложении только линейных членов. Различия между МУП и табличными методами заключаются в выборе исходных топологических уравнений и вектора базисных координат. Для получения ММС в нормальной форме наиболее приемлем метод переменных, характеризующих состояние системы, называемой обычно более коротко – метод переменных состояния (МПС). Метод узловых потенциалов. Исходные топологические уравнения в МУП – уравнения закона токов Кирхгофа (ЗТК) или аналогичные им уравнения, выражающие равновесие переменных типа потока во всех узлах эквивалентной схемы, за исключением лишь одного узла, принимаемого за базовый: АI = 0, (34) где А – матрица инциденций; I – вектор токов. Матрица инциденций характеризует связи узлов и ветвей эквивалентной схемы. В матрице инциденций i-я строка соответствует i-му узлу, а j-й столбец – j-й ветви дерева. Всего в матрице a столбцов и p строк, где а и p – число ветвей и узлов в эквивалентной схеме. Элемент матрицы 34 ai,j = + 1, если i-й узел инцидентен j-й ветви и положительное направление тока в этой ветви выбрано от i-гo узла; аij = –1 при тех же условиях инцидентности, но при противоположном направлении тока, иначе aij = 0. Исходные компонентные уравнения в классическом варианте МУП должны иметь вид I n  G( dU ,U n , t n ), (4.43) dt n (35) где In и Un – векторы токов и напряжений ветвей на n-м шаге. Алгебраизация с помощью (33) приводит к системе уравнений In = G[nUn + n), Un, tn]. После линеаризации относительно вектора Un компонентные уравнения имеют вид In = YnUn + Qn, (36) где Yn = ∂In/∂Un – матрица проводимостей; Qn – вектор, зависящий от значений вектора U, полученных на предыдущих шагах интегрирования. В качестве основных базисных координат в МУП используют узловые потенциалы, вектор которых на n-м шаге обозначим n. Отметим, что связь между векторами Un и n выражается с помощью матрицы инциденций: Un = – Atn (37) Подставив (37) в (36) и (36) в (34), имеем – AYnAtn + AQn, т. е. получена алгебраизованная и линеаризованная форма ММС в виде Яnn = Bn (38) t где Яn = AYnA – матрица Якоби; Bn = AQ, j – вектор правых частей. В классическом варианте МУП имеются ограничения на вид компонентных уравнений. Применительно к схемной форме представления моделей эти ограничения выражаются в недопустимости таких ветвей, как идеальные источники напряжения и любые ветви, параметры которых зависят от каких-либо токов. В модифицированном варианте МУП эти ограничения снимаются благодаря расширению вектора базисных координат – дополнительно к узловым потенциалам к базисным координатам относят также токи особых ветвей. Особыми ветвями при этом называют: 1) ветви источников напряжения; 2) ветви, токи которых являются управляющими (аргументами в выражениях для параметров зависимых ветвей); 3) индуктивные ветви. 35 Обозначим через I1 и I2 векторы токов, а через A1 и А2 – подматрицы инциденций узлов с неособыми и особыми ветвями соответственно. Тогда уравнение (34) перепишем в виде A1I1 + A2I2 = 0. (39) Компонентные уравнения неособых ветвей имеют вид I1 = Y1U + Y2I2 + Q, где в дополнение к уравнению (36) в правой части фигурирует слагаемое, определяемое токами особых ветвей. Здесь Y1 = ∂I1/∂U, Y2 = ∂I1/∂I2; n – индекс, указывающий номер шага, опущен, но подразумевается, что все переменные величины относятся к п-му шагу интегрирования. Подстановка U в (39) дает подсистему из (β – 1) - го уравнения: – A1Y1A1t + (A1Y2 + A2) I2 + A1Q = 0, в которую входят P + Y – 1 неизвестных (Y – число особых ветвей). Эта подсистема доопределяется с помощью у компонентных уравнений особых ветвей. Алгоритм вычисления матрицы Якоби в методе узловых потенциалов. Вычисление матрицы Якоби как матричного произведения AYAt без учета разреженности матриц А и Y нерационально, так как приводит к излишне большим затратам машинных времени и памяти. Например, в схеме средней сложности, включающей β = 51 и  = 80, матрица А имеет размер 50  80, а матрица Y – размер 80  80, т.е. только эти две матрицы для хранения всех их элементов требуют около 10 000 ячеек памяти. В то же время статистические исследования показывают, что ненулевыми в этих матрицах оказываются лишь около 240 элементов. Поэтому на практике используют алгоритмы формирования матрицы Якоби, учитывающие сильную разреженность матриц А и Y. Рассмотрим алгоритм вычисления матрицы Якоби в классическом варианте МУП. Этот же алгоритм используют в модифицированном варианте МУП для вычисления подматрицы A1YA1t, входящей в матрицу Якоби. Исходная информация о матрицах А и Y в этом алгоритме задается в виде списка, где k-я строка списка соответствует k – n ветви эквивалентной схемы. В этой строке указываются номера i и j узлов, инцидентных k-й ветви, местонахождение (адрес) проводимости ∂Ik/∂Uk (где Ik и Uk – ток и напряжение k-й ветви) и, если Ik зависит также от потенциала φp некоторого узла р, отличного от i и j, местонахождение производной Ik по потенциалу 36 этого узла. Каждая k-я строка приводит к добавлению проводимостей в некоторые клетки матрицы Я. Строки списка просматриваются поочередно. Просмотр k-й строки приводит к добавлению проводимости ∂Ik/∂Uk в диагональные клетки яii и яij и к вычитанию этой проводимости из клеток яii и яjj. Если имеется ∂Ik/∂φp  0, то эта производная добавляется в клетки яip и яjp, причем если ток Ik направлен от узла j к узлу i, то добавление в яip происходит со знаком «плюс», а в яjp – со знаком «минус». Пример формирования матрицы Якоби приведен ниже. Табличный метод. В качестве базисных координат используют токи и напряжения всех ветвей схемы, а в качестве исходных топологических уравнений – уравнения Кирхгофа. Эти уравнения записывают для системы контуров и сечений, выбранной в схеме так, чтобы получить  топологических линейно независимых уравнений, где  – число ветвей в схеме. В этих уравнениях фигурируют 2 неизвестных токов и напряжений, поэтому система уравнений доопределяется с помощью  компонентных уравнений. Применяемый способ выбора системы независимых контуров и сечений основан на построении фундаментального дерева в графе схемы. Используется полюсный граф, повторяющий структуру эквивалентной схемы. Фундаментальное дерево связного графа есть связный подграф, включающий β – 1 ребро и не имеющий циклов. Ребра, вошедшие в дерево, образуют множество ветвей дерева (ВД), а остальные ребра – множество ветвей, называемых хордами (ВХ). Контуром k-й хорды называют подмножество ребер графа (ветвей схемы), входящих в замкнутый контур, образуемый при подключении k-й хорды к дереву. Сечения образуются следующим образом: отделим часть вершин графа от остальных с помощью замкнутой линии сечения, проведя ее так, чтобы ни одно ребро не пересекалось более одного раза, и при этом пересекалась одна и только одна ветвь дерева. Следовательно, каждому сечению соответствует определенная ветвь дерева. На рис. 10, (а) для примера приведена некоторая схема, а на рис. 10, (б) – ее граф с выделенным жирными линиями фундаментальным деревом. Штрихом показаны линии сечения. Уравнения токов Кирхгофа для сечений ветвей дерева и напряжений Кирхгофа для контуров хорд образуют систему независимых топологических уравнений 37 I в.д.  М t I x  0,   (4.48) (40) U x  MU в.д  0, где Iв.д., Ix, Uв.д., Ux – векторы токов и напряжений ветвей и хорд соответственно; М – матрица контуров и сечений. 2 C2 R2 1 I вх R1 4 R6 С1 5 I вх 3 2 С1 C4 R4 R2 R3 С3 R5 C5 4 R6 R4 C2 1 3 C4 С3 R3 R5 C5 5 R7 R1 6 R7 6 а б Рис. 9. Эквивалентная схема (а) и ее граф (б) Деление ветвей на хорды и ветви дерева может выполняться по различным правилам. В табличном методе рекомендуется выбирать фундаментальное дерево так, чтобы минимизировалось количество ненулевых элементов в матрице М. Элемент mij матрицы М равен ±1, если в контур i-й хорды входит j-я ветвь дерева, иначе mij = 0. Знаки у mij  0 выбирают по специальному правилу: mij = +1, если положительные направления токов в i-й хорде и j-й ветви дерева при обходе контура совпадают; mij = –1, если эти направления не совпадают. В столбцах источников напряжения записываются +1, если ток подключенной хорды втекает в выбранный за положительный полюс источника. На рис. 9 матрица М имеет вид табл. 2. Таблица 2 С1 С3 С4 С5 К3 С2 +1 –1 R1 –1 R2 –1 +1 R4 +1 –1 R5 –1 +1 R6 –1 +1 +1 R7 –1 Iвх +1 Компонентные уравнения F( V , V, t), где V = (Ux, Iв.д., Uв.д., Ix), V = dV/dt 38 предварительно алгебраизуются относительно вектора Vn: с помощью WVn = Q, (33) и линеаризуются (41) где W = F/V + F/V; Q – вектор правых частей. Окончательно ММС представляет собой систему из 2  уравнений (40) и (41). Сравнительно высокий порядок получающейся системы уравнений – недостаток табличного метода. Метод переменных состояния. Метод ориентирован на получение ММС в виде системы обыкновенных дифференциальных уравнений в нормальной форме с последующим ее интегрированием явными методами. В качестве базисных координат (переменных состояния) выбирают независимые емкостные напряжения и индуктивные токи, а в качестве исходных топологических уравнений – уравнения типа (40). Подмножество ветвей дерева БД образуется путем включения в него ветвей схемы в следующем порядке – ветви источников напряжения, емкостные, резистивные, индуктивные. Именно такой приоритет ветвей обусловливает последующее получение ММС в форме либо нормальной, либо легко преобразуемой в нормальную. После выбора дерева и построения матрицы М систему (40) преобразуют – из нее исключают все токи и напряжения, не относящиеся к переменным состояниям. Исключения производят с помощью компонентных уравнений. Наиболее просто такие исключения выполняют в том случае, если в схеме отсутствуют топологические вырождения (рис. 10), под которыми понимаются емкостные контуры (а) и индуктивные звезды (б). К топологическим вырождениям относят также такие ситуации, при которых в контур резистивной хорды оказывается включенной хотя бы одна резистивная ветвь дерева. Наличие топологических вырождений усложняет процедуру получения ММС в нормальной форме и ее использование – требуется либо решение систем линейных алгебраических уравнений на каждом шаге численного интегрирования системы обыкновенных дифференциальных уравнений, либо предварительное устранение топологических вырождений с помощью изменений схемы. Очевидно, что последний прием может привести к дополнительным погрешностям решения или к ухудшению обусловленности системы уравнений. 39 а б Рис. 10. Топологические вырождения Пример получения ММС различными методами. Для примера рассмотрим схему, показанную на рис. 9. Предполагается, что численное интегрирование уравнений будет выполняться с помощью метода первого порядка точности. Неявная формула такого метода Vn = (Vn - Vn -1)/hn, (42) где hn – величина n-го шага интегрирования В методе узловых потенциалов компонентные уравнения In = yiUn; I = CjdU/dt. Последнее уравнение, алгебраизованное на n-м шаге интегрирования с помощью (42), принимает вид In – yjUn = Qn. (43) Здесь yi = 1/R; yj = Cj/hn; Qn = –CjUn-1/hn. Использование вышеописанного алгоритма формирования матрицы Якоби непосредственно приводит к ММС вида (44). yR1 + yR2 + yC1 – yR2 yR2 + yR3 + yR4 – yR2 – yR4– yC2 + yC2 yR4 + yC2 + – yR4– yC2 – yC4  yC3 + yC4 – yC4 yR5 + yR6 + yC4 – yR6 – yR6 yR6 + yR7 + yC5 1 1 1 1 1 Iвх + yC11 yC2(2 – 3) (44) = yC2(3 – 2) + yC33 + yC4(3 – 4)  yC4(4 – 3) yC55 Здесь уp – проводимость ветви р; φk – потенциал k-го узла на данном шаге интегрирования; φk –то же на предыдущем шаге интегрирования. В табличном методе полученная матрица М дана в табл. 3.  40 С1 –1 –1 –1 С2 +1 –1 –1 С3 +1 –1 –1 –1 С4 +1 +1 Таблица 3 C5 +1 –1 R1 R2 R3 R4 R5 R6 R7 Iвх Итоговая ММС есть следующая система линейных алгебраических уравнений: 1 - Yx 1 1 M - Yвд Mt 1 Ux I вд  U вд Q1 Ix Q2 где 1 – единичная подматрица; Yвд и Yx,–диагональные подматрицы с проводимостями ветвей дерева и хорд соответственно на диагонали; Q1 и Q2 – подвекторы правых частей компонентных уравнений ветвей (43). Здесь первые две строки – топологические уравнения (40), а две последующие – компонентные уравнения. В методе переменных состояния граф и дерево, выбранное в соответствии с приоритетами ветвей, показаны на рис. 11. C2 R2 С1 R6 С3 R5 C 5 R4 R3 I вх C4 R7 R1 Рис. 11. Граф и выбранное дерево в методе переменных состояния Матрица М получается в виде табл. 3. Вектор правых частей (V, I) в системе уравнений V = φ(V, T) с помощью полученной матрицы М находят по следующему алгоритму: 1. Вычисление вектора резистивных напряжений с помощью второго из уравнений (40): Ux = –MUв.д.. 41 Это матричное уравнение удобно развернуть в следующую совокупность уравнений путем сканирования матрицы М по строкам: UR1 = UC1; UR2 = UC1 – UC2 – UC3; UR3 = UC2 + UC3; UR4 = UC2; UR5 = UC3 – UC4; UR6 = UC3 – UC4 – UC5; UR7 = UC5. Отметим, что значения фигурирующих в правой части переменных состояния известны – на первом шаге это начальные условия, на каждом из последующих шагов – это значения, полученные на предыдущем шаге. 2. Вычисление вектора резистивных токов с помощью компонентных уравнений: IRi = URi / Ri 3. Вычисление вектора емкостных токов с помощью первого из топологических уравнений (40): Iв.д. = M tIx Это матричное уравнение при вычислениях разворачивается в системе уравнений путем сканирования матрицы М по столбцам: IC1 = –IR1 – IR2 – Iвх; IC2 = –IR2 – IR3 – I4; IC3 = –IR2 – IR3 – IR5 – IR6; IC4 = IR5 + IR6; IC5 = IR6 – IR7. 4. Вычисление вектора производных переменных состояния с помощью компонентных уравнений: dUCj / dt = ICj / Cj. Далее применяют ту или иную формулу численного интегрирования, преобразующую вектор производных в вектор переменных состояния для очередного момента моделируемого времени, после чего переходят к новому шагу интегрирования. 1.8. МАТЕМАТИЧЕСКИЕ МОДЕЛИ НА МЕТАУРОВНЕ Математические модели аналоговой РЭА. Использование рассмотренных положений схемотехнического моделирования для 42 проектирования сложной аналоговой РЭА на метауровне оказывается затруднительным из-за чрезмерно больших размерностей задач. Необходимы упрощения. Основой снижения размерности задач является макромоделирование. Часто используют ряд дополнительных упрощений и допущений. Главные из них формулируются следующим образом: 1. Однонаправленность в передаче сигналов, т.е. использование макромоделей, в которых отсутствует влияние выходных переменных на состояние входных цепей. 2. Отсутствие влияния нагрузки на параметры и состояние моделируемых систем. 3. Использование вместо фазовых переменных двух типов (напряжение и ток) переменных одного типа, называемых сигналами. При этом компонентные уравнения элемента представляют собой уравнения связи сигналов на входах и выходах этого элемента. 4. Линейность моделей инерционных элементов. Перечисленные допущения характерны для функционального моделирования, широко используемого для анализа систем автоматического управления. Элементы (звенья) систем при функциональном моделировании делят на три группы: 1) линейные безинерционные звенья для отображения таких функций, как повторение, инвертирование, чистое запаздывание, идеальное усиление, суммирование сигналов; 2) нелинейные безинерционные звенья для отображения различных нелинейных преобразований сигналов (ограничение, детектирование, модуляция и т.п.); 3) линейные инерционные звенья для выполнения дифференцирования, интегрирования, фильтрации сигналов. Инерционные элементы представлены отношениями преобразованных по Лапласу или Фурье выходных и входных фазовых переменных. При анализе во временной области применяют преобразование Лапласа, модель инерционного элемента с одним входом и одним выходом есть передаточная функция, а при анализе в частотной области – преобразование Фурье, модель элемента есть выражения амплитудно-частотной и частотно-фазовой характеристик. При наличии нескольких входов и выходов ММ элемента представляется матрицей передаточных функций или частотных характеристик. 43 Рассмотрим примеры моделей элементов аналоговой РЭА. Анализ однокаскадного RC-усилителя, проводимый в курсах основ электроники, позволяет получить следующее приближенное выражение для передаточной функции каскада: h( p )  K0 , (1  p в )(1  1/(p н )) где K0 – коэффициент усиления на средних частотах; τв и τн – постоянные времени каскада на высоких и низких частотах соответственно. Выражение частотной характеристики того же усилителя h( j )  Для однокаскадного функция имеет вид K0 . (1  j в )(1  1/(j н )) резонансного h( p )  S RС-усилителя pL  r , (4.53) LCp 2  Crp  1 передаточная (45) где S – крутизна усиления на резонансной частоте; L, С, r – индуктивность, емкость и последовательное сопротивление потерь параллельного колебательного контура. Примером нелинейного элемента является амплитудный модулятор, для которого при функциональном моделировании используют модель в виде Uвых(t) = Um(1 + mUвх(t))sin(t + ), где Uвых и Uвх – напряжения модулированных и модулирующих колебаний; т – коэффициент модуляции; ω и φ – частота и фаза несущей; Uт – амплитуда несущей на выходе модулятора. Для автогенераторов в качестве модели используется описание формы генерируемых колебаний, так, в случае генератора гармонических колебаний имеем U (t )  U m sin(t   ), Допущения, принимаемые при функциональном моделировании, существенно упрощают алгоритмы получения математических моделей систем (ММС) из математических моделей элементов (ММЭ). Математическая модель системы представляет собой совокупность ММЭ, входящих в систему, при отождествлении переменных, относящихся к соединяемым входам и выходам. Пример 1. Электронный усилитель работает в малосигнальном 44 режиме и состоит из двух каскадов и цепи обратной связи, имеющих передаточные функции K1(р), K2(р) и K3(р) соответственно. Математическая модель может быть получена непосредственно по схеме усилителя, представленной на рис. 12: Uвых(p) = K*(Uвх(p) + K3(p)Uвых(p)), т.е. передаточная функция усилителя h(p) = Uвых(p) / Uвх(p) – K* / [1 – K*K3(p)], (46) * где K = K1(p)K2(p). K3  p  Uвх  p  K1  p  K2  p  Uвых  p  Рис. 12. Схема усилителя Для исследования линейных систем во временной области на основе модели типа (46) можно использовать два подхода. Первый подход связан с применением правил операционного исчисления и требует выполнения прямого преобразования Лапласа над входными сигналами и обратного преобразования над выходными. Второй подход связан с представлением модели (46) в виде системы обыкновенных дифференциальных уравнений с ее последующим численным интегрированием. Пример 2. Модель резонансного усилителя в виде системы обыкновенных дифференциальных уравнений получается из (45) с учетом того, что Uвых(p) = h(p)Uвх(p) и p = d / dt: LC d 2U вых dU вых  dU вх   Cr  U вых  S  L  rU вх  2 dt dt  dt  или, вводя обозначение I = CdUвых / dt, образуем систему уравнений в нормальной форме Коши I = CdUвых / dt; LdI / dt + rI + Uвых = S(LdUвх / dt + rUвх), где Uвх, Uвых – заданная и рассчитанная функция времени соответственно. При наличии в схеме нелинейных звеньев применяют второй из вышеназванных подходов. При проектировании систем, в которых информация представлена в виде огибающей высокочастотных колебаний, возможны два способа введения переменных в модели. При первом способе несущей переменные изображают высокочастотные модулированные колебания. При анализе приходится имитировать поведение объекта в течение большого числа 45 периодов несущей, что зачастую делает неприемлемо крупными затраты машинного времени. При втором способе огибающей переменные отображают огибающие высокочастотных колебаний. Отражение только низкочастотной огибающей существенно ускоряет вычисления, однако построение моделей может оказаться затруднительным. Математические модели логических схем цифровой РЭА На функционально-логическом уровне необходим ряд положений, упрощающих модели устройств и тем самым позволяющих анализировать более сложные объекты по сравнению с объектами, анализируемыми на схемотехническом уровне. Часть используемых положений аналогична положениям, принимаемым для моделирования аналоговой РЭА. Вопервых, это положение о представлении состояний объектов с помощью однотипных фазовых переменных (обычно напряжений), называемых сигналами. Во-вторых, не учитывается влияние нагрузки на функционирование элементов–источников. В-третьих, принимается допущение об однонаправленности, т.е. о возможности передачи сигналов через элемент только в одном направлении – от входов к выходам. Дополнительно к этим положениям при моделировании цифровой РЭА принимается положение о дискретизации переменных, их значения могут принадлежать только заданному конечному множеству – алфавиту, например двоичному алфавиту {0, 1}. Моделирование цифровой РЭА возможно с различной степенью детализации. На логическом (вентильном) подуровне функциональнологического проектирования в качестве элементов аппаратуры рассматривают простые схемы типа вентилей, на регистровом подуровне элементами могут быть как отдельные вентили, так и любые более сложные сочетания простых схем, например регистры, счетчики, дешифраторы, сумматоры, арифметико-логические устройства и т.п. Рассмотрим математические модели элементов на логическом подуровне. Для одновыходных комбинационных элементов математическая модель представляет собой выражение (в общем случае алгоритм), позволяющее по значениям входных переменных (значениям входов) в заданный момент времени t вычислить значение выходной переменной (значение выхода) в момент времени t + t3, где t3 – задержка сигнала в элементе. Такую модель элемента называют асинхронной. При t3 = 0 модель элемента называют синхронной. Модель многовыходного элемента должна 46 включать в себя алгоритм вычисления задержек и значений всех выходных сигналов. Пример 3. Модель элемента 3И – НЕ с фиксированной задержкой t3, входами u1, u2, u3 и выходом y(t  t З )  u1 (t ) & u2 ((t ) & u3 (t ) , где & – символ операции конъюнкции. Для элементов последовательностных схем (элементов с памятью) используют модели, в которых аргументами выходных переменных уj могут быть как входные ui, так и внутренние uk переменные. Вектор внутренних переменных V отражает состояние элемента (состояние его памяти). Пример 4. Модель J – K триггера с входами J и K и состояниями V0 в предыдущем такте и V – в последующем такте: V  V0 & J & K  V0 & J & K  J & K , где «+» – символ операции дизъюнкции. Объединение моделей элементов в общую математическую модель системы выполняется на основе вышеперечисленных допущений отождествлением переменных на соединяемых входах и выходах элементов. Пример 5. На рис. 13 представлен фрагмент схемы цифрового устройства. Соответствующая ему модель есть следующая система логических уравнений: a  u4 & (u1  u2 ); b  u4 & u3 ; y1  a & y2 ; y2  b & y1. U4 U1 U2 U3 & 1 a & y1 & y2 & & b Рис. 13. Фрагмент схемы цифрового устройства Двузначные и многозначные модели. Применение двузначного алфавита приводит к наиболее экономичным алгоритмам моделирования, однако двузначный алфавит ограничивает возможности анализа работоспособности схем. Поэтому чаще используют многозначное моделирование. В трехзначном алфавите сигналы могут принимать значения из множества {0, 1, X}, где 0 и 1 соответствуют низкому и высокому уровням сигнала, а значение X интерпретируется как неопределенное (неизвестное). 47 В моделях элементов с трехзначным (иначе – троичным) алфавитом, называемых трехзначными (троичными) моделями, сигналы могут принимать значение X во время переходных процессов, например во время своего переключения из единичного в нулевое состояние и обратно, а также в установившихся режимах при неопределенных значениях входных и внутренних переменных. Трехзначный алфавит используют для выявления ситуаций, в которых возможны сбои аппаратуры из-за рассогласования времен прохождения сигналов по различным цепям (из-за состязаний сигналов). Дополнительные возможности обнаружения сбойных ситуаций предоставляет использование пятизначного алфавита {0, 1, X, Д, Е}, в котором с помощью значений Д и Е различают переходы сигналов из состояния 1 в состояние 0, и из состояния 0 в состояние 1 соответственно. Иногда алфавит расширяют в еще большей степени, добавляя в него символы, соответствующие таким ситуациям, как импульсные помехи, воздействие некоторых внешних факторов и т.п. Результаты выполнения основных логических операций над аргументами а и b в двух-, трех- и пятизначном алфавитах представлены в табл. 4 (здесь аргументом для операции НЕ служит переменная b). При этом к двух- и трехзначному моделированию относятся только те столбцы, в которых оба аргумента имеют значения соответственно из множества {0, 1} и {0, 1, X}. Таблица 4 a 1 X E Д b 01X E Д 01X E Д 01 X E Д 01 X E Д 01X E Д Результат операции И 0 0 0 0 01X E Д 0 X X X X 0 E X E X 0 Д X X Д Результат операции ИЛИ 01X E Д 1 1 11 1 X1X X X E1X E X Д1X X Д Результат операции НЕ 10 X E Д                     Синхронные и асинхронные модели. Синхронные модели применяют для анализа установившихся состояний логических схем. Они представляют собой системы логических уравнений вида V = F(V, U), (47) где V – вектор, состоящий из выходных уj и внутренних переменных; U – вектор входных переменных. Примером синхронной модели может служить вышеприведенная система уравнений для схемы, показанной на рис. 13. 48 Система уравнений (47) для последовательных схем имеет столько решений, сколько устойчивых состояний при заданном U имеет моделируемая схема. Как правило, для анализа синхронных моделей используют методы, позволяющие получить то решение, которое соответствует исходному значению вектора V и заданному входному набору U. Получение такого решения называют синхронным моделированием. Синхронное моделирование на основе двузначного алфавита позволяет проверить схему на отсутствие грубых ошибок типа неправильных соединений элементов. Дополнительную информацию о наличии в схеме рисков сбоя получают при применении трех- и пятизначного алфавита. Трехзначное синхронное моделирование позволяет обнаружить статические риски сбоя. Статический риск сбоя выражается в появлении ложных сигналов на выходе схемы при неблагоприятном рассогласовании времен переключения входных сигналов. Пример 6. На рис. 14 приведены возможные временные диаграммы для входов u1, u2 и выхода у схемы И. Из рисунка видно появление ложного сигнала на выходе схемы, что и является проявлением статического риска сбоя. Рис. 14. Диаграммы статического риска сбоя Для обнаружения статического риска сбоя требуется на каждом такте синхросигналов двукратное решение уравнений синхронной модели. Первое решение проводится при промежуточных значениях входных переменных: все изменяющиеся из состояний 1 или 0 входные переменные получают значение X, не изменяющиеся сохраняют свои исходные значения. Второе решение проводится при итоговых значениях входных переменных. Если у какой-либо переменной в схеме исходное, промежуточное и итоговое значения имеют последовательности 0 – X – 0 или 1 – X – 1, то данная переменная изображает ложный сигнал, т.е. указывает на наличие статического риска сбоя. Пример 7. Для временной диаграммы (рис. 14) входная переменная u1 49 имеет последовательность исходного, промежуточного и итогового значений 1 – X – 0; переменная u2 – последовательность 0 – Х – 1 (значение X используется для отображения неопределенного состояния переменных при их переключении). В соответствии с правилами выполнения операции И в трехзначном алфавите (табл. 4) имеем для переменной у последовательность 0 – X – 0, что отождествляется со статическим риском сбоя. Пятизначное синхронное моделирование позволяет дополнительно обнаруживать динамические риски сбоя. Динамический риск сбоя выражается в возможности многократного изменения некоторой переменной вместо правильного однократного изменения в течение одного такта синхронизации схемы. Для выявления динамического риска сбоя выполняют двукратное решение системы логических уравнений при промежуточных и итоговых значениях входных переменных. Если у какой-либо изменяющейся переменной последовательность исходного, промежуточного и итогового значений отличается от возможных корректных последовательностей (корректными являются последовательности 0 – Е – 1 и 1 – Д – 0), то в схеме имеет место динамический риск сбоя. Пример 8. На рис. 15, б приведена временная диаграмма, соответствующая такому сочетанию входных сигналов, при котором в схеме (рис. 15, а) проявляется динамический риск сбоя. Этой временной диаграмме соответствуют следующие последовательности значений переменных: для u1 – 1 – Д – 0, для u2 – 0 – Е – 1, для u3 – 1 – Д – 0. Для переменной у по правилам выполнения логических операций в пятизначном алфавите получаем 0 – X – 1, что отождествляется с динамическим риском сбоя. U1 U2 U3 & y & а б Рис. 15. Схема (а) и диаграммы (б) для иллюстрации динамического риска сбоя Следует, однако, заметить, что синхронное моделирование указывает 50 на возможность сбоев, которые в действительности происходят лишь при неблагоприятных рассогласованиях моментов переключения входных сигналов. Так, если (рис. 14) переключение переменной и1 произойдет раньше переключения переменной u2, то в действительности ложного сигнала на выходе не будет, хотя синхронное моделирование указывает на риск сбоя. В этом отношении результаты синхронного моделирования аналогичны результатам анализа на наихудший случай – они могут привести к отрицательным заключениям о работоспособности вполне корректных схем. Асинхронные модели схем составляют из асинхронных моделей элементов и применяют для анализа переходных процессов в цифровой РЭА. Время t в асинхронных моделях дискретизируется и измеряется в количестве тактов. Продолжительность такта достаточно малая – не должна превышать допустимую погрешность расчета временных параметров. Асинхронную модель схемы можно представить в виде V’ = F(V, U), (48) где U – вектор входных переменных; V – вектор текущих значений внутренних и выходных переменных; V’ – вектор будущих значений тех же переменных. Примером логического выражения, входящего в систему (48), является выражение (46). Алгоритм асинхронного моделирования заключается в подстановке в текущий момент времени t известных значений V и U в правую часть выражений (48), в вычислении новых значений переменных и их задержек, в увеличении модельного времени на длительность такта и корректировке вектора V с учетом результатов расчета задвижек и значений вектора V. Далее такты асинхронного моделирования повторяются до исчерпания заданного временного интервала анализа. Асинхронные модели обычно используют с двузначным или трехзначным представлением переменных. Трехзначное асинхронное моделирование позволяет учесть разбросы задержек распространения сигналов в элементах. Пусть в момент времени t1 на вход элемента приходит сигнал, изменяющий состояние элемента с 0 на 1 с задержкой t3, лежащей в интервале [t3min, t3max]. Тогда в асинхронной модели элемента значение выходной переменной 0 при t1  t  t1  t З min ;  y   X при t1  t З min  t  t1  t З max ; 1 при t  t  t 1 З max .  51 С помощью асинхронных моделей можно проанализировать прохождение сигналов во времени в цифровой РЭА, с учетом реальных задержек в элементах, при различных последовательностях входных сигналов. Однако асинхронное моделирование обычно требует заметно больших затрат машинного времени по сравнению с синхронным. Математические модели функциональных схем цифровой РЭА на регистровом подуровне. Первая особенность ММ на регистровом подуровне связана с разнообразием типов функциональных узлов, рассматриваемых в качестве элементарных при моделировании. Разнообразие типов элементов влечет за собой разнообразие их математических моделей. В ММ элементов могут использоваться различные типы данных, в частности величины булевы, целые, вещественные. Эти величины могут быть скалярными и векторными. Введение векторных переменных позволяет лаконично описывать многоразрядные счетчики, регистры, их входные и выходные сигналы. С помощью вещественных величин и операций над ними, которые присущи алгоритмическим языкам общего назначения, можно описать разнообразные алгоритмы, реализуемые в функциональных узлах различной сложности. Математические модели на регистровом подуровне могут быть алгоритмического и схемного типов. Модели алгоритмического типа описывают алгоритмы функционирования устройств без привязки к их схемной реализации. Модели схемного типа отражают связи между переменными на входах и выходах функциональных узлов, составляющих анализируемую схему. Возможны смешанные модели, состоящие из алгоритмических и схемных описаний. В большинстве случаев модели схемного типа представляют в виде конечных автоматов KA = {V0, U, Y, V, , } и называют автоматными моделями. Здесь V0 – исходное состояние автомата; U, Y, V – векторы входных, выходных и внутренних переменных, принимающих значения в конечных множествах;  – функция переходов;  – функция выходов. Функция переходов задает преобразование наборов входных и внутренних переменных на предыдущем такте в набор внутренних переменных на последующем такте: 52 V(t + tЗ) = (V(t), U(t)), где t3 – задержка, а функция выходов задает преобразование векторов V(t) и U(t) в вектор Y(t). Функции  и  могут выражаться с помощью формул, таблиц (матриц), описаний на входных языках прикладных программ. Задержки tЗ могут быть различными для разных путей прохождения сигналов от входов к выходам и в общем случае функционального узла с п входами и т выходами задаются в виде матрицы [ij], i = 1,, n , j = 1,, m , где ij – задержка распространения сигнала от i–гo входа к j–му выходу. В отдельных случаях модель функционального узла может быть представлена в виде алгоритма, в котором действия выполняются над переменными U и Y вещественного типа. В таком виде удобно представлять сложные устройства, например арифметико-логические, выполняющие действия над числами с плавающей запятой. Пример модели функционального узла: Y:= X1*X2, где Y – идентификатор выхода; X1 и Х2 – идентификаторы входов; * – символ выполняемой в функциональном узле операции (алгоритма) над операндами X1 и Х2. Переменные Y, X1 и Х2 могут быть скалярами или векторами булевого или вещественного типа. Подобным образом описываются шифраторы, сумматоры, схемы контроля четности и т.п. 1.9. МАТЕМАТИЧЕСКИЕ КОНСТРУИРОВАНИЯ ЭВА МОДЕЛИ ДЛЯ ЗАДАЧ При проектировании конструкций пользователю удобнее иметь дело с моделями, которые легко образуются, если элементы конструкций принять за точки, а связи между ними – за линии. Такое представление объекта отличается высокой наглядностью, позволяет сосредоточить внимание на наиболее существенных связях, находить оптимальное решение задач проектирования. Использование аппарата теории графов для разработки алгоритмов конструкторского проектирования приводит нас к введению лишь некоторых определений, правил, теорем и положений из общей теории графов, которые будут представлять интерес в дальнейшем изложении. При решении задач автоматизации основные свойства и характеристики объектов описывают с помощью формальных 53 математических объектов, обеспечивающих адекватность и сохраняющих наглядность и необходимую содержательность. При решении задач проектирования возникает необходимость построения различных ММ и выбора из них наиболее приемлемой. Сформулируем основные требования, предъявляемые к ММ: адекватность и простота представления исходного объекта; информационная сложность, т.е. возможность перехода от одной ММ к другой, от объекта к модели и обратно; разумный объем памяти ЭВМ, отводимый для хранения информации о модели; степень разработанности математического аппарата для оперирования с ММ; простота обработки. Очевидно, что самой лучшей моделью является сам объект. Построение его математической модели вызвано попытками формализовать и алгоритмизировать процесс проектирования, а также необходимостью минимизации объема памяти ЭВМ для представления модели. В этой связи рассматриваемые модели не могут полностью отвечать всем требованиям описания объекта, а могут лишь описывать его части. Для этой цели используют теоретико-множественные подходы. Для правильного моделирования объекта необходимо определить все свойства объекта, и каждый раз подбирать абстрактную формальную модель, чтобы между объектом и его моделью можно было устанавливать различные виды изоморфизма. В конструкторском проектировании выделяют ММ схем (структурных, функциональных, электрических), монтажного пространства, самих конструкций. Использование графотеоретических моделей объектов сохраняет всю наглядность и содержательность описываемых объектов и позволяет строить формальные алгоритмы конструирования, которые легко обрабатываются на ЭВМ. Любая функциональная или принципиальная схема объекта состоит из набора элементов, соединенных между собой заданным образом. Тогда схему ЭВА можно рассматривать как некоторое множество элементов x1, x2, …, xn, соединенных между собой цепями из множества Е. Такое представление схемы обычно называют схемой соединений или коммутационной схемой. На рис. 16 показан условный фрагмент схемы. Каждый элемент схемы имеет некоторое множество соединительных выводов, которые будем называть множеством контактов и обозначать 54 С = {с1, c2, ..., сp}. Кроме элементов и контактов в схеме имеются внешние контакты С0, которые осуществляют связь рассматриваемой схемы с другими. e2 C12 e1 C01 C13 C11 e3 C31 e5 x2 C32 e4 e6 2 C41 C33 x3 C03 2 C23 C21 x1 C02 C22 1 C42 C04 x4 Рис. 16. Условный фрагмент схемы Два контакта считаются связанными, если объединяются одной электрической цепью. Под электрической цепью понимается некоторое множество эквипотенциальных контактов. Рассмотрим несколько способов задания схем РЭА и ЭВА графами, гиперграфами и их матричными и списковыми эквивалентами. Зададим схему в виде графа G = (X  E  C, U), где X – вершины графа, соответствующие элементам схемы; Е– вершины, соответствующие цепям схемы; С – вершины, соответствующие контактам элементов. Множество ребер U состоит из элементных F и сигнальных W ребер, причем U = F  W. Ребра-подмножества F определяют принадлежность контактов из С элементам X и задаются парами (xi, ck). Ребра-подмножества W задают вхождение контактов из С в цепи Е и описываются парами (ck, lj). На рис. 17 показан граф схемы рис. 16. C32 C11 l3 C13 l1 C21 C33 x1 C01 l5 x3 C12 C03 C41 l2 C31 C22 x4 l4 x2 C13 C42 C02 l5 C04 x0 Рис. 17. Граф G = (X  E  C, U) фрагмента схемы Обычно граф G задается в виде двух матриц A1  aij1 E C 1, если контакт с i  l j ; и a 1ij  0 в противном случае; 55 A2  aij2 ХС 1, если с j  X i ; и aij2  0 в противном случае. Граф также задают в виде матрицы цепей T = |tij|nk, строки которой соответствуют элементам схемы, а столбцы – контактам. Если элементы имеют различное число выводов, то в качестве k принимается max Ki, i = 0, 1, … Элемент tij – номер цепи, связанный с контактом сj элемента xi. Заметим, что если в матрице Т элементы tij = tpq = … = t, то это означает, что контакты j, q,.... , расположенные на элементах с индексом i, p, ..., , принадлежат одной электрической цепи (одному узлу). Для фрагмента схемы рис. 16 матрица имеет вид x1 x T 2 x3 x4 x0 c1 c2 l1 l2 l3 l2 l4 l3 l5 l6 l1 l4 c3 l3 l5 l5  l5 c4   .   l6 Например, элемент t33 матрицы Т означает, что цепь l5 соединена с контактом с3 в элементе схемы х3. При решении задач компоновки и покрытия на конструкторском этапе проектирования между входами и выходами логических элементов схем устанавливаются различия. Они реализуются путем приписывания ребрам графа схемы направления. Входной сигнал логического элемента исходит из соответствующей вершины, а выходной сигнал направлен к вершине. Каждое ребро имеет вес, равный номеру контакта, что позволяет полностью идентифицировать схему коммутации. Тогда фрагмент схемы рис. 16 можно представить в виде двудольного орграфа (рис. 18, а). x1 C11 C13 C22 C12 x2 C 23 x3 x4 C33 C21 C41 C32 C31 C01 C02 x0 x0 C03 C04 x1 x2 x4 x3 C42 l1 l2 l3 l4 l5 2 l6 1 а б Рис. 18. Двудольный орграф (а), гиперграф (б) фрагмента схемы, изображенной на рис. 16 56 x1 l1 x0 x2 l2 l4 x3 l3 2 x1 x4 l6 1 l5 x2 2 1 1 1 x0 1 x3 1 1 x4 в г Рис. 18. Окончание граф G = (X  V, U) (в) и взвешенный граф (г) фрагмента схемы, изображенной на рис. 16 Так как любой двудольный граф может быть представлен гиперграфом, то иногда коммутационные схемы удобно задавать гиперграфами (рис. 18, б). Основное преимущество такого задания – «плавающая» информация о цепях, каждая из которых может быть представлена любым из покрывающих деревьев. Иногда удобно задать схему в виде графа G = (X  V, U), где V – узловые точки схемы (рис. 18, в). Существует много модификаций описанных методов. Одним из них является построение взвешенного графа схемы. Для этого подсчитывают для каждой пары элементов число соединяющих их цепей. Затем строят граф G = (X, U), в котором вершины соответствуют элементам, а ребра uij с весами – числу цепей между хi и xj. Для фрагмента схемы взвешенный граф показан на рис. 18, г. Существуют и другие модели, производные от описанных, причем основные критерии использования любой модели – это затраты памяти ЭВМ, расчет вычислительных процедур и полиномиальная оценка алгоритмов. 57 ГЛАВА 2. АНАЛИЗ ТЕХНИЧЕСКИХ ОБЪЕКТОВ В САПР 2.1. ТРЕБОВАНИЯ К МЕТОДАМ АНАЛИЗА Многообразие применяемых методов. Анализ технических объектов в САПР основан на математическом моделировании, т.е. на исследовании проектируемых объектов путем оперирования их математическими моделями. В предыдущей главе показано, что функциональными моделями проектируемых объектов на макроуровне являются системы ОДУ, которые могут быть представлены в общем виде, либо предварительно приведены линеаризацией, либо алгебраизацией и линеаризацией к виду системы линейных алгебраических уравнений. К таким же формам уравнений с помощью методов конечных разностей или конечных элементов приводятся ММ объектов на микроуровне. Таким образом, математические модели объектов проектирования на микро- и макроуровнях сводятся к системам обыкновенных дифференциальных и конечных уравнений (под конечными уравнениями понимаются алгебраические и трансцендентные уравнения). Оперирование такими моделями в процедурах одновариантного анализа означает решение соответствующих уравнений. Поэтому методы одновариантного анализа на этих уровнях суть численные методы решения систем дифференциальных и конечных уравнений. То же относится к моделям и методам анализа аналоговой РЭА на метауровне. Многовариантный анализ заключается в многократном повторении решения систем названных уравнений при варьировании внутренними и (или) внешними параметрами. Типовыми процедурами многовариантного анализа, реализуемыми в САПР, являются процедуры анализа чувствительности и статистического анализа. Особенностью цифровой РЭА является использование на функционально-логическом уровне проектирования в качестве моделей систем логических уравнений, а в качестве методов анализа – методов решения этих систем. В вычислительной математике известно большое количество методов численного решения систем уравнений. Однако применение большинства из них в САПР РЭА оказывается неэффективным, что объясняется особенностями ММ проектируемых объектов. Поэтому при создании 58 математического обеспечения САПР усилия направлены не только на разработку математических моделей, но и в не меньшей мере на развитие численных методов и алгоритмов анализа. Поскольку эффективность метода зависит от особенностей решаемой задачи, целесообразна реализация в САПР более чем одного метода для каждого класса решаемых уравнений. Выбор метода в большинстве случаев возлагается на пользователя, что требует от него соответствующей подготовки в области численных методов анализа. При неудачном выборе моделей или методов анализа пользователь САПР может столкнуться с рядом неприятностей: чрезмерной продолжительностью вычислений, несходимостью или неустойчивостью вычислительного процесса, малой точностью получаемых результатов. В САПР целесообразно применять методы, исключающие возможность возникновения подобных ситуаций, т.е. методы, обладающие свойствами высокой экономичности, надежности и точности. Однако эти требования противоречивы и не всегда удается их одновременное удовлетворение в должной мере, поэтому важно уметь распознавать неблагоприятные ситуации и знать факторы, изменение которых может привести к исправлению положения. Экономичность. Экономичность метода характеризуется затратами вычислительных ресурсов (машинного времени ТМ и памяти ПМ) на его применение в некоторых заранее оговоренных условиях (например, в тестовых задачах, в среднем по группе задач определенного класса и т.п.). На показатели ТМ и ПМ обычно оказывают влияние многие факторы и в первую очередь размерность решаемой задачи N. В качестве N принимают порядок решаемой системы уравнений, число элементов, из которых состоит моделируемый объект, и т.п. При сравнении методов по экономичности часто не интересуются абсолютными показателями ТМ и ПМ в конкретной ситуации, а исследуют характер зависимости ТМ и ПМ от N. Наиболее эффективные методы имеют линейную или близкую к линейной зависимость показателей экономичности от сложности задачи. Для многих численных методов характерна полиномиальная зависимость ТМ от N: ТМ  cN, где с – коэффициент пропорциональности. При   2 метод имеет заметное 59 ограничение по сложности решаемых задач. Надежность. Надежность метода оценивается как вероятность получения правильных результатов при использовании метода для решения задач заданного класса. Обычно условия применимости метода связаны с такими характеристиками ММ анализируемых объектов, которые пользователь не может оценить заранее имеющимися в его распоряжении средствами, поэтому возможны ситуации, когда вычислительный процесс оказывается неустойчивым или отсутствует сходимость, что может выражаться в зацикливании или останове ЭВМ из-за переполнения разрядной сетки. В САПР стараются применять надежные методы. Однако высоконадежные методы часто характеризуются недостаточной экономичностью. В этом случае целесообразно комбинирование методов с переходом к трудоемким, но надежным методам только в результате автоматического распознавания ситуаций несходимости или неустойчивости вычислений. Точность. Погрешности решения задачи определяются особенностями используемых моделей, численных методов, ограниченностью разрядной сетки ЭВМ. Каждый источник погрешности должен контролироваться, с тем, чтобы погрешности не превысили предельно допустимые. Обычно точность результатов, получаемых с помощью численного метода, зависит от некоторых параметров, выбираемых «по умолчанию» или задаваемых среди исходных данных. С помощью этих параметров можно управлять погрешностями решения, но необходимо помнить, что снижение погрешностей возможно лишь до некоторого отличного от нуля предела и, как правило, сопровождается увеличением затрат машинного времени. Целесообразно в математическом обеспечении САПР иметь не один, а несколько методов одинакового целевого назначения, но с различными возможностями компромиссного удовлетворения противоречивых требований точности и экономичности. Пользователь САПР должен также знать, что явления зацикливания вычислений или переполнения разрядной сетки могут происходить не только из-за недостатков выбранного численного метода, но и из–за ошибок в задании исходных данных. Некоторые ошибки, связанные с нарушением формальных правил грамматики входного языка, распознаются автоматически. Однако ряд ошибок не может быть выявлен формальными средствами без участия пользователя. Примерами таких ошибок могут быть 60 ошибки в задании численных значений параметров или в задании соединений в анализируемой схеме. Если эти ошибки приводят к получению модели самовозбуждающейся схемы, то возможны явления зацикливания и переполнения разрядной сетки. Направления повышения эффективности методов анализа. Высокие размерности задач проектирования, необходимость выполнения многих вариантов решения систем уравнений при проектировании ЭВМ и других сложных технических объектов обусловливают большие затраты вычислительных ресурсов. Поэтому повышение экономичности методов анализа при соблюдении требований точности является актуальной задачей создания и совершенствования математического обеспечения САПР. Эта задача решается на основе идей и методов, группируемых в несколько направлений. Декомпозиция – деление модели проектируемого объекта на части и раздельный анализ получающихся частей. Если ТМ = cN, то после деления модели на т равных частей затраты машинного времени приближенно оцениваются величиной cm(N/m), т.е. уменьшаются в q ≈ m–1 раз. Однако раздельный анализ происходит в условиях принятия упрощающих предположений о взаимном влиянии частей, т.е. сопровождается увеличением погрешностей расчетов. Декомпозиция составляет основу блочно-иерархического подхода к проектированию. Это направление широко используют как в автоматизированных, так и в неавтоматизированных методах проектирования. Диакоптика – направление исследования сложных систем по частям, отличающееся от декомпозиции тем, что раздельный анализ осуществляется без упрощающих предположений о влиянии частей друг на друга. Экономичность диакоптических методов соизмерима с экономичностью обычных декомпозиционных методов, а точность выше. Учет разреженности матриц – направление экономичной организации операций над разреженными матрицами. Матрицу называют разреженной, если в ней преобладают нулевые элементы. Отказ от хранения нулевых элементов и реализация алгоритмов, в которых игнорируются арифметические действия над нулевыми элементами, могут дать значительную экономию ТМ и ПМ. Учет событийности – направление, называемое также учетом 61 временной разреженности моделей и основанное на исключении из вычислительного процесса действий над неактивными переменными. Неактивной на интервале [t, t + ∆t] переменной называют величину, изменения которой на этом интервале не превышают достаточно малого заранее заданного значения. В моделях сложных систем в каждый момент модельного времени большинство переменных неактивно. Моделирование, основанное на учете событийности, принято называть событийным моделированием. В алгоритмах событийного моделирования необходимо реализовать критерии своевременного включения переменных и соответствующих им частей моделей в группу неактивных (латентных) и своевременного их исключения из этой группы. Комбинирование моделей и методов – одновременное использование при решении конкретной задачи нескольких разнотипных моделей или методов анализа одинакового целевого назначения. Комбинирование может быть пространственным, если разнотипные модели или методы применяют в разных частях общей модели, или временным, если их применяют на разных этапах вычислительного процесса. Пространственное комбинирование является частным случаем диакоптического подхода, так как подразумевает разделение модели на части (фрагменты). Повышение эффективности при комбинировании моделей и методов основано на использовании наиболее подходящих моделей и методов для данного фрагмента и данного этапа вычислений. Пространственное комбинирование моделей, относящихся к разным иерархическим уровням, называют многоуровневым (или смешанным) моделированием. 2.2. МЕТОДЫ АНАЛИЗА СТАТИЧЕСКИХ РЕЖИМОВ Численные методы решения систем конечных уравнений. Большинство методов численного решения конечных уравнений F(V) = 0 (49) относится к итерационным. В этих методах последовательные приближения к корню V* системы (49) производятся по итерационной формуле Vi+1 = Vi + Vi, начиная с выбранного исходного приближения V0, где ∆Vi – поправка на (i + 1)-й итерации. Сходящийся вычислительный процесс заканчивается при выполнении условия ||Vi+1 – Vi|| < ε1 или ||F(Vi)|| < ε2, 62 где ε1 > 0, ε2 > 0. В случае несходимости вычисления следует прекращать, если число итераций И превысит заданное М. Численные методы решения (49) различаются способом вычисления поправки ΔVi. В методе простой итерации Vi = hF(Vi), где h – положительная константа (шаг); F(Vi) – вектор невязок. Условие сходимости метода при решении системы линейных уравнений |1 + hj| < 1 (50) должно выполняться для всех j = 1, 2, ..., п. Здесь λj и п – j–е собственное значение и порядок матрицы коэффициентов системы (49). Если система (49) нелинейна, то условие (50) используется как приближенное, фигурирующие в нем собственные значения относятся к матрице Якоби системы (49): Я = F(V) / V. Метод простой итерации характеризуется медленной сходимостью. В большинстве задач скорость сходимости удается увеличить, применяя релаксационные методы, когда требуется такое упорядочение уравнений, чтобы диагональные элементы матрицы Якоби были отличны от нуля. На очередной (i + 1)-й итерации вектор неизвестных ~ Vi 1  Vi  (1   )Vi , где  – параметр релаксации, выбирается в пределах 0  2 (при  > 1 имеем метод последовательной верхней релаксации (ПВР), при  = 1 – метод Зейделя, при  < l – метод последовательной нижней релаксации (ПНР). Однако не во всех случаях релаксационные методы оказываются эффективнее метода простой итерации); ~i  i  Hf (~1i ,~2i ~ 1,i ,i ni ) – элементы вектора Vi; H = h / akk, fk(V) = 0 – k-е уравнение упорядоченной системы (49); akk – k-й диагональный элемент матрицы Якоби. В методе Ньютона поправка Vi = –Яi–1F(Vi), где Яi–1 – обратная матрица Якоби, вычисленная на i-й итерации. Метод Ньютона, характеризуемый высокой скоростью сходимости, широко распространен в процедурах автоматизированного проектирования. Однако по сравнению с предыдущими методами реализация метода Ньютона связана с увеличенными затратами памяти, требующимися для 63 размещения матрицы Якоби. Кроме того, увеличивается трудоемкость вычислений на одной итерации. В алгоритмах решения системы конечных уравнений по методу Ньютона для вычисления поправки ∆Vi вместо обращения матрицы Якоби используют решение системы линейных алгебраических уравнений ЯiVi = –F(Vi), т.е. метод Ньютона по своей сути есть метод линеаризации решаемой системы уравнений. Важным фактором, управляя которым можно добиться выполнения условий сходимости метода Ньютона, является близость точки начального приближения V0 к точке корня V*. Это обстоятельство привело к появлению метода, повышающего вероятность сходимости метода Ньютона и называемого методом продолжения решения по параметру. В этом методе в решаемой системе уравнений выделяют параметр, влияющий на положение точки корня в пространстве фазовых переменных. Например, при анализе электронной схемы таким параметром может быть напряжение источника питания. Система (49) решается методом Ньютона многократно при ступенчатом изменении параметра. Пусть параметр Е выбран так, что при E  0 имеем V*  0. Тогда при первом решении выбираем V0 = 0 и находим значение корня V*, соответствующее начальному значению параметра E1. Далее увеличиваем Е и решаем систему уравнений при начальном приближении V0 = V0* и т.д. вплоть до решения уравнений при истинном значении параметра Е. Следует отметить, что анализ статических состояний можно рассматривать как частный случай анализа переходных процессов, при котором определяется установившееся состояние объекта. Метод анализа статических состояний с помощью интегрирования системы дифференциальных уравнений, описывающей переходные процессы, называемый методом установления, широко используют в программах анализа проектируемых объектов. Численные методы решения систем линейных алгебраических уравнений (ЛАУ). Процедура решения системы ЛАУ AV = B (51) занимает важное место в математическом обеспечении САПР. Решение системы ЛАУ требуется при анализе статических режимов на каждой 64 итерации метода Ньютона. Эта же процедура многократно выполняется при анализе переходных процессов после приведения математической модели объекта к алгебраизованной и линеаризованной форме (32). При этом эффективность решения названных задач анализа в значительной мере зависит от эффективности используемых методов решения ЛАУ. Для решения систем ЛАУ в большинстве проектных процедур анализа используют метод Гаусса или его разновидности. Вычисления по методу Гаусса состоят из прямого и обратного ходов. При прямом ходе из уравнений последовательно исключают неизвестные, т.е. исходную систему приводят к виду, в котором матрица коэффициентов становится треугольной. Такое приведение основано на k-кратном применении формулы пересчета коэффициентов aij := aij – aik akj / akk (52) ко всем элементам матрицы коэффициентов с i > k и j >k, k = 1, 2, ..., (п – 1), где п – порядок системы уравнений. При пересчете столбец свободных членов рассматривают как (n + 1)-й столбец преобразуемой матрицы. Обратный ход заключается в вычислении неизвестных, начиная с определения из n-го уравнения единственной фигурирующей в нем неизвестной vn. Далее из (п – 1)-го уравнения находят неизвестную vп–1, так как в этом уравнении кроме vn–1 фигурирует только неизвестная vn, значение которой к данному моменту вычислено, и т.д. вплоть до нахождения v1. Методы разреженных матриц. Если выполнять вычисления, пользуясь (52), для всех элементов матрицы коэффициентов, то экономичность метода Гаусса характеризуется кубической зависимостью затрат машинного времени Тm от порядка системы уравнений п. Это приводит к ограничению области целесообразного применения метода Гаусса значениями п в несколько десятков. Однако во многих практических задачах п имеет порядок сотен или тысяч. Применение метода Гаусса к таким задачам оказывается эффективным, если учитывать свойство разреженности матрицы коэффициентов в системе решаемых уравнений (51). Разреженной называют ту матрицу, в которой преобладают элементы, равные нулю. Разреженность S оценивается отношением числа нулевых элементов к общему числу элементов матрицы. Анализ показывает, что в математических моделях большинства проектируемых объектов число ненулевых элементов пропорционально первой степени п. Поэтому если 65 учитывать разреженность матрицы, то ТM можно сделать линейной функцией п и существенно расширить пределы эффективного применения метода Гаусса. Учет разреженности при этом заключается в том, что арифметические действия по (52) не производят, если выполняется хотя бы одно из условий aik = 0 или акj = 0. Следует различать исходную и итоговую разреженность матрицы. Исходная разреженность определяется числом ненулевых элементов, имеющихся в матрице до начала выполнения исключений неизвестных по методу Гаусса. Такие нулевые элементы называют первичными. Однако в процессе вычислений по (52) некоторые коэффициенты aij, бывшие нулевыми, могут стать ненулевыми. Такие коэффициенты называют вторичными ненулевыми элементами. Итоговая разреженность определяется суммарным числом первичных и вторичных ненулевых элементов, и эффективность учета разреженности матрицы тем выше, чем больше итоговая разреженность. Итоговая разреженность SИТ матрицы зависит от ее исходной разреженности и от порядка, в каком расположены уравнения и неизвестные (строки и столбцы матрицы). Существует некоторый оптимальный порядок, при котором SИТ максимальна. Процедура упорядочения строк и столбцов матрицы с целью максимизации SИТ имеет важное значение для эффективности вычислений и включается в алгоритм метода Гаусса. В программах анализа используют приближенные алгоритмы максимизации SИT. Простейший из них – алгоритм, по которому строки матрицы коэффициентов располагаются в порядке увеличения числа ненулевых элементов в них. Эффективность учета разреженности матриц проявляется и в экономии оперативной памяти, так как в алгоритмах, учитывающих разреженность, требуется хранение в памяти ЭВМ только ненулевых элементов. Метод прогонки. Примерами сильно разреженных матриц являются матрицы Якоби в системах конечных уравнений, получаемых по методам конечных разностей или конечных элементов из дифференциальных уравнений в частных производных. Если алгебраизация дифференциального уравнения производится на основе регулярной сетки, то разреженная матрица Якоби оказывается ленточной, т.е. матрицей, у которой ненулевые элементы располагаются только на k главных диагоналях. Специфические особенности структуры ленточных матриц можно использовать для 66 упрощения алгоритмов учета разреженности. Для решения систем ЛАУ с трехдиагональными матрицами коэффициентов используют разновидность метода Гаусса, называемую методом прогонки. Нетрудно заметить, что в трехдиагональных матрицах при исключении очередной неизвестной vi–1 из системы уравнений пересчет по (52) следует производить только в отношении диагонального элемента аii и свободного члена i-го уравнения bi. Обозначим преобразованные по (52) значения аii и bi через ri и qi соответственно. Тогда прямой ход по методу Гаусса сводится к расчету коэффициентов ri и qi, i = 2, 3, ..., п, по рекуррентным формулам, являющимся частными случаями (52): ri = aii – ai,i–1 ai–1,i / ri–1; qi = bi – ai,i–1 qi–1 / ri–1. Начальные значения коэффициентов r1 = a11 и q1 = b1. Обратный ход в методе прогонки также выполняют по очевидной рекуррентной формуле vi = (qi – aijvi+1) / ri, при этом принимают vn+1 = 0. Коэффициенты ri и qi называют прогоночными. В качестве прогоночных коэффициентов часто используют li = –ai–1,i / ri–1 и pi = qi–1 / ri–1, тогда обратный ход осуществляют по формуле vi = pi+1 + li+1vi+1. Алгоритмы метода прогонки в отличие от более общих алгоритмов учета разреженности матриц с нерегулярной структурой характеризуются большей простотой программной реализации. Выбор метода решения системы алгебраических уравнений. Решение систем алгебраических уравнений (АУ) имеет место во многих проектных процедурах и прежде всего в процедурах функционального проектирования. Эффективность решения этих задач вносит существенный вклад в общую эффективность выполнения проектных процедур, поэтому необходимо правильно выбрать метод решения системы АУ. Такой выбор приходится осуществлять разработчику пакета прикладных программ (ППП) для подсистем функционального проектирования. Если же пакет выполнен открытым по отношению к численным методам решения систем АУ и, следовательно, содержит ряд модулей, реализующих альтернативные методы, то выбор метода возлагается на пользователя. На эффективность применения метода оказывают влияние не только особенности самого метода, но и в не меньшей мере особенности решаемой 67 задачи и используемой ЭВМ. Среди наиболее существенных особенностей задач, называемых ниже факторами, отметим размерность п (порядок системы уравнений), число обусловленности Ц и разреженность S матрицы Якоби, а среди особенностей ЭВМ – быстродействие Б, определенное для класса научно-технических задач, емкость оперативной памяти и разрядность машинного слова. Разработчик ППП должен ориентироваться на некоторые диапазоны значений этих факторов, характерные для моделей проектируемых объектов в соответствующей предметной области. Эти диапазоны должны быть либо указаны в техническом задании на разработку ППП, либо спрогнозированы самим разработчиком на основе исследования статистических данных об имеющихся или предполагаемых моделях. Собственно выбор метода основан на установлении взаимосвязи между факторами и характеристиками используемых ЭВМ, с одной стороны, и такими показателями эффективности применения метода, как время решения, вероятность получения правильного результата и его точность, с другой стороны. Рассмотрим, как устанавливаются подобные взаимосвязи для основных методов решения АУ. Для этих методов TM = nИ / Б, где ТM – затраты машинного времени; п – порядок решаемой системы АУ, принимаемый за оценку сложности задачи;  – среднее число арифметических операций, приходящихся на единицу сложности задачи, на одной итерации; И – среднее число итераций. Для решения системы линейных алгебраических уравнений (ЛАУ) вида AV = B выбирают либо метод Гаусса, либо итерационные методы. Для метода Гаусса И = 1, и если не учитывать разреженность матрицы коэффициентов А, то   2(n2 / 3 + 2n). Неучет разреженности ограничивает целесообразность применения метода Гаусса решением задач только невысокой размерности. При n > 50 учет разреженности становится необходимым. Для метода Гаусса при учете разреженности и оптимальном упорядочении строк и столбцов матрицы А в задачах проектирования технических объектов имеем γ = const. Так, для моделей переключательных электронных схем γ ≈ 25, а для распределенных моделей с трехдиагональной матрицей коэффициентов при применении метода прогонки γ ≈ 8. Для решения систем ЛАУ итерационными методами с учетом 68 разреженности матрицы коэффициентов имеем И > 1, а γ = Qn, где q = 1 – S – насыщенность матрицы. Так как Q = K / n, где K – среднее арифметическое для числа ненулевых элементов в одной строке матрицы А, то  = 2K. Так, для моделей переключательных электронных схем получаем по результатам статистических исследований γ ≈ 7,8, т.е. одна итерация выполняется быстрее, чем по методу Гаусса. Однако из–за того, что И >> 1, итерационные методы по показателю ТM практически всегда проигрывают методу Гаусса. Решение систем нелинейных АУ выполняется итерационными методами, при этом на требуемое число итераций И в методе Ньютона решающее влияние оказывает выбор начального приближения, а в остальных итерационных методах – число обусловленности Ц матрицы Якоби решаемой системы уравнений. В методе Ньютона, применяемом в рамках методов установления или продолжения решения по параметру, обычно И не превышает трех. В случаях, если И превышает некоторый порог Ипр (например, Ипр = 7), лучше уменьшать значения коэффициентов, управляющих процессом установления, чем продолжать итерации при И > Ипр. Следует отметить, что при решении нелинейных АУ величина у растет, так как при ее подсчете должны быть учтены затраты на вычисление элементов матрицы Якоби. В методе простых итераций И может достигать неприемлемо больших значений, поэтому целесообразно ввести на И ограничение Игр сверху. Если принять Игр= 1,5*104, то из соотношения Игр = –0,5 Ц lgε при ε = 10–3 получаем, что метод простых итераций можно применять только к решению системы уравнений, у которых матрица Якоби имеет Ц < 104. Методы Зейделя, Якоби, последовательной верхней релаксации (ПВР) имеют аналогичный характер зависимости И от Ц, хотя скорость сходимости у них часто оказывается несколько выше, чем в методе простых итераций. Экономичность метода решения систем АУ определяется также затратами оперативной памяти. При неучете разреженности только на хранение матрицы Якоби нужно n2 ячеек памяти. Поэтому если для одного слова используется 8 байт, то при n = 100 для хранения требуется 80 Кбайт, а при n = 500 – уже 2 Мбайт. Итак, подтверждается вывод о необходимости учета разреженности при решении задач с п > ппр где nпр, зависит от характеристик используемой ЭВМ и, как правило, составляет несколько 69 десятков. В задачах анализа распределенных моделей, в которых п может превышать 104, экономичность метода по затратам машинной памяти становится одной из важнейших характеристик. В таких случаях применяют либо релаксационные методы, либо метод Ньютона с использованием на каждой итерации метода Гаусса, но в рамках рассматриваемого ниже диакоптического подхода. На точность решения задачи оказывают влияние задаваемые пользователем в исходных данных значения допустимых погрешностей ε1 или ε2, а также обусловленность модели. Однако задаваемые значения ε1 или ε2 могут вообще оказаться недостижимыми или из-за несходимости, или изза слишком медленной сходимости вычислительного процесса. Надежность применения метода определяется не только фактом принципиальной сходимости к корню, но и тем, каковы затраты времени Тм на получение решения с требуемой точностью. Ненадежность итерационных методов проявляется либо при неудачном выборе начального приближения к корню (метод Ньютона), либо при плохой обусловленности задачи (методы релаксационные и простых итераций), либо при повышенных требованиях к точности решения (метод простых итераций), либо при высокой размерности задач (метод Гаусса при неучете разреженности). 2.3. МЕТОДЫ АНАЛИЗА ПЕРЕХОДНЫХ ПРОЦЕССОВ Классификация методов численного интегрирования систем обыкновенных дифференциальных уравнений (ОДУ) Методы численного интегрирования ОДУ являются методами преобразования дифференциальных уравнений в алгебраические. После дискретизации независимой переменной t система ОДУ (V*, V, t) = 0 (в каждой точке tk этой переменной представляется в виде системы п алгебраических уравнений (Vk*, Vk, tk) = 0 (53) * * с 2п неизвестными Vk = V (tk) и Vk = V(tk). Система (53) доопределяется уравнениями Vk* = f(Vk) (54) задаваемыми выбранным методом численного интегрирования. Система алгебраических уравнений (53), (54) решается в каждой точке tk, k = 1, 2, .... Ш, где Ш – число точек дискретизации (шагов интегрирования). 70 Формулу численного интегрирования (54), в которой в качестве неизвестных величин фигурируют V*k и Vk, и соответствующие этой формуле методы интегрирования называют неявными. В неявных формулах кроме V*k–1 и Vk могут присутствовать значения переменных V* и (или) V в р предыдущих точках дискретизации tk–i, i = 1, 2, ..., р. При p  2 метод интегрирования называют многошаговым. Следует отметить, что к моменту решения системы (53), (54) значения V*k–1 и Vk–1 для i ≥ l, фигурирующие в (54), уже вычислены на предыдущих шагах. Название метода «многошаговый» происходит из-за использования в формуле интегрирования результатов нескольких предыдущих шагов. Величину р при этом называют порядком многошагового метода. Вместо V*k–1 или Vk–1, i  2, в формуле интегрирования могут присутствовать производные V по t порядка выше первого или заменяющие их результаты некоторых дополнительных вычислений на данном шаге. В этом случае метод называется одношаговым, а порядок одношагового метода совпадает с порядком старшей из используемых производных. Систему алгебраических уравнений, решаемых на каждом шаге численного интегрирования, можно записать также в следующем виде: (V*k–1, Vk–1, tk–1) = 0, (55) * V k–1 = f(Vk), (56) * где V k–1 и Vk – неизвестные величины; Vk–1 вычислены на предыдущем шаге. Формулу численного интегрирования (56) и соответствующие ей методы интегрирования называют явными. Явные методы по аналогии с неявными могут быть одно- и многошаговыми, аналогично определяются порядки явных методов. Очевидно, что необязательно на каждом шаге интегрирования численно решать систему из 2n конечных уравнений. В большинстве случаев выполняют предварительное исключение неизвестного вектора V*k из (53) или V*k–1 из (55) с помощью формул интегрирования (54) или (56) в общем виде и на каждом шаге численно решают систему п уравнений с неизвестным вектором Vk. Методы численного интегрирования ОДУ, применяемые в САПР. В практике машинных вычислений наиболее распространены для решения ОДУ методы Гира, Адамса и Рунге – Кутта. Общий вид формул интегрирования в неявных методах Гира: 71  p  Vk    aiVk i  /h k ,  i 0  где hk = tk – tk–1 – величина k-го шага интегрирования; ai – коэффициенты, значения которых зависят от порядка р метода и величин порядков последних шагов. Формулы Гира называют также формулами дифференцирования назад (ФДН) по той причине, что в них аппроксимация производных в точке tk производится с помощью значений функций, относящихся к данному и предыдущим моментам времени. Формула Гира при р = 1 совпадает с неявной формулой Эйлера: Vk  (Vk  Vk 1 )/h k . Часто применяют формулу Гира второго порядка, называемую также формулой Шихмана, которая при h = const имеет вид 1 3  Vk   Vk  2Vk 1  Vk 2 /h. 2 2  Общий вид формул интегрирования в явных методах Адамса при р ≥ 2: p Vk 1  a1 (Vk  Vk 1 ) /h k   aiVk i ; i 2 в неявных методах Адамса при р ≥ 2: p 1 Vk  a0 (Vk  Vk 1 ) /h k   aiVk i . i 1 Явная формула Адамса при р = 1 называется также явной формулой Эйлера: Vk 1  (Vk  Vk 1 )/h k ; явная формула Адамса при p = 2: 2 Vk 1  Vk  Vk 1 /h k  Vk 2 /3. 3 Неявный метод Адамса второго порядка точности называют также методом трапеций, ему соответствует формула интегрирования Vk  2(Vk  Vk 1 )/h k  Vk 1 (5.9) (57) Рассмотренные методы при p ≥ 2 являются многошаговыми. К одношаговым методам относится метод Рунге – Кутта. В САПР распространены неявные методы трапеций и Гира, а в отдельных случаях применяют явный метод Эйлера. Использование методов возможно, если порождаемый ими вычислительный процесс является устойчивым. Неустойчивость вычислений может возникнуть в связи с катастрофическим ростом 72 погрешностей. Различают локальную погрешность интегрирования, допущенную на данном шаге интегрирования, и погрешность, накопленную к моменту tk за все предыдущие шаги. В неустойчивых методах погрешность решения увеличивается от шага к шагу, что приводит к полному искажению результатов и, возможно, к переполнению разрядной сетки. Среди рассмотренных методов интегрирования имеются А-устойчивые и ограниченно устойчивые методы. А-устойчивым называют метод, при применении которого к интегрированию системы линейных ОДУ V  AV (5.10) (58) с начальными условиями V0 ≠ 0 погрешность решения стремится к нулю при любом значении постоянного шага h > 0 и при tk  . Здесь А – постоянная матрица, собственные значения которой имеют отрицательные действительные части. К А-устойчивым относятся неявные методы Гира и Адамса первого и второго порядков точности. Ограниченно устойчивыми являются остальные из рассмотренных методов, для них характерно сохранение устойчивости вычислений только при выполнении ограничений, накладываемых на значение шага интегрирования. Так, для явного метода Эйлера при h = const в задаче (58) условие устойчивости имеет вид неравенства |1 + hj| < 1, (59) которое должно выполняться для всех собственных значений λj матрицы А. Для матрицы А с отрицательными вещественными собственными значениями K условие (59) можно представить в виде 0 < h < –2 / j. Если при этом система уравнений (58) есть модель динамической системы (например, электронной схемы), то величины – 1 /λj принято называть постоянными времени xt. Тогда условие устойчивости явного метода Эйлера приводится к виду 0 < h < 2min, (60) где τmin – минимальная среди постоянных времени моделируемого объекта. Условия (58) или (59) устойчивости методов интегрирования в применении к нелинейным системам ОДУ можно рассматривать как приближенные, при этом под λi понимают собственные значения матрицы Якоби Я = V* / V. Так как в нелинейных задачах элементы матрицы Якоби непостоянны, то непостоянны и ее собственные значения. Поэтому 73 априорный выбор значения постоянного шага h, удовлетворяющего условиям устойчивости на всем интервале интегрирования [0, Tкон], оказывается практически невозможным (случай гарантированного выполнения условий устойчивости за счет выбора h << τmin неприемлем, так как приводит к большим затратам машинного времени). Интегрирование с постоянным шагом нецелесообразно и в Аустойчивых методах, так как h влияет на точность и время решения. Влияние h на точность решения по-разному проявляется на различных участках моделируемого переходного процесса. Поэтому минимизация затрат машинного времени при соблюдении точностных ограничений возможна только в условиях интегрирования с переменным шагом. Большинство алгоритмов автоматического выбора шага основано на контроле локальных погрешностей интегрирования. Локальные погрешности включают в себя погрешности методические, обусловливаемые приближенностью формул интегрирования, и округления, обусловливаемые представлением чисел с помощью ограниченного количества разрядов. Локальная методическая погрешность многошагового метода порядка р, допущенная на k-м шаге интегрирования, зависит от значения шага hk и оценивается по формуле  k  c V ( p1) ( ) hkp1 , (5.13) (61) где с – коэффициент, зависящий от характера и порядка метода; ||V(p+1)(t)|| – норма (p + 1)-х производных V по t, вычисленных в точке   [tk–p, tk]. В алгоритмах интегрирования задаются значением допустимой локальной погрешности , а значение шага hk ограничивают сверху k < , (62) при этом вектор (р + 1)-х производных приближенно оценивают с помощью вектора (p + 1)-х конечных разностей, вычисленных по уже найденным значениям Vk–i, i = l, 2, ..., (р + 2). Однако выполнение условия (62) с большим запасом приводит к заниженным значениям hk и увеличивает ТM. Поэтому наряду с (62) вводят также ограничение на допущенную погрешность εk снизу и значение шага выбирают автоматически, например, по следующему правилу: m1hk при  k   2 ;  hk 1  h k при  2   k  1 ; m h при    , k 1  2 k 74 где m1 < l; m2 > l; δ2 / δ1 = 2. (Значения m1, m2, δ2 выбирают, исходя из имеющегося опыта вычислений по данному алгоритму.) Для того чтобы сделать несущественным влияние погрешностей округления, предусматривают в соответствующих частях алгоритма вычисления с удвоенным числом разрядов, а также не допускают снижения значения шага ниже некоторого минимального уровня hmin. Сравнение методов и обоснование их выбора для конкретных задач автоматизированного проектирования. Эффективность метода численного интегрирования оценивается его влиянием на экономичность и точность вычислений. Качественно характер влияния значения шага на погрешность интегрирования для различных методов можно представить с помощью графиков, изображенных на рис. 19. Рис. 19. Зависимость погрешности интегрирования от величины и порядка метода В области малых h < hmin значительны погрешности округления. В явных методах наблюдается резкий рост погрешностей при h > hкр, где hкр – максимально допустимая по условиям устойчивости величина шага интегрирования. С ростом порядка метода снижаются погрешности интегрирования, но только при умеренных значениях h < h'. Поскольку величина h' зависит не только от порядков сравниваемых методов, но и от особенностей переходного процесса, что подтверждается формулой (61), то заранее определить, метод какого порядка точности даст оптимальный компромисс между показателями точности и экономичности, затруднительно. В связи с этим в программах анализа электронных схем распространен метод ФДН, основанный на автоматическом выборе не только значения шага, но и порядка р неявной формулы Гира. Если на предшествующем (k–1)-м шаге использовалась формула порядка pk–1, то на следующем шаге выбор производился между формулами порядков pk–1 – 1, или pk–1 + 1. Выбирается та формула, использование которой при заданной погрешности обеспечивает наибольшую величину шага hk. При этом 75 порядок pk, формул Гира не должен выходить за пределы диапазона [1, 6]. Среди неявных методов интегрирования при p = const применяют методы Эйлера, трапеций, Шихмана. Их положительными особенностями являются А-устойчивость и сравнительно малый объем памяти, требующийся для хранения результатов интегрирования, полученных на предыдущих шагах. Однако метод Эйлера не обеспечивает необходимой точности при анализе переходных процессов в слабодемпфированных системах. Метод трапеций в его первоначальном виде имеет недостаток, заключающийся в появлении в численном решении ложной колебательной составляющей уже при сравнительно умеренных значениях шагов, поэтому метод трапеций удобен только при принятии мер; устраняющих ложные колебания. Значительное уменьшение ложных колебаний, но при несколько больших погрешностях, дает формула Шихмана. В целом затраты машинного времени на анализ переходных процессов неявными методами существенно зависят от экономичности алгоритмов численного решения конечных уравнений, применяемых на каждом шаге интегрирования. Обычно для решения конечных уравнений используют метод Ньютона, тогда TM  bNИШ, (63) где b – коэффициент пропорциональности, зависящий в основном от быстродействия используемой ЭВМ; N – показатель сложности анализируемого объекта; И, Ш – число ньютоновских итераций на одном шаге и шагов интегрирования;   [1, 3] и зависит от свойств выбранного алгоритма решения системы линейных алгебраических уравнений (ЛАУ) на каждой ньютоновской итерации. Если разреженность матрицы Якоби не учитывается, то  = 3 и возможности применения неявных методов ограничиваются задачами сравнительно малой размерности. Поэтому в САПР сложных объектов (таких, как БИС) необходим учет разреженности матриц. При этом  в (63) оказывается в интервале [1, 2] и существенно повышает эффективность неявных методов. Величины И и Ш зависят от особенностей ММ схемы и характера анализируемых процессов. По классу задач анализа реакции электронных схем на одиночное импульсное воздействие получены средние типичные значения И = 13 и значения Ш, равные нескольким десяткам – сотням. Явные методы интегрирования целесообразно применять к решению 76 систем ОДУ, представленных в нормальной форме Коши: V* = (V, t) (64) Тогда отпадает необходимость решения систем конечных уравнений на каждом шаге. Например, сопоставляя формулу Эйлера Vk–1* = (Vk – Vk–1) / hk и Vk–1* = (Vk–1, tk–1) получаем явное относительно искомого вектора Vk выражение Vk = Vk–1 + hk(Vk–1, tk–1) (65) (явность этого выражения и послужила причиной названия явных методов интегрирования). Вычисления по (65) фактически сводятся к расчету правых частей системы ОДУ (64) и соответствуют значению а, близкому к единице, в формуле (63). Малый объем вычислений на одном шаге и малый объем требующейся оперативной памяти – положительная черта явных методов. Таким образом, сравнение явных и неявных методов интегрирования ОДУ свидетельствует о большей универсальности последних. Поэтому неявные методы являются основными методами анализа переходных процессов в подсистемах схемотехнического проектирования современных САПР БИС и РЭА. Явные методы могут давать лучшие результаты только в отдельных случаях анализа объектов с хорошо обусловленными ММ и играют в современных САПР вспомогательную роль. 2.4. МЕТОДЫ АНАЛИЗА ПОВЫШЕННОЙ ЭФФЕКТИВНОСТИ Диакоптические методы. Диакоптические методы основаны на фрагментации модели сложного объекта, организации раздельных вычислений по фрагментам с периодическим согласованием результатов, получаемых в отдельных фрагментах. Диакоптические методы применяют для решения систем различных уравнений совместно с традиционными численными методами. Для решения систем линейных алгебраических уравнений (ЛАУ) AV = B применяют диакоптический вариант метода Гаусса, основанный на приведении матрицы коэффициентов к блочно-диагональному виду с окаймлением (БДО). При анализе электронных схем этот вариант называют методом подсхем. В методе подсхем исходную схему разбивают на фрагменты (подсхемы). Фазовые переменные (например, узловые потенциалы) делят на внутренние переменные фрагментов и граничные переменные. Вектор фазовых переменных V = (V1, V2, …, Vi, …, VL), 77 где Vi – подвектор внутренних переменных i-го фрагмента; l – число фрагментов; VL – подвектор граничных переменных (L = l + 1). Матрица коэффициентов имеет структуру БДО A11  AL1 A22  A33 A1L A2 L A3 L 0  All AL 3  ALl AlL ALL  AL 2      Прямой ход метода Гаусса состоит из L этапов. На i-м этапе исключаются переменные Vi, при этом пересчет коэффициентов по формуле Гаусса производится только в подсистеме уравнений Aii ALi AiL ALL Vi B  i , VL BL где Вi и BL – подвекторы правых частей соответствующих уравнений. Следовательно, решение одной сложной задачи заменяется решением L задач меньшей размерности, что при  > 1 в формуле (71) дает снижение ТM. Кроме того, при нехватке оперативной памяти в методе подсхем можно минимизировать число обменов информацией между оперативной и внешней памятью. Для решения систем нелинейных конечных уравнений используют диакоптический вариант метода Ньютона с контролем сходимости итерационного процесса отдельно по выделенным фрагментам. Выполнение условий сходимости в i-м фрагменте является основанием для прекращения вычислений по уравнениям этого фрагмента. Очевидно, что раздельное интегрирование означает и раздельное решение подсистем ЛАУ, относящихся к отдельным фрагментам. Раздельное по фрагментам интегрирование дифференциальных уравнений довольно просто организуется лишь при использовании явных методов. Покажем это на примере решения методом Эйлера системы ОДУ, представленной в нормальной форме Коши и разделенной на две подсистемы: V* = 1(V, Z, t); Z* = 2(V, Z, t), где V и Z – подвекторы фазовых переменных, относящиеся к разным фрагментам (граничные переменные предполагаются включенными в один из подвекторов V или Z). Обозначая величины шагов в подсистемах h и Н и 78 применяя формулу Эйлера для момента времени tk имеем V (t k  h)  V (t k ) Z (t k  H )  Z (t k ) V (t k )  ; Z (t k )  , h H имеем V(t k  h)V(t V+ (t kh) ) =hV(t 1 (V (t k ),),t kZ(t ), ), t ), ) +(t kh), Z(V(t k k 1 k k Z (t k  H )  Z (t k )  H 2 (V (t k ), Z (t k ), t k ), k Z(tk + H) = Z(tk) + H2(V(tk), Z(tk), tk), Пусть h  H. откуда вычисляем V(tk + h) и Z(tk + H). Пусть h < H. Для выполнения очередного шага по переменным V нужно вычислить Z(tk + h), что осуществляют с помощью интерполяционного полинома степени р построенного по p + 1 последним рассчитанным точкам кривой Z(t). Здесь р – порядок метода интегрирования. При применении метода Эйлера p = 1, и значение Z(tk + h), вычисляют на основе линейной интерполяции по известным точкам Z(tk) и Z(tk + h). После осуществления шага по переменным V шаги повторяют – очередной шаг делают по переменным V или Z в зависимости от того, в каком фрагменте моделируемое время имеет меньшее значение. Интегрирование подсистем ОДУ с оптимальным для каждого фрагмента значением шага может привести к существенной экономии затрат машинного времени, особенно при применении неявных методов интегрирования. Однако организация неявного пофрагментного интегрирования оказывается более сложной, чем явного. Примеры методов пофрагментного неявного интегрирования – методы однонаправленных моделей и релаксации формы сигнала (РФС). Методы однонаправленных моделей и релаксации формы сигнала. Модели многих сложных элементов являются однонаправленными. В них могут быть выделены входные и выходные фазовые переменные, причем выходные не влияют на входные. Примерами однонаправленных моделей служат большинство моделей логических элементов. При представлении макромодели логического элемента в схемной форме ее типичная структура имеет вид, показанный на рис. 20. откуда вычисляем V(t k  h)и Z(tk  H ). Рис. 20. Структура модели логического элемента Входной блок I отражает входные характеристики схемы, выходной блок О – выходное сопротивление, блок F – функциональные 79 преобразования сигналов и задержку их распространения. Передача сигналов через линии а–а' и b–b' происходит только в направлении слева направо. Метод однонаправленных моделей применяется при выделении фрагментов в схеме в соответствии со следующими правилами:  границы фрагментов проходят через схемы однонаправленных моделей по линиям а–а' или b–b';  контур любой обратной связи должен полностью находиться в пределах одного фрагмента. После выделения фрагментов необходимо их ранжирование – упорядочение в соответствии с последовательностью прохождения сигналов. Ранг 1 получают фрагменты, на входы которых поступают только внешние возбуждения. Ранг 2 получают фрагменты, на входы которых поступают сигналы с выходов фрагментов ранга 1 и, возможно, также внешние возбуждения, и т.д. Другими словами, ранг г присваивают фрагменту, на входы которого поступают сигналы с выходов только ранжированных фрагментов (внешние возбуждения рассматриваются как сигналы с выходов фрагментов ранга 0), причем старший из рангов предшествующих фрагментов равен r – 1. После фрагментации и ранжирования выполняют раздельное численное интегрирование подсистем дифференциальных уравнений, относящихся к различным фрагментам в порядке увеличения их рангов. Интегрирование выполняют на всем заданном отрезке интегрирования Tкон. При интегрировании уравнений фрагмента с рангом r в качестве входных воздействий используют результаты интегрирования уравнений фрагментов с более низкими рангами. Раздельное интегрирование позволяет организовать вычисления в каждом фрагменте с оптимальным для фрагмента значением шага, что может привести к значительной экономии вычислительных затрат. Однако метод однонаправленных моделей имеет ограниченное применение из-за необходимости соблюдения указанных правил фрагментации. Эти ограничения устраняются в методе РФС. Метод РФС является итерационным методом раздельного интегрирования дифференциальных уравнений. Условие однонаправленности моделей снимается благодаря введению фрагментации 80 схем с перекрытием, поясняемой рис. 21. Заштрихованный участок соответствует подсхеме, включаемой при раздельном интегрировании и в фрагмент А, и в фрагмент В. Чем шире зона перекрытия, тем точнее учитывается нагрузка для фрагмента А и точнее рассчитываются входные сигналы для фрагмента В. Если в схеме нет межфрагментных обратных связей, то достаточно ранжирования фрагментов и выполнения одной итерации пофрагментного интегрирования на интервале [0, Ткон]. При наличии межфрагментных обратных связей требуется выполнение нескольких итераций, с помощью которых форма сигналов на выходах фрагментов приближается к истинной (происходит релаксация формы сигналов). Рис. 21. Фрагментация схемы с перекрытием Очевидно, что выполнение пофрагментного интегрирования на всем заданном интервале [0, Tкон] может привести к излишне большому числу итераций и к колебаниям фазовых переменных в процессе релаксации с излишне большой амплитудой, что потребует значительных затрат машинного времени. Во избежание этих неприятностей можно использовать метод прогнозируемой реакции. В отличие от метода РФС здесь пофрагментное интегрирование производится на отрезках [tk, tk + Hj], где Hj – шаг прогноза, соизмеримый с длительностью задержки прохождения сигнала по j-му контуру обратной связи. Значение Нj определяют автоматически аналогично тому, как рассчитываются значения шагов интегрирования. Комбинированные методы и алгоритмы анализа. При решении задач анализа в САПР получило достаточно широкое распространение временное комбинирование численных методов. Наиболее известны рассмотренные выше алгоритмы ФНД для численного интегрирования ОДУ, являющиеся алгоритмами комбинирования формул Гира. Другим примером временного комбинирования методов служат циклические алгоритмы неявно–явного интегрирования ОДУ. В этих алгоритмах циклически меняется формула интегрирования – следом за шагом неявного интегрирования следует шаг явного интегрирования. В базовом алгоритме неявно–явного интегрирования используют формулы первого порядка 81 точности – формулы Эйлера. Такой комбинированный алгоритм оказывается реализацией А-устойчивого метода второго порядка точности, повышение точности объясняется взаимной компенсацией локальных методических погрешностей, допущенных на последовательных неявном и явном шагах. Следует отметить, что в качестве результатов интегрирования принимаются только результаты неявных шагов, поэтому в алгоритме комбинированного неявно–явного интегрирования устраняются ложные колебания, присущие наиболее известному методу второго порядка точности – методу трапеций. Комбинирование неявных и явных формул интегрирования успешно применяют для повышения эффективности решения нестационарных двумерных задач на микроуровне в рамках метода, называемого методом переменных направлений. Если на очередном временном шаге использовать неявные формулы аппроксимации производных по обеим пространственным координатам x1 и х2, то потребуется решать систему конечных уравнений порядка pq, где р и q – уменьшенное на единицу число участков дискретизации по осям x1 и x2. Однако можно циклически менять характер аппроксимации, например на г–м временном шаге для производных по X, использовать неявную, для производных по х2 – явную формулы, а на (i+1)-м шаге, наоборот, для производных по x1 – явную, для производных по х2 – неявную формулы. Тогда на каждом шаге по времени нужно решать q систем уравнений порядка р или р систем уравнений порядка q. Учет латентности фрагментов. Локальные погрешности интегрирования зависят от значения шага интегрирования h и от характера переходных процессов. Если фазовые переменные претерпевают быстрые изменения, то погрешность не выше заданной обеспечивается при малых h. Если же фазовые переменные меняются медленно, то значения h при тех же погрешностях могут быть существенно больше. В сложных схемах ЭВА, как правило, большинство фрагментов в любой момент времени относится к неактивным (латентным), т.е. к таким, в которых не происходит изменений фазовых переменных, причем отрезки латентности ТЛАТ могут быть довольно продолжительными. В латентных фрагментах допустимо увеличивать шаг интегрирования вплоть до значения Тлат, что эквивалентно исключению уравнений фрагментов из процесса интегрирования на период их латентности. Такое исключение выполняется в алгоритмах учета латентности, относящихся к алгоритмам событийного моделирования. 82 Основу этих алгоритмов составляет проверка условий латентности. Примером таких условий может служить ||∆V || < ε, (66) где ∆V – вектор изменений фазовых переменных; ε – малая положительная константа. При проверке (66) у нелатентных фрагментов в вектор AV входят изменения как внешних, так и внутренних переменных фрагмента, и если условие (66) выполняется на протяжении нескольких шагов интегрирования подряд, то фрагмент включается в группу латентных. При проверке (66) у латентных фрагментов вектор AV состоит из изменений только внешних фазовых переменных, происшедших за время с начала латентности, и фрагмент исключается из числа латентных, если условие (66) нарушается. Следует отметить, что в таких методах, как РФС, латентность учитывается естественным образом, в них не требуется специальных проверок статуса латентности фрагментов. Адаптивное моделирование. Адаптивное моделирование – метод автоматического выбора подходящих моделей для фрагментов в процессе анализа сложной системы. Метод адаптивного моделирования составляют способы решения следующих основных вопросов: 1. Фрагментация. 2. Критерий смены моделей. 3. Определение начальных значений переменных во включаемых моделях. Рассмотрим один из возможных вариантов метода адаптивного моделирования. В нем используются модели трех уровней сложности: M1 – простейшая макромодель фрагмента; М2 – промежуточная по сложности макромодель, в которой блоки I и О (см. рис. 20) адекватно отражают входные и выходные цепи схемы, М3 – полная модель фрагмента. Модели M1 и М2 имеют блочную структуру (см. рис. 21). Фрагментацию используют двоякую: во-первых, функциональную фрагментацию, выделяющую фрагменты, которым соответствуют сменяемые модели; во-вторых, диакоптическую фрагментацию, выделяющую раздельно интегрируемые фрагменты. Диакоптическая фрагментация выполняется с перекрытием, зону перекрытия составляют блоки I макромоделей. При функциональной фрагментации нужно стремиться к получению максимально возможного числа типовых фрагментов, под которыми понимаются фрагменты, для 83 которых заранее разработаны и включены в библиотеку макромодели типов M1 и М2 вместе с рассчитанными для макромоделей типа M1 областями адекватности (ОА). Диакоптический фрагмент есть функциональный фрагмент вместе с зоной перекрытия, которой соответствует блок I из модели схемы-нагрузки. Предполагается, что анализ начинается при представлении всех фрагментов моделями M1. При дальнейшем анализе возможны замены моделей M1 моделями М3. Такие замены выполняются для тех фрагментов, в которых происходит нарушение условий адекватности, обнаруживаемое при выходе вектора внешних переменных фрагмента за пределы ОА. Одновременно с включением полной модели М3 для некоторого фрагмента А происходит замена блоков О в моделях схем, являющихся источниками сигналов для фрагмента А, и замена блоков I в моделях схем, являющихся нагрузками для фрагмента А. Эта замена – включение блоков из моделей M2 вместо блоков из моделей M1. Включение новой модели сопровождается появлением новых фазовых переменных. Очевидно, что для продолжения анализа должны быть определены значения этих переменных в момент tвкл включения модели. Расчет этих значений целесообразно выполнять повторным интегрированием подсистемы дифференциальных уравнений включаемой модели на интервале [0, tвкл]. 2.5. МЕТОДЫ АНАЛИЗА ЛОГИЧЕСКИХ И ФУНКЦИОНАЛЬНЫХ СХЕМ ЭЛЕКТРОННОЙ ВЫЧИСЛИТЕЛЬНОЙ АППАРАТУРЫ Методы решения логических уравнений. Анализ переходных процессов в логических схемах выполняют с помощью асинхронных моделей (47), т.е. на основе асинхронного моделирования. К началу очередного такта U известны значения векторов внутренних Vi = (v1i, v2i, …, vni) и входных Ui - переменных. Подставляя Vi и Ui в правую часть выражения (48), получаем новые значения Vi, которые примут внутренние переменные в моменты времени ti + τk, где τk – внутренняя задержка распространения сигнала vk в соответствующем элементе схемы. Далее переходим к следующему такту, в котором вычисления по (48) повторяются со значениями векторов V и U, соответствующими новому моменту времени ti+1 = ti + 1 (напомним, что время измеряется в количестве тактов). Асинхронное моделирование называют потактовым. 84 Если длительность такта превышает задержку в некоторых элементах, т.е. для некоторых k имеем τk ≈ 0, то в модели (48) последовательностной схемы появляются отдельные неявные относительно vk выражения, а это приводит к необходимости решать подсистемы логических уравнений в пределах каждого такта. Если задержки не учитывать во всех элементах, то имеем синхронную модель (47), с помощью которой анализируются установившиеся состояния в схеме и могут определяться статические и динамические риски сбоя. Синхронная модель – это система логических уравнений V = F(V, U) (67) Для синхронного моделирования (решения систем логических уравнений) используются итерационные методы простой итерации и Зейделя. Алгоритм метода простой итерации при решении (67) совпадает с алгоритмом асинхронного моделирования при τk = 1. На первой итерации (такте) выбирается начальное приближение V0 и подставляется в правую часть (67), при этом определяется новое приближение V1. На второй итерации рассчитывается V2 при подстановке V1 в правую часть (67) и т.д. Условием прекращения итераций является или совпадение результатов двух последних итераций Vi = Vi–1, где i – номер итерации, или выполнение равенства i = imax, где imax – предельно допустимое число итераций. В последнем случае делается заключение о несходимости вычислительного процесса из-за генерирования колебаний в моделируемой схеме. Итерации по методу Зейделя отличаются от простых итераций тем, что новое значение k-го элемента вектора Vi сразу же после его вычисления на i-й итерации заменяет старое значение и используется для вычисления новых значений следующих элементов вектора Vi на той же i-й итерации, т.е. vki = F(v1i, v2i, …, vk–1,i, vk,i–1, ..., vn,i–1). В правильно спроектированной комбинационной схеме каждому входному вектору U соответствует одно определенное решение V* системы логических уравнений (67). Корень V* может быть найден при произвольном начальном приближении как методом простой итерации, так и методом Зейделя. Однако в последовательностной схеме система уравнений (67) имеет несколько решений, соответствующих различным устойчивым состояниям схемы. Обычно требуется найти то состояние, в 85 которое переходит схема из заданного исходного состояния при заданном входном воздействии. Пусть заданному исходному состоянию соответствует значение VИCX вектора V. Тогда для решения сформулированной задачи рекомендуется выбрать начальное приближение V0 = VИСХ и, кроме того, в методе Зейделя необходимо привести в соответствие последовательности вычислений по уравнениям модели и прохождения сигналов в схеме. Последнее условие выполняется с помощью алгоритмов ранжирования, что способствует определению искомого корня за малое число итераций. Например, для комбинационной схемы решение получается за одну итерацию. Ранжирование логических схем и уравнений. Алгоритмы ранжирования комбинационных схем аналогичны алгоритмам ранжирования электронных схем в методе однонаправленных моделей. Ранг 0 присваивается входным цепям, по которым поступают сигналы на схему извне. Ранг r присваивается элементам, все входы которых ранжированы и старший из рангов равен r – 1. Ранг r присваивается также выходам элементов ранга r. После этого логические уравнения упорядочиваются по значениям рангов соответствующих элементов схемы. В последовательностных схемах ранжирование может быть выполнено таким же образом после преобразования схемы в комбинационную с помощью разрывов обратных связей. Ранжирование последовательностей схемы начинается так же, как и комбинационной, но проранжировать таким образом удается только часть элементов. Далее нужно выполнить разрывы контуров обратной связи. Контуры обратной связи ищутся среди элементов, которые проранжировать не удалось. Выбирается любой из этих элементов, например А, к нему добавляется элемент В, выход которого связан с не имеющим ранга входом элемента А. Очевидно, что такой элемент В имеется, иначе элемент А был бы проранжирован. Далее таким же образом присоединяется элемент С к элементу В и т.д., пока в цепочке из элементов А, В, С... не встретится повторно какой–либо из элементов, уже включенных в цепочку. Эта ситуация означает обнаружение контура обратной связи. Далее одна из связей в контуре обратной связи обрывается. Входы элементов от разорванных связей, называемые псевдовходами, получают ранг 0. После этого ранжирование возобновляется. Операции выделения контуров и ранжирования повторяются до тех пор, пока в схеме не останется непроранжированных элементов. 86 Повышение эффективности моделирования логических и функциональных схем. Для повышения эффективности решения уравнений методом Зейделя целесообразно использовать диакоптический подход, в рамках которого итерации выполняются отдельно по фрагментам логической схемы. Введем следующие понятия: составной элемент – множество контуров обратной связи, имеющих попарно общие связи; фрагмент логической схемы – составной элемент или комбинационная схема, состоящая из взаимосвязанных логических элементов, не вошедших в составные элементы. Очевидно, что после обнаружения контуров обратной связи нетрудно найти составные элементы и, следовательно, представить моделируемую схему, состоящую из фрагментов. В такой схеме нет межфрагментных обратных связей и, следовательно, система логических уравнений может быть решена путем однократного обращения к проранжированным подсистемам уравнений. При этом решение подсистем комбинационных фрагментов будет выполнено за одну итерацию Зейделя, а решение подсистем составных элементов – за число итераций, определяемое особенностями данного составного элемента и способом разрыва в нем обратных связей. Это число итераций можно уменьшить, если осуществить разрывы всех контуров обратной связи в минимальном количестве межэлементных связей. Пофрагментное выполнение итераций существенно сокращает общие затраты машинного времени. Другим способом повышения эффективности является параллельное моделирование, основанное на том, что для представления логической переменной достаточно k разрядов, где k = 1 в двузначном алфавите и k = 2 в трехзначном. Тогда моделирование одной и той же схемы можно выполнять одновременно для m = s/k различных наборов входных сигналов, где s – количество разрядов в разрядной сетке ЭВМ. Подобное параллельное моделирование эффективно используется при синтезе тестов для проверки логических схем, где требуется определить реакцию схемы на большое количество входных тестовых наборов. Наиболее общим направлением повышения эффективности математического обеспечения как синхронного, так и асинхронного моделирования является учет событийности. При анализе логических и функциональных схем событием называют изменение состояния любого элемента или, что то же самое, изменение значения любой переменной 87 состояния. В процессе событийного моделирования вычисления производят только по уравнениям активных элементов, т.е. таких элементов, на входах которых на данном такте или итерации произошли события. Рассмотрим алгоритм асинхронного событийного моделирования. В алгоритме используются списки текущих и будущих событий. Все события привязаны к моментам дискретного модельного времени. Ссылки на события, происходящие в текущий момент, находятся в списке текущих событий, а ссылки на те события, наступление которых можно предвидеть, помещаются в список будущих событий. Моделирование текущих событий означает обращение к моделям соответствующих элементов, после чего становятся известными изменяющиеся выходные переменные и задержки в их изменении по отношению к текущему моменту модельного времени. Эти изменения суть события, и сведения о них заносятся в список будущих событий. После того как список текущих событий полностью отработан, модельное время увеличивается, и соответствующая часть списка будущих событий становится списком текущих событий. Так как в сложных схемах на каждом такте активными являются не более нескольких процентов элементов, то учет событийности позволяет существенно уменьшить затраты машинного времени. В маршрутах проектирования БИС и СБИС к числу основных проектных процедур относятся верификация логических и функциональных схем, синтез и анализ тестов. В этих процедурах требуется многократное выполнение моделирования логических схем. Однако высокая размерность задач логического моделирования (СБИС насчитывают десятки–сотни тысяч вентилей) существенно ограничивает возможности многовариантного анализа. Так, современные программы анализа логических схем на универсальных ЭВМ могут обеспечить скорость моделирования приблизительно 103 вентилей в секунду (т.е. на анализ реакции схемы из 103 вентилей на один набор входных воздействий затрачивается 1 с машинного времени), что значительно ниже требуемого уровня. Преодоление затруднений, обусловливаемых чрезмерной трудоемкостью вычислений, происходит в двух направлениях. Первое из них основано на использовании общих положений блочно-иерархического подхода и выражается в переходе к представлениям подуровня регистровых передач. Второе направление основано на применении специализированных вычислительных средств логического моделирования, называемых спецпроцессорами или машинами 88 логического моделирования (МЛМ). Важно отметить, что появление СБИС не только порождает потребности в таких спецпроцессорах, но и обусловливает возможности их создания с приемлемыми затратами. Разработанные к настоящему времени МЛМ функционируют совместно с универсальными ЭВМ и обеспечивают скорость моделирования 106– 109 вентилей в секунду. Логико-электрическое моделирование. Логико-электрическое моделирование представляет собой разновидность многоуровневого моделирования, при котором в анализируемой схеме одновременно используются модели и методы, относящиеся как к схемотехническому, так и к функционально-логическому иерархическому уровню. В анализируемой схеме выделяются подсхемы, подлежащие анализу с помощью логических и электрических моделей. Сопряжение моделей подсхем осуществляется с помощью специальных переходных моделей элементов и алгоритмов синхронизации событий в логической и электрической частях. Переходные модели служат для отображения процессов в элементах с преобразованием аналоговых переменных в логические и наоборот. Логико-электрическое моделирование является незаменимым способом анализа сложных цифроаналоговых схем, поскольку позволяет резко сократить размерность задач, благодаря использованию логических моделей для цифровой части при сохранении необходимой точности анализа, благодаря использованию электрических моделей для аналоговой части схемы. Логико-электрическое моделирование позволяет повысить эффективность решения также ряда задач проектирования цифровой аппаратуры, если в анализируемой логической схеме имеются отдельные фрагменты, требующие для своего адекватного представления моделей схемотехнического уровня. 2.6. МЕТОДЫ МНОГОВАРИАНТНОГО АНАЛИЗА Анализ чувствительности. Анализ чувствительности входит составной частью в алгоритмы решения многих задач, в частности в алгоритмы оптимизации градиентными методами. Для анализа чувствительности задаются ММ объекта и вектор тех внутренних и внешних параметров X, влияние которых на вектор выходных параметров Y требуется определить. 89 В большинстве случаев анализ чувствительности выполняют на основе численного дифференцирования, при котором поочередно задают приращения ∆xi элементам вектора X и определяют получающиеся при этом изменения ∆уj выходных параметров. Тогда абсолютный коэффициент влияния (коэффициент чувствительности) i-го элемента xi вектора X на j-й выходной параметр уj определяется по формуле aji = yj / xi  yj / xi; относительный коэффициент влияния bji = ajixiном / yjном, где xiном и yjном – номинальные значения параметров xi и yj соответственно. Такой метод анализа чувствительности называют методом приращений. Если п – размерность вектора X, то в методе приращений требуется n + 1 раз выполнить одновариантный анализ. Сравнительно большие трудоемкость и погрешности вычислений, присущие численному дифференцированию, относятся к недостаткам этого метода, а универсальность метода – к его преимуществам. Повысить точность вычислений можно, если, определяя aji, производить одновариантный анализ при значениях параметра xi, равных xiном + ∆xi и xiном – ∆xi (остальные внутренние параметры при этом сохраняют номинальные значения). Тогда аij есть отношение разности значений уj, полученных в этих двух вариантах, к 2∆xi. Однако в этом методе увеличивается трудоемкость – нужно 2п раз выполнить одновариантный анализ. Применяют также методы, имеющие менее универсальный характер, например методы анализа чувствительности функционалов зависимостей V(t), получающихся при интегрировании систем ОДУ. В этих методах или сокращается число вариантов интегрирования уравнений, или упрощается система интегрируемых уравнений. Статистический анализ. Статистический анализ имеет целью получение информации о распределении вектора выходных параметров Y при заданном законе распределения случайного вектора X внутренних параметров объекта. Основным методом статического анализа в САПР является метод статистических испытаний (метод Монте–Карло). Каждое k-e статистическое испытание заключается в присвоении элементам Xi вектора 90 X случайных значений Xik и расчете вектора выходных параметров Yk с помощью одновариантного анализа. После выполнения запланированного числа N статистических испытаний их результаты Yk обрабатываются с целью оценки числовых характеристик распределений выходных параметров. В основе алгоритма задания случайных значений параметрам лежит формула  xik  p( x )dx i i  с помощью которой случайное значение ξ величины, равномерно распределенной в интервале [0, 1], преобразуется в случайное значение xik величины xi, имеющей плотность распределения p(xi). Для выработки значений ξ используют стандартные подпрограммы, имеющиеся в программном обеспечении любой универсальной ЭВМ. Выработанное значение ξ интерпретируется как значение F(xik) функции распределения величины xi, так как F ( xik )  xik  p( x )dx . i i  Функция распределения F(xi) монотонная, обычно представлена в памяти ЭВМ в табличной форме, и, следовательно, преобразование ξ в хik, сводится к поиску нужного интервала таблицы и определению результата с помощью интерполяции. Вместо значений ξ равномерно распределенной случайной величины можно использовать значения uk нормированного нормально распределенного вектора U = (u1, u2, …, un) и в алгоритме осуществлять преобразование uik и xik. Но элементы вектора U некоррелированы, поэтому если требуется выработка значений коррелированных случайных величин xi, то вводятся промежуточный вектор коррелированных нормальных величин Z и матрица А преобразования U в Z. Тогда алгоритм задания случайных значений параметрам Xi сначала вырабатывает п случайных значений uik, затем преобразует их в вектор Zk = AUk и далее значения этого вектора в искомые значения xik коррелированных параметров элементов. Матрица преобразования А, так же как и таблица преобразования Z (или U) в X, определяется на основе обработки результатов предварительно выполненных измерений параметров на партии тестовых образцов либо 91 изделий данного или аналогичного типа. Результаты статистических испытаний Yk используются для построения гистограмм, подсчета математических ожиданий и дисперсий выходных параметров. Можно рассчитать также коэффициенты корреляции между выходными уj и внутренними xi параметрами, которые используются для определения коэффициентов регрессии yj на xj. Поскольку относительные коэффициенты регрессии являются аналогами коэффициентов влияния xi на уi, регрессионный анализ, совмещаемый со статистическим анализом, следует рассматривать как возможный подход к анализу чувствительности. Точность и трудоемкость статистических испытаний зависят от их числа N. Для получения большинства интересующих разработчика результатов статистического анализа с приемлемой погрешностью требуется выбирать N = 50–200. Однако получение некоторых результатов (таких, как вероятность выхода годных изделий при значениях этой вероятности, близких к единице или нулю) с приемлемой точностью требует значительно большего числа испытаний. Отсюда следует вывод о значительной трудоемкости статистического анализа. Именно по этой причине статистический анализ проводят лишь на заключительных итерациях процесса проектирования изделий. Разработка тестов цифровой аппаратуры. Проектирование цифровой аппаратуры включает в себя разработку тестов, с помощью которых проверяется отсутствие неисправностей в изготовляемой или эксплуатируемой аппаратуре. Назовем тестируемый объект блоком. Входное воздействие обозначим вектором X = (x1, x2, …, xn), а выходную реакцию – вектором Y = (y1, y2, …, ym), где xi – булева переменная на i-м входе; yj – то же на j-м выходе. Определенному значению Хk вектора X в исправном блоке соответствует значение Yk вектора Y. Пару {Xk, Yk} называют элементарной проверкой. С помощью одной элементарной проверки нельзя выявить все неисправности даже ограниченного класса и тем более нельзя локализовать неисправность. Поэтому тест представляет собой множество элементарных проверок. Число элементарных проверок в тесте называют его длиной. Тест, предназначенный только для установления факта неисправности блока, называют контролирующим (проверяющим), а тест, с помощью которого дополнительно устанавливается место неисправности в блоке, – 92 диагностическим. Тест, который выявляет все неисправности заданного класса, называют полным, а тест, из которого нельзя исключить ни одну элементарную проверку без изменения его полноты, – неизбыточным. Очевидно, что тест, в который включены всевозможные комбинации значений входных сигналов, является полным. Однако громадное число таких комбинаций, равное 2n в комбинационном блоке и 2n+q в последовательностном, делает нереальным применение подобных полных тестов (здесь q – число элементов памяти в блоке). Поэтому возникают задачи синтеза и анализа тестов. При синтезе генерируются входные наборы Xk, а при анализе устанавливаются характеристики теста. Обычно тесты рассчитывают на выявление только одиночных константных неисправностей. Это класс устойчивых неисправностей, возникающих поодиночке в единственном элементе блока. К ним относятся неисправности «константный нуль», «константная единица», связанные с наличием постоянного низкого или высокого уровня напряжения, и «инверсная неисправность», связанная с появлением непредусмотренного инвертирования сигнала. Для синтеза тестов применяют вероятностные и детерминированные методы. В вероятностных методах наборы генерируются с помощью датчиков случайных чисел. Основные затраты машинного времени приходятся при этом на анализ проверяющих возможностей генерируемых наборов. Анализ каждого набора состоит в расчете реакции на воздействие Xk как исправного блока, так и всех его возможных разновидностей. Если блок состоит из N элементов, то имеем 3N таких разновидностей, и общее число вариантов моделирования блока окажется пропорциональным произведению sN, где s – число проверяемых входных наборов. Практика показывает, что при заданной полноте теста s зависит от N и в результате затраты машинного времени ТM оказываются пропорциональными N, где  = 2–3. Большие значения ТM обусловливают применение для анализа тестов наиболее экономичных методов моделирования логических и функциональных схем. Обычно используют параллельное синхронное трехзначное моделирование. Трехзначный алфавит целесообразен для отбраковки входных векторов Хk, приводящих к состязаниям сигналов в блоке, из-за которых результаты применения теста могут стать 93 неопределенными. В детерминированных методах синтеза тестов для каждой неисправности из заданного списка подбирают свой входной набор. Одним из наиболее распространенных алгоритмов, реализующих детерминированный подход к синтезу тестов, является алгоритм Рота. В соответствии с этим алгоритмом для очередной неисправности, связанной с элементом Э тестируемой схемы, подбираются входные для Э воздействия, такие, что выходы исправного и неисправного элементов Э будут иметь неодинаковые значения. После этого ищется путь транспортировки неисправности (т.е. путь передачи информации о наличии неисправности) к выходам схемы, которые доступны для наблюдения. Поиск пути заключается в установлении цепочки элементов между элементом Э и выходом схемы и в определении таких значений входов элементов цепочки, при которых неисправность транспортируется к выходу. Очевидно, что не любая цепочка может служить путем транспортировки, поскольку требуемые для транспортировки значения входов могут оказаться противоречивыми. Поэтому в алгоритме Рота осуществляется последовательное наращивание цепочки с проверкой непротиворечивости условий транспортировки для ранее и вновь включаемых элементов. После установления пути требуемые значения входов у элементов цепочки пересчитываются в значения тех входов тестируемой схемы, которые доступны для управления. Эти значения входов и представляют собой очередной искомый набор теста. В программах синтеза тестов обычно используют комбинирование вероятностных и детерминированных методов. С помощью первых генерируется и анализируется некоторое заранее заданное число наборов. Если полнота теста после этого недостаточна, то выбирают неисправности, не выявляемые уже полученными наборами, и для этих неисправностей ищут входные наборы детерминированным методом. 94 БИБЛИОГРАФИЧЕСКИЙ СПИСОК 1. Норенков И.П. Основы автоматизированного проектирования: Учебник для вузов. – М.: Изд-во МГТУ им. Н.Э. Баумана, 2010. 2. Корячко В.П., Курейчик В.M., Норенков И.П. Теоретические основы САПР: Учебник для вузов. – M.: Энергоатомиздат, 1987. 3. Казеннов Г.Г. Основы проектирования интегральных схем и систем. – М.: БИНОМ. Лаборатория знаний, 2005. 4. Автоматизация проектирования радиоэлектронных средств: Учеб. пособие для вузов / Под ред. О.В. Алексеева. – М.: Высшая школа, 2000. 5. Аветисян Д.А. Автоматизация проектирования электротехнических систем и устройств. – М.: Высшая школа, 2005. 6. Конструкторско-технологическое проектирование электронной аппаратуры / Под ред. В.А. Шахнова. – М.: Изд-во МГТУ им. Баумана, 2005. 7. Евгенев Г.Б. Интеллектуальные системы проектирования. – М.: Издво МГТУ им. Баумана, 2009. 8. Кунву Ли. Основы САПР (CAD/CAM/CAE). – СПб.: Питер, 2004. 9. Норенков И.П. Введение в автоматизированное проектирование технических устройств и систем: Учебное пособие для студ. втузов. – М.: Высшая школа, 1986. 95 СОДЕРЖАНИЕ ВВЕДЕНИЕ ............................................................................................................ 3 ГЛАВА 1. МАТЕМАТИЧЕСКИЕ МОДЕЛИ ПРОЕКТИРУЕМЫХ ОБЪЕКТОВ ........................................................................................................... 4 1.1. ОБЩИЕ СВЕДЕНИЯ О ПРОЕКТИРОВАНИИ ........................................ 4 1.2. КЛАССИФИКАЦИЯ МАТЕМАТИЧЕСКИХ МОДЕЛЕЙ ...................... 5 1.3. ТРЕБОВАНИЯ К МАТЕМАТИЧЕСКИМ МОДЕЛЯМ.......................... 10 1.4. МЕТОДЫ ПОЛУЧЕНИЯ МОДЕЛЕЙ ЭЛЕМЕНТОВ ............................ 13 1.5. МАТЕМАТИЧЕСКИЕ МОДЕЛИ ОБЪЕКТОВ ПРОЕКТИРОВАНИЯ, ИСПОЛЬЗУЕМЫЕ НА МИКРОУРОВНЕ ...................................................... 17 1.6. МАТЕМАТИЧЕСКИЕ МОДЕЛИ ОБЪЕКТОВ ПРОЕКТИРОВАНИЯ, ИСПОЛЬЗУЕМЫЕ НА МАКРОУРОВНЕ ...................................................... 25 1.7. ПОЛУЧЕНИЕ МАТЕМАТИЧЕСКИХ МОДЕЛЕЙ СИСТЕМ НА МАКРОУРОВНЕ ............................................................................................... 32 1.8. МАТЕМАТИЧЕСКИЕ МОДЕЛИ НА МЕТАУРОВНЕ .......................... 41 ГЛАВА 2. АНАЛИЗ ТЕХНИЧЕСКИХ ОБЪЕКТОВ В САПР ................... 57 2.1. ТРЕБОВАНИЯ К МЕТОДАМ АНАЛИЗА .............................................. 57 2.2. МЕТОДЫ АНАЛИЗА СТАТИЧЕСКИХ РЕЖИМОВ............................. 61 2.3. МЕТОДЫ АНАЛИЗА ПЕРЕХОДНЫХ ПРОЦЕССОВ .......................... 69 2.4. МЕТОДЫ АНАЛИЗА ПОВЫШЕННОЙ ЭФФЕКТИВНОСТИ ............ 76 2.5. МЕТОДЫ АНАЛИЗА ЛОГИЧЕСКИХ И ФУНКЦИОНАЛЬНЫХ СХЕМ ЭЛЕКТРОННОЙ ВЫЧИСЛИТЕЛЬНОЙ АППАРАТУРЫ............... 83 2.6. МЕТОДЫ МНОГОВАРИАНТНОГО АНАЛИЗА ................................... 88 БИБЛИОГРАФИЧЕСКИЙ СПИСОК ........................................................... 94 Учебное издание Гладков Леонид Анатольевич Гладкова Надежда Викторовна Модели и методы анализа проектных решений Конспект лекций Ответственный за выпуск Редактор Корректор Заказ № Формат 60 х 84 1/16. Гладков Л.А. Надточий З.И. Надточий З.И. Подписано в печать Тираж 300 экз. Усл. п.л. – 6,1. Уч. - изд.л. – 6,0. Издательство Южного федерального университета 344091, г. Ростов-на- Дону, пр. Стачки, 200/1. тел. (863) 2478051 Отпечатано в Секторе обеспечения полиграфической продукцией кампуса в г. Таганроге отдела полиграфической, корпоративной и сувенирной продукции ИПК КИБИ МЕДИА ЦЕНТРА ЮФУ ГСП 17А, Таганрог, 28, Энгельса, 1. Тел. (8634) 371717, 371655