БИБЛИОТЕКА ПО АВТОМАТИКЕ
Выпуск 95
А. В. ШИЛ ЕЙ КО
ЦИФРОВЫЕ
МОДЕЛИ
ИЗДАТЕЛЬСТВО «ЭНЕРГИЯ»
МОСКВА 1964 ЛЕНИНГРАД

Р е д а к ц и о и и а я к о л л е г и я:
И. В. Антик, А. И. Бертинов, С. Н. Вешеневский, Л. М. Закс,
Н. Е. Кобринский, В. С. Кулебакин, В. Э. Низе, В. С. Малов,
А. Д. Смирнов, Б. С. Сотсков, А. С. Шаталов.
ЭЭ-5(4)-3
УДК 681.142.324
Ш 55
В книге рассматривается ряд вопросов теории и
некоторые принципы построения цифровых моделей.
Описываются существующие конструкции цифровых
моделей, а именно цифровых дифференциальных ана-
лизаторов и сходных с ними машин. Рассматриваются
различные численн?.1е методы решения дифференциаль-
ных уравнений и различные варианты цифровых струк-
тур, реализующих эти численные методы. Выводится
ряд оценок «качества» алгоритмов и структур цифро-
вых моделей. В заключение предлагается методика
'Синтеза оптимальных структур цифровых моделей.
Книга предназначается для широкого круга инженеров
и техников, интересующихся вопросами вычислитель-
ной техники и вопросами применения вычислительной
техники для задач управления. Она может быть полез-
ной также студентам старших курсов вузов.
Шилеико Алексей Вольдемарович.
Цифровые модели, М.—Л., Издательство „Энергия", 1964, 112 с.
с черт. („Библиотека по автоматике", вып. 95).
Редактор И. А. Б ран Техн. редактор Н. И. Борунов
Сдано в набор П/Х 1963 г. Подписано к печати 14/1 1964 г.
Т-12291 Бумага 84Х1087з2. 5,74 п. л. Уч.-изд. л. 7,6
Тираж 18 000 экз. Цена 38 коп. Зак. 561
Московская типография № 10 Главполиграфпрома
Государственного комитета Совета Министров СССР по печати
Шлюзовая наб., 10.

ВВЕДЕНИЕ
Цифровое кмодслироваи'ие — это одна из наиболее молодых обла-
стей 'вычислительпой техники, |начало которой было положено
в 1950 г., когда была разработана первая 'машина, названная циф-
ровым дифференциальным анализатором. В СССР первые работы
в этой области проводились Ф. В. Майоровым [Л. 29] и Е. А. Дроз-
довым [Л. 16]. Разработки цифровых дифференциальных анализа-
торов проводились также К. С. Неслуховским [Л. 40], Л. М. Голь-
денбергом, А. А. Гураковым и А. Г. Шевелевым [Л. 15] и А. В. Ка-
ляевым [Л. 19]. Большой цикл работ по созданию и использованию
цифровых моделей в системах управления станками был выполнен
в V нституте электромеханики АН СССР под руководством А. А. Во-
ронова. Ему же принадлежит одна из первых теоретических работ
в рассматриваемой области [Л. 9]. Разработкой и использованием
цифровых моделей занимались также В. А. Брик [Л. 7] и В. В. Ка-
рибский [Л. 21].
История развития рассматриваемой области вычислительной тех-
ники достаточно показательна. Первой машиной, которую можно
отнести к классу цифровых моделей, был цифровой дифференциаль-
ный анализатор (ЦДА) MADDIDA, разработанный фирмой Northrop
Aviation (США). Машина была предназначена для решения обык-
новенных линейных и ислииейных дифференциальных уравнений.
В ней использовалась простейшая формула численного интегрирова-
ния, так называемая формула прямоугольников, в соответствии с ко-
торой процесс интегрирования иа каждо;м шаге сводился по суще-
ству к сложению двух чисел. Машина состояла из ряда интеграто-
ров. Информация, передаваемая между интеграторами, кодирова-
лась ио методу дельта-модуляции, что предельно упрощало способ
связи между отдельными элементами. Все перечислсеные особежто-
сти привели к тому, что машина получилась исключительно про-
стой, а способ набора задач па ней в основном соответствовал спо-
собу набора задач иа аналоговых вычислительных машинах. Созда-
лось впечатление, что появился совершенно новый тип машин, со-
вмещающий в себе все достоинства аналоговых вычислительных
машин в части простоты программирования с потенциальной точ-
ностью цифровых вычислительных машин. Ыа основании ошибочного
отожгтествления большого количества отдельных ншгоп интсгрирочза-
ния, выполняемых в машине в одну секунду, с быстродействием
в литературе неоднократно высказывалось мнение о высоком бы-
стродействии ЦДА [Л. 1, 30]. В создании такого впечатления ие-
3

малую роль сыграли, несомненно, рекламные сообщения фирм, вы-
пускающих ЦДА.
Вслед за разработкой ЦДА MADDIDA появился целый ряд
аналогичных разработок ,в США (фирмы Computer Research Corp.,
Bendix, Litton), Австралии (Австралийский инженерный институт)
и Японии (фирма Токио Шибаура Электрик). Однако очень скоро,
как 1в результае накопившегося опыта работы, так и 'В результате
более глубокого исследования, мнение о ЦДА последовательного ти-
па резко изменилось. Прежде всего выяснилось, что использование
формулы прямоугольников, а также метода дельта-модуляции для
связи .между блоками резко ограничивает как точность получаемых
решений, так и быстродействие. Так, например, решение простейшей
задачи внешней баллистики, сводящейся к системе двух дифферен-
циальных урашений .второго порядка, при точности в три верных
десятичных знака занимало 30 мин чистого машинного времени
[Л. 26]. Ограничения, свойственные ЦДА, в большом ряде случаев
приводили к тому, что решения отдельных уравнений вообще не
удавалось получить с приемлемой точностью. Первоначальное ^пред-
ставление о простоте программирования также в большой степени
оказалось иллюзорным, так как здесь в полной мере сохранялись все
трудносгги, связанные с подготовкой уравнений и отладкой состав-
ленных програм.м. Известны случаи, когда попытки устранить ука-
занные недостатки приводили к тому, что в процессе разработки,
первоначально имевшей целью построение ЦДА, ряд вводимых
усложнений приводил в конечном итоге ж построению универсаль-
ной цифровой вычислительной 'Машины. Наконец, в результате раз-
вития методов автоматического программирования вопрос о просто-
те программирования окончательно потерял значение.
Некоторое оживление внесли сообщения американских фирм
Paocard Bell Co. и Bendix и английской фирмы А. V. Roe о разра-
ботке ЦДА параллельного типа, в которых относительно высокое
быстродействие достигалось ценой построения структур с чрезвычай-
но большим количеством элементов. Однако и это полностью не
решило проблемы, поскольку оставались в силе ограничения, касаю-
щиеся точности.
iB настоящее время вопрос о ЦДА можно считать окончательно
решенным. Сам принцип, заложенный в этих машинах, является
плодотворным и может быть с успехом использован при разработке
специализированных вычислительных устройств, предназначенных
для решения одного какого-либо вида уравнений, и особенно в слу-
чаях, когда не требуется большого быстродействия. Однако всякая
попытка построения универсальной машины этого типа, т. е. маши-
ны, предназначенной для решения дифференциальных уравнений
вообще, IB лучшем случае должна приводить к получению конструк-
ции того же порядка сложности, что и универсальная цифровая вы-
числительная машина, обладающей такими же внешними парамет-
рами, но лишенной при этом преимущества универсальных вычис-
лительных машин, а именно гибкости. Это вполне объяснимо, так
как преимущества цифровых моделей —это преимущества, давае-
ем ы е СП ед и а л из а ци ей.
У читателя может возникнуть законный .вопрос: почему в свете
подобных высказываний ЦДА уделяется так много внимания в дан-
ной книге? Однако это вполне естественно. Прежде всего разработ-
ки и в отдельных случаях серийный выпуск ЦДА продолжаются
4

в настоящее время как в нашей стране, так и за рубежом. В этой
связи ощущается настоятельная необходимость в работе, в которой
бы достаточно подробно и по возможности строго был исследован
вопрос о вычислительных возможностях ЦДА и была бы очерчена
область целесообразных их применений. Именно такую цель и ста-
вил перед собой автор при написании главы, посвященной ЦДА.
Кроме того, при всех своих недостатках ЦДА сыграли значи-
тельную роль, так как именно в них впервые были реализованы
принципы, составляющие в настоящее время одно из двух, в боль-
шой степени иротивоиоложных друг другу направлений развития
цифровой вычислительной техники. Одно из этих направлений сво-
дится 1к непрерывным усовершенствованиям способов автоматическо-
го программирования. Лри этом (возникают и развиваются специаль-
ные алгоритмические языки (Алгол, Фортран и т. in.). В предельном
случае программа, записанная иа одном из алгоритмических языков,
оказывается ^полностью абстрагированной от реальной вычислитель-
ной машины, иа которой она может быть реализована.
Второе направление в иротивоиоложность первому предусматри-
вает построение цифровой вычислительной системы, структура кото-
рой в наибольшей возможной степени оказывается опециализиро-
ванной на реализации одной данной программы или даже отдель-
ных входящих в состав этой программы операторов. Такое направ-
ление является направлением структурной специализации. Если
в первом случае достигается простота подготовки задачи и возмож-
ность выполнять вычисления по однажды составленной программе
практически иа любой универсальной цифровой вычислительной
машине, то во втором случае достигаются наилучшие динамические
характеристики получаемых вычислительных систем. Цифровое мо-
делирование и представляет собой одну из конкретных задач, ре-
шаемых в рамках этого второго направления.
Первые разработки цифровых моделей, не относящихся к типу
ЦДА, были описаны в работах Мейера [Л. 33], Харриса [Л. 18],
Муррея [Л. 36] и Е. А. Дроздова [Л. 16]. В настоящее время коли-
чество работ по разработке и применению цифровых моделей не-
прерывно растет, что является лучшим свидетельством перспектив-
ности этого направления. Правильное использование принципов циф-
рового моделирования позволяет получать существенную экономию
в количестве оборудования и выигрыш в быстродействии н точности.
При создании цифровых моделей возникают две самостоятельные
задачи: задача построения алгоритма и задача построения струк-
туры.
Одна из первых работ по синтезу алгоритма решения обыкно-
венных дифференциальных уравнений была выполнена Р. Боксером
и С. Талером [Л. 6]. Дальнейшее развитие этот вопрос получил
в работах Ю. А. трейдера, Н. Я. Матюхина [Л. 31] и автора [Л. 50].
Задача построения структуры требует для своего решения наличии
численных оценок «качества» таких структур. При цифровом моде-
лировании в первую очередь представляют интерес динамические
характеристики. Понятие частотной характеристики цифровой си-
стемы впервые было сформулировано в работах Дж. Солзера [Л. 43]
и Я. 3. Цыпкина [Л. 46].
'Предлагаемая читателю книга представляет собой одну из пер-
вых попыток систематического изложения ряда основных проблем
цифрового моделирования. Книга разбивается на две самостоятель-
5

ные части. В первой из них после введения ряда определений и фор-
мулировки основных задач проводится критический обзор цифро-
вых дифференциальных анализаторов и ряда других машин, близ-
ких к ЦДА ПО структуре и принципу действия. Во второй части рас-
сматриваются задачи синтеза алгоритмов и структур цифровых мо-
делей. При написании книги автор стремился ограничить себя кру-
гом вопросов, имеющих непосредственное отношение к проблемати-
ке цифрового моделирования, избегая при этом как пересказа ма-
териалов, имеющихся в других источниках, так и чрезмерного
усложнения, необходимого подчас для придания 'изложению соот-
ветствующей строгости. В частнооти, автор предполагает хорошее
знакомство читателя с основами цифровой вычислительной техники.
Поэтому такие, например, понятия, как регистр, вентиль, запоми-
нающее устройство, сумматор, двоичная система счисления и т. п.,
вводятся в тексте без каких-либо ком'ментариев. Читателей, не
имеющих соответствующей подготовки, мы отсылаем к работе
А. И. Китова и Н. А. Криницкого [Л. 22] или любому другому си-
стем атичеокому курсу основ цифровой вычислительной техники. При
рассмотрении некоторых вопросов, как, например, вопроса о накоп-
лении ошибок при численном интегрировании обыкновенных диффе-
ренциальных уравнений, было сочтено целесообразным ограничиться
рядом иллюстративных примеров, дающих общее представление
о сущности явлений. Строгое и последовательное изложение этого
материала читатель может найти в работах М. Урабе [Л. 38],
М. Р. Шура-Бура [Л. 51] и Ю. В. Ракитского [Л. 42].
Автор считает приятным долгом выразить глубокую благодар-
ность доктору техн. наук Б. Я. Когану, в беседах с которым роди-
лась основная идея этой книги. Главы 4 и б написаны совместно
с Г. М. Козыревой. Автор сердечно благодарит доктора техн. наук
Я. 3. Цыпкина и канд. физ.-математич. наук И. А. Брина за внима-
ние, проявленное к работе, и ряд ценных указаний.
Л. В. ШИЛЕИКО

ГЛАВА ПЕРВАЯ
ЦИФРОВОЕ МОДЕЛИРОВАНИЕ
1. ИСХОДНЫЕ СООБРАЖЕНИЯ
Эту главу мы посвуТтим основным определенияим и попытаемся
сформулировать ряд проблем новой области вычислительной техни-
ки— цифрового моделирования. С этой точки зрения она является
как бы введением к остальным главам книги.
Резкое усложнение систем автоматического регулирования про-
изводствешых процессов, систем управления сложными технодоги-
ческими и военными ком1Илекса|ми, а также появление комплексных
систем, составным звеном которых является человек-ошератор, при-
вели к необходимости специальных средств воспроизведения процес-
сов в этих системах. Такие средства получили o6u.iee название мо-
делей физических систем. Развитие техники моделирования привело
в свою очередь к появлению так называемых аналоговых модели-
рующих вычислительных мащш.
Благодаря своей простоте, удобству эксплуатации и легкости
программирования аналоговые вычислительные машины являются
идеальным средством моделирования. С их помош.ью удалось ре-
шить огромное количество задач, связанных с разработкой всевоз-
можных систем автоматического управления. Своего значения ана-
логовые вычислительные машины ие потеряли и в настоящее вре-
мя. С другой стороны, развитие техники выдвигает сейчас требова-
ния к построению систем, процессы в которых требуется исследо-
вать с точностью, значительно превышающей предельную точность
аналоговых вычислительных машин. В настоящее время известно
большое количество задач, которые принципиально не могут быть
решены на аналоговых вычислительных машинах из-за ограничен-
ной точности последних. Таким образом возникла проблема оты-
скания средств моделирования, обеспечивающих практически не-
ограниченную точность вычислений.
Появление в 40-х годах универсальных цифровых вычислитель-
ных машин, казалось бы, должно было разрешить эту проблему.
Однако здесь возникла вторая трудность, связанная с тем, что бы-
стродействие универсальных цифровых вычислительных машин,
обладающих классической структурой, принципиально ограничено,
в то время как для современной техники характерно постоянное
убыстрение процессов в системах. Таким образом, была поставлена
7

вторая ороблема — построение вычислительных средств, работаю-
щих в натуральном масштабе времени.
Многочисленные неудачные попытки использовать для моделиро-
вания в натуральном масштабе времени универсальные цифровые
вычислительные машины или комбинированные .вычислителыные ком-
плексы, состоящие из уни'версальных цифро'вых и аналоговых вы-
числительных машин, /привели к достаточно прочно установившему-
ся в литературе мнению, что цифровая техника вообще непригодна
для создания систем, работающих в натуральном масштабе време-
ни. Существующие средства вычислительной техники принято делить
на быстродействующие и неточные, т. е. аналоговые, и медленные,
но точные, т. е. цифровые.
Подобное мнение справедливо применительно к универсальным
цифровым вычислительным машинам, но неверно в общем случае.
Одна из основных задач данной к}гиги состоит в том, чтобы пока-
зать, что п'ри правильном подходе к построению структуры цифро-
вой вычислительной машины можно получить модель некоторого
процесса, граничные возможности которой в смысле быстродействия
и точности близки к пределам, определяемым теоремой Котельни-
кова. Практически это означает, что, используя современные цифро-
вые элементы, работающие на частотах следования такто'вых сиг-
налов лорядка нескольких мегагерц, можно получить модель, вос-
производящую с точностью 5—6 десятичных знаков процессы, соб-
ственные частоты которых лежат в пределах от нуля до неакольких
сотен и даже тысяч герц, что вполне достаточно для удовлетворе-
ния требований современной практики.
2. АНАЛОГОВЫЕ И ЦИФРОВЫЕ МОДЕЛИ
Метод моделирования, широко используемый в инженерной
практике, основывается на том, что вместо исследования процессов
в некоторой исходной динамической системе исследуются процессы
в другой системе, называемой моделью и подобной в некотором
определенном смысле исходной системе, подлежащей исследованию.
Таким образом, процесс моделирования предусматривает решение
двух самостоятельных задач. Первая из них — это п ос т р о е н и е
модели, т. е. создание физической системы, удовлетворяющей тем
или иным критериям подобия [Л. 44]. Вторая задача состоит в по-
лучении от модели необходимой информации. Обычно эта задача
сводится к измерению физических величин, определяющих состоя-
ния модели.
Способ решения каждой из этих задач и определяет тип полу-
чаемой модели. В случае динамических систем, а именно эти систе-
мы будут нас интересовать в дальнейшем, происходящие в них
процессы описываются рядом переменных величин, изменяющихся
во времени, т. е. функциями времени. Сама исходная система при
этом может быть описана дифференциальным уравнением. Системы
L сосредоточенными .параметрами описываются обыкновенными диф-
ференциальными уравнениями, а системы с распределенными пара-
метрами—дифференциальными уравнениями в частных производ-
ных. В этой книге мы ограничимся только рассмотрением способов
цифрового моделирования систем, описываемых обьшювенными
дифференциальными уравнениями. Мы ограничимся также расомо-
8

трением опособо.1з так называемого матемапгчесжого (Моделирования,
* г. е. будем предполагать, что уравнения, описывающие работу ис-
ходной моделируемой системы, всегда известны.
Попытаемся теперь установить, ib чем состоит основное разли-
чие между аналоговыми и цифровьгми средствами моделирования.
При аналоговом моделировании физическому процессу в исходной
моделируемой системе или в какой-либо ее части ставится в соот-
ветствие также физический процесс, происходящий в модели.
Другими словами, аналоговая модель также ошисывается диффе-
ренциальным уравнением. Это дифференциальное уравнение может
отличаться от исходного дифференциального уравнения, описываю-
щего моделируемую систему, только масштабами представления со-
ответствующих переменных. Наиболее существенным с нашей точки
зрения здесь является то обстоятельство, что функции времени, опи-
сывающие процессы в модели, с той или иной степенью точности
совпадают с функциями времени, описывающими соответствующие
процессы в исходной системе во все без исключения моменты
времени.
В общем случае аналоговая моделирующая машина состоит из
ряда отдельных элементов, работа каждого из которых может быть
описана некоторой математической зависимостью между его вход-
ными величинами, выходньгми величинами и временем или другой
независимой переменной К В этом смысле говорят, что аналоговый
элемент выполняет некоторую определенную математическую опера-
цию. При этом также наиболее характерным является то обстоя-
тельство, что значения входных и выходных величин в любой без
исключения момент времени в пределах рассматриваемого интерва-
ла представляются значениями соответствующих физических вели-
чин, определяющих состояние элемента.
Итак, основу работы аналоговых вычислительных хМашин состав-
ляет использование физических (Процессов, отображающих процессы
в моделируемой системе во все без исключения (моменты времени
на интервале моделирования. В большинстве случаев оказывается,
что дифференциальные уравнения, описывающие работу модели,
в точности совпадают с дифференциальными уравнениями, описы-
вающими работу моделируемой системы. Погрешность же при мо-
делировании определяется как неизбежным пренебрежением некото-
рыми факторами при составлении дифференциальных уравнений
исходной системы, так и несовершенством элементов самой модели
(дрейф операционных усилителей, нестабильность параметров дета-
лей и т. п.). Результатом работы аналоговой модели являются зна-
чения физических величин, изменяющихся иногда достаточно быст-
ро. Для получения необходимой информации эти значения должны
быть измерены. Таким образом возникает третья причина погреш-
ности, связанная с конечной точностью измерений.
Цифровые вычислительные машины, так же как и аналоговые,
состоят из отдельных элементов. В каждом таком элементе под
влиянием внешних воздействий, представляющих собой входные сиг-
налы этого элемента, также возбуждается некоторый физический
процесс. Однако значения входных и выходных переменных каждо-
го элемента рассматриваются только в отдельные дискретные мо-
* у аналоговых вычислительных машин, построенных на основе электрон-
ных операционных усилителей, независимой переменной может служить толь-
ко время.
9

моиты времени, отстоящие друг от друга на конечные интервалы.
Длительности этих интервалов выбираются таким образом, чтобы
переходные .процессы, возникающие в элементе под .влиянием внеш-
него воздействия, в основном успевали затухнуть. Таким образом,
состояние каждого цифрового элемента определяется значением фи-
зической величины в установившемся режиме. С этой точки зрения
можно сказать, что состояние элемента отражает не сам тфотекав-
ший в нем физический процесс, а только результат этого процесса.
Обычно цифровому элементу ставится в соответствие некоторое ко-
нечное количество состояний и соответственно область возможных
изменений его выходного сигнала разбивается на конечное количе-
ство конечных интервалов. Это практически полностью исключает
погрешность измерений, поскольку считается, что элемент находится
в данном С0СТ0ЯН1П1 в том случае, если значение его выходного сиг-
нала оказывается заключенным в пределах некоторого достаточно
большого интервала.
Таким образом, если аналоговая вычислительная машина или
аналоговая модель устанавливает требуемую зависимость между не-
которыми функциями, определенными на всем множестве то'Чек дан-
ного отрезка оси независимой переменной, то цифровая вычисли-
тельная машина или цифровая модель устанавливает зависимость
между функциями, определенными на дискретном множестве точек
оси независимой .переменной, разделенных конечными интервалами.
На всем .протяжении каждого такого интервала, за исключением
его границ, функция не определена, а цифровая модель не дает ни-
каких сведений как о своем состоянии, так, следовательно, и о со-
стоянии исследуемой системы. Подобные функции получили назва-
ние решетчатых функций [Л. 46].
Одной из возможных форм установления зависимости между
решетчатыми функциями являются разностные уравнения [Л. 10].
Таким образом, если работа аналоговой вычислительной машины
в общем случае описывается дифференциальными уравнениями, то
работа цифровой вычислительной хмашины может быть описана
только разностным уравнением или, как иначе говорят, уравнением
в конечных разностях. Из сказанного следует, что основное разли-
чие между аналоговыми и цифровыми вычислительными машинами
заключается в способе представления соответствующих величин.
Если в случае аналоговой машины величины представляются непре-
рывными функциями независимой переменной, то в случае цифро-
вой вычислительной машины эти величины могут быть представле-
ны только решетчатыми функциями.
Присвоение каждому цифровому элементу конечного числа со-
стояний позволяет практически полностью исключить аппаратурную
погрешность и погрешность измерения. Поэтому цифровым вычисли-
тельным машинам и, в частности, цифровым моделям свойственны
только два основных вида погрешностей. Это так называемые ме-
тодическая погрешность и погрешность округления (квантования).
Причиной появления методической погрешности является то об-
стоятельство, что в общем случае при заданных начальных и гра-
ничных условиях решения разностного уравнения, описывающего ра-
боту цифровой модели, могут только приближенно совладать
с соответствующими решениями дифференциального уравнения, опи-
сывающего работу исследуемой системы. Причина появления погреш-
ности округления состоит в том, что с .помощью конечного числа со-
10

стоянии одного элемента или совокупности элементов можно пред-
ставить значение любой величи'ны только с конечной точностью.
Это значение будет всегда представляться числом с конечным коли-
чеством цифр. К рассмотрению этих видов погрешностей мы будем
неоднократно возвращаться на страницах этой книги.
Теперь мы достаточно подготовлены к тому, чтобы дать опре-
деление понятию «цифровая модель». Пусть имеется некоторая фи-
зическая система, исследование которой предполагается проводить
методом ^моделирования. Пусть далее состояние (Этой системы в лю-
бой момент времени может быть полностью определено заданием г
переменных величин, представляющих собой функции времени или
другой независимой переменной. В этом смысле говорят, что иссле-
дуемая система является г-мерной или обладает г степенями сво-
боды.
Цифровой моделью г-мериой системы мы будем называть в об-
щем случае некоторую другую физическую систему, обладающую
следующими характерными особенностями:
1. В пределах этой системы можно выделить по меньшей ме-
ре г величин, представляющих собой решетчатые функции времени
или другой независимой переменной. Эти величины мы будем
в дальнейшем называть координатами модели.
2. Между координатами модели существует зависимость, опи-
сываемая разностным уравнением. Решение этого разностного урав-
нения при соответствующих начальных и граничных условиях с не-
которой наперед заданной степенью точности совпадает с решения-
ми дифференциального уравнения, описывающего работу моделируе-
мой системы.
3. Все переменные в цифровой модели представляются цифро-
вым способом в указанном выше смысле этого слова.
Подобное определение является достаточно широким. Ему бу-
дет отвечать, в частности, и универсальная цифровая вычислитель-
ная машина, запрограммированная для решения уравнений, описы-
вающих работу исследуемой системы. С другой стороны, всякую
цифровую модель можно рассматривать .как вычислительную маши-
ну, решающую данное уравнение или класс уравнений. Однако по-
лучить требуемые характеристики цифровой 1М0дели, в первую оче-
редь в смысле точности и быстродействия, в большинстве случаев
удается только в том случае, если полностью используются преиму-
щества, даваемые сп е ц и а л и з а ц lie й цифровой модели. При
этом возникают две самостоятельные задачи, к рассмотрению ко-
торых мы и переходим.
3. ЗАДАЧИ ЦИФРОВОГО МОДЕЛИРОВАНИЯ
Общую проблему создания цифровой модели данной системы
можно в свою очередь разделить на две самостоятельные задачи.
Первая из них состоит в построении разностного уравнения, обла-
дающего свойства^ми, отмеченными в .данном выше определении.
При этом чрезвычайно существенным является то обстоятельство,
что в общем случае одному и тому же дифференциальному уравне-
нию можно поставить в соответствие несколько разностных уравне-
ний, причем решения этих уравнений при заданных начальных и
граничных условиях будут по-разному совпадать с решениями ис-
11

ход.ного дифференциального уравнения, отвечая при это-м выбран-
ному критерию (например, критерию минимума модуля абсолютно-
го отклонения, критерию минимума среднеквадратичной ногрешно-
сти и т. д.). Таким образом, перед конструктором цифровой модели
возникает проблема выбора одного из этих разностных уравнений,
в наилучшей степени отвечающего его целям. Получение требуемых
характеристик модели в большой степени зависит от правильности
такого выбора.
Эта задача, названная в дальнейшем задачей синтеза алгоритма
цифровой модели, будет рассмотрена в гл. 4 настоящей работы.
Естественно, что выбор алгоритма или разностного уравнения воз-
можен только в том случае, если известны некоторые его характе-
ристики, имеющие смысл качества. Вопросам оценки качества алго-
ритмов цифровых моделей посвящена гл. 5.
Разностное уравнение, описывающее работу цифровой модели,
устанавливает в свою очередь связь между работой отдельных ее
элементов, т. е. в известной степени определяет схему соединений
между этими элементами, или структуру цифровой модели. Однако
в общем случае одно и то же разностное уравнение может быть
реализовано несколькими различными структурами. Таким образом,
возникает вторая задача, а именно задача синтеза структуры циф-
ровой модели. Для успешного решения этой задачи, так же как и
в предыдущем случае, необходимо выделять некоторые характери-
стики структуры, имеющие смысл ее качества.
Совместное решение обеих рассмотренных задач, как это будет
показано ниже,' позволяет получить структуру цифровой модели,
наилучшую в некотором определенном смысле. В этом и состоит
преимущество, даваемое специализацией. Ограничение функций циф-
ровой модели решением одного уравнения или небольшого количе-
ства сходных между собой уравнений позволяет ставить саму зада-
чу синтеза наилучшего алгоритма и наилучшей структуры. Боль-
шую роль при этом играет конечно выбор величины, принимаемой
за критерий качества модели в целом.
Рассмотренные соображения позволяют нам сформулировать
еще одну характерную особенность цифровых моделей. Эта особен-
ность состоит в том, что цифровые мюдели должны быть специали-
зированы по структуре, т. е. обладать структурой, яри данных усло-
виях наилучшей или близкой к наилучшей в смысле одного из кри-
териев. Конечно ие все цифровые модели всегда обладают этим
свойством. Однако caiMa возможность постановки подобной задачи
позволяет выделить их в самостоятельный класс машин, отличный
от класса универсальных цифровых вычислительных машин.
В заключение этой главы отметим одно чрезвычайно важное
обстоятельство. Как уже говорилось, переменные на входе и вы-
ходе цифровой модели представляют собой решетчатые функции или
совокупности значений, взятых в дискретных точках оси независи-
мой переменной К Каждая решетчатая функция в свою очередь
определяет бесконечное множество непрерывных функций, а имен-
но все функции, значения которых в узловых точках совпадают со
значениями данной решетчатой функции, а на интервалах между
этими узловыми точками могут сколь угодно сильно отличаться
друг от друга. В частности, такими функциями могут служить лю-
' Подробнее о решетчатых функциях см. [Л. 46].
12

бые интерполяционные полиномы, йостроенные на имеющихся зна-
чениях решетчатой функции. Рхтественно, что поскольку цифровая
модель моделирует непрерывный процесс, то среди всех рассмотрен-
ных фун?кцнй каждый раз необходимо выбирать одну и считать ее
выходным сигналом модели или описанием процесса в моделируе-
мой системе. В качестве такой функции можно выбрать любую
функцию, представляющую собой результат интерполяции между
значениями решетчатой функции. Однако в этом случае выбранная
функция 'будет зависеть не только от вида разностного урашения,
граничных и начальных условий, т. е. от той информации о моде-
лируемой системе, которая введена в модель, но и от способа ин-
терполяции.
С другой стороны, как показано в гл. 5 , среди всего множе-
ства функций, совпадающих в узловых точках с данной решетчатой
функцией, существует одна, обладающая тем свойством, что она
полностью определяется значениями решетчатой функции. Послед-
нее свойство делает совершенно естественным выбор именной этой
функции — функции Котельникова в качестве выходного сигнала
цифровой модели. Именно она является той единственной функцией,
которая зависит только от разностного уравнения, граничных и
начальных условий и, следовательно, несет в себе только 'информа-
цию о моделируемой системе, первоначально введенную в модель.
Как известно из [Л. 45], функция Котельникова всегда имеет
ограниченный спектр. Отсюда мы можем сделать вывод, имеющий
существенное значение не только для .цифровой модели, но и для
всякой вычислительной машины. Сущно'сть этого вывода состоит
в том, что входным и выходным сигналом любой цифровой систе-
мы может быть только функция с ограничеиньш спектром, причем
ширина этого спектра определяется величиной интервала между
дискретными точками, который в свою очередь зависит от конструк-
тивных особенностей машины.
На первый взгляд последнее обстоятельство резко ограничивает
класс уравнений, которые можно решать на цифровых вычислитель-
ных машинах и, в чаотнюсти, на цифровых моделях. Однако при-
менительно к моделированию реальных физических систем, содер-
жащих инерционные элементы, подобное ограничение не является
существенным, поскольку спектры функций, описывающих прО'Цессы
в таких системах, также всегда оказываются ограниченными или, во
всяком случае, ограниченной оказывается та область частот, в ко-
торой амплитуда этих спектров существенно отличается от нули.
Поэтому вопрос может стоять только о согласовании граничных ча-
стот опектр'ов модели и моделируемой системы. Это и есть другая
формулировка задачи обеспечент1я возможности работы в реальном
масштабе времени.

ГЛАВА ВТОРАЯ
ЦИФРОВЫЕ ДИФФЕРЕНЦИАЛЬНЫЕ
АНАЛИЗАТОРЫ
4. ИСХОДНЫЕ СООБРАЖЕНИЯ
Цифровые диффереидиальиые анализаторы (ЦДА) занимают
особое место .среди прочих средств вычисл.ител1^ной техники. К на-
стояи^ему (Времени существует уже довольно обширная литература
[Л. 1, 15, 17, 19, 29], посвященная описанию конструкций и особен-
ностям экоплуатации машин этого типа. К сожалению, в большин-
стве известных aiBTopy (Публикаций по вопросу об ослювных характе-
ристиках и, чго особенно .важно, о предельных вычислительных воз-
можностях ЦДА в?.1сказывается м'ного лрогиворечиных и шдчас не-
верных М'нений. В отдельных работах допускаются также неточности
при onncaiHHH (Принципа действия ЦДА.
В связи с этим при написании данной главы автор ставил перед
собой следующие задачи:
а) Раз'обраться в вопросе о том, что такое ЦДА и в какой ме-
ре к машинам этого типа применимы основные положения теории
обыкновенных (аналоговых) дифференциальных анализаторов.
б) Дать оценки таким характеристикам ЦДА, как быстродей-
ствие и точность (Вычислений, и по возможности выявить связь меж-
ду этими пар.аметра.ми.
в) Оценить предельные вычислительные возможности ЦДА и хо-
тя бы приближенно очертить области -возможных и целесообразных
их применений.
При этом сознательно был опущен ряд второстепенных техни-
ческих .подробностей (системы счисления, способы образования при-
ращений и т. п.), которые достаточно хорошо освещены в литера-
туре.
Важность перечисленных выше вопросов обусловлена тем, что
количество ЦДА, уже находящихся (В эксплуатации или готовящих-
ся к выпуску, в настоящее .время непрерывно растет. Кроме того,
анализ перечисленных .проблем позволит нам на примере ЦДА вы-
явить некоторые характерные особенности машин, объединяемых
здесь под общим названием дифровых .моделей.
Во избежание недоразумений под ЦДА мы будем понимать
цифровую вычислительную машину, обладающую следующими отли-
чИТельными признака ми:
1. В составе машины можно выделить ряд самостоятельных ре-
шающих блоков, выполняющих вычисления по жесткому, неизме-
няемому алгоритму.
2. Единственным алгоритмом .для всех без исключения решаю-
щих блоков ЦДА является разностное уравнение одного определен-
ного вида, позволяющее реализовать одношаговую операцию при-
ближенного численного решения обыкновенных дифференциальных
уравнений и операцию суммирования переменных. Подобные блоки
получили название «цифровых интеграторов».
3. Для решения задачи цифровые интеграторы соединяются
между собой до одределенной схеме, называемой схемой набора.
14

4. Для кодирования данных, передаваемых по каналам связи
между интеграторами, используется метод разностно-дискретной мо-
дуляции (РДМ [Л. 3]).
Последнее обстоятельство является особенно 1важ,ным, посколь-
ку, как будет показано ниже, оно в большой степени определяет
основные характеристики ЦДЛ.
В литературе был описан также ряд машин [Л. 27, 33], близких
к ЦДЛ по своей структуре и назначению, но не полностью отве-
чающих приведенному определению. Поскольку эти отличия резко
из'меняют некоторые их параметры, нецелесообразно рассматривать
подобные машины в классе ЦДА. Описанию двух машин этого типа
будет посвящена следующая глава.
5. ПРИНЦИП ДЕЙСТВИЯ ЦИФРОВОГО ИНТЕГРАТОРА
Несмотря на большое разнообразие конструкций современных
ЦДА, схемы их основных блоков — цифровых интеграторов — отли-
чаются друг от друга только во второстепенных деталях. Поэтому
,1 I
Рис. 1. Типовая схема цифрового интегратора.
для рассмотрения принципа действия цифрового интегратора вос-
пользуемся типовой схемой, показанной на рис. 1. Цифровой инте-
гратор устанавливает некоторую определенную зависимость между
выходной переменной величиной, которую мы в дальнейшем будем
обозначать буквой г, и рядом входных переменных величин, kotoj
рые мы будем обозначать буквами у\, У2 .. ., Ус и х. Вид этой
зависимости таков, что она позволяет, кроме основной математиче-
ской операции — нриближенного интегрирования пли приближенного
решения — обыкновенного дифференциального уравнения первого
порядка, реализовать с помощью цифрового интегратора (прибли-
15

же1М1о или точно) также ряд других математических операций: сум-
мирование, умножение на постоянную и т. п. Это первая особен-
ность цифрового интегратора ЦДА, в известной степеии затр}'1Д1Няю-
щая по'ндмание принципа его действия.
Вторая особенность цифровото интегратора состоит в том, что
его работа вылолняетоя в виде последовательности одинаковых по-
вторяющихся циклов. Каждый такой цикл, называемый лтерацл-
е й, совершается во времени и имеет некоторую длительность Т,
обычно называемую длительностью, или временем, одной итерации.
Длительность итерации зависит от большого количества причин и,
в частности, в одном и том же ЦДА может изменяться в зависи-
мости от вида решаемой задачи. Однако она не может быть мень-
ше некоторой величины Гмин, определяемой конструктивными осо-
бенностями данного цифрового интегратора. Моменты начала каж-
дой итерации мы будем обозначать как to, tu . • •, tn - - • При
этом момент to совпадает с моментом начала решения задачи,
на ЦДА.
Наконец, третья особенность цифрового интегратора, уже от-
амечавшаяся выше, состоит в том, что как выходная леременная г,
так и входные переменные уи У2, . • Ус их кодируются по ме-
тоду РДМ. Практически это означает, что сигналы, передающие
значения перечисленных переменных, имеют смысл приращений, или
первых разностей, этих переменных и могут действовать только
в олределенные фиксированные хмоменты времени, а именно момен-
ты tu /2, .. tny ... При модуляции истинные значения разностей
любой переменной ^ округляются до одного из ближайших значе-
ний —б^, О или +^^, поэтому передающий это значение сигнал мо-
жет принимать только одно из трех фиксированных значений. Та-
ким образом, любая входная или выходная переменная величина
в цифровом интеграторе за время одной итерации или вообще не
изменяется (разность равна нулю), или изменяется на ±6^, причем
эти изменения могут пооисходить только в моменты времени tu
t2, ^n, ...
Конечно, все переменные величины, связанные с цифровым ин-
тегратором, изменяются во времени и, следовательно, являются
функциями времени. Однако из того обстоятельства, что измене-
ния этих величин могут происходить только в моменты времени
tu ^2, ..^п, ..., а интервалы Г, разделяющие эти моменты, как
уже отмечалось, могут изменяться и притом в достаточно широких
пределах, вытекает, что сам вид функциональных зависимостей пе-
ременных величин от времени также может изменяться при изме-
нениях режима работы ЦДА. Поэтому для лас более удобно (и бо-
лее правильно) считать^ что все переменные величины, связанные
с цифровым интегратором, представляют собой функции некоторой
независимой переменной я, представляющей собой порядковый но-
мер момента времени /п и принимающей только целочисленные
значения О, 1, 2, . . ., п, ... Другими словами, мы будем считать,
что выходной и входными переменными цифрового интегратора яв-
ляются решетчатые функции z(n), yi{n), У2(п), .. ., Ус(Щ, х{п)
(определение решетчатой функции см. в предыдущей главе). Нако-
нец, разность любой переменной | мы будем определять как
a5(AZ)^$(AZ-f l)-S(A2), (1)

а округленные значения разностей, передаваемые по каналам связи
между интеграторами, будем обозначать через А|(п):
М(п) = К^(п) + г^{п)у (2)
где (/г) —остаток, получающийся после округления.
Теперь можно перейти непосредственно к рассмотрению типо-
вой схемы цифрового интегратора, показанной на рис. 1. Эта схема
состоит из двух основных частей: сумматора входных величин СВ
и собственно цифрового интегратора ЦИ. Работу схемы мы про-
следим в течение одной итерации, начинающейся в момент вре-
мени tn-
Буквой Y на схеме, показанной на рис. 1, обозначен регистр,
хранящий последовательные 31начения некоторой переменной вели-
чины у. Непосредственно перед моментом tn в регистре Y хранится
очередное значение этой функции, образованное в момент tn-\,
т. е. у{п—\).
В момент на входы г/,, У21 - - -у Уо поступают входные сиг-
налы, имеющие смысл округленных разностей Ayi(n—\),
^y2[f^— 1),..., LyQ{n—\). Эти сигналы суммируются в сумматоре S
блока СВ. Результатом этой операции является чйЬло
G
^y{n-\) = Y,^yi(n^\). (3)
Число Аг/ (п—1) с помощью сумматора 2у блока ЦИ сумми-
руется с числом у(п—1), в результате чего образуется новое зна-
чение _
r/(Az) = r/(n-l)+ ^y(n-^\). (4)
Переменная величина у(п) является одной из основных величин
в цифровом интеграторе. Операции вида (4) выполняются в течение
каждой итерации. Кроме того, перед 'началом решения задачи в ре-
гистр Y вводится «начальное» значение величины (/ — число «/(0).
Следовательно, мы можем написать:
п
y[n) = y(0) + Y,'^y(i^^y
/=1
(5)
Одновременно с рассмотренной операцией число (/2—1) умно-
жается с помощью узла А блока ЦИ на некоторое постоянное
число а. Результат умножения с помощью сумматора Sa сумми-
руется с числом у(п—1). Таким образом, на выходе сумматора Sa
образуется число
уЧп)^у(п^\) + а'Еу(п^\), (6)
С выхода сумматора Sa число у"^ (п) поступает на входы венти-
лей ^1 и В2. Эти вентили управляются входным сигналом Lx(n — 1)
таким образом, что, если ~Кх(п— =Ь^, вентиль 5, открыт
2 А. в. Шилейко. 17

в течение всей итерации, начинаюш;ейся в момент ^д, а вентиль В2
закрыт; если Ах (п—1) = — дх, то открыт вентиль 82 и закрыт
вентиль Bi, а если (п — I) = 0, то закрыты оба вентиля.
С выхода вентиля Bi число непосредственно поступает
на вход сумматора 2г. С .выхода же вентиля В2 число г/* (/г) посту-
пает на вход узла ОД, образующего дополнительный код этого
числа или, что то же самое, изменяющего знак этого числа на
обратный. Таким образом, когда открыт вентиль В2, на вход сум-
матора 2г поступает число —Наконец, когда оба вентиля
закрыты, на в.ход сумматора 'Lr поступает число 0.
Из сказанного следует, что работа вентилей Bi, В2 и узла ОД
может быть описана соотношением
^ Д-д. _ ].) y-i^ (д) =
= sgn Jx (/2 - 1) [у {п - 1) + аКу {п — 1)], (7)
где _
1 при Ах = + дх\
sgn Ах= о при Ах = 0;
— 1 при Дл: = — дх.
Буквой R на рис. I обозначен числовой регистр, аналогичный
регистру Y. В этом регистре хранятся последовательные значения
некоторой переменное величины . В частности, в момент tn в
этом регистре хранится очередное значение этой переменной
г^^{п — ]). С помощью сумматора число г^^(/2 —1) суммируется
с числом z/*"^ (п), и в результате этого образуется новое число,
которое мы обозначим через
ДгЛп) = г^„(я-1) + г/**(/г). (8)
Число Аго(/^;) по существу представляет собой окончательный
результат всей совокупности операций, выполнявшихся в цифровом
интеграторе в течение данной итерации. Однако в соответствии
с принятым в ЦДА способом кодирования входных • и выходных
сигналов [см. соотношение (2)] число Azo{n) разделяется на две
самостоятельные части
Az,{n) = ~Az,{n) + r^^(n), (9)
причем часть г^^{п) помещается в регистр R на место числа
f^^{n — 1) и хранится там вплоть до следующей итерации. Число
AzQ(n) передается на выход цифрового интегратора. Рассмотренная
операция округления выполняется специальным узлом округления
(узел О на рис. 1).
Комбинируя теперь соотношения (3), (5), (9), выведем основ-
ные уравнения, описывающие работу цифрового интегратора. Эти
уравнения имеют вид
а G
У{п) = у{0) + "^УАу^(1-\) (10)
/=1/=1

где 5 — разряд з.нака, а X — значения разрядов числа, представлен-
ного в некоторой В-ричной системе счисления {В — основание си-
стемы счисления). Вес разряда ,в каждой данеой позиции может
быть определен только после того, как будет введено условие о ме-
сте размещения В-ричной залятой.
Предположим, что 5-рнчная запятая в числах у(п) расположена
так, что вес старшего разряда равен Л^^, а вес младшего разряда—
QM-N-v\ Здесь Л^ —общее количество разрядов в числе у{п) (не
считая знака), а М — некоторое, пока произвольное, целое число.
Размещая таким образом В-ричную запятую, ;мы определяем базу
для отсчета осталыных 1масштабов. Заметргм также, что число раз-
рядов в цифровых интеграторах ЦДА не является фиксирован-
ным. Лри решении каждой данной задачи величина 1выбирается
для каждого интегратора. Ограничено только максимальное значе-
ние числа Л^. Это (Максимальное значение является общим для всех
цифровых илтеграторов данного ЦДА.
При образовании if^{n) числа [у{п—^\) + a^y — умно-
жается на ± 1 или 0. Следовательно, вес старшего разряда числа
у-^-^{п) также равен В^^. При суммировании чисел у'''* (п) и г^^(п—\)
можно складывать только разряды с одинаковыми весами. Следо-
вательно, старший разряд числа г^^(п — 1) также должен иметь
вес В^,
Из сказанного следует, что значения переменных у и г^^ пред-
ставля10тся в цифровом интеграторе iV-разрядными числами, причем
веса старших разрядов чисел у и г^^ одинаковы и равны В'^К Зна-
чения переменной г/*^ также представляются yV-разрядными чис-
лами с тем же весом старшего разряда. При суммировании чисел
г^^ и f/** в общем случае образовывается N + 1-разрядное число,
вес старшего разряда которого равен В^'^^'^^\
2^= 19
Уравнения (Ш) и (;\\\) полностью и однозначно описывают ра-
боту цифрового И'нтегратора. Это линейные разностные ура.В1нения,
заиисаипые относительно переменлых z, уи у2, ..., Уо, ^
До сих пор мы не вводили никаких условий относительно пере-
менных Z, уи У2, ..., Ус, X, считая про"сто, что это решетчатые
функции независимой переменной п, 1причем последовательные зна-
чения inepeMeiHHbix у и г хранятся в регистрах Y я цифрового
интегратора. Расомотри.м теперь более .подробно вопрос о масшта-
бах представления этих переменных.
Числа, хранящиеся в регистрах Y и iR цифровых интеграто-
ров ЦДА, представляются в системе счисления с фиксироБанной за-
пятой. Другими словами, каждое такое число может быть записано
в виде

Способ кЁантования выходного сигнала АгоН, наиболее часто
используемый в ЦДА, можно описать следующими соотношениями^:
sgnA^o (п) =
+ I при Azo(/2)>L;
О при 0< Azo(n)<L;
- 1 при Аго (п) < 0;
\Azoin)\= (12)
где L = max | 1 — максимально возможное значение остатка г^^.
Согласно сказанному выше
Из соотношений (9), (12) и (13) можно сделать также вывод,
что числа г^^ всегда положительны. Следовательно, разряд знака
в регистре R не используется.
Заметим далее, что величину М нельзя выбирать произвольно.
Она определяется полным числом разрядов в регистрах R и Y
"и весом сигналов Af/j(/i —1). Представим каждый сигнал Af/j (/i—1)
в виде
Ayj (п — \) = sgnlyj (Д — 1).| lyj I,
причем, очевидно,
I I =
где Byj —принятый шаг квантования величины yj.
Для цифровых интеграторов ЦДА устанавливается условие
5у1 = 5у2 = ... = (14)
где My — некоторое постоянное число, называемое масштабом ве-
личины Ai/j.
При таких условиях вес младшего разряда числа Af/(n —1) =
G
= ^ Af/j (/г — 1) будет равен Б Поскольку этот младший раз-
ряд в схеме цифрового интегратора суммируется с младшим раз-
рядом числа у(п — 1), имеем:
M-N + \ ==-Му
или
\Iz,{n)\ = B^^^'^=d''-''y\ (15)
' Это так называемая тернарная система образования приращений.
В некоторых конструкциях ЦДА использовалась также другая система — би-
нарная, однако она приводит к большим погрешностям и поэтому здесь не
рассматривается.
20

Подставляя равенство (15) в уравнение (И) и умножая обе часТН
этого уравнения на |Ах|, получим:
sgn Azo (п) \Az,(n)\ i А ^ i + \~Ех \ г^^(п) = \ Ах \ г^^{п - \) +
+ [у{п-'1) + аЛу (п _ 1)] sgn Ajc (А2 - 1) 1 Ал: |. (16)
Обычно в ЦДА величина | Алг 1 выбирается равной
|AJC|=J5~^^ (17)
где Мх — масштаб величины Ал". Подставляя равенства (15) и (17)
в уравнение (16) и вводя обозначения
Ai (/7) - sgn Azjn) i А"го I i Ал: i;
гЛп) = \'Кх\г^^(п)у (18)
получим окончательно:
Аг(п) + г,[п) = Гг(п-\) + [у(п-\) + аЛу(п-\)]1х{п-\у (19)
Уравнение (19) совместно с уравнением (10) полностью описывает
работу цифрового интегратора с учетом масштабов представления
переменных величин. Полагая | Az (л) | = 5 получим из (18) еще
одно важное соотношение, называемое соотношением масштабов\
Мг = М^ + Му — N, (20)
где Мх — масштаб представления величин Az.
Сделаем теперь одно замечание, которое представляется нам
существенным. Суммируя по аргументу обе части уравнения (19) и
п
обозначая z(n+\) = z (0)+ ^Аг (/), где Az (/) =а2 (/) +Гг (1), получим:
i=\
п
2{n+\)=z(0) + rx(n-\)+'£ (/)А"л: (/). (21)
i = l
Из сопоставления выражения (21) с определением интеграла Стнль-
тьеса
b b
\ У {X) dg (X) = у(а)+ Urn Уу{х,) ^g (л:,) (22)
а а
некоторые авторы [Л. 36] делают вывод, что цифровой интегратор
вычисляет определенный интеграл в смысле Стильтьеса. Это утверж-
дение швершенно неверно, поскольку никакой предельный переход
здесь невозможен (/г+1—/2=1), а х — произвольная переменная ве-
личина, поступающая в интегратор извне.
^ в ЦДЛ можно устанавливать определенное значение величины N
в каждом интеграторе.
21

Кроме того, подобное отождествление только запутывает и без
того сложный BOinpoc. Работа цифрового интегратора с точностью до
погрешностей округления ^ описывается раз,носгны,м!и уравнениями
(1Q) и (19i). |Прн определенных дополнительных предположениях
решения этих уравнений могут приближенно пли точно совпадать
с некоторыми функциями одной или более переменных. Только
в это,м смысле мы и говорим, что цифровой интегратор выполняет
те или иные 'математические операции, к рассмотрению которых мы
и переходим.
6. РАБОТА ЦИФРОВОГО ИНТЕГРАТОРА В РЕЖИМЕ
РЕШЕНИЯ ОБЫКНОВЕННЫХ ДИФФЕРЕНЦИАЛЬНЫХ
УРАВНЕНИЙ
Этот режим является основным режи.мом работы цифрового
интегратора, поскольку ЦДА предназначены для решения обыкно-
венных диффереициальных уравнений. Для того чтобы получить воз-
можность исследовать работу цифрового интегратора в этом режи-
ме, рассмотрим общий подход к численному решению обыкновенного
дифференциального уравнения. Такой подход заключается в сле-
дующем.
Пусть необходимо решить дифференциальное уравнение по-
рядка г вида
zl^] = F (zlr - Ч, z[^-^],..., 2', z, X) (23)
при начальных условиях
zlr-Ч {0) = Zr-г{0);
z[r-4(0) = Zr.2m
ziO) = z, (0);
где z[i] — f-я производная функции z (х) по х.
Вводя подстановки
Z (х) = Z, {X)]
z'{x)^z, (х);
z[r-4 = Zr^i{x),
преобразуем уравнение (23) в систему уравнений:
2% - 1 = f г-1 (^г - 1» Zr -2f - - • у Zq, х)\
Zr-2 — fr-2{Zr~-lj Zr-2'> ' ' ' ■> Zq, X)\
z'o = fo{^r-u 2:.-2,...,2:o, x). (24)
Замечая, наконец, что все переменные Zr-u Zr-2...-» Zq в свою
очередь являются функциями переменной х, введем сокращенные
обозначения
.fi{Zr.-U Zr.2,'-',Zo, x) = yi{Zi, х)
' Вопрос о погрешностях округления будет подробно рассмотрен в гл. 5.
22

и запишем систему (24) в более простом виде
z'i=yi(Zi, х)\ i = 0, 1, 2,...,г-1. (25)
Функция Zi выделяется здесь потому, что относительно нее каждое
■из уравнений (24) оказывается неявным.
Таким образом, мы видим, что задача решения уравнения (23) сво-
дится к задаче одновременного совместного решения г уравнений
вида (25). Рассмотрим одно из уравнений (25) и запишем его (опу-
ская индексы) в виде
Z' = y(z, х)\
^(О)-^о. (26)
Введем теперь в рассмотрение последовательность точек х^,
ЛГ!,..., Хд,..., разбивающих ось независимой переменной х на ряд
равных интервалов:
Lx = Xi — Xq = Х2 — Xi = . . . = Xn + i — Xfi = . . .
предположим, что каким-либо способом нам удалось вычислить п
значении z(Xi), 2 (Xg),..., ^ (x„) функции z{x) в точках Хи Ха,...,
Хп- Тогда точное значение z{Xn+i) можно определить из выраже-
ния
2:(Xn + i) = ^(x„) + J y[z(s), s]ds. (27)
Вычислить точно интеграл >в правой части выражения (27) ib общем
случае невозможно, поскольку нам еще не известно ни одно зна-
чение переменной z(x) на интервале Хп+и Хп, кроме значения
в точке Хп. Поэтому все приближенные численные методы решения
обыкновенных дифференциальных уравнений представляют собой
попытки приближенно оценить величину этого интеграла. Здесь воз-
можны два способа. Первый из них, экстраполяционный, состоит
в том, что интеграл
j У [Z (S), S] ds
приближается линейной комбинацией предыдущих значений функ-
ции
Az(n)= J (/ [г (s), s] ds::^AxY^ b^y (n - v), (28)
где
Lz{n) = z(Xn+\) — z{Xn)\
6 (v)—постоянные коэффициенты, выбираемые из соображений по-
лучения наименьшей погрешности.
23

Второй слюсоб, лнтерполяционлый, лредусматривает замелу
фулщил y[z(x), х] на интервале Хп+ь Хп некоторым приближением
y[z{x), х]. Значения интеграла уточняются затем в результате по-
ел едювате льны х приближений. Второй способ иелригоден для лас,
поскольку, как было показано выше, цифровой интегратор выпол-
няет в течение одной операции только один цикл вычислений по
формуле (19). Поэтому мы будем рассматривать здесь только пер-
вый способ.
Ограничимся двумя членами суммы в правой части выражения
(28) и перепишем это выражение в виде
А2 (/г) = [у1п-\) + аАу {п - 1)] Ал:, (29)
т. е положим bi = \—a\ Ь^ = а.
Полагая в выражении (19)
~Kz{n) + r,(n) = bz{n)
и
Ах (I) = Ах = const,
мы видим, что выражение (29) совпадает с выражением (19) с точ-
ностью до слагаемого Гг(п — \)у имеющего, как было показано
выше, порядок величины | Ал: |.
Отсюда следует вывод, что при Ал: (/) = Ах = const с помощью
цифрового интегратора ЦДА можно реализовать один из прибли-
женных экстраполяционных численных методов решения обыкно-
венных дифференциальных уравнений порядка не выше 1 ^ 1
в формуле (28)J. Получаемая при этом точность вычисления зави-
сит от многих факторов и, в частности, от выбора величины а.
В современных ЦДА, описанных в литературе, для а выбирают
одно из трех значений: л = 1, а = v2 и а = ^4- Этому соответ-
ствуют три различные формулы численного интегрирования:
Az(n)=^r,(n-\) + y (п) Ах; (30а)
Az (п) = г, (A2 - 1) + у Ал: [у (п) + у(п-\)]; (306)
Az(n) = r,(n-l)+^Ax [Зу (п)-у{п-\)]. (ЗОв)
Различные свойства этих формул, а также соображения о целе-
сообразности выбора каждой из лих будут рассмотрены ниже.
7. СТРУКТУРЫ ЦДА
Из сказанного в предыдущем параграфе можно сделать вывод,
что единственной арифметической операцией, выполняемой в ЦДА,
является операция суммирования. На каждом шаге вычислений вы-
полняются две самостоятельные операции суммирования, [одна из
которых имеет целью образовать новое значение у(п) путем сумми-
24

рования предыдущего значения у(п — \) с „приращением^ Д^/(/2—1),
а вторая имеет целью образование нового значения AZo(n) в ре-
зультате суммирования чисел f^^in — 1) и i/**(/2). В соответствии
с этим можно представить себе схему ЦДА, аналогичную показан-
ной на рис. 2. Здесь Д — динамический регистр (например, дорожка
магнитного барабана) общей емкостью 2NQ, где — количество
разрядов в числах у и г^^; Q — общее количество интеграторов.
Числа записываются в регистре Д, как показано на рис. 2. В пер-
вую ячейку дорожки помещается младший разряд числа Уи отно-
Рис. 2. Структурная схема последовательного ЦДА.
сящегося К интегратору № 1. Во вторую ячейку записывается млад-
ший разряд числа r^i, также относящегося к интегратору № 1, и
т. д. В ячейку с номером 2Л^ + 1 записывается младший разряд числа
У2, о'Т'носящегося к интегратору № 2, затем младший разряд числа
Г2о2 и т. д.
Выход динамического регистра соединен с входом а однораз-
рядного двоичного сумматора. Выход переноса этого сумматора со-
единяется с входом с через линию задержки, осуществляющую за-
держ1К]у иа два элементарных тактовых интервала, т. е. на время
прохождения двух нчеек регистра Д. Выход суммы н вход b
сумматора соединены с коммутатором /С, в котором выполняются
следующие операции:
1. При поступлении на вход а сумматора очередного разряда
числа yi на вход b сумматора подается соответствующий разряд
числа Ar/i. Результат суммирования поступает в коммутатор К и за-
поминается в нем на время одного тактового интервала. Сигнал
переноса задерживается на два тактовых интервала. Он будет
действовать на входе с сумматора к моменту поступления следую-
щего разряда числа Уг.
2 При поступлении на вход а сумматора очередного разряда
числа г^^^ коммутатор посылает на вход b задержанный на один
тактовый интервал соответствующий разряд числа г/**,-. Результат
суммирования направляется коммутатором К на вход регистра Д,
а сигнал переноса снова задерживается на два тактовых интер-
вала,
25

3. После суммирования старших разрядов очередной пары чи-
сел r^^i^ и г/**г цепь переноса разрывается, а сигнал переноса,
имеющий здесь смысл сигнала переполнения или выходного сиг-
нала интегратора, посылается коммутатором К в блок набора про-
граммы БНП. В блоке БИП этот сигнал запоминается и хранится
до тех пор, пока он не будет использован в качестве сигнала Ау
в интеграторе с номером k.
Описанная схема относится к структурам чисто последователь-
ного типа. Она наиболее проста, так как, кроме необходимых во
всех случаях запоминающих ячеек, хранящих числа yi и r^oi' в ней
используется один единственный одноразрядный сумматор. В то же
время эта схема обладает наименьшим быстродействием. Длитель-
ность выполнения одного вычислительного цикла во всех интегра-
торах, соответствующая длительности одной итерации, 'очевидно,
равна 2NQx.
Величина т определяется типом .использованных элементов.
Если типы элементов заданы, то единственный способ увеличения
быстродействия состоит в увеличении числа параллельных каналов.
Обычно сумматор обладает большим быстродействием, чем дина-
мический регистр, поэтому вместо одного динамического регистра
мож,но использовать два вдвое меньшей длины: один для чисел //г,
другой —для чисел г^^^ . Единственный сумматор будет выполнять
при этом две операции суммирования в течение одного тактового
интервала. Быстродействие, естественно, повысится вдвое. Подоб-
ная структура осуществлена в ряде выпускаемых промышленностью
конструкций ЦДА.
В ЦДА, работающих с десятичной системой счисления, исполь-
зуются четыре параллельно включенных одноразрядных сумматора,
выполняющих операцию сложения четырехразрядных двоичных
чисел. Количество динамических регистров увеличивается при этом
до четырех или до восьми (на.пример, четыре регистра для числа
У1 и четыре—для числа r^J.
Увеличивая дальше число параллельных каналов, мы приходим
к схеме, содержащей N параллельно включенных одноразрядных
сумматоров. При этом] операция сложения одного числа yi'^
и одного числа г^^^ выполняется в течение одного тактового ин-
тервала. Число динам.ических регистров может быть равно здесь
или 2N. Это так называемая структура последовательно-параллель-
ного типа [Л. 19]. Подобная структура также применена в несколь-
ких конструкциях, выпускаемых промышленностью ЦДА {Л. 50].
Однако возможности увеличения числа параллельных каналов
еще далеко не исчерпаны. Следующий тип структуры — параллельно-
последовательный— будет содержать 2Q независимо работающих
одноразрядных сумматоров, каждый из которых замкнут на дина-
мический регистр емкостью разрядов. Цифровой интегратор будет
состоять при этом из двух сумматоров и двух динамических регист-
раторов (рис. 3). Интересно заметить, что в случае параллельно-по-
следовательной структуры каждый интегратор представляет собой
автоно1Мный блок. Здесь отпадает необходимость в использовании
блока набора программы, так как соединения между интеграторами
можно осуществлять шнурами, как это делается в электронных ана-
логовых вычис/гательных машинах. Структуры параллельно-последо-
26

вательного типа реализованы в ЦДЛ TRICE фирмы Paccard iBell Co.
(США) и в ЦДА AD(DAM-2 фирмы А. V. Roe (Англия) [Л. 49].
Наконец, схема, содержащая Q .групп по N параллельно вклю-
ченных сум'маторов в каждой, соответствует параллельно-параллель-
ной структуре. Такая схема обладает наибольшим быстродействием
и использует наибольшее количество оборудования. Структуры по-
добного типа исследовались экспериментально [Л. 48], однако ввиду
с
42,
Р\ Hi
с
■i/(n-0
Входной]
блок
Рис. 3. Структурная схема цифрового интегратора ЦДА
параллельно-последовательного типа.
чрезвычайно большого объема оборудования цромышленный выпуск
параллельно-параллельных ЦДЛ, насколько известно, не предпри-
нимался.
Кроме рассмотренных основных типов структур, естественно,
возможны любые промежуточные типы. Для каждото ЦДА, состоя-
щего из Q интеграторов и предусматривающего выполнение опера-
ций над Л^-разрядными числами, существует всего NQ различных
физически реализуемых структур, каждой из которых соответствует
определенное значение числа параллельных каналов F. При этом
точность ЦДА не зависит от типа структуры, а определяется только
принципом действия и количеством разрЯ1Дов N. С другой стороны,
такие важные характеристики ЦДА, как сложность и быстродейст-
вие, пеликом определяются числом параллельных каналов, т. е.
структурой. Отсюда возникает идея выбора некоторого оптимального
значения для числа параллельных каналов. Один из подходов к ре-
шению подобной задачи описан в работе [Л." 49].
8. РАБОТА ЦИФРОВОГО ИНТЕГРАТОРА
В СПЕЦИАЛЬНЫХ РЕЖИМАХ
В этом параграфе мы рассмотрим работу цифрового интеграто-
ра в режимах сумматора, нуль-органа, устройства умножения на
постоянную и устройства для вычисления определенных интегралов
Эти и ряд подобных режимов достаточно подробно описаны в ли
тературе [Л. 1, 29, 40, 49]. Здесь это рассмотрение проводится толь
ко в той мере, в какой это будет необходимо для дальнейшего ана
лиза. Для удобства мы будем пользоваться условным обозначением
цифрового интегратора, показанным па рис. 4.
27

Сумматор. При работе цифрового интегратора в режимах сум-
матора и нуль-юргаиа используется свойство чисел, хранящихся в ре-
гистре У, называемое цикличностью. Пусть в какой-либо момент вре-
мени IB регистре Y оказывается записанным число г/макс, раннее сво-
ему максимально воз^можному лоложительному значению. Это число,
очевидно, имеет вид 0,il 11 ... 1 (здесь и в дальнейшем мы будем
предполагать, что числа в цифровом интеграторе представляются
0
С иг "ала Jz
'Начальное значение
Рис. 4. Условное обозначение цифрового интегратора.
В двоичной системе счисления. Однако все последующее рассмотре-
ние без каких-либо ограничений распростраляется ла случай любой
системы счисления. В частности, максимальное положительное зна-
чение числа у{п) в системе счисления с основанием В будет иметь
вид г/(^г)=0, XXX ... X, где Х = В—1. Если теперь прибавить к^шслу
^макс некоторое положительное и отличное от нуля число Ау{п),
имеющее вид О, 00 .., О YYY, где Y — разряды, могущие иметь значе-
ние либо О, либо 1, то в результате получится число '1, 0O...OYYY,
представляющее собой, очевидно, дополнительный шд отрицательно-
го числа — Аг/(л). Легко проверить и обратное утверждение: если
в регистре Y хранится максимально возможное отрицательное чис-
ло, имеющее вид у(п\) = 1, О.. . О, и к нему прибавляется некоторое
отрицательное число — At/(/zj) =1, U ...И YYY, то в результате по-
лучается положительное число +Ау(п). Обычно подобное явление
называется «переполнением» регистра У и при его возникновении ге-
нерируется специальный сигнал, останавливающий машину. Однако
при работе цифровых интеграторов в режимах сум|матора и нуль-
органа аигнал переполнения блокируется, а свойство лереполнения
используется для 1ВЫ1полнен1Ия операции.
При работе в режиме сумматора осуществляется схема соеди-
нений, показанная на рис. 5. Выходной сигнал цифрового интеграто-
ра передается на его вход. На остальные входы поступают слагае-
мые yi, У2> ■ >. , Уо-1. Перед началом решения задачи в регистр У
S
£=:М
Ушке
ВыХоб CUMMbi
28
Рис. 5. Схема включения для работы в режиме сумма-
тора.

вводитоя максимальное положительное число ^(0)= г/макс. На вход
1л: цифрового интегратора в каждую итерацию поступают ангналы,
имеющие смысл +Лл:.
Предположим теперь, что в течение некоторого количества
итераций значения входных сигналов уи У2У'--^Уо—х остаются рав-
ными нулю, а начальное значение (0) переменной (п) было вы-
брано отличным от нуля*. Тогда в первой же итерации образуется
выходной сигнал + Аг. В течение следующей итерации этот сигнал
прибавится к числу t/(0) и благодаря свойству цикличности преоб-
разует его в число г/(1) = —г/макс. В течение этой же итерации
образуется выходной сигнал —~Кг [действительно, согласно (19)
''г(0)+f/макс = + А2+Гг(1); (1 )= у Гамаке И (1) — Г/макс < 0]
и в течение следующей, второй итерации снова получим у{2) =
= + У макс-
Таким образом, при равенстве нулю входных сигналов yi =
= ^2. на выходе сумматора мы будем иметь последова-
тельность +А2, — Az, + Az, — Аг,..., в сумме также равную нулю.
Если теперь в начале некоторой /.й итерации сумма входных
сигналов SAf/i(/—1) окажется равной a\Ayi\ (а — целое), то в ре-
гистре У образуется число —a\Ayi\ и понадобится а итераций
(т. е. а передач сигнала + Az на вход интегратора), чтобы восста-
новить первоначальное значение г/макс в регистре Y. Следователь-
но, поступлению на вход цифрового интегратора, работающего
в режиме сумматора, группы входных сигналов, ^определяющих
в сумме положительное число а\Ау{\, будет соответствовать по-
следовательность из а сигналов + Az на его выходе. Так выпол-
няется операция суммирования. Легко проверить, что аналогичный
результат будет иметь место в том случае, если сумма входных
сигналов определяет отрицательное число — а \Ауг |.
Нуль-орган. Этот режим известен также в литературе под
названием режима „следящего интегратора" или „сервоинтегратора".
При работе в режиме нуль-органа в регистр У цифрового интегра-
тора перед началом работы также вводится максимальное положи-
тельное число f/макс, но цспь обратной связи с выхода на вх^д
здесь отсутствует. Предполагается, что нормально на входах Ау
нуль-органа действует нулевая последовательность входных сигна-
лов вида —- А(/, +~Дг/, — Аг/, —"Дг/,... Тогда 'последовательность
выходных сигналов согласно только что рассмотренному также
будет иметь вид + Az, — Az, + Az, — Az,... Однако достаточно
хотя бы одного нарушения входной последовательности, чтобы на
вы_ходе установилась либо максимальная положительная (+Az, +Az
+Az, 4-Az,...), либо максимальная отрицательная (—Az, —Az,
—'Az, —- Az,...) последовательность. Так будет происходить до тех
^ Обычно перед началом работы в регистры R всех интеграторов ЦДА
вводятся «начальные значения», равные 0,5 от максимально возможных зна-
чений соответствующих чисел.
29

пор, пока указанное нарушение входной последовательности не
будет компенсировано. Работа в таком рел^нме будет эквивалентна
работе в качестве релейного усилителя.
Применение в схемах набора цифровых интеграторов, работа-
ющих IB режиме сумматора или нуль-органа, может слулсить причи-
ной двух весьма серьезных видов погрешности. Первый нз них со-
стоит в том, что сумматор н нуль-орган вносят задержку длитель-
ностью от 1 до G—1 итераций. Действительно, как ,мы видели при
рассмотрении работы сумматора, нри ноступлении на вход группы
сигналов, определяющих число a\Ayi\, последний нз сигналов Дг, пе-
редающих значение этой суммы; П'ОЯ!Внтся на выходе цифрового ин-
тегратора только через а итераций. Как будет показано ниже, на-
личие задержки в за'мкнутой вычислительной цепи приводит к по-
вышению порядка соответствующего разностного уравнения, а это
в свою ючередь может полностью исказить решение.
Другой вид погрешности обусловлен необходимостью представ-
ления нулевых значений сигналов чередующимися последователь-
ностями вида + Дг, — Дг, + Дг, — Дг, ... Это так называемая
«фазовая погрешность», которая в «отдельных случаях может дости-
гать весьма больших значений. Подробно вопрос о фазовой ошибке
рассмотрен в работе [Л. 49]. Из сказанного следует, что всегда,
когда это только возможно, следует избегать (Включения в схемы на-
бора цифровых интеграторов, работающих в режимах сум!матора
или нуль-ограна. Последнее обстоятельство совершенно игнорирует-
ся рядом авторов [Л. 1].
Умножение на постоянную. При работе в этом режиме,
часто называемом в литературе также режимом масштабного ин-
тегратора, в регистр Y цифрового интегратора помещается неко-
торое „начальное значение" у(0) = А, а входы Дг/i, ..,AyQ ни с чем
не соединяются. На вход Ах подаются сигналы AS (//), имеющие
смысл приращений некоторой величины ? (п) (рис. 6).
0—-
+> 0
Рис. б. Схема включении для работы
в режиме умножения на постоянную.
Предположим, что на некотором интервале изменения перемен-
ной п все сигналы А? (п) имеют положительный знак, а А = ^^^1г^ >
/72
где г/макс по-прежнему максимально возмол<ное положительное зна-
чение переменной у{п), а т — произвольное число, большее еди-
ницы. Тогда очевидно, что для образования каждого очередного
сигнала + Дг необходимо т раз передать число А из регистра У
в регистр R (см. выше). Непосредственно из того обстоятельства,
30

что каждая такая передача возможна только при поступлении от-
личного от нуля сигнала Л^, делаем заключение, что
b b b
V Д2 (/) = 4г У! (О = Ум (О, (31)
1=а
где а я b — границы рассматриваемого интервала изменения пере-
ме1пюи п. Соотношение (31) легко распространить на случаи отри-
цательных значений величин А и (/г).
В случае, если в регистре Y хранятся значения некоторой пере-
менной величины у{х) с достаточно хорошей степенью приближе-
ния, можно считать, что
Iz(x)^y{x)J^(x), (32)
где Az{x) и АS (х)—„средние" значения переменных Az{n) и А^ (/2)
на некотором интервале изменения п.
Описанный режим работы интегратора используется достаточно
часто (в ЦДА это единственная возможность осуществить умноже-
ние на постоянный коэффициент), однако при .работе ib этом режи-
ме также вносится запаздывание длительностью .миннмум в одну
итерацию.
Вычисление определенных интегралов. В случае, когда сиг-
налы Ayi (п) и Ах (п) имеют смысл приращений двух произвольных
переменных у(п) и х{п), можно считать, что цифровой интегратор
приближенно выполняет операцию вычисления определенного ин-
теграла
b b
^ydx^Yj^^^{l)Ax{l). (33)
а а
Интеграл в правой части выражения (33) можно определить как
в смысле Римана, так и в смысле Стильтьеса. Это зависит от кон-
кретного вида функции х(п\). ,В обоих случаях речь идет о прибли-
жении этото интеграла конечной суммой, причем для совершенно
произвольных переменных у и х ничего нельзя__сказать о качестве
такого приближения [вспомним, что величина Ах ограничена соот-
ношением масштабов (20)]. Этот последний режим, как правило, не
используется при нормальной работе ЦДА.
9. ОСНОВНЫЕ ХАРАКТЕРИСТИКИ ЦДА
В этом параграфе мы попытаемся подвести некоторые итоги
всему сказанному выше. Итак, ЦДА — это цифровая вычислитель-
ная машина, состоящая из отдельных решающих блоков — цифро-
вых интеграторов. Каждый цифровой интегратор выполняет вычисле-
ния в соответствии с жестко фиксированным алгоритмом, определяе-
мым разностными уравнениями (10) и (19). При решении каждой
31

шикретной задачи необходимое количество цифровых интеграторов
соединяетоя между собой по схеме набора. При составлении схемы
набора, как правило, исходят из предположения, что цифровые инте-
граторы точно (или во всяком случае достаточно точно) выполняют
операцию интегрирования некоторой функции у ло независимой ле-
ре.менной х. Именно это лредлоложелие и позволяет свести методи-
ку решения задач на ЦДА к хорошо разработанной методике
решения задач на обычных (аналоговых) дифференциальны,х анали-
3aTqpaix и обусловливает «простоту nporpaiM'MHpoBaHHH», выдвигае-
мую в качестве одного из основных преимуществ ЦДА.
Однако на самом деле подобное Л|редположение далеко не всег-
да оказывается уместны1М. Каждому, кто работал с ЦДА, Х)Орошс
известно, что после набора задачи и пробного ее решения сплошь
и рядом приходится изменять первоначально введенные значения
коэффициентов и начальных условий, другими словам1И, лроизводить
«подгонку» получаемых решений.
Не говоря уже о том, что подобный процесс «отладки програм-
мы» не всегда можно провести, так как он Т1ребует хотя бы мини-
мума первоначальных сведений о характере ожвдаемых решений, он
может оказаться значительно более трудоемким, чем само составле-
ние схемы набора. Конечно, огромную роль здесь играет наличие
опыта в решении однотипных задач. Однако даже при наличии та-
кого опыта далеко не всегда удается получить решение с желаемой
точностью. Причина этого состоит в следующем.
Схема набора, состоящая из некоторого количества цифровых
интеграторов, олисывается разностным уравнением, в общем случае
нелинейньвм, порядок которого равен количеств1у интеграторов, ра-
ботающих в режиме приближенного решения дифференциального
уравнения (см. выше), плюс суммарная задержка (измеренная в ко-
личестве итераций), вносимая масштабными интегратора1Ми, оум,ма-
торами, нуль-органами и т. п. Как травило, порядок этого разно-
стного уравнения оказьввается выше лорядка исходного дифференци-
ального уравнения. Будем считать, что в общем случае
r=P+q. (34)
где г —порядок раз1ностного уравнения;
р — порядок исходного дифференциального уравнения;
q — суммарная задержка.
Как известно из [Л. 10], раз1ностное уравнение порядка г при
соответствующих началыных условиях имеет г независимых реше-
ний, лредставляющих собой функции незашшмой переменной п. Вы-
делим среди них р решений, наиболее близко совпадающих с соот-
ветств1ующими решениями исходного дифференциальнош уравнения,
и будем называть их оонювньими, а остальные q решений — побочны-
ми. Часто считают, что побочные решения можно исключить соот-
ветствующ.им выбором начальных условий. Заметим, однако, что
сами начальные условия можно задавать только с некоторой опре-
деленной точностью (это числа с ограниченным количеством раз-
рядов), и, следовательно, полностью исключить побочные решения
таким методом в общем случае нельзя,
32

Основные решения разностного уравнения благодаря .наличию
двух видов погрешности, а именно методической погреш-ности и по-
греш1Ност1И округления, также лишь прибиженно совпадают с ре-
шениями исходного дифференциального уравнения. Из сказанного
следует, что при решении на ЦДА обыкновенных дифференциальных
уравнений имеют место следующие три основных вида погрешности:
1) иметодическая погрешность;
2) погрешность округления;
3) наличие побочных решений разностного уравнения, обуслов-
ленных задержками в вычислительной схеме.
Рассмотрим теперь более подробно каждую из этих погреш-
ностей.
Методическая погрешность. Выше было показано, что на каж-
дом шаге решения с помощью цифрового интегратора дифференци-
ального уравнения первого порядка очередное приращение функции
Аг(Хп)= { y[z{s), s]ds (35
вычисляется лишь приближенно по одной из формул (30а), (ЗОб')
или (30в|). Получаемая при этом погрешность называется методиче-
ской погрешностью (или погрешностью аппроксимации) одного шага.
Мы будем обозначать ее Jo п.
Для того чтобы оценить величину методической погреш,ности
одного шага, можно разложить функцию у [z(x), х] в ряд Тейлора
на интервале аГп+ь Хп, проинтегрировать этот ряд и сравнить полу-
ченное выражение с каждым из выражений (ЗОа), v-^б) и (30в|).
Раз1ность этих выражений снова будет представлять собой степен-
ной ряд. Используя для оценки величины Еп первый отличный от
нуля член этого ряда, получаем непосредственно:
Y i^xf max I z ''I для формулы (30a);
\Еп\<
(AxY max | z'' | для формулы (306);
5
Y^i^xy max |2'" I для формулы (ЗОв).
Заметим теперь, что погрешность Еп представляет собой абсо-
лютную погрешность в вычислении интеграла (35):
En=Az(Xn) —Az(n),
где Аг(Аг)—величина, вычисляемая в течение данной итерации по
формулам (30а1), (306) или (ЗОв).
Определяем теперь относительную методическую погрешность
одного шага как
"~ max I Z Г
3 А. в. Шилейко. 33

Эти оценки позволяют сделать вывод, что наибольшую величи-
ну методической погрешности одного шага дает формула (306).
Формула (306,) сложнее формулы (30а) в том смысле, что она тре-
бует дополнительного оборудования для умножения на a = V2. По-
этому применять формулу (306) для выполнения операций в циф-
ровом интеграторе coBepuiein-io нецелесообразно. Использование ее
в некоторых конструкциях ЦДЛ [Л. 1, 30] связано, видимо, с до-
садным ыедоразу.меиием. В далы1ейи1ем мы будем рассматривать
только формулы (30а) и (ЗОз).
Методическая погрешность одного шага накапливается от шага
к шагу в процессе вычислений. Вопрос о накоплении этих погреш-
ностей мы рассмотрим ниже.
Погрешность округления. Основная причина возникновения
в цифровых интеграторах ЦДА погрешностей округления состоит
в том, что в соответствии с методом РДМ вместо полной величины
l\z{n), вычисляемой по формулам (30а) или (306), по каналам
связи между решающими блоками передаются ее округленные зна-
чения Az{ii). Погрешность округления мы будем обозначать буквой
Гп. Она, очевидно, равна «остатку» округления
Часто высказывается мнение, что как методическая погрешность
одного шага, так и погрешность округления могут быть сделаны
в ЦДА сколь угодно малыми при соответствующем выборе величи-
ны интервала Ах или, что то же самое, масштаба Мх. Из выраже-
ния (37) видно, что относительная погрешность округления не зави-
сит от Ах и полностью определяется числом разрядов N. Следова-
тельно, в данном случае это мнение является явно ошибочным. По-
кажем теперь, что оно является ошибочным также и для случая
относительной методической погрешности одного шага.
Заметим прежде всего, что результатом работы данного цифро-
вого интегратора является последовательность чисел Дг(п), которые
34

для получения послело'вателыюст,и значений переменной г{п\} сум-
мируются в регистре Y другого 'Интеграто.ра
zin) = z(0)+yiz(i).
ы\
Обозначим индексами 1 величины, относящиеся к интегратору, по-
грешность которого исследуется, и индексами 2 — 'величины, относя-
щиеся к интегратору, в регистре Y которого накапливаются значе-
ния z{n). с учетом соотношения (20) получим:
Подставляя выражения (38) и (39) в формулы (36), получим
окончательно
для формулы (ЗОв).
Полагая для простоты Myi — My^ nG=l, получаем оконча-
тельно
для формулы (ЗОв).
Соотношения (40а) и (406) дают возможность сделать вывод,
что относительная методическая погрешность одного шага также
полностью определяется количеством разрядов в соответствую-
щих регистрах. Заметим также, что методическая погрешность фор-
мулы (ЗОв!) на три порядка меньше погрешности округления. Это
обстоятельство делает использование формулы (ЗОв) совершенно
нецелесообразным, поскольку обеспечиваемая этой формулой точ-
ность вычислений не может быть использована благодаря наличию
погрешности округления. Поэтому все дальнейшие рассуждения бу
дут проводиться только применительно к формуле (30ai).
3* 35

Как уже отмечалось выше, оба рассмотренных виаа погрешно-
сти накапливаются в процессе вычислений. Процесс накопления по-
грешности описывается разностным уравнением для погрешности.
Вьпвод этого уравнения будет рассмотрен в гл. 5. Здесь мы ограни-
чимся тем, что проиллюстрируем процесс накопления погрешности на
ряде простых примеров.
Рассмотрим простейшее дифференциальное уравнение вида
у' = ау. (41)
При начальном условии
^(0) = ^о.
При наборе на ЦДА уравнение (41)
преобразуется к виду
Ау{п) = аАху{п) (42)
или после введения масштабов
Рис. 7. Схема набора для а / \ hU ( \.
решения уравнения у'-^ау. аи[П) ^йПи[П);
где
h =
и{0)-
(43)
Уб .
Му J.\ly
Ь = аМ:с-
Соответствующая схема набора показана на рис. 7.
Если пренебречь задержкой, вносимой масштабным интеграто-
ром, то работа схемы, показанной на рис. 7, может быть описана
следующим разностным уравнением:
и(п+1) — (:\—Ьк)и(п)=0 (44)
Решение уравнения (44) при начальном условии u{Q)=Uq имеет вид
(45)
где
Заметим, прежде всего, что
lim u,(l+bh) ^ ^и,е^\
/г->0
(46
т. е. при стремлении величины шага h к нулю решение уравнения (44
совпадает с решением исходного уравнения (41).
36

Для определения погрешности при конечных h полол<им
(\+bhf^ =Л (47)
где
«=;^1п(1+6Л).
Раскладывая \n{\ + bh) в ряд Тейлора и ограничиваясь двумя пер-
выми членами, получим:
bh
a^l--^-. (48)
На основании выражений (47) и (48) имеем:
2
е
Таким образом, при конечном h решение уравнения (44) отли-
чается от истинного решения уравнения (41) множителем е . ко-
торый изменяется с изменением независимой переменной по экспо-
ненциальному закону.
Учтем теперь влияние задержки, вносимой масштабным интегра-
тором. При наличии этой задержки разностное уравнение, описы-
вающее работу схемы (рис. 7), приобретает вид
u{n + 2) — a(n+\) — bha(n) = 0. (49)
Заметим прежде всего, что это уравнение второго порядка. Одним
из наиболее серьезных последствий запаздываний, вводимых блока-
ми масштабных коэ'ффициентов ЦДА, является повышением порядка
разностного уравнения по сравнению с порядком исходного диффе-
ренциального уравнения.
Общее решение уравнения (49) имеет вид
а(%) = С,\х^ +€21^2 , (50)
где
— корни характеристического уравнения, а С\ и Сг — постоянные,
определяемые из начальных условий.
Исследование выражения (51) показывает, что при отрицатель
ном b и при условии
37

решение (50) будет иметь периодический характер, т. е. даже качест-
венно будет отличаться от решения исходного дифференциального
уравнения. Однако такой случай практически исключен, поскольку h
обычно выбирается достаточно малым
|5|<1.
В случае, когда
0<6Й
4 '
имеем:
Л,>1; -^<Х2<0.
С,
а)
Рис. 8. Решения уравнения у'=ау с учетом задержки в масштабном интегра-
торе.
±
При этом решение С,Х, ^ представляет собой функцию, монотонно
возрастающую с ростом независимой переменной а решение
€2^2 — знакопеременную функцию, убывающую с ростом ^ (рис. 8,а)
В случае, когда
—^<bh<0,
оба корня Xi и ^2 положительны и по абсолютной величине меньше
единицы. Решения Cik\^^ и Cz'^^J^ имеют вид, показанный на рис. 8,6,
Здесь, как и в предыдущем случае, решение С^2^ затухает с ро-
стом независимой переменной
Для определения постоянных Ci и Cg составим систему урав-
нений
С, + С2 = и(0).
(52)
38

Определяя Ам(0) из уравнения (49), получаем:
u[Q) — \2ii {Щ-ЬНа(—\)
.у-
+ bh
С2=
bhu (—1) —^2"(l)
-+bh
(53)
Поскольку при достаточно малых h решение С>^2 быстро затухает
с ростом 5, можно положить приближенно
(\-\,)u(Q)-bhu(-\) ,
и (£) ^ л ■ ^1
. + bh
Или, раскладывая + bh в ряд Тейлора и производя необхо-
димые преобразования, получим окончательно:
(1 —X,)u{0) + bhu (—1) /l+2bh\T
+ bh
(г
+ bhj
(54)
Здесь снова при стремлении h к нулю решение разностного уравне-
ния (49) стремится к решению исходного дифференциального урав-
нения (51), однако при конечных h вносится погрешность, которая
для каждого конкретного случая может быть определена из выра-
жения (54).
Рассмотрим теперь дифференциальное уравнение второго по-
рядка:
У(0) = Уо\ У'(^) = У'^ )
(55)
Схема соединений решающих блоков для этого уравнения пока-
зана на рис. 9. После введения масштабов и выполнения необхо-
димых преобразований получаем разностное уравнение
u(n+'d)'-2a(n + 2) + a(n+\) — bh''u{n)=^0.
(56)
Здесь снова благодаря наличию в схеме масштабного интегра-
тора порядок разностного уравнения оказался на единицу выше
порядка исходного дифференциального уравнения. Соответствующее
характеристическое уравнение имеет вид
(57)
39
?l(?l— 1)2 —6/i2 = 0.

Уравнение (57) может быть решено в общем виде. Однако для
проводимого здесь упрощенного анализа достаточно качественного
рассмотрения характера получаемых решений.
Положим
Л, = 1 + е.
Подстановка в уравнение (57) дает
ф >
Рис. 9. Схема набора для решения
уравнения у"^ау.
Поскольку обычно очень малая величина, можно положить при-
ближенно
z^hV~b для 6>0
t = bh'^—\ для 6 < 0.
Деля уравнение (57) на (А. —А,), получим:
Х2-(2-Л,)А+(Л,-1)а = 0,
откуда
(59)
(60)
На основании проведенного анализа можно сделать следующие
выводы:
1. Корню Xi соответствует решение уравнения (56) вида CJ^i^ ,
которое представляет собой чистую погрешность. Эта погрешность
быстро затухает с ростом ^ при 6 <^ О и возрастает с увеличением
S по экспоненциальному закону при 6 > 0.
2. Корням ^2 и Яз соответствуют решения уравнения (56), в пер-
вом приближении при достаточно малых h совпадающие с реше-
ниями исходного дифференциального уравнения. В частности, при
6<0 решение уравнения (56) может быть записано в виде
и (S)::::: е^^^^ (А sin со5 + 5 COS соЕ),
(61)
где а>;^ "|/6-а.6^ ; Л и
начальных условий.
40
постоянные, определяемые из

Частота этого решения при малых h очень близко совпадает
с частотой истинного решения, в то время как амплитуда отличается
от амплитуды истинного решения множителем увеличиваю-
щимся с ростом S по экспоненциальному закону.
Итак, мы видим, что точность решения задач на ЦДА опреде-
ляется в основном количеством разрядов в регистрах цифровых
интеграторов, причем получаемая при этом погрешность, а также
побочные решения соответствующего разностного уравнения могут
непрерывно расти в ходе решения и в отдельных случаях полностью
исказить получаемые результаты.
В большинстве случаев оказывается возможным значительно
повысить точность путем подбора начальных значений и некоторого
изменения схе.мы набора. Так, например, для образования функций
riZH
7
Рис. 10. Схема набора для образования функций sin сол:
и cos (ИХ.
sin сол: или cos сох можно воспользоваться схемой набора, показан-
ной на рис. 10. Здесь А и В — дополнительные интеграторы, рабо-
тающие в режиме умножения на постоянные а и 6.
Работа схемы, показанной на /рис. 10, описывается системой
разностных уравнений ^
у, (А2 + 1) = (\+Ь) у,{п) + аАху2 Щ
У2 (/^ + I) = У2 (п) + аАху, (/г). J
При условии, что
у г (0) = 0;
1 — о =.-COS
ahy^ (0) = sinp,
решением системы (63) будет функция
yi{x) = sinj^ х.
(62)
(63)
(64)
Из соотношений (63) и (64) получаем условия для выбора начальных
условий и коэффициентов:
У. (0) = 0;
Без учета задержек в масштабных интеграторах.
41

/а" (hxf .
1- 4 '
sin (О
b = — a^ (Axf.
Выбор этих значений и составляет процесс отладки программы.
Из рассмотренного примера видно, что даже в тако.м простом слу-
чае этот процесс может быть достаточно трудоемким. При всех ус-
ловиях минимальная погрешность, очевидно, не может быть меньше
погрешности округления одного шага, т. е. величины В-^. Именно
эту величину мы будем принимать в дальнейшем в качестве оценки
точности ЦДА, хорошо понимая при этом, что такая оценка являет-
ся оценкой снизу. Другими словами, реальная погрешность всегда
будет выше этой оценки.
Рассмотрим теперь вопрос о быстродействии ЦДА. Пусть не-
которая переменная величина z (п) представляет собой одно из
приближенных решений исходного дифференциального уравнения.
Соответствующим точным решением этого уравнения является функ-
ция Z (х) (х ==Хо-\- пАх). Если мы хотим получить приближенное
решение с минимальной относительной погрешностью е^, мы должны
положить B^^^Bz. Будем считать, что B~'^=zi.
Введем теперь в рассмотрение приближенное соотношение
max
dz (X)
dx
Отсюда ясно, что максимум модуля производной- в полученном ре-
шении определяется количеством разрядов N и масштабом My.
Предположим теперь, что результатом вычислений на ЦДА яв-
ляется гармоническая функция
z(x) = А sin (ИХ.
Максимум производной этой функции можно получить из соотно-
шения
max
dz (X)
dx
или окончательно, с учетом соотношения (;15),
ЛсОмакс =В У СОмакс =В ^ ;
СОма
;= 1.
(65)
Предположим, теперь, что длительность одной итерации в данном
ЦДА равна Т. За время этой итерации независимая переменная х
увеличивается на величину Ал:. С учетом этих соображений получим:
/макс — «макс 27C^f/ ~2^T
(66)
42

Соотношение (бб) дает оценку Для максимальной частоты синусо-
идального сигнала, который может быть воспроизведен с заданной
точностью на ЦДА, в котором длительность одной итерации состав-
ляет Т сек. Таким образом, мы видим, что точность и быстродейст-
вие ЦДА оказываются связанными между собой соотношением (66).
На первый взгляд подобный прием для оценки быстродействия
кажется несколько искусственным, однако он находится в полном
соответствии с общепринятым методом оценки динамических свойств
непрерывных систем. Действительно, рассматривая синусоидальный
сигнал неременной частоты, мы как бы снимаем частотную характе-
ристику системы, моделируемой схемой набора на ЦДА. Более по-
дробно вопрос о частотных характеристиках цифровых моделей будет
■рассмотрен в гл. 5 этой книги. Там же будут выведены более стро-
гие оценки. Вычисления по формуле (66) показывают, что для
ЦДА, длительность итерации в котором составляет 20 мсек (это
стандартная длительность итерации для ЦДА с запоминающим
устройством на магнитном барабане), при точности в пять десятич-
ных знаков граничная частота синусоидального сигнала составляет
всего 0,8-10-4 гц. Для ЦДА, имеющего длительность одной итера-
ции в 10 мксек ,(это наименьшая известная в настоящее время дли-
тельность итерации), граничная частота синусоидального сигнала со-
ставляет 0,16 гц. Из проведенных выше выкладок видно, что быстро-
действие ЦДА целиком и полностью определяется принятым мето-
дом квантования входных и выходных сигналов, т. е. методом РДМ.
При этом наблюдается явное несоответствие между динамическими
свойствами использованных разностных уравнений, описывающих
работу цифрового интегратора, и динамическими свойствами спосо-
ба передачи данных. Вычисления по формулам, выведенным в гл. 5,
показывают, что при прочих равных условиях граничная частота
синусоидального сигнала, обеспечиваемая разностным уравнением,
для первого примера [длительность итерации 20 мсек, разностное
уравнение (ЗОа)] составляет 6,2-10-2 гц, а граничная частота сину-
соидального сигнала во втором примере [длительность итерации
10 мксек, разностное уравнение (ЗОв)] составляет 730 гц.
Это несоответствие представляет собой один из самых главных
недостатков ЦДА.
Теперь можно сделать окончательный вывод. Среди прочих ма-
шин, предназначенных для решения обыкновенных дифференциаль-
ных уравнений, ЦДА безусловно являются наиболее простыми, по-
скольку они состоят из решающих блоков одинаковой структуры,
а в случае ЦДА чисто последовательной структуры физически реа-
лизуется только один такой блок, что позволяет использовать мини-
мум элементов. Если применять ЦДА для решения одной и той же
задачи или достаточно узкого класса однотипных задач (что обыч-
но имеет место при включении вычислительной машины в состав
системы автоматического управления), то тем самым можно избе-
жать трудностей, связанных с отладкой программы. Однако при этом
остается основное ограничение, связанное с чрезвычайно малым
быстродействием. Таким образом, ЦДА можно рекомендовать для
использования в составе систем автоматического управления весьма
медленно протекающими процессами (некоторые химические про-
цессы, тепловые процессы и т. п.). С другой стороны, совершенно
нецелесообразно использовать ЦДА в вычислительных центрах и
в лабораториях для решения уравнений. Действительно, сложность
43

программирования ЦДА с учетом последующей отладки программы
имеет, во всяком случае, тот же порядок, что и сложность програм-
мирования современных УЦВМ, предусматривающих возможность
использования любой из современных систем автоматического про-
граммирования. С другой стороны, наличие .жестких алгоритмов ра-
боты решающих блоков полностью исключает возможность исполь-
зования существующих методов контроля точности (работа с пере-
менным шагом, повторное решение другим численным методом
и т. д.), поэтому в большинстве случаев при работе на ЦДА оказы-
вается невозможным даже приближенно оценить точность получае-
мых решений. Это обстоятельство достаточно ярко иллюстрируется
рассмотренными выше простыми примерами. Несколько более по-
дробно этот вопрос будет рассмотрен в гл. 4 данной книги. Наконец,
использовать ЦДА в составе систем управления быстро протекаю-
щими процессами практически невозможно ввиду ограниченного их
быстродействия.
Таки.м образом, единственной возможной областью применения
ЦДА является управление медленно протекающими процессами, где
может быть реализовано единственное их преимущество, а именно
их простота. Однако и здесь необходимо проявлять известную осто-
рожность. Во всяком случае перед включением ЦДА в контур управ-
ления каждый раз необходимо убеждаться, что набранная на нем
задача решается с необходимой точностью. Это можно сделать, на-
пример, проводя контрольное решение другим численным методом
на универсальной цифровой вычислительной машине и сравнивая
полученные результаты.
Как уже от.мечалось в предисловии, ЦДА сыграли значительную
роль в развитии техники цифрового моделирования. Однако разра-
ботку машин этого типа в настоящее время вряд ли можно счи-
тать целесообразной,
ГЛАВА ТРЕТЬЯ
ЦИФРОВЫЕ МОДЕЛИ С ПЕРЕМЕННЫМИ
ПРИРАЩЕНИЯМИ
10. ИСХОДНЫЕ СООБРАЖЕНИЯ
В этой главе мы рассмотрим конструкции двух машин, близких
к ЦДА по основной структуре, но имеющих принципиальные отличия
как в части реализуемых алгорит.мов, так и в части способов коди-
рования выходных сигналов каждого функционального блока. Об-
щим для этих конструкций является то, что обе они появились
в результате попыток преодолеть основные недостатки ЦДА, связан-
ные с использованием для кодирования выходных сигналов метода
РДМ.
Первая из них —машина GEVIC (General Electric Variable
Increment Computer), разработана фирмой General Electric (США)
в качестве бортовой машины и предназначена для решения навига-
ционных задач, задач управления огнем, бомбометания и т. п.
[Л. 32]. Реализуемые в машине алгоритмы предусматривают в ос-
новном решение алгебраических и тригонометрических уравнений.
Машина может быть использована также и для решения дифферен-
44

Ц'Иаль'Ных уравнений, однако получаемая при этом точность в об-
щем случае должна быть даже хуже, чем у обычных ЦДА. Поэто-
му машину GEVIC следует отнести к классу цифровых моделей,
предназначенных для моделирования объектов, работа которых
описывается алгебраическими и тригонометрическими уравнениями.
Конструкция второй машины была предложена в 1958 г. Л. А. Ко-
жарским [Л. 26]. Это машина для решения обыкновенных диф-
ференциальных уравнений одним из численных методов интерполя-
ционного типа, так называемым методом трапеций. Значения
переменных на каждом шаге вычисляются в ней в результате несколь-
ких последовательных итеративных циклов. Насколько известно ав-
тору, машина никогда не была построена, однако она представ-
ляет несомненный интерес прежде всего как одна из попыток
получения рациональной структуры цифровой модели, предназначен-
ной для решения обыкновенных дифференциальных уравнений.
11. МАШИНА GEVIC
Основной отличительной особенностью машины GEVIC является
использование в ней способа передачи выходных сигналов функцио-
нальных блоков, получившего название способа «переменных при-
ращений» или способа «передачи но частям». Основная идея этого
способа состоит в следующем. Пусть результатом работы некоторого
функционального блока является неременная величина t/(n), прини-
мающая ряд последовательных значений г/(1), у(2), ... , у(\п) для
значений независимой переменной 1, 2, ... , п. Предположим, что
в некоторый момент времени вычисленное значение этой величины
равно y(i). Это значение представлено Л^-разрядным двоичным чис-
лом, хранящи.мся в специально.м регистре, называемом регистром Р.
Возьмем старший двоичный разряд числа y{i) и будем рассматри-
вать его в качестве выходного сигнала данного функционального
блока Ау (/—1), соответствующего значению i независимой пере-
менной. Очевидно, что этот сигнал может рассматриваться как при-
ближенное значение числа y(i) с точностью до его младших двоич-
ных разрядов. Следующему значению независимой переменной со-
ответствует значение y(i+l) рассматриваемой функции.
Образуем величину
p(i+\) = p (О + Ау (1) -.Ay{i-ll
где
^y{i) = y{l+^)r-y(^)^
т. е. вычтем из содержимого регистра Р_величину Ау (i) и возьмем
в качестве нового выходного сигнала Af/(/+l) старший разряд
этой разности. Продолжая описанный процесс, получим последо-
вательность выходных сигналов Ау {i), Ау (I + \), ..., Ау{п).При
этом будет справедливо очевидное соотношение
п
Ау(у) + р{п) = у{п), (67)
v = l
где р{п) — ^остаток*, хранящийся в регистре Р. При наличии опре-
деленных ограничений на возрастание абсолютной величины функ-
ции у[п] величина остатка р(п) будет стремиться к нулю, а сумма
45

сигрталов ^ Ay (v) будет стремиться к точному значению функц11й
1
у(п). Последнее обстоятельство иллюстрируется кривыми, по-
казанными на рис. 11. Учитывая сказанное, величины At/(v) можно
считать „приращениями" функции у(п). Заметим также, что по-
скольку сигнал Af/(v) всегда представляет собой старший двоичный
разряд числа /7 (v), его численное значение точно равно 2^ т. е.
целой степени 2. Поэтому, для того чтобы передать такой сигнал
Рис. и. к расчету последова-
тельности выходных сигналов
машины GEVIC.
1-^1 1*2 i+3 i+Ч 14-5
по каналу связи, достаточно передать только его знак и значение
показателя степени k. Так, например, сигнал Af/(v) =-|-2^ будет
передан в форме 0,011, а сигнал Ai/(v) = — 2^ будет передан в фор-
ме 1 101. Это позволяет значительно уменьшить загрузку каналов
связи по сравнению с тем случаем, когда значения у{п) передают-
ся полностью.
Заметим сразу же, что основной характерной особенностью
рассмотренного способа передачи является то обстоятельство, что
п
сумма приран^ений 'Va^(v) сходится к правильному значению фуик-
1
ции у{п) к концу вычислительного процесса, в то время как на
промежуточных этапах этого процесса погрешность смол<ет дости-
гать в обп;см случае у{^). Для большого класса задач и, в ча-
cthoctii, перечисленных выше задач 'навигации, управления огнем
и т. п. подобное положение является .вполне допустимьгм, поскольку
интерес здесь представляет именно конечная точка решения (напри-
мер, вывод на цель). Однако в общем случае подобный способ пере-
дачи может привести к появлению погрешности округления, во
много раз прсвьш1ающей истинные значения решений. Подробнее
этот вопрос будет рассмотрен в следующих главах.
Полная структурная схема машины GEViIC показана на рис. 12.
Машина coctoiit из четырех регистров К, 5 и Р, представляю-
щих собой, как и в случае ЦДА последовататыного типа, дорожки
мапштного барабана. Емкость каждой такой дорожки равна MN,
где Л' — количество разрядов в числах, над которыми выполняется
операция, а М — полное количество различных функциональных
блоков. В регистрах U w V хранятся (последовательные значения
4G

чисел Ui(n) и Vi{n) Эти регистры относятся к циркуляционному
типу и замкнуты через сумматоры 2^ и 2у. На вторые входы этих
сумматоров поступают приращения Aui\(n) и AVi'{n). Таки!м образом,
после завершения каждого полного оборота барабана все значения
Ui{n—1)и Vi(n—1)заменяются значениями Ui(n) = Ui{n—\)+Aui{n)
и Vi(n) = Vi(n—\) + Avi(n). В регистре S хранятся постоянные чис-
ла S. Наконец регистр Р замкнут через сумматоры и :^з. В этом
регистре хранятся последовательные значения чисел р{п). Выходы
регистров и, V и S связаны с тремя аналогичными схемами умно-
жения иа приращение. Каждая такая схема состоит из элемента,
образующего дополнение, и элемента переменной задержки. Вели-
чина задержки в этом элементе устанавливается каждьп"! раз в со-
ответствии с весом данного сигнала приращения. Таким образом
умножение, например, числа Ui(n) на приращение Ati{n) =—2^ свя-
зано с изменением знака числа Ui{n) на обратный, т. е. образова-
нием дополнительного кода числа Ui(n) и сдвигом этого числа вле-
во иа k двоичных разрядов.
Выходные сигналы схем умиожения поступают на три после-
довательно включенных одноразрядных двоичных сумматора 2i, '^2
и 2з. Выходной сигнал су.м1матора Z3 представляет собой число
р{(п), поступающее в соответствующий регистр Р. Непосредственно
из схе.мы видно, что в течение каждого вычислительного цикла об-
разуется новое значение каждого числа Рг(п), равное
Рг н = Рг{П—1)— sKpi (а2 — I) + Ui (п) At i (п) +
+ Vi{n-\)/^i{n)-^A,.{n). (68)
Выходной сигнал сумматора 2з поступает также на так назы-
ваемый блок выбора величины приращения {БП на рис. 12), задачей
которого является формирование приращения Api{n), т. е. выделе-
ние старшего разряда чиата в соответствии со сказанным выше.
Блок БП состоит из двух двоичных счетчиков, системы триггеров
и переключательных элементов. Рассмотрим подробнее работу этого
блока.
|Каждое очередное число Рг {п) поступает на вход блока БП по-
следовательно, начиная с младшего разряда. Перед прохождением
каждого очередного разряда содержимое первого счетчика (назовем
его С, рис. 13) увеличивается на единицу. Таким образом, к мо-
менту прохождения первого, самого младшего разряда числа Р\{п),
в счетчике оказывается заннсанным число 1, к моменту прохожде-
ния второго разряда—число 2 и т. д. Если значение очередного
разряда числа Pi{n) равно 1, то число, накопленное к этому мо-
менту в счетчике, передается в промежуточный регистр Q, состоя-
щий из триггеров Т\—Ti. После этого процесс продолжается. При
поступлении некоторого следующего разряда того же числа Pi{n),
значение которого также равно единице, предыдущее содержимое
регистра Q стирается и на его место записывается новое число,
накопленное к этому 1Моменту в счетчике С. iK концу описанного
процесса в регистре Q окажется записанным число, равное весу
* Здесь и в дальнейшем при описании машины GEVIC индекс i будет
обозначать номер вычислительного блока, к которому относятся числа и^(а2),
v^{n) и т. д.
47

старшего разряда числа pi{n) при условии, что этот вес не превы-
шает максимальной емкости счетчика, т. е. 2К Разряды числа pi{n)
с весом, большим 2^ просто не рассматриваются. Значение / раз-
лично у различных модификаций машины GBViIC (см. ,ниже).
Последний разряд числа Рг{п), поступающий па вход блока БП,
представляет собой разряд знака. Если значение этого разряда рав-
но нулю [число Рг{п) положитсльно], ТО В регистре Q оказывается
Рис. 12. Полная схема машины GEVIC.
помещенным правильное значение приращения. Если значение раз-
ряда знака равно единице, то число, находящееся в регистре Q, не
соответствует истинной величине нриращения. Для_ образования от-
рицательных приращений служит второй регистр Q (рис. 13), кото-
рый, так же как и регистр Q, связан со счетчиком С системой вен-
тилей. Однако эти вентили управляются значениями разрядов числа
Pi(n), представляющего собой дополнение к числу Pi{n). Таким
образом, в регистре Q образуется правильное значение приращений
Api(n) в том случае, если число Рг(п) отрицательно. В зависимости
от знака числа р{(п) открывается один из вентилей В-\- или В—
и правильное значение приращения Арг{п) выдается на выход бло-
ка БП,
(Кроме описанных операций, блок БП осуществляет также деле-
ние последовательности выходных сигналов Api'{n) на величину
постоянного коэффициента s (ключ Ki на рис. 12 в верхнем поло-
жении) либо на величину числа vi(n) (ключ Ki в нижнем поло-
жении). Процесс выполнения операции деления состоит в следую-
щем. Числа 5 или Vi{n), поступающие на вход блока БП, также
аппроксимируются одним старшим разрядом, причем в данном слу-
чае аппроксимация производится «'с избытком». Например, число 3
48

4 А. В. Шилейко.

аппроксимируется 'ближайшим значением 2^, т. е. 4. Элементы, вы-
полняющие эту аппроксимацию, для простоты не показаны на
рис. 13.
В момент поступления сигнала, представляющего собой резуль-
тат такой аппроксимации числа 5 (или Vi{n)), счетчик С сбрасы-
вается на нуль и затем начинает считать сначала. Например, если
число 5 = 3, оно аипроксимируется значением As = 4 и последова-
тельность чисел в счетчике имеет следующий вид: 0,1,0,1,2,3,...
Рассмотрим подробнее процесс выполнения операции деления .на
следующем примере. Будем предполагать, что в исходном состоя-
нии содержимое регистров U и V равно нулю, содержимое регистра
.S равно 3, а содержимое регистра Р=19, т. е. 00 ... 010011 в дво-
ичной системе счисления. (Предположим также, что выход Ар(п)
(см. рис. 12) соединен со входом Ар{п—1). В первом цикле на вход
блока БП поступают числа s и рг(п). Если бы число 5 было равно
единице, то соответствующее /приращение, образованное на выходе
блока БП, равнялось (бы 2\ т. е. выходной сигнал был бы равен
О 100 (^=4). Однако в результате сброса счетчика сигнал, образо-
ванный на выходе блока БП, будет иметь вид 0 010 (k=2). В сле-
дующем цикле этот сигнал поступит на вход Ар(п^1)_, и в резуль-
тате этого на 1вх0д сумматора 2з поступит число sApi{n), равное
00 ... 01100 {sApi(n) =А' 3=\2. Это число вычтется сумматором 2з
из числа рг{п), и в результате в регистр Р и на вход блока БП
поступит число pi(n+l)=00 ... 001М [р^ (a2-|-1) =49-^12=7]. В ре-
зультате уже описанного процесса новый сигнал Лрг(^+'1) будет
равен Арг(^ + 1) =0011. В данном случае при выборе величины
приращения аппроксимация при образовании Api{v) выполнена
«с избытком», т. е. вместо р{{п4-1) =7 взято рг(п-\-\) =8. (Выбор спо-
соба аппроксимации осуществляется специальными элементами бло-
ка БП, также не показанными на рис. 13.
В следующем {п + 2)-м цикле на вход сумматора S3 поступает
число sApi{n -j- \) = о, а на вход блока БП поступает число
Pi (п + 2)=^ О 0.. .00100 (/2 + 2) = 7 — 3 = 4]. Соответствующее
значение А/?, (я + 2) = 0,001 (^ = 1). Таким образом, в результате
трех вычислительных циклов мы получили значение частного,
равное
2
Y^Api (n + v) = 4+l + l=6.
v=0
Это значение соответствует величине ^^/з с точностью до младшего
разряда.
Из сказанного следует, что вычислительный процесс, выполняе-
мый в схеме, показанной на рис. 12, на каждом шаге, может быть
описан соотношением
- Рг{П—\)- slpi (n—\) + Ui (/2) Iti(n) +
Api (n) + Гг(п) = --. .
+ Vi(n — \)wi(n) +Kxi{n)
^ ... (69)
50

Это основной алгоритм, реализуемый во всех функщтопальных бло-
ках машины GEVIiC. Формуле (69) соответствует условное обозна-
чение функционального блока машины, показанное на рис. 14,а.
В соответствии с формулой (70) и рнс. 12 такой блок имеет пять
входных каналов Аи, Av, At, Aw и Ах и один выходной канал Ар(п),
нормально соединенный с входом Ар(п—\) (см. рис. 12). В зависи-
мости от различных способов соединения ме>аду собой этих каналов
получаются различные частные алгоритмы функциональных блоков.
Например, алгоритм сложения и вычитания имеет вид
Lpi (/?) + ri (/?)
р.(п-\)- sApi (п-\) + u,Mi (п) +
S
+ v,Iwi{n)+AXi(n)
(70)
Здесь величины Uq и Vq представляют собой постоянные весовые
коэффициенты. Подобный функциональный блок выполняет сумми-
рование_ трех_ переменных, заданных своими приращениями
Ati{n)y Awi{n), Axi(n), причем переменные ti и Wi суммируются
с весовыми коэффициентами. Результат суммирования делится на
масштабный коэффициент s. Соответствующая схема соединений
приведена на рис. 14,6.
Алгоритм умножения имеет вид
- Pi (/2 — 1) — sApi {n — \) + Ui {n)Avi {n)+
^Pi(n) + ri{n) = --
+ Vi(n-l)AUi{n)
(71)
Здесь используется формула умножения
Д (uv) (п) = и (п) Av (п) + у (/2 — \)Аи (п).
Результат умножения делится на масштабный коэффициент s. Соот-
ветствуюишя схема соединений показана на рис. 14,е.
Алгоритм деления имеет вид
— Vi(n — \ )AUi(n ^ \ ) — Ui (п — \ ) Av i (п)
(72)
Здесь на каждом шаге из остатка р{{п—А\) 1вычитается очередное
значение произведения UiVi и прибавляется очередное значение де-
лимого. Операция выполняется в соответствии с процессом, рас-
смотренным выше. Заметим здесь, что сходимость процесса деления
обеспечивается тем обстоятельством, что при аппроксимации дели-
теля каждый раз берутся значения «с избытком». Это соответствует
следящей системе с избыточным демпфированием. Соответствующая
схема соединений показана на рис. 14,г.
4* 51

Алгоритм извлечения квадратного корня имеет вид
Pi{n-\)--[Ui{n-\)+Ui{n-2)]Kui(n-\)+sAni{n)
Ui (п—\)
(73)
Здесь также вводится ib действие сходящийся итеративный процесс.
Соответствующая схема соединений показана на рис. 14,^.
Наконец формула ириближениого дифференцирования имеет
вид
dt
^ti (n)
(74)
Соответствующая схема соединений показана иа рис. 14,в.
Кроме перечисленных, может вьшолияться ряд более сложных
операций, также представляющих собой различные модификации ос-
COSX.
Stnx
Рис. 15. Схема набора для образования функций sin х и cos х.
новиого алгоритма. Для образования функции sin исполь-
зуется метод разложения этой функции в степенной ряд, причем
для образования последовательных степеней последовательно приме-
няется алгоритм умножения. Схема соединения функциональных бло-
ков для образования функций sin л; и cos л: показана на рис. 15.
Функция sin X образуется по формуле
sin х = С^х + С,х* + С,х' + С^х\
а функция cos X —по формуле
(75)
cos х = У\— sin'^x.
(76)
В качестве постоянных коэффициентов используются значения S».
При вычислении этих коэффициентов с помощью полиномов Че-
бышева обеспечивается точность порядка 10-^ в диапазоне изме-
нения X — -2<,х<,-2'
Машина GEVIC выполнена на магнитных элементах, изготов-
ленных на основе миниатюрных сердечников из пермаллоевой ленты.
Каждый элемент состоит из одного сердечника и двух диодов. Эле-
менты питаются от источника переменного тока высокой частоты.
52

Всего в машине GEVIC используются .600 таких элементов и запо-
минающее устройство иа магнитном барабане. Основная тактовая
частота равна 157,5 кгц. Этому соответствует возможность выполне-
ния 25 полных вычислительных циклов в секунду в 300 функцио-
нальных блоках. ^Полная длина каждого слова составляет 21 дво-
ичный разряд (19 разрядов мантиссы, 1 разряд знака и 1 разряд
для контроля но паритету). В одной из модификаций (машины, опи-
санной в работе [Л. 32], приращения могут принимать только сле-
дующие четыре |фиксированных значения О, +'1, -f 4 и -1-32. Для со-
единений между различными блоками используются адресные дорож-
ки м.агнитного барабана, так (же как и в случае обычных ЦДА.
Весьма знаменательным является то обстоятельство, что при раз-
работке машины GEVLC весьма широко использовался метод моде-
лирования отдельных ее узлов и всей машины в целом на универ-
сальной цифровой вычислительной машине. Для этой цели были спе-
циально разработаны две системы автоматического программирова-
ния, одна из которых, GiLOP, предусматривала автоматическую про-
верку правильности схем различных узлов и машины в целом, авто-
рая, GASP, предусматривала автоматическую проверку алгоритмов
машины.
12. МАШИНА ДЛЯ РЕШЕНИЯ ДИФФЕРЕНЦИАЛЬНЫХ
УРАВНЕНИЙ ИНТЕРПОЛЯЦИОННЫМ МЕТОДОМ
Схема машины, предложенная Л. А. Кожарским [Л. 26], пред-
усматривает решение обыкновенных дифференциальных уравнений
с помощью интерполяционной формулы численного интегрирования
вида
Ахц^,){п), (77
где Az(i)(/2)— приращение интеграла на выходе данного функцио-
нального блока на п-ы шаге интегрирования в j-й
итерации;
у{п — 1) — значение подынтегральной функции, вычисленное
в течение (п — 1)-го шага;
Af/(i_i)(/2) — приращение подынтегральной функции, вычислен-
ное в течение (i —1)-й итерации п-го шага ин-
тегрирования;
A^(i-i) (/2) —приращение независимой" переменной в течение
(/ —1)-й итерации п-го шага интегрирования.
Блок-схема машины показана на рис. 16.
Машина состоит из семи групп регистров, представляющих со-
бой, как и в предыдущем случае, дорожки :магнитного барабана.
Емкость каждой такой группы составляет MN двоичных единиц, где
М — общее количество функциональных блоков-интеграторов, а N—
количество разрядов в числах, над которыми выполняется операция.
Числа представляются в естественной двоичной системе счисления.
В региетре Y хранятся значения подынтегральной функции у (п),
В регистре АУ хранятся значения приращений Ау(п). В регистре К
хранятся значения вычисленных сумм А^/(л — 1)-f(az). Ре-
53

гистр АХ cлyжиt ДЛИ хранения значений Ах{п), Регистр Р слу-
жит для распределения вычисленных значений Аг мелсду ин-
теграторами. Эти значения по мере их вычисления поступают
в группу регистров AZ. На.конед, в регистре А хранятся специальные
коды пли адреса, ib соответствии с которыми осуществляется рас-
пределение вычисленных значений Az.
Все аперации, предусматриваемые формулой (78), выполняются
с помощью единственного одноразрядного двоичного сумматора 2.
Наибольший интерес здесь представляет способ 'Выполнения опера-
ции умножени.я на величину Ах, которая в данном случае представ-
ляет собой нолноразрядное двоичное число. В течение каждого
к
д
блон
\управлениА
Рис. 16. Схема машины для решения дифференциальных уравнений
интерполяционным методом.
оборота барабана (Производится умножение всех чисел, хранящихся
в группе регистров /С, на один из разрядов, соответствующих чисел
Ал:. К началу следующего оборота барабана числа Ал: и уже образо-
ванные частичные произведения сдвигаются на один разряд и снова
выполняется /умножение на один из разрядов Ах. Таким образом,,
для выполнения операции умножения на Ал: во всех интеграторах
требуется всего оборотов барабана.
(Е'оли вычисления данной группы значений Аг могут быть выпол-
нены с требуемой точностью за / итераций, то полная длительность
одного шага вычислений во всех интеграторах составит N1 оборотов
барабана. ,В остальном процесс решения дифференциальных уравне-
ний аналогичен соответствующему процессу в обычных ЦДА.
Несмотря на значительное увеличение (в N1 раз) длительности
одного шага вычислений по сравнению с обычным ЦДА, описанная
машина дает существенный выигрыш как в точности, так и в бы-
стродействии. На рис. 17 показаны кривые зависимости ;10лучаемой
точности в случае решения дифференциального уравнения у—у'=^0
от полного времени решения, измеренного в количестве оборотов ба-
рабана. Кривая / соответствует обычному ЦДА, а кривая 2—описы-
ваемой машине.
54

Подобный выигрыш объясняется двумя факторами. С одной
стороны, в описываемой машине увеличена пропускная способность
каналов связи между отдельными интеграторами (по этим каналам
передаются полноразрядные числа), а во-вторых, используется зна-
чительно более точная интерполяционная формула численного инте-
грирования. Более подробно этот вопрос будет рассмотрен в сле-
дующих главах.
Общим для двух описанных здесь машин является иапользо-
вание интерполяционных или итеративных методов вычислений. Дей-
ствительно, если при рассмотрении работы машины GEVLC условно
/.V
hO
0,8
0,6
0,4
о,г
о
N
\
S
s
\
ш т то 2560320038W тот57606т
Скорость 6ара6ана,од/мия
Рис. 17. Зависимость точности решения
от скорости барабана.
считать, ЧТО в течение некоторого фиксированного интервала вре-
мени числа, над которыми выполняются операции, остаются постоян-
ны.ми, то описанные выше алгоритмы превращаются в обычные ите-
ративные формулы выполнения соответствующих операций. Выпол-
нение такими методами алгебраических и тригонометрических опера-
ций связано с рядом преимуществ и поэтому представляет большой
интерес. С другой стороны, как это будет показано в последующих
главах, операцию интегрирования лучше выполнять по формулам
так называемого экстраполяционного типа.
ГЛАВА ЧЕТВЕРТАЯ
АЛГОРИТМЫ ЦИФРОВЫХ МОДЕЛЕЙ
13. ИСХОДНЫЕ СООБРАЖЕНИЯ
В этой главе мы рассмотрим несколько способов построения
алгоритмов, или программ для цифровых моделей. Как было отме-
чено в гл. 1, всякий такой алгоритм вне зависимости от способа
его получения должен представлять собой систему разностных
уравнений или рекуррентных формул. Одна из отличительных осо-
бенностей цифровых моделей состоит в том, что, если не считать
так называемые погрешности округлени-я, которые всегда можно
55

сделать асоль угодно малыми, разностное уравнение, описывающее
работу модели, решается точно, т. е. сама аппаратура не вносит ка-
ких-либо дополнительных погрешностей. Это значит, что после на-
хождения непротиворечивого алгоритма, отвечающего поставленным
требованиям, задача построения цифровой модели в основном может
считаться решенной. Рассмотренные свойства являются одним из
значительных преимуществ цифровых (моделей по сравнению с ана-
логовыми и аналого-цифровыми, 'у которых аппаратурная погреш-
ность играет существенную роль и притом не всегда может быть
учтена.
С другой стороны, выбор алгоритма в большой степени опре-
деляет все основные хара(ктеристики цифровой .модели, включая та-
кие, как точность, быстродействие, сложность, удобство эксплуата-
ции и т. п. Поэтом(у оказывается возможной постановка вопроса
о «качестве» того или иного алгоритма. Подробнее этот вопрос бу-
дет рассмотрен в следующей главе.
Задача выбора алгоритма или (программы для цифровой модели
сходна с аналогичной задачей выбора численного метода и построе-
ния программы для универсальной (цифровой вычислительной ма-
шины. Однако здесь имеется своя специфика. Цифровая модель
является специализированной цифровой вычислительной машиной,
поэтому, как уже отмечалось, ее характеристики в большой степени
будут зависеть от конкретного вида реализуемого в ней алгоритма.
При программировании для универсальной цифровой вычислитель-
ной машины, обладающей фиксированной емкостью запоминающего
устройства, количество запоминающих ячеек, нужное для решения
той или иной задачи данным способОхМ, практически не имеет ни-
какого значения, если только оно не превышает имеющейся емкости.
В случае же цифровой модели это количество запоминающих .ячеек
является одним из основных факторов, определяющих сложность
машины. Поэтом|у часто для цифровых -моделей выбирают численные
методы, обладающие меньшей точностью, но требующие -меньшего
количества аппаратуры.
>В данной работе мы ограничимся рассмотрением нескольких
способов построения алгоритмов для моделирования объектов, ра-
бота которых описывается обыкновенными дифференциальными
уравнениями. Задача построения такого алгоритма может быть сфор-
мулирована следующим образом. Найти разностное уравнение, ре-
шение которого при данных начальных или граничных условиях
совпадало бы в узловых точках с решением исходного дифферен-
циального уравнения. Разностное уравнение, отвечающее этому тре-
бованию, называется парным к исходному дифференциальному.
В теории разностных уравнений [Л. 20] доказывается существование
парных уравнений для достаточно широкого класса дифференциаль-
ных уравнений. Однако практичеокие методы получения парных
уравнений разработаны еще очень слабо. Кроме того, в большом
ряде случаев парное разностное уравнение, даже если оно найдено,
оказывается чрезвычайно неудобным для реализации. Поэтому
обычно вместо парного разностного уравнения используют какое-ли-
бо другое разностное уравнение, решение которого приближенно
соответ^ствует решению исходного дифференциального уравнения.
Методы численного решения обыкновенных дифференциальных
уравнений составляют самостоятельный раздел математики и описа-
56

ны в большом количестве источников [Л. 3, 28, 34, 8]. Количество
практически используемых методов также достаточно велико, одна-
ко подавляющее их большинство можно свести к небольшому ко-
личеств|у основных типов. Предлагаемая ниже классификация не
претендует иа исчерпывающую полноту. Однако, по нашему мне-
нию, она может облегчить понимание основных идей, используемых
при построении каждого опдельного метода.
По способу преобразования исходного диф-
ференциального уравнения 1мы будем различать разно-
стные и квадратурные численные методы. Сущность разностных ме-
тодов состоит в том, что производные, входящие в состав исходного
дифференциального уравнения, приближенно заменяются конечными
разностя'Ми, или, точнее, линейными комбина^циями значений зави-
симой переменной. |Подо«бная замена и приводит к получению ис-
комого разностного уравнения или, как говорят, разностной схемы.
При построении квадратурных методов исходное дифференциальное
уравнение сначала преобразуется к системе интегральных уравнений
первого порядка. Один из примеров такого преобразования мы рас-
сматривали в гл. 2. Затем каждый из интегралов приближенно
заменяется линейно комбинацией значений подынтегральной функ-
ции. Подобный прием приводит к получению системы разностных
уравнений. В следующем параграфе мы рассмотрим способы построе-
ния наиболее часто употребляемых разностных схем и на несколь-
ких примерах постараемся проиллюстрировать их характерные осо-
бенности. Остальная часть главы будет посвящена рассмотрению
различных квадратурных методов.
По способу вычислений различают экстраполяционные
и интерполяционные численные методы. Методы экстраполяционного
типа предусматривают вычисление каждого «нового» значения ис-
комого решения за один вычислительный цикл, состоящий из конеч-
ного, строго определенного количества отдельных операций. Методы
интерполяционного типа вводят в действие процесс последователь-
ных приближений. |Количество операций при вычислении каждого
отдельного значения искомого решения зависит здесь от требуе-
мой точности и сходимости процесса. В данной книге мы почти
исключительно будем рассматривать методы экстраполяцион-
ного типа.
Наконец, третьим признаком для классификации может служить
способ приближения, используемый при построении прибли-
женных формул. (По этому признаку мы будем различать формулы,
получаемые путем приближения в смысле, наименьшего отклонения
на данном отрезке оси независимой переменной, и формулы, полу-
чаемые путем приближения в смысле минимума среднего квадрата
погрешности в частотной области.
При программировании универсальных цифровых вычислитель-
ных машин почти исключительно используются приближенные мето-
ды. В случае же цифровых моделей благодаря преимуществам, да-
ваемым специализацией, в отдельных случаях, особенно при реали-
зации заданных функциональных зависимостей, оказывается целе-
сообразным использовать парные разностные уравнения. Один из
методов нахождения таких уравнений для весьма ограниченного кру-
га задач будет описан ниже,
57

14. РАЗНОСТНЫЕ СХЕМЫ
Здесь мы рассмотрим практические приемы составления разно-
стных схем для решения обыкновенных дифференциальных уравне-
ний, 'разрешенных относительно старшей производной
Основная идея этих лриемов состоит в том, что производные, вхо-
дящие в исходное уравнение (78), приближенно заменяются линей-
ными комбинациями конечого числа значений зависимой переменной,
взятых в |узловых точках
где Xi
I и k
точки разбиения оси независимой переменной;
положительные целые числа;
постоянные коэффициенты.
Выражения (79) мы будем называть разностными формулами,
а величину ^ — порядком разностной формулы. В результате такой
замены получается разностное уравнение, или, как иногда говорят,
разностная схема.
Для составления разностной схемы необходимо провести сле-
дующие преобразования:
1. В исходном дифференциальном уравнении (78) производится
замена независимой переменной по формуле
(80)
где /i=Xi+i—;^i — выбранный интервал разбиения оси независимой
переменной. В результате такой подстановки для s-й производной
в уравнении (78) имеем:
(81)
58

Подставляя (81) в уравнение (78), получаем:
(82)
2. Производные, входящие в состав уравнения (82), приближен-
но заменяются по формулам (79). Результатом такой замены яв-
ляется разностное уравнение, устанавливающее зависимость -между
решетчатыми функциями R^, Д(^~^), ..R независимой неременной
Заметим теперь, что переменная g в соответствии с выражением
(80) принимает в узловых точках целочисленные значения. Полагая
^п = п, запишем полученное разностное уравнение в окончательном
виде:
1
Rr{y (п)} = hrf Rr-^^y Щ ^—i ,...,/? Щ - , у (п), hn
(83)
Рассмотрим теперь две наиболее употребительные разностные
формулы:
1) разностная формула первого порядка (1 = п, k = l)
Ri {у(п)} = у{п+\)-у{п);
R'l {у Щ = у(п + 2)-2у(п + 1) + у («);
^r = |](-i)"(:),(n + /n-a):
0=0
(84)
2) разностная формула второго порядка (/ = п+ 1, fe = 2)
Rt{y{n+l)} = -2[У('г + ^)-УШ'
т
а = 0
(85)
Каждая из этих формул обладает своими характерными осо-
бенностями. Чтобы выявить их, мы рассмотрим несколько простых
дифференциальных уравнений и составим для них разностные схе-
мы, пользуясь обеими фор.мулами, а зате.м проанализируем решения
полученных разностных уравнений, сравнивая их с точным реше-
нием исходного дифференциального уравнения.
Пример 1. Рассмотрим неоднородное дифференциальное
уравнение первого порядка
_dy_
dx
= а
(86)
59

при начальном условии у{0)=уо. Точное решение этого уравнения
имеет вид
уМ = Уо + С1Х. (87)
Теперь составим две разностные схемы, предварительно лреобразо
вав уравнение '(86) (вводя подстановки x=nh) к виду
y' = ah, (88
что, очевидно, совпадает с формой записи дифференциального урав
нения (83).
1-й способ. Пользуясь приближенными равенствами (84), за
меним уравнение (86) следующим разностным уравнением:
у {п + \)—у (п) = ah. (89;
Это неоднородное разностное уравнение первого порядка [Л. 10],
Его решение будем искать в виде суммы
y(n)==Y{n) + y*(n). (90;
где Y(n) — общее решение однородного уравнения у{п+\)—'У(п)=0,
а y*in)—частное решение неоднородного уравнения (89). Для оты
екания Y(n) составим характеристическое уравнение
X—1 =0.
Оно имеет единственный корень ^i = l. Тогда
У(я) = СЛ^ = С,
где С — произвольная постоянная. Полагая /г=0, получим:
(^=Уо-
Частное решение у* {п) будем искать в виде
y*(n)==A,nh. (91)
Подставив выражение (91) в уравнение (90J, определим
Ло = а.
Тогда
y''{n)==anh, (92)
а
У (п) = Уо + anh.
С учетом соотношения x—nh получаем окончате-1ьно:
у{х) = Уо + ах. (93)
2'й способ. Теперь составим разностное уравнение, пользуясь
равенствами (86),
у(п + 2}-у(п)^2аН. (94)
60

При этом .необходимо задать уже два начальных условия, напри-
мер у(0)=Уо, ^(—1) = ^-1. Выражение (94) представляет со'бой не-
однородное разностное уравнение второго порядка. Как и прежде,
решение его ищем в виде
y(n) = Y(n) + y^{n),
где
У(п) = С,\^+С2Ц (95)
— решение однородного уравнения. Из характеристического урав-
нения
Х2 - 1 = О
определяем == 1; == —1 и, подставляя эти значения в выра-
жение (95), получим:
Y(n) = C, + C2(-\)-. (96)
Значения постоянных С, и Cz найдем из уравнения (96), полагая
п = 0 и п = — и
Ci = У^ + У-^ .
2
Уо — У-1
Тогда
y(„)==i^^±^ + ^i^(_l)n (97)
Частное решение, как и прежде, будет:
у* (п) = anh.
Таким образом,
у in) = + (- 1)" + anh. (98)
Или с учетом равенства ^ х = nh путем несложного преобразова-
ния окончательно имеем:
x
у [х) = {у о + ах) —
(99
Сравним теперь полученные результаты. Сопоставляя решеьгие
(93), полученное с применением первой разностной формулы, и
решение (99), полученное с применением второй разностной форму-
лы, с точным решением исходного дифференциального уравнения
(87), мы может сделать следующие выводы:
а) Выражение (93) в точности совпадает с выражением (87).
Это общее свойство первой разностной формулы, состоящее в том,
что она дает точный результат во всех случаях, когда правая часть
исходного дифференциального уравнения—постоянная величина
(полином нулевого порядка).
61

б) Выражение (99) состоит из двух членов. Один из них
УоЛ-ах представляет собой основное решение (см. гл. 2), которое
в данном случае совпадает с выражением (87). Второй из них
X
___ Уо — У-1 I Уо — У-1 ,
2 2
представляет собой побочное решение. Заметим, что значение y_i
можно выбирать произвольно, поскольку для исходного дифферен-
циального уравнения (86) задано только одно начальное условие
Уо. iB данном случае, полагая У-\ = Уо, можно сделать побочное ре-
решение тождественно равньш шулю. Однако в общем случае, ког-
да аналитическое решение разностного уравнения неизвестно, возни-
кает очень серьезная задача выбора второго начального условия
У-\. Как видно из (99), неправильный выбор значения у-\ может
привести к появлению очень большой погрешности. Например, если
^/-1 = 0, то побочное решение приобретает вид
-f+f (-1)"^.
т. е. вводится погрешность, максимальное значение которой равно
\уо\. Наличие побочных решений является общим свойством второй
разностной формулы.
Пример 2. (В качестве второго примера рассмотрим опять не-
однородное дифференциальное уравнение первого порядка, правой
частью которого является полином первого порядка,
dy
-^ = ах (100)
при начально.м условии у{0) = уо. Точное решение этого дифферен-
циального уравнения имеет вид:
y{x) = y(0) + Y^'' (1^^)
Как и в прошлый раз, преобразуем исходное уравнение, вводя за-
мену x = nh,
y'=:anh^ (102)
и снова составим две разностные схемы.
J-й способ. Используя приближенные равенства (84), полу-
чим следующее неоднородное разностное уравнение первого порядка:
у{п+1)-у{п) = апк2 (103)
Очевидно, что уравнения (103) и (89) отличаются только правыми
частями, и поэтому общие решения однородных уравнений для них
совпадают. Частное решение неоднородного уравнения (103) будем
искать в виде
уЦп)=ЛопН + Л,(пк\2 (104)
62

Сопоставляя решение (106), полученное с применением первой раз-
ностной формулы, и решение (109), полученное с применением вто-
рой разностной формулы, с точным решением (101) исходного диф-
ференциального .уравнения, приходим к следующим выводам:
а) Выражение (105) отличается от точного решения (101)
а
членом — Здесь, так же как и в любом другом случае,
когда правая часть исходного дифференциального уравнения не яв-
ляется постоянной, первая разностная формула дает погрешность,
причем, как видно, из выражения (106), эта погрешность может не-
ограниченно возрастать с ростом х. Заметим, однако, что при стрем-
лении /г к О выражение (106) сходится к выражению (101).
63

6) Решение (109), так же как и в предыдущем случае, состоит
из двух членов. Один из них
а
представляет собой основное решение и в точности совпадает с вы-
ражением (101). Это общее свойство второй разностной формулы,
состоящее в том, что ее основное решение всегда оказывается точ-
ным, если правая часть исходного дифференциального уравнения
представляет собой полином не выше первого порядка. Второй член
x
_ Уо — У-1 , Уо — У-1 .
2 2 ^""^^
представляет собой побочное решение, относительно которого спра-
ведливо все сказанное при рассмотрении предыдущего примера.
Пример 3. Рассмотрим теперь линейное однородное диффе-
ренциальное уравнение второго порядка
с начальными условиями у (0) = у' (0) = у'\.
1-й способ. После необходимых преобразований получаем разно-
стное уравнение
y(n + 2)-2y{n+\) + (\ + h'ii>^)y{n)==Q. (Ill
Оно представляет собой линейное разностное уравнение с постоян-
ными коэффициентами. Характеристическое уравнение для него за-
пишется в виде
Х2 — 2Х + (1 + h^<^^) = 0. (112)
Корнями этого характеристического уравнения будут:
Л» = 1 + jhi^ =
X, = 1 - /to =
(113)
где, очевидно.
(114)
Р = arctg Лео,
Общее решение разностного уравнения (111) имеет вид
у (п) = ^'^''(Cicos -f Cgsin f л). (115)
где Cl и Сг — произвольные постоянные.
64

Сравним теперь оба полученных решения с точным решением
уравнения (110), которое, очевидно, имеет вид
у (дг) = COS сох + ^ sin (ЛХ. (126)
Сопоставляя выражения (126) и (117), мы видим, что в данном
случае использование разностной схемы, составленной по первому
способу, приводит к появлению погрешностей трех самостоятель-
ных видов. Первая из них—амплитудная обусловливается тем, что
а
в вмражении (117) имеется множитель е'^ . Этот множитель, а
следовательно, н связанная с ним погрешность неограниченно воз-
растают с ростом X по экспоненциальному закону. Второй вид по-
грешности — частотный. Действительно, круговая частота в гармони-
ческих членах выражения (117) равна
arctg^/ico
F h
и, следовательно, отличается от истинной круговой частоты о.
Наконец, третий вид погрешности, фазовый, зависит от выбора
начального условия уи Этот вид 'погрешности может быть полно-
стью исключен, если положить
У1 = Уо + hy\.
Однако в общем случае, когда аналитическое уравнение для ре-
шения разностного уравнения «ам неиз'вестно, величина фазовой
погрешности будет зависеть от способа задания начального усло-
вия Уи
Сравним теперь выражение (126) с решением (125), полученным
с помощью разностной схемы, составленной по второму способу.
Заметим прежде всего, что амплитудная ошибка здесь отсутствует.
Частотная погрешность имеет место, поскольку величина ""^^
arcsin to
= д отличается от истинного значения круговой частоты со.
Фазовая погрешность также имеет место и также может быть ис-
ключена, если положить
^, = г/о cos р + hy\.
В заключение этого примера заметим, что как амплитудная, так
п частотная погрешности стремятся к нулю при стремлении к нулю
величины h. Фазовая же погрешность определяется выбором началь-
ных условий.
Сделаем теперь несколько общих замечаний относительно ме-
тода решения обыкновенных дифференциальных уравнений с по-
мощью разностных схем. Из рассмотренных примеров видно, что
в обптем случае подобному методу решения свойственны два вида
погрешностей. Один из них обусловливается тем, что производные
66

в исходном дифференциальном уравнении приближенно заменяются
разностными формулами. Эта погрешность, как правило, уменьшает-
ся, стремясь к еулю ири стремлении к нулю интервала разбиения h.
Поэтому 1может быть предложен эффективный способ контроля точ-
ности получаемых решений. Согласно этому способу решение сна-
чала проводится при некоторой фиксированной величине h. Затем
эта величина уменьшается вдвое и решение повторяется. Затем оба
полученных решения сравниваются. Если они отличаются друг от
друга на величину, большую некоторой наперед заданной величин?.!
погрешности, интерн ал разбиения h снова уменьшается вдвое и так
до тех пор, пока не будут получены два последовательных решения,
достаточно мало отличающихся друг от друга.
Второй вид погрешности обусловлен наличием побочных реше-
ний, которые в свою очередь возникают благодаря тому, что при
использовании разностных формул лорядка выше первого получае-
мое разностное уравнение имеет порядок, больший, чем порядок
исходного дифференциального уравнения. Погрешности этого вида
в общем случае не уменьшаются при стремлении h к нулю и зависят
только от выбора начальных значений. 'Поэтому в тех случаях, ко-
гда в наличии имеется мало данных о характере решений ис-
ходного дифференциального уравнения, пользоваться разностными
схемами следует с большой осторожностью. Заметим далее, что ес-
ли сами побочные решения возрастают достаточно быстро с ростом
независимой переменной л, то даже при правильном выборе на-
чальных значений достаточно одного олучайного сбоя в вычислитель-
ном устройстве или просто достаточно большой погрешности из-за
округления, чтобы соответствующее побочное решение стало отлич-
ным от нуля и в дальнейшем возрастало. В таких случаях говорят,
что построенная разностная схема является неустойчивой. Подробно
вопрос об устойчивости численных методов решения дифференци-
альных уравнений рассмотрен -в [Л. 42 и 51].
15. КВАДРАТУРНЫЕ ФОРМУЛЫ
Большинство применяемых в настоящее время приближенных
численных методов решения обыкновенных дифференциальных урав-
нений основано на использовании так называемых квадратурных
формул. Эти чиоленные методы в свою очередь могут быть подраз-
делены на три основных типа: экстраполяционные, интерполяционные
и формулы Рунге-Кутта.
Формулы экстраполяционного типа. Предположим, что исходное
дифференциальное уравнение порядка г задано в виде системы г
дифференциальных уравнений первого порядка
j|^=fi(^b ^2...., ^»,..., Уг. X), (127)
/= 1, 2, 3,..., г
при начальных условиях
yi {Хо) == ^го»
/=1, 2, 3,..., Г, (128)
где X — независимая переменная.
5* 67

Задача состоит в том, чтобы вычислить приближенные значе-
ния переменных yt в точках Xi, Xg,... оси независимой переменной
= -^0 + rih,
где h — постоянный шаг разбиения оси х.
На основании (127) можно написать:
yi(Xh + i) = yi(Xh) +
Будем предполагать, что нам известны k + 1 значений переменных
Уг(Хо), yi(xi),..., yi(Xk). Тогда могут быть вычислены также зна-
чения
^-^^ = fi[yA^j)^''-^ Уг{х^). хл, (130)
7 = 0, 1, 2,..., k.
Для упрощения дальнейших выкладок введем обозначения:
Y(Xj)=\\yA^j)^-"^ Уг{х^)\\ и fi(xj) = fi[Y(Xj), хЛ.
Аппроксимируем теперь функцию fi полиномом степени ky исполь-
зуя при этом известные ее значения fi{Xj). Для такой аппрокси-
мации можно использовать, например, экстраполяционную формулу
Ньютона [ЯА]:
+-^-+'^-,r^-'4^f(x.). (.3,)
где
Подставляя теперь выражение (131) в формулу (129), получим не-
посредственно:
yiiXk + i)=^yi{Xk) + h \ fi{Xk+8)ds^
о
k
^yi(Xk) + h Y^ajvni(Xk), (132)
где
68
,^=j£ii±ib_^£Jiizii),, (.33)
0

произведя вычисления по формуле (132), получим:
,3 251 95 \ ^
(134)
Для определения погрешности такой аппроксимации воспользуемся
величиной остаточного члена экстраполяционной формулы Ньютона
С учетом выражения (132) и очевидного соотношения
(135)
получим:
Е, = I ^ + Ц + у, 1^+^] (^) ds, (136)
О
где Ei — погрешность, получаемая при вычислении одного значе-
ния yi(Xn + i) по формуле (134). Учитывая далее, что коэффициент^
при 1/1^не меняет знака на интервале (0,1), получим:
(137)
где
£i = aA+:/i(''+=)(/.-I*+4(=),
Хп + 1 ^ 5 _ h'
в более общем случае вместо выражения (136) можно восполь-
зоваться выражением
1
yi (Xk + i) = yi (Xk - p) + Л j fi (Xh + s) ds,
(138)
где /? —положительное целое число. Выкладки, подобные преды-
дущим, приводят к группе формул
Уг (Xk + г) ^ yi (Xk - о + /1 ^2 + Ov +
1 1 29 14 \
+у V= + 1" V'+90 ^+45 V= + • • • jf < (-^О
Выражение для ошибки аппроксимации имеет при этом вид
1
(139)
= + yi^^4{yi)ds. (140)
-Я
69

Формулы, полученные при целом и нечетном /?, обладают тем
интересным свойством, что в них коэффициент при р-й разности
равен нулю. Следовательно, они дают возможность, удерживая
только (р—1)-ю разность, 'получать точность, которая нормально
достигалась бы при удержании р разностей в разложениях вида
(139).
Пользуясь описанными методами, можно получить целый класс
квадратурных формул, отличающихся друг от друга значениями
величин k я р. Некоторые из них приведены в табл. 1. 'Все эти фор-
мулы будут относиться к формулам экстраноляционного, или «откры-
того», типа, поакольку они предусматривают вычисление нового зна-
чения yi(xh + i) по известным значениям yi(Xk-j). Полная система
разностных урашений, позволяющая получать приближенное реше-
ние системы 'уравнений ((127) в случае использования квадратурных
формул экстраполяционного типа, будет иметь вид
k
yi (Xn+l) ^ yi (Хп-р)+ h^djfi {Хп ~ i),
/= 1, 2,..., г.
(141)
(Из выражения (141) видно, что для начала вычислений необхо-
димо знать, кроме начальных условий (1'28), еще k «первых» зна-
чений yi(xi), ..., yi(Xk). В дальнейшем каждый (шаг вычислений
по формулам (141) разбивается на два этапа: вычисление значений
fi(xn-j) по полученным в течение предыдущих шагов значениям
y(xn-j) и вычисление значений y(Xn+i). Некоторые способы вычис-
ления функций fi будут рассмотрены ниже. (Величину k мы будем
называть в дальнейшем порядком квадратурной форхмулы.
Формулы интерполяционного типа. Формулы интерполяционного
типа отличаются от рассмотренных выше тем, что для аппроксима-
ции функции fi используется интерполяционный полином
f{^n + s)'
f iXn + l) + (5 - 1) Vf (^n+ ,) + V'f (^n-b,) + . . .
...+
2!
(s^l)s{s+\)...{s + k-2)
v4{xn^iy (142)
Подводя выкладки, аналогичные рассмотренным выше, получаем
формулы:
yi
29^^"~ 720^*'" 160^^"'•••j fi(-^n + i)\
yi (Xn + i) ^ yi (Xn - i)+h(2 - 2v + 4 +
(143)
70

Таблица 1
Некоторые формулы численного интегрирования
о tt
ts ?!
ad.
о о
I. Экстраполяционные формулы
Примечание
0
Уп + 1 = Уп + hy'n
Формула Эйлера I
1
Уп + 1 = yn-i + 2hy'n
Формула Эйлера II
1
1
Уп + 1— Уп-1+ 2 ^(^У'п—y'u'i)
Формула трапеций
2
3
3
Уп + 1-=Уп-2+ h{y'n-2 + ^y'n)
4
Уп + 1—Уп-3+ 3 — г/'п-1 +
Формулы Нью-
тона—Котесса
2
2
yn+i = yn + h(l,8S\77Sy'n-
- l,294418r/',_i+0,41218V,.2)
yn + i = yn-i + h(2,237S32y'n-
Метод приближе-
ния в среднем в ча-
стотной области (по-
лучен автором)
— 0,629006f/'^_i +0,327549f/'^_2)
II. Интерполяционные формулы
0
+ i == Уп + hy'n + i
Уп + 1—Уп+2 ^(y'n + y'n + i)
Формула прямо-
угольников
1
Формула трапеций^
2
1
Уп + \'= Уп-1+ 3 /i(^n-l + Vn +
Формула Симп-
+ ^'n + l)
сона
3
3
УпЛ-1 — Уп-2+ 8 /г(^'п-2 + 3(/'„>1 +
3 •
Формула g Симп-
сона
71

Все получаемые таким образом формулы относятся к формулам ин-
терполяционного, или «закрытого», типа, так как они предусматри-
вают вычисление значений yi{Xn+i) на основании k-\-\ значения
функции /г, включая сюда fi(Xn+i). Некоторые из таких формул
также приведены в табл. .1.
72

Только что отмеченное обстоятельство не дает возможности не-
посредственно использовать квадратурные формулы интерполяцион-
ного типа для решения обыкновенных дифференциальных уравне-
ний, так как значения /г((л:„ + 1) согласно (427) могут быть вычислены
только в том случае, если уже известны значения Уг(Хп + 1). 1По-
этому при использовании квадратурных формул интерполяционного
Т1ша процесс вычисления разбивается на следующие этапы:
1) Приближенные значения yi(Xn + i) вычисляются по одной из
формул типа (141). Эти значения мы будем называть нулевым при-
ближением и обозначать как yf^ (Xn + i)- Используемая при этом кон-
кретная формула типа (141) называется предсказывающей.
2) Значения f/(°)(Xn + i) подставляются в уравнения (127), вычис-
ляются значения ff^{Xn + i), которые мы также будем рассматри-
вать как нулевые приближения значений fiiXn + i)-
3) Полученные таким образом значения fP(Xn-ti) подставляют
в одну из формул интерполяционного типа, которая в этом случае
рассматривается как „исправляющая". С помощью этой формулы
вычисляется первое приближение значений y[^\xn + i).
4) Значения у1{Хп + г) подставляются в уравнение (127) для
вычисления значений j\^^Xn+i).
5) С помощью исправляющей формулы вычисляют значения
y^i'Hxn + x).
6) Этапы 4 и 5 повторяют до тех пор, пока разность
y^^Hxn + i)— tfi^^Hxn + i) не станет меньше некоторой наперед за-
данной величины.
Сходимость подобного процесса может быть значительно улуч-
шена, если воспользоваться одним специальным приемо.м, который
мы рассмотрим на примере построения так называемых формул
Милна [Л. 34].
Пусть в качестве предсказывающей используется формула
yf^ {Хп ^i) = yi {Хп - з) + 4 ^ (^гг) - fi (Хп - i) + 2fi {Хп -s)]- (145)
Здесь под ^г(^п-з) и fi{Xj) мы понимаем „точные" значения, вы-
численные на предыдущем шаге. Введем также обозначение
Pi{Xn + i)=y^^\xn + i). Погрешность, с которой вычисляется Pi{Xn+i),
28
согласно (140), равна gg ('^)- Пусть в качестве исправляющей
используется формула
y^/f^{Хп+х) = y^^-'Hxn-i) + J h [f^-'^{Xn+x) +
+ ^f^r'^(Xn)+f^r'^{Xn..)l (И6)
^=1,2,.,,
73

Тогда при условии, что выполнено 'достаточное количество циклов
последовательных приближений, погрешность вычисления величины
yi{Xn-i) будет согласно (144) равна — qq/I'^^pl ($). Если предполо-
жить теперь, что в пределах рассматриваемого интервала измене-
ния X величина у\^^ остается практически постоянной [у^^^ (?) =
29
= i'P'(v) = #']. то Pi(x„^r)-yi(Xn.,) = QQll'yfy
На практике можно предположить, что рг(Хп)—yi(Xn) очень 1мало
изменяется от шага к шагу. Тогда (МОЖно ввести в рассмотрение
«модифицированное» значение
28
mi{Xn + i) = Pi {Xn+i) — 2^[Pi (Хп) — У г {Хп)]> (147)
для которого погрешность предсказывающей формулы (145) почти
полностью окажется скомпенсированной.
Значение Шг(;Сп+1) подставляется затем в формулу '(146), при-
чем, поскольку погрешность этого значения достаточно мала, в боль-
шинстве случаев оказывается достаточным привести один цикл вы-
числений по формуле |(146). Таким образом, мы приходим к построе-
нию следующего алгоритма:
1. Предсказание
4
Pi (Xn+i) = yi {Хп~г) + ^ h [2/i (Хп) — fi (Хг,.г) + 2fi (Хп~2)]'
2. Модификация
28
mi {Xn+i)=Pi {Xn+i)— 29 [fi (^n) — Ci (Xn)J.
3. Исправление
Ci(Xn+x) = yi (Xn-,) + у /г [m'i (x„+,) + 4f»(Xn) + fi (x.,, -1)].
Здесь m'i{Xn+i) = fi Hi (x„+,),.,., (Xn+i), ^n+i]-
4. Вычисление окончательного значения
yi {Xn+i) = с i (Xn+i) + ^ [Pi (Xn-^-i) — Ci (x„+,)j.
Величины Pi(Xn)—Ci(Xn) могут использоваться в качестве критерия
правильности вычислений. Если значения одной или нескольких из
этих величин, медленно изменяясь, выходят за пределы некоторой
заданной области, то это указывает на неправильный выбор вели-
чины h или порядка k квадратурных формул. Резкое внезапное уве-
личение значения Pi(x7i)—Ci{Xn) почти однозначно свидетельствует
об ошибке в вычислениях.
74

Некоторые формулы РунгемКутта приведены 'В табл. 1. Погрешность
формул Рунге-Кутта имеет порядок величины первого отброшенного
члена разложения Тейлора.
16. КВАДРАТУРНЫЕ ФОРМУЛЫ, ОСНОВАННЫЕ
НА ПРИБЛИЖЕНИИ В СРЕДНЕМ
ТРИГОНОМЕТРИЧЕСКИМИ РЯДАМИ
Рассмотренные выше способы построения квадратурных формул
основывались на аппроксимации функции fi степенньими полиномами.
Однако такая аппроксимация не всегда оказывается наилучшей. Мо-
жет быть предложен другой способ построения формул приближен-
ного решения обыкновенных дифферендиальных уравнений, основан-
ный на аппроксимации тригонометрическим рядом.
:Будем считать, что нам известно значение какого-либо решения
у(х) дифференциального уравнения (\27) в некоторой произвольной
точке Хг. Для вычисления следующего значения этого решения в точ-
ке Xi + h воспользуемся формулой
k
y(x + h)^y{x) + hj^ а J {X ~ v/z), (150)
v=0
где — неизвестные пока постоянные коэффициенты. Точное зна-
чение yix-^-h) может быть получено из выражения
x+h
y(x + h) = y(x)+ J f(s)ds. (151)
x
Применяя преобразование Фурье к выражениям (150) и (151) и
обозначая со/г = со получим:
k
(е^^-1)^(7;;)^/г/(«) ^ а,^-/- (152)
75
Формулы Рунге-Кутта. Формулы Рунге-Кутта основываются на
соотношении
причем величины а, \1 и Ъ выбираются таким образом, чтобы разло-
жение правой части формулы (148) полностью совпадало с k пер-
выми членами правой части разложения Тейлора,

{ei'"-l)y{w)=^h}(co) ^ _ ' , (153)
A A
где y{(o) и f (со) — преобразования Фурье функций у ух) и f (s) со-
ответственно.
Из сравнения выражений (152) и (153) видно, что в такой по-
становке задача сводится к отысканию совокупности значений
коэффициентов а^, обеспечивающей наилучшее приближение в
среднем на некотором заданном отрезке |(л> | <^ [а функции
е^'""- 1
/(О
суммой
Рассмотрим более общую задачу нахождения наилучшего при-
ближения в среднем на заданном отрезке функции F(ы) суммой
k
где (о)) — некоторые функции переменной со, ортогональные на
отрезке — р., [х.
Как известно, искомое наилучшее приближение в среднем бу-
дет получено в том случае, если постоянные Л^ будут представ-
лять собой коэффициенты Фурье, определяемые из вырал<ения
и-
j F (^) ф% (©) flfw
j (*^) ^^*v (^)
где Ф*у(со) — функция, сопряженная с функцией ф^(со).
В качестве Ф^(со) молено выбрать систему функций:
Фо(«)=1;
Ф1(ю) = ^-^"^ _Л,^ф^с;^;
Ф, (^) = е-'^'^ - М - Ко^о (^); \ (^
ф, (^) = ^-^^"^-Х,.,.,ф,.1(^;)^..,-Хл^о(с7),
где А./,^ — постоянные коэффициенты.
76

Значения этих коэффициентов получаются непосредственно из ус-
ловия ортогональности
j Ф/ Ш'ш Й сйа=^0- I ф т. (157)
Для практического вычисления значений коэффициентов можно
воспользоваться следующими рекуррентными формулами:
т
J] «m,pI^pSin(/-p)lx
'-^^—<—■■
1 при т = р;
т—Р
2j Р(— ^m.m-i) при т^р,
где ^f?= j Фг((о)ф*,-Й^с7.
-jj.
Соответствующая рекуррентная формула для вычисления зна
чения имеет вид
(158)
(159)
р=0
Несколько первых значений Xi^rn и М- вычисленных по форму-
■к п п
лам (158) и (159) для jj. = -g-, ^ "З » сведены в табл. 2.
Подставляя значения Л/^^ и Л1? в выражение (155) и полагая
получаем формулу для вычисления коэффициентов Л^. Эта формула
имеет вид:
(160)
где значения ^ вычисляются по формуле (158). Несколько пер-
вых значений коэффициентов А^, вычисленных по формуле (160) для
тс тс 71
[k = '^j и у , также сведены в табл. 2.
77

Таблица 2
тс
тс
т
т
т
0,9744952
0,900316
0,826993
0,9003162
0,636619
0,413497
X,.
1,9289768
1,727056
1,534517
Ml
0,0395519
0,297557
0,662002
Ml
0,0016081
0,046632
0,177453
Ло
0,991475
0,966258
0,941044
Ai
—0,499317
—0,479792
—0,465374
Л,
0,412189
0,399850
0,385856
Подставляя значения коэффициентов А^, вычисленных по фор-
муле (160) в выражение (154), получаем разложение
— h
v=0
(161)
Наконец, подставляя в выражение (161) значения функций Ф^(со),
взятые из (156), получаем разложение вида (152), в котором
Это и есть искомое разложение, обеспечивающее наилучшее при-
ближение в среднем на отрезке —^[г, |л. Выражения (150) и (162) опре-
деляют класс формул численного интегрирования. (Некоторые из этих
формул приведены в табл. 1. Все приведенные выше соображения
легко обобщить на случай системы дифференциальных уравнений.
17. СПОСОБЫ ВЫЧИСЛЕНИЯ ФУНКЦИЙ
Все рассмотренные выше способы позволяют построить разно
стное уравнение, решения которого с определенной степенью точно
сти приближают решения исходного дифференциального уравнения
В то же время, как было показано Р. Калманом [Л. 20], для до
статочно широкого класса обыкновенных дифференциальных урав
нений существуют так называемые парные разностные уравнения
Парным называется разностное уравнение, решения которого в. оТ'
дельных дискретных точках оси независимой переменной в точности
совпадают с решениями соответствующего дифференциального урав-
нения. Однако способы построения парных разностных уравнений
для общего случая в настоящее время еще не разработаны. Более
того, имеющиеся способы нахождения парных разностных уравне-
ний для линейных дифференциальных уравнений приводят к чрезвы-
78

чайно громоздким расчетам, и использование их ири цифровом мо-
делировании вряд ли )Можио признать целесообразным.
|Г1ри цифровом 'Моделировании объектов, работа кото1рых опи-
сывается нелинейными дифференциальными уравнениями, возникает
самостоятельная задача вычисления значений функций, входящих
в состав этих уравнений. При этом существенное значение приобре-
тают .следующие два обстоятельства. iBo-первых, вычисление функций
в общем случае должно производиться с точностью большей, чем
требуемая точность самого решения. Подробнее этот вопрос будет
рассмотрен в следующей главе. iBo-вторых, класс функций, встре-
чающихся иа практике, достаточно сильно ограничен.
Одним из используемых в настоящее 1время способов вычисле-
ния функций является так называемый способ «производящих» диф-
ференциальных уравнений, согласно котороцу для получения зна-
чений данной функции решается дифференциальное уравнение, од-
ним из решений которого является эта фун^щия. Обычно в цифро-
1вых моделях производящие дифференциальные уравнения решаются
та'Кже приближенными методами. Однако здесь имеет прямой смысл
использовать парные разностные уравнения. Один из способов по-
строения разностного уравнения, точным решением которого являет-
ся заданная функция, состоит вг следующем. К функции независи-
мой переменной F{x), заданной аналитически, применяется дискрет-
ное преобразование Лапласа, а затем по полученному изображению
строится разностное уравнение. Поскольку не по всякому изображе-
нию можно построить разностное уравнение, допускающее реализа-
цию, способ не является универсальным. Если функция задана своим
дифференциальным уравнением, то к нему сначала применяют обыч-
ное преобразование Лапласа, а затем с помощью ^-преобразования
[Л. 46] отыскивается соответствующее дискретное преобразование
Лапласа. В рез1ультате 1может быть построено разностное уравнение,
парное к исходному дифференциальному уравнению.
Закачивая главу, посвященную рассмотрению алгоритмов цифро-
вых моделей, необходимо отметить, что в настоящее время не суще-
ствует методов получения алгоритмов, позволяющих находить точ-
ные решения любого обыкновенного дифференциального уравнения.
Можно предположить также, что воз1можность получения таких ме-
тодов в будущем является сомнительной. Это совпадает с основными
положениями численного анализа. Из сказанного следует, что одним
из преимуществ, даваемых специализацией, и, следовательно, пре-
имуществ цифровых моделей 1является возможность подбора алго-
ритма для каждого конкретного устройства и связанная с этим воз-
можность достижения максимальной эффективности.
ГЛАВА ПЯТАЯ
ХАРА1КТЕРИСТИКИ ЦИФРОВЫХ МОДЕЛЕЙ
18. ИСХОДНЫЕ СООБРАЖЕНИЯ
Вопрос об оценке качества реальной конструкции цифровой
модели в большой степени связан с условиями ее эксплуатации
и с теми конкретными техническими условиями, которые были поло-
жены в основу ее разработки. Из двух конструкций, моделирующих
79

один и тот же объект, лучшей в каждом данном случае будет счи-
таться та, которая в наибольшей степени отвечает конкретным
условиям эксплуатации и требует при этом наименьших затрат. Сре-
ди ряда критериев, которые можно принять для оценки качества
цифровой модели, в первую очередь нас будут интересовать крите-
рии, допускающие четкое количественное выражение. С их помощью,
как это будет показано ниже, может быть решена задача синтеза
наилучшей в некотором определенном смысле структуры устройства.
К числу таких критериев прежде всего следует отнести точность,
быстродействие и относительную сложность.
Точность является основной характеристикой, определяющей
применимость данной конструкции. Требования к повышению точно-
сти и привели к появлению цифровых вычислительных машин и,
в частности, цифровых моделей. С другой стороны, чрезмерная точ-
ность, как правило, приводит к неоправданным затратам оборудо-
вания или времени вычисления. Поэтому точность обычно прини-
мается в качестве основного пункта технических условий.
Быстродействие определяет динамические характеристики моде-
ли. Проблема повышения быстродействия ставится особенно остро
при работе в истинном масштабе времени, т. е. при использовании
моделей в составе систем автоматического управления. Наконец,
относительная сложность характеризует такие важные свойства мо-
делей, как удобство эксплуатации, первоначальная стоимость, экс-
плуатационная стоимость, габариты, вес и надежность.
В зависимости от реальных условий эксплуатации каждая из
рассмотренных характеристик может приобретать первостепенное
значение. Например, при использовании модели в составе бортового
оборудования самолета первостепенное значение приобретают вес
и габариты, т. е. в конечном итоге сложность. Можно указать так-
же ряд применений, где первостепенное значение играет надеж-
ность, например при использовании модели в составе системы управ-
ления атомным реактором. В этих случаях увеличение надежности
часто покупается ценой резервирования, т. е. увеличения сложно-
сти. При необходимости управления быстро протекающими процес-
сами основное значение приобретает быстродействие и т. д.
С другой стороны, рассмотренные характеристики связаны м(^)к-
ду собой и улучшение одной из них всегда может быть достигнуто
за счет остальных. В гл. 2 мы уже видели, как при увеличении чис-
ла параллельных каналов ЦДА можно было достичь увеличения
быстродействия при фиксированной точности. При фиксированном
количестве оборудования увеличить быстродействие можно за счет
снижения точности и наоборот. Поэтому вполне естественным яв-
ляется желание получить некоторый единый критерий оценки каче-
ства цифровой модели. Такой критерий может быть получен на
основе информационного подхода, и мы будем называть его эффек-
тивностью. Понятие эффективности имеет много общего с понятием
коэффициента полезного действия для машин — преобразователей
энергии.
При использовании цифровых моделей для анализа различных
объектов, особенно при разработке и проектировании этих объекто-в,
существенную важность приобретают такие вопросы, «ак удобство
подготовки задач, ввод исходных данных и интерпретация получен-
ных результатов. По существу в'се эти вопросы могут быть сведены
■к единому критерию, который 'МОжно назвать эффективностью си-
80

схемы связи между машиной и работающим с ней человеком-опера-
тором. К сожалению, -в настоящее время эта проблема еще недо-
статочно разработана и эффективности системы связи не может
быть дано четкое количественное выражение.
Каждая из рассмотренных характеристик однозначно опреде-
ляется алгоритмом, реализуемым в цифровой модели, характери-
стиками отдельных элементов, используемых при ее построении, и
способом соединения между собой этих элементов, т. е. структурой.
Улучшение характеристик элементов, естественно, влечет за собой
улучшение общих характеристик модели. Однако в дальнейшем мы
будем предполагать, что в распоряжении конструктора всегда
имеется некоторый набор стандартизованных нормальных элементов
и поэтому требуемые характеристики он может получать только
путем выбора алгоритма и изменения структуры. Такой подход
к проблеме является общим для теории конечных автоматов [Л. 23]
и, в частности, для структурной теории цифровых вычислительных
машин.
19. ТОЧНОСТЬ
Здесь, как и 'в предыдущей главе, мы ограничимся только рас-
смотрением погрешностей, возникающих при цифровом моделирова-
нии объектов, работа которых описывается обыкновенными диффе-
ренциальными уравнениями. При этом
всегда необходимо иметь в виду, что iu(hT)
цифровое моделирование объектов, иод- ^
вергающихся нроизвольньгм внешним
воздействиями, всегда может быть толь-
ко нриближенньим. Это положение вы-
текает из самой природы цифрового
представления величин. Действительно,
как отмечалось в гл. 1, цифровое пред-
ставление предусматривает замену ис-
ходной (функции f{x) непрерывно'го ар-
гумента X соответствующей ей решетча-
той функцией 1{п). Следовательно, в
любом случае разностное уравнение,
описывающее работу модели, будет
содержать в своей правой части эту
решетчатую функцию. С другой стороны, решетчатой функции
соответствует бесконечное (множество различных огибающих,
совпадающих между собой в узловых точках. Таким образом,
даже тогда, когда разностное уравнение, описывающее работу мо-
дели, является парным к исходному дифференциальному уравне-
нию, решение будет действительно точным только для одной из
функций, соответствующих данной решетчатой функции.
Сказанное можно проиллюстрировать следующим простым при-
мером. Пусть решетчатая функция, представляющая внешнее воз-
действие и поступающая на вход цифровой модели, имеет вид, по-
казанный на рис. 18. На этом же рисунке нанесены две возможные
непрерывные функции, соответствующие исходной решетчатой. Одна
из них постоянна, а другая изменяется периодически с периодом Т.
Очевидно, что решения простейшего дифференциального уравнения
6 А. в. Шилейко. 81
Рис. 18. Решетчатая функ-
ция, поступающая на вход
модели.

Длй каждой из этих функций могут отли11аться одно от другого Нг1
сколь угодно большую величину.
Разностное уравнение для погрешности. При описании способа
построения квадратурных формул каждого класса проводилась так-
же оценка погрешности аппроксимации [см. выражения (137) и
(140)]. Наличие этой погрешности, а также ряда других факторов
приводит к тому, что решения дифференциальных уравнений, полу-
чаемые с помощью квадратурных формул, вообще говоря, отли-
чаются от истинных решений этих уравнений.
Рассмотрим наиболее общий вид квадратурной формулы экс-»
траполяционного типа:
k k
v=0 v=0
где и 6^ — постоянные коэффициенты. Величину k, как и раньше,
будем называть порядком формулы. Пусть формула (163) исполь-
зуется для решения дифференциального уравнения первого по-
рядка
Предположим, что это уравнение при некоторых заданных началь-
ных условиях имеет решение у(х). Тогда мы сможем написать:
y{^n+i) = Yi b,y(x^^^) +
k
+ hY,ciJ у (x^^,)] + E (x^^i), (164)
v=0
где Xi=xo-{-ih\ £"(jCn+i)—погрешность аппроксимации, получаемая
на (n+l)-M шаге вычислений.
Кроме погрешности аппроксимации, при численном решении
дифференциальных уравнений имеет место также погрешность округ-
ления, которая получается за счет того, что реальные вычисления
по формулам вида (163) выполняются с числами, имеющими ограни*
ченное количество разрядов. При вычислении значений функции
f[xi, y(Xi)] также допускается определенная погрешность. Наконец,
начальное значение г/о тоже может быть задано только с опреде-
ленной точностью. Обозначим приближенные значения функции у(х)^
полученные в результате выполнения / шагов вычислений по фор-
муле (163), через y(i). Тогда можно написать:
yln+n^Y^ b,y[n^v] +
v=0
k
+ aj{n--v, y[n-y])+T{Xn+,). (165)
v=0
82

где T(Xn+i) — погрешность, получаемая на (п+1)-м шаге вычисле-
ний по (163) в результате перечисленных выше причин.
Разность
г[1]^у{Хг)-у[1] (166)
будет представлять собой, очевидно, полную величину ошибки, по-
лученной после выполнения t-ro шага вычислений. Вычитая урав-
нение (165) из уравнения (164), получим:
k
v=0
k
+ £ «V {f k-v^ у (-^,.-v)l -f(n-^, y[n~- V])} +
+ Е{Хп + г)-Т(Хп + ,). (167)
Ha основании теоремы о среднем запишем:
flxu у (Xi)] - f (/, у [i]) = г и,
где 5i —некоторая точка, лежащая на интервале между y(Xi) и
y[i]. Теперь окончательно можно написать:
k
^[п+Ц = ^Ь,г [/2-^v] +
v=0
k
Уравнение (168) представляет собой разностное уравнение с пере-
менными коэффициентами, решением которого является функция
е[/г]. Повторяя аналогичные выкладки, можно построить разностные
уравнения для ошибки, получаемой при использовании квадратур-
ных формул интерполяционного типа и..формул Рунге-Кутта. Эти
уравнения будут иметь в основном тот же нид, что и уравне-
ние (168).
Из рассмотрения уравнения (168) можно сделать вывод, что
поведение функции е [п] зависит от самой квадратурной формулы
(коэффициенты и 6J, от характера истинного решения исходного
f df(xub)\
дифференциального уравнения
величины
и от точно-
сти выполнения вычислений [величины T(Xi)\ Наконец, в уравнение
(168) входит член E(Xi), который, как это следует из сказанного
выше, зависит как от порядка квадратурной формулы, так и от ви-
да исходного дифференциального уравнения. Возможны случаи не-
ограниченного роста функции г[п\ при увеличении п. Тогда говорят,
б* 83

что данная квадратурная формула оказывается «неустойчивой» по
отношению к данному дифференциальному уравнению [Л. 42].
Исследованию устойчивости и точности квадратурных формул по-
свяшено больпюе количество работ. Наиболее полно .этот вопрос
рассмотрен в работах М. Р. Шура-Бура [Л. 51] и М. Урабе [Л. 38],
где выводятся критерии устойчивости квадратурных формул и
даются оценки максимальной величины г[п] как для одномерного,
так и для многомерного случаев.
Однако эти оценки основываются на исследовании свойств ре-
шешш разностных уравнений вида (168) и, следовательно, предпо-
лагают наличие определенного объема априорных сведений о харак-
тере решения исходного дифференциального уравнения. Другими
словами, с их помощью можно судить о пригодности данной квад-
ратурной формулы для численного решения некоторого конкретного
дифс()еренцнального уравнения.
Для цифровой модели характерно наличие жесткого, раз на-
всегда выбранного алгоритма работы системы в целом или во вся-
ком случае ее отдельных решающих блоков. Поэтому для рацио-
нального выбора такого алгоритма было бы желательно иметь си-
стему оценок, основанных на некоторых более общих характери-
стиках класса уравнений, которые предполагается решать с по-
мощью этой модели. Кроме того, при работе цифровой модели
совместно с другими устройствами и, в частности, с реальной аппа-
ратурой выдвигается требование к тому, чтобы при заданной точ-
ности вычислений модель «успевала» перерабатывать полностью
весь объем поступающей к ней информации. Это требование изве-
стно в литературе под названием работы в «истинном», или «на-
туральном», масштабе времени. До последнего времени условия ра-
боты всякой цифровой вычислительной машины в натуральном мас-
штабе времени формулировались в терминах быстродействия этой
машины, исчисляемого средним количеством операций, выполняемых
в единицу времени. Однако подобную оценку нельзя считать
исчерпывающей, поскольку истинные динамические свойства систе-
мы, описываемой разностным уравнением, зависят прежде всего
от вида этого уравнения, а также от организации вычислительного
процесса, выбранного для его решения.
Другой подход к решению проблемы может состоять в том,
чтобы оценивать динамические свойства цифровой модели так же,
как оцениваются динамические свойства моделируемых объектов не-
прерывного действия, т. е. объектов, работа которых описывается
обыкновенными дифференциальными уравнениями. Одним из таких
общих методов анализа и синтеза систем, содержащих цифровые
устройства, является метод частотных характеристик [Л. 12]. Что-
бы иметь возможность пользоваться этим методом, необходимо
сначала определить понятие частотной характеристики цифровой
модели.
Частотные характеристики цифровых моделей. Пусть имеется
некоторая цифровая модель, обладающая для простоты одним вход-
ным и одним выходным каналами. Предположим, что на вход этой
модели в момент времени t=0 подано число, равное единице,
а во все предыдущие и во все последующие моменты времени зна-
чения входных сигналов равны нулю. На выходе модели при этом
будет иметь место некоторая последовательность выходных сигна-
лов или чисел, действующих в моменты времени tn = nbt. Такую по-
84

следовательыость будем рассмагривать как решетчатую функцию
независим.ой переменной п.
По аналогии с непрерывными системами решетчатую функцию
h[n\ полученную на выходе цифровой модели при описанных выше
условиях, можно рассматривать как реакцию модели на воздействие
типа единичного импульса. Функцию h[n\ будем называть, следуя
Я. 3. Цыпкину [Л. 46], импульсной переходной функцией цифровой
модели. Продолжая аналогию дальше, определим частотную харак-
теристику цифровой модели как спектр ее импульсной переход{юй
функции S{h[n]}.
При таких условиях задача определения частотной характери-
стики сводится к задаче определения понятия спектра решетчатой
функции. При этом необходимо, чтобы такое определение правиль-
но учитывало реальные условия взаимодействия цифровой модели
и работающей совместно с ней непрерывной системы.
Будем считать, что каждому значению произвольной решетча-
той функции /[/г] ставится в соответствие бесконечно короткий им-
пульс, «площадь» которого равна /М- Другими словами, мы пред-
полагаем, что выполняется формальное равенство
00
nn\=Yi^[n]4t-n$t), (169)
—00
где 5 —— дельта-функция, равная единице в момент времени
n^t и равная нулю во все остальные моменты времени.
Теперь можно воспользоваться преобразованием Фурье и за-
писать:
со
5*{fm}=j ^ fMSC^-nSO^-'^'rf^. (170)
—00 /г=—00
На основании свойства линейности преобразования Фурье и тео-
ремы запаздывания получим непосредственно:
оо оо
S*{f[n]}= v f[n] j 8(/-nS,)^-''°'rf< =
/г=—ОО •—оо
00
= S П«]. (171)
П =—00
и, наконец, вводя подстановку
запишем окончательно:
00 _
5*{П«]}= S fWe-'""- (172)
«=—00
85

Из рассмотрения выражения (172) видно, что спектр решетча-
той функции представляет собой ряд Фурье, коэффициентами кото-
рого являются значения функции /[/г].
Заметим далее, что функция 5* |f [/ijj является периодической
функцией переменной со с периодом 2 = 2л. Отсюда можно сделать
вывод, что, для того чтобы получить полное представление о харак-
тере импульсной переходной функции, а следовательно, и о дина-
мических свойствах цифровой модели, достаточно знать участок
спектра этой функции на интервале шириной Q. Физически это
означает, что если на выходе цифровой модели установлен идеаль-
ный фильтр с шириной полосы пропускания, равной Q, то все ко-
личество информации, содержащейся в выходных сигналах модели,
будет передано через этот фильтр и поступит на вход соединенного
с моделью устройства.
С другой стороны, если представлять реальный выходной сиг-
нал цифровой модели в виде последовательности очень коротких
импульсов, как это было сделано выше, и воздействовать таким
сигналом непосредственно на вход непрерывной системы, то в этой
системе может возникнуть ряд процессов, обусловленных самой
формой представления выходного сигнала модели. В частности, это
может служить причиной неустойчивости системы в целом.
Отсюда следует вывод, что в качестве истинного выходного
сигнала цифровой модели следует принимать сигнал, получаемый
в результате фильтрации пocлeдoвaтevIьнocти импульсов вида (169)
идеальным фильтром с полосой пропускания Q. Такой сигнал дол-
жен содержать всю информацию о происходящих в модели процес-
сах и не содержать информации об условиях, принятых относитель-
но формы представления чисел, образующихся на ее выходе.
Спектр такого сигнала 5|/[/i]} будет совпадать с полученным
выше спектром в полосе частот от —до и равняться нулю
повсюду вне этого интервала. Будем называть спектр S {f [п]}
истинным спектром выходного сигнала цифровой модели.
Применим теперь к спектру S{f[n]} обратное преобразование
Фурье. В результате очевидных преобразований получим:
^ sin-^(/-nS,)
= У, f^"^— • (173)
Функция
86
7Z

стоящая в правой части выражения (173), используется при выводе
известной теоремы Котельникова [Л. 45] и получила название функ-
ции отсчета. Функция /к (О» определенная ио формуле (173), пред-
ставляет собой одну из огибающих исходной решетчатой функции.
Эта функция обладает тем свойством, что_ее спектр ограничен и
полностью заключен в пределах изменения со от —тс до я.
Функцию
00 sin^ (t — ndt)
n=—oo
будем называть истинной импульсной переходной функцией цифро-
вой модели, а спектр этой функции
S*{h,(t)} = S{h[n]} (175)
и будет являться искомой частотной характеристикой.
Если исходное разностное уравнение линейно, .мы можем при-
менить к нему дискретное преобразование Лапласа [Л. 47] и полу-
чить дискретную передаточную функцию модели в виде
U7*(,)^b2 . (176)
|х=0
Учитывая основные свойства дискретного преобразования Лап-
ласа, можно сделать вывод о том, что для получения спектра им-
пульсной переходной функции цифровой модели достаточно подста-
вить в выражение для ее дискретной передаточной функции /со
вместо q. Наконец, рассматривая этот спектр в пределах от —л;
до -fn, мы получим частотную характеристику цифровой модели,
описываемой передаточной функцией (176).
Практически интересующий диапазон частот заключен в пределах
О < (О < 7Г.
Таким образом, для того чтобы получить частотную характери-
стику цифровой модели по заданному разностному уравнению, до-
статочно получить ее дискретную передаточную функцию, заменить
в этой функции всюду q на ](»} п затем образовать выражение
S {/г[п]} = П(со)Т17*(/со), (177)
где
ПН = / ' приО<1со|<.;
\ О при i СО i ^ Л.
87

20. ЧАСТОТНЫЕ ОЦЕНКИ «КАЧЕСТВА»
КВАДРАТУРНЫХ ФОРМУЛ
Определив понятие частотной характеристики цифровой моде-
ли, мы можем оценивать теперь качество квадратурной формулы
по расхождению между «идеальной» частотной характеристикой,
или частотной характеристикой устройства, точно решающего
исходное дифференциальное уравнение, и частотной характеристи-
кой системы, работа которой описывается данной квадратурной
формулой. Применяя преобразование Фурье к уравнению (129), по-
лучим выражение для идеальной частотной характеристики в виде
5//z(0} = 4-. (178)
Таким образом, для каждой квадратурной формулы мол^ет быть
построена функция
k
—/cov
h ^ a^e-i"-
— h v=o
М")=дГ k ('79)
2 b^e-i-"-
V = —1
представляющая собой функцию погрешности в частотной обла-
сти. Кривые функции |8(а)) для некоторых квадратурных формул
из табл. 1 показаны на рис. 19.
Выражение (179) определяет комплексную функцию действи-
тельного аргумента са. Поэтому непосредственное использование
этого вырал^ения для оценки погрешности вычислений в отдельных
случаях может быть связано со значительными трудностями. Для
получения оценки, свободной от этого недостатка, мы воспользуем-
ся следующими соображениями.
Предположим, что исходное дифференциальное уравнение вида
(127) при некоторых заданных начальных условиях имеет решением
функцию у{х), непрерывную на некотором отрезке 0<д:<Л' и
имеющую на этом отрезке по меньщей мере k + 2 непрерывных
производных. Будем считать также, что нам известны для указан-
ного отрезка максимальные значения модулей этих производных
max I у' |, max inax | ^/1^ + ^^ |. Примерный вид функции у {х)
показан на рис. 20.
Предпотол<им теперь, что тем или иным способом было вычис-
лено п первых значений функции у(х,). Перед началом вычислений
(^-f!l)-ro шага на основании имеющихся данных люжно утверждать,
что значение ординаты в точке Хп+\ находится в пределах отрезка
цщриной А, показанного на рис. 20, причем А, очевидно, равно
А = 2/г max I i/'|. (180)
Если для вычисления используется одна из рассмотренных выше
формул численного интегрирования, то после выполнения вычисле-

дб
20
10
-20
L
(5 (Си
п
1
1
/
/
/
I
/
/
у
/
1
LffL
С
V
0,5
1
\
\
\
\
1
1
1
Метод прямоигольника I
JS/Z"
Метод ЗалераЕ
Метод опрапеций
Метод Симпсона З/в
Рис. 19. Частотные характеристики для
некоторых формул численного интегри-
рования.
>
——-
Рис. 20. Примерный вид функ-
ции у{х).
89

НИИ (/2+1)-го шага значение ординаты в точке Xn+i снова окажет-
ся заключенным в пределах некоторого отрезка Еп, такл<е показан-
ного на рис. 20. Ширина этого отрезка определяется суммарной по-
грешностью, получающейся при выполнении вычислений одного
шага. Учитывая только погрешность аппроксимации и имея в виду,
что все остальные составляющие могут быть сделаны сколь угодно
малыми путСхМ соответствующего выбора количества разрядов, в со-
ответствии с формулой (137), можно положить
En = 2Rh(^ + ') max|f/l^+2]|^
где R — постоянный коэффициент.
iB качестве оценки формулы численного интегрирования прини-
мается величина
_A__J__ hmax\y'\
Е "■■^/,(^ + 2)„iax|r/l^+2l|- (181)
Величина с может служить своеобразной оценкой качества квадра-
турной формулы, поскольку она показывает, до какой степени зна-
чение очередной ординаты оказалось уточненным в результате вы-
числений по сравнению с теми сведениями, которые имелись об этом
значении перед началом очередного шага вычислений.
Предположим теперь, что функция у{х) имеет преобразование
Фурье
у(х) ^уН>
Тогда для преобразования Фурье 1-й. производной функции у{х)
имеем непосредственно:
yV](x)^{icoyl(co).
Значения самой /-й производной функции у [х] могут быть получены
в результате применения обратного преобразования Фурье:
00
^1^1 (-^) i J ii^Yy («) ^^'^"^ dc,, (182)
—00
Будем считать, что спектр функции у{х) ограничен по оси частот.
Другими словами,
У (со) при |(о|<(Ос;']
а I Y((o) при |(о|<(Ос;]
^(")=\ О при |(0|>0)с.
где сос — граничная частота спектра. При этом выражение (182)
может быть переписано в виде
с
90
(/со)^У((о)^^'^^^а). (183)
1

Для оценки максимального значения 1-н производной функции
у(х) воспользуемся очевидным соотношением
max ul^iwk 2^
с
271
SCO
/+1
1(/со)ч \S\\e^'-^ld<o = щ^^^
с
5 = max 1 г(сдз)|.
(184)
Подставляя теперь соотношение (184) в выражение (181), получим
окончательно:
с^-2^(Лсос)-(^ + ^).
(185)
Формула (185) выражает оценку качества квадратурной фор-
мулы через значения порядка этой формулы граничной частоты
спектра решения Ос и шага разбиения оси независимой перемен-
ной h. Из рассмотрения выражения (185) видно, что в зависимо-
сти от конкретных значений /кос величина с может ли'бо возра-
стать, либо убывать при возрастании порядка квадратурной фор-
мулы k. В частности, увеличение точности при повышении поряд-
ка формулы k будет иметь место в том и только в том случае,
если удовлетворяется неравенство
/гсос<1.
(186)
Если /гсос^1, k = Oy то величинам, вычисленная из выражения
(185), не будет превышать 1. Ясно, что выполнять вычисления
по формуле численного интегрирования при таких условиях бес-
смысленно. Следовательно, при заданном h неравенство (186) мож-
но рассматривать как основное соотношение, ограничивающее класс
дифференциальных уравнений, которые могут решаться методами
численного интегрирования в реальном масштабе времени. Этот
класс включает в себя только такие уравнения, решения которых
имеют спектр, ограниченный по (186).
Обычно при конструировании цифровой модели величина h за-
дается, исходя из возможностей используемой аппаратуры, а вели-
чина 0с определяется характеристиками моделируемого объекта.
При таких условиях соотношение (il85) может быть использовано
для выбора порядка квадратурной формулы k.
Аналогичным образом может быть получено выражение для от-
носительной погрешности вычислений одного шага. Это выражение
имеет вид
R
(187)
91

Логарифмируя'выражение (187), получим для фиксированного h:
log i в к log + (^ + 2) log /гсос - log (k + 3).
Кривые изменения 1о^ |8| нри измененни граничной частоты
спектра сОс для значений /г^Ю-^ сек и k = 0, 1, 2, я показаны на
рис. 21. Каждая кривая на
^2
20
IS
этом рисунке ограничивает об-
ласть, определяющую возмож-
ности квадратурной формулы
данного порядка, .в смысле по-
лучаемой точности и полосы
пропускаемых частот.
Выведенные соотношения
(185) н (186) мо'гут быть при-
няты в основу при выборе од-
ной из известных формул чи-
сленного интегрирования, при
заданной точности и заданных
частотных свойствах модели-
руемого объекта. KpoiMe того,
они позволяют в общем виде
решать вопрос о применимости
данного цифрового вычисли-
тельного устройства для реше-
ния определенного класса за-
дач.
|Вьгведем теперь выраже-
ние для оценки погрешности
аппроксимации, получаемой при использовании формул численного
интегрирования, ооно!ванных на приближении в среднем. Перепи-
шем выражение ч'152) в виде
гл=гч
О
Рис. 21. Области возможных значений
быстродействия и точности для некото-
рых квадратурных формул.
(^/-(со) (СО) У] АЦ^^[^):
v=0
(188)
где £/* (о)) — преобразование Фурье функции г/[п], представляющей
собой приближенное решение исходного дифферен-
циального уравнения, полученное с помощью рассмат-
риваемых формул численного интегрирования.
Вычитая соотношение (188) из соотношения (153), получим:
(^^■"-1)е(со)=/г/(со)
0^"^ —\
v=0
(с)
(189)
где e(w) = (й>) — f/* ((о). Или, учитывая, что f (дг, у) = у\
k
(^^"' —1)£((о)=: hj(£ty (со)
_£^1
v=0
(со)
. (190
92

Будем предполагать, что спектр функции у{х) ограничен или,
другими словами.
Y (со) при i со к JJ-;
О при i со i > [x.
Применим теперь к выражению (190) обратное преобразование
Фурье. С учетом сделанного предположения имеем:
в(х + h) — e(x)= jay у (со)^
J(i>X
v = 0
da.
(191)
Разность в левой части выражения (191) представляет собой
приращение погрешности аппроксимации за один шаг вычислений и,
следовательно, имеет тот же смысл, что и величины, получаемые
выше. Для оценки максимального значения этой погрешности вос-
пользуемся неравенством Коши — Буняковского
E=\e(x+h)-^{x)\nai.o<Yi
X
doi
(192)
или в результате очевидных преобразований
АТ7
|/4Si(.)-4bipiI^:^
(193)
Выражение (193) и дает искомую оценку максимума погрешности
аппроксимации, получаемой за один шаг вычислений. Подставляя
выражение (193) в формулу (181), получим соответственно:
4 Si (jx) —4
1 — cos (J.
v=0
(194)
21. ИНФОРМАЦИОННАЯ ПРОИЗВОДИТЕЛЬНОСТЬ
При выводе оценки информационной производительности ре-
шающего блока цифровой модели, выполняющего операцию прибли-
женного интегрирования дифференциального уравнения первого по-
93

рядка по одной из формул численного интегрирования, рассмотрен-
ных выше, заметим прежде всего, что величина
А
/ = log2C = log2-^, (195)
где с определяется из выражений (185) или i(194), может рассма-
триваться как количество информации, перерабатываемой в течение
одного шага. Действительно, величина 1/А в выражении (195) опре-
деляет количество сведений о значении данной ординаты функции
у(х), имеющихся перед началом данного цикла вычислений. Вели-
чина \1е определяет количество сведений, которые будут иметься
о значении данной ординаты функции у{х) после того, как вычис-
ления данного шага уже произведены. При таких условиях выра-
жение (195) совладает по смыслу с определением количества инфор-
мации, введенным Хартли [Л. 45].
Информационной производительностью решающего блока циф-
ровой модели будем называть величину
9 = -^, (196)
где Т — количество времени, потребное для выполнения всех вычис-
лительных операций одного шага.
Поскольку оценки для величины с были нами уже получены
в предыдущей главе, задача получения оценки для информацион-
ной производительности сводится теперь к определению величины Т.
При этом будем исходить из следующих соображений. Пусть для
вычислений используется одна из квадратурных формул вида (163).
Будем предполагать, что все числа, участвующие в вычислениях,
представляются в двоичной системе счисления с фиксированной за-
пятой, а вычисление значений функции f{x, у) производится не
в данном решающем блоке, а в каком-либо другом устройстве
цифровой модели. При таких условиях вычисления по формуле сво-
дятся к выполнению операций суммирования и умножения на по-
стоянные числа. Предположим далее, что операция умножения вы-
полняется в форме суммирования сдвинутых частичных произведе-
ний, причем каждое такое частичное произведение получается
в результате умножения множимого на данный, отличный от нуля
разряд множителя с последующим сдвигом на требуемое количество
разрядов. Тогда вычисления по формуле (163) сводятся к сложе-
нию W чисел, где величину w можно определить из выражения
k ^
(.197)
v=0 v=0
где (6J и (aJ — количества отличных от нуля двоичных разрядов
в двоичном представлении коэффициентов и а^. Полученную
сумму необходимо умножить на величину h. Однако эту величину
можно всегда выбрать равной целой степени 2, и тогда операция
умножения сведется к операции сдвига.
94

Будем предполагать, наконец, что длительность одного шага
вычислений определяется возможностями арифметических элемен-
тов решающего блока, в то время как запоминающие элементы,
используемые для хранения чисел, над которыми выполняются опе-
рации, могут обладать значительно большим быстродействием.
Естественно, что перечисленные предположения довольно силь-
но ограничивают класс рассматриваемых структур. В общем виде
5
6)
0—^
б)
Рис. 22. Различные структуры арифметического
устройства цифровой модели.
задача синтеза структуры арифметического устройства, обладаю-
щего заданными временными характеристиками, может быть реше-
на средствами теории конечных автоматов.
(При рассмотрении различных структур арифметического
устройства мы будем исходить из следующих соображений. Пред-
положим прежде всего, что единственным элементом арифметиче-
ского устройства является одноразрядный двоичный функциональ-
ный сумматор, который мы будем характеризовать двумя величи-
нами. Первая из них ti представляет собой длительность переход-
ного процесса в схеме сумматора. Если на входы сумматора одно-
временно в некоторый момент поступают сигналы, имеющие вид
прямоугольных импульсов, то правильное значение результата бу-
дет действовать на обоих выходах сумматора, начиная с момента
^o+ti. Вторая 'величина т определяет длительность тактового интео-
1?ала. Если на входы сумматора в момент /р были поданы сигналы,
95

то следующая группа одновременных входных сигналов может быть
подана на входы того же сумматора, начиная с момента to+r.
Время X складывается из времени ti — времени, в течение которого
на выходе сумматора должно существовать правильное значение
результата, и длительности переходного процесса после окончания
действия входных сигналов.
Пусть на одном шаге вычислений необходимо сложить w чи-
сел. Будем считать, что обращение к любой ячейке запоминающего
устройства может производиться в любом порядке и в любой по-
следовательности. Заметим теперь, что в любом случае начинать
операцию сложения разрядов веса 2« можно только после того, как
получено правильное значение результата сложения разрядов веса
2^~К Таким образом, к моменту начала вычислений данного шага
мы имеем всего w двоичных разрядов, над которыми можно выпол-
нять независимые операции суммирования. Эти операции могут
быть выполнены либо с помощью .Q = wl2 (если один из входов
сумматора используется для сигнала переноса), либо с помощью
Q = wl3 (если все три входа сумматора свободны) независимых
сумматоров. Величину Q мы будем называть числом независимых
каналов вычислений. Однако можно показать, что значениям Q>1
в рассматриваемом случае соответствуют нерациональные структу-
ры арифметических устройств. Поэтому будем полагать в дальней-
шем, что Q=l.
Пусть в момент на входы сумматора Ец (рис. 23) поступают
младшие разряды двух чисел. Результаты операции над этими дву-
мя разрядами, а именно сигнал суммы и сигнал переноса, дейст-
вуют на выходе этого сумматора, начиная с момента to + Xi. В даль-
нейшем возможны следующие варианты:
1. Оба результата посылаются в оперативное запоминающее
устройство. В этом случае получается структура арифметического
устройства, состоящая из одного единственного сумматора. Дли-
тельность вычислений одного шага при этом равна:
T = N(w—\)z.
2. Сигнал переноса посылается в оперативное запоминающее
устройство, а сигнал суммы используется немедленно на входе сле-
дующего сумматора. При этом получается структура арифметиче-
ского устройства, показанная на рис. 22,а. Количество сумматоров,
включенных подобным образом, мы будем называть числом после-
довательных каналов вычислений F^. Очевидно, что в общем случае
величина Гз может изменяться от 1 до w—\. Длительность вычис-
лений одного шага при этом равна:
т = м'^^[^ + (Г,-\).,1
3. Сигнал суммы посылается в оперативное запоминающее
устройство, а сигнал переноса используется немедленно на входе
следующего сумматора. При этом получается структура арифмети-
ческого устройства, показанная на рис. 22,6. Количество суммато-
ров, включенных подобным образом, мы будем называть числом
параллельных каналов вычислений F2. Очевидно, что в общем слу-
96

чае величина F2 может изменяться от 1 до N, Длительность вычис-
лений одного шага при этом равна:
4. Оба выходных сигнала немедленно используются на входах
следующих сумматоров. При этом получается структура арифмеги-
чсского устройства, показанная на рис. 22,е. Общее количество сум-
маторов в схеме равно 8 = р2ръ- Длительность вычислений одного
шага при этом равна:
N w—X
T-J^ -77- Ь + {F2 + т,]. (198)
Формула (198) является общей формулой для расчета дли-
тельности вычислений одного шага.
Окончательно оценка информационной производительности бу-
дет получена, если подставить в формулу (196) значения w и Т,
вычисленные по формулам (197) и (198).
22. ОТНОСИТЕЛЬНАЯ СЛОЖНОСТЬ
В (результате проведенного выше анализа мы получили количе-
ственную оценку информациониой производительности цифровой
системы. В конечном итоге информационная производительность
является мерой того полезного эффекта, который мы получаем при
работе всякой вычислительной машины. Здесь также вполне умест-
на аналогия с машинами — преобразователями энергии, мерой вы-
ходного эффекта которых служит мощность, т. е. количество энер-
гии, преобразуемой в единицу времени. Теперь необходимо рассмо-
треть вопрос о затратах, ценой которых достигается этот полезный
эффект. К этим затратам относятся затраты энергии для питания
машины, затраты труда на ее обслуживание и ремонт и т. п,
Большую роль при оценке затрат играет надежность, определяю-
щая как полезное рабочее время машины, так и затраты на оты-
скание и устранение неисправностей. Можно утверждать, что все
эти факторы так или иначе зависят от количества использованных
в машине элементов или от ее сложности. В настоящее время от-
сутствует возможность проанализировать -конкретный вид всех этих
зависимостей и потому количественная оценка затрат, связанных
с эксплуатацией, встречает значительные трудности.
В то же время из сказанного можно сделать вывод, что до-
статочно хорошей приближенной оценкой всей суммы затрат, це-
ной которых достигается данная величина информационной произво-
дительности, может служить так называемая относительная слож-
ность машины, представляющая собой некоторую функцию от ко-
личеств использованных в ней элементов разных типов. Определим
относительную слолсность X как
р
X=Y^ainu (199)
/=;1
7 А. В. Шилейко. 97

где rii — количество элементов гипа г, ai— весовой коэффициент,
соответствующий элементам данного типа. 'Конкретные значения ве-
совых 'коэффициентов можно получить 'на основе статистического
анализа данных об эксплуатации машин одного и того же вида.
Один из примеров определения этих коэффициентов для достаточ-
но простого случая приводится в работе [Л. 48].
В пользу того, что за меру относительной сложности берется
именно линейная комбинация количеств элементов, говорит следую-
щее, хотя 'И не строгое, но достаточно убедительное рассуждение.
Пусть имеются две одинаковые вычислительные машины, работаю-
щие независимо и перерабатывающие информацию, поступающую
от различных независимых источников по одному и тому же алго-
ритму. Суммарная информационная производительность этих ма-
шин, очевидно, будет вдвое больше информационной производитель-
ности каждой отдельно взятой машины. Общая сумма затрат на
две машины будет также вдвое больше затрат на одну машину.
Следовательно, при одинаковых элементах и одинаковых 'условиях
эксплуатации относительная сложность, определенная по (199),
инвариантна по отношению к структуре вычислительной системы.
23. ЭФФЕКТИВНОСТЬ
Все соображения, развитые в этой главе, приводят непосред-
ственно к определению эффективности ri цифровой системы как
отношения информационной производительности к относительной
сложности
■П = -Т' (200)
Такое определение полностью соответствует общему определению
эффективности, используемому в современной теории исследования
операций, где эффективность определяется как отношение получае-
мой пользы к понесенным затратам. Частным случаем эффективно-
сти является коэффициент полезного действия машин — преобразо-
вателей энергии. Понятие эффективности даст нам возможность
проводить независимую оценку различных структур цифровых мо-
делей. Этому вопросу и посвящена следующая глава.
ГЛАВА ШЕСТАЯ
СИНТЕЗ СТРУКТУР ЦИФРОВЫХ МОДЕЛЕЙ
24. ИСХОДНЫЕ СООБРАЖЕНИЯ
Задача разработки цифровой модели в общем случае сводится
к составлению схемы соединений различных элементов или цифро-
вой структуры, соответствующей основному назначению и удовлет-
воряющей при этом определенным техническим условиям. Техни-
ческие условия формулируются в виде системы ограничений, накла-
дываемых на основные характеристики модели, рассмотренные
в предыдущей главе. Говоря другими словами, задание технически^
98

условий эквийалентио выделению в иространсгве характеристик
ограниченной области.
Структура цифровой модели состоит из конечного числа эле-
ментов и, следовательно, имеет конечное число различных вариан-
тов соединения этих элементов между собой. Каждой структуре
можно поставить в соответствие число или ряд чисел, которые
в силу конечности общего количества различных структур будут
однозначно определять каждую данную структуру. Такие числа мы
будем называть показателями структуры. С другой стороны, из ма-
териала предыдущих глав следует, что характеристики цифровой
модели зависят от конкретного вида структуры, т. е. являются
функциями ее показателей. Эти функции совместно с системой огра-
ничений, представляющей собой технические условия, определяют
в пространстве показателей некоторую ограниченную область. Дру-
гими словами, если выбран ряд показателей структуры fi, f 2,
то всегда может быть составлена система неравенств:
h(Fu Л)<Ф\-; ) (201)
р,(Л. F,) = Vi, ]
определяющая область допустимых значений параметров. Сами
функции у, ф, р в неравенствах (201) отображают ограниченную об-
ласть в пространстве характеристик на соответствующую область
в пространстве показателей. Функции (р, р могут быть получены
на основе соотношений, рассмотренных в предыдущей главе.
|Внутри рассмотренной области значения показателей структу-
ры могут изменяться. Отсюда следует вывод, что одним и тем же
техническим условиям, вообще говоря, будет удовлетворять несколь-
ко различных структур, причем число их будет тем больше, чем
шире область допустимых изменений параметров. Если каждой
из таких структур поставить в соответствие величину, принимаемую
за критерий качества, то может быть поставлена задача выбора
отруктуры, наилучшей в Схмысле этого критерия.
За критерий качества структуры цифровой модели мы будем
принимать ее эффективность, определяемую по (200). Это дает воз-
можность сформулировать в общем виде задачу синтеза структуры
цифровой модели следующим образом: ^1айти совокупность значе-
ний показателей структуры, удовлетворяющих системе ограничений
(201) и дающих при этом максимальное значение эффективносги.
Подобная задача находится в полном соответствии с общей задачей
синтеза оптимальных систем, рассматриваемой в рамках теории
исследования операций. Она может быть решена развиваемыми
в настоящее время методами линейного и в общем случае динами-
ческого программирования. Кроме того, сказанное выше может слу-
жить отправной точкой для создания единой структурной теории
цифровых моделей и вообще специализированных цифровых вы-
числительных машин. В настоящее время подобная теория только
начинает разрабатываться. Описываемая ниже методика синтеза
структуры цифровой модели представляет собой одну из первых
попыток в этом направлении.
7* 99

25. ПОКАЗАТЕЛИ СТРУКТУРЫ
Здесь мы ограничимся рассмотрением задачи синтеза структуры
решающего блока цифровой модели, выполняющего вычисления
по одной из формул численного интегрирования, описанных в гл. 4.
Некоторые рекомендации по выбору конкретного вида такой фор-
мулы будут приведены ниже.
Будем считать, что в техническом задании на конструирование
решающего блока оговорены значения граничной частоты спектра
выходного сигнала сос и точности, с которой должны выполняться
вычисления одного шага с. Тогда на основании формул (185) или
(194) может быть получена зависимость между значениями шага
интегрирования h и порядка формулы численного интегрирования k.
Величины h \1 k мы будем принимать в качестве основных показа-
телей структуры решающего блока.
В дальнейшем мы будем предполагать, что все числа, над ко-
торыми выполняются операции, представлены в двоичной системе
счисления. Количество разрядов в этих числах выбирается по сле-
дующим соображениям. Из рассмотрения разностного уравнения
для погрешности (168) можно сделать вывод, что характер накоп-
ления ошибки при интегрировании дифференциального уравнения
в равной степени определяется как погрешностью аппроксимации
Е{Хп), так н погрешностью округления Т{Хп). Если погрешность
округления оказывается значительно ниже погрешности аппрокси-
мации, то это не увеличивает общей точности решения, но влечет
за собой неоправданное увеличение либо количества аппаратуры,
либо длительности 'вычислений. В противоположном случае не бу-
дут полностью использованы возможности данной формулы числен-
ного интегрирования. Следовательно, целесообразно выбирать
Т(х)^Е(х). (202)
В свою очередь погрешность округления определяется количеством
разрядов в числах, над которыми выполняются операции, а именно
T(jc)<:W'2-^^-^^\ (203)
где — количество двоичных разрядов;
W — общее количество чисел, участвующих в одной операции.
Отсюда
N-=\-\og,^. (204)
Как видно из выражений (185), (193) н (197), значения величин
Е{х) WW определяются порядком формулы численного интегрирова-
ния k. Величину N мы также будем принимать в качестве показа-
теля структуры решающего блока. Из сказанного выше следует, что
количество разрядов представляет собой некоторую функцию
величины k.
Для хранения чисел, участвующих в операциях одного шага,
а именно чисел y[n—v] и )[п—v] {см. (163)], в состав решающего
блока должно входить местное запоминающее устройство, Полная
100

емкость этого запоминающего устройства Q также является показа-
телем структуры. В соответствии со сказанным выше она равна
wN двоичным ячейкам:
Q = wN. (205)
В предыдущей главе было показано, что структура арифмети-
ческого устройства решающего блока при сделанных предположе-
ниях полностью определяется двумя показателями: числом парал-
лельных каналов вычислений F2 и числом последовательных кана-
лов вычислений Fs. Длительность вычислений одного шага Т опре-
деляется по формуле ,(!198) через значения показателей F2 и Fs, ве-
личин N и W, рассмотренных выше, и значения постоянных т и Ть
представляющих собой характеристики используемых при конструи-
ровании элементов.
Будем предполагать, что управление последовательностью вы-
полнения отдельных операций в решающем блоке осуществляется
поступаюшт1ми извне управляющими сигналами, а количество эле-
ментов, осуществляющих функции управления, мало по сравнению
с общим количеством элементов^ занятых в блоке. При таких усло-
виях и сделанных выше предположениях показатели /г, k, N, Q, F2
и Fs полностью определяют структуру решающего блока. Другими
словами, при заданных значениях сос, т, ti и с каждой совокуп-
ности конкретных значений показателей соответствует единственная
схема соединений элементов решающего блока.
26. МЕТОДИКА СИНТЕЗА СТРУКТУРЫ
РЕШАЮЩЕГО БЛОКА
Задача синтеза наилучшей структуры решающего блока сводит-
ся к нахождению совокупности значений показате^тей структуры /г,
k, N, Q, F2 и ^3, при которых iB пределах области возможных изме-
нений этих параметров величина эффективности г\ будет достигать
своего максимального значения. Как уже отмечалось выше, эта за-
дача может быть решена методами линейного программирования.
При построении функциональной зависимости величины т], опре-
деленной по формуле (200), от значений показателей структуры
необходимо учитывать следующее:
il. Величины N я Q полностью определяются значениями huh
в соответствии с выражениями (204) и (205).
2. 'Величина с полностью определяется значениями ^ и Я и за-
данным значением (Ос в соответствии с выражениями (165) или
(194).
3. Величина длительности одного шага вычислений Т полностью
определяется значениями Л^, F2, /^з, ^ и заданными значениями t
и Ть Таким образом, функциональная зависимость величины 6, вхо-
дящей в состав выражения (200), от значений показателей струк-
туры может быть получена на основании выражений (195), (196),
(204),,(205) и (185) или (194).
4. Функциональная зависимость величины А,, также входящей
в состав вырал<ения (200), от значений показателей структуры
в рассматриваемом случае имеет вид
X = a,Q + a2F2F,, (206)
101

гДе ai и U2 — Постоянные коэффициенты, определяющие относитель-
ные «веса» запоминающих элементов и одноразряд-
ных сумматоров.
Из сказанного можно сделать вывод, что в рассматриваемом
случае определения структуры решающего блока цифровой модели
имеются в наличии все необходимые данные для построения функ-
циональной завнсимости
ri = G(h,kyNyQ,F,,F,). (207)
Рассмотрим теперь условия, ограничивающие область возмож-
ных изменений показателей структуры. В том случае, когда неза-
висимой переменной задачи является время, необходимо, чтобы
удовлетвор^ялось условие /г < Г.
Полагая h = T и подставляя выражение (198) в одну из формул
(185) или (194), получим первое основное ограничение
(iV, w, т, F,, F„ k, o)c) > c, (208)
где T, Tj и coc — заданные параметры; и w — как уже отмеча-
лось выше, представляют собой известные функции величины k.
Второе основное ограничение может быть получено в том слу-
чае, если 'В техническом задании оговорены максимальные значения
габаритов или веса решающего блока. Это ограничение имеет вид
g2 = M + hF2Fz<B, (209)
где постоянные Pi и в данном случае имеют смысл удельных
весов или удельных объемов выбранных элементов, а величина
В—максимально допустимые вес или объем.
Кроме двух рассмотренных основных ограничений, можно
ввести также ряд дополнительных неравенств. Легко, например,
показать справедливость неравенства
F,F^<Nw. (210)
Предлагаемая методика синтеза наилучшей структуры решаю-
щего блока цифровой модели предполагает, таким образом, выпол-
нение следующего ряда этапов:
1. Выбор типа формулы численного интегрирования на основа-
нии имеющихся в наличии сведений о характере изменения соответ-
ствующих переменных. Здесь могут быть предложены следующие
рекомендации:
а) Как известно [Л. 4], погрешность аппроксимации квадратур-
ной формулы экстраполяционного типа порядка k равна нулю
^в том случае, если функция f(y, х) представляет собой алгебраиче-
*ский полином порядка не выше ^—1. Поэтому применение квадра-
турных формул экстраполяционного типа можно рекомендовать во
всех случаях, когда закон изменения входных переменных решаю-
щего блока хорошо аппроксимируется алгебраическим полиномом
невысокого порядка (^^4—6). Эти же формулы выгодны в тех
случаях, когда граничная частота спектра функции у(х) относи-
тельно невелика (порядка сотых долей герца и ниже).
б) В случаях, когда граничная частота спектра функции у(х)
велика, выгодно применять формулы численного интегрирования,
основанные на приближении в среднем.
102

в) При необходимости непрерывного контроля за накоплением
погрешности (такая необходимость возникает, иапример, когда
функция f(y, х) испытывает резкие скачки), целесообразно исполь-
зовать формулы Мил1на.
Существует также ряд дополнительных соображений, которые
можно найти в специальной литературе.
2. На основании имеющегося технического задания строится
2/1
г)
Рис. 23. Различные структуры цифрового интегра-
тора.
функциональная зависимость вида (207) и система неравенства
вида (208), (210).
3. Методом линейного ирограммирования или любым другим
подходящим методом отыскивается совокупность значений показа-
телей структуры, обеспечивающая максимальное значение величины
эффективности т] в заданной области допустимых изменений пока-
зателей структуры. Этот этап обычно требует большого объема вы-
числений, и при его выполнении целесообразно использовать цифро-
вую .вычислительную машину.
В [Л. 27] был приведен пример построения структуры цифро-
вого интегратора с использованием описанной методики. Расчеты
проводились применительно к случаю использования квадратурных
103

формул экстраполяционного типа для 0<^<3. Результаты расче«
тов сведены в табл. 3.
Из таблицы видно, что максимальная эффективность достигает-
ся при значениях k=\, F2=\, Гз=\, т. е. наилучшей оказывается
схема цифрового интегратора последовательного типа, работающе-
го по квадратурной формуле типа Эйлера И (см. табл. 1). Блок-
схемы цифровых интеграторов для различных значений k показаны
на рис. 23.
Таблица 3
k
f2
w
Г
5
X
0
1
1
2
43
1
93
64
0
43
1
2
15
43
387
45
1
1
1
2
34
1
109
2 700
1
34
1
2
12
34
340
2 500
2
1
1
3
66
1
205
967
2
33
2
3
23
33
429
1 330
2
1
2
3
22
2
212
1 390
2
33
2
3
12
66
658
1 650
3
i
1
4
93
1
224
770
3
31
1
4
33
31
434
1 092
3
1
3
4
52
3
238
1 270
3
31
3
4
12
93
868
1 535
С увеличением числа последовательных и параллельных каналов
значение эффективности уменьшается при k==0 и k=\ и увеличи-
вается при k^2. Значения относительной сложности X при k=\
имеют тот же порядок, что при k = 0. Следовательно, даже при на-
личии жестких ограничений на сложность (габариты, вес) следует
выбирать структуру, характеризуемую значением k=\, обладающую
при этом наивысшей эффективностью. При наличии ограничений по
быстродействию можно переходить к методам более высоких поряд-
ков, хотя получаемый при этом выигрыш, как видно из табл. 3,
относительно невелик.
В качестве второго примера рассмотрим применение описанной
методики к структуре цифрового дифференциального анализатора
(см. гл. 2). В этом" случае алгоритм является заданным. Для оцен-
ки информационной производительности используется то обстоя-
тельство, что .выходной сигнал каждого интегратора ЦДА коди-
руется по способу РДМ при трех возможных значениях для каж-
дого сигнала (—^1, О, 1). Выражение для пропускной способности
такого канала было получено В. М. Байковским [Л. 3]. Оно имеет
вид
1 / ти \
0 = — log2 ^1 + 2^cos ^ _|_ i j . дв. ед/сек, (211)
где Л —число уровней квантования сигнала, которое в данном слу-
чае равно
Л = 2^-2; (212)
104

N — полное количество разрядов в регистрах Y и R (см. гл. 2);
Т — интервал времени между двумя последовательными выход-
ными сигналами или длительность одной итерации.
Для достаточно больших значений А можно приближенно счи-
тать
o^yriog^a, (213)
и в окончательном виде выражение для критерия качества ЦДА,
состоящего из М независимых интеграторов, может быть записано
следующим образом:
Л4е М loggS
^ = -r=-jr-' (2^4)
Операции, выполняемые в каждом интеграторе, сводятся
к сложению двух чисел на каждом шаге. Всего на каждо-м шаге
независимо складывается 2М чисел, где М — общее количество
интеграторов. Емкость олеративного запоминающего устройства
равна 2MN, где — полное количество разрядов в числах, над
которыми выполняются операции. Таким образом, единственным
изменяемым показателем здесь является показатель Fz или число
параллельных каналов в арифметическом устройстве.
При расчетах за единицу количества элементов принимался
условный элемент, способный хранить одну двоичную единицу
информации (например, триггер) или перерабатывать одну двоич-
ную единицу информации за один такт. Величина коэффициента а
для такого элемента принималась равной единице. Рассматривались
также запоминающие устройства двух различных типов. К первому
типу были отнесены устройства, использующие для хранения одной
двоичной единицы информации один реальный физический элемент
(магнитный сердечник, электронный триггер и т. д.). Относительная
сложность такого устройства оценивалась количеством элементов,
т, е. емкостью.
К запоминающим устройствам второго типа были отнесены
устройства типа линий задержки (магнитный барабан, магнито-
стрикционные линрш задержки и т. д.). При этом для относитель-
ной сложности была принята оценка в 10 элементов на один канал
хранения информации (одна дорожка магнитного барабана).
При расчете быстродействия предполагалось, что элементы
арифметического устройства и запоминающих устройств первого ти-
па могут работать при тактовой частоте "1 Мгц. Элементы запоми-
нающих устройств второго типа могут работать при тактовой ча-
стоте порядка 100 кгц. Задержка, вносимая элементами арифмети-
ческого устройства, предполагалась равной ti = 0,l мксек.
Результаты расчетов, проведенных при указанных исходных со-
ображениях, представлены в виде графиков на рис. 24. Эти графики
построены для случая ЦДА, состоящего из 60 интеграторов
(M=60), при количестве разрядов в числах, равном 20 (Л^=20).
Кривая / соответствует запоминающим устройствам первого типа,
а кривая 2 — запоминающим устройствам второго типа. На кривых
отмечены точки, соответствующие чисто последовательной, после-
довательно-параллельной, параллельно-последовательной и чисто
параллельной структурам. Из рис. 24 видно, что при использовании
105

запоминающих устройств первого типа наилучшей оказывается чи-
сто параллельная структура, а при использовании запоминающих
устройств второго типа — параллельно-последовательная.
Два рассмотренных простых примера ни в коей мере не исчер-
пывают проблемы. Они могут только служить указанием на важ-
ность постановки самой задачи синтеза оптимальной структуры
цифровой модели, поскольку в обоих случаях эффективность резко
3000
2000
ж
600
400
200
ЮО
60
40
20
10
6
3
Z
'А
2
Z3 4 6 810 20
Рис. 24. Кривые зависимости эффективности
от числа параллельных каналов вычислений.
меняется при изменении тех или иных параметров структуры. До-
стигаемый при этом .выигрыш в эффективности будет тем больше,
чем уже класс задач, предполагаемых к решению на модели.
Второй вывод, который можно сделать нз приведенного выше
рассмотрения, состоит в том, что при синтезе структуры цифровой
модели наряду с чисто конструктивными параметрами необходимо
учитывать также и показатели алгоритма. Например, численный
метод нулевого порядка, или так называемый метод прямоугольни-
ков, достаточно часто используемый в ЦДА, оказывается наименее
целесообразным, так как, не будучи менее сложным и не давая
выигрыша в быстродействии, он приводит к резкому снижению
эффективности по сравнению с другими численньши методами.
Совместно с существующими методами динамического програм-
мирования описанная методика синтеза оптимальной структуры
цифровой модели дает возможность запрограммировать процесс
синтеза для универсальной цифровой вычислительной машины и,
следовательно, автоматизировать один из этапов разработки модели.
Поскольку в выражение для эффективности входит оценка относи-
тельной сложности, зависящая в свою очередь от оценок относи-
тельной сложности отдельных элементов, успех работы в большой
степени будет зависеть от наличия конкретных данных по соответ-
ствующим элементам.
106

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

ЛИТЕРАТУРА
1. Цифровые дифференциальные анализаторы, об. переводов п'од
ред. Б. Я. 'Когана, Изд. ино,сгранн'ОЙ литературы, 1959.
2. Б а 3 и л е в с к и й Ю. Я- и Ш р е й д е р Ю. А., Методы оцен-
ки про'изводитатьности уииверсальных цифровых вычислительных
машин 1с программным управлением, сб. «Вопросы теории матема-
тических машин», Физ1матиздат, '1958.
3. Банковский В. М., О пропускной способности канала
с разностно-дискретной модуляцией, «Радиотехника», 1961, т. 16,
№ 7.
4. Б ере ЗИП И. С п Жидков Н. П., Методы вычислений,
Фпзматиздат, 1959.
5. Boxer iR., 'Analysis of sampled-data systems land digital com-
puters in the frequency domain, IiRE Convention iRecord, p. 10, 1955.
6. Boxer R., Thaler S., A simplified method of solving li-
near an-d nonlinear systems, Proc. LRE, -1956, № 4.
7. Б p и к В. A., Цифровое вычислительное ;устройство для со-
ставления программы обработки деталей на фрезерном станке, сб.
«Автоматическое управление п регулирование», Изд. АН СОСР, 1962.
8. Бут Э. Д., Численные методы, Физматиздат, 1959.
9. iBopoHOB А. А. п др.. Цифровые аналоги для систем авто-
'Матичеокого управления. Изд. АН ОСОР, 1960.
10.'Г ель фонд А. О., Исчисление конечных разностей, Физ-
матиздат, 1959.
И. Genichi Н. and А о к 1 К., Interpretive differential equation
analyzer (IDEA) on USSC, 'Codes iReactor Computat., Vienna,
1961.
12. Гетманов A. Г., Динамические особенности реализации
некоторых линейных операций иа цифровой вычислительной машине,
сб. «Автоматическое (управление н вычислительная техника», вьвп.З,
Машгиз, 1960.
13. 'Г ИТ'И с Э. И., Преобразователи информации для электрон-
ных цифровых вычислительных устройств, Госэнергоиздат, 1961.
14. G г е е п s t е i п J. L., Leger R. M., Simulate digitally, or
by combining analog and digital computing facilities, Control Engi-
neering, Sept., 1956.
15. Г у p a к 0 в A. A. и Шевелев A. Г., Цифровой интегратор
с троичной системой кодирования прираи1еиий, сб. «Комбинирован-
ные вычислительные машины», Изд. АН СССР, 1962.
108

16. Дроздов Е. А., Цифровые аналоги, «Приборостроение»,
1957, № 5.
17. Жидков Е. П., Математические основы цифрового диффе-
ренциального анализатора, сб. «Автоматическое управление и вычи-
слительная техника», вып. -1, Машгиз, 1958.
18. Harris J. N., А programmed variable-rate counter for ge-
nerating the sine-function. Trans. IRE, PGEC EC-5, № 1, March,
1956.
19. Каляев A. В., О путях повышения быстродействия и
о расширении логических возможностей цифровых дифференциаль-
ных анализаторов, сб. «Комбинированные вычислительные .машины».
Изд. АН СССР, 1962.
20. К а 1 m а п R. Е., iNonlinear aspects of sampled-data control
systems, Proceedings of the symposium on nonlinear circuit analysis,
Politechnic Ins.titute of Broocklyn, April 25, 26, 27, 1956.
21. К a p и б СКИ Й В. В., Специализированное вычислительное
устройство для задания движения объекта по прямой, параболе и
окружности, сб. «Автоматическое управление п регулиронанне».
Изд. АН СССР, 1962.
22. Китов А. И. и К'риницкий Н. А., Электронные цифро-
вые машины и программирование, Физматиздат, 1959.
23. Кобринский Н. Е. и Т р а х т е н б р о т Е. А., Введение
в теорию конечных автоматов, Физматиздат, 1962.
24. К о г а н Б. Я., Электронные моделирующие устройства и их
применение для исследования систем автоматического (регулирования,
Физматиздат, (1959.
25. Коган Б. Я., О принципах построения комбинированных
аналого-цифровых вычислительных машин, сб. «аКомбинированные
вычис^тительные машины». Изд. АН COOP, 1962.
26. К о ж а р с к и й Л. А., Об одном способе построения цифро-
вых дифференциальных анализаторов, сб. «Вопросы теории матема-
тических машин», Физматиздат, 1958.
27. Козырева Г. М. и Шилейко А. В., О структурах спе-
циализированных цифровых вычислительных машин, сб. «Комбини-
рованные (Вычислительные машины», :Изд. АН ОООР, 1962.
28. К о л л а т ц Л., Численные методы решения дифференциаль-
ных уравнений. Изд. иностранной литературы, 1953.
29. М а й о р о в Ф. В., Электронные цифровые интегрирующие
машины, Машгиз, 1962.
30. М а й о р о в Ф. В., Принципы построения и применение вы-
числительных машин, работающих на приращениях, сб. «(Комбини-
рованные 'Вычислительные 1мапгины», Р1зд. АН ОООР, *1962.
31. Матюхин Н. Я., Дискретные линейные фильтры. Известия
высших учебных заведений. Радиофизика, т. Л, 1959.
32. Merz D. М., GEVIO—а real time variable increment digi-
tal computer. Automatic Control, v. ill, № 3, September, '1959.
33. Meyer M. A., iDigital techniques in analog systems. Trans.
IRE, PGEiC :E0-3, № 2, June, 1954.
34. Милн В. Э., Численное решение дифференциальных урав-
нений, Изд. иностранной литературы, 1955.
35. М о п г о е А. J., Digital processes for sampled-data systems,
John Willey and iSons, iNew York, 1962.
36. M u г r a у D. v., A variable reate binary scaler, Trans. IRE,
POEC EC-4, № 2, June, 1955.
109

37. Мухин И. С, К накоплению ошибок нрн численном инте-
грировании уравнений, «Прикладная математика н механика», 1952,
т. XVI, вып. 6.
38. iM i п о г U и г а b е. Theory of errors in numerical integration
of ordinary differential equations. Journal of Science of the Hiroshima
University, series A-1, v. 25, n9 1, July, 1961.
39. H e с л у X 0 в с к и Й К. С, Некоторые вопросы построения
цифровых дифференциальных анализаторов последовательного типа,
об. «Комбинированные вычислительные машины», Изд. АН СССР,
11962.
40. Н е с л у X о в 'С к и й К. С, Цифровые дифференциальные
анализаторы, Физматиздат, 1963.
41. Полляк Ю. Г., Об экстраполяции нредсказаиия функции
с ограниченным спектром но ее дискретным значениям, «Радиотех-
ника», 1962, т. 47, № 11.
42. Ракитский Ю. В., О некоторых свойствах решений си-
стем обьжновенных дифференциальных уравнений одношаговыми
методами чишенного интегрирования. Журнал вычислительной ма-
тематики и математической физики, 1961, т. I, № 6.
43. S а 1 Z е г J. М., iFrequency analysis of digital computers ope-
rating in real (time, Proc. of the iIiRE, 1954.
44. Тетельбаум И. М., Электрическое моделирование, Физ-
матиздат, 11959.
45. X а р к е в и ч А. А., Очерки общей теории связи, Гостехиз-
дат, 1955.
46. Цыпки н Я. 3., Теория импульсных систем, Физматиздат,
1958.
47. Шил ей ко А. В., Цифровые модели (обзор), «Автоматика
и телемеханика», 1959, т. iXX, OSfe /12.
48. Шил ей к о А. 'В., Методика выбора оптимальной структуры
цифровой модели, «^Автоматика и телемеханика», 1,1961, т. XXII, № Ь
49. Ш и л е й к о А. В., Цифровые дифференциальные анализато-
ры, изд. ВИНИТИ, 1961.
50. Шил ей к о А. В., Об одном классе оценок качества фор-
мул численного интегрирования дифференциальных уравнений, сб.
«Вычислительная техника в управлении», Изд. АН СССР, 1963.
51. Шур а-Бур а М. Р., Оценки ошибок численного интегри-
рования обыкновенных дифференциальных уравнений, «Прикладная
математика и механика», 4952, т. XVI, вып. 5.
52. Western S/c ODA Council meeting of 5 May on problems
unique to computers of a kind, Instrum. and Automation, v. 31,
July, 4958.
53. Western S/c joint meeting of 8 January in interrelationship
of analog and digital computers, Instrum. and Automation, v. 31,
№ 3, 1958.

СОДЕРЖАНИЕ
Введение . . . • • 3
Глава первая. Цифровое моделирование 7
1. Исходные соображения 7
2. Аналоговые и цифровые модели 8
3. Задачи цифрового моделирования И
Глава вторая. Цифровые дифференциальные анали-
заторы • 1^
4. Исходные соображения 14
5. Принцип действия цифрового интегратора 15
6. Работа цифрового интегратора в режиме решения обык-
новенных дифференциальных уравнений 22
7. Структуры ЦДА 24
8. Работа цифрового интегратора в специальных режимах 27
9. Основные характеристики ЦДА 31
Глава третья. Цифровые модели с переменными при-
ращениями 44
10. Исходные соображения 44
11. Машина GEVIG 45
12. Машина для решения дифференциальных уравнений
интерполяционным методом 53
Глава четвертая. Алгоритмы цифровых моделей . . 55
13. Исходные соображения 55
14. Разностные схемы • . . . . 58
15. Квадратурные формулы 67
16. Квадратурные формулы, основанные на приближении в
среднем тригонометрическими рядами 75
17. Способы вычисления функций 78
Глава пятая. Характеристики цифровых моделей ... 79
18. Исходные соображения 79
19. Точность 81
20. Частотные оценки „качества" квадратурных формул . . 88
111

21. Информационная производительность 93
22. Относительная сложность 97
23. Эффективность • 98
Глава шестая. Синтез структур цифровых моделей. 98
24. Исходные соображения 98
25. Показатели структуры
26. Методика синтеза структуры решающего блока .... 101
Литература 108