CC BY
31
Изв. вузов «ПНД», т. 20, № 6, 2012 УДК 519.6
ДВУМЕРНЫЕ САМООРГАНИЗОВАННО КРИТИЧЕСКИЕ МОДЕЛИ ТИПА КУЧИ ПЕСКА С АНИЗОТРОПНОЙ ДИНАМИКОЙ РАСПРОСТРАНЕНИЯ АКТИВНОСТИ
А. В. Подлазов
Работа посвящена численному и аналитическому исследованию двух самоорганизо-ванно критических моделей типа кучи песка, имеющих анизотропную динамику распространения активности, - модели Дхара-Рамасвами и дискретной модели Федеров. Теоретически определен полный набор критических показателей для этих моделей.
Дается систематическое изложение метода конечно-размерного скейлинга и его применения для решения самоорганизованно критических систем.
При изучении дискретной модели Федеров обнаружен и объяснен ряд нетривиальных явлений - таких как спонтанная анизотропия, аномальная диффузия и возникновение срединного рва заполнения.
Ключевые слова: Самоорганизованная критичность, куча песка, масштабная инвариантность, степенные распределения, конечно-размерный скейлинг, аномальная диффузия, спонтанная анизотропия.
Введение
Четвертьвековое развитие теории самоорганизованной критичности [1, 2], объясняющей возникновение степенных распределений вероятности в нелинейных далеких от равновесия системах, привело к появлению ряда базовых моделей, которые при очень простых правилах обладают нетривиальной масштабно-инвариантной динамикой. Одной из основных задач при исследовании таких моделей является аналитическое и компьютерное определение показателей наблюдаемых распределений.
Обыкновенно данная задача решается численно, тогда как теоретически показатели определяются лишь частично, то есть удается выразить все их через один-два ключевых показателя. Последние всё равно приходится либо брать из эксперимента, либо определять с помощью приближенных методов (см., например, [3]). Полное аналитическое решение самоорганизованно критических моделей легко находится лишь тогда, когда их динамика может быть сведена к ветвящемуся процессу с независимыми частицами (например, в случаях анизотропных правил [4], случайного
соседства [5] или бесконечномерного пространства [6]). А в менее тривиальных ситуациях поиск решения оказывается очень сложной задачей, на протяжении продолжительного времени не поддающейся усилиям исследователей.
Следует особо подчеркнуть, что методы ренормгруппы, хорошо зарекомендовавшие себя при изучении обычных критических явлений [7], не всегда оказываются применимы в случае самоорганизованной критичности. Характерным примером могут служить двумерные модели Бака-Танга-Визенфельда [1, 2] и Манна [8, 9]. Данные модели, правила которых обладают одинаковыми симметриями и законами сохранения, должны, с ренормгрупповой точки зрения, характеризоваться одними и теми же показателями. Однако в действительности это не так, причина чего заключается в различии симметрий, возникающем на макроскопическом уровне при развитии лавин [10, 11]. Таким образом, в общем случае не приходится рассчитывать на решение моделей без учета динамики, обеспечивающей их самоорганизацию в критическое состояние и функционирование в нем как целостных систем.
В настоящей статье предлагается принципиально новый метод, позволяющий получить полное точное решение двух классических моделей типа кучи песка, существенно различных по правилам.
Статья имеет следующую структуру. Часть 1 носит обзорный характер - в ней описываются математические свойства масштабно-инвариантного состояния и на примере концептуальной модели кучи песка поясняется его самоорганизованно-критическая природа. Часть 2 демонстрирует применение этих общих идей для решения модели Дхара-Рамасвами - одной из простейших моделей теории самоорганизованной критичности. Часть 3 посвящена детальному изучению свойств и решению дискретной модели Федеров, которая имеет значительно более сложное поведение, связанное с аномальной диффузией. Это явление получает качественное объяснение и количественное описание в части 4. Завершается статья выводами, в которых суммируются и обобщаются ее результаты.
1. Масштабная инвариантность и самоорганизованная критичность
Скейлинговые свойства масштабно-инвариантных систем. Одним из ключевых признаков сложного поведения является масштабная инвариантность происходящих процессов, то есть отсутствие у описывающих их переменных собственных характерных значений, обуславливаемое степенным видом плотности распределения
Эта формула, разумеется, справедлива лишь в промежуточной асимптотике, протяженность которой ограничена с обеих сторон: сверху - в силу конечности размеров рассматриваемой системы, а снизу - в силу конечности размеров слагающих ее элементов. И если отклонение и (х) от степенного вида при малых х обычно не представляет особого интереса, то при больших оно становится принципиальным. Нестепенное поведение плотности учитывается с помощью метода конечно-размерного скейлинга [12], предполагающего замену формулы (1) более общей записью вида
и (х) - х-(1+а
(1)
и (х) = х (1+а)Л, (х/х1) ,
(2)
где х\ - правая граница области промежуточной асимптотики, а Н (Е) - скейлинговая функция, примерно постоянная при Е ^ 1 и убывающая быстрее любой степени аргумента при Е ^
В силу масштабной инвариантности, величина события, крупного настолько, что оно уже не помещается в систему размера Ь, растет как некоторая его степень:
Х1 ~ V. (3)
Сочетание формул (2) и (3) позволяет определить скейлинговое поведение и для среднего события
(х) = ! хи (х) йх ~ Ьа, где о = V (1 — а)
(4)
Скейлинговые показатели о и V, связывающие характерные значения (х) и х1 с обусловившим их появление конечным размером системы Ь, часто удается установить из общих соображений, что существенно облегчает теоретический анализ. А для практического анализа удобно представить плотность в виде
и(х) = Ь-вд(хЬ-^). (5)
Эквивалентность записей (2) и (5) предполагает, что д (Е) ~ Е-(1+а) при Е ^ 1, а кроме того, выполнено скейлинговое соотношение
в = V (1 + а). (6)
Графики зависимостей Ьви (х) от Ь-Чх, полученные в результате компьютерного моделирования при различных размерах Ь, можно совместить, подобрав правильные значения скейлинговых показателей в и V. Это дает возможность найти характеристический показатель а, не зная точный вид функций Н и д (примеры см. далее на рис. 3 и 17).
Наряду с рассмотрением отдельных переменных интерес представляет анализ взаимосвязи между разными переменными, описывающими одно и то же событие. Пусть имеются две величины х и у с характеристическими показателями ах и ау, соответственно. Зависимость между этими величинами будем считать функциональной, понимая ее как условное среднее
У (х) = (У)х ~ х"!ух ■
Записав формулу преобразования вероятностей их (х) йх = иу (у) йу, выражаем показатель связи через характеристические показатели
У ух = ах/ау. (7)
С другой стороны, сопоставив размеры крупных событий у1 ~ х\ух, получаем соотношение и на скейлинговые показатели
Уух = Vy/Vx■ (8)
В комбинации формулы (7) и (8) дают инвариант
е = утат = V,, а.
у^у,
не зависящий от рассматриваемой переменной. Универсальный показатель е определяет скейлинговое поведение доли сверхкрупных событий
тете те
Ргоь {,>„(=/„ (,) * ^ ь = (!)« -
Если событие не помещается в систему данного размера, то оно не помещается сразу по всем описывающим его переменным. Поэтому не существенно то, какая из них используется для расчета показателя е.
«Песочная» парадигма. Механизм возникновения масштабно-инвариантных свойств у динамических систем дает теория самоорганизованной критичности [1, 2, 13]. Ее базовой моделью является куча песка.
Рассмотрим уголок с песком, изображенный на рис. 1. Будем считать, что возможно лишь поверхностное перемещение песка, причем инерцией его движения можно пренебречь. Тогда состояние системы вполне определяется средним наклоном поверхности г. Если он невелик, то песок в целом неподвижен. А при превышении наклоном некоторого порогового значения гс возникает спонтанный ток песка 3 по поверхности, непрерывно возрастающий по мере увеличения г (см. врезку к рис. 1). То есть налицо непрерывный фазовый переход, в котором управляющим параметром является наклон г, а параметром порядка - ток песка 3.
Критическое состояние соответствует моменту отрыва параметра порядка от нуля. Обычные критические системы попадают в это состояние благодаря тонкой подстройке управляющего параметра к заранее не известному критическому
значению. Однако в ряде систем оказывается возможным, устанавливая параметр порядка в +0, вынудить управляющий параметр самостоятельно отыскать критическую точку, что и называется самоорганизованной критичностью.
Иначе говоря, вместо того, чтобы крутить ручку прибора, можно начать сдвигать с нулевой отметки стрелку на его шкале, вынуждая ручку повернуться до нужного положения [14]. Такое управление параметром порядка обычно достигается при помощи разделения временных масштабов [15], при котором время релаксации системы много меньше времени между последовательными возмущениями, то есть события едва происходят.
Самоорганизация кучи песка в критическое состояние происходит при
Рис. 1. Уголок с песком. Состояние песка определяется углом наклона поверхности г. При его изменении происходит непрерывный фазовый переход (зависимость параметра порядка от управляющего параметра приведена на врезке) от неподвижного состояния (7 = 0) к состоянию непрерывного тока песка (7 > 0). При токе 7 = +0, соответствующем добавлению одной песчинки за один шаг, система самоорганизуется в состояние с критическим наклоном г = гс
токе ■] = +0. Чтобы обеспечить такую величину параметра порядка, будем добавлять песчинки по одной на вершину кучи (см. рис. 1), дожидаясь завершения процесса релаксации. При этом ток песка через систему, очевидно, имеет минимально возможное значение - в среднем одна песчинка за время одного события.
Если наклон поверхности мал, то осыпание, вызванное добавленной песчинкой, скорее всего, не достигнет края кучи и наклон увеличится. При большом наклоне возможно возникновение глобального осыпания, в результате которого большое количества песка покинет систему и наклон уменьшится. Равновесие между количеством песка, добавляемого в систему, и количеством песка, покидающего ее, достигается при критическом наклоне поверхности, когда возмущение может распространяться по куче сколь угодно далеко, не затухая и не разрастаясь. При этом описывающие его переменные не имеют собственных характерных размеров, подчиняясь степенным распределениям вероятностей.
Таким образом, имеет место отрицательная обратная связь, вынуждающая наклон принять со временем значение г = гс вне зависимости от начального профиля поверхности. При этом куча песка, состоящая из локально взаимодействующих песчинок, начинает вести себя как единое целое. То есть в результате самоорганизации в критическое состояние система приобретает свойства, которых не было у ее элементов, демонстрируя сложное целостное поведение.
Простейшей моделью, описывающей рассмотренную выше динамику кучи песка, является DR-модель, предложенная Д. Дхаром и Р. Рамасвами в 1989 году [4].
Правила модели. Правила модели формулируются как клеточный автомат на двумерной гексагональной решетке размерами Ь х Ь (рис. 2). В ячейках расположены целые числа , характеризующие локальный наклон поверхности кучи. Горизонтальные слои решетки (с заданным г) условно соответствуют линиям уровня поверхности. Если превышает пороговое значение гс = 1, ячейка объявляется неустойчивой и опрокидывается. Опрокидывание ячейки заключается в уменьшении на 2 стоящего в ней числа с одновременным увеличением на 1 чисел в двух ячейках, примыкающих к данной снизу (см. рис. 2):
что символизирует пересыпание песчинок вниз по склону. При наличии нескольких неустойчивых ячеек они опрокидываются одновременно.
Элементарное событие состоит из возмущения и релаксации. Возмущение устойчивого состояния производится путем увеличения на единицу значения в случайно выбранной ячейке верхнего слоя, что соответствует добавлению одной песчинки на вершину кучи. Если в результате возмущения ячейка теряет устойчивость, то она опрокидывается и начинается процесс релаксации. Опрокидывание ячейки приводит к увеличению наклона в нижележащих ячейках, что, в свою очередь, способно нарушить их устойчивость и т.д. по принципу цепной реакции. Таким образом, потеря устойчивости одной ячейкой может вызвать лавину опрокидываний, продолжа-
2. DR-модель
Рис. 2. БЯ-модель: правила и пример лавины. Устойчивыми считаются ячейки с нулевым или единичным наклоном. При потере ячейкой устойчивости из нее изымаются две песчинки и передаются в пару нижележащих ячеек. Лавина, инициируемая добавлением одной песчинки в случайно выбранную ячейку верхнего слоя, распространяется строго вниз. Слева приведено состояние системы до лавины опрокидываний, справа - после. Темной заливкой показана область лавины, светлой - ячейки на ее границе, которые, получив песчинку, сохранили устойчивость
ющуюся до тех пор, пока все ячейки вновь не обретут устойчивость. После этого релаксационный процесс считается завершенным и дается старт следующему событию.
Нижний край решетки является открытым, так что при опрокидывании ячейки из нижнего слоя две песчинки покидают систему. Это обеспечивает существование стационарного состояния и возможность самоорганизации. Ради простоты левый и правый края решетки отождествляются, то есть она свернута в вертикальный цилиндр (периодические граничные условия).
В БЯ-модели активность не возвращается, то есть лавины распространяются строго сверху вниз, не затрагивая дважды один слой. Поэтому их вполне можно охарактеризовать всего двумя числами: длительностью Т, определяемой как число шагов моделирования, в ходе которых происходили опрокидывания, и площадью Б, определяемой как число опрокинувшихся ячеек. Длительность лавины для данной модели совпадает с числом затронутых ею слоев решетки, а площадь - с размером лавины N, определяемым как число опрокидываний.
Как показывает компьютерное моделирование, обработка результатов которого приведена на рис. 3, распределение лавин в БЯ-модели по длительности характеризуется скейлинговыми показателями Vт = 1 и |3у = 3/2, а по площади - Vs = 3/2 и = 2. Характеристические показатели, в соответствии с формулой (6), равны ат = 1/2 и аs = 1/3 (см. врезки на рис. 3), а универсальный показатель, даваемый формулой (9), е = 1/2.
Эти значения были теоретически рассчитаны еще авторами модели [4]. Однако их решение опирается на ряд нетривиальных математических результатов, которые,
Ю-3 10-2 ю-1 Ь'1-Т 10"5 10"4 10"3 10"2 10"1 1"3/2.5
Рис. 3. Распределение лавин в DR-модели по длительности Т и площади Я. При правильном подборе скейлинговых показателей графики плотности распределения, полученные для систем различного размера, совмещаются. Правильность определения соответствующих характеристических показателей подтверждается врезками, на которых скомпенсирована степенная компонента плотности, в результате чего график в промежуточной асимптотике становится горизонтальным. Пики (для длительности) и горбы (для площади) в правой части графиков соответствуют лавинам, которые не поместились в решетку, оборвавшись из-за достижения ее нижнего края
в действительности, являются излишними и только маскируют суть явления. Кроме того, решение Д. Дхара и Р. Рамасвами не может быть распространено на рассматриваемую далее DFF-модель, родственную по свойствам. Поэтому мы здесь предлагаем альтернативное решение, дающее общее представление о том, как могут быть найдены показатели для самоорганизованно критических систем.
Решение модели. Длительность лавин, крупных настолько, что они уже не помещаются в систему, определяется, очевидно, числом имеющихся в решетке слоев
Т = Ь, или ут = 1. (10)
На каждом шаге в систему добавляется одна песчинка, поэтому в стационарном режиме в среднем одна песчинка будет ее и покидать. Единичное опрокидывание перемещает две песчинки на один слой, а для прохождения песчинки через всю систему необходимо перемещение размером в Ь песчинко-слоев. Таким образом, средняя площадь области лавины
(5) = Ь/2, или о^ = 1. (11)
Чтобы завершить решение, необходимо установить связь между площадью и длительностью лавин. Для этого вспомним, что мы имеем дело с самоорганизованно критической системой. На каждом слое I, который проходит лавина, она имеет некоторую ширину ш (I), определяемую как число опрокинувшихся ячеек данного слоя. Среднее изменение ширины лавины от слоя к слою
(Аш) = 0. (12)
В самом деле, будь оно отрицательно, вероятность достижения лавиной слоя I экспоненциально убывала бы с его номером, а будь оно положительно, лавина распространялась бы неограниченно с ненулевой вероятностью. Поскольку в первом случае
Рис. 4. Зависимость площади лавины от ее длительности для БЯ-модели. Зависимость понимается в смысле условного среднего, то есть типичной площади лавины при заданной ее длительности. Для наглядности обе переменные отмасштабированы на размер системы со скейлинговыми показателями Vт = 1 и Vs = 3/2, в результате чего графики, полученные при разных Ь, совместились
количество песка в системе увеличивается, а во втором - уменьшается, возникает отрицательная обратная связь, подстраивающая систему в критическое состояние. В нем при прохождении лавины по слоям ее ширина изменяется лишь диффузионным образом -как координата частицы, совершающей несмещенные случайные блуждания.
Следовательно, типичная ширина лавины, достигшей слоя I, определяется формулой Эйнштейна-Смолуховс-кого (см. также формулу (30) далее)
w (I) -VI. (13)
Суммирование ширин по слоям дает типичную площадь лавины, достигшей слоя Т,
т
Б (Т) - I w (I) (I - Т3/2, о
что подтверждается рис. 4.
Таким образом, показатель связи
Ysт = 3/2, (14)
что позволяет легко найти приведенные выше значения показателей как решение системы уравнений {(4), (7), (9), (10), (11), (14)}.
3. DFF-модель
Сравнительная простота теоретического исследования БЯ-модели связана с тем, что в ней лавина не затрагивает повторно одну и ту же ячейку в силу существенной анизотропии правил. Оказывается, существует самоорганизованно критическая модель с очень похожими свойствами, но при этом с изотропными правилами. Это дискретный вариант модели, предложенной Х. и Е. Федерами в 1991 году [16], -DFF-модель.
Правила модели. ББР-модель представляет собой клеточный автомат с двумерной ортогональной решеткой размерами Ь х Ь. Находящиеся в ячейках целые числа интерпретируются как количество песчинок, способных участвовать в процессах пересыпания. Выделенного направления склона нет. Если ^ превышает пороговое значение хс = 3, ячейка объявляется неустойчивой и опрокидывается.
Опрокидывание заключается в обнулении стоящего в ячейке числа с одновременным увеличением на 1 значений в четырех ячейках, имеющих с данной общую сторону (см. рис. 5, в центре):
г^,] ^ 0,
г±\]±\ ^ + 1.
Если есть несколько неустойчивых ячеек, они опрокидываются одновременно.
Принципиально, что правила опрокидывания неконсервативны, то есть не предполагают сохранения количества песчинок при их раздаче соседям. Если теряет устойчивость ячейка, содержащая более 4 песчинок, то их излишек необратимо теряется, то есть происходит диссипация. Пример развития лавины, в ходе которой происходит событие диссипации, приведен на рис. 5.
Возмущение устойчивого состояния производится путем добавления одной песчинки в случайно выбранную ячейку. Важна равномерность возмущения по пространству, то есть все ячейки решетки должны выбираться для вброса песчинок с равной вероятностью.
Может показаться, что в силу наличия в системе диссипации для достижения стационарного состояния открытый край не нужен. Однако, если отождествить противоположные края решетки, свернув ее в тор, то система в какой-то момент впадает в бесконечный цикл. То есть количество песка, достаточное для того, чтобы сколь угодно долго поддерживать развитие лавины, оказывается недостаточным для возникновения диссипации.
Примечательно, что во время бесконечной лавины среднее заполнение ячеек (г) = х^ = 2 точно (без дисперсии). Объясняется данное значение устройством бесконечной лавины, которая представляет собой один или несколько фронтов опрокидывания, циклически обходящих решетку. При этом каждая ячейка, потерявшая устойчивость, должна при опрокидывании в среднем вернуть назад (туда, откуда фронт уже ушел) столько же песчинок, сколько она передает вперед (туда, куда идет фронт). А вперед передается столько песчинок, сколько в среднем не хватает ячейке до потери устойчивости. В сочетании с отсутствием диссипации два указанных обстоятельства и приводят к значению (г) = (гс + 1)/2 = 2.
Рис. 5. DFF-модель: правила и пример развития лавины. Устойчивыми считаются ячейки со значениями 0, 1, 2 или 3. При опрокидывании неустойчивой ячейки из нее изымаются четыре песчинки и передаются соседям, а ее значение обнуляется. Лавина инициируется добавлением одной песчинки в случайно выбранную ячейку решетки. Светлой заливкой выделены ячейки, только что получившие песчинки, темной - только что опрокинувшиеся. На третьем шаге лавины опрокидывается ячейка с 5 песчинками, одна из которых утрачивается
В дальнейшем мы будем полагать верхний и нижний края решетки открытыми, а левый и правый отожествим, свернув ее в вертикальный цилиндр. Такая топология решетки наиболее удобна для анализа, поскольку система становится квазиодномерной, распадаясь на слои ячеек, находящихся на одинаковом удалении от края.
В качестве основных характеристик лавины рассматриваются ее длительность Т - число шагов моделирования, в ходе которых происходили опрокидывания, размер N - число событий опрокидывания и площадь Б - число опрокинувшихся ячеек (без учета кратности опрокидывания). Кроме указанных величин интерес представляют также диссипация В - число песчинок, утраченных при опрокидывании, и падение ¥ - число песчинок, вывалившихся за край.
Роль открытого края. У каж-
<Р> 0.6 0.5 0.4 0.3 0.2 0.1
0 10°
« ехр([1.1-/0]/4.39)
\ 10-2-
\ 10"4-
Ь= 28^15
\ 10-6- усреднение ^^
\ ю-8-
¿=4096 \ °...............20.............40.............60..............80"/0
101
10
ю3 1а
Рис. 6. Среднее число падений за край в зависимости от слоя вброса. По мере удаления точки вброса от края число падений стремительно (в 10 раз на 10 слоях) уменьшается. На врезке - аналогичные данные (собранные для решеток различного размера) в полулогарифмическом масштабе
Рис. 7. Слой центра масс области лавины в зависимости от слоя вброса. Среднее положение центра масс почти совпадает с точкой вброса инициирующей песчинки, так как лавины в своем большинстве невелики по размеру. Наибольшее зарегистрированное смещение центра масс от его среднего положения к краю решетки весьма незначительно (не превышает нескольких десятков слоев). Однако смещение центра масс от края макроскопически велико (может достигать порядка трети линейного размера решетки)
дой песчинки, попавшей в систему, есть два пути ее покинуть: диссипация при опрокидывании ячейки с избытком песчинок или падение за край при опрокидывании краевой ячейки. Чтобы понять соотношение этих двух путей, рассмотрим число падений в ходе одной лавины, инициированной вбросом песчинки на заданном слое 10, номер которого вычисляется как расстояние до ближайшего края решетки (1 < 10 < Ь/2).
Рис. 6 дает наглядное представление о том, что падение за край как механизм выбытия песчинок из системы значимо только для событий, начинающихся в непосредственной близости от края. Это обусловлено неконсервативностью правил, ограничивающей длину свободного пробега песчинок по решетке, в результате чего они покидают ее на сравнительно небольшом расстоянии от места вброса. Поэтому почти все песчинки выбывают в результате диссипации, а за край падают только те немногие, которым удается его достичь.
Вместе с тем, край решетки очень важен, поскольку он играет роль источника активности. Как можно видеть из рис. 7, лавина легко распространяется вглубь решетки, но весьма неохотно идет в обратную сторону.
Самые крупные лавины, происходящие в системе, начинаются вблизи ее края и завершаются в глубине. На рис. 8 показано распределение
средних размеров лавины по слоям вброса. Чем ближе к середине решетки вброшена инициирующая песчинка, тем меньше получается лавина, что связано с сокращением открывающегося перед ней пространства.
Рис. 9 демонстрирует типичный пример развития крупной лавины, которое, несмотря на изотропию правил модели, оказывается отчетливо анизотропным. Такое поведение модели в самоорганизованно критическом состоянии уместно характеризовать как спонтанную анизотропию. Она связана с распространением краевых эффектов на всю систему (точнее говоря, на ближайшую к данному краю ее половину).
Следует отметить, что спонтанная анизотропия DFF-модели, хоть и является очень похожей на анизотропию, изначально заложенную в правила DR-модели, все-таки отличается тем, что имеет место лишь в глобальном масштабе. При этом локальная динамика остается изотропной, что приводит к такому важному явлению, отсутствующему в DR-модели, как возвращение активности на уже пройденные лавиной слои. Как будет показано далее, именно оно обуславливает различия в значениях показателей распределений для этих моделей.
Кратность опрокидывания. Показанный на рис. 9 вариант развития лавины, идущей одним фронтом, является хотя и типичным, но не единственно возможным. Рис. 10 демонстрирует альтернативную ситуацию, когда фронтов оказывается более одного. Для их возникновения необходимо, чтобы лавина, инициированная вбросом песчинки на некотором отдалении от края, на начальном этапе успешно развивалась сразу во все стороны, а не только в глубину. При этом возможен ее распад на набор
Рис. 8. Средний размер лавины в зависимости от слоя вброса. Положение максимума не зависит от размера решетки. Он отстоит от края на такое расстояние, на котором пренебрежимо малы потери от падений, препятствующие развитию крупной лавины. После прохождения максимума график спадает линейно практически до нуля в середине решетки. На врезке - те же данные, но в полулогарифмическом масштабе
Рис. 9. Типичное развитие крупной лавины. Линией обозначены границы области лавины на четыре различных момента времени и после ее завершения, точками - места опрокидываний (за 8 предшествующих шагов). Лавина идет одним фронтом от края к середине решетки, где обрывается. Поле: Ь = 2048. Лавина: Т = 2166, Я = 407092, N = 407148, Б = 2764
Рис. 10. Развитие крупной лавины со множественными опрокидываниями. Лавина представляет собой последовательный набор приблизительно эквидистантных фронтов, идущих от края решетки к ее середине, преодолеть которую они не могут. Развитие лавины показано для трех различных моментов времени. Поле: Ь = 2048. Лавина: Т = 2022, Я = 284777, N = 1287730, Б = 1325
не связанных между собой лавин, каждая из которых в дальнейшем начинает движение в глубину решетки. Если они отстают одна от другой, то и получается несколько фронтов.
Реализация описанного сценария маловероятна, поэтому последовательные фронты опрокидывания наблюдаются сравнительно редко. И в силу анизотропного характера развития лавин редким явлением должны быть вообще повторные опрокидывания ячеек, что подтверждается видом распределения кратности опрокидывания, показанным на рис. 11.
Отсюда следует, что можно не различать размер лавин N и их площадь 5 (см. врезку к рис. 11). И в самом деле, как показывают результаты моделирования, эти характеристики лавины описываются одними и теми же показателями распределения.
Кроме того, поскольку пространственное развитие лавины направлено к середине решетки и ею же и ограничено, а повторных опрокидываний практически не происходит, длительность самых крупных лавин определяется числом затрагиваемых слоев
Т\ ~ Ь, или Уу = 1, (15) что полностью аналогично формуле (10).
Рис. 11. Роль множественных опрокидываний. Львиная доля лавин характеризуется отношением размера к площади, очень близким к 1. Лавинам, в ходе которых возник второй фронт, соответствует диапазон значений отношения размера к площади от 1 до 2 (чем дольше существовал второй фронт и чем большую часть площади он прошел, тем ближе отношение к 2), лавинам, в ходе которых возник третий фронт, - от 2 до 3 и т.д. На врезке - связь размера лавины с ее площадью. График условного среднего неотличим от тождественной зависимости. Условный максимум размера лишь незначительно превосходит его среднее значение за счет редких повторных опрокидываний
Распределение опрокидываний по слоям. Песчинки не могут далеко перемещаться по системе, поэтому каждая из них покидает ее примерно там же, где и была вброшена. Процессы добавления песчинок и их диссипационного выбытия уравновешивают друг друга всюду за исключением узких областей, находящихся в середине решетки и вблизи ее краев (рис. 12). Соответственно, должны уравновешиваться и процессы пере- Рис. 12. Распределение активности по слоям ре-
дачи/получения песчинок в ходе опро- шетки.Средняя диссипация равна одной выбывшей
песчинки на одну вброшенную на каждый слой
кидываний. (учитывается домножением на Ь) по всей решетке
Пусть N (¡) - число опрокидыва- за исключением краев и середины. На врезке - за-
ний, происшедших на слое I в ходе ла- висимость среднего числа отртквдьтшш на сл°е
при разных размерах решетки. Чем больше размер
вины. В их результате этот слой переда- системы, тем меньшую его долю составляют при-ет соседним по N (¡) песчинок, получая краевые и срединная области и тем точнее выполот них N (I - 1) + N (I + 1). Поскольку няется ф°рмула (17)
песчинки не могут накапливаться на каком-либо слое, выполнено балансовое соотношение
—2
ж* V'¡»
N (I -1)) - 2 N (¡)) + N (I +1)) = 0,
где угловые скобки означают усреднение по всем лавинам. Интегрируя его, получаем поток песчинок, направленный к краю решетки
й
— N (¡)) = ео^.
(16)
Возникновение потока обусловлено тем, что часть песчинок покидает систему через открытый край, подтягивая на свое место песчинки из более глубоких слоев. Примечательно, что, хотя активность распространяется от краев к середине, песчинки перемещаются навстречу ей - от середины к краям.
Учитывая, что при I = 0 (то есть сразу за краем) опрокидываний нет, интегрированием уравнения (16) получаем зависимость
N (¡,Ь))~ ¡¡Ь (17)
с коэффициентом пропорциональности, не зависящим от Ь в пределе бесконечного размера, в чем дает возможность убедиться врезка к рис. 12. Такой вид коэффициента объясняется тем, что среднее число опрокидываний, происходящих на слое за время Ь событий, соответствующее одному вбросу инициирующей песчинки на него, не зависит от размера системы. Интегрированием формулы (17) находим средний размер лавины
Ь/2
N) = 2/ N (¡,Ь)) —1 ~ Ь, или ом = 1. (18)
1
Как уже было сказано, хотя для БРР-модели площадь лавин и не равна их размеру, принципиальной разницы между этими величинами нет. То есть формула (18) выступает аналогом формулы (11) для БЯ-модели. Тем не менее, природа этих формул существенно различается.
Поскольку лавины распространяются от края в глубину, то среднее число опрокидываний в слое ¡, вызванных вбросом инициирующей песчинки на слой ¡о, есть
N (1,1о,ь)) = е (I - 1о) I (1,ь),
где е - функция Хевисайда. Интегрирование этой формулы по слоям вброса должно дать формулу (17), откуда немедленно получается
I (1,Ь) ~ 1/Ь,
что дает окончательное соотношение
N (¡,1о,Ь))~ ^^ (19)
для среднего числа опрокидываний, произошедших на слое I в результате вброса инициирующей песчинки на слой ¡о при размере системы Ь.
Полученное соотношение позволяет найти среднее число опрокидываний в системе в целом, вызванных вбросом на слое ¡о:
и /2
N (¡о, Ь)) = 2 у N (¡, ¡о, Ь))
Ь/2
77 1 ¡о
^ 2 - Ь
что объясняет линейное спадание правой части графика на рис. 8.
Распределение диссипации по слоям. По аналогии с формулой (19) для N (¡, ¡о, Ь)) найдем выражение и для {Б (¡, ¡о, Ь)) - среднего числа песчинок, дис-сипировавших в слое ¡ в результате вброса инициирующей песчинки на слое ¡о при размере системы Ь.
Средняя диссипация в слое, происходящая за время Ь событий, соответствующая одному вбросу инициирующей песчинки на этот слой, не зависит от Ь и равна 1. Соответственно, средняя диссипация в слое ¡ за один шаг
{Б (¡,Ь)) = 2 I {Б (¡, ¡о, Ь)) Мо = Ь■ (20)
Доля опрокидываний, сопровождающихся диссипацией, определяется шансами встретить локальный избыток песка, которые зависят только от номера слоя, то есть
{Б (¡, ¡о, Ь)) = N (¡, ¡о, Ь)) с (¡) ■
В сочетании с формулами (19) и (20) это дает
с (¡) ~ (21)
и, следовательно,
{Б (¡, ¡о, Ь))
е ^ - ¡о) ¡ь '
Полученная общая формула позволяет найти среднюю диссипацию, вызванную вбросом песчинки на слой ¡о
{Б (¡о, Ь)) =
Ь/2 1 Ь
= 2 1 {Б (¡,¡о,Ь)) —а ~ Ь ы—.
Рис. 13. Диссипация в зависимости от слоя вброса при разных размерах решетки. Графики чем-то схожи с изображенным на рис. 8 - выраженный пик на небольшом удалении от края с дальнейшим убыванием до нуля. Однако в отличие от среднего размера лавины, убывающего линейно с номером слоя вброса, средняя диссипация убывает лишь логарифмически. На врезке - те же зависимости, но с масштабированием по оси абсцисс на размер системы
Логарифмическое поведение этой величины подтверждается рис. 13.
Заполнение решетки. Наличие выделенного направления, в котором распространяются лавины, приводит к неравномерному распределению песка по слоям решетки, представленному на рис. 14. Распределение имеет антиинтуитивный вид. Казалось бы, меньше всего песчинок, способных участвовать в процессах пересыпания, должно быть у открытого края, а больше всего - в глубине. Однако возрастание их количества по мере удаления от края происходит всего лишь на нескольких слоях решетки, сменяясь затем плавным убыванием, которое продолжается почти до самой ее середины, где становится скачкообразным. В середине отрыт своеобразный ров - область шириной в несколько
десятков слоев со значительным дефицитом песка (см. врезку к рис. 14).
Срединный ров заполнения может служить яркой и наглядной иллюстрацией самоорганизационного возникновения целостных свойств, то есть свойств, имеющихся у системы, но в принципе отсутствовавших у ее составных частей. Правила ВБР-модели вполне локальны - при опрокидывании неустойчивой ячейки песчинки из нее передаются лишь ее ближайшим соседям. Более того, дистанция, проходимая каждой отдельной песчинкой за время ее пребывания в системе, ограничена несколькими десятками ячеек. Тем не менее система умеет найти свою середину, отстоящую от краев на многие тысячи ячеек. Это воистину удивительно!
Рис. 14. Среднее заполнение ячеек по слоям решетки. Заполнение максимально вблизи края, а по мере удаления от него уменьшается, стремясь к предельному значению = 2. Вблизи середины решетки график имеет резкий провал, увеличенный на врезке
Структура заполнения решетки вполне определяется особенностью распределения активности по слоям. Так, ров отрывается потоками песка, направленными от середины решетки к ее краям. При этом следует заметить, что за счет флуктуаций заполнения лавины могут немного перехлестывать за середину решетки, вынося там песчинки на другую сторону рва, то есть перемещая их сразу по направлению к обоим краям. В противном случае ров не образовывался бы, поскольку правила модели изотропны и суммарный перенос песка любой лавиной равен нулю.
Несколько менее очевидной оказывается ситуация с уменьшением заполнения по мере удаления от края решетки, что находится в кажущемся противоречии с преимущественным распространением лавин в ее глубину. На первый взгляд, при опрокидывании ячейки те из ее соседей должны чаще терять устойчивость, которые имеют в среднем большее заполнение, то есть лавина, скорее, должна идти туда, где песка больше, а не меньше.
Однако, как уже было отмечено ранее, поддержание лавины не столь требовательно к заполнению ячеек, как диссипация. Поэтому именно последняя управляет структурой заполнения. По мере удаления от края увеличение числа опрокидываний, происходящих на слое, ведет к уменьшению вероятности диссипации при единичном опрокидывании, что возможно лишь при соответствующем уменьшении заполнения ячеек.
Естественным представляется предположение, что доля опрокидываний, сопровождающихся диссипацией, прямо пропорциональна избыточному заполнению ячеек, то есть
с (0 ~ г-
В силу формулы (21) левая часть обратно пропорциональна ¡. Поведение правой, переданное на рис. 15, подтверждает сделанное предположение.
Решение модели. Для определения набора показателей модели осталось найти связь между размером лавины и ее длительностью.
На каждом слое ¡, затронутом лавиной, инициированной на слое ¡о, определим ее ширину и> (¡, ¡о) как число опрокинувшихся ячеек данного слоя. Как и в случае БЯ-модели, пребывание системы в самоорганизованно критическом состоянии влечет диффузионное изменение ширины лавины по мере ее продвижения по слоям. Однако, в силу изотропии правил БРР-модели, в ней возможно возвращение активности на уже пройденные слои, в силу чего формула (13) оказывается неприменимой.
Рис. 15. Избыточное среднее заполнение ячеек. В глубине решетки среднее заполнение приближается к предельному по гиперболическому закону. Область построения графиков ограничена серединой решетки
Далее будет показано, что
1-5/3
ш (¡,к) - ^ - ¡о)2/3 ■ (22)
ш
Поскольку, как уже отмечалось, различия между площадью лавины и ее размером непринципиальны, суммированием ширин по слоям находим
г0+т
10"8
N (Т) - Б (Т) - / ш (¡^о) <И - Т5/3,
/
10"4 10"3 10"2 ю-1 ьл-т
10
что подтверждается рис. 16.
Таким образом, показатель связи
Рис. 16. Зависимость размера лавины от ее длительности для БРР-модели. Переменные отмасштаби-рованы на размер системы со скейлинговыми показателями ут = 1 и = 5/3. Ср. с рис. 4
Умт = 5/3,
(23)
что замыкает систему уравнений {(4), (7), (9), (15), (18), (23)}, из которой находятся скейлинговый показатель = 5/3 и характеристические показатели ат = 2/3 и аN = 2/5. Рис. 17 дает возможность убедиться в правильности полученных значений.
Универсальный показатель е = 2/3 оказывается для БРР-модели больше, чем для БЯ-модели. Соответственно, доля лавин, не помещающихся в решетку, - меньше, что понятно, поскольку развитие лавин ограничивается не только конечностью размеров системы, но и диссипацией.
Рис. 17. Распределение лавин в БРР-модели по длительности и размеру. Совмещение графиков при конечно-размерном скейлинге демонстрирует правильность определения скейлинговых показателей, а появление горизонтальных участков на графиках со скомпенсированной степенной компонентой плотности (врезки) - правильность определения характеристических показателей. В отличие от графиков на рис. 3 здесь плотность является монотонно убывающей функцией
Если в случае БЯ-модели можно строго показать, что область лавины является односвязной и не имеет дырок, то в случае БРР-модели односвязность, вообще говоря, ниоткуда не следует. Однако возникновение дырок в этой модели является,
4. Аномальная диффузия
по-видимому, событием весьма маловероятным. В ходе компьютерного исследования модели автору ни разу не встретились дырявые лавины, зато неоднократно наблюдалось зарастание уже, казалось бы, вполне оформившейся дырки. Поэтому будем считать, что дырок нет, и изменение ширины области лавины обусловлено исключительно сдвигом ее границ при прохождении лавины по слоям.
Уширение области лавины. Без потери общности можно следить лишь за одной границей области лавины, поскольку полностью аналогичные эволюции второй границы могут сказаться лишь на коэффициенте, но не на показателе зависимости (22). Соответственно, нет необходимости различать ширину области лавины и положение ее границы. Обозначим последнюю через Wl.
Сдвиг границы при переходе лавины со слоя ¡ на слой ¡ +1 представляет собой шаг случайного блуждания
"^+1 = Wl + 81, (24)
которое является несмещенным
{81) = ^I) = 0 (25)
в силу природы самоорганизованно критического состояния. Кроме того, оно должно быть и некоррелированным
{wl8l) = {81), (26)
так как лавина «не знает» своей ширины.
Принципиальным отличием ЭРР-модели от ЭЯ-модели является возможность многократного сдвига границы области лавины вовне из-за возвращения активности на уже посещенные слои. При этом блуждание оказывается несимметричным. Так, если происходит уменьшение ширины области лавины, ее граница сдвигается внутрь один раз на некоторое расстояние, которое для простоты примем за единицу. А если происходит увеличение ширины, то граница может сдвигаться несколько раз, фактически совершая скачок на некоторое расстояние щ. Иными словами,
81 Л (27)
1 -1 1 - Р1,
где вероятность совершения скачка однозначно определяется из условия несмещенности блужданий (25)
VI (П1 + 1) = 1. (28)
Возведем формулу (24) в квадрат и усредним по реализациям
К+1> = + 2 ^81) + (82).
Как >ледует из формул (25) и (26), второй член в правой части равен нулю, а (82) = п1 в силу формул (27) и (28).
Таким образом, изменение среднего квадрата ширины при переходе к следующему слою совпадает с длиной скачка, который может совершить граница на данном слое ( )
А = п.
Если предположить степенное возрастание длины скачка по мере удаления от
стартового слоя
П1 ~ ^ - ¡о)
(29)
то {w
=
п
(¡ - ¡о)
1+а
что ведет к степенной зависимости ширины от
номера слоя
1 + а
w (¡) ~ ^ - ¡о) 2
(30)
В тривиальном случае а = 0 получается формула (13), а при нетривиальных значениях а возникает аномальная диффузия, характеризуемая значением показателя в формуле (30), отличным от 1/2.
Однако и нетривиальные значения а, и даже сама зависимость (29), явным образом предполагающая, что блуждатель умеет измерять сделанное им число шагов, нуждаются в объяснении, которое в общем случае предложить не удается.
Самосогласованный рост ширины. Для ЭРР-модели величина скачка определяется, конечно, не числом пройденных слоев, а локальными условиями. Возможность совершения границей лавины длинного скачка вовне ее области обуславливается возвращениями активности со следующих слоев. Чем больше раз активность вернется, принося на слой дополнительные песчинки, тем длиннее может быть скачок.
Успешному, то есть приводящему к новым опрокидываниям, возвращению активности способствует медленный рост ширины области лавины по слоям, поскольку вдоль менее крутой границы активность может возвратиться на большее число слоев. Схематично эта ситуация представлена на рис. 18.
Рис. 18. Пример возвращения активности. Стрелками показана передача песчинок при опрокидывании ячеек. Светлой заливкой отмечены ячейки, опрокинувшиеся при первом же прохождении лавины по слоям, промежуточной - ячейка, опрокидывание которой происходит при возвращении активности на один слой, а темной - ячейка, опрокидывание которой происходит (вверху) или не происходит (внизу) при возвращении активности на два слоя. Жирной линией показана граница области лавины. Чем быстрее она смещается вовне при движении по слоям, тем ниже шансы на успешное возвращение активности
В некотором приближении можно считать, что длина скачка пропорциональна обратному локальному градиенту границы (или, что то же самое, числу слоев, на которых активность может вернуться, сдвигая границу, за один раз):
gH ~ (l - lo)~ . (31)
Сопоставление полученного выражения с формулой (29) приводит к нетривиальному значению a = 1/3. Соответственно, формула (30), определяющая зависимость ширины лавины от номера слоя, принимает вид (22).
Выводы
Модель Дхара-Рамасвами и дискретная модель Федеров имеют анизотропную динамику распространения активности, что делает их свойства во многом схожими и позволяет исследовать эти модели одними и теми же методами. Однако, если в случае первой модели анизотропия явным образом заложена в правила, то в случае второй она возникает как результат самоорганизации в критическое состояние в условиях диссипации. При этом лавины распространяются от открытого края в глубину решетки, в середине которой возникает непреодолимый для них ров заполнения.
Методика решения самоорганизованно критических моделей основана на совместном применении конечно-размерного скейлинга для двух различных характеристик происходящих событий. Анализ правил моделей дает возможность установить для каждой из этих характеристик по одному скейлинговому показателю, а учет природы самоорганизованно критического состояния дает показатель связи между ними.
Анизотропное распространение лавин позволяет рассматривать формирование их области как процесс блуждания границы при продвижении лавины по слоям решетки, которое в критическом состоянии оказывается несмещенным. Изменение ширины области лавины представляет собой диффузию - обычную для модели Дхара-Рамасвами и аномальную для дискретной модели Федеров. Аномальный характер диффузии в последнем случае объясняется возможностью возвращения активности на уже пройденные слои, вследствие чего блуждание границы становится несимметричным.
Полагая длину скачка, совершаемого границей вовне области лавины, обратно пропорциональной градиенту границы, удается замкнуть уравнения для аномальной диффузии, получая нетривиальное значение ее показателя без предположений о лежащих в ее основе каких-либо степенных закономерностях.
Работа выполнена при поддержке РФФИ (проекты 10-01-00786-а и 11-01-00887-а).
Библиографический список
1. BakP., Tang C., Wiesenfeld K. Self-organized criticality: An explanation of 1/f-noise
// Phys. Rev. Lett. 1987. Vol. 59, № 4. P. 381.
2. Bak P., Tang C., Wiesenfeld K. Self-organized criticality // Phys. Rev. A. 1988. Vol. 38, № 1. P. 364.
3. Ivashkevich E.V. Critical behavior of the sandpile model as a self-organizing branching process// Phys. Rev. Lett. 1996. Vol. 76, № 18. P. 3368.
4. DharD., Ramaswamy R. Exactly solved model of self-organized critical phenomena// Phys. Rev. Lett. 1989. Vol. 63, № 16. P. 1659.
5. Christensen K., Olami Z. Sandpile models with and without an underlying spatial structure// Phys. Rev. E. 1993. Vol. 48, № 5. P. 3361.
6. Dhar D., Majumdar S.N.Abelian sandpile model on the Bethe lattice // J. Phys. A: Math. Gen. 1990. Vol. 23. P. 4333.
7. Ма Ш. Современная теория критических явлений. М.: Мир, 1980. 298 с.
8. Manna S.S. Two-state model of self-organized criticality// L. Phys. A: Math. Gen. 1991. Vol. 24. P. L363.
9. Подлазов А.В. Двумерная самоорганизованно-критическая модель Манна. М.: ИПМ им. М.В. Келдыша, 2012. Препринт № 42. 20 с. http://library.keldysh.ru/preprint.asp?id=2012-42
10. Подлазов А.В. Сравнение двумерных изотропных консервативных самоорга-низованно-критических моделей типа кучи песка. М.: ИПМ им. М.В. Келдыша, 2012. Препринт № 43. http://library.keldysh.ru/preprint.asp?id=2012-43
11. Малинецкий Г.Г., Подлазов А.В. Сравнение двумерных изотропных консервативных самоорганизованно-критических моделей типа кучи песка// Вестник МГТУ им. Н.Э. Баумана. Естественные науки. 2012. Спец. выпуск № 2 «Математическое моделирование в технике». C. 119.
12. Kadanoff L.P., Nagel S.R., Wu L., Zhou S. Scaling and universality in avalanches// Phys. Rev. A. 1989. Vol. 39, № 12. P. 6524.
13. Bak P. How nature works: The science of self-organized criticality. Springer-Verlag, New York, Inc. 1996.
14. Sornette D., Johansen A., Dornic I. Mapping self-organized criticality onto criticality // J. Phys. I (France). 1995. Vol. 5. P. 325.
15. Clar S., Drossel B., Schwabl F. Forest fires and other examples of self-organized criticality// J. Phys.: Cond. Mat. 1996. Vol. 8. P. 6803.
16. Feder H.J.S., Feder J.Self-organized criticality in a stick-slip process// Phys. Rev. Lett. 1991. Vol. 66, № 20. P. 2669.
ИПМ им. М.В. Келдыша РАН, Поступила в редакцию 13.04.2012
Москва После доработки 28.11.2012
TWO-DIMENSIONAL SELF-ORGANIZED CRITICAL SANDPILE MODELS
WITH ANISOTROPIC DYNAMICS OF THE ACTIVITY PROPAGATION
A. V. Podlazov
We numerically and analytically investigate two self-organized critical sandpile models with anisotropic dynamics of the activity propagation - Dhar-Ramaswamy and discrete Feder-Feder models. The full set of critical indices for these models is theoretically determined.
We also give systematical statement of the finite-size scaling ansatz and of its use for the solving of self-organized critical systems.
Studying the discrete Feder-Feder model we find and explain a number of nontrivial phenomena, such as spontaneous anisotropy, anomalous diffusion and the appearance of midline ditch of filling.
Keywords: Self-organized criticality, sandpile models, scale invariance, power laws, finite-size scaling, anomalous diffusion, spontaneous anisotropy.
Подлазов Андрей Викторович - родился в Москве (1973), окончил Московский физико-технический институт (1996). После окончания МФТИ работает в Институте прикладной математики им. М.В. Келдыша РАН старшим научным сотрудником. Защитил диссертацию на соискание ученой степени кандидата физико-математических наук в ИПМ им. М.В. Келдыша РАН (2001) в области теории самоорганизованной критичности и теории масштабно-инвариантных процессов. Автор монографий «Управление риском. Риск, устойчивое развитие, синергетика» [М.: Наука, 2000. 432 с.] (в соавторстве с В.А. Владимировым, Ю.Л. Воробьевым, Г.Г. Малинецким и др.) и «Нелинейная динамика: Подходы, результаты, надежды» [Изд.3/ Синергетика: от прошлого к будущему. М.: ЛИБРОКОМ, 2011. 280 с.] (в соавторстве с Г.Г. Малинецким и А.Б. Потаповым). Опубликовал 60 научных статей по направлениям, указанным выше, а также по демографии, математической истории и вопросам развития системы образования.
125047 Москва, Миусская пл., 4
Институт прикладной математики им. М.В. Келдыша РАН E-mail: Tiger@Keldysh.ru
