Текст
                    В.В.Амелькин
ДИФФЕРЕНЦИАЛЬНЫЕ
УРАВНЕНИЯ
В ПРИЛОЖЕНИЯХ


В.В.Амелькин ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ В ПРИЛОЖЕНИЯХ МОСКВА «НАУКА» ГЛАВНАЯ РЕДАКЦИЯ ФИЗИКО-МАТЕМАТИЧЕСКОЙ ЛИТЕРАТУРЫ 19 87
ББК 22.161.6 А61 УДК 517.91 @23) А м е л ь к и н В. В. Дифференциальные уравнения в приложе- приложениях.— М.: Наука. Главная редакция физико-математической литера- литературы, 1987.—160 с. Книга популярно знакомит с возможностями использования обык- обыкновенных дифференциальных уравнений при изучении реальных явле- явлений и процессов. Приемы составления дифференциальных уравнений, а также некоторые методы их качественного исследования иллюстри- иллюстрируются задачами, возникающими в различных областях знаний. Для школьников старших классов, преподавателей, студентов, для специалистов нематематических профессий, использующих мате- математику в своей работе. Табл. 15. Ил. 97. Библиогр. 8 назв. Рецензент доктор физико-математических наук М. В. Федорюк 1702050000-031 по „ А —Псо /аО\ Q7 38"87 © Издательство «Наука». UOo [yZ)-o/ Главная редакция физико-математической литературы, 1987
СОДЕРЖАНИЕ Предисловие . , , ПОСТРОЕНИЕ ДИФФЕРЕНЦИАЛЬНЫХ МОДЕЛЕЙ И ИХ РЕШЕНИЙ Чей кофе более горячий? . . . . 9 Стационарный тепловой поток 11 Случай в заповеднике 13 Истечение жидкости из сосудов. Водяные часы 18 Эффективность рекламы 20 Спрос и предложение 21 Химические реакции . 22 Дифференциальные модели в экологии .... 24 Одна задача математической теории эпидемий 28 Кривая погони 33 Модели боевых действий 35 Почему маятниковые часы не являются точными? 44 Циклоидальные часы 46 Задача о брахистохроне , 51 Среднее арифметическое, среднее геометрическое и дифферен- дифференциальное уравнение . 54 О полете тела, брошенного под углом к горизонту . . 57 Невесомость 59 Законы Кеплера движения планет . 61 Прогиб балок . . G8 Транспортировка леса . . . 71 КАЧЕСТВЕННЫЕ МЕТОДЫ ИССЛЕДОВАНИЯ ДИФФЕРЕНЦИАЛЬНЫХ МОДЕЛЕЙ Кривые с постоянным направлением магнитной стрелки . . , . 80 Зачем инженеру знать теоремы существования и единственности? 84 Динамическая интерпретация дифференциальных уравнений второго порядка 91 Консервативные системы в механике 95 Устойчивость точек равновесия и периодических движений . . 103 Энергетические функции .... . 107 Простые состояния равновесия 111 3
Движение тела единичной массы под действием линейных пружин в среде с линейным трением 114 Адиабатический поток идеального газа в канале переменного диаметра 119 Точки равновесия высшего порядка 123 Преобразование обратными радиусами и однородные коорди- координаты 127 Поток идеального газа во вращающемся канале постоянного диаметра 130 Изолированные замкнутые траектории 138 Периодические режимы в электрических цепях .... . . 145 Кривые без контакта . 150 Список литературы , , , 153
ПРЕДИСЛОВИЕ Дифференциальное уравнение является одним из ос- основных математических понятий. Дифференциальное урав- уравнение — это уравнение для отыскания функций, производ- производные которых (или дифференциалы) удовлетворяют некото- некоторым наперед заданным условиям. Дифференциальное урав- уравнение, полученное в результате исследования какого-либо реального явления или процесса, называют дифференциаль- дифференциальной моделью этого явления или процесса. Понятно, чго дифференциальные модели — это частный случай того мно- множества математических моделей, которые могут быть по- построены при изучении окружающего нас мира. При этом необходимо отметить, что существуют и различные типы самих дифференциальных моделей. А^ы будем рассматри- рассматривать лишь модели, описываемые так называемыми обыкно- обыкновенными дифференциальными уравнениями, одной из ха- характерных особенностей которых является то, что неизве- неизвестные функции в этих уравнениях зависят только от одной переменной. В процессе построения обыкновенных дифференциаль- дифференциальных моделей (да и не только их) важное, а подчас и первен- первенствующее значение имеет знание законов той области науки, с которой связана природа изучаемой задачи. Так, напри- например, в механике это могут быть законы Ньютона, в теории электрических цепей — законы Кирхгофа, в теории скоро- скоростей химических реакций — закон действия масс и т. д. 5
Конечно, на практике приходится иметь дело и с такими случаями, когда неизвестны законы, позволяющие соста- составить дифференциальное уравнение, и поэтому необходимо прибегать к различным предположениям (гипотезам), -ка- -касающимся протекания процесса при малых изменениях параметров — переменных. К дифференциальному урав- уравнению тогда приводит предельный переход. При этом, если окажется, что результаты исследования полученного диффе- дифференциального уравнения как математической модели согла- согласуются с опытными данными, то это и будет означать, что высказанная гипотеза правильно отражает истинное поло- положение вещей *). Работая над книгой, автор ставил перед собой две цели. Первая из них заключалась в том, чтобы на примерах (в основном содержательных, а не чисто иллюстративных) из различных областей знаний показать возможности ис- использования обыкновенных дифференциальных уравнений в процессе познания окружающей нас действительности. Ко- Конечно, рассмотренные примеры далеко не охватывают тот круг вопросов, которые могут быть решены с помощью обыкновенных дифференциальных уравнений. Но, во-пер- во-первых, «никто не обнимет необъятного», а во-вторых, уже и приведенные примеры дают представление о той роли, ко- которую играют обыкновенные дифференциальные уравнения при решении практических задач. Вторая цель — познакомить читателя с простейшими приемами и методами исследования обыкновенных диффе- дифференциальных уравнений, характерными для качественной теории дифференциальных уравнений. Дело в том, что лишь в редких случаях удается решить дифференциальное урав- уравнение в так называемой замкнутой форме, т. е. представить решение в виде аналитической формулы, использующей конечное число простейших операций над элементарными функциями. И это тогда, когда известно, что дифферен- *) Более подробно о математических моделях можно прочитать, например, в увлекательных книгах А. Н. Тихонова, Д. П. Костомарова «Рассказы о прикладной математике» (М.: Наука, 1979) и Н. Н. Мои- Моисеева «Математика ставит эксперимент» (М.: Наука, 1979).
циальное уравнение решение имеет! Другими словами, оказывается, что решения дифференциальных уравнений в своем многообразии таковы, что для их представления в замкнутой форме конечного числа аналитических операций недостаточно. Такая ситуация схожа с имеющей место в теории алгебраических уравнений: в случае алгебраических уравнений первой и второй степеней их решения могут быть легко получены в радикалах; если обратиться к уравнениям третьей и четвертой степеней, то решения в радикалах еще могут быть получены, но формулы становятся весьма слож- сложными; что же касается алгебраического уравнения общего вида степени выше четвертой, то решение такого уравне- уравнения в радикалах, вообще говоря, уже не может быть получено. Возвращаясь к дифференциальным уравнениям, отме- отметим, что если для представления их решений пользоваться бесконечными рядами того или иного вида, то удастся ре- решить значительно больше уравнений, чем в замкнутой фор- форме. Но, к сожалению, часто бывает так, что наиболее суще- существенные и интересные свойства решений никак нельзя выявить из вида полученных рядов. Более того, даже если удается решить дифференциальное уравнение и в замкну- замкнутой форме, то далеко не всегда такое решение можно проана- проанализировать, ибо полученная зависимость между различны- различными параметрами часто оказывается весьма и весьма сложной. Таким образом, становится очевидной необходимость в приемах и методах, которые позволяли бы, не решая самих дифференциальных уравнений, все же получать необходимые сведения о тех или иных свойствах решений. Так вот, такие приемы и методы существуют, и они и составляют содер- содержание качественной теории дифференциальных уравнений, в основе которой лежат общие теоремы о существовании и единственности решений, о непрерывной зависимости реше- решений от начальных данных и параметров. Частичное обсуж- обсуждение роли теорем существования и единственности реше- решений проводится в параграфе «Зачем инженеру знать теоремы
существования и единственности?». Что же касается каче- качественной теории обыкновенных дифференциальных урав- уравнений вообще, то, начиная с работ А. Пуанкаре и Л. М. Ля- Ляпунова (конец XIX-го века), в которых были заложены ее основы, она интенсивно развивается и ее методы широко используются в процессе познания окружающей нас дей- действительности. Автор благодарен профессорам Ю. С. Богданову и М. В. Федорюку за полезные советы и замечания, выска- высказанные в процессе работы над книгой. В. В. Амелькин
ПОСТРОЕНИЕ ДИФФЕРЕНЦИАЛЬНЫХ МОДЕЛЕЙ И ИХ РЕШЕНИЙ Чей кофе более горячий? Анатолий и Владимир заказали в кафе кофе и сливки. Когда им одновременно подали по чашке одинаково горя- горячего кофе и сливки, они поступили следующим образом. Анатолий добавил в кофе немного сливок, накрыл чашку бумажной салфеткой и вышел позвонить по телефону. Вла- Владимир сразу же накрыл чашку бумажной салфеткой, а до- добавил то же количество сливок только через 10 мин, когда вернулся Анатолий, и они начали пить кофе вместе. Кто же пил более горячий кофе? Задачу будем решать с учетом естественных предполо- предположений, которые отражают физическое содержание проис- происходящих процессов и заключаются в следующем. Считаем, что теплообмен через поверхность стола и салфетки намного меньше теплообмена через боковые стенки чашек; темпера- температура пара в чашке над поверхностью жидкости равна тем- температуре жидкости. Выведем сначала соотношение, показывающее, как с течением времени изменялась температура кофе в чашке Владимира до смешивания кофе со сливками. В соответствии с принятыми допущениями на основе известного закона физики количество теплоты, полученное воздухом от чашки Владимира, определяется соотношением dQ = r\^p^sdt, A) где Т — температура кофе в момент времени t, 6 — тем- температура воздуха в кафе, т\ — теплопроводность материала чашки, / — толщина стенок чашки, s — площадь боковой поверхности стенок чашки. С другой стороны, количество
теплоты, отданное кофе, находим из равенства B) где с — удельная теплоемкость кофе, т — масса кофе в чашке. Рассматривая теперь вместе уравнения A) и B), приходим к уравнению ■П —J— sdt = — cm dTy которое, разделяя переменные, можно переписать в виде *Г 2Edt. C) Г —0 1ст к ' Обозначая начальную температуру кофе через То и интегри- интегрируя дифференциальное уравнение C), находим, что I!L/ Г = е + (Г0—в)е ьт . D) Формула D) и есть аналитическое описание закона, по которому изменялась температура кофе в чашке Владимира до смешивания кофе со сливками. Посмотрим теперь, какой будет закон изменения темпе- температуры кофе после того, как Владимир добавил в чашку сливки. Для этого воспользуемся уравнением теплового баланса, которое в нашем случае запишется в виде сго(Г--вв) = с1т1(ев--Г1), E) где 6В — температура смеси в момент времени t, T± — температура сливок, с* — удельная теплоемкость сливок, 1Щ — масса сливок, добавленная в кофе. Из уравнения E) находим, что cm Принимая во внимание равенство D), формулу F) можно переписать в виде ев= , ЦОТУ Ti+ct CmC le + (To—G)<? lcm J • G) Равенство G) и задает закон изменения температуры кофе после добавления в чашку Владимира сливок. Для вывода закона изменения температуры кофе в чашке Анатолия снова воспользуемся уравнением теплового ба- баланса, которое в данном случае принимает вид cm (То—90)= схтх F0—7\), (8) 10
где 0О — температура смеси. Из равенства (8) получаем, что 9 -_, 0 ст-\- cm - Т А тогда, воспользовавшись уравнением D), где роль началь- начальной температуры играет уже 60, а произведение cm заме- заменяется суммой cm -j- c-jrii, окончательно получаем, что за- закон изменения температуры 6А кофе в чашке Анатолия ана- аналитически задается формулой О cm cm -\- Citnx (9) Таким образом, для ответа на поставленный в задаче вопрос остается лишь обратиться к формулам G) и (9) и провести численные расчеты, имея в виду, что сх ж 3,9 х х108 Дж/(кг-К), ^^4,1-10;1Дж/(кг..К),'ПЯ=!0,6 В/(м-К), и полагая для определенности nii~2-10~2 кг, т=8-10~? кг, Ti~2Q C, 0=20 vJ, T0=80 C, s=lblO 3 м2, /=2« 10 3 м. Вычисления показывают, что более горячий кофе пил Ана- Анатолий. Стационарный тепловой поток Прежде всего напомним, что о стационарном тепловом потоке говорят в «том случае, когда температура тела в каждой точке со временем не меняется. При решении задач, физическое содержание которых связано с влиянием тепловых потоков, существенную роль играют так называемые изотермические поверхности. „ 20 ■ Ряс. 1 Рис. 2 Для пояснения рассмотрим, например, теплопроводную трубу (рис. 1) диаметром 20 см, сделанную из однородного материала и защищенную покрытием из магнезии толщиной 10 см. Предположим, что температура трубы равна 160°С, 11
а внешнее покрытие имеет температуру, равную 30 °С. Тогда интуитивно ясно, что существует поверхность, се- сечение которой на рис. 2 показано пунктиром, в каждой точке которой температура будет одной и той же, например, рапной 95 °С. Пунктирная кривая на рис. 2 называется изотермической кривой, соответствующая же ей поверх- поверхность называется изотермической поверхностью. В общем случае изотермические кривые могут иметь самый разнооб- разнообразный вид, что, в частности, связано с нестационарностью теплового потока и неоднородностью материала. В рассмат- рассматриваемом нами случае изотермическими кривыми (поверх- (поверхностями) будут концентрические окружности (цилиндры). Выведем закон распределения температуры внутри по- покрытия и найдем количество теплоты, выделенное трубой на участке длиной 1 м в течение суток, если коэффициент теплопроводности /г=1,7-10~4. Для этого воспользуемся законом теплопроводности Фурье, согласно которому количество теплоты, излучаемое в единицу времени телом, находящимся в неизменном тепло- тепловом состоянии, температура Т которого в каждой точке есть функция только одной координаты х, находится по формуле ^ = const, A0) где F(x) — площадь сечения, перпендилулярного направле- направлению распространения тепла, k — коэффициент теплопро- теплопроводности. Из условий задачи следует, что в рассматриваемом слу- случае F(x) — 2тсх1, где / — длина трубы, см, к — радиус основания цилиндрической поверхности, расположенной внутри внешнего цилиндра. А тогда на основании формулы A0) приходим к равенствам 30 20 J dTz=z~~ 0,00017.2л/J T' (П) F0 10 [ Лт Q f dl J ' 0,00017- 2л/ J I' 160 10 Интегрируя соотношения A1) и A2)', получаем, что 160 — T_JnO,lY__lgO,iA: 130 ~~ In2 lg 2 " Отсюда Т = 591,8 — 431,8 \gx. 12
Последней формулой и задается закон распределения тем- температуры внутри покрытия. Как видим, длина трубы здесь никакой роли не играет. Чтобы ответить на второй вопрос, обратимся к уравне- уравнению A1). Тогда при 1 — 100 см получаем, что п 130.0,00017.2я-100 __200я-130-0,00017 У~ In 2 ~ 0,69315 а поэтому количество теплоты, выделенное в течение суток, равно 24-60.60Q=726852 Дж. Случай в заповеднике При обходе заповедника два егеря обнаружили тушу уби- убитого дикого кабана. Ее осмотр показал, что выстрел браконь- браконьера был точным и кабан убит наповал. Рассудив, далее, что браконьер должен вернуться за добычей, егеря решили дож- дождаться его, укрывшись недалеко от того места, где лежала туша. Вскоре показались два человека, прямо направляв- направлявшиеся к убитому животному. Задержанные неизвестные вся- всячески отрицали свою причастность к браконьерству. Однако у егерей уже были косвенные улики их виновности, но для ее полного доказательства следовало еще уточнить время, когда был убит кабан. Это удалось сделать с помощью закона излучения тепла. Покажем какими соображениями можно было при этом поль- пользоваться. Согласно закону излучения тепла скорость охлаждения тела в воздухе пропорциональна разности между темпера- температурой тела и температурой воздуха, т. е. ^^~k(x-a), A6) где х — температура тела в момент времени t; a — темпе- температура воздуха; k — положительный коэффициент пропор- пропорциональности. Решение задачи связано с исследованием соотношения, получающегося в результате интегрирования дифференци- дифференциального уравнения A3). При этом следует учитывать, что после того, как кабан был убит, температура воздуха могла оставаться неизменной, а могла и меняться с течением вре- времени. В первом случае интегрирование дифференциального уравнения A3) с разделяющимися переменными приводит к равенству \n^- = -kt, х^=а, A4) ■Л 0— " 13
где х0 — температура тела в момент времени *=0. А тогда если в момент задержания неизвестных температура туши кабана к была равна 31 °С, а спустя час составляла 29 °С, то, сбитая, что В момент выстрела в кабана его температура была х = 37 °С, а температура воздуха а = 21 °С, можно, полагая t—О временем задержания неизвестных, опреде- определить и время выстрела. Так, воспользовавшись имеющимися данными, из соотношения A4) получим, что /г=1п|Ь^ = 1п 1,25 = 0,22314. A5) Подставляя теперь в формулу A4) значение k из равенства A5) и значение х=37, находим Ш3121 022314 m ]'°" Z, '"" 0,22314 Ш31-21 0,22314 Иначе говоря, между моментом выстрела и тем моментом, когда неизвестные были задержаны, прошло 2 часа и 6 ми- минут. В том случае, когда температура воздуха меняется со временем, закон охлаждения тела запишется в виде линей- линейного неоднородного дифференциального уравнения ~ + fcc=foz(*), A6) где a (t) — температура воздуха в момент времени /. Для иллюстрации одного из методов определения момен- момента времени, когда был убит кабан, предположим, что в мо- момент задержания неизвестных температура туши кабана была равна 30 °С. Пусть известно также, что в день случив- случившегося температура воздуха падала в течение каждого часа после полудня на 1 °С и в момент обнаружения, туши была равна 0°С. Предположим, далее, что через час после обна- обнаружения температура туши стала равной 25 °С, а темпера- температура воздуха понизилась до —1°С. Если теперь принять за момент выстрела браконьера /=0 и считать, что в этот мо- момент л:0=37°С, то, полагая время обнаружения убитого ка- кабана £=/*, получим a{i)^t*—t. Интегрируя теперь уравнение A6), придем к соотноше- соотношению 14
Далее, имея в виду, что л:=30 при t=t* и #=25 при t=t*-\- -Ы, из последнего равенства получаем соотношения C7—/*—!)<>-*<*+1 = 30, которые позволяют вывести уравнение относительно /г, а именно уравнение 1)-26+1 = 0. A7) К уравнению A7) можно прийти, исходя и из других начальных предпосылок. Действительно, примем за t=0 время обнаружения убитого кабана, тогда а@=—t, и мы приходим к дифференциальному уравнению . A8) (с начальным условием лто=ЗО при £=0), из которого требу- требуется найти л: как явную функцию L Решая уравнение A8), получаем, что Полагая в последнем соотношении t—l и #=25, мы и при- приходим к уравнению A7), позволяющему численно решить исходную задачу. Действительно, как известно, уравнение A7) не может быть алгебраически разрешено относительно к. Вместе с тем оно легко решается численными методами нахождения корней трансцендентных уравнений, в частности, методом последовательных приближений Ньютона- Метод Нью- Ньютона, как и другие методы последовательных приближений, является способом, посредством которого грубая оценка истинного значения корня используется для получения более точных его оценок. Причем процесс продолжается до тех пор, пока не достигается желаемая точность. Чтобы показать, как пользоваться методом Ньютона, приведем уравнение A7) к виду 30/г—1-f A— 26/г)<?*=;0, B0) 15
а уравнение A9), полагая в нем х~37, — к виду м — 30/е+ 1 -0. B1) Последние два уравнения — это уравнения вида (ах 4 Ь) еКх + cx + d = Q. B2) Если теперь левую часть уравнения B2) обозначить через ср (х), то дифференцирование по х дает следующие равенства: <р' (х) = (box 4- М> 4- а) ср" (х) = (Ь* А тогда метод Ньютона нахождения корня уравнения B2) состоит в том, что если для /-го приближения xt выполня- выполняется неравенство то следующее приближение xi+i находится по формуле х =>х -_-!ЕМ. Для непосредственного вычисления корня (с точностью, например, до 10~6) с помощью микрокалькулятора «Элект- «Электроника БЗ-34» составим программу: ПЗ; ИПЗ; ИП6; X; Fe*; П4; ИПА; ИПЗ; X; ИПВ; 'П5;'ИПС; ИПЗ; X; +; ИГО; +; ПО; С/П; ИПА; ИП4; X; П4; ИП5; ИП6; X; ИП4; +; П5; ИПС; +; П1; С/П; ИП5; ИП4; +; ИП6; X; П2; С/П; ИПЗ; ИПО; ИП1; Н-; —; ПЗ; С/П; В/О. Эта программа использует регистры памяти 0—6 и А, В, С, D, назначение которых видно из следующей таблицы: Регистры памяти Содержа- А а В ъ с с D d 6 X 0 ф 1 ф' 2 Ф" 3 */ 4 екх 5 (ах+Ь) еКх Нахождение корня осуществляется в следующей после- последовательности: 1) занести в регистры памяти А, В, С, D и 6 коэффи- коэффициенты уравнения B2) в соответствии с последней таблицей; J6
2) набрать приближенное (начальное) значение корня Xi(x0) и нажать клавишу «В/О»; и далее: 3) С/П; Выписать значение ip (xf); 4) С/П; Выписать значение cp'(Xj); 5) С/П; Выписать значение ср"(л;г); 6) С/П; Выписать следующее приближение x{+i\ 7) Повторить процедуру, начиная с пункта 3). Отметим здесь, что вычисление заканчивается в том слу- случае, когда приближения xt nxi+i содержат требуемое коли- количество одинаковых значащих цифр. Используя эту общую процедуру, обратимся к уравне- уравнению B0). Для его левой части ср(/г) дифференцирование по k приводит к равенству А тогда нетрудно проверить, что <р@)-^0, фA)<1, q/@)>(L Таким образом, функция ф возрастает в малой окрестности начала координат, а затем убывает до отрицательного значения при k—\. Отсюда следует, что на интервале @, 1) существует корень уравнения ц>{к)=0. Беря в качестве начального приближения /го=О,5 и учитывая, что в пашем случае а~—26, 6=1, е=30, d——1, К=\, в соответствии с написанной выше программой приходим к таблице п 0 1 2 3 4 5 6 К 0,5 0,322835 0,228162 0,188497 0,179453 0,178953 0,178952 —5,784655 —1,525956 —0,351424 —0,055194 —0,002747 —0,000008 0 —32,651408 — 16,118043 —8,859807 —6,103385 —5,497019 —5,463736 ч>"(кп) — 105,5 —82,0 —71,5 —67,5 —66,6 —66,5 Окончательный шаг в решении задачи заключается в подстановке вычисленного значения kc^k^Q, 178952 в урав- уравнение B1) и решении последнего относительно t (времени когда был убит кабан). Чтобы воспользоваться описанной выше схемой, обозначим левую часть уравнения B1) через g(t). А тогда, выбирая в качестве t0 значение £0——1 и имея в виду, что в данном случае a=k, b~37k—1, с=0, d——30/г-Ь + X—k, получаем таблицу 17
n 0 1 2 3 4 *n — 1,0 — 1,188775 — 1,192556 — 1,192558 —1,192558 g (tn) 0,181972 0,003505 0,000001 0,0000002 —0,0000003 0,963962 0,927054 0,926329 0,926329 g"(tn) 0,19 0,19 0,19 0,19 Из полученных результатов следует, что кабан был убит приблизительно за 1 час и 12 минут до обнаружения его егерями. Истечение жидкости из сосудов. Водяные часы Приводимые ниже две задачи иллюстрируют связь между их физическим содержанием и геометрией. Предварительно остановимся на некоторых общих теоре- теоретических выводах. Итак, рассмотрим сосуд (рис. 3), пло- площадь горизонтального сечения которого является произ- произвольной функцией расстояния сечения от дна сосуда. Пусть Рис. 3 высота уровня жидкости в сосуде в начальный момент вре- времени t = 0 равна h метров. Пусть, далее, площадь сечения на высоте х равна S(x), а площадь отверстия на дне сосуда есть s. Известно, что скорость истечения жидкости v в тот момент, когда высота ее уровня равна х, определяется ра- равенством v — k V%gx, где g=9,8 м/с2, k — коэффициент ско- скорости истечения жидкости из отверстия. На бесконечно малом промежутке времени dt истечение жидкости можно считать равномерным, а потому за время dt вытечет столбик жидкости, высота которого v dt и пло- 18
щадь сечения s, что в свою очередь вызовет понижение уровня жидкости в сосуде на —dx. В результате этих рассуждений приходим к дифферен- дифференциальному уравнению ks V2gx dt^ — S (x) dx, B 2) которое можно переписать в виде dt = Щ=дх. B3) Решим теперь следующую задачу. Цилиндрический резервуар с вертикальной осью высотой б м и диаметром 4 м имеет па дне круглое отверстие радиусом 1/12 м. Требуется установить зависимость уровня воды в резервуаре от вре- времени I, а также определить время, в течение которого выте- вытечет вся вода. По условиям задачи S(x)=4:rc, s= 1/144. Так как для годы /г=0,6, то уравнение B3) примет вид .. 217,152 . У х Интегрируя это дифференциальное уравнение, приходим к соотношению (= 434,304 Vx> которое и дает искомую зависимость уровня воды от вре- времени t. Если теперь в последнем равенстве положить х—6, то получим, что вся вода вытечет из резервуара приблизи- приблизительно через 18 минут. Вторая задача состоит в следующем. Известно, что древ- этие водяные часы представляли собой чашу (рис. 4), из которой через небольшое отверстие на дне вытекала вода. Такие часы использовались в греческих и римских судах для хронометрирования речей адвокатов, чтобы не допускать слишком долгих выступлений. Требуется найти форму водя- водяных часов, при которой уровень воды убывал бы в чаюе с постоянной скоростью. Задача легко решается с помощью выведенного выше уравнения B3), которое мы только перепишем в виде Именно, учитывая, что чашу можно рассматривать как поверхность вращения, в соответствии с обозначениями на 19
рис. 4 из уравнения B4) получаем, что V~x = ^=а, B5) где a = vx= ~ — проекция скорости свободной поверхности жидкости на ось х, которая по условию задачи есть величи- величина постоянная. Возведя обе части уравнения B5) в квадрат, приходим к уравнению х=сг*, B6) где c=a?n2/Bgk2s2). Последнее означает, что форма поверх- поверхности водяных часов получается вращением кривой B6) Еокруг оси х. Эффективность рекламы Предположим, что торговыми учреждениями реализуется продукция В, о которой в момент времени t из числа потен- потенциальных покупателей N знает лишь х покупателей. Пред- Предположим далее, что для ус- ускорения сбыта продукции В были даны рекламные объяв- объявления по радио и телевидению. Последующая информация о продукции распространяется среди покупателей посредст- г вом общения друг с другом. рис> 5 С большой степенью досто- достоверности можно считать, что после рекламных объявлений скорость изменения числа знающих о продукции В пропорциональна как числу знаю- знающих о товаре покупателей, так и числу покупателей, о нем еще не знающих. Если условиться, что время отсчитывается после рек- рекламных объявлений, когда о товаре узнало N/y человек, то приходим к дифференциальному уравнению dx _Ц±. — Uy /Д/ Y\ /О7\ I. — /ел yi\ —Л/) V / с начальными условиями x=Nly при t=Q. В уравнении B7) коэффициент k — это положительный коэффициент пропорциональности. Интегрируя уравнение B7), находим, что 1 . х 20
Полагая NC=Ci, приходим к равенству xl(N — х) = A eNkt, где А= ес». Если последнее уравнение разрешить относительно х, то получим соотношение где Р=\/А. В экономической литературе уравнение B8) обычно называют уравнением логистической кривой. Если учесть теперь начальные условия, то уравнение B8) /V перепишется в виде х Y На рис. 5 схематически изображена логистическая кри- кривая при Y=2. В заключение отметим, что к уравнению B7) сводится, в частности, задача о распространении технологических новшеств. Спрос и предложение Как известно, спрос и предложение — экономические категории товарного производства, возникающие и функ- функционирующие на рынке, в сфере товарного обмена. При этом спрос — представленная на рынке потребность в товарах, а предложение — продукт, который есть на рынке или мо- может быть доставлен на пего. Одним из экономических за- законов товарного производства является закон спроса и предложения, который заключается в единстве спроса и предложения и их объективном стремлении к соответствию. Рассмотрим следующую задачу. Пусть в течение неко- некоторого (достаточно продолжительного) времени крестьянин продает на рынке фрукты (например, яблоки), причем про- продает их после уборки урожая, с недельными перерывами. Тогда при имеющихся у крестьянина запасах фруктов не- недельное предложение будет зависеть как от ожидаемой цены в наступающей неделе, так и от предполагаемого изменения цены в последующие недели. Если в наступающей неделе предполагается, что цена упадет, а в последующие недели повысится, то предложение будет сдерживаться при усло- условии превышения ожидаемого повышения цен над издерж- издержками хранения. При этом предложение товара в ближай- ближайшую неделю будет тем меньшим, чем большим предпола- предполагается в дальнейшем повышение цены. И наоборот, если в наступающей неделе цена будет высокой, а затем ожи- 21
дастся ее падение, то предложение увеличится тем больше, чем большим предполагается понижение цены в дальней- дальнейшем. Если обозначить через р цену на фрукты в наступаю- наступающей неделе, а через р' — так называемую тенденцию фор- формирования цены (производную цены по времени), то как спрос, так и предложение будут функциями указанных величин. При этом, как показывает практика, в зависи- зависимости от разных факторов спрос и предложение могут быть различными функциями цены и тенденции формирования цены. В частности, одна из таких функций задается ли- линейной зависимостью, математически описываемой соот- соотношением у=ар'-\-Ьр-\-с, где а, Ъ, с — некоторые вещест- вещественные постоянные. А тогда если, например, в рассматри- рассматриваемой задаче цена на фрукты вначале составляла 1 р. за 1 кг, через t недель она была уже р(t) p. за 1 кг, а спрос q и предложение s определялись соответственно соотноше- соотношениями q = 4p'— 2р + 39, s= 44//+ 2/7—1, то для того чтобы спрос соответствовал предложению, не- необходимо выполнение равенства 4//—2/?+ 39 = 44//+ 2/? — 1. Отсюда приходим к дифференциальному уравнению dP ._ р—10 Интегрируя, находим, что p=Ce~lot + 10. Если же учесть начальные условия р— 1 при <=0, то окончательно полу- получаем р*=— 9<?-ш+10. B9) Таким образом, если требовать, чтобы между спросом и предложением все время сохранялось равновесие, необ- необходимо, чтобы цена изменялась в соответствии с фор- формулой B9). Химические реакции Химическое уравнение показывает, как в процессе вза- взаимодействия одних веществ образуется другое вещество. Так, например, уравнение 2Н2 + О2-*2Н2О 22
показывает, что в результате взаимодействия двух молекул водорода и одной молекулы кислорода получаются 2 моле- молекулы воды. В общем случае химическое уравнение записывается в виде аА-\-ЬВ + сС-\-... —+ тМ-\-nN + рР-\-..., где А, В, С, . . . — молекулы взаимодействующих веществ, М, N, Р, ... — молекулы веществ, полученных в ре- результате химической реакции, а постоянные а, Ь, с, . . . ..., т, п, р, . . . — положительные целые, указывающие на число молекул, принимающих участие в реакции. Скорость, с которой образуется новое вещество, назы- называется скоростью реакции. Действующая же масса или концентрация реагирующего вещества описывается коли- количеством молей этого вещества в единице объема. Одним из основных законов теории скоростей хими- химических реакций является закон действующих масс, соглас- согласно которому скорость химической реакции при постоянной температуре пропорциональна произведению концентраций веществ, участвующих в данный момент в реакции. Решим следующую задачу. Два жидких химических вещества А и В объемом 10 и 20 литров соответственно в процессе химической реакции образуют новое жидкое химическое вещество С. Считая, что температура в процессе реакции не меняется, а также что из каждых двух объемов вещества А и одного объема вещества В образуется три объема вещества С, определить количество вещества С в произвольный момент времени /, если за 20 мин его обра- образуется б л. Обозначим через х объем (в литрах) вещества С, образо- образовавшегося к моменту времени t (в часах). Тогда из условий задачи следует, что к этому моменту времени в химическую реакцию вступило 2л;/3 литров вещества А и х/3 литров вещества В. Последнее означает, что к указанному моменту осталось 10—^ литров вещества А и 20—-^ литров веще- вещества В. Таким образом, в соответствии с законом действую- действующих масс приходим к дифференциальному уравнению которое можно переписать в виде It''
где k — постоянная пропорциональности (&=2/(/9). При этом следует иметь в виду, что так как в начальный момент времени /=0 вещества С еще не было, то можно считать, что в этот момент времени х—0. Что же касается момента времени /=1/3, то здесь уже х=6. Итак, решение исходной задачи свелось к решению так называемой краевой задачи i^ = fcA5—л;) F0—х), л;@) = 0, хA/3) = 6. Для ее решения проинтегрируем сначала последнее дифференциальное уравнение с учетом начального условия лг(О):=О. В результате получим соотношение F0—х)/{15—х) = №ш. Теперь, так как х=в при /=1/3, то, подставляя эти значения переменных в последнее равенство, находим, что elbk=3/2. А тогда F0— х)/A5—х) =-- 4 (е1^)** = 4 C/2K*, т. е. х=15A—B/3K1)/A—A/4) B/3K'). Последнее равенство и определяет количество вещества С, образовавшегося в результате реакции к моменту времени /. Сделаем следующее замечание. Из практических сооб- соображений понятно, что в процессе химической реакции 10 л вещества Л и 20 л вещества В может образоваться веще- вещество С лишь конечного объема. Вместе с тем формальное рассмотрение полученной зависимости х от / показывает, что при конечном значении /, а именно в случае, когда B/3)зг=4, переменная х обращается в бесконечность. Этот факт, однако, не противоречит практическим соображе- соображениям, ибо последнее равенство возможно только при отри- отрицательном значении /, а процесс химической реакции рас- рассматривается лишь при Дифференциальные модели в экологии Экология изучает взаимоотношение человека и вообще живых организмов с окружающей средой. Основным объек- объектом экологии является эволюция популяций. Ниже описываются дифференциальные модели популя- популяций, которые связаны с размножением или вымиранием последних, а также с сосуществованием различных видов животных в ситуации «хищник — жертва»*). *) Murray J. D. Some simple mathematical models in ecology// Math. Spectrum.— 1983—1984.— V. 16, № 2.— P. 48—54. 24
Пусть x(t)— число особей в популяции в момент вре- времени t. Тогда если А — число особей в популяции, рож- рождающихся в единицу времени, а В — число особей, уми- умирающих в единицу времени, то с достаточным основанием можно утверждать, что скорость изменения х со временем задается формулой ~ = Л-В. C0) Задача теперь состоит в том, чтобы описать зависимость А и В от х. Простейшим случаем является ситуация, когда А=ах, B=^bxy C1) где а и Ь — коэффициенты рождения и смерти особей в единицу времени соответственно. С учетом равенств C1) дифференциальное уравнение C0) перепишется в виде *L = {u-b)x. C2) Полагая, что в момент времени t=t0 число особей в популя- популяции есть х=х0, из уравнения C2) находим л: (*) = -»> <t-fo> Из полученного равенства следует, что если а > Ь, то при t-^oo число особей х-> оо. С другой стороны, если a<L <С Ь, то х ->■ 0 при /->-оои популяция становится вымираю- вымирающей. Хотя приведенная модель является упрощенной, ока все-таки в ряде случаев соответствует действительности. Практически же все модели, которые описывают реальные явления и процессы, нелинейны, и вместо дифференциаль- дифференциального уравнения C2) следует рассматривать уравнение вида dxг / ч где f(x) — нелинейная функция, например, уравнение вида —- = f(x) = ах—Ъх%, где а>0, 6>0. Полагая, что х=х0 при t—t0, из последнего уравнения 25
находим, что ^r«-w C3) Отсюда видно, что при /->-оо число особей в популяции x(t) ->alb. При этом возможны два случая: alb~>xu и а/Ь< <Сх0. Различие между этими Л{fслучаями хорошо видно из рис. 6. Отметим, что соотношение C3) описывает, в частности, попу- популяции фруктовых вредителей и некоторых видов бактерий. Если рассматривается не- сколько сосуществующих видов, t например, больших и малых Рис, 6 рыб, где малые рыбы являются кормом для больших, то, состав- составляя дифференциальные уравнения для каждого вида, полу- получим систему дифференциальных уравнений Рассмотрим более подробно двухвидовую модель «хищ- «хищник— жертва», которая впервые была построена Вольтер- ра для объяснения колебаний рыбных уловов в Адриати- Адриатическом море, имеющих один и тот же период, но отличаю- отличающихся по фазе. Пусть х — число больших рыб-хищников, которые пи- питаются малыми рыбами-жертвами, число которых обозна- обозначим через у. Тогда число рыб-хищников будет расти до тех пор, пока у них будет достаточно пищи, т. е. малых рыб- жертв, но в конце концов наступит ситуация, когда корма не будет хватать и в результате число больших рыб начнет уменьшаться. Это приведет к тому, что с некоторого мо- момента число малых рыб снова начнет увеличиваться. Это будет способствовать новому росту числа больших особей, и цикл снова повторится. Модель, построенная Вольтерра, имеет вид dx ax + by C4) -JL^ ex—d xy, C5) где a, b, c, d — положительные константы. В уравнении C4) для больших рыб слагаемое Ьху вы- выражает зависимость прироста больших рыб от численности 26
малых рыб. В уравнении же C5) слагаемое —dxy выражает уменьшение числа малых особей в зависимости от числен- численности больших. Для удобства исследования последних двух уравнений введем в рассмотрение безразмерные переменные / ^ d /46 4. LL \ t J — ^^ Л| С I 11 — у j L — isV | \л> —— Ы'1 к)• С (X В результате дифференциальные уравнения C4) и C5) при- примут вид где а > 0, а штрих обозначает дифференцирование по т. Предположим, что в некоторый момент времени i:=To число особей обоих видов известно, т. е. и (То) _ и v (То) _ v ^37) Заметим, что в дальнейшем мы интересуемся только поло- положительными решениями. Выявим связь между и и v. Для этого, разделив первое уравнение из системы C6) на второе и затем проинтегрировав полученное дифференциальное уравнение, найдем, что av + и — In Фи = av0 + и0—In v$u0 = Я, где Н — постоянная, определяемая начальными условиями C7) и параметром а. На рис. 7 показан вид графиков и как функции v при различных значениях Н. Как видно из этого рисунка, в плоскости (и, v) имеются только Vk замкнутые кривые. Предположим теперь, что начальные значения щ и 1>0 задаются точкой А на траек- траектории, соответствующей значению 1 _ Н=Н3- Поскольку ио>\, со<1, то первое уравнение из системы C6) показывает, что переменная и вна- вначале убывает. Аналогичный факт имеет место и для переменной v. Далее, когда переменная и дости- достигает значения, равного единице, то v'—О и затем в течение длительного времени т перемен- переменная v будет возрастать. Когда же v=\, то и'=0 и затем уже возрастать начинает переменная и. Таким образом, как переменная и, так и переменная v пробегают замкнутую траекторию. А это означает, что решения являются функ- функциями, периодическими по времени При этом максимум и 27 Рис. 7
й£ попадает на максимум v, т, е. колебания в популяции происходят в разных фазах. Типичный график зависимости и и v от времени т показан на рис. 8 (в случае v£>l, uo<\). Рис В заключение отметим, что изучение сообществ, взаи- взаимодействующих более сложным образом, дает более инте- интересные практические результаты. Так, например, если две популяции конкурируют в борьбе за один и тот же источ- источник питания (третья популяция), то можно показать, что один из видов начнет вымирать. При этом понятно, что если этим видом окажется источник питания, то такая же участь постигнет и два других вида. Одна задача математической теории эпидемий Рассмотрим одну из дифференциальных моделей, кото- которая встречается в теории эпидемий. Предположим, что некая популяция, состоящая из N особей, подразделяется на три группы. В первую из них включаются особи, которые вос- восприимчивы к некоторой конкретно имеющейся в виду бо- болезни, по здоровы. Число таких особей в момент времени / будем обозначать через S(t). Во вторую группу объединя- объединяются особи, которые являются инфекционными — они сами больны и являются источником распространения болезни. Число таких особей в популяции в момент времени / обоз- обозначим через /(/). Наконец, третья группа — это особи, которые здоровы и обладают иммунитетом к данной болезни. Число таких особей в момент времени / обозначается через R(t). Таким образом, S(/) + /@-!-Я (О-Л- C8) Предположим далее, что в случае, когда число инфек- инфекционных особей превосходит некоторое фиксированное 28
число /*, скорость изменения числа восприимчивых к болезни особей будет пропорциональна числу самих вос- восприимчивых особей. Что же касается скорости изменения числа инфекционных, но выздоравливающих особей, то ее будем считать пропорциональной числу инфекционных особей. Понятно, что эти предположения упрощают реаль- реальную ситуацию, но в ряде случаев они отражают существо дела. В связи с первым предположением будем считать, что когда число инфекционных особей /(/)>/*, то они способны заражать восприимчивых к болезни особей. Последнее оз- означает, что принимается во внимание факт изоляции (до некоторого момента времени) инфекционных особей (каран- (карантин или нахождение вдали от восприимчивых к болезни особей). Таким образом, приходим к дифференциальному уравнению dS___(—aS, если /(/)>/*, qq dt " \ 0, если /(/)</*. Теперь, поскольку каждая восприимчивая к болезни особь, которая в конце концов заболевает, сама становится инфекционной, то скорость изменения числа инфекционных особей представляет разность за единицу времени между вновь заболевшими особями и теми, которые уже выздорав- выздоравливают. Итак, S—W> ссли 7 (')>/*. —р/, если /(/)</*. Постоянные пропорциональности аир будем называть коэффициентами заболеваемости и выздоровления соответ- соответственно. Наконец, скорость изменения числа выздоравливающих ^ „ dR aT особей задается уравнением -тг — \>'- Для того чтобы решения соответствующих уравнений определялись однозначно, необходимо задать начальные условия. Для простоты предположим, что в момент времени £=0 в популяции нет особей с иммунитетом к болезни, т. е. R@)=-0, и что первоначально число инфекционных особей равно /@). Далее предположим, что коэффициенты забо- заболеваемости и выздоровления равны, т. е. а—{3 (случай неравенства этих коэффициентов предлагается для рассмот- рассмотрения читателю). В результате приходим к необходимости рассмотрения двух случаев. Случай 1. Число /@)^ /*. В этом случае с ростом времени особи в популяции не будут подвергаться 29
заражению болезнью, поскольку в этом случае dS/dt=O и, значит, в соответствии с уравнением C8) и условием R @)—О, для всех / справедливо равенство S(/) = S@) = N — /@). Рассматриваемый случай соответствует той ситуации, когда довольно много инфекционных особей оказываются в изо- изоляции. В этом случае из уравнения D0) приходим к диф- дифференциальному уравнению dl . Отсюда I (l)—I @)e~at и, значит, R(t) = N—S(l) — /@ = /@)[l— e~at\. На рис. 9 графически показано изменение числа особей с ростом t в каждой из трех групп. Рис. 9 Случай 2. Число /@)>/*. В этом случае дол- должен существовать интервал 0 ^ / < 7\ для всех значений t которого справедливо неравенство /(/)>/*, ибо по смыслу задачи / как функция t должна быть функцией непрерыв- непрерывной. Отсюда следует, что для всех t из промежутка [0, Т) болезнь будет распространяться на восприимчивых к ней особей. Таким образом, из уравнения C9) следует, что S{t) = S{0)e-°* для 0 < t < Т. Подставляя значение 5@ из последнего равенства в уравнение D0), приходим к дифференциальному уравнению dt D1) 30
Если теперь умножить обе части уравнения D1) на eat, то это уравнение примет вид dt Отсюда /еа*=а5@)/-ЬС и, значит, множество всех реше- решений уравнения D1) задается соотношением / (t) = Ce~at + a5 @) te~at. D2) Полагая здесь /=0, получаем С=/@), и, таким образом, уравнение D2) принимает вид /(/) = [/ @) -f а5 @) /] e~at D3) для 0 < t < Т. Дальнейшие исследования мы свяжем с нахождением конкретного значения Т и нахождением того момента вре- времени /тах, при котором число инфекционных особей ока- оказывается максимальным. Ответ на первый вопрос важен с той точки зрения, что в момент времени Т заражаемость восприимчивых к болезни особей прекращается. Если обратиться к уравнению D3), то из предыдущего следует, что при t=T его правая часть принимает значение /*, т. е. /* ;= [/ @) + а5 @) Т] е~ат. D4) Но 5 (Г) = lim S(t) = S (oo) есть число восприимчивых особей, которые избегают'заболе- вания и для которых справедлива такая цепочка равенств: Отсюда Г=-1п|Д. D5) a S(oo) v ' Таким образом, если мы сможем указать явное значение для 5(оо), то тем самым мы сможем использовать уравнение D5) для предсказания времени прекращения эпидемии. Подставляя тогда значение Т из соотношения D5) в урав- уравнение D4), получим равенство или с, . S(oo) '* =-Шг+Ш
которое можно переписать в виде D6) Поскольку /* и все члены в правой части уравнения D6) известны, мы можем использовать это уравнение для опре- определения 5(оо). Рис. 10 Чтобы ответить на второй вопрос, обратимся к уравне- уравнению D3), из которого, в соответствии с поставленным вопро- вопросом, приходим к равенству -^ = [aS @)—а/ @)—aaS @) /] е~ш = 0. Отсюда время tfmax, за которое / достигает максимального значения, задается соотношением 'max ~ 2 S(O) Если теперь подставить значение tmM из последнего ра- равенства в уравнение D3), то получаем, что Полученное равенство показывает, в частности, что в мо- момент времени tfmax число восприимчивых к болезни особей совпадает с числом инфекционных особей. Когда же t> T, то восприимчивые особи уже не ста- становится инфекционными и 32
На рис. 10 схематически показано, как с течением вре- времени изменяется число особей в каждой из трех рассматри- рассматриваемых групп. Кривая погони Приведем один из примеров использования дифферен- дифференциальных уравнений для выбора правильной стратегии при решении задач поиска. Пусть, например, миноносец охотится за подводной лодкой в густом тумане. В какой-то момент времени туман поднимается и подводная лод- лодка оказывается обнаруженной на поверхности воды на рас- расстоянии 3 миль от минонос- миноносца. Скорость миноносца вдвое больше скорости подводной лодки. Требуется определить траекторию (кривую погони), по которой должен следовать __ миноносец, чтобы он прошел ~ ~ точно над подводной лодкой, если последняя сразу же по- погрузилась после ее обнаружения и ушла на полной скорости прямым курсом в неизвестном направлении. Для решения сформулированной задачи введем полярные координаты г, G таким образом, чтобы полюс О находился в точке обнаружения подводной лодки, а полярная ось г проходила через точку, в которой в момент обнаружения подводной лодки был миноносец (рис. И). Дальнейшие рас- рассуждения основаны на следующих соображениях. Прежде всего миноносцу надо занять такую позицию, чтобы он и подводная лодка находились на одном расстоянии от по- полюса О. Затем миноносец должен двигаться вокруг полюса О по такой траектории, чтобы оба движущихся объекта все время находились на одинаковом расстоянии от точки О. Только в этом случае миноносец, обходя вокруг полюса О, пройдет над подводной лодкой. Из вышесказанного следу- следует, что сначала миноносец должен идти прямым курсом к точке О до тех пор, пока он не окажется на том же расстоя- расстоянии х от полюса О, что и подводная лодка. Очевидно, что расстояние х можно найти либо из урав- уравнения х _3—х v ~ Ъ ' 2 В. В Амелькив 33
либо из уравнения х _3+х v ~~ 2v ' где v — скорость подводной лодки, a 2v — скорость ми- миноносца. Решая последние уравнения, находим, что либо расстояние х равно одной, либо трем милям. Теперь, если «встречи» не произошло, то миноносец должен в дальнейшем двигаться вокруг полюса О (по направлению движения часовой стрелки или против), удаляясь от последнего со скоростью подводной лодки v. Разложим скорость миноносца 2v на две составляющие: ра- радиальную vr и тангенциальную vT (рис. 11). Радиальная составляющая — это скорость, с которой миноносец удаляется от полюса О, т. е. dr vr—df Тангенциальная составляющая — это линейная ско- скорость вращения миноносца относительно полюса. Она, dG как известно, равна произведению угловой скорости -п- на радиус г, т. е. VT r dt' Но так как vr=v, то vT = VBv)*—v2 = КЗ v. Итак, решение исходной задачи сводится к решению системы двух дифференциальных уравнений dr rfO v Г которая в свою очередь может быть сведена к одному урав- dr dO пению — = —тггг исключением переменной t, г V 3 F Решая последнее дифференциальное уравнение, полу- получаем, что где С — произвольная постоянная. Учитывая теперь, что миноносец начинает движение вокруг полюса О с полярной оси г на расстоянии х миль от точки О, т. е. учитывая, что г=1 при 6=0 и г=3 при 0=—я, приходим к выводу, что в первом случае С—\, а во втором С — Зел^ 8. Таким образом, чтобы выполнить свою 34
задачу, миноносец должен пройти две или шесть миль прямым курсом по направлению к месту обнаружения под- подводной лодки, а затем двигаться либо по спирали r^=e°^v 8» либо по спирали /■= Зе<0+я^у я. Модели боевых действий Во время первой мировой войны английский инженер и математик Ф. У. Ланчестер построил несколько математи- математических моделей ведения воздушных сражений. Затем эти модели были обобщены и распространены на случаи боевых действий регулярных войск, партизанских соединений, а также как тех, так и других одновременно. Эти три модели мы и рассмотрим ниже. Итак, пусть в боевых действиях участвуют две противо- противоборствующие стороны х и у. Их численный состав в мо- момент времени t, где / измеряется в днях, начиная с первого дня боевых действий, обозначим через x(i) и y(t) соответст- соответственно. В дальнейшем именно численность сторон будет играть определяющую роль в построении упомянутых выше моделей. Дело в том, что практически трудно указать критерии, которые учитывали бы при сравнении противо- противоборствующих сторон еще и степень боевой готовности и вооруженности, уровень и опыт командирского состава, моральный дух и многие, многие другие факторы. Мы предположим также, что x(t) и y(t) изменяются не- непрерывно и, более того, что они дифференцируемы как функции времени. Конечно, эти предположения являются упрощением реальной ситуации, поскольку x{t) и y(t) — это целые числа. Но вместе с тем понятно, что при достаточ- достаточно большом численном составе каждой из противоборствую- противоборствующих сторон увеличение численности на одного или двух человек дает с практической точки зрения бесконечно ма- малую величину по сравнению с уже имеющимся наличным составом. Поэтому можно считать, что за малые промежутки времени численный состав также изменяется па малые количества (не целые). Этих соглашений, конечно, недо- недостаточно для того, чтобы выписать конкретные формулы для x(t) и y(t) как функций t. Однако можно указать ряд факторов, которые позволяют описать скорость изменения численности противоборствующих сторон. Именно, обозна- обозначим через OLR величину, выражающую скорость,4 с которой сторона х несет потери от болезней и других факторов, не связанных с непосредственными боевыми действиями. Далее, пусть CLR — скорость, с которой сторона л; несет 2* 35
потери от непосредственных столкновений в процессе бое- боевых действий со стороной у. Наконец, через RR обозначим скорость подхода подкреплений к силам стороны х. Тогда ясно, что общая скорость изменения x(t) задается уравне- уравнением *LJ2 = _ (OLR + CLR) + RR. D7) Аналогичное уравнение имеет место и для y{t). Задача те- теперь состоит в том, чтобы указать соответствующие формулы для введенных выше величин OLR, CLR и RR и затем ис- исследовать построенные дифференциальные уравнения. Полученные выводы позволят ответить на вопрос о вероят- вероятном победителе. В дальнейшем используем следующие обозначения- а, Ь, с, d, g, h — неотрицательные постоянные, характеризую- характеризующие степень влияния различных факторов на потери в живой силе обеих сторон х и у; P(t) и Q(t) — члены, учитывающие возможность подхода подкреплений силам л: иг/ в течение дня; Хо, у0 — численный состав сил х и у перед началом боевых операций. Выпишем три модели *), построенные Ланчестером. Первая из них относится к описанию боевых действий между регулярными войсками, и она имеет вид В дальнейшем эту систему будем называть дифференциаль- дифференциальной системой типа (А). Вторая модель описывает боевые действия между партизанскими соеди- соединениями. Ее мы будем называть дифференциальной систе- системой типа (В). Наконец, третья модель, которую будем называть дифференциальной системой типа (С) и которая *) В статье Coleman С. S. «Combat models» (Differential equa- equation models.— New York e.a., 1983.-— P. 109—131) приведены конкрет- конкретные примеры боевых действий, причем показано, как они согласуются с рассматриваемыми дифференциальными моделями. 36
имеет вид dJLM = _ ах @ -gx (t) у (t) + P (t), dy (t) ^ (L\—d (t\ + O(t\ описывает смешанный тип боевых действий, в которых при- принимают участие как регулярные части, так и партизанские соединения. Каждое из выписанных выше дифференциальных урав- уравнений выражает скорость изменения численного состава противоборствующих сторон в зависимости от действия различных факторов и имеет вид D7). Потери в живой силе, которые не связаны с непосредственными боевыми действия- действиями и которые определяются членами —ax{f) и —dy(t), дают возможность выписать постоянные относительные скорости потерь (в отсутствие боевых действий и подкреплений) по- посредством уравнений ~хЖ = ~~а' ~у~Ш = Если в моделях Ланчестера присутствуют лишь члены, соответствующие только подкреплениям и потерям, не связанным с боевыми действиями, то это означает, что пос- последние вообще отсутствуют. Наличие же членов —by(t), —cx(t), —gx(f)y(l) и —hx{t) y(t) означает уже, что боевые действия имеют место. Рассматривая дифференциальную систему типа (А), предположим, во-первых, что каждая из противоборствую- противоборствующих сторон находится в зоне действия огневых средств другой стороны, а во-вторых, что огонь ведется только по живой силе, непосредственно участвующей в боевых действи- действиях. При таких предположениях Ланчестер предлагает для подразделений регулярных войск стороны х ввести член —by(t), который отражал бы боевые потери. Коэффициент b будет тогда указывать на эффективность боевых действий стороны у. Таким образом, уравнение показывает, что постоянная b — это единица измерения средней эффективности каждой единицы боевых сил сторо- стороны у. Подобное объяснение можно дать и члену —cx(t). Понятно, что совсем непросто вычислить коэффициенты эф- эффективности b и с. Один из путей здесь — рассматривать 37
эти коэффициенты в виде Ь = гуру, с=гхрх, D8) где гv и гх — коэффициенты огневой мощи сторон у и х соответственно, а ру и ру — это вероятности того, что каж- каждый из выстрелов со стороны у и х соответственно окажется метким. Отметим далее, что члены, которые соответствуют бое- пым потерям в дифференциальной системе типа (А), яв- являются линейными. Что же касается систем типа (В), то здесь аналогичные члены уже являются нелинейными и объясняется это следующим. Пусть боевые силы партизан численностью x{t) человек занимают некоторую террито- территорию R, .оставаясь невидимыми для противника. И хотя противник держит под огнем территорию R, он не может знать эффективности своих действий. При этом весьма прав- правдоподобно, что потери партизанских подразделений х пропорциональны, с одной стороны, числу х{1) на R, а с другой — числу y(t) боевых сил противника. Таким обра- образом, член, который соответствует потерям партизанских соединений х, имеет вид —gx(t)y(t), где коэффициент эффек- эффективности боевых действий стороны у, вообще говоря, более труден для оценки, чем коэффициент Ъ в первом равенстве D8). Однако для его определения мы можем использовать коэффициент огневой мощи гу, а также принять во внима- внимание рассуждения Ланчестера, согласно которым вероят- вероятность меткого выстрела со стороны у прямо пропорциональ- пропорциональна так называемой территориальной эффективности Arlt одного выстрела со стороны у и обратно пропорциональна площади Лх территории R, занимаемой силами х. Через Агу здесь обозначается площадь, которую занимает один из партизан. Таким образом, вероятные формулы для опре- определения g и h таковы- *-'»Зг- h-r*%- <49> Ниже рассмотрим более подробно каждую из приведен- приведенных выше трех дифференциальных моделей. Случай А. Дифференциальные системы типа (А). Квадратичный закон. Предположим, что регулярные вой- войска двух противостоящих сил ведут боевые действия в той упрощенной ситуации, когда потери, не связанные с такими действиями, отсутствуют. А тогда если обе стороны не полу- получают еще и подкреплений,то математическая модель сводит- 38
ся к дифференциальной системе dx dt = — by, dy_ dt = — ex. E0) Разделив второе уравнение на первое, получим, что dy сх dx by E1) Интегрируя дифференциальное уравнение E1), приходим к равенству b[y4t)-yl] = c[x*(t)-xl]. E2) Соотношение E2) и объясняет, почему система E0) соответ- соответствует модели с квадратичным законом. Если обозначить через К постоянную Ьу\— сх\, то уравнение Ьу2—сх2-^К, E3) полученное из равенства E2), задает гиперболу (пару пря- прямых, если К—О), и мы можем более точно классифицировать систему E0). Именно такую систему можно назвать диффе- дифференциальной системой с гиперболическим законом. Л" >О: и побеждает равиобесие побеждает so(t) Рис. 12 На рис. 12 изображены гиперболы для различных значе- значений /С, при этом по очевидным соображениям рассматри- рассматривается лишь первый квадрант (х^О, у^О). Стрелки на кривых указывают направление изменения численности сил с течением времени. Чтобы ответить на вопрос, кто побеждает в построенной модели E0), условимся, прежде всего, говорить, что побеж- побеждает сторона у (соответственно х), если она первой уничто- уничтожает боевые силы стороны х (соответственно у). Так, в 39
нашем случае побеждает сторона у, если /С>0, ибо в соот- соответствии с уравнением E3) переменная у никогда не может обратиться в нуль, в то время как при значении y(t)— VKlb переменная х обращается в нуль. Таким обра- зОхМ, чтобы победили силы у, им необходимо стремиться достичь такой ситуации, при которой К > 0, т. е. когда Ьу\ > сх\. E4) Из равенств D8) следует, что неравенство E4) можно пере- переписать в виде (МЛ* >!*£*. E5) Левая часть неравенства E5) показывает, что изменения в отношении сил yjx» дают преимущество одной из сторон в соответствии с квадратичным законом! Так, например, изменение в отношении сил от yjxu=\ до yjxu=2 дает четырехкратное преимущество силам у. Отметим далее, что уравнение E3) определяет соотношение между силами про- противоборствующих сторон, не учитывая явно время. Чтобы вывести формулы, которые учитывали бы явно и время, поступим следующим образом. Продифференцируем по t первое из уравнений системы E0), а затем воспользуемся вторым уравнением этой же системы. В результате придем к дифференциальному уравнению Ьсх = 0. E6) Используя тогда в качестве начальных условий соотношения решение уравнения E6) получим в виде #@ = *ocos^—Y^oSinp/, E7) где р = Vbc, у = Vbjc. Подобным образом получаем, что ipf. E8) yQ$sinpf. На рис. 13 показаны графики функций, заданных урав- уравнениями E7) и E8) в случае, когда К > 0 (т. е. когда Ьу£> >cxl или когда ууо>хо). В заключение отметим, что для победы сил стороны у не обязательно, чтобм число #0 было больше числа х0. Тре- Требуется лишь выполнение неравенства ууо>хо. 40
Случай В. Дифференциальные системы типа (В). Линейный закон. Уравнения динамики, моделирующие боевые действия двух противоборствующих сторон, могут быть легко решены, если, как и в предыдущем случае, х ос, Рис. 13 исключаются потери, не связанные с боевыми действиями, и ни одна из сторон не получает подкреплений. При таких ограничениях дифференциальная система типа (В) прини- принимает вид ^■^—gxy, ^ = —hxy. E9) Разделив второе уравнение из системы E9) на первое, полу- получим простое уравнение dx~ g • интегрирование которого приводит к равенству yo] = h[x(t)-xo]. F0) Линейная зависимость, выраженная соотношением F0), и объясняет, почему нелинейная система E9) соответствует модели с линейным законом ведения боя. Равенство F0) можно переписать в виде gy-hx = L, F1) где L=gy0 — hx0. Отсюда, в частности, следует, что если L>0, то в результате боевых действий побеждает сторона у. Если же L < 0, то победителем будет сторона х 41
На рис. 14 дана геометрическая интерпретация линейной функциональной зависимости F1) при различных значе- значениях L. Рассмотрим более подробно ситуацию, при которой побеждает одна из сторон. Пусть таковой будет сторона у(Щ у побеждает рабновесие О x(t) у. Тогда, как мы уже знаем, должно выполняться нера- неравенство gy0—/глго>-О, которое перепишем в виде Если же теперь обратиться к формулам D9), то условие победы стороны у выразится в виде неравенства М *0 г„А, Таким образом, стратегия сил у состоит в том, чтобы сделать максимальным отношение yjxu и минимальным отношение AJAy. С практической точки зрения неравенство F2) лучше рассматривать в виде > Итак, мы видим, что произведения Ауу0 и Ахх0 являются в определенном смысле критическими величинами. Отметигм, наконец, что на основании равенства F1) из системы E9) нетрудно вывести формулы, показывающие, как силы обеих сторон изменяются с течением времени. Случай С. Дифференциальные системы типа (С). Параболический закон. В модели (С) партизанские силы противостоят регулярным частям. Мы снова сделаем упро- упрощающие предположения, заключающиеся в том, что обе 42
противоборствующие стороны не обеспечиваются подкрепле- подкреплениями и не несут потерь, не связанных с непосредствен- непосредственными боевыми действиями. В этом случае имеем дифферен- дифференциальную систему вида dx dy_ dt = — ex, F3) где x(t) — силы партизан, a y(t) — силы регулярных войск. Разделив второе уравнение системы F3) на первое, получим дифференциальное уравнение dy с Интегрируя его в соответствующих пределах, приходим к соотношению gtf{t) = 2cx(l) + M, F4) где M=gyl — 2сх0. Таким образом, дифференциальная система F3) соответствует модели с параболическим зако- законом ведения боевых действий. При этом, если УИ<СО, то побеждают партизаны, а если М ~> 0, то они терпят пора- поражение. у(Щ s М>о: у побеждает равнодеоие -fi/Zc x(t) Рис. 15 На рис. 15 схематически изображены параболы, опреде- определяемые уравнением F4) при различных значениях М. Опыт показывает, что регулярные части могут нанести пора- поражение партизанским соединениям только в случае, когда отношение yjxu значительно больше единицы. Основываясь же на параболическом законе ведения боевых действий и учитывая условие М > 0, приходим к выводу, что победа регулярных сил будет гарантирована, если выполняется 43
/ у \2 2с 1 неравенство ( — } > , которое с учетом формул D8) \ хо / 8 хо и D9) можно переписать в виде Уо_У- >> о Тх АхРх 1 хо ) ry rv х° Почему маятниковые часы не являются точными? Чтобы ответить на поставленный вопрос, рассмотрим идеализированную модель маятниковых часов, состоящую из стержня длиной / и гири массой m на его конце (масса стержня предполагается такой, что ее можно не принимать в расчет по срав- сравнению с массой гири) (рис. 16). Если гирю отклонить на угол а и затем отпустить, то в соответствии с зако- законом сохранения энергии —2~ = mg(lcos В—/cosa), F5) где v — скорость движения гири, a g— ускорение силы тяжести. Рассматривая только малые отклонения гири от поло- положения равновесия, всегда можно считать, что длина дуги, по которой гиря отклоняется от положения равновесия на угол 0, определяется равенством s=/0. А в этом случае v==~7f — l-jr и соотношение F5) приводит к дифференци- дифференциальному уравнению fJ F6) Учитывая теперь, что 0 убывает с возрастанием / (для ма- малых t), уравнение F6) можно переписать в виде —cosa А тогда если Т— период колебаний маятника, то 2g J l^cosO —cos a a ИЛИ a m o У СОа 44
Как видно из последней формулы, период колебаний маят- маятника зависит от угла а. Этот факт и является основной причиной того, что маятниковые часы не точные, ибо прак- практически гиря всякий раз отклоняется в крайнее положение на угол, отличный от а. Обращаясь к формуле F7), заметим, что ее можно пере- переписать и в более простом виде. Действительно, так как cos 0=1—2 sin2 у, cosct—1— TO F8) где k — sin -j • Заменим теперь переменную 6 на переменную ф, полагая sin -s- = k sin ф . Из последнего равенства следует, что когда 6 возрастает от 0 до а, то ср возрастает от 0 до л/2, причем Y cos -j d6 = k cos ф dy, или cos! 1^ 1-й sin? q> Последнее соотношение дает возможность формулу F8) переписать в виде - Г -. d(p =4 I/' - F(k, я/2), g J /1-/2? sin? ф У S ' где функцию о ' " ф называют эллиптическим интегралом первого рода, в отличие от эллиптических интегралов второго рода, имею- имеющих вид ф Е (k, Ч>) = \ Vl— /^sin^q) dq>. о 45
Эллиптические интегралы не могут быть вычислены в элементарных функциях, и поэтому дальнейшее обсуждение задачи о маятнике мы свяжем с другим подходом, который будет рассматриваться при исследовании консервативных систем в механике. Здесь мы лишь отметим, что исходным пунктом в дальнейших исследованиях будет дифференциаль- дифференциальное уравнение которое получается из уравнения F6) дифференцированием по t. Циклоидальные часы Мы выяснили, что часы с обыкновенным (круговым) маятником не могут идти точно. Поэтому возникает есте- естественный вопрос: существует ли какой-либо другой маят- маятник, время качания которого не зависело бы от размаха? Впервые эта задача была поставлена и решена еще в XVII в. Ниже приводится ее решение, но предварительно обратимся к выводу уравнения одной замечательной кривой, которую Галилей назвал циклоидой (от греческого kykloeides — кругообразный, круглый). Под циклоидой понимается плоская кривая, представляю- представляющая собой траекторию точки, лежащей на окружности круга (называемого производящим кругом), катящегося без скольжения по прямой линии. Пусть этой прямой линией является ось абсцисс (рис. 17), а радиус производящего круга равен г. Предположим, что в начальном положении вычерчивающая циклоиду точка на- находится в начале координат, а после того как круг повер- повернется на угол 0, она займет положение М. Тогда, исходя из геометрических построений, получаем, что x^OS^OP—SP, y^MS-^CP-CN. Но Поэтому параметрически циклоида задается уравнениями x=r(G—sin G), r/=r(l— cos 9). F9) Если из уравнений F9) исключить параметр 6, то в пря- прямоугольной системе координат Оху уравнение циклоиды примет вид #= г arccos —-—]/г2гу—уй. 46
Из самого способа образования циклоиды следует, что ока должна состоять из конгруэнтных арок, каждая из кото- которых соответствует полному обороту производящего круга *). Отдельные арки соединяются в точках, в которых они имеют общую вертикальную касательную. Эти точки, называе- называемые точками возврата циклоиды, соответствуют самым Рис. 17 низким положениям той точки на окружности катящегося круга, которая описывает циклоиду. Самые высокие поло- положения этой точки находятся посредине между точками воз- возврата и называются вершинами циклоиды. Отрезок же пря- прямой линии между двумя соседними точками возврата, длина которого равна 2пг, называется основанием арки циклоиды. Для циклоиды имеют место, например, такие свойства! а) плошддь, ограниченная аркой циклоиды и ее основа- основанием, равна утроенной площади производящего круга {тео- {теорема Галилея); б) длина одной арки циклоиды равна четырем диамет- диаметрам производящего круга (теорема Рена). Последний результат особенно неожиданный: ведь даже для вычисления длины такой простой кривой, как окруж- окружность, пришлось вводить иррациональное число зх, вычис- вычислить которое не так просто, а длина арки циклоиды выра- выражается через радиус (диаметр) целым числом! Циклоида обладает и многими другими интересными свойствами, имею- имеющими исключительное значение для физики и техники. В частности, профили зубьев шестерен, очертание многих типов эксцентриков, кулачков и иных деталей машин имеют форму именно циклоиды. Обратимся теперь к задаче, решение которой позволило нидерландскому ученому X. Гюйгенсу сконструировать в *) Много интересного о циклоиде и родственных ей кривых можно почерпнуть, например, из книги Г. Н. Бермана «Циклоида» (М.: Наука, 1980). 47
1673 г. точные часы. Эта задача заключается в построении ъ вертикальной плоскости такой кривой, чтобы время, необ- необходимое для спуска по ней до фиксированного горизонта тяжелой материальной точки, находящейся в начальный момент времени t~-t0 в состоянии покоя, не зависело от исходного положения точки на этой кривой. Как показал Гюйгенс, такой изохронной (от греческого isos — равный, одинаковый и chronos — время), или таутохронной, кривой оказалась циклоида. Рис. 18 Практически решение этой задачи можно провести сле- следующим образом. Предположим, что на торцовой стороне вертикально поставленной доски вырезан желоб в форме циклоиды (рис. 18). Не учитывая трение, попытаемся опре- определить время, в течение которого металлический шарик ска- скатится из точки М до низшей точки К вырезанного желоба. Пусть х0, у о — координаты исходного положения ша- шарика, т. е. точки М, а б0 — соответствующее ей значение параметра. Когда шарик скатится из положения М в не- некоторое положение N @), то он снизится по вертикали на расстояние /г, которое в силу второго уравнения F9) опре- определится следующим образом: h = y—у0 = г A — cos 6)—г(\ — cos 0О) = г (cos Go—cos 0). Мы знаем, что скорость падающего тела находится по формуле где g — ускорение силы тяжести. В нашем случае последняя формула примет вид v = (cos 00—-cos 0). С другой стороны, так как скорость — это производная пути s по времени Т, то приходим к соотношению ds_ dT = J/r2gr(cos0o—cosO). 48
Для циклоиды ds = 2r sin-g-dB, и поэтому последнее диф- дифференциальное уравнение эквивалентно такому: 2r sin -5- № Интегрируя это уравнение в соответствующих пределах, получаем, что " 2rsin~du г—* dcos-|- ^_ Г z 2 1/ — \ J /gr(cosB0—cosG) ' g J -1/ 2 6о ,6 о 6 о 1/ cos? -у— cos2 д- Таким образом, время Т, в течение которого шарик скатится из положения М в положение К, определяется по формуле показывающей, что период Т не зависит от 0О, т. е. пе- период не зависит от исходного положения М шарика, и поэтому очевидно, что два шарика, начавшие одновременно скатываться по желобу из точек М и N, окажутся в точке /( в один и тот же момент времени. Так как мы условились трение не учитывать, то шарик, скатившись в точку /С, будет по инерции продолжать дви- движение и через указанный последней формулой промежуток времени поднимется до точки Mi, находящейся на одной высоте с точкой М. Проделав затем обратный путь, шарик и осуществит движение так называемого циклоидального ма- маятника с периодом колебания Гв=4я1/7/ё. G0) Отличительным свойством циклоидального маятника но сравнению с обычным (круговым) является то, что период его колебаний не зависит от амплитуды. Покажем теперь, как можно заставить обычный круговой маятник двигаться таутохронно, не прибегая к желобу и подобным ему приспособлениям с большим трением. Для этого достаточно изготовить шаблон (например, из дерева), состоящий из двух одинаковых полуарок циклоиды, имею- имеющих общую точку возврата (рис. 19). Шаблон укрепляется вертикально, а в точке возврата О привязывается шарик 49
на нити, длина которой равна удвоенному диаметру произ- производящего круга циклоиды. Шарик, отклоненный в произвольную точку М, начнет совершать колебания, период которых не будет зависеть от выбора точки М. Даже если под влиянием трения и сопро- сопротивления воздуха размах колебаний будет уменьшаться, Рис. 20 время колебания маятника останется неизменным! Для кругового маятника, движущегося по дуге окружности, свойство изохронности практически приближенно удовлет- удовлетворяется для небольших амплитуд, когда дуга окружности незначительно отклоняется от дуги циклоиды. Рассмотрим в качестве примера малые колебания маят- маятника по дуге А В циклоиды (рис. 20). Если эти колебания очень малы, то влияние направляющего шаблона практи- практически не ощущается и маятник раскачивается почти как обыкновенный круговой маятник с длиной нити, равной 4г. Путь Л В циклоидального маятника практически не бу- будет отличаться от пути СЕ кругового маятника с длиной нити, равной 4г. Поэтому период малых колебаний кругового маятника с длиной нити /=4г практически одинаков с пе- периодом колебаний циклоидального маятника. Если теперь в формуле G0) положить /-=//4, то период малых колебаний кругового маятника можно выразить через длину его нити /: Последняя формула выводится (исходя из других сообра- соображений) в школьном курсе физики. В заключение заметим, что циклоида является единст- единственной кривой, двигаясь по которой тяжелая материальная точка совершает изохронные колебания. 50
Задача о брахистохроне Задача о брахистохроне (от греческого brachistos—крат- brachistos—кратчайший и chronos — время), или о кривой скорейшего спус- спуска, была поставлена швейцарским математиком И. Бер- нулли в 1696 г. и заключалась в следующем. Ё вертикальной плоскости даны две точки Л и В {рис. 21), не лежащие на одной вертикали. Требуется среди всех кривых, проходящих через указанные точки, найти такую, спускаясь по которой под действием силы тяжести, материальная точка скатится из точки А в точку В за кратчайшее время. Рис, 21 Рис. 22 Эта задача решалась лучшими математиками разных времен и была решена как самим И. Бернулли, так и Г. В. Лейбницем, И. Ньютоном, Г. Лопиталем, Я. Бернул- Бернулли. Задача знаменита тем, что, помимо значимости се реше- решения с точки зрения естественнонаучной, она стала источни- источником идей совершенно новой области математики — гариа- ционного исчисления. Решение задачи о брахистохроне можно связать с реше- решением другой задачи, возникающей ё оптике. Обратимся к рис. 22, на котором схематически изображен луч света, падающий из точки А в точку Р со скоростью vJ} а затем проходящий в более плотной среде от точки Р до точки В с меньшей скоростью vz. Полное время Т, требуемое для прохождения луча света от точки А до точки В, определя- определяется, очевидно, из равенства -1 Если предположить, что луч света проходит от точки А до точки В по указанному пути за минимальное время Т, то 51
AT производная ^- = 0. А тогда x с—х ^ ^""" VST\—a:)? ' или ^^ = ^^ . Последнее равенство выражает извест-' ный закон преломления Снеллиуса, который был первона- , sin «i чально открыт экспериментально в форме ginct ==a* гДе а~ const. Высказанное предположение о том, что луч света про- проходит от точки Л до точки В за минимальное время, назы- называется принципом Ферма наименьшего времени. Значение этого принципа заключается не только в том, что он дает рациональную основу для вы- вода закона Снеллиуса, но, в ^|частности, и в том, что он мо- vz Г^2Ч<^_ жет быть применен для нахож- ~ л Лй"~ ~™ дения траектории луча све- s^Ята> когДа последний проходит через среду переменной плот- ности, вообще говоря» не по отрезкам прямых. Рис- 2 Для примера рассмотрим рис. 23, на котором изобра- изображена слоистая среда. В каждом отдельном слое скорость света постоянна, но она убывает при переходе от вышележа- вышележащего слоя к нижележащему. Падающий луч света при пере- переходе от слоя к слою преломляется все больше и больше по направлению к вертикали. Применив закон Снеллиуса к границам между слоями, получим sin at sincc2 sinoc3 slna4 Щ v2 vs v4 Допустим теперь, что толщина слоев неограниченно умень- уменьшается, а число слоев неограниченно растет. Тогда в пре- пределе скорость света убывает непрерывно, и мы заключаем (рис. 24), что Й52-О. G1) где a=const. Подобная ситуация наблюдается (с определен- определенными оговорками) в случае падающего на землю луча солн- солнца, когда он замедляет скорость при прохождении через атмосферу с возрастающей плотностью. 52
Обращаясь теперь к задаче о брахистохроне, введем си- систему координат в вертикальной плоскости так, как это по- показано на рис. 25, и представим себе, что шарик (подобно лучу света) способен выбрать себе такую траекторию спуска из точки А в точку В, чтобы время спуска было наименьшим. Рис. 24 Рис. 25 Тогда, как это следуег из предыдущих рассуждений, имеет место формула G1). Исходя из принципа сохранения энергии, получаем, что скорость, достигнутая шариком на заданном уровне, зави- зависит только от потери потенциальной энергии при достижении этого уровня, а не от вида траектории, по которой скаты- скатывается шарик. Это означает, что v = ]^2gy. С другой стороны, геометрические построения показы- показывают, что sin a = cos р — 1 1 seep Комбинируя тогда равенство G1) с последними двумя соот- соотношениями, получаем У{1+(У'П = С. G2) Уравнение G2) и есть дифференциальное уравнение брахи- брахистохроны. Покажем сейчас, что брахистохроной может быть циклоида, и только она. Действительно, так каку' — -^, ах то, разделив переменные в уравнении G2), придем к урав- уравнению вида 53
Введем теперь новую переменную ср, полагая Тогда у = С $>m2q>, dy — 2С sin ф cos (р tftp, dx = tg ф dy = 2С sin2 ф dtp = С A —cos 2ф) dcp. Интегрирование последнего уравнения приводит к соотно- соотношению Q х=~2 Bф—sin 2ф) -f С{, где, в силу начальных условий, х~у~0 при ф=0 и постоян- постоянная Ci=0. Итак, С С я Bф81п2ф) у С з5п2ф A Полагая здесь С/2=г, 2ф=0, придем к стандартным пара- параметрическим уравнениям циклоиды F9). Удивительная кривая циклоида — она является не только изохроной, но и брахистохроной! Среднее арифметическое, среднее геометрическое и дифференциальное уравнение Рассмотрим следующую любопытную задачу, которая впервые была решена немецким ученым К. Ф. Гауссом. Пусть Шо и п0 — два произвольных положительных числа (то>/го). Составим из них два других числа т1 и пи являющихся для т0 и п0 соответственно средним арифмети- арифметическим и средним геометрическим. Другими словами, поло- положим Поступая с mi и «i так же, как с т0 и щ, полагаем Продолжая этот процесс и дальше, получаем две числовые последовательности {mft}, {nk} (k = 0, 1, 2, . . .), которые, как нетрудно показать, являются последовательностями сходящимися. Ставится вопрос: чему равна разность пре- пределов этих последовательностей? 54
Приведем принадлежащее немецкому математику К. В. Борхардту изящное решение этой задачи, связанное с построением линейного дифференциального уравнения вто- второго порядка. Итак, пусть а — искомая разность. Очевидно, что а зависит от т0 и па. Аналитически этот факт выра- выразим записью a—f(m0, л0), где f — некоторая функция. Из определения числа а вытекает также, что a=f(m.i, fit). Теперь, если умножить т0 и п0 на одно и то же число k, то каждое из введенных чисел ти пи т2, п2, . . ., не исключая и а, также приобретет множителем число k. Последнее означает, что а является однородной функцией первой степени относительно т0 и п0 и, следовательно, a = mof(\, no//7zo) = mj(l, ihlmt). Полагая по/то=х, nJmt—Xi, ... и обозначая 1//A, п0/т0) через у, а 1//A, njm^ через t/i, ..., находим Так как х± связано с х уравнением то получаем, что dxj __ 1-х __ (xi—xf)(l+x)*. dx ~~A+х)*У"~х 2 (х-xs) Уравнение же G3) приводит к следующему соотношению: dy __ 2 2 dy± dxj dx~~ A+д:)?^1'' l+ж dxi dx ' Заменяя здесь производную -~ ее значением из предыду- . .-JL (IX щего равенства и освобождаясь от знаменателя х—Xs, получаем соотношение Продиф(}зеренцировав обе части этого равенства по х, най- найдем, что d\x{x~ dx\ _п а\ l+x dx *yi dx 2x (x— 1) 1+* Л 1 (x dxi L dyi dx± \ — dxf , dx "T" dx ' 55
Элементарными преобразованиями последнее равенство приводится к такому: Если в последнем уравнении заменить х на хх, то хг перей- перейдет в х2. Если затем заменить хх на xit то #2 перейдет п х3 и т. д. Поэтому, полагая придем к равенству а* (и) — ^~~х — ки) A+х)\Гх Если теперь устремить п к оо, то 1 — хп будет стремиться к 0, и, таким образом, придем к тому, что а* (у) = 0. А это означает, что у удовлетворяет дифференциальному уравнению Заметив теперь, что нетрудно указать и значение этого числа. Действительно, исходя из того, что у обязано быть постоянным решением уравнения G4), находим, что таковым является только #=0. Таким образом, разность пределов последовательностей {mh} и {пк} равна 0. Дифференциальное уравнение G4) интересно не только тем, что оно позволило исходную задачу свести к очевидной, но и тем, что это уравнение имеет непосредственную связь с решением задачи о вычислении периода малых колебаний кругового маятника. Как было показано, период малых колебаний кругового маятника находится по формуле (k, я/2), 56
где л/2 F(k, я/2)= Так вот, оказывается, что если О^/г<С1, то Г ^Ф где „2 2?-4?-62 ... B/г)? есть решение дифференциального уравнения G4). О полете тела, брошенного под углом к горизонту Пусть тело брошено под углом а к горизонту с началь- начальной скоростью v0. Требуется вывести уравнения движения тела, пренебрегая силами сопротивления. У Рис. 26 Выберем оси координат так, как показано на рис. 26. В произвольном положении М на тело массой т действует лишь одна сила — его вес P=mg. Поэтому в соответствии со вторым законом Ньютона дифференциальные уравнения движения в проекциях на оси хну запишутся в виде Сокращая на т, получаем уравнения dt*. Начальные же условия движения тела таковы: х = 0, у —0, -~- = i>0 sin a при ^ = G5) ). G6) 57
Интегрируя уравнения G5) с учетом начальных условий G6), приходим к выводу, что уравнения движения тела за- задаются формулами х = (и0 cos a) t, y — (vos\na) t—gt2j2. G7) Из уравнений G7) можно сделать ряд выводов о харак- характере движения брошенного тела. Например, можно отве- ответить на вопросы: каково время полета тела до его падения на землю, каково расстояние полета по горизонтали, какова максимальная высота летящего тела, какова траектория полета тела? На первый вопрос можно ответить, найдя значение t, при котором у—О. Из второго уравнения G7) видно, что это будет тогда, когда т. е. когда либо t—О, либо t= By0sin a)/g. Второе значение и дает ответ на первый вопрос. Чтобы ответить на второй вопрос, вычислим значение х при значении t, равном времени полета. Из первого урав- уравнения G7) получаем, что расстояние полета по горизонтали задается формулой (t'o cos a) Bi>0 sin cc) v% sin 2a _ _ . Из последнего равенства, в частности, следует, что расстоя- расстояние будет наибольшим, когда 2а=90°, т. е. а=45°. В этом случае расстояние будет равно v\lg. Ответ на третий вопрос мы сразу же получим, если ука- укажем условие максимальности у. Но это означает, что в той точке, где у максимально, производная -тг—®- Учитывая же, что ~ = —g/ + aosina, приходим к равенству —gt+v0 sin a=0, откуда t = {v0 sin a)/y. Подставляя теперь полученное значение t во второе урав- уравнение G7), находим, что максимальная высота равна v\ sin2a/Bg). На четвертый вопрос ответ уже получен выше. А именно, траекторией полета будет парабола, так как уравнения G7) — это параметрическое задание параболы, которая в 58
декартовых координатах запишется в виде nv2 Невесомость Состояние невесомости может быть достигнуто различ- йыми способами, хотя оно (вольно или невольно) и ассоци- ассоциируется с «плаванием» космонавтов в кабине космического корабля. Рассмотрим следующую задачу. Пусть человек весом Р находится в кабине пассажирского лифта, движущегося вниз с ускорением (a=ag, где 0<Ca<C\,ag — ускорение силы тяжести. Требуется определить давление человека на дно кабины, а также уско- ускорение лифта, при котором это давление отсутст- отсутствует. На человека в лифте действуют две силы (рис. 27): вес человека Р и сила реакции дна кабины лифта Q (равная давлению человека на дно кабины). Дифференциальное уравнение дви- движения человека запишется в виде Так как ~тк" — g> —a£> m=P/g, то из уравнения G8) полу- получаем соотношение Q = P-m^--P(l-a). G9) Принимая, далее, во внимание, что 0<Са<С1, приходим к выводу, что Q<CP. Итак, давление человека на дно кабины движущегося вниз лифта задается силой <Э = ЯA—а). В том же случае, когда кабина лифта поднимается вверх с ускорением (o=ag, где 0<а<С1, давление человека на дно кабины определяется силой Q=P A+a). Выясним теперь, при каком ускорении лифта давление человека на дно каби- кабины отсутствует. Для этого достаточно обратиться к равенст- равенству G9), положив Qt=0. В результате приходим к выводу, что <х=1, т. е. для того, чтобы Q=0, ускорение лифта обязано равняться ускорению силы тяжести. Таким образом, при свободном падении лифта с ускоре- ускорением силы тяжести g давление человека на дно кабины лифта 69
отсутствует. Именно это состояние и называют состоянием невесомости. При невесомости отсутствуют взаимные давле- давления отдельных частей тела человека. Это вызывает у него необычные ощущения. При состоя- состоянии невесомости все точки тела име- имеют равные ускорения. Конечно, состояние невесомости может иметь место и не только при свободном падении. Для иллюстра- иллюстрации рассмотрим такую задачу. Какова должна быть скорость космического корабля, движущего- движущегося вокруг Земли как искусствен- искусственный спутник, чтобы человек нахо- находился в кабине в состоянии невесо- невесомости? Задачу будем решать в предпо- предположении, что космический корабль движется по круговой орбите ра- радиуса r+Л, где г — радиус Земли, a h — высота полета космического корабля над поверхностью Земли. Из преды- предыдущей задачи следует, что в состоянии невесомости давле- давление на дно кабины космического корабля, а следователь- следовательно, и реакция Q дна кабины равны нулю. Итак, Q=0. Об- Обратимся теперь к рис. 28. Ось х направлена вдоль главной нормали п круговой траектории космического корабля. Воспользовавшись дифференциальным уравнением движе- движения материальной точки в проекции на главную нормаль, которое имеет вид п tnv2 Рис. 28 и гдер = г + /1, 2 Fkn — F' причем сила F направлена по fei главной нормали к траектории движения, придем к урав- уравнению ri2 _ //2 у mv 7+h или к уравнению v2 (I2 К Подставив значение -—• из последней формулы в соотно- 60
шеиие G8), получим равенство Здесь сила Р равна силе F притяжения к Земле, которая в соответствии с всемирным законом тяготения обратно про- пропорциональна квадрату расстояния r-\rh до центра Земли, т. е. Р km где т — масса космического корабля, а постоянная k опре- определяется из следующих соображений. На поверхности Зем- Земли, т. е. при /г=0, сила притяжения F=mg. Из предыдущей формулы тогда получаем, что k=gr2 и, таким образом, где g — ускорение силы тяжести на поверхности Земли. Подспавив теперь полученное значение Р в формулу (80) и учтя тот факт, что Q=0, приходим к выводу, что тре- требуемая в задаче скорость находится из равенства Законы . Кеплера движения планет В соответствии с законом всемирного тяготения лю- любые два тела, находящиеся на расстоянии г друг от друга и имеющие массы т и М соответственно, притягиваются с силой (8D У\ где у — универсальная пос- \ |< ^%т тоянная притяясения. Основываясь на этом законе, опишем движение . ^ . планет, считая, что т — I /У \У масса планеты, движущей- . . ся вокруг Солнца, М—мае- ^ уР са Солнца. Влияние других ^0 х ~Р планет при этом учитывать Рис. 29 не будем. Пусть (рис. 29) Солнце находится в начале координатной системы Оху, а планета в момент времени / находится в точке 61
А о текущими координатами х, у. Сила притяжения F, действующая на планету, раскладывается на две составляю- составляющие: параллельную оси х, равную F cos ср, и параллельную оси у, равную F sin ф. Используя формулу (81) и второй закон динамики, получаем уравнения тх — — F cos ф = — yt\ cos ф, (82) ту = — F sin ф = — ^р- sin ср. (83) Принимая же во внимание, что sin ц>=у/г, a cos q>—xlr, уравнения (82) и (83) можно переписать в виде kx •• ky л гз » У — гз » где постоянная k=yM Наконец,учитывая,что r = Vx2-{-y2t приходим к диффе- дифференциальным уравнениям " kx •' ky Х~~ У Далее, не умаляя общности рассуждений, можно считать, что выполняются следующие условия: х~а, у = 0, х = 0, y — v0 при £ = 0. (85) Таким образом, задача свелась к исследованию уравне- уравнений (84) при начальных условиях (85). Исходя из специфики уравнений (84), удобно воспользоваться полярными коор- координатами x—r cos ф, у=г sin ф. А тогда х—г cosy—(гзтф)ф, у = г sin ф + (г cos ф) ф, (86) х — г cos ф—2 (г sin ф) ф—(г sin ср) ср—(г cos ф) ф3 у = г sin ф + 2 {г cos ф) ф -\- (г cos ср) ср—(г sin ср) Отсюда 3 / у = (г—rep2) sin ф + Bгф + гф) cos ф Используя последние два равенства, дифференциальные уравнения (84) можно записать в виде (г_гф2)со5Ф_B/:ф + гфMтф = -^^, (87) (г—гф2) sin ф + Bг ф + гц) cos ср = — ~|^. (88) 62
Умножая обе части уравнения (87) на cos ф, уравнения (88) на sin ф и складывая результаты, находим, что ;__/-фа = — £/,.?. (89) Если же обе части уравнения (87) умножить на sin ф, а уравнения (88) — на cos ф и из первого полученного соот- соотношения вычесть второе, то получим уравнение 2гф + гф==0. (90) Что же касается начальных условий (85), то в полярных координатах они принимают вид г = а, ф = 0, г = 0, q> — vja при t — 0. (91) Итак, исследование уравнений (84) при начальных усло- условиях (85) свелось к исследованию уравнений (89) и (90) при начальных условиях (91). При этом можно заметить, что уравне- ние (90) можно переписать в виде •^-ИО-О. (92) Но уравнение (92) дает нам равен- о ство re«p = Clf (93) Рис- 30 где Сг — постоянная, которая имеет интересную геометри- геометрическую интерпретацию. Именно, предположим, что тело (рис. 30) движется из точки Р в точку Q по дуге PQ. Пусть S — площадь сектора, ограниченного отрезками OP, OQ и дугой PQ. Тогда из курса математического анализа из- известно, что или dS = -2 г2 dq>. Отсюда dS 1 dcp _ 1 2 dcp 1 2 Производная же -^- представляет собой так называемую векториальную скорость. А поскольку в силу равенства (93) величина г2ф является постоянной, то отсюда следует, 63
что и секториальная скорость также постоянна. Но это в свою очередь означает, что тело движется таким образом, что радиус-вектор «заметает» равные площади за равные промежутки времени. Этот закон площадей и есть один из трех законов Кеплера. Полностью он формулируется сле- следующим образом: каждая из планет движется по плоской кривой относительно Солнца таким образом, что радиус- вектор, связывающий Солнце и планету, «заметает» равные площади за равные промежутки времени. Чтобы вывести следующий закон Кеплера, касающийся типа траекторий движения планет, вернемся к уравнениям (89) и (90), рассматриваемым при условиях (91). Из началь- начальных условий (91), в частности, следует, что г—а, q>=i>0/a при /—0. А тогда из равенства (93) получаем Ci—avQ. Итак, г2ф = аи0, или ф = аио/г?. (95) Отсюда уравнение (89) принимает вид - _ a2vo 8 г2 Полагая здесь г—р, последнее уравнение можно переписать так: dp dp dr dp a2vl k r2 » dp a2t>o k „ или p -—• =—§ s-. Разделяя переменные в последнем дифференциальном уравнении и интегрируя, получаем со- соотношение Поскольку р—г=0 при г—а, то из последнего равенства имеем 2 , 2 2 а Таким образом, приходим к уравнению г* _ k a*pg vt k_ 2 ~~ r 2r* ~*~ 2 a * или, если рассматривать только положительный квадрат- квадратный корень, к дифференциальному уравнению 64
Разделив теперь уравнение (96) на уравнение (95), находим, что где l 2k Последнее уравнение интегрируется подстановкой г—Ни. В результате получаем где е = V~a + Pe/P = {аюЦЩ —1.11остоянная Сз определяется из условия г—а при <р=0. Нетрудно проверить, что С3—0. Таким образом, окончательно получаем, что г — l+ecosrp* Из аналитической геометрии известно, что (97) — это урав- уравнение в полярных координатах конического сечения, имею- имеющего эксцентриситет е. При этом возможны следующие случаи: 1) эллипс, если е<С1, т. е. v\<2kla\ 2) гипербола, если £>1, т. е. v^>>2kla; 3) парабола, если е—\% т. с. vl—2kla; 4) окружность, если е=0, т. е. vl=kla. Из астрономических наблюдений следует, что для всех планет Солнечной системы величина vl всегда меньше, чем 2k/a. Таким образом, приходим еще к одному закону Кеп- Кеплера: траектории планет являются эллипсами, в одном из фокусов которых находится Солнце. Отметим, что орбиты Луны и искусственных спутников Земли также являются эллипсами, но в большинстве слу- случаев эти эллипсы близки к окружности, т. е. эксцентриситет е мало отличается от нуля. Что же касается возвращающихся комет, таких как, на- например, комета Галлея, то они имеют орбиты в виде «вытя- «вытянутых» эллипсов, у которых эксцентриситет хотя и меньше единицы, но очень близок к ней. В частности, комета Галлея появляется в зоне видимости Земли приблизительно через каждые 76 лет. Следующее ее появление произошло в конце 1985 — начале 1986 года. Небесные тела, которые имеют параболические и гипер- гиперболические орбиты, могут наблюдаться лишь однажды, они никогда не возвращаются. 3 В. В. Амелькин 65
Выясним теперь физический смысл эксцентриситета е. Предварительно отметим, что компоненты к и у вектора скорости планеты V по осям х и у соответственно и вели- величина v вектора V удовлетворяют равенству которое с учетом формул (86) можно записать в виде Отсюда следует, что кинетическая энергия планеты мас- массой т определяется по формуле 1 mv2 = |m (rV 4- h. (98) Потенциальная же энергия системы есть взятая со знаком «минус» величина работы, требуемой для передвижения планеты в бесконечность (где потенциальная энергия равна нулю), и, следовательно, Г km , km — \ -я- йГ = — —. (99) Если теперь обозначить через Е полную энергию системы, которая в силу закона сохранения энергии является вели- величиной постоянной, то из формул (98) и (99) получаем, что -^■т[г2(р2-\-'г2] -==£. A00) Полагая ф=0, из соотношений (97) и A00) имеем о mr2a2i'e km i- r Исключая из последних двух равенств г, находим, что эксцентриситет = ]/ В результате уравнение орбиты планеты (97) принимает вид аЧ'о/k Г == + Е Ba4%/mk2)cos ц> Из последней формулы следует, что орбита будет эллипсом, гиперболой, параболой, окружностью, если соответственно 66
пт £<0, £>0, Е=0, Е=—mk2/Bazvl). Таким образом, харак- характер орбиты планеты полностью определяется полной энер- энергией Е. В частности, если бы планета, например Земля, могла получить извне толчок, позволивший бы увеличить полную энергию Е до положительной величины, то Земля перешла бы на гиперболическую орби- орбиту и покинула нашу Солнеч- Солнечную систему. Обратимся теперь к треть- третьему закону Кеплера, который касается периодов обращения планет вокруг Солнца. Исхо- Исходя из предыдущего закона Кеплера, мы, естественно, or- p раиичимся лишь случаем эл- эллиптической орбиты, уравне- уравнение которой в декартовых координатах записывается, как известно, в виде где (рис. 31) эксцентриситет с=С/|, а С2=£2—г]2, так что пли A01) Рассматривая теперь совместно равенства (97) и A01) и учитывая свойства эллипса, приходим к выводу, что k\f. ' т. е. A02) Обозначим через Т период обращения планеты, т. е. время, необходимое для одного полного оборота планеты по своей орбите. Тогда, поскольку площадь области, ограниченной эллипсом, равна я|т), на основании формул (94) и (95) при- приходим к выводу, что п\ц=сюцТ12. Наконец, принимая 67
во внимание равенство A02), находим, что А это и есть формализованная запись третьего закона Кеп- Кеплера: квадраты периодов обращения планет пропорциональ- пропорциональны кубам больших осей их орбит. Прогиб балок Рассмотрим горизонтально расположенную балку АВ (рис. 32) постоянного поперечного сечения, сделанную из однородного материала. Ось симметрии балки указана на рис. 32 пунктирной линией. Предположим, что под влиянием сил, которые действуют на балку в вертикальной плоскости, содержащей ось симметрии, балка прогибается (рис. 33). Рис. 32 Рис. 33 Действующие силы могут быть обусловлены несом балки, внешне приложенной нагрузкой или как той, так и другой силами вместе. I [онятно, что под действием сил ось симмет- симметрии будет искривляться. Обычно искривленную ось сим- симметрии называют упругой линией. Определение формы этой линии играет важную роль в теории упругости. д Р • Q Рис. 34 Рис. 35 Отметим, что существуют различные типы балок в зави- зависимости от способов их крепления или опоры. Например, на рис. 34 изображена балка, у которой конец Л жестко за- закреплен, а конец В свободен. Такая балка называется кон- консольной балкой. На рис. 35 показана балка, лежащая сво- свободно на опорах А и В. Еще один тип балок с опорами показан на рис. 36. Существуют] и различные способы приложения внешних нагрузок. Например, на рис. 34 пока- показана равномерно распределенная нагрузка. Конечно, иа- 68
грузка может быть и переменной вдоль всей длины балки или некоторой ее части (рис. 35). На рис. 36 указан случай сосредоточенной нагрузки. \в Рис. 36 Рис. 37 О X Рассмотрим горизонтальную балку О А (рис. 37). Пусть ее ось симметрии (показанная на рисунке пунктиром) лежит на оси х, где за положительное направление выбирается направление вправо от точки О, являющейся началом коор- координат. За положительное направление на оси у выберем направление вниз от точки О. Под действием внешних сил Fu F2, . . . (и веса балки, если он большой) ось симметрии искривляется в упругую линию, которая показана на рис. 38 пуьктиром. Смещение супругой линии от оси х называется прогибом балки в положении х. Таким образом, если изве- известно уравнение упругой ли- линии, то всегда можно ука- указать и прогиб балки. Ниже мы покажем, как это может быть сделано практически. Обозначим через М(х) из- изгибающий момент в вертикаль- вертикальном поперечном сечении балки с координатой х. Изгибаю- Изгибающий момент определяется как алгебраическая сумма момен- моментов тех сил, которые действуют с одной стороны балки в положении х. При подсчете моментов будем считать, что силы, которые действуют на балку снизу вверх, дают отри- отрицательные моменты, а силы, действующие сверху вниз, дают положительные моменты. В сопротивлении материалов доказывается, что изгибаю- изгибающий момент в положении х связан с радиусом кривизны упругой линии соотношением у" рис> EJ (ЮЗ) где Е — модуль упругости Юнга, который зависит от мате- материала, / — момент инерции поперечного сечения балки в 6Э
положении х относительно горизонтальной прямой, прохо- проходящей через центр тяжести этого поперечного сечения. Поо- изведеиие EJ обычно называют жесткостью' при изгибе■; ее величину в дальнейшем будем считать постоянной. Теперь, если предположить, что балка лишь слегка про- прогибается, что часто бывает на практике, то угловой коэффи- коэффициент у' упругой линии будет очень мал, и поэтому вместо уравнения A03) можно рассматривать приближенное урав- уравнение = M{x). A04) Чтобы показать, как на практике используется уравне- уравнение A04), рассмотрим следующую задачу. Горизонтальная однородная стальная балка длины /, свободно лежащая на Pl/2 Pi/2 со 1-х | Y J7 -^ч-i D т Q Рис. 39 двух опорах, прогибается под действием собственного веса, равного р кгс на единицу длины. Требуется найти уравнение упругой линии и максимальный прогиб балки. На рис. 39 упругая линия показана пунктиром. Посколь- Поскольку балка является двухопорной, то каждая из опор создает направленную вверх реакцию, равную половине веса балки (равную р112). Изгибающий момент М(х) есть алгебраиче- алгебраическая сумма моментов этих сил» действующих на балку с одной стороны от точки Q (рис. 39). Рассмотрим сначала действие сил слева от точки Q. На расстоянии х от точки Q сила pt/2 действует на балку снизу вверх и создает отри- отрицательный момент. Сила же рх, которая действует на балку сверху вниз на расстоянии х/2 от точки Q, создает положи- положительный момент. Таким образом, суммарный изгибающий момент в точке Q задается формулой рх plX ... /->г\ = ~2 g". A00) Если же рассмотреть действие сил справа от точки Q, то в этом случае на расстоянии (/—х)/2 от точки Q на балку действует сверху вниз сила рA—л;), которая создает поло- положительный момент. Отрицательный же момент создает сила 70
pl/2, которая действует на балку снизу вверх на расстоянии /—х от точки Q. Суммарный изгибающий момент подсчиты- вается в данном случае по формуле М(х)^рA-х)^-^-{1-х) = -^—^. A06) Как показывают формулы A05) и A06), изгибающие момен- моменты в обоих случаях оказываются равными. Теперь, зная, как находится изгибающий мохмент, легко выписать и ос- основное уравнение A04), которое в нашем случае принимает вид J?£-J*L. A07) Учитывая же, что на концах О и Л балка не прогибается, для нахождения у из уравнения A07) воспользуемся усло- условиями на концах балки: у=0 при х=0 и у=0 при х=1. А тогда интегрирование уравнения A07) с учетом последних условий дает x*-2lx* + l*x). A08) Уравнение A08) является уравнением упругой линии. Фор- Формула A08) используется на практике для определения мак- максимального прогиба. Так, в нашем конкретном случае, ос- основываясь на соображениях симметрии (это можно сделать и прямыми вычислениями), находим, что максимальный прогиб будет при х=112 и равен он 5/?/4/C84£7), где Е= =2Ы05кгс/см2, /=3-104 см4. Транспортировка леса При транспортировке спиленных деревьев на деревооб- деревообрабатывающие предприятия лесовозы определенную часть пути обязательно проезжают по лесным дорогам. Обычно ширина лесной дороги такова, что по ней может проехать только один лесовоз. Для того чтобы встречные машины могли разъехаться, на дороге устраивают несколько разъ- разъездных участков. Не касаясь вопроса об оптимальном режи- режиме движения, при котором нагруженные и пустые лесовозы должны встречаться именно на разъездных участках, выяс- выясним, насколько широким должен быть поворот дороги и по какой траектории водитель должен вести автомобиль на повороте, чтобы лесовоз мог перевозить, скажем, тридца- 71
тиметровые стволы. При этом предполагается, что лесовоз имеет большую маневроспособиость, позволяющую манев- маневрировать на довольно ограниченном участке дороги. С Q 'а X Рис. 40 Обычно лесовоз состоит из тягача и трейлера (прицепа), свободно сцепленных друг с другом. Тягач имеет один пе- передний (ведущий) мост и два задних моста, над которыми ус- устанавливается так называе- называемый коник, состоящий из по- поворотной балки со стойками и круглой платформы, свободно вращающийся "В горизонталь- горизонтальной плоскости относительно симметрично расположенной вертикальной оси. На коник одним из концов и укладыва- укладываются стволы деревьев. Трейлер имеет только два задних мос- моста, но также с присоединенным к ним коником, к которому крепятся свободные концы стволов деревьев. Шасси трей- трейлера состоит из двух метал'ли- ческих цилиндров, один из которых может входить внутрь другого. Шасси соединяет коник с осью, скрепляющей трей- трейлер с тягачом. Таким образом, длина шасси может изме- изменяться во время движения, а это позволяет тягачу и трей- трейлеру двигаться независимо. Схематически лесовоз изобра- изображен на рис. 40. При этом точки А и В соответствуют осям коников, расположенным друг от друга на расстоянии h. Через XY обозначен ствол дерева, для которого AX-=%h. Точка С соответствует небольшой оси, которая скрепляет тягач с трейлером, причем AC—all. Обычно а=0, 3; что же касается простейшего случая буксировки, то а=0. Далее, FF — передняя ось тягача, РР и QQ — его задние оси, RR и SS — оси трейлера. Все оси имеют одну и ту же длину 2L. Будем считать, что ширина лесовоза также равна 2L. Что же касается ширины груза в хвостовой его части, то ее обозначим через 2W. В дальнейшем нам придется поль- 72
зоваться таким понятием, как размах груза лесовоза. Обыч- Обычно под размахом груза понимают максимальное отклонение хвостовой его части (для простоты точки X) от траектории движения лесовоза. Будем считать, что дорога имеет ширину 2р/г и что обычно поворот на дороге имеет вид дуги окруж- окружности радиуса hi а с центром в точке О (рис. 41). Для прос- простоты предположим, что лесовоз входит в поворот таким об- образом, что тягач и трейлер располагаются на одной прямой линии; водитель при этом управляет лесовозом так, что точ- точка Л, соответствующая оси переднего коника, находится точно над средней линией дороги. Точка А на рис. 41 опре- определяется углом %, который тягач АС составляет с началь- начальным направлением. При этом удобно задать систему коор- координат Оху так, чтобы ось абсцисс была параллельна выб- выбранному начальному направлению, а ось ординат была ему перпендикулярна. В общей ситуации перевозимый лесово- лесовозом груз будет составлять с начальным направлением неко- некоторый угол 0. Что же касается угла ВАС на рис. 41, то, обозначая его через и, находим, что и=%—6. Обычно угол и называют углом запаздывания лесовоза. Требуемая полу- полуширина дороги Л, которая определяет размах груза лесово- лесовоза на повороте и которая называется полушириной дороги на внешней стороне поворота дороги, определяется алгеб- алгебраической суммой ОХ—ОА 4- W. Полуширина дороги на внутренней стороне поворота определяется уже алгебраи- алгебраической суммой ОА -f- L—ОР, где ОР — перпендикуляр из точки О на" Л В. Условимся, что при движении лесовоза его колеса либо вовсе не испытывают бокового скольжения, либо если оно есть, то оно мало. Это требование, в частности, означает, что средняя линия тягача АС является касательной к дуге ок- окружности в точке /4, так что О А есть перпендикуляр к АС, а угол % задается движением точки А по дуге окружности. Заметим далее, что при строительстве дорог кривизна пово- поворота дороги определяется углом №, который соответствует приблизительно длине дуги поворота в 30 м. В наших обоз- обозначениях П h * V ' где h измеряется в метрах. Таким образом, для й=9 м и а=0, 1 получаем, что N&W. Если же h—\2 м, а а=1, 0, то W«142°. Сточки зрения практики необходимо рассмат- рассматривать лишь такие а, которые удовлетворяют условию 73
О < a < 1, Причем для увеличения маневроспособности лесовоза а должно быть по возможности большим. Длина ствола транспортируемого дерева Ш будет боль- больше, чем h, по опять-таки, исходя из практических сообра- соображений, она не должна превосходить ЗЛ. Таким образом, значения X изменяются в пределах 1 <Я<3. Относительно константы а предполагается, что 0^д<0,5. Наконец, отметим, что величина h выбирается в каждом конкрет- конкретном случае по-разному, но обычно ее эначение колеблется в пределах 9—12 м. Так как колеса тягача не испытывают бокового сколь- скольжения, то координаты точки А на рис. 41 опрел°ляются так: h . h x = -smx, # = -cosx. Координаты точки В определяются так: £cos + uie A10) При этом, так как колеса трейлера также не скользят, то точка В движется в направлении ВС и где \\>—это угол, который ВС составляет с начальным направлением. Далее, воспользовавшись равенством у, = — и + 6, из треугольника ABC получаем цепочку равенств sin(%—тр) _ sin @—-ф) _ sin» П19 h ~~ ah ~~ bh * lllZ> где 0<Ь<1, а я)), 6, и является функциями %. Если вос- воспользоваться. равенством A11), то на основании формул (НО) получаем, что 1 V smyv + h -5— cos C cos ib + a dx J sine) sin sine) Упрощая последнее равенство, приходим к соотношению sin(x—^ = o^-cos@—ф). A13) Если исключить из уравнения (ИЗ) переменную-ty при за- заданном а, воспользовавшись равенствами A12) и соотноше- соотношением % — и + в, то, так как 0 = 0 при % = 0, придем *) к диф- *) ТауТег А. В. The Sweep of а Logging Truck//Math. Spect- Spectrum.—1974—1975.—V. 7, Nil.— P. 19—26. 74
ференциальному уравнению du — 1 sin и A14) a(l— a cos u) с начальным условием u@)=0, где роль искомой функции играет угол запаздывания. Подстановкой tg — = v дифференциальное уравнение A14) может быть проинтегрировано в замкнутой форме. Однако получаемая при этом зависимость между перемен- переменными и и х оказывается очень сложной, что, естественно, затрудняет исследование полученного соотношения. Вме- Вместе с тем уравнение A1.4) может быть легко проинтегриро- проинтегрировано и исследовано численно. Воспользуемся, например, численным методом Рунге — Кутта второго порядка, суть которого состоит в следующем. Для интегральной кривой и=<р(%) дифференциального уравнения du/d%—f(%t и), проходящей через точку (%0, иа)г п равноотстоящих друг от друга точках %о, %ь %2, ... G/+1—Хг^&У^О) берутся значения и0, щ, и2, . . такие, что u^(pi(%), причем для последовательного определения иг, и2, . . . пользуются формулами где Для численного нахождения решения уравнения A14) с начальным условием и@)=0 при помощи микрокалькуля- микрокалькулятора «Электроника БЗ-34» составим следующую программу: В/О; ГЮ; С/П; ПП; 11; +; ГТО; ИПВ; ИПА; -|-; ПВ; ПП; 22; ИПА; X; 2; -=-; f; ИПС; +; ПС; В/О; НПО; F sin; ИП9; НПО; F cos; x; ИП7; —; -^; 1; +; В/О. Эта программа использует регистры памяти 7—9 и А, В, С, Д, назначение которых видно из следующей таблицы: Регистры памяти Содержание 7 а 8 а 9 аа А Дх В С щ D Отметим, что для составления таблицы значений реше- решения, для конкретных значений а и а следует вычислить пре- предельное значение С згого решения. Число С находится из 75
условия аA—a cos C)=sin С, которое приводит к формуле При этом, если а=1, то С=(я/2)—2 arctg а; если же ^, то С~аA—а). Предельное значение С аппроксимируется экспонентом, а именно, С—ы«е-*«, a (cos С—а) _, ^ где 7 = —gin2C—-. Для малых а тогда у довольно большое и оно приблизительно таково: ТЛ^ _______ ~ аA-а)' Если же а=1, то у 1 В качестве шага Ах изменения аргумента % следует вы- выбирать число, не превосходящее С. Это требование особенно а-/ с-0, Рис. 42 важно при малых значениях а, а также в окрестности зна- значения %—0. Непосредственные вычисления проводятся в следующем порядке: а; П7; а\ П8; х; П9; Ах; ПА; 0; ПВ; ПС; В/О; С/П; щ\ С/П; щ; С/П; и_2; ...; С/П; щ\ ..., где черта указывает на показание табло после очередного нажатия клавиши С/П. Другими словами, черта указывает па очередное значение иг~и(%^)у где %t=i Д%. 76
Результаты численного счета при а=0,3 графически изображены на рис. 42, и они показывают, как с увеличе- увеличением угла х изменяется угол запаздывания лесовоза Для большеп наглядности масштаб на осях и и % выбран разным. Н 18*0,3 Рис, 43 Определим теперь размах груза лесовоза, воспользовав- воспользовавшись полушириной дороги на внешней стороне ее поворота, которая, как уже указывалось выше, задается алгебраичес- алгебраической суммой ОХ—ОА-i-W (рис. 41). Прежде всего находим, что Отсюда понятно, что размах груза убывает с возрастанием X, поскольку возрастает угол запаздывания и. Таким обра- образом, максимальная полуширина дороги p/i определяется равенством На рис. 43 сплошные линии соответствуют графикам, a W связывающим переменные р—у-и а при различных зна- значениях Я. Пунктирные же линии показывают, каким долж- должно быть значение р —-^ для необходимого «зазора» на внутренней стороне поворота дороги. При этом для любого 77
места нахождения лесовоза на повороте должно выполнять- выполняться соотношение A15) Далее, так как значение С убывает с уменьшением а, то в простейшем случае буксировки, когда а—0, угол запазды- запаздывания оказывается наибольшим. В этом случае из соотно- соотношения A15) находим, что Остановимся теперь на тех выводах, которые следуют из проведенных выше рассуждений. Прежде всего приведем типичный пример, иллюстрирующий полученные резуль- результаты. Если лесовоз с расстоянием между кониками в 12 м про- проходит поворот по дуге окружности радиуса 60 м, тоа=0,2, и тогда в соответствии с формулой A09) поворот составляет приблизительно 28°. Если при этом ширина лесовоза равна 2,4 м, а ширина связки груза в хвостовой части равна 1,2 м, то в случае длины стволов 24 м (от переднего коника) при- приходим к выводу, что А,=2, W7/i=0,05, a L//i=0,l. Таким об- образом, из рис. 43 следует, что для любого значения а для внешней стороны поворота дороги {5=0,45, а для внутренней стороны поворота |3=О,2. Теоретически необходимая полу- полуширина дороги на ее внешней стороне поворота тогда равна 5,4 м, а на внутренней стороне поворота — 2,4 м. Если лесо- лесовоз перевозит спиленные деревья длиной 14,4 м от перед- переднего коника и шириной связки в хвостовой части, равной 1,8 м, то в этом случае Я,= 1,2. Значение Р в этом случае оказывается одним и тем же как для внешней стороны пово- поворота, так и для его внутренней стороны, и оно равно 0,22. Поэтому, как нетрудно видеть, необходимая полуширина дороги на ее внешней стороне поворота равна 2,64 м, а на внутренней стороне поворота, как и в предыдущем случае, она равна 2,4 м. Из проведенных рассуждений следует, что чем длиннее перевозимый груз, тем шире должна быть до- дорога на повороте. В частности, если сравнивать рассмотрен- рассмотренные только что два случая, то увеличение длины перевози- перевозимых деревьев на 9,6 м требует расширения дороги на пово- повороте на 2,76 м для того, чтобы водитель мог вести машину по кривой, длина которой на повороте приблизительно равна 78
длине средней линии дороги. Практика показывает, что водитель с малым опытом не может вести машину по такой кривой и для него полная ширина дороги на повороте долж- должна быть по крайней мере равна 10,8 м (если он перевозит груз длиной 24 м) при ширине лесовоза в 2,4 м. Развитая выше теория показывает, что размах груза наибольший, когда лесовоз въезжает в поворот, ибо в этом случае угол запаздывания увеличивается. Этот вывод имеет силу и в топ ситуации, когда лесовоз въезжает с одного участка зигзагообразного поворота на другой его участок через точку перегиба. Результаты, которые графически изображены на рис. 43, соответствуют тому случаю, когда до въезда в поворот тягач и трейлер находятся на одной прямой. Если же имеется начальный угол запаздывания Со, обусловленный зигзагообразностью поворота, то в диффе- дифференциальном уравнении A14) начальное условие надо вы- выбирать в виде и@)~—Со. В этом случае необходимая шири- ширина дороги определяется формулой sinCe—-+~. ° а ' h Обратим теперь внимание на то, что в случае простой буксировки, т. е. когда а—О, невозможно преодолеть пово- поворот при сО1. Однако для сравнительно больших значений а значения а>1 уже возможны и они должны удовлетво- удовлетворять равенству аA—a cos C)=sin С. Таким образом, а имеет максимальное значение 1/Kl—а2- и для а=0,5 как экстремального практического значения а=1,25. В заключение отметим, что для значений сС>0,5 значи- значительная экономия ширины дороги достигается за счет уве- увеличения значения а (рис. 43). При этом, если груз такой, что Я не намного больше, чем Л,= 1, то значение а выбирается таким, чтобы необходимая полуширина дороги па внутрен- внутренней стороне любого поворота всегда была меньше, чем на его внешней стороне.
КАЧЕСТВЕННЫЕ МЕТОДЫ ИССЛЕДОВАНИЯ ДИФФЕРЕНЦИАЛЬНЫХ МОДЕЛЕЙ При решении приведенных выше задач мы строили диф- дифференциальные модели и ответы на поставленные вопросы получали после интегрирования дифференциальных урав- уравнений. Однако, как уже указывалось в предисловии, по- подавляющее большинство дифференциальных уравнений не может быть проинтегрировано в замкнутой форме. Поэтому для исследования дифференциальных моделей реальных явлений и процессов приходится изыскивать методы, кото- которые позволяли бы получать необходимую информацию, исходя из свойств самого дифференциального уравнения. Ниже на конкретных, примерах показывается, как при решении практических задач используются простейшие приемы и методы качественной теории обыкновенных диффе- дифференциальных уравнений. Кривые с постоянным направлением магнитной стрелки Покажем, как при качественном интегрировании, т. е. при выяснении общей природы решений обыкновенных диф- дифференциальных уравнений, можно воспользоваться одним их общим свойством, аналогом которого является, напри- например, свойство магнитного поля па поверхности Земли, за- заключающееся в том, что на земной поверхности можно ука- указать кривые, вдоль которых направление магнитной стрелки постоянно. Итак, рассмотрим обыкновенное дифференциальное урав- уравнение первого порядка где функция f предполагаается однозначной и непрерывной по совокупности переменных х, у в некоторой области D 80
плоскости (л;, у). Каждой точке М(х, у) из области D опре- определения функции f дифференциальное уравнение A16) ста- ставит в соответствие значение j-—угловой коэффициент К касательной к интегральной кривой в точке М(х, у). Имея это в виду, говорят, что в каждой точке М (х, у) обла- области D дифференциальное уравнение A16) задает направле- направление или линейный элемент. Совокупность же всех линейных элементов в D называют полем направлений или полем ли- линейных элементов. Графически линейный элемент изобража- изображается отрезком, для которого М(х, у) является внутренней точкой и который образует с положительным направлением оси х угол 0 такой, что /C=tg Q=f(x, у). Отсюда следует, что геометрически дифференциальное уравнение A16) вы- выражает тот факт, что направление касательной в каждой точке интегральной кривой совпадает с направлением поля в этой точке. При построении поля направлений полезно использо- использовать изоклины (от греческого isos — равный, одинаковый и kllno — наклоняю), т. е. множества точек в плоскости (х, у), в которых направление поля, задаваемое дифференциаль- дифференциальным уравнением A16), одно и то же. Если говорить об изоклинах магнитного поля на по- поверхности Земли, то ими являются кривые, в каждой точке которых направление магнитной стрелки одно и то же. Что же касается рассматриваемого дифференциального уравне- уравнения A16), то его изоклины задаются уравнением f(x, y) = v, где v — переменный действительный параметр. Знание изоклин дает возможность приближенно выяс- выяснить поведение интегральных кривых заданного дифферен- дифференциального уравнения. Так, например, рассмотрим диффе- дифференциальное уравнение не итерируемое в замкнутой форме. Из вида дифферен- дифференциального уравнения следует^ что семейство его изоклин задастся уравнением x* + tf = v, v>0, т. е, изоклинами являются концентрические окружности радиуса У v с центром в начале координат, расположенные 81
в плоскости (лг, у). В данном случае в каждой точке изоклин угловой коэффициент касательной к интегральной кривой, проходящей через эгу точку,, равен квадрату длины радиуса соответствующей окружности. Уже этой информации у к достаточно, чтобы полу- получить представление о по- поведении интегральных кри- кривых рассматриваемого диф- дифференциального уравнения (рис. 44). Мы быстро пришли к конечной цели, так как приведенный пример срав- сравнительно прост, но и в слу- случае более сложных уравне- уравнений знание изоклин может оказаться полезным при решении поставленной за- задачи. Рассмотрим геометрический метод интегрирования диф- дифференциальных уравнений вида A16), основанный на ис- использовании геометрических свойств кривых, заданных уравнениями Рис. 44 Уравнение (/п), т. е. уравнение изоклин нуля, задает кри- кривые, в точках которых -г, —О- Это означает, что точки указанных кривых могут оказаться для интегральных кри- кривых исходного дифференциального уравнения точками мак- максимума или минимума. Именно этим обстоятельством и ру- руководствуются, выделяя из всего множества изоклин изо- изоклину нуля. Для большей точности построения интегральных кривых находят также множество их точек перегиба (если таковые существуют). Точки перегиба, как известно, следует искать среди тех точек, где у" — 0. Используя уравнение A16), получаем, что У — у) дх 1 У) У) Тогда кривые, заданные уравнением (L), и будут возмож- 82
ными линиями точек перегиба *). Заметим, в частности, что точка перегиба интегральной кривой — это точка касания интегральной кривой с изоклиной. Кривые экстремумов (точек максимума или минимума) и точек перегиба интегральных кривых разбивают область определения функции / на такие подобласти Su 52, . . . ..., Sm, в которых первая и вторая производные решения дифференциального уравнения имеют определенные знаки. В каждом конкретном случае эти области и необходимо уста- установить, что позволит дать схематическую картину поведе- поведения интегральных кривых. Рис. 45 В качестве примера рассмотрим дифференциальное урав- уравнение у' — х-\-у. Уравнение кривой (/0) в данном случае имеет вид к -\-у— =0, или */=—х. Непосредственной проверкой убеждаемся, что кривая G0) не является интегральной. Обращаясь же к кривой (L), уравнение которой у + х + 1=0, находим, что она уже оказывается интегральной и, таким образом, не есть линия точек перегиба. Прямые (/0) и (L) разбивают (рис. 45) плоскость (х, у) на три подобласти: Si(y'>0, //*>0)—справа от прямой (/о); 52(//'<0, #">0) — между прямыми (/„) и {L)\ S3(t/<0, *) Предполагаем при этом, что интегральные кривые, заполняющие некоторую область, обладают тем свойством, что ч^рез каждую точку этой области проходит только одна интегральная кривая. 83
*/"<0) — слева от прямой (L). На прямой (/0) расположены точки минимумов интегральных кривых. Справа от прямой (/0) интегральные кривые поднимаются вверх, слева опус- опускаются вниз (на рис. 45 слева направо). Точек перегиба пет. Справа от прямой (L) кривые выпуклы вниз, слева — вверх. Поведение интегральных кривых в целом показано на рис. 45. Заметим, что в данном случае интегральная прямая (L) является своего рода «разделительной» кривой — она отде- отделяет одно семейство интегральных кривых от другого. Такие кривые обычно называют сепаратрисами (от латинского se- separator — отделитель). Зачем инженеру знать теоремы существования \ и единственности? Выше, говоря об изоклинах и линиях точек перегиба- мы молчаливо предполагали, что соответствующее диффе- дифференциальное уравнение имеет решение. Вопрос же о том, когда решение существует, когда оно единственно, решается ■так называемыми теоремами существования и единственно- единственности. Эти теоремы очень важны как для самой теории, так и для практики. Теоремы существования и единственности имеют прин- принципиальное значение, гарантируя законность применения качественных методов теории дифференциальных уравнений для решения задач естествознания и техники. Они являются обоснованием для создания новых методов и теорий. Часто доказательства самих теорем существования и единственно- единственности являются конструктивными, т. е. методы доказательст- доказательства дают и методы приближенного отыскания решений с лю- любой степенью точности. Таким образом, теоремы существо- существования и единственности лежат в основе не только упоми- упоминавшейся выше качественной теории дифференциальных уравнений, но и в основе методов численного интегриро- интегрирования. К настоящему времени разработаны многочисленные ме- методы численного решения дифференциальных уравнений. Хотя эти методы обладают тем недостатком, что всегда дают лишь какое-то конкретное решение, что сужает возможности их использования, они тем не менее широко используются на практике. Следует, однако, отметить, что численному интегрированию дифференциального уравнения обязатель- обязательно должно предшествовать обращение к теоремам существо- существования и единственности. И это необходимо делать для того, 84
чтобы избежать недоразумений или вообще неправильных выводов. Для иллюстрации сказанного рассмотрим два простых примера *), но предварительно приведем формулировки одного из вариантов теорем существования и единственно- единственности. Теорема существования. Если в уравнении A16) функция f определена и непрерывна в некоторой огра- ограниченной области D плоскости (х, у), то для любой точки (*o. Уо)£В существует решение у (х) начальной задачи **) % = fix> У), У{хо) = Уо, (П7) определенное на некотором интервале, содержащем точку х0. Теорема существования и единст- единственности. Если в у равнении'A 16) функция f определена и непрерывна в некоторой ограниченной области D плоскости (х, у), причем она удовлетворяет в области D условию Лип- Липшица по переменной у, т. е. \f{x, У*)—fix, yi)\<L\y2—yi\, где L — положительная постоянная, то для любой точки ix0, У о) € D существует единственное решение у (х) начальной задачи A17), определенное на некотором интервале, содержа- содержащем точку х0. Теорема о продолжении. При выполнении условий теоремы существования или теоремы существования и единственности всякое решение уравнения A16) с началь- начальными данными (xOt y0) £D может быть продолжено до точ- точки, сколь угодно близкой к границе области D. При этом в первом случае продолжение, вообще говоря, будет не обяза- обязательно единственным, во втором же случае оно единственно. Рассмотрим следующую задачу. Требуется, используя численный метод интегрирования Эйлера с итерационной схемой yl+1=yt +hf(xit yt) и шагом Л=0,1, решить на- начальную задачу y' = -xly, у (-1) = 0,21 A18) на отрезке [—1,3]. *) Roberts С. Е., Jr. Why teach existence and uniqueness theo- theorems in the first course in ordinary differential equations?// Int. J. Math. Educ. Sci. Technol.—1976.—V. 7, № 1.—P. 41—44. **) Если требуется найти решение дифференциального уравне- уравнения, удовлетворяющее некоторому начальному условию (в данном слу- случае условию у(хо)—уо), то такая задача называется начальной за- задачей, 85
Заметим зцесь? что к исследованию уравнения #'=—x/'j приводит, например, рассматриваемая на с. 96 конеео- вативная система, состоящая из тела, которое совершает горизонтальные движения в вакууме под действием линей- линейных пружин. Решая начальную задачу A18) с помощью микрокаль- микрокалькулятора «Электроника БЗ-34», составляем сначала про- программу, которую запишем так: ПП; 15; ИГГО; X; ИПВ; +; ПВ; ИПА; ИГО; +; ПА; XY; С/П; БП; 00; ИПА; ИПВ; ~; /—/; В/О. Инструкции к программе такова: х0; ПА; у0', ПВ; h\ ПД; В/О; С/П; уг\ С/П; уг\ . . . Результаты численного интегрирования приводят к сле- следующей таблице: X У(х) ~1 0,21 -0,9 0,686 —0,8 0,817 —0,7 0,915 —0.G 0.992 —0,5 1,052 —0,4 1,10.0 —0.3 1,136 —0,2 1,163 У(х) —0,1 1.1S0 0,0 1,188 0,1 1,188 0,2 1,180 0,3 1,163 0,4 1,137 0,5 1,102 0,6 1,056 0,7 1,000 0,8 0,930 X У(х) 0,9 0,844 1,0 0,737 1,1 0,601 1,2 0,418 1,3 0,131 1,4 —0,859 1,5 —0,696 1,6 —0,480 X Vi*) 1,7 —0,146 1,8 1,014 1,9 0,837 2,0 0,609 2,1 0,281 2,2 —0,465 2.3 0,007 2,4 —31,625 Графическая интерпретация полученных результатов показана на рис. 46. Обратимся теперь к теореме существования. Для иссле- исследуемой начальной задачи A18) функция f> определяемая равенством f(x, у)=—х/у, определена и непрерывна во
всей плоскости (jc, у), за исключением точек оси абсцисс. Таким образом, в соответствии с теоремой существования существует решение у(х) начальной задачи A18), опреде- определенное на некотором интервале, содержащем точку xv=—1, и это решение по теореме о продолжении может быть про- продолжено до значения у(х), близкого к значению у{х)~0. -4-4—I I I ! 1 1 I I- ■/ -0.6 -дг -0,2 -0,6- -1 i i i 2,6 Я Рис. 46 В результате численного интегрирования получаем решение начальной задачи A18) на некотором интервале (а, Ь), где я<1, 1,3<СЬ<С1,4. Однако, учитывая конкретный вид диф- дифференциального уравнения, .можно установить истинный промежуток существования решения начальной задачи A18). Действительно, так как в исходном уравнении пере- переменные разделяются, то в,-21 Интегрируя, получаем, что t/= К 1,0441—хг. Отсюда ре- решение начальной задачи A18) существует только для М<1Л,0441 «1,0218. Итак, обращение к теореме существования (и к теореме о продолжении) позволило «отсечь» отрезок, на котором решение исходной начальной задачи заведомо не существует. Одно же только численное интегрирование приводит к оши- 67
бочному результату. Дело здесь в том, что при приближе- приближении решения у=у(х) к оси х угол наклона кривой прибли- приближается к 180°. Поэтому пока аргумент х изменяется на вели- величину 0,1, значение у успевает «перескочить» ось х, и мы попадаем на интегральную кривую, отличную от исходной. А это происходит потому, что метод Эйлера учитывает угол наклона только в текущей точке. Еще более поучительным является следующий пример. Требуется, используя сначала метод Эйлера, а затем улуч- улучшенный метод Эйлера, с шагом /г=0,1 и итерационной схе- схемой y^i^yi+/г/(^+1/2, yi+l/2), ще yi+l/2=yi+hf(Xi, yt)l2, решить начальную задачу на промежутке [—1, 1]. Решая начальную задачу A19) с помощью микрокальку- микрокалькулятора «Электроника БЗ-34», составляем сначала программу численного интегрирования методом Эйлера. Программа такова: ПИ; 15; ШТО; X; ИПВ; +; ПВ; ИПА; ИПБ; +; ПА; <-» X Y* С/П; БП; 00; ИПВ; Fx=Q; 19; В/О; Fxa; ИПЗ; XY; Fxv; ИПВ; X; ИПЭ; X; 3; X; В/О. Инструкция же к программе следующая: л:0; ПА; у0; ПВ; /г; ПО; 3; /—/; F\lx\ ПЗ; В/О; С/П; ух\ С/П; у2; . . . Результаты численного счета сведем в таблицу X ^ 1 — 1,000 —0,9 —0,700 -0,8 —0,460 -0,7 —0,275 —0,6 —0,138 —0,5 —0,045 X У1(х) —0,4 —0,008 —0,3 —0,016 —0,2 0,007 -0,1 —0,005 0,0 0,000 0,1 0,000 0,2 0,003 X У\ (х) 0,3 0,011 0,4 0,031 0,5 0,068 0,6 0,129 0,7 0,220 0,8 0,347 0,9 0,516 1,0 0,732
Рис. 47 Графически имеем картину, представленную на рис. 47 Что же касается улучшенного метода Эйлера, то здесь программа такая: ПС; С/П; ПП; 29; ИГО; X; ИПВ; +; ПС; ИПА; ИГО; +; ПА; ПП; 29; ИПЭ; X; 2; X; ИПВ; +; ПВ; ИПА; ИГО; +; ПА; XY; БП; 00; ИПС; Fx=0; 33; В/О; Fx2; ИГО; XY; Fxy; ЙПС; X; ИПА; X; 3; X; В/О. Инструкция к программе следующая: 3; /—/; F\lx; ПЗ; Л/2; ПЭ; х0; ПА; у0; ПВ; В/О; С/П; С/П; у,; С/П; г/2; . . . Полученные числовые значения заносим в таблицу: X X 1 — 1,000 —0,4 —0,068 —0,9 —0,730 —0,3 —0,031 —0,8 —0,514 —0,2 —0,012 -0,7 —0,346 —0,1 —0,004 —0,6 —0,219 0,0 —0,002 —0,5 —0,129 0,1 —0,004
X Г/2 (А) 0,2 —0,013 0,3 —0,033 0,4 —0,071 0,5 —0J33 0,6 —0,225 X 0,8 —0,522 0,9 —0,739 1,0 —1,010 0,7 —0,352 Графически получаем картину (рис. 48), отличную от картины, изображенной на рис. 47. 0,6- 0,1- -1 -0,6 -0,2 Чтобы лучше разобраться в причине столь разительного расхождения в результатах, проинтегрируем исходную на- начальную задачу. Разделяя переменные, имеем или, окончательно, у=±хэ. Отсюда уже можно заметить, что решение по методу Эйлера приближает функцию Ух{х) = =ха, а по улучшенному методу Эйлера — функцию У 2 М л:3, если -л:3, если х > 0. 90
при этом как i/i, так и у2 являются решениями начальной задачи A19), а значит, для рассматриваемой на промежутке [—1, 1] начальной задачи имеет место неединственность. Обращаясь теперь к теореме существования и единст- единственности, отметим прежде всего, что так как функция f» заданная равенством f(x, у) = 3х\/ у, непрерывна во всей плоскости (л:, у), то из теоремы существования следует, что существует решение начальной задачи A19), определенное на некотором промежутке, содержащем точку л:0——1, и это решение по теореме о продолжении может быть продолжено на любой промежуток. Далее, поскольку *' — ху~2'-\ то функция / (х, у) = 3х\/ у удовлетворяет условию Липшица по переменной у в любой области, не содержащей точки оси л:. Если же область содержит точки оси х, то нетрудно показать, что в ней указанная выше функция условию Лип- Липшица уже не удовлетворяет. Поэтому из теоремы существо- существования и единственности (и теоремы о продолжении) следует, что в данном случае решение начальной задачи может быть продолжено единственным образом по крайней мере до оси х. Но поскольку прямая у—0 является особой интеграль- интегральной прямой для дифференциального уравнения у' — Зхр'у, то мы уже знаем, что как только у станет равным нулю, реше- решение начальной задачи A19) не может быть единственным образом продолжено за точку 0@, 0). Итак, обращение в данном случае к теореме существо- существования и единственности (и теореме о продолжении) позволи- позволило разобраться в результатах численного интегрирования. Именно, если речь идет о единственном на промежутке [—1, 11 решении начальной задачи A19), то оно существует и определено лишь на отрезке [—1, 0]. В общем же случае таких решений несколько. Динамическая интерпретация дифференциальных, уравнений второго порядка Рассмотрим нелинейное дифференциальное уравнение частным случаем которого является дифференциальное урав- уравнение второго порядка, полученное на с. 46 при рас- рассмотрении маятниковых часов, и представим простую ди- динамическую систему, состоящую из частицы единичной 91
массы, которая движется по оси х (рис. 49) и на которую действует сила f(x,-£). Тогда дифференциальное урав- \ at I пение A20) — это уравнение движения частицы. Значениям х и tj- , в любой момент времени характеризующим со- / йх\ стояние системы, соответствует точка на плоскости ( х, -тг I (рис. 50), называемой плоскостью состояний или фазовой dt ' dt) О sb Q Рис. 49 Рис. 50 плоскостью. Фазовая плоскость изображает совокупность всех возможных состояний рассматриваемой динамической системы. Каждому новому состоянию системы соответству- соответствуют различные точки фазовой плоскости. Таким образом, изменению состояний системы можно поставить в соответ- соответствие движение некоторой точки на фазовой плоскости. Такую точку называют изображающей точкой. Траектория изображающей точки называется фазовой траекторией, а скорость этой точки — фазовой скоростью. Если ввести переменную y=jr> то уравнение A20) можно свести к системе двух дифференциальных урав- уравнений ~ = #, 37 = /(*. У)- A21) При этом, если t рассматривать как параметр, то реше- решением системы A21) является пара функций х{1) и y(t), опре- определяющих в фазовой плоскости (х, у) некоторую кривую (фазовую траекторию). Можно показать, что система A21), как и более общая система вида где функции X и У непрерывны вместе со своими частными производными в некоторой области D, обладает тем свойст- 92
вом, что если л:(/), y(t)— решение дифференциальной си- системы, то A23) где С — произвольная вещественная постоянная, также есть решение рассматриваемой дифференциальной системы. Решениям же A23) при всевозможных значениях С соответ- соответствует на фазовой плоскости (х, у) одна и та же фазовая траектория. Далее, если две фазовые траектории имеют общую точку, то они совпадают. При этом возрастанию или убыванию параметра t соответствует определенное движение изображающей точки по траектории. Иначе говоря, фазо- фазовая траектория является направленной или ориентирован- ориентированной кривой. В том случае, когда мы будем обращать вни- внимание на этот факт, на соответствующих рисунках направ- направление движения изображающих точек по траектории будет указываться стрелками. Отметим далее, что системы вида A22) относятся к так называемым автономным или стационарным дифферен- дифференциальным системам, т. е. системам обыкновенных диффе- дифференциальных уравнений, правые части которых не зависят явно от времени /. Если же хотя бы в одном из уравнений, входящих в систему, правая часть зависит явно от времени t, то такая система уже называется неавтономной или неста- нестационарной. В связи с приведенной классификацией дифференци- дифференциальных уравнений обратим внимание на следующий факт. Если решение x(f) уравнения A20) является непостоянным периодическим решением, то на фазовой плоскости (х, у) ему соответствует простая замкнутая кривая, т. е. замкну- замкнутая кривая без самопересечений. Справедливо и обратное утверждение. Далее, если дифференциальные системы вида A22) за- заданы во всей плоскости (х, //), то, вообще говоря, фазовые траектории полностью покроют фазовую плоскость, не пересекаясь друг с другом. Если при этом окажется, что в некоторой точке М0(х0, у0) выполняются равенства то траектория вырождается в точку. Такие точки называют- называются особыми точками. В дальнейшем будут рассматриваться в основном лишь изолированные особые точки. Особая точка М0(х0, у о) называется изолированной, если можно указать ее окрестность, в которой нет других особых точек, кроме самой точки М0(х0, у0), 93
С точки зрения физической интерпретации для уравне- уравнения A20) особой точкой будет точка М0(х0, 0), в которой */—0 и f (х0, 0)=0. Таким образом, в данном случае особая точка соответствует такому состоянию частицы единичной массы, в котором как скорость ~, так и ускорение — =* d2x = -Т7% одновременно равны нулю. Это означает, что час- частица находится в состоянии покоя или равновесия. В связи с последним обстоятельством особые точки называют еще точками покоя или точками равновесия. Рис. 51 Состояния равновесия физической системы характери- характеризуют особенные ее состояния, и поэтому изучение типов особых точек занимает видное место в теории дифференци- дифференциальных уравнений. Впервые вопрос о классификации особых точек диффе- дифференциальных систем вида A22) был подробно рассмотрен известным русским ученым Н. Е. Жуковским в его магис- магистерской диссертации «Кинематика жидкого тела» A876). Этот вопрос возник в связи с изучением теории скоростей и ускорений жидкости. Названия же различных типов осо- особых точек впервые были предложены французским матема- математиком Л. Пуанкаре. Попытаемся теперь ответить на вопрос о том, какой фи- физический смысл можно придать фазовым траекториям и особым точкам дифференциальных систем вида A22). Для наглядности Еведем двумерное векторное поле (рис. 51), определяемое функцией V{x, y) = X{x, y)l + Y{x, y)j, где i, J— орты осей декартовых координат х и у соответст* венно. Это поле в любой точке Р(х, у) имеет компоненты: 94
горизонтальную Х(х, у) и вертикальную У (х, у). Поскольку -~ = Х(х, у), а -^• = К(л:1 у), то вектор, связанный с каж- каждой неособой точкой Р(х, у), будет касаться в этой точке фазовой траектории. Если переменную t интерпретировать как время, то век- вектор V можно представить как вектор скорости изображаю- изображающей точки, движущейся вдоль траектории. Таким образом, можно считать, что вся фазовая плоскость заполнена изо- изображающими точками и что каждая фазовая траектория представляет собой след движущейся изображающей точки. В результате приходим к аналогии с плоским движением несжимаемой жидкости. При этом, так как система A22) является автономной, то вектор V в каждой фиксированной точке Р (х, у) не изменяется со временем, и поэтому движение жидкости является стационарным, т. е. установившимся. Фазовые траектории тогда являются траекториями движу- движущихся частиц жидкости, а особые точки О, О', О" (рис. 51) представляют собой неподвижные частицы. Наиболее характерными особенностями движения жидко- жидкости, показанного на рис. 51, являются: 1) наличие особых точек; 2) различное расположение траекторий вблизи осо- особых точек; 3) устойчивость и неустойчивость сособых точек, т. е. реализация двух возможностей: остаются ли с течением времени частицы жидкости в окрестности особой точки или они удаляются от особой точки в другую часть плоскости; 4) наличие замкнутых траекторий, которые в данном случае соответствуют периодическим движениям. Указанные характерные особенности и составляют глав- главную часть фазового портрета или полной качественной кар- картины поведения фазовых траекторий системы общего вида A22). Поскольку, как уже указывалось выше, в общем слу- случае дифференциальные, уравнения не могут быть решены в замкнутой форме, то целью качественной теории обыкновен- обыкновенных дифференциальных уравнений вида A22) является по возможности наиболее полное построение фазового портрета непосредственно по функциям Х(х, у) и У (х, у). Консервативные системы в механике Из практики хорошо известно, что в любой реальной динамической системе энергия рассеивается. Рассеяние {диссипация) энергии обычно происходит в связи с наличи- наличием той или иной формы (вида) трения. Вместе с тем в неко- некоторых конкретных случаях рассеяние энергии бывает на- 95
столько медленным, что им можно пренебречь, если ограни- ограничиться относительно непродолжительным отрезком време- времени. В таких случаях можно считать, что при рассмотрении конкретной физической системы имеет место закон сохране- сохранения энергии, а именно: сумма кинетической и потенциаль- потенциальной энергий постоянна. Системы такого типа называют N^ т ос ^ L- Рис. 52 консервативными. Таким образом, вращающийся земной шар может рассматриваться как консервативная система, если взять промежуток времени в несколько столетий. Если же изучать движение земного шара за несколько миллионов лет, то следует учитывать рассеяние энергии, связанное с приливно-отливиым движением воды в океанах и морях. Простейшим примером консервативной системы являет- является система, состоящая из тела, которое совершает горизон- горизонтальные движения в вакууме под действием двух пружин (рис. 52). Если х обозначает смещение тела массой га из состояния равновесия, а сила, с которой пружины дейст- действуют на тело (восстанавливающая сила), пропорциональна смещению х, то уравнение движения имеет вид m^-\-kx = 0, k>0. Пружины этого типа называются линейными, так как их восстанавливающая сила является линейной функцией х. Если тело массой т движется в среде с сопротивлением и сопротивление (демпфирующая сила), действующее на это тело, пропорционально скорости движения, то уравне- уравнение движения такой кеконсервативной системы будет m^ + cc£ + kx = 0, c>0. A24) Здесь мы имеем линейное затухание, так как демпфирую- демпфирующая сила представляет собой линейную функцию скоро- dx сти Ж' Если / и g являются такими произвольными функциями, что /@)—О, g@)^0, то более общее уравнение ^Ш + 8§ + Пх) = 0 A25) 96
можно интерпретировать как уравнение движения тела массой т под действием восстанавливающей—f(x) и демп- демпфирующей —g-^- сил. В общем случае указанные силы нелинейны, и поэтому уравнение A25) можно рассматри- рассматривать как основное уравнение нелинейной механики. Кратко рассмотрим специальный случай нелинейной консервативной системы, которая описывается уравнением mc§ + f(x) = 0, A26) где демпфирующая сила равна нулю и, следовательно, нет рассеяния энергии. От уравнения A26) можно перейти к автономной системе вида теперь в системе A27) исключить время t, то полу- получим дифференциальное уравнение траекторий системы на фазовой плоскости dx ту ' Последнее уравнение можно переписать так: mydy=—f{x)dx. A29) Тогда, полагая, что х = х0 при t = to,a у = у0, после интегри- интегрирования уравнения A29) в пределах от t0 до t получаем равенство которое можно переписать так: f{l)dl = \myt + ^f{x)dx. A30) и 1 1 /dx\2 Заметим, что -^m\f — -^m \-jA есть формула кинетиче- кинетической энергии динамической системы, а A31) 4 В. Б. Амелькин 87
— формула ее потенциальной энергии. Таким образом, ура- уравнение A30) выражает закон сохранения энергии; ^my* + V{x) = E, A32) где Е = у туо -[- V (х0) — полная энергия системы. Ясно, что уравнение A32) — это уравнение фазовых траекторий системы A27), поскольку оно получено в результате интег- интегрирования уравнения z~*V(x)/ A28). Таким образом, у различным значениям Е на фазовой плоскости соответствуют различ- различные кривые постоянной энергии. Особыми точка- точками системы A27) явля- являются точки Mv (xv, 0), где xv — корни уравнения f(x)=O. Как уже отме- отмечалось, особые точки яв- являются точками равно- равновесия динамической си- системы, описываемой уравнением A26). Из уравнения же A28) сле- следует, что фазовые траек- траектории пересекают ось х под прямым углом, а прямые x=xv —горизон- —горизонтально. Кроме того, уравнение A32) показывает, что фазовые траектории сим- симметричны относительно оси х. В таком случае, если переписать уравнение A32) в виде Рис 53 A33) то можно легко построить фазовые траектории. Действи- Действительно, введем в рассмотрение плоскость (х, z) — «плос- «плоскость баланса энергии» (рис. 53), у которой ось z лежит на той же вертикальной прямой, что и ось у фазовой плоско- плоскости. Затем изобразим график функции z~V(x) и несколько горизонтальных прямых z=E в плоскости (х, z) (одна такая прямая указана на рис. 53). Отметим на рисунке значение разности Е—У(х). После этого для каждого х умножим раз- 98
ность Е—V (х) на 2/т, что с учетом формулы A33) дает воз- возможность нанести соответствующие значения у на фазовую плоскость. Заметим, что поскольку тг — У* то положи- положительное направление вдоль любой траектории определяется движением изображающей точки слепа направо выше оси х и справа налево ниже оси х. Проведенные рассуждения общего характера дают воз- возможность исследовать уравнение движения маятника в среде без сопротивления, которое имеет вид (см. с. 46) сРх Ь / sin л: = 0, A34) где k — положительная постоянная. Поскольку уравнение A34) является частным случаем уравнения A26), то его можно интерпретировать и как урав- уравнение, описывающее прямолинейное движение без трения тела единичной массы под действием нелинейной пружины, где восстанавливающая сила равна —k sin x. В этом случае автономная система, соответствующая уравнению A34), запишется в виде Особыми точками здесь будут @, 0), (±зт, 0), (±2зт, 0), . . ., а дифференциальное уравнение фазовых траекторий систе- системы A35) примет вид dy k sin x dx~~ у Разделяя переменные в последнем уравнении и интегрируя, получим уравнение фазовых траекторий Последнее уравнение есть частный случай уравнения A32), где т—1, а потенциальная энергия, определяемая по фор- формуле A31), задается соотношением V {х) =\f (I) dl = k{\ —cos x). Изобразим теперь в плоскости (х, z) график функции 2— = V(x), а также несколько прямых z=E (на рис. 54 указана только прямая г=£=2/г). Определив значения Е—V(^), 4* G9
мы можем схематически набросать картину поведения тра- траекторий на фазовой плоскости, если воспользоваться соот- соотношением y=±V2[E-V(x)]. Полученный фазовый портрет показывает (рис. 54), что если энергия Е изменяется от 0 до 2k, то соответствующие фазо- фазовые траектории оказываются замкнутыми и уравнение A34) Рис. 54 имеет периодические решения. С другой стороны, если E>2k, то соответствующие фазовые траектории не являются замкнутыми и уравнение A34) периодических решений не имеет. Значению же E=2k на фазовой плоскости соответст- соответствует фазовая траектория, которая разделяет два различных типа движений, т. е. является сепаратрисой. Волнистые фазовые траектории, расположенные вне сепаратрис, соот- соответствуют вращательным движениям маятника, а замкнутые траектории, расположенные в областях, ограниченных сепа- сепаратрисами,— его колебательным движениям. Из рис. 54 видно, что в окрестности особых точек (±2ят, G), гдет=0,1,2,..., поведение фазовых траекторий отличается от поведения фазовых траекторий в окрестности особых точек (±лл, 0), где п=\, 2, ... 100
Существуют различные типы особых точек, и с некоторы- некоторыми из них мы познакомимся ниже. Что же касается послед- последнего примера, то точки (±2зхт, 0), где т=0, 1,2, . . ., от- относятся к особым точкам типа центр. Особые точки (±я/г, 0), где /г=1, 2, . . .,— это точки типа седло. Особая точка автономной дифференциальной системы вида A22) называ- называется центром, если некоторая ее окрестность сплошь запол- заполнена непересекающимися замкнутыми фазовыми траекто- траекториями, окружающими эту точку. Седлом же называют такую особую точку, к которой примыкает конечное число фазо- фазовых траекторий («усов»), разбивающих некоторую окрест- окрестность рассматриваемой точки на области, где траектории ведут себя подобно семейству гипербол, заданных уравнени- уравнением ху=const. Посмотрим теперь, какое влияние оказывает на поведе- поведение фазовых траекторий консервативной системы линейное трение. В этом случае уравнение примет вид d2x , dx . i . г, . Л 3_ + c_. + *sin*=a0l c>o. А тогда если трение достаточно мало, т. е. возможны коле- колебания маятника относительно положения равновесия, то можно показать, что фазовые траектории таковы, как это схематически показано на рис. 55. Если трение не допуска- допускает никаких колебаний маятника относительно положения равновесия, то картина фазовых траекторий будет иметь вид, показанный на рис. 56. Если теперь сравнить фазовый портрет консервативной системы с последними двумя фазовыми портретами некон- неконсервативных систем, то можно заметить, что седла остались седлами (рассматриваем только достаточно малую окрест- окрестность особых точек), а вот в окрестности точек (±2лт, 0), где т=0, 1,2, . ., замкнутые фазовые траектории при слабом трении перетили в спирали, а при сильном трении — в траектории, которые «входят» в особые точки в определен- определенных направлениях. В случае спирали приходим к особой точке типа фокус; в последнем же случае особая точка назы- называется узлом. Особая точка двумерной автономной дифференциальной системы общего вида A22) (если таковая существует) назы- называется фокусом, если существует окрестность этой точки, сплошь заполненная непересекающимися фазовыми траек- траекториями, подобными спиралям, которые «наматываются» на особую точку либо при £->+оо^ либо при £-*—оо. Узлом же называют такую особую точку, в некоторой окрестности 101
которой каждая из фазовых траекторий ведет себя подобно ветви параболы или полупрямой, примыкающей к иссле- исследуемой точке по определенному направлению. Следует отметить, что если какая-либо из консерватив- консервативных систем имеет периодическое решение, то последнее не Рис. 56 может быть изолированным. Более того, если Г — замкну- замкнутая фазовая траектория, соответствующая периодическому решению консервативной системы, то некоторая окрестность траектории Г всегда сплошь заполнена замкнутыми фазо- фазовыми траекториями. Заметим далее, что приведенные выше определения не- некоторых типов особых точек носят чисто качественный, 102
описательный характер. Что же касается аналитических признаков их различия, то, к сожалению, в общем случае систем вида A22) таких критериев не существует, однако для некоторых классов дифференциальных уравнений их удается получить. Простейшим примером является линей- линейная система вида dx 7 1 du — /У v I }> II •— п V I /т it где au bi, «2, b2 — вещественные постоянные. Если матрица коэффициентов этой системы невырожде- невырождена, т. е. если определитель at bi а2 I то начало координат О@, 0) фазовой плоскости — единст- единственная особая точка дифференциальной системы. Считач последнее неравенство выполненным, обозначим через %i и К2 собственные значения матрицы коэффициентов. Тогда можно показать, что 1) если Ai и к2 вещественны и одного знака, то особая точка — узел; 2) если %i и А2 вещественны и разных знаков, то особая точка — седло; 3) если /if и Х2 невещественны и не чисто мнимы, то особая точка — фокус; 4) если %i и А2 чисто мнимы, то особая точка — центр. Отметим, что первые три типа особых точек относят к так называемым «грубым» особым точкам: их характер не меняется при малых возмущениях правых частей исходной дифференциальной системы. Особая же точка центр — «не- «негрубая» точка: ее характер меняется уже при малых возму- возмущениях правых частей дифференциальной системы. Устойчивость точек равновесия и периодических движений Как мы уже знаем, различные по своему типу особые точки характеризуются различным расположением фазовых траекторий в достаточно малой окрестности этих точек. Вместе с тем существует еще одна характеристика — ус- устойчивость точки равновесия, которая позволяет получить дополнительную информацию о поведении фазовых траек- траекторий в окрестности особых точек. Рассмотрим маятник, изображенный на рис. 57. На рисунке показаны два состоя- 103
ния равновесия: а) тело массой т находится в состоянии равновесия в верхней точке; б) тело массой т находится в состоянии равновесия в нижней точке. Первое состояние неустойчиво, а второе устойчиво и вот почему. Если тело массой т находится в состоянии равновесия в верхней точ- точке, то достаточно слегка толкнуть его, чтобы оно начало с возрастающей скоростью отклоняться от положения равно- равновесия и тем самым удаляться от исходного положения; если же тело находится в состоянии равновесия в нижней точке, У Рис. 57 Рис. 58 то после толчка оно будет двигаться с уменьшающейся ско- скоростью, причем чем слабее толчок, тем на меньшее расстоя- расстояние произойдет отклонение от исходного положения. Состояние равновесия физической системы соответству- соответствует особой точке на фазовой плоскости. Малые возмущения неустойчивой точки равновесия приводят к большим откло- отклонениям от этой точки; в случае же устойчивой точки равно- равновесия малые возмущения приводят к малым отклонениям. Отправляясь от таких наглядных интуитивных соображе- соображений, рассмотрим изолированную особую точку системы A22), предполагая для простоты, что она находится в начале ко- координат 0@, 0) фазовой плоскости. Будем говорить, что эта особая точка устойчива, если для любого положитель- положительного числа R существует положительное число r^.R такое, что любая фазовая траектория, выходящая в начальный момент времени t—t0 из точки Р, лежащей в круге х*-\-у'2= =л2, при всех f>t0 будет лежать в круге x2+t/2=JR2 (рис. 58). Не придерживаясь строгой формулировки, можно сказать, что особая точка является устойчивой, если все фазовые траектории, которые в начальный момент времени находят- находятся вблизи особой точки, с течением времени также остаются вблизи этой точки. Далее, особая точка называется асимп- 104
тотшески устойчивой, если она устойчива и если существу- существует круг хг-\гу2=г1 такой, что каждая траектория, которая в момент времени t—U находится в этом круге, при t—>—{-co стремится к началу координат. Наконец, если особая точка не является устойчивой, то ее называют неустойчивой. Особая точка типа центра всегда устойчива (но не асимп- асимптотически). Седло всегда является неустойчивой особой точкой. На рис. 55, где показано поведение фазовых траек- траекторий в случае колебаний маятника в среде с малым трени- трением, особые точки — фокусы — асимптотически устойчивы; особые точки, являющиеся узлами (рис. 56), также асимп- асимптотически устойчивы. Введенное понятие устойчивости точки равновесия явля- является понятием чисто качественным, так как ни о каких свой- свойствах, касающихся характера поведения фазовых траекто- траекторий, здесь не говорится. Что же касается понятия асимпто- асимптотической устойчивости, то по сравнению с понятием прос- простой устойчивости здесь дополнительно требуется, чтобы лю- любая фазовая траектория с течением времени стремилась к началу координат. Однако и в этом случае никаких условий на характер приближения к точке О@, 0) также не наклады- накладывается. Понятие устойчивости (асимптотической устойчивости) играет важную роль в приложениях. Дело в том, что если какой-то агрегат сконструирован без учета устойчивости, то он будет чувствителен даже к самым незначительным внешним воздействиям, а это в конечном счете может при- привести к нежелательным последствиям. Подчеркивая важ- важность понятия устойчивости, известный советский специа- специалист в области математики и механики Н. Г. Четаев писал *): «...если конструируется пассажирский самолет, то его проектным движениям нужно обеспечить известного рода устойчивость, чтобы тем самым получить машину, спокой- спокойную в полете и безаварийную на взлете и посадке. Коленча- Коленчатый вал нужно рассчитать так, чтобы он не поломался от вибраций, какие могут возникнуть в реальных условиях работы двигателя. Чтобы обеспечить артиллерийскому ору- орудию наибольшую меткость и кучность боя, надо строить орудия, снаряды и мины так, чтобы была известного рода устойчивость траекторий и правильность полета снарядов. Можно привести много примеров, но они ничего не при- прибавят к тому, что при решении вопроса о действительных *) Четаев Н. Г. Устойчивость движения.— М.; Наука, 1965.— С. 8—9. 105
< Ж 'Jo \ г > 1»/j 6B / /(£) / Рис. 59 движениях необходимо из возможных решений уравнений останавливаться на решениях, отвечающих устойчивым состояниям, и что в тех случаях, когда желательно избе- избежать в действительности какого-либо решения, разумно путем некоторого изменения в конструкции системы делать отвечающее этому решению состояние движения неустой- неустойчивым». Возвращаясь к маятнику, изображенному на рис. 57, отметим следующий любопытный и в определенной степени неожиданный факт. Оказывается, что верхнее неустойчивое положение равновесия маят- маятника можно сделать устой- устойчивым с помощью вертикаль- вертикальных колебаний точки подвеса А. Более того, не только верх- верхнее положение маятника, но и любое другое его положе- положение (в частности, горизонталь- горизонтальное) можно сделать устойчи- устойчивым за счет вибрации точки подвеса *). Обратимся теперь к не менее важному понятию, чем понятие устойчивости точки равновесия,— к понятию устой- устойчивости периодических движений (решений). Предполо- Предположим, что рассматривается консервативная система, имею- имеющая периодические решения. На фазовой плоскости этим решениям соответствуют замкнутые траектории, сплошь заполняющие некоторую область. Таким образом, всякому периодическому движению консервативной системы соот- яетсгвует движение изображающей точки по некоторой фиксированной замкнутой траектории в фазовой плоскости. В общем случае период обхода изображающих точек по разным траекториям различен. Иначе говоря, период коле- колебаний в консервативной системе будет зависеть от началь- начальных данных. В геометрической интерпретации это означает, что две близкие друг к другу изображающие точки, начав- начавшие в некоторый момент времени t=t0 свое движение, на- например, с оси х, с течением времени разойдутся на некото- некоторое конечное расстояние. Однако может случиться и так, что с течением времени указанные точки и не разойдутся. Чтобы различать эти две возможности, вводится понятие *) Много примеров стабилизации различных видов маятников при- приведено в книге Т. Г. Стрижак «Методы исследования динамических сис- систем типа «маятник» (Алма-Ата: Наука, 1981.—254 с), i 106
устойчивости по Ляпунову периодических решений. Суть его в следующем. Если по заданной сколь угодно малой е-окрестиости движущейся по замкнутой траектории Г точки М *) (рис. 59) можно указать такую движущуюся б (е)-окрестность той же точки М, что всякая изображаю- изображающая точка, лежащая в начальный момент времени в 6(е)-ок- рестности, с течением времени никогда не выйдет из £-ок- рестности, то периодическое решение, соответствующее траектории Г, называют устойчивым по Ляпунову. Если периодическое решение не является устойчивым по Ляпуно- Ляпунову, то его называют неустойчивым по Ляпунову. Говоря о неустойчивых по Ляпунову периодических решениях, следует заметить, что они все-таки обладают определенным видом устойчивости: так называемой орби- орбитальной устойчивостью, заключающейся в том, что при малых изменениях начальных данных изображающая точка переходит с одной фазовой траектории на другую, но сколь угодно близко лежащую к первоначально рассматриваемой. Примером устойчивых по Ляпунову периодических реше- решений могут быть решения, которые возникают, например, при рассмотрении дифференциального уравнения, описы- описывающего горизонтальные движения тела массой т под действием двух линейных пружин в вакууме (рис. 52). Примером же неустойчивых по Ляпунову, но орбитально устойчивых периодических решений являются решения дифференциального уравнения A34), описывающего дви- движение кругового маятника в среде без сопротивления. В первом случае период колебаний не зависит от на- начальных данных и вычисляется по формуле Т = 2n]/rm/k. Во втором случае период колебаний зависит от начальных дан- данных и выражается, как мы знаем, через эллиптический ин- интеграл первого рода, взятый в пределах от 0 до я/2. Отметим, наконец, что вопрос об устойчивости по Ляпу- Ляпунову периодических движений непосредственно связан с изучением вопроса об изохронных колебаниях **). Энергетические функции Интуитивно ясно, что если полная энергия некоторой физической системы имеет минимум в точке равновесия, *) Под е-окрестностыо точки М понимается круг радиусом & g цент- центром d точке М. **) См., например, книгу В. В. Амелькина, Н. А. Лукашевича, А. П. Садовского «Нелинейные колебания в системах второго порядка» (Минск: Изд-во БГУ, 1982.—208 с). 107
то эта точка является точкой устойчивого равновесия. Вы- Высказанная идея лежит в основе одного из двух методов изу- изучения различных задач устойчивости, предложенных из- известным русским математиком и механиком А. М. Ляпуно- Ляпуновым. Этот метод обычно называют прямым или вторым методом Ляпунова *). Прямой метод Ляпунова проиллюстрируем на системах пида A22) в случае, когда начало координат для этой си- системы является особой точкой. Итак, пусть Г — фазовая траектория системы A22). Рассмотрим функцию V=V(x, у), которая является непре- dV рывпои вместе со своими частными производными у- и -г- в области фазовой плоскости, содержащей траекто- траекторию Г. Если изображающая точка (x{i), y(l)) движется вдоль кривой Г, то вдоль этой кривой функцию V(x, у) мо- можно рассматривать как функцию /. А тогда скорость изме- изменения функции V(x, у) вдоль траектории Г определяется равенством dV dVdx . dVdy dV v, , . dV „, ч где Х(х, у) и Y(x,-y)—это правые части системы A22). Формула A36) играет существенную роль в реализации прямого метода А. М. Ляпунова, для практического ис- использования которого важны следующие понятия. Пусть V=V(x, у) непрерывна вместе со своими частными dV dV „ л п производными -Q- и -5- в некоторон области О, содер- содержащей начало координат, причем V@, 0)—0. Тогда эту функцию называют определенно положительной (отрица- (отрицательной), если во всех точках области G, исключая начало координат, выполняется неравенство V(x, у)~>0 (<0). Если же в точках области G имеет место нестрогое неравенство V(x, y)^0 (<Ю), то функцию V(x, у) называют знакополо- знакоположительной (знакоотрицательной). Так, например, функция V, определяемая равенством V{x, y)—x2-\-y2 и рассматри- рассматриваемая на плоскости (л:, у), будет функцией определенно положительной. Функция же V(x, y)=x2 будет на плоскости *) Много интересных примеров исследования дифференциальных моделей на устойчивость можно найти в книге Н. Руша, П. Абетса, М. Лалуа «Прямой метод Ляпунова в теории устойчивости» (М.: Мир, 1980). 108
(*» У) функцией знакоположительной, ибо она обращается в нуль на всей оси у. Если функция У(х, у) является определенно положи- положительной, то, потребовав, например, чтобы во всех точках области и\О выполнялось неравенство (з—) + (#"") Ф®, уравнение z=V(x, у) можно интерпретировать как уравне- уравнение поверхности, которая похожа на параболоид, касаю- касающийся плоскости (х, у) в точке О@, 0) (рис. 60). Вообще же Рис. 60 Рис. 61 уравнение z=V(x, у) при определенно положительной функ- функции V можег задавать поверхности и более сложной струк- структуры. Одна из таких поверхностей показана на рис. 61, где сечение ее плоскостью z=C дает не кривую, а целое кольцо. Если положительно определенная функция У{х, у) обладает тем свойством, что функция у) = у) дх v . dV (х, у) у) A37) является знакоотрицательной, то V называют функцией Ляпунова или энергетической функцией для системы A22). Заметим здесь, что в силу равенства A36) требование того, чтобы функция W была знакоотрицательной, означает, что и, следовательно, функция V не возрастает вдоль траекто- траектории Г в окрестности начала координат. Приведем один из результатов А. М. Ляпунова, который состоит в следующем: если для системы A22) существует энергетическая функция V {х, у), то начало координат, яв- являющееся особой точкой, устойчиво. Если определенно поло- 109
оюительная функция такова, что функция W, определяемая равенством A37), определенно отрицательна, то начало координат устойчиво асимптотически. Покажем на примере, как практически применяется только что приведенный результат. Рассмотрим уравнение движения тела единичной массы под действием пружины, которое на основании соотношения A24) записывается в виде fg + c^ + fo^O, c>0. A38) Напомним, что в уравнении A38) постоянная с>0 харак- характеризует вязкость среды, в которой движется тело, а посто- постоянная £>0 характеризует свойства пружины. Автономная система, соответствующая уравнению A38), имеет вид зН*' %t=-kx-cy. A39) Для этой системы начало координат фазовой плоскости (х, у) является единственной особой точкой. Кинетическая энергия тела единичной массы в данном случае равна #2/2, а потенциальная энергия (энергия, накопленная пружиной) определяется равенством о Отсюда полная энергия системы определяется равенством * \^> У)== "о* У~ ~* о" RX-» Легко видеть, что функция V, заданная соотношением A40), является функцией определенно положительной. А тогда поскольку в данном случае ^Х{х, У)-\гщ-У(х, y) = kxy+y{—kx—q/) = — а/2<0, то эта функция будет для системы A39) энергетической функ- функцией, и поэтому особая точка О@, 0) устойчива. В рассмотренном примере мы быстро получили резз^ль- тат. Однако это не всегда задается. Дело в том, что приведен- приведенный критерий А. М. Ляпунова носит чисто качественный характер. Исходя из него, не всегда можно построить энер- энергетическую функцию, даже если мы и знаем, что она су- 110
шествует. Это существенно усложняет выяснение вопроса об устойчивости в каждом конкретном случае. Следует иметь в виду, что приведенный критерий А. М. Ляпунова необходимо рассматривать как принцип получения эффективных признаков устойчивости. Этому вопросу посвящено много исследований, и в настоящее время получен ряд интересных результатов *). Простые состояния равновесия Уже из динамической интерпретации дифференциаль- дифференциальных уравнений второго порядка ясно, что исследование характера состояний равновесия или, что то же самое, осо- особых точек дает ключ к выяснению поведения интегральных кривых. Понятно также, что так как, вообще говоря, дифферен- дифференциальные уравнения не интегрируются в замкнутой форме, то необходимо иметь критерии, которые дапали бы возмож- возможность по виду исходного дифференциального уравнения ответить на вопрос о типе особой точки. К сожалению, в общем случае этот вопрос решить, как правило, очень слож- сложно, но все-таки всегда можно выделить классы дифферен- дифференциальных уравнений, когда эта задача решается довольно просто. Ниже на примере исследования движения тела еди- единичной массы под действием линейных пружин в среде с линейным трением мы покажем, как используются некото- некоторые факты качественной теории дифференциальных урав- уравнений. Но предварительно остановимся на рассмотрении системы вида A22). Оказывается, что при выяснении типа особой точки такой системы простейшим является случай, когда определитель Якоби (или якобиан) J(x, */) = отличен от нуля в рассматриваемой точке. Если (л:*, у*) — особая точка системы A22) и если «/* = —J (х*, у*)Ф0, то тип особой точки, называемой в данном случае простой особой точкой, зависит в основном от знака константы J*. Так, если У*<0, то особая точка (х*, у*) является седлом. В случае же J*>0 особая точка может дХ дх 3Y дх дХ ду dY ду *) См., например, книгу Е. Л. Барбашина «Функции Ляпунова» (М.: Наука, 1970.—240 с).
быть либо центром, либо узлом, либо фокусом. Центр мо- может быть лишь в случае, когда дивергенция (от позднела- тинского divergentia — расхождение) г», ч дХ , дУ обращается в нуль в особой точке, т. е. когда D*—D(x*t у*)=0. Следует, однако, заметить, что выполнение условия D*=0 недостаточно в общем случае для того, чтобы в осо- особой точке (#*, у*) был центр. Для наличия центра требует- требуется выполнение дополнительных условий, где участвуют частные производные высших порядков. И этих условий, вообще говоря, бесчисленное множество. Вместе с тем, если функции X и Y оказываются линейными по переменным х и у, то условие D*=0 становится уже достаточным для того, чтобы особая точка (л;*, у*) была центром. Если J*X), но особая точка (л;*, у*) не является цент- центром, то ее достаточно малая окрестность сплошь заполнена траекториями, либо стремящимися к этой особой точке по спиралям, либо входящими в нее в определенном направ- направлении. При этом, если D*>0, то особая точка достигается при £->—оо и является неустойчивой. Если же D*<0, то особая точка достигается при £->+<х> и оказывается устой- устойчивой. Если фазовые траектории, достигающие особой точ- точки, являются спиралями, то эта точка есть фокус. Если же все интегральные кривые «входят» в особую точку с определенной касательной, то такая точка есть узел (рис. 62). Независимо от знака якобиана J* касательные к траек- траекториям дифференциальной системы A22) в особой точке (л:*, у*) определяются из так называемого характеристи- И2
ческого уравнения dY (**, у") - ■ BY (**, у») ~ У_--Ах дЛ 1 пл\\ где х = л-—х*, у=у—у*. A42) Если X и Y содержат линейные члены, то частные про- производные в уравнении A41) представляют собой коэффициен- коэффициенты при хи у в системе, получающейся из системы A22) пос- после замены A42). Уравнение A41) — это однородное уравнение. Поэтому если ввести в рассмотрение угловой коэффициент К—у/х так называемых исключительных направлений, то для на- нахождения Я имеем квадратное уравнение XlK + iXl-YDK-Yl-O. A43) Дискриминант этого уравнения Д = (XI + KJ)?_4J* = D*2 — AJ\ Поэтому если J*<0, что соответствует случаю седла, ть уравнение A43) всегда дает два вещесгвенных исключитель- исключительных направления. Если же У*>0, то либо нет веществен- вещественных исключительных направлений, либо их два или одно кратное. В первом из указанных случаев особая точка будет либо центром, либо фокусом. Существование вещественных исключительных направ- направлений означает (при условии J*>0) наличие особой точки типа узла. В частности, если существуют два вещественных исключительных направления, то можно доказать, что су- существуют точно две траектории (по одной с каждой сторо- стороны), касающиеся в особой точке одной из исключительных прямых, в то время как все другие траектории «входят» в особую точку, касаясь другой такой прямой (рис. 62, а). Далее, если Л—0, а уравнение A43) не есть тождество, то мы имеем одну исключительную прямую. Картина пове- поведения траекторий в данном случае показана на рис. 62, б, и она получается из предыдущего случая, когда два исклю- исключительных направления совпадают. Особая точка разделя- разделяет исключительную прямую на две полупрямые U и 12. Окрестность же особой точки разбивается на два сектора, один из которых сплошь заполнен траекториями, «входя- «входящими» в особую точку, касаясь полупрямой lit а второй 113
сплошь заполнен траекториями, «входящими» в особую точку, касаясь полупрямой /2. Граница между двумя сек- секторами состоит из двух траекторий, одна из которых каса- касается в особой точке полупрямой lit а вторая — полупря- полупрямой /2. Если в уравнении A43) все коэффициенты обращаются п нуль, то мы имеем тождество, и тогда все прямые, прохо- проходящие через особую точку, являются исключительными и при этом существуют точно две траектории (по одной с каж- каждой стороны), касающиеся каждой из этих прямых в особой точке. Эта особая точка (рис. 62, в) подобна точке с одной двойной вещественной исключительной прямой. Движение тела единичной массы под действием линейных пружин в среде с линейным трением Выше было показано, что дифференциальное уравнение движения тела единичной массы под действием линейных пружин в среде с линейным трением имеет вид S+4<+fa=°- A44) Для полноты исследования дифференциальной модели A44) мы не будем фиксировать направление действия сил dx c-jj- и kx. Уравнению A44), как уже указывалось выше, можно поставить в соответствие автономную систему вида dx dy , ЪГ^У ■Т="~~кх~~сУ А тогда если исключить тривиальный случай k—О, что мы пока и предполагаем, то дифференциальная система A45) имеет в начале координат изолированную особую точку. Система A45) является специальным видом общей системы A22). В нашем конкретном примере Далее, якобиан J (х, y)=k, а дивергенция D (х, у) ——с. Характеристическое уравнение принимает вид где Д=с2—4k v— дискриминант этого уравнения. Отсюда, в соответствии с результатами, изложенными на с. 111—114^ приходим к следующим возможным случаям: 114
1. Если k<CO, то особая точка является седлом с одним положительным и одним отрицательным исключительными направлениями. Поведение фазовых траекторий показано на рис. 63, где можно выделить три различных типа движе- движений. Когда начальные условия соответствуют точке а, Рис. 63 Рис. 64 в которой вектор скорости направлен к началу координат и скорость достаточно большая, изображающая точка будет двигаться по траектории по направлению к особой точке с убывающей скоростью; затем, пройдя мимо начала коор- координат, изображающая точка будет удаляться от него уже с возрастающей скоростью. Если начальная скорость убы- убывает до некоторого критического значения, соответствую- соответствующего точке Ь, то изображающая точка будет приближаться к особой точке с убывающей скоростью и «достигнет» на- начала координат за бесконечный промежуток времени. На- Наконец, если начальная скорость меньше, чем критическое значение, и соответствует, например, точке С, то изобра- изображающая точка будет двигаться к началу координат с убы- убывающей скоростью, обращающейся в нуль на некотором расстоянии Xi от особой точки. В точке (л*, 0) вектор ско- скорости меняет направление, и изображающая точка удаляет- удаляется от начала координат. Если фазовая точка, которая соответствует начальному состоянию динамической системы, находится в любом из трех оставшихся квадрантов, то интерпретация движений уже очевидна. 2. Если /г>0, то J*>0 и тип особой точки зависит от значения коэффициента с. Рассмотрим представляющиеся здесь возможности. 2а) Если с=0, т. е. если трение отсутствует, то особая точка является центром (рис. 64). Движения являются 115
периодическими, и их амплитуда зависит от начальных ус- условий. 26) Если С>0, т. е. если имеет место случай положитель- г, дХ . dY иого затухания, то дивергенция & = *—Ьз— отрицатель- отрицательна и, следовательно, изображающая точка движется по Рис. 65 Рис. 66 траекториям по направлению к началу координат и дости- достигает его за бесконечный промежуток времени. Более точно: 2бх) Если Д<0, т. е. если c*<.4k, то особая точка оказы- оказывается фокусом (рис. 65) и, значит, динамическая система совершает затухающие колебания относительно состояния равновесия с убывающей амплитудой. 2б2) Если Д=0, т. е. если c2=4k, то особая точка —узел с одним исключительным направлением отрицательного наклона (рис. 66). Движение в данном случае является апериодическим, и оно соответствует так называемому критическому затуханию. 26j) Если Д>0, т. е. если cs;>4&, то особая точка явля- является узлом с двумя исключительными направлениями от- отрицательного наклона (рис. 67). Качественно движение динамической системы таково, как и в предыдущем слу- случае, и оно соответствует затуханию колебаний. Из приведенных результатов следует, что в случае, когда £>0 и fc>0, т. е. когда трение положительно, а восстанав- восстанавливающая сила является притягивающей, динамическая система стремится к состоянию равновесия и ее движение устойчиво. llfi
2в) Если с<0, т. е. если имеет место случай отрицатель- отрицательного затухания, то качественная картина поведения фазо- фазовых траекторий подобна той, которая имеет место в слу- случае 26). Единственное различие состоит лишь в том, чтс *х\ Рис. 67 динамическая система в данном случае перестает быть ус- устойчивой. Диаграмма на рис. 68 содержит сводку полученных ре- результатов и указывает характер (тип) особой точки в за- Рис 68 висимости от значения параметров с и k. Заметим, что при- приведенная диаграмма может быть использована и как сводка результатов исследования типов особых точек системы A22), 417
когда якобиан J*=?M) при значениях с——D*, ft—J*. Од- Однако в общем случае равенство с=0 еще не означает, что система A22) имеет центр, а равенство k—О еще не означает, с>0 Рис. 69 что система общего вида не имеет особой точки. Указанные случаи — случаи сложных особых точек, рассматриваемых ниже. Возвращаясь к рассматриваемой в этом параграфе динамической системе, отметим, что если /г=0 {сФО), то Рис. 70 соответствующая уравнению A44) автономная система A45) примет вщ dx __ dy It~y' dT~~cy' Отсюда следует, что прямая у=0 сплошь заполнена особы- особыми точками. Поведение же фазовых траекторий показано на рис 69. Наконец, если k—c=Q, то уравнение A44) таково; ^1 йп 118
Соответствующая ему система записывается так! dx dy р. Здесь, как и в предыдущем случае, ось х сплошь заполнена особыми точками. Картина поведения фазовых траекторий показана па рис. 70. Адиабатический поток идеального газа в канале переменного диаметра Изучение потоков сжимаемых вязких сред имеет боль- большое значение для практики. В частности, такие потоки воз- возникают вокруг крыла и фю- фюзеляжа самолета; с ними свя- связана работа паровых и газо- газовых турбин, воздушно-реак- воздушно-реактивных двигателей и ядерных реакторов. Ниже мы остановимся на Рис. 71 рассмотрении потока идеаль- идеального газа, проходящего через сопло с изменяющейся пло- площадью поперечного сечения А (рис. 71) и имеющего удель- удельную теплоемкость ср. Поток будет трактоваться как одно- одномерный, т. е. все свойства потока предполагаются единооб- единообразными на одном и том же поперечном сечении сопла. Трение в пограничном слое обуславливается тангенциаль- тангенциальным напряжением т, которое задается равенством т=--<7ри2/2, A46) где q — коэффициент трения, зависящий в основном от чис- числа Рейнольдса, ко предполагаемый постоянным вдоль соп- сопла, р — плотность массы потока, v — скорость потока. На- Наконец, мы предположим, что выполняются условия адиаба- тичности, которые исключают лобовое сопротивление внут- внутренних тел, горение, химические изменения, испарение или конденсацию и т. д. Одним из основных уравнений, описывающих поток, является известное уравнение неразрывности, которое в рассматриваемом случае записывается в виде w = pAv, A47) где скорость изменения массы потока w постоянна. Из равенства A47) следует, что ~~!-^Ч-- = 0. A48) р ' А ' v v ' ]!9
Обращаясь к уравнению для энергии установившегося потока, отметим, что в общем случае такое уравнение свя- связывает действие внешней работы и внешних тепловых воз- воздействий с возрастанием в потоке энтальпии (теплосодержа- (теплосодержания), кинетической и потенциальной энергий. В рассмат- рассматриваемом случае поток является адиабатическим, поэтому уравнение для энергии записывается в виде О = w {h + dh)—wh + w [vV2 + d (u?/2)] —wv2/2, или /2) = 0, A49) где h — энтальпия потока (термодинамический потенциал) при абсолютной температуре Т. Но в уравнении A49) dh= =cpdT, и поэтому уравнение энергии потока можно записать в виде l2) = 0. A50) Обратимся к выводу уравнения (количества) движения потока. Заметим здесь, что для задач, касающихся устано- установившихся потоков, обычно используется второй закон Нью- Ньютона. Считая, что угол расходимости стенок сопла являет- является малым, уравнение движения можно записать в 'виде или — Adp—dAdp—TdA = wdv, A51) где р — статическое давление. Слагаемое dA dp в уравнении A51) имеет более высокий порядок малости, чем остальные члены, и поэтому всегда можно считать, что уравнение движения потока имеет вид — A dp—х dA = wdv. A52) Обозначая через D гидравлический диаметр, заметим, что его изменение вдоль оси сопла определяется функцией F такой, что D=F(x)~, где л: ■— координата вдоль оси сопла. Из определения гидравлического диаметра следует, что dA 4dx Отмечая, что ри2/2 = ^рМ2, где v — коэф/фиииеТ1т удельной теплоемкости, ?» ЛЛ—число Маха, формулу A46) можно пе^сп.ис; ч, ь и.де т-- .'?ГМЧ. . A54)
А тогда, учитывая равенства A47), A53) и A54), для урав- уравнения движения потока A52) получаем предегавление -^-|--^-j = 0. A55) Обозначив квадрат числа Маха через у и воспользовавшись соотношениями A48), A50) и A55), после алгебраических комбинаций и преобразований придем *) к дифференциаль- дифференциальному уравнению аи *y(l+lt:r-y)iviy-F'(x)) d~ )F) Ее. dx (y)() где штрих обозначает дифференцирование по х. Знаменатель правой части уравнения A56) обращается в нуль при у—I, т. е. когда число Маха М=1. Это означает, что интегральные кривые последнего уравнения пересека- пересекают так называемую звуковую линию с вертикальными каса- касательными. Так как правая часть уравнения A56) меняет знак при переходе звуковой линии, то интегральные кри- кривые «поворачиваются» и точки перегиба исключаются. Из физического смысла явления следует, что вдоль интеграль- интегральных кривых значение х должно непрерывно возрастать. Следовательно, сечение, на котором интегральные кривые пересекают звуковую линию с вертикальными касатель- касательными, должно быть сечением выхода из сопла. Таким обра- образом, переход от дозвукового потока к сверхзвуковому, а также обратный процесс, может происходить внутри сопла только через особую точку с вещественными исключитель- исключительными направлениями, т. е. через седловую или узловую точку. Координаты особых точек уравнения A56) задаются равенствами которые означают, что указанные точки располагаются в расходящейся части сопла. Седловая точка появляется в том случае, когда якобиан «/*<С0, т. е. когда F"(x*y>0. Поскольку q — достаточно малая константа, то седловая точка появляется возле горловины сопла. Узловая же точ- *) К е s t i п ]., Zaremba S. К. One-dimensional high-speed flows. Flow patterns derived for thp flow of gases through nozzles, inclu- including compressibility and viscosity effects// Aircraft Engin.—1953.—V. 25, № 292.—P. 172—175, 179. 121
ка возникает лишь при выполнении условия F"(x*)<0. Та- Таким образом, узловая точка образуется в той части сопла, которая следует за точкой перегиба в его профиле или, на практике, на некотором расстоянии от горловины сопла при условии, что профиль содержит точку перегиба. Из характеристического уравнения видно, что два исключительных направления имеют угло- угловые коэффициенты, противоположные по знаку в случае седла и совпадающие (отрицательные) по знаку в случае узла. Последнее означает, что только седловая точка допус- допускает переход как от сверхзвуковых к дозвуковым, так и от Рис. 72 дозвуковых к сверхзвуковым скоростям (рис. 72). Случай узловой точки (рис. 73) говорит о том, что здесь допускает- допускается непрерывный переход только от сверхзвукового к до- дозвуковому потоку. Далее, поскольку уравнение A56) не может быть про- проинтегрировано в замкнутой форме, то для дальнейших исследований необходимо применять методы численного интегрирования. В связи с этим отметим, что построение четырех сепаратрис седла как интегральных кривых лучше начинать с учетом того, что сама особая точка является как бы точкой, из которой эги интегральные кривые вы- выходят. И такое построение действительно возможно, ибо из характеристического уравнения мы знаем направление двух касательных в особой точке 5 (л;*, 1). Если это заме- замечание не принять во внимание и начать следить за движе- движением, скажем, с начальными точками а и Ъ на рис. 72, кото- которые расположены по разные стороны от сепаратрисы, то соответствующие точки будут двигаться вдоль кривых а и 122
р, которые сильно расходятся и, следовательно, не дают никакой информации об интегральной кривой, которая «входит» в особую точку S. С другой стороны, если по ин- интегральной кривой двигаться, «выходя» из особой точки S, считая начальным сегментом интегральной кривой сег- сегмент исключительной прямой, то погрешность может быть минимизирована, если учесть сходимость интегральных кривых в направлении убывания переменной х. На рис. 72 показана картина расположения интеграль- интегральных кривых в окрестности особой точки. Прямая, которая проходит через точку xt горловины сопла, соответствует значениям, при которых числитель в правой части уравне- уравнения A56) обращается в нуль, указывая на наличие точек экстремума. Точки равновесия высшего порядка В предыдущих параграфах были исследованы типы осо- особых точек, которые возникают в случае, когда якобиан J*=f£=O. В случае же, когда, например, все частные производ- производные функций X и Y, фигурирующих в представлении пра- правых частей системы A22), обращаются в нуль до порядка п включительно, в окрестности особой точки возможно бес- бесконечно много картин поведения фазовых траекторий. Вместе с тем, если исключить из рассмотрения точки рав- равновесия типа центра и фокуса, то оказывается, что окрест- окрестность особой точки может быть разбита на конечное число секторов, принадлежащих трем стандартным типам. Это гиперболические, параболические и эллиптические секторы. К описанию этих секторов мы и приступим ниже, но пред- предварительно сделаем ряд упрощающих исследование пред- предположений. Прежде всего будем считать, что начало координат пере- перенесено в особую точку, т. е. х*=у*=0; правые части сис- системы A22) могут быть представлены в виде * (х, У) = Хп (д, у) + Ф (х, у), Y Y Ч{х, у), где Хп и Уп — однородные относительно переменных хну полиномы степени п (один из которых может обращаться в нуль тождественно), а функции Ф и W имеют в окрестности начала координат непрерывные частные производные пер- первого порядка. Кроме того, предполагается, что в той же окрестности начала координат ограничены функции, J23
определяемые выражениями Ф (*. у) Фх(х, у) Фу (х, у) (*, у) ^ (*, у) Уу (х, у) При этих предположениях могут быть доказаны следующие утверждения: 1. Любая траектория системы уравнений A22) с правы- правыми частями вида A57), «входящая» в начало координат с определенной касательной, касается одной из исключитель- исключительных прямых, задаваемых уравнением хУ„{х, у)—уХп(х, у) = 0. A58) Поскольку функции Хп и Yn являются однородными, то уравнение A58) может быть переписано как уравнение от- относительно углового коэффициента %~у!х. При этом исклю- исключительная прямая называется особой, если на этой пря- прямой выполняются равенства хп (х> У) = Уп (я, у) = 0. Примером таких прямых являются прямые, показанные на рис. 62. Прямые же, которые определяются уравнением A58) и которые не являются особыми, называются регуляр- регулярными. 2. Поведение фазовых траекторий системы A22) в окрест- окрестности одного из двух лучей, «выходящих» из начала коор- координат и образующих вместе исключительную прямую, может быть исследовано путем рассмотрения достаточно ма- малого круга (с центром в начале координат), в котором вы- выбирается сектор, ограниченный двумя радиусами, располо- расположенными по обе стороны от полупрямой и достаточно близ- близкими к ней. Такой сектор обычно называют нормальной областью. Более точно, в случае регулярной исключительной пря- прямой, которая соответствует линейному множителю в урав- уравнении A58), нормальная область, рассматриваемая в круге достаточно малого радиуса, принадлежит одному из двух типов: притягивающему или отталкивающему. 2а) Притягивающая нормальная область характеризует- характеризуется тем, что каждая траектория, проходящая через эту об- область, достигает начала координат в направлении касатель- касательной, совпадающей с исключительной прямой (рис. 74). 26) Отталкивающая нормальная область характеризует- характеризуется тем," что только одна фазовая траектория, проходящая 124
через эту область, достигает особой точки в направлении касательной, совпадающей с исключительной прямой. Все остальные фазовые траектории системы A22), которые вхо- входят в нормальную область через границу круга, покидают эту область, проходя через один из радиусов, ее ограничи- ограничивающих (рис. 75). Обратим внимание на следующее обстоятельство. Если круг с центром в начале координат достаточно мал, то рас- рассмотренные два типа нормальных областей могут быть Рис. 74 Рис. 75 классифицированы по поведению вектора (X, Y) на грани- границе рассматриваемой области. При этом поведение вектора (X, Y) можно отождествить с поведением вектора (Хп, Yn). Более того, можно показать, что если (как мы и предпола- предполагаем) заданное исключительное направление не соответствует кратному корню характеристического уравнения, то вектор, рассматриваемый на одном из радиусов, ограничивающих нормальную область, направлен либо внутрь области, либо наружу. А тогда если в первом случае вектор, рассматривае- рассматриваемый на той части границы нормальной области, которая со- состоит из дуги окружности, направлен также внутрь области, а во втором случае направлен наружу, то нормальная область является притягивающей. Если же имеет место противо- противоположная ситуация, то нормальная область является от- отталкивающей. Следует отметить, что в любом случае век- вектор, рассматриваемый на той части границы нормальной области, которая состоит из дуги окружности, всегда на- направлен либо внутрь области, либо наружу, ибо он почти параллелен радиусу. Нормальные области, соответствующие особым исключи- исключительным направлениям или кратным корням характеристи- характеристического уравнения, имеют более сложную природу, но вви- ввиду того, что они в какой-то мере представляют собой явле- явление исключительное, мы их здесь описывать не будем. Если теперь обратиться, например, к седловой особой точке, то можно отметить, что она допускает четыре оттал- отталкивающие нормальные области. В окрестности же особой 125
точки типа узла имеются две притягивающие и две оттал- отталкивающие нормальные области. 3. Если существуют вещественные исключительные пря- прямые, то окрестность особой точки может быть разбита на конечное число секторов, каждый из которых ограничен двумя фазовыми траекториями си- системы A22), которые «входят» п на- начало координат с определенными касательными. При этом каждый из таких секторов принадлежит од- одному из следующих трех типов: За) Эллиптический сектор (рис. 76) содержит бесконечно много фа- фазовых кривых в виде петель, про- проходящих через начало координат и касающихся с каждой стороны границы сектора. 36) Параболический сектор (рис. 77) заполнен фазовыми кривыми, которые связывают (соединяют) особую точку с границей окрестности. Рис. 76 Рис. 77 Рис. 78 Зв) Гиперболический сектор (рис. 78) заполнен траекто- траекториями, которые подходят к границе окрестности в обоих направлениях. Более точно: 4а) Эллиптические секторы образуются между двумя траекториями, принадлежащими двум последовательным нормальным областям, которые обе являются притягиваю- притягивающими. 46) Параболические секторы образуются между двумя фазовыми кривыми, принадлежащими двум последователь- последовательным нормальным областям, одна из которых притягиваю- притягивающая, а другая отталкивающая. Все фазовые кривые, кото- которые проходят через эту область; касаются в особой точке 126
исключительной прямой, которая определяет притягиваю- притягивающую область. 4в) Гиперболические секторы образуются между двумя траекториями, принадлежащими двум последовательным отталкивающим нормальным областям. Легко, например, отличить четыре гиперболических сектора в седловой точке и четыре параболических сектора в узловой точке. Эллиптические секторы не появляются в случае простых особых точек, когда якобиан /*=£0. Если особая точка не допускает существования вещест- вещественных исключительных направлений, то фазовые кривые в ее окрестности обязательно имеют структуру центра или фокуса *). Преобразование обратными радиусами и однородные координаты Выше были описаны методы выяснения локального поведения фазовых траекторий дифференциальных систем вида A22) в окрестности особых точек. И хотя во многих случаях вся требуемая информация может быть получена указанным выше образом, может, однако, возникнуть не- необходимость исследовать поведение траекторий в бесконеч- бесконечно удаленных частях фазовой плоскости, когда x^-hy^-t-oo. Простейший путь изучения асимптотического поведения траекторий дифференциальной системы A22) состоит в том, чтобы соответствующим преобразованием первоначальной дифференциальной системы, таким, например, как инвер- инверсия, определяемая подстановками ввести в рассмотрение точку в бесконечности. Геометричес- Геометрически это преобразование представляет собой так называемое преобразование обратными радиусами, которое начало координат переводит в точку на бесконечности и наоборот. Любую же конечную точку М (х, у) фазовой плоскости пре- преобразование A59) переводит в точку М'(£, т)) той же плос- *) Методы различения центра от фокуса рассматриваются, напри- например, в книге В. В. Амелькина, Н. А. Лукашевича, А. П. Садовского «Нелинейные колебания в системах второго порядка» (Минск: Изд-во БГУ, 1982). 127
кости, причем точки М и Мг лежат на одном луче, выходя- выходящем из начала координат, и для них выполняется равенст- равенство ОМ-ОМ'—г* (рис.79). Хорошо известно, что при таком преобразовании окружности переходят в окружности (пря- (прямые линии рассматриваются как окружности, проходящие через точку в бесконечности). В частности, прямые линии, проходящие через начало координат, инвариантны относи- относительно преобразования A59). Отсюда угловые козффици- J4 Рис. 79 Рис. 80 енты асимптотических направлений являются угловыми коэффициентами касательных в новом начале координат £=т]=0. Отметим, что в большинстве случаев новое начало координат будет особой точкой. Причины этого обстоя- обстоятельства будут обсуждены ниже. Что же касается того, как построить в первоначальной плоскости (х, у) кривую, которая имеет определенное асимп- асимптотическое направление, т. е. которая имеет определенную касательную в начале координат новой плоскости (£, т]), то построение такой кривой можно начать с плоскости (£» К\)г а точнее, с рассмотрения этой кривой, скажем, в единичном круге в плоскости (£, ц). Дело в том, что посколь- поскольку единичная окружность подстановками- A59) преобразует- преобразуется в себя, то всегда можно найти точку пересечения кривой с единичной окружностью в плоскости {х, у); дальнейшее же исследование уже проводится обычным образом. Отметим далее, что пополнение плоскости (#, у) одной точкой в бесконечности топологически эквивалентно инвер- инверсии стереографической проекции (рис. 80), при которой точки сферы отображаются на плоскость, касательную к сфере в точке S. Центр проекции N диаметрально противо- противоположен точке S. Ясно, что центр проекции N соответству- соответствует бесконечно удаленной точке плоскости (х, у). Обратно, 128
если мы отооражаем плоскость па сферу, то векторное поле па плоскости переходит в векторное поле на сфере и точ- точка в бесконечности может оказаться особой точкой ня сфере. Следует, однако, заметить, что хотя метод обратных радиусов и является полезным, он все-таки становится гро- громоздким и малоудобным, когда особая точка в бесконечности имеет сложную структуру. В таких случаях используегся другое, более удобное преобразование плоскости (х, у) путем введения так называемых однородных координат х = \lz% у = ц/z. При таком преобразовании каждой точке (х, у) соответст- соответствуют тройки действительных чисел (|, ц, z), не равных одно- одновременно нулю, причем не делается различия между трой- тройками (£, г), z) и (kl, kit], kz) при любом действительномk^Q. Рис. 61 Если точка (х, у) конечная, то z=^=0. В случае же г—0 полу- получаем прямую в бесконечности. Плоскость (х, у), пополнен- пополненная прямой в бесконечности, называется проективной плоскостью. Ма такой прямой может существовать несколь- несколько особых точек, и их природа будет, как правило, проще чем природа одной особой точки, вводимой посредством пре- преобразования обратными радиусами. Если теперь рассмотреть пучок прямых и описать во- вокруг его центра сферу, например, радиуса 1,то каждая пря- прямая пучка пересечет эту сферу в двух диаметрально проти- противоположных точках. Отсюда следует, что любая точка проективной плоскости отображается взаимно однозначно и иепре;.»ы:«1'о на пару диаметрально противоположных то- точек едини 4яой сферы. Таким образом, проективную плос- 5 В. В л*а,1„кин 12Э
кость можно рассматривать как совокупность всех пар диа- диаметрально противоположных точек единичной сферы. Что- Чтобы представить себе проективную плоскость, достаточно поэтому рассмотреть, например, нижнюю половину сферы и считать ее точки точками проективной плоскости. Если спроектировать ортогонально нижнюю полусферу на каса- касательную к ней в полюсе S плоскость а (рис. 81), то проек- проективная плоскость отображается на единичный круг, у ко- которого отождествлены диаметрально противоположные точ- точки границ. Каждая пара диаметрально противоположных точек границы соответствует при этом несобственной (бес- (бесконечно удаленной) прямой, присоединение которой к ев- евклидовой плоскости превращает ее в замкнутую поверх- поверхность — проективную плоскость. Поток идеального газа во вращающемся канале постоянного диаметра В некоторых системах турбовинтовых вертолетов и са- самолетов, а также воздушно-реактивных турбин газообраз- газообразная рабочая смесь пропускается через вращающиеся кана- каналы постоянной площади поперечного сечения, которые раз- размещаются в лопатках компрессора и которые в свою оче- очередь связаны посредством пустотелой вертикальной оси. Рис. 82 Чтобы установить оптимальные условия вращения, необ- необходимо проанализировать прохождение смеси через враща- вращающийся канал и связать решение с граничными условиями, обусловленными конструкцией канала. В лопатке газ при- принимает участие во вращательном движении относительно ее оси с постоянной угловой скоростью со и движется относи- относись . тельно канала с ускорением v-r, где v — скорость частицы 130
газа относительно канала, а г — координата, измеряемая вдоль вращающейся лопатки компрессора. На рис. 82 схематически изображен один вращающий- вращающийся канал лопатки компрессора*). Предполагается, чго рабочая смесь, начальное состояние которой известно, подается вдоль пустотелой оси к полости на оси, в которой скорость потока может считаться весьма незначительной. Граница нахождения газа в этой полости будет обозначать- обозначаться символом а; сам же газ предполагается идеальным с по- постоянной удельной теплоемкостью. Показатель изэнтро- пийного, т. е. обратимого адиабатического процесса будем обозначать через у. Далее предполагается, что газ расширяется через соп- сопло, заканчивающееся сечением Ь, которое является вход- входным сечением канала постоянного диаметра. Расширение газа из состояния а до состояния b предполагается изэнтро- пийным; скорость газа после расширения будем обозначать через v±. Расстояние от оси вращения О* до сечения Ъ обоз- обозначается через гх. При прохождении через канал, площадь постоянного поперечного сечения которого обозначается через А, а гидравлический диаметр — через D, газ ускоряется блгго- даря комбинированному действию градиента давления и динамического ускорения во вращающейся лопатке комп- компрессора. Влиянием изменения уровня давления, если оно вообще имеет место, а также изменением градиентов давле- давления, действующих на плоскость поперечного сечения, ко- которые являются следствием ускорений Кориолиса, будем пренебрегать. Это последнее предположение, вообще гово- говоря, требует экспериментальной проверки, так как сущест- существование поперечных градиентов давления может быть при- причиной вторичных потоков. Однако если диаметр канала мал по сравнению с его длиной, то такое предположение оправ- оправдано. Теперь ясно, что уравнения движения и энергии сжи- сжимаемой смеси, проходящей через канал постоянного дна- метра, находящийся в состоянии покоя, должны видоизме- видоизменяться с учетом силы инерции. Что же касается уравнения неразрывности, то оно остается неизменным. Далее предполагается, что, начиная от сечения с в конце канала, газ сжимается изэнтропийно, проходя через рас- *) Kestin J., Zaremba S. К- Adiabatic one-dimensional flow of a perfect, gas through a rotating tube of uniform cross-section// Aeronaut. Quart.—1954.—V. 4.- P. 373-399. 5» 131
ширяющееся сопло. При этом он приходит в состояние по- покоя относительно лопатки компрессора во второй полости на расстоянии г3 от оси вращения и достигает давления Pd и температуры 7Y Из второй полости газ расширяется изэнтропийно через сужающееся или сужающе-расширяющееся сопло так, что он покидает полость под прямым углом к оси канала. Та- Таким образом, возникает сила тяги, обуславливаемая на- наличием крутящего момента. Ниже для простоты изложения считаем, что выходное сопло является сужающимся и имеет выходное поперечное сечение площади Л*. Обозначая внешнее (атмосферное) давление через Ра, рассмотрим два случая прохождения смеси через сопло. Это, во-первых, случай, когда отно- отношение PJPa превосходит критическое значение, т. е. когда В этом случае поток на выходе из сопла будет дозвуковым и, следовательно, давление Р3 на выходе из сопла будег равно атмосферному давлению, т. е. Р — Р Во-вторых, это случай, когда отношение PaIP<i меньше кр1ь тического значения, т. е. когда PjPd < B/G+1 ))^>. Здесь давление Я3 на выходе из сопла имеет фиксированное значение, которое зависит от давления Pd и не зависит от давления Ри. Таким образом, В последнем случае поток на выходе из сопла будет иметь скорость звука где значение величины ad зависит только от температу- температуры Td. При дальнейшем анализе поток всюду считается адиа- адиабатическим, причем он рассматривается как изэнтропий- иый, за исключением части канала между сечениями Ьн с. Как и в случае вывода дифференциального уравнения, описывающего адиабатический поток идеального газа в соп- 132
ле переменного диаметра, остановимся теперь на уравне- уравнениях неразрывности, движения и энергии. Так, уравнение неразрывности в рассматриваемом случае имеет вид v v3 A* где V — удельный объем, г{; — плотность массы потока. Иначе, где т — масса потока. Чтобы вывести уравне- уравнение движения, обратимся к рис. 83. При этом заме- заметим, что динамическое действие вращательного движения лопатки компрес- компрессора может быть использовано для описания движения потока относительно движущегося канала, в котором в со- соответствии с принципом Даламбера предполагается, что сила инерции dl &rdr Рис. 83 действует в положительном направлении г. Следовательно, А , dv элемент массы dm = ydr движется с ускорением v -^ под действием силы инерции dl, силы давления A dP и силы л л л л трения dF = -^dr. Здесь т = ^ -, где X зависит от числа Рейнольдса R. В первом приближении можно считать, что % — величина постоянная вдоль всего канала. С учетом этого предположения уравнение движения записываемся в виде Тш udr~ нли, после упрощения, в виде 21 ., A61) Что же касается уравнения энергии, то оно может быть лег- легко выведено, если воспользоваться первым законом термо- термодинамики для незамкнутой системы и имгть п Риду, что количество работы, совершаемой системой, есть w2rdr. 133
Таким образом, dh + vdo—w*rdr = О, где h — энтальпия. Определяя скорость звука а с помощью равенства Л—^т. получаем, что a da + ^-77- и <&>—^- и*8/" dr = 0. Отсюда уравнение для энергии записывается в виде al = a* + '^vi---'2^-w*rh A62) Если теперь ввести безразмерные величины M0 = v/a0, x=^r/D, G* = w*D*jal, A63) то приходим к равенству йо , V- а* A64) Исключая давление и удельный объем из основных уравнений A60), A61) и A64), можно вывести уравнение, связывающее безразмерную величину Мо с безразмерным расстоянием, которое является основным дифференциаль- дифференциальным уравнением для решения рассматриваемой задачи. Из уравнения неразрывности A60) получаем, что A65) А тогда, учитывая равенства из уравнения A62) приходим к соотношениям J Y-I у 2у dr L y u \ ^y Yy^ ^y ь; Подставляя значения V из равенства A65) и -^т— из ра- равенства A66) в уравнение движения A61) и учитывая соот- 134
ношения A63), приходим к дифференциальному уравнению Полагая здесь = у(? —1). уравнение A67) можно переписать в виде dy _ 2y(my—G2x) e Дифференциальное уравнение A68) и является предме- предметом дальнейших исследований. Подстановкой A59) оно приводится к виду Ясно, что начало координат является особой точкой для уравнения A69). Члены наименьшей степени как в числите- числителе, так и в знаменателе являются квадратичными. Если отбросить члены высшего порядка, как это описано на с. 123 при рассмотрении точек равновесия высшего по- порядка, то получим, что 3r] + 2т]{if—I*) (т-ц-G*l), и, таким образом, начало координат является особой точ- точкой высшего порядка. Характеристическое уравнение дифференциального ураз- нения A69) после упрощений принимает вид lr\№ + rf){(q + 2)G*l-2mr]] = 0. A70) Отсюда имеем следующие три вещественные исключитель- исключительные прямы а а) 1 = 0, б) ч = 0, в) (<7-f2)G?£—2тт} = 0, A71) каждая из которых является регулярной и соответствует одному из множителей в уравнении A70). Легко видеть, что для положительных | и г\ значение Х4 (Н, !]) является положительным в окрестности всех трех 133
исключительных прямых, и поэтому в первом квадранте и в окрестности указанных прямых выражение У* (£. И) И ^ ft (Е' + Ч1) [(9 + 2) 6^-2/т]} имеет тот же знак, что и левая часть равенства A71 в), и, следовательно, отрицательный знак ниже прямой A71в) и положительный выше нее. Последнее объясняется тем, что выражение A72) может менять знак только при переходе Ч\. ч Рис. 84 Рис. 85 через исключительную прямую. Геометрически это озна- означает, что в первом квадранте на радиусах, расположенных между прямыми A716) и A71 в), правая часть дифферен- дифференциального уравнения dx\ __ Y4 (Е, ц) dl ~ X, (I, Ц) f которое определяет поведение интегральных кривых в кру- круге с центром в начале координат и с достаточно малым ра- радиусом, задает угол, который больше угла наклона радиу- радиусов. Что же касается радиусов, расположенных между пря- прямыми A71в) и A71а), то на них правая часть последнего диф- дифференциального уравнения задает угол, который меньше угла наклона радиусов (рис. .84). Таким образом, нормаль- нормальная область, содержащая прямую A71 в), должна быть от- отталкивающей, в то время как области, содержащие прямые A71а) и A716), являются притягивающими. В силу сим- симметрии поля, задаваемого вектором (Х4, К4), приведенные факты справедливы и для нормальных областей, получен- полученных из упомянутых выше поворотом вокруг особой точки на угол 180°. 136
Таким образом, существуют точно две интегральные кри- кривые, которые «входят» в особую точку по касательной A71 в), и бесконечно много таких кривых, касающихся ко- координатных осей A71а) и A716) в точке покоя. Таким образом, видим, что второй и четвертый квадран- квадранты содержат эллиптические секторы, так как они располо- расположены между двумя последовательными притягивающими нормальными областями (рис. 85). Каждый из первого и третьего квадрантов разделяется на два сектора интеграль- интегральными кривыми, которые касаются в начале координат ис- исключительной прямой A71 в). Эти секторы являются пара- параболическими, так как они рас- расположены между двумя после- yi довательными нормальными областями, одна из которых является притягивающей, а другая — отталкивающей. Бо- Более того, все интегральные кривые, за исключением тех, которые касаются прямой A71 в), касаются осей коорди- координат A71а) и A716) в особой точке. РиСт 86 Исходя из физических со- соображений, дифференциальное уравнение A68) имеет смысл рассматривать только в нер- нервом квадранте плоскости (х, у). Возвращаясь к этой плос- плоскости, видим, что существует точно одна интегральная кри- кривая, имеющая асимптотическое направление с угловым коэф- коэффициентом (g+2) G2/B га). Все же другие интегральные кри- кривые допускают в качестве асимптотического направления одну из осей координат. Действительно, нетрудно доказать, что все эти кривые приближаются асимптотически к одной из осей координат, т. е. что на каждой из этих кривых при стремлении к бесконечности не только у/х-^0 или х/у->0, но в действительности //->0 или х->0. Этот факт иллюстри- иллюстрируется на рис. 86, где начерчены некоторые из интеграль- интегральных кривых при больших значениях к и у. Обратим вни- внимание на то, что картина поведения интегральных кривых па конечном расстоянии от начала координат зависит преж- прежде всего от расположения особых точек и требует4 специаль- специального исследования. Можно также доказать, что интеграль- интегральные кривые, которые асимптотически приближаются к осям координат (рис. 86), покидают первый квадрант на конечном расстоянии от начала координат. 137
Изолированные замкнутые траектории Мы уже знаем, что в случае особой точки типа центра некоторая область фазовой плоскости сплошь заполнена замкнутыми траекториями. Вместе с тем возможна и более сложная ситуация, когда имеется изолированная замкну- замкнутая траектория, т. е. траектория, в некоторой окрестности которой нет других замкнутых траекторий. Последний случай непосредственно связан с решением вопроса о су- существовании изолированных периодических решений. При этом интересно, что изолированные замкнутые траектории могут иметь только нелинейные дифференциальные урав- уравнения и системы. Изолированные периодические решения соответствуют самым разнообразным свойствам явлений' и процессов, происходящих в биологии и радиофизике, в теории коле- колебаний и астрономии, в медицине и теории конструирова- конструирования приборов. Такие решения возникают при изучении диф- дифференциальных моделей в экономике, при рассмотрении различных вопросов автоматического регулирования, са- самолетостроения и т. д. Далее будет изучаться возмож- возможность возникновения изолированных периодических реше- решений при исследовании процессов, происходящих в электри- электрических цепях; здесь же мы рассмотрим как модель нелиней- нелинейную дифференциальную систему $L^-y + x(\-x*-tf), % = x + y(\-x*-if). A73) Чтобы решить ее, введем полярные координаты г, 8, где х=гсо$в, у=г sin 6. Тогда, продифференцировав соотноше- соотношения хг~\-у2—гг и 0=arctg (ylx) no t, получим равенства dx , dy dr dy dx odQ /лтлх +yt=rJt> xii~y-di-r"-di' <174> Умножая первое уравнение из системы A73) на х, а вто- второе — на у и складывая, с учетом первого равенства из A74) находим, что Если же умножить второе уравнение из системы A73) на х, а первое — на у и вычесть, с учетом второго равенства из A74) получим соотношение rzft=-r\ A76) 138
Система A73) имеет единственную особую точку О (О, 0). Так как мы интересуемся сейчас только построением траек- траекторий, то молаю считать С>0. А тогда уравнения A75) и A76) означают, что система A73) приводится к виду at = /-(! ■о. £- A77) Каждое из полученных уравнений системы A77) легко ин- интегрируется и все семейство решений, как нетрудно ви- видеть, задается формулами =-. е=/ + /0. A78) или, в старых переменных х и у, формулами cos (I -f- to) ,, sin X ——" * ?-2* V 1 + Се-х Если теперь в первом уравнении системы A78) положить С=0, то получим r=\, G—tf-J-/n. Эти равенства определяют замкнутую траекторию — окруж- окружность х2-{-у*= 1. Если С<С0, то ясно, что г>1 и г-*-\ при ^-^+°°. Если же О-О, то г<1 и снова г—^-1 при /->-|-оо. Это означает, что сущест- существует единственная замкнутая траек- траектория г=\, к которой все осталь- остальные траектории с течением време- времени приближаются по спиралям (рис. 87). Замкнутые фазовые траектории, обладающие таким свойством, назы- называют предельными циклами или, точнее (орбитально), устойчивы- устойчивыми предельными циклами. Дело в том, что различают еще два типа предельных циклов. Предельный цикл называется (орбитально) неустойчи- неустойчивым, если все соседние к нему траектории спиралевидно от него удаляются при £->+°°. Предельный цикл называется (орбитально) полуустойчивым, если все траектории с од- одной стороны (например, изнутри) наматываются на него, а с другой стороны (извне) разматываются с него при В рассмотренном выше примере мы смогли в явном ви- виде найти уравнение замкнутой фазовой траектории, однако в общем случае, конечно, этого сделать нельзя. Поэтому 139 Рис, 87
в теории обыкновенных дифференциальных уравнений иг- играют большую роль критерии, которые позволяют хотя бы выделить те области, где может содержаться предельный цикл. Заметим, что замкнутая траектория системы A22), если таковая существует, обязательно содержит внутри себя по крайней мере одну особую точку этой системы. Отсюда, в частности, следует, что если в некоторой области фазовой плоскости нет особых точек дифференциальной системы, то в этой обласги нет и замк- замкнутых траекторий. Пусть D — ограниченная область вместе со своей гра- границей, лежащая в фазовой плоскости и не содержащая осо- особых точек системы A22). Тогда доказывается (к р и т е- рий Пуанкаре — Бендиксона), что если Г— траектория системы A22), которая в начальный момент времени t~U выходит из точки, леоюаюей в области D, и остается в D при всех t^t0, то траектория Г либо сама яв- является замкнутой траекторией, либо с течением времени она по спирали приближается к замкнутой траектории. Проиллюстрируем этот факт с помощью рис. 88. Здесь D состоит из двух замкнутых кривых 1\ и Г2 и кольцевой области между ними. Свяжем с каждой граничной точкой области D вектор V(x, у) = Х{х, y)l + Y{x, Тогда если траектория Г, выходящая в начальный момент времени i=U из точки границы, входит в область D и оста- остается там при всех Г>£ь то» согласно сформулированному выше утверждению, траектория Г будет по спирали при- приближаться к некоторой замкнутой траектории Го, лежащей в области D. При эгом кривая Го должна окружать особую точку дифференциальной системы (точку Р), не лежащую в области D. Дифференциальная система A73) дает простой пример применения приведенного критерия для отыскания предель- предельных циклов. Действительно, система A73) имеет единствен- единственную особую точку 0@, 0), и поэтому в области D между окружностями радиусов г=1/2 и г=2 нет особых точек. Исходя из первого уравнения системы A77), видим, что -jt!>0 на внутренней окружности и ~л"<0 на внешней 140
окружности. Вектор V, связанный с точками на границе области D, всегда направлен внутрь D. Эго и означает, что в кольцевой области между окружностями радиусов г=1/2 и г=2 должна быть замкнутая траектория дифференциаль- дифференциальной системы A73). Такая замкнутая траектория действи- действительно существует, и ею является окружность радиуса г—1. Следует, однако, отметить, что в общем случае системы вида A22) практическая реализация критерия Пуанкаре — Бендиксона встречает большие трудности, ибо нет общих методов построения соответствующих областей, и поэтому успех зависит как от вида системы, так и от опыта исследо- исследователя. Вместе с тем необходимо иметь в виду, что отыска- отыскание признаков отсутствия предельных циклов не менее важ- важная задача, чем нахождение критериев их существования. Наиболее распространенным в этом плане является признак Дюлака: если существует 'непрерывная вместе с непрерывными частными производными функция В (х, у) такая, что в некоторой односвязной области D (ра- (разовой плоскости функция д(ВХ) d(DY) дх ду является функцией знакоопределенной *), то в области D нет предельных циклов дифференциальной системы A22). При В (х, г/)н==1 сформулированный признак называют при- признаком Бендиксона. Если обратиться к рассмотренному выше дифференци- дифференциальному уравнению A56), представляющему собой диффе- дифференциальную модель адиабатического одномерного потока идеального газа с постоянной удельной теплоемкостью в канале с трением, то для этого уравнения Х(х, y) = {\-y)F{xI Y(x, y) = 4 А тогда если положить В{х. y)= то оказывается, что в данном случае д(ХВ) , д{ХВ)_ дх ' ду ~- уд дх ' ду ~- F (х) ^ *) То есть определенно положительной или определенно отрица- отрицательной. 141
и, значит, уравнение A56) не имеет замкнутых интеграль- интегральных кривых. Остановимся еще на одном понятии, которое может быть использовано при выяснении вопроса о существовании предельных циклов. Это понятие индекса особой точки. Пусть Г — простая замкнутая кривая (т. е. кривая без самопересечений), которая не обязательно является тра- траекторией дифференциальной системы A22), лежит в фазо- фазовой плоскости и не проходит через особые точки этой сис- системы. Тогда если Р(х, у) — точка кривой Г, то вектор V(x9 y)~X{xy y)l+Y{x, y)J, где I, J — единичные векторы осей декартовых координат, будет ненулевым вектором и, следовательно, будет иметь определенное направление, заданное некоторым углом Э Рис. 89 Рис. 90 (рис. 89). Если точка Р{х, у), двигаясь по кривой Г, на- например, против хода часовой стрелки, совершит один обо- оборот, то вектор V сделает при этом некоторое целое число обо- оборотов, т. е. угол G получит приращение Д6=2зт, где п — целое положительное, нуль или целое отрицательное число. Это число и называют индексом замкнутой кривой Г (цик- (цикла Г). Если цикл Г начать непрерывно стягивать так, чтобы он не проходил при деформации через особые точки заданного векторного поля, то индекс цикла должен, с одной стороны, меняться непрерывно, а с другой — всегда быть равным целому числу. Это означает, что при непрерывной дефор- деформации индекс цикла не меняется. Исходя из этого свойства, под индексом особой точки понимают индекс простой замк- замкнутой кривой, окружающей эту особую точку. Отметим следующие свойства индекса: 1) индекс замкнутой траектории дифференциальной системы A22) равен +1; 142
2) индекс замкнутой кривой, окружающей несколько особых точек, равен сумме индексов этих точек; 3) индекс замкнутой кривой, не окружающей ни одной особой точки, равен нулю. Отсюда, в частности, следует, что так как индекс замк- замкнутой траектории системы A22) всегда равен +Ь то замкну- замкнутая траектория должна окружать либо одну особую точку с индексом +1, либо несколько таких точек, но обяза- обязательно с суммарным индексом, равным +1. Это часто ис- используется при доказательстве отсутствия предельных циклов. Индекс особой точки вычисляется по формуле A79) где е — число эллиптических секторов, a h — число ги- гиперболических секторов. Практически при вычислении ин- индекса особой точки можно пользоваться следующим прос- простым приемом. Пусть L — цикл, не проходящий через осо- особые точки дифференциальной системы A22) и гакой, что любая из траекторий дифференциальной системы имеет с кривой L не более конечного числа общих точек. При этом траектории могут пересекать кривую L или касаться ее. В случае касания учитываются (рис. 90) только внешние (типа А) и внутренние (типа JB) точки касания. Касания же типа С, т. е. когда точка С является точкой перегиба, во внимание не принимаются. Для вычисления индекса осо- особой точки используют формулу A79), где уже е — число точек внутреннего касания, a h — число точек внешнего касания траекторий системы A22) с циклом L. На рис. 91 показаны особые точки с индексами, равными соответственно 0, +2, +3, +1, —2. Выше уже отмечалось, что построение полной картины поведения фазовых траекторий дифференциальной системы A22) облегчается введением точки на бесконечности по- посредством преобразований A59). Топология дает здесь весь- весьма общую теорему, которая говорит о том, что когда не- непрерывное векторное поле с конечным числом особых точек задается на сфере, то сумма их индексов равна +2. Таким образом, если сумма индексов всех расположенных в ко- конечной части фазовой плоскости особых точек дифферен- дифференциальной системы, имеющей их конечное число, отлична от +2, то точка в бесконечности будет особой точкой с не- ненулевым индексом. 143
Если же вместо инверсии использовать однородные ко- координаты, то сумма индексов всех особых точек будет уже равна +1. То, что это так, можно заметить из того, что если плоскость проектируется на сферу с центром проекции в центре сферы, то две точки на сфере соответствуют одной точке проективной плоскости. При этом окружность боль- большого круга, параллельного плоскости, соответствует пря- прямой в бесконечности. * Ц/7 Рис. 91 Еслч обратиться к уравнению A68), описывающему ади- адиабатический одномерный поток идеального газа, проходя- проходящего через вращающийся канал постоянного диаметра, то, как видно из рис. 85, для особой точки в бесконечности е=-2, h~0. Отсюда следует, что индекс этой точки равен •'-2. При зтоад индекс не зависит от зна ;еиий констант в уравнении A08). Из предыдущего, в частности, следует, что сумма индексов конечных особых точек равна нулю. Можно показать, что r зависимости от того, ni^eei ли пря- прямая, заданная уравнением ту — G2x—О, две, одну двой- двойную или ни одной точки пересечения с параболой, заданной уравнением 1 — py-\-qG2x*~Q (что в свою очепедь экъи- галентно исполнению соотьетстеенно соотношений 144
G~G0, Cj<G0, где G0—2mql/s/p), имеюг место следующие комбинации конечных особых точек: а) G>G0. Седло и узел. б) G~G0. Особая точка является особой точкой высшего порядка с двумя гиперболическими (ft—2) и двумя парабо- параболическими {е—0) секторами. в) G<Xjq. Особые точки отсутствуют. Мы видим, чго во всех случаях сумма индексов равна нулю, как это и должно быть. Периодические режимы в электрлческих цепях Покажем, как появляются предельные циклы при рас- рассмотрении динатронного осциллятора (рис. 92), анализ работы которого приводит к так называемому уравнению Ваи-дср-Поля. Хотя явление, связанное с возникновением Рис. 92 предельных циклов, можно п|ч)иллюстрировать, напри- например, на задачах из механики, биологии, экономики, пока- покажем, как они появляются при изучении электрических цепей. На рис. 92 схематически изображен динатронный осцил- осциллятор, чьи характеристики ia, va на рис. 93 показаны сплош- сплошной линией. Здесь ia — ток, a va — напряжение в элек- электронной лампе. Цепь содержит сопротивление г, индук- индуктивность L и емкость С, соединенные параллельно и после- последовательно с динатроном. Реальная цепь может быть в данном случае заменена цепью, показанной схематически на рис. 94. При этом характеристика электронной лампы может быть аппроксимирована полиномом третьей степени 1=аН-уД что и показано пунктирной линией на рис. 93. Здесь i и v обозначают координаты в системе с началом, перенесенным в точку перегиба О. Как видно из рис. 93, 145
справедливы неравенства а > 0. у > 0. Далее, в соответствии с одним из законов Кирхгофа полу- получаем, что где ir — —, L -~- = v, ic — Cv. А отсюда в результате Рис. D4 простых операций приходим к дифференциальному урав- уравнению а , 1 . з7. Л • . 1 1 LC Если теперь положить а . 1 то предыдущее уравнение можно переписать в виде v + (а = 0. А это и есть так называемое дифференциальное уравнение Ван-дер-Поля. Если воспользоваться преобразованием вида то уравнению Ван-дср-Поля можно поставить в соответствие дифференциальное уравнение первого порядка A80) u, у) у V (v, у) Единственной конечной особой точкой здесь является начало координат, при этом j*@, O) = co^ D(v, //)== —(я 4- bv£). 146
Поскольку #>0, то, полагая а>0, приходим к выводу, что дивергенция D не меняет знака и, следовательно, замкнутые интегральные кривые у последнего дифференциального уравнения появиться не могут. Поэтому будем рассматри- рассматривать только случай, когда а<0, т. е. когда а<—Mr. Отсюда следует, что D @,0)=—а>0 и, значит, особая точка являет- является либо узлом, либо фокусом. Если теперь рассмотреть диф- дифференциальную систему, соответствующую дифференциаль- дифференциальному уравнению A80), т. е. рассмотреть дифференциальную систему 1-Й*. *). w=v<°-0. то видно, что при возрастании t изображающая точка дви- движься по траектории в направлении от особой точки. Та- Таким образом, траектория, которая выходит из точки в бес- бесконечности, не может достигнуть особой точки в начале координат ни. при каком значении t, включая /—+оо. От- Отсюда следует, что если доказать что траектория, которая выходит из бесконечно удаленной точки, имеет вид спира- спирали, навивающейся на начало координат при /->—f-oo, то тем самым можно гарантировать существование по край- крайней мере одного предельного цикла. Доказательство существования такого цикла может быть получено или численно, или аналитически. Численный метод доказательства, предложенный впервые голланд- голландским физиком и математиком Ван-дер-Полем, состоит в построении траектории, выходящей из точки, достаточно удаленной от начала координат, и из проверки того, обла- обладает ли она указанным выше свойством. Такая процедура дает приближение к предельному циклу, но она может ис- использоваться только в случае конкретных численных зна- значений. Ниже приведем схему доказательства *), основанную на аналитических соображениях, использующую исследо- исследование особенностей в бесконечности. При этом, в отличие от способа рассмотрения уравнения A68), воспользуемся более удобным в данном случае преобразованием однород- однородных координат v^t/z, tj^-ц/г. A81) Прямая в бесконечности здесь задается уравнением z=0. Для того чтобы числа переменных свести к двум, исключим *) Kestin J., Zaremba S. К. Geometrical methods in the analysis of ordinary differential equations//App!. Sci. Res., sect. В.— 1953. V. 3.—P. 144—189. 147
сначала точку £=z=0. Это можно сделать, полагая £=1. А тогда V— \/z, у^-ц/z. В этом случае дифференциальная система, соответствующая уравнению (ISO), преобразуется к виду dz (az1 \b) 11Ч- ©о г2 ■II3. dt га Здесь удобно ввести новый параметр 0, полагая А тогда предыдущая система переписывается так: A82) A83) Отметим прежде всего, что прямая в бесконечности z=0 является траекторией дифференциальной системы A83) и при возрастании t (а значит, и 0) изображающая точка движется по ней по направлению к единственной особой точке z=t) — 0. Характеристическое уравнение, которое в данном случае имеет вид Ш —Ьцг = 0, определяет регулярное исклю- исключительное направление 2=0, которому соответствуют две отталкивающие области, что может быть установлено рас- рассмотрением векторного поля (рис. 95). Второе исключитель- исключительное направление г)=0 явля- является особым, и поэтому здесь необходимы дополнительные рассуждения. Геометрическое место то- точек ?У=0 — это кривая, кото- которая касается прямой rj=O в начале координат и проходит через второй и третий квадранты (рис. 95), позволяя опре- определить три направления I, II, III по разные стороны от оси симметрии z=0. Область между кривой, заданной уравие- 148 Ш Рис. 95
нием 2/—0, и осью rj=U топологически эквивалентна двум отталкивающим областям. Таким образом, существует по крайней мере по одной траектории с каждой стороны от прямой i]=0, которые касаются последней в начале коор- координат. Если рассмотреть теперь дифференциальное уравнение — = - + -т+— + — то можно заметить, что во втором квадранте между кривой 2/=0 и осью т]=0 дц для малых значений |tj|. Следовательно, если взять две фа- фазовые траектории с одним и тем же значением tj, to нетруд- нетрудно видеть, что изображающие точки, движущиеся по этим траекториям, при убывании z будут расходиться. Это озна- означает, что существует только по одной фазовой траектории с каждой стороны от прямой rj=O, которые ее касаются в начале координат и принадлежат рассматриваемой облас- области. Ясно также, что качественная картина симметрична от- относительно оси z—0. Более того, поскольку для малых значений z и положительных т], то нет кривых, касающихся оси т]=0 в начале координат и проходящих через первый или четвертый квадрант. Такие кривые также отсутствуют и в области слева от кривой, заданной уравне- уравнением 2/=0, ибо в этом случае 6С>0, а % имеет тот же знак, что и г (рис. 95, стрелки IV). Проведенные рассуждения говорят о том, что особая точка т]— г=0— седло. Две фа- фазовые траектории, достигающие этой точки при /->-foo, являются сегментами прямой в бесконечности B=0), свя- связывающими предыдущую точку с точкой £>=z-0. Две же другие фазовые траектории достигают ссдловой точки при t—*—сю. В проведенных рассуждениях мы не касались точки £--z=0. Чтобы завершить исследование, положим ц=\ в формулах A81). Тогда дифференциальная система, соот- соответствующая дифференциальному уравнению A80), запи- запишется в виде -149
где переменная 9 задается формулой A82). Рассматривае- Рассматриваемая точка £,=z=0 оказывается особой. Характеристичес- Характеристическое уравнение здесь z3=0. Таким образом, исключительное направление z—О оказывается особым. Кривая, заданная уравнением Р(£, z)=0 (рис. 96), касается оси 2=0 в начале координат и имеет в нем точку возврата. Кривая же, за- заданная уравнением Q(%, z)=0, имеет в начале координат кратную точку. Исследуя знаки функций Р и Q, можно Рис. 96 Рис. 97 выяснить характер векторного поля (рис. 96), а также кар- картину поведения фазовых траекторий в окрестности особой точки £=z=0 (рис. 97). В частности, отметим, что вдали от оси |=0 пет фазовых траекторий, которые бы входили в особую точку справа. Это следует из того, что в первом и четвертом квадрантах выполняется неравенство Проведенные рассуждения показывают, что уравнение Ван-дер-Поля не имеет фазовых траекторий, которые стре- стремятся к бесконечности с ростом t, но существует бесконеч- бесконечно много фазовых кривых, которые удаляются от бесконеч- бесконечности при возрастании L Это и доказывает существование по крайней мере одного предельного цикла у дифференци- дифференциального уравнения Ван-дер-Поля. Кривые без контакта В сравнительно простых случаях полная картина пове- поведения интегральных кривых заданного дифференциаль- дифференциального уравнения или, что то же самое, фазовых траекторий соответствующей ему дифференциальной системы опреде- определяется типом особых точек и замкнутыми интегральными 150
кривыми (фазовыми траекториями), если таковые имеются. Иногда качественная картина может быть построена, если удается, кроме выяснения типов особых точек, найти кри- кривые (сепаратрисы), которые «связывают» особые точки. Однако, к сожалению, нет общих методов, которые позво- позволяли бы решать эффективно последнюю задачу. В связи с этим полезно при качественном интегрировании использо- использовать так называемые кривые без контакта. Напомним, что кривая или дуга кривой с непрерывной касательной назы- называется кривой {дугой) без контакта, если она нигде не ка- касается вектора (X, Y), заданного дифференциальной систе- системой A22). Из определения ясно, что вектор (X, Y) должен быть направлен на кривой всегда в одну сторону. Таким образом, кривая без контакта может пересекаться фазовыми кривыми системы A22) в одном направлении только при воз- возрастании t, а в другом — при убывании /. Поэтому знание соответствующих кривых без контакта может дать требуе- требуемую информацию о ходе выбранной частной фазовой траек- траектории. При качественном интегрировании дифференциальных уравнений могут применяться также различные вспомога- вспомогательные неравенства. Так, если заданы два дифференциаль- дифференциальных уравнения f-/(*.»). £=£'(*•</) и известно, что в некоторой области D выполняется нера- неравенство / {х, y)^g(x, у), то, обозначая через уг(х) решение первого уравнения такое, что yi(xo)=yo, где (xo,yo)^D, a через у»(х) — решение второго уравнения с теми же на- начальными условиями, можно доказать, что yi{x)^y.2(x) для х^х0 в области D. Если же в области D выполняется строгое неравенство f (x, y)<g(x, у), то yi(x)<:yz(x) для jv>a'o в области D и кривая у=Уг{х) является кривой без контакта. Для примера рассмотрим дифференциальное уравнение A68). Выше было показано, что если (j>G0, to это уравне- уравнение допускает две конечные особые точки — седло и узел, которые являются точками пересечения прямой ту—G'2x=0 и параболы 1 — py-j-qG2 х2=0. Сегменты прямой и пара- параболы, связывающие эти две конечные особые точки, явля- являются кривыми без контакта, и они ограничивают область плоскости, которую обозначим через А. Если возьмем . X {х, у) = 1 —ру + qG*x\ У (*, у) = 2у \tny—G*x)t 151
то легко видеть, что вектор (X, Y) направлен па границе области А наружу, исключая особые точки. Следовательно, если изображающая точка выходит из любой внутренней точки области А и следует вдоль интегральной кривой, ког- когда t убывает, то изображающая точка не сможет покинуть область А, не проходя через одну из особых точек. Но так как внутри области А выполняется неравенство Х(х, #)<0, то особая точка, которая является притягивающей, будет узлом. Найдя угловые коэффициенты исключительных направ- направлений для седла, можно увидеть, что одна из исключитель- исключительных прямых проходит через область А. Отсюда следует, что интегральная кривая, касательная к этой прямой в сед- ловой точке, входит внутрь области А и отсюда должна до- достигнуть узловой точки. В заключение отметим, что при исследовании конкрет- конкретных дифференциальных моделей часто возникает необходи- необходимость применять методы, о которых не шла речь в настоя- настоящей книге. Все зависит от степени сложности дифференци- дифференциальной модели, от того, насколько глубоко развит соответст- соответствующий математический аппарат, и, конечно, от эрудиции и опыта исследователя.
СПИСОК ЛИТЕРАТУРЫ 1. Амелькин В. В., Садовский Л. П. Математические модели и дифференциальные уравнения.— Минск: Вышэйшая шко- школа, 1982.— 272 с. 2. Андроноъ Л. А., ВиттЛ. А., X а й к и н С. Э. Теория ко- колебаний.— М.: Физматгиз, 1959.— 916 с. 3. D e r r i с k W. R., Grossman S. 1. Elementary differential equations with applications.— 2-nd ed.— Reading. Mass.: Addison- Wesley, 1981,—532 p. 4. Differential equation models / Ed.: Braun M.— New York etc.: Sprin- Springer, 1983.—380 p. 5. E p у г и ii H. П. Книга для чтения по общему курсу дифференци- дифференциальных уравнений.— Минск: Наука и техника, 1979.—744 с. 6. Пономарев К- К- Составление дифференциальных уравне- уравнений.— Минск: Вышэйшая школа, 1973.— 560 с. 7. S i m m о n s Q. F. Differential equations with applications and historical notes.— New York, N. Y.: McGraw-Hill Book Co., 1972.— 465 p. 8. S p i e g e 1 M. R. Applied differential equations.— Englewcod Cliffs, N. J.: Prcnticc-Ilall, 1981.—654 p.
Владимир Васильевич Амелькин ДИФФЕРЕНЦИАЛЬНЫЕ УРАВНЕНИЯ В ПРИЛОЖЕНИЯХ Рсдакюр М. М. Горячая Художественный редактор Т. И. Кольченко Технический редактор С. fl. Шкляр Корректор Н. Б. Румянцева ИБ .№ 12918 Сдано в набор 07.05.86. Подписано к печати 11.12.86. Формат 84X108/32. Бумага тип. ЛйЗ. Гарнитура литературная. Печать высокая. Усл. неч. л. 8,4. Усл. кр.-отт. 8,61. Уч.-изд. л. 7,65, Тираж 100 000 экз. Заказ № 2528. Цена 20 коп. Ордена Трудового Красного Знамени издательство «Наука» Главная редакция физико-математической литературы 117071 Москва В-71, Ленинский проспект, \Ъ Ордена Октябрьской Революции и ордена Трудового Красного Знамени МПО «Первая Образцовая типография» имени А. А. Жданова Союзполиграфпрома при Государственном комитете СССР по делам издательств, полиграфин и книжной торговли. 113054 Москва М-54, Валовая, 28 Отпечатано во 2-й типографии ИЗд-ва «Наука». 121099 Москва i -09s Шубинский пер., С.дак 54
ПРОИЗВОДНЫЕ ЭЛЕМЕНТАРНЫХ ФУНКЦИЙ Функция С (шллоянная) X хп 1 X 1 X" V* nf— у х ех ах In* loge х \gx sin л COS X tgx Производная 0 1 пхп~* 1 n Xn+X 1 2 J/T 1 «/ ex a^ In a 1 I , 1 x Oga° x\na 1 , 0,4343 x ь л: COS Jt — sin x: cos2 л; .—5—= — esc2 x 155
Продолжение Функция Произподная SCX esc X arcsin х arccos x arctgx arcctgx arcsc x arccsc x shx ch x tlu Arch .v Arch x Arth x Aiclhx sm x COS2 X COS-Г Sib8* -=t£ XSCX = — dgx esc x V'\~- ch a' sh x \ Y'x* - i 1 1 - x* I
ОСНОВНЫЕ ИНТЕГРАЛЫ СТЕПЕННЫЕ ФУНКЦИИ X1' п+1 — = lnUl X ТРИГОНОМЕТРИЧЕСКИЕ ФУНКЦИИ \ sin х dx = — cos x \ cos x dx — sin л; V tg xdx — — Iti|cosx| \ ctgxdx= In | I ь sin2 x :— ctg.r ДРОБ НО-Р АЦИОН А Л ЬН Ы Е ФУ Н КЦИИ С dx 1 л: \ T~i—.,-=— arctg — Г-^-5- = — Arth — = J-ln^iii <| л: f < a) J a2 — *3 a a 2a a~x v 1 dx I , х \ . х - а =— Arcth 1 x2-— a2 a a 2a i 157
ПОКАЗАТЕЛЬНЫЕ ФУНКЦИИ С ех dx=ex 1 In a ГИПЕРБОЛИЧЕСКИЕ ФУНКЦИИ \ sh х dx = ch я \ ch xdx — shx f th xdx = \n\ сЬл;[ \ cth xdx = \n\sh x[ С dx J chax \ г"я —- Ct П ИРРАЦИОНАЛЬНЫЕ ФУНКЦИИ К а%-\- I dX =Arrh - = ln |(ж4-}A 22 «
ИЗДАТЕЛЬСТВО «НАУКА» ГЛАВНАЯ РЕДАКЦИЯ ФИЗИКО-МАТЕМАТИЧЕСКОЙ ЛИТЕРАТУРЫ 117071 Москва В-71, Ленинский проспект, 15 ГОТОВИТСЯ К ПЕЧАТИ: Фаддеев Д. К., Никулин М. С, Соколов- Соколовский И. Ф. Элементы высшей математики для школьников.— 18 л. (Темплан 1987 г., № 70) В книге излагаются основные понятия дифференциального п интегрального исчислений, их приложения к исследованию эле- элементарных функций, применения к приближенным вычислениям, решению некоторых задач механики и физики. Имеются главы, посвященные изучению тригонометрических функций, комплекс- комплексных чисел, элементов теории вероятностей. Каждая глава снаб- снабжена упражнениями. Для учащихся старших классов школ и ПТУ, студентш техникумов и вузов, а также преподавателей математики, ин- инженеров и техников» Заказы на книгу принимаются без ограничения всеми книжными магазинами Книготорга и Академкниги.
ИЗДАТЕЛЬСТВО «НАУКА» ГЛАВНАЯ РЕДАКЦИЯ ФИЗИКО-МАТЕМАТИЧЕСКОЙ ЛИТЕРАТУРЫ 117071 Москва В-71, Ленинский проспект, 15 ГОТОВИТСЯ К ПЕЧАТИ: Горелов И. Н. Разговор с компьютером. Пси- холиигвистические аспекты проблемы.— 12 л. (Темплан 1987 г., № 157) Ученые и инженеры разных стран заняты разработкой ЭВМ пятого поколения. Наступит черед и следующих поколений компьютеров. Эти машины будут способны общаться с челове- человеком не только на своем, машинном, но и на естественном язы- языке, уметь ориентироваться в среде и ситуации общения. Чтобы их создать, надо знать особенности мышления человека и сся^и его мышления с речью, знать, как понимается текст. Обо всем этом, а также о многом другом (например, о рассудочной дея- деятельности животных, об опытах обучения обезьян языкам же- жесте ii и знаков) рассказывается в книге, написанной строго научно и вместе с тем доступно и увлекательно. Для всех, кто интересуется современной наукой о языке и мышлении, проблемами искусственного интеллект* и и»-1 ^- матики. Заказы на книгу принимаются оез ограничь»,. isvx,«n ма- магазинами Книготорга и Лкя1<'*'кпиги.