CC BY
13
УДК 53.072+550.34.013.4
ИСПОЛЬЗОВАНИЕ САМОАФИННОГО ПРЕОБРАЗОВАНИЯ В МОДЕЛИ ОЛАМИ-ФЕДЕРА-КРИСТИАНСЕНА ПРИ РАСЧЕТЕ КРИТИЧЕСКОГО СОСТОЯНИЯ
© 2008 г. А.С. Черепанцев
The opportunity of use self-affine transformation for spatial distribution of energy for acceleration of achievement of a critical state in Olami-Feder-Christensen model is considered. The given problem is actual for small values of dissipation parameter and modelling calculations on the big grids. Efficiency of offered method for calculation a critical state has power law growth with increase the size of a grid.
В соответствии с традиционным классическим взглядом на геофизическую среду литосферы как сплошную с конечной прочностью по отношению к деформации, возникновение землетрясения обусловлено образованием плоского разрыва. Вместе с тем в настоящее время активно развивается представление о среде как системе различных неоднородных объектов - блоков (литосферная плита, отдельности горной породы, трещиноватость, дефекты кристаллической решетки). Континуальная модель позволяет получить адекватное представление об отдельных, единичных процессах в среде (образование и распространение разрыва, виды сейсмических волн, их спектральный состав и т.д.). Но вместе с тем при данном подходе возникает проблема при описании таких основных закономерностей сейсмического режима, как локализация землетрясений в пространстве, стационарность режима землетрясений, автомодельный характер закона повторяемости землетрясений, неоднородное, самоподобное распределение землетрясений в пространстве и во времени [1, 2].
Для данного подхода характерны следующие положения:
- реальная геофизическая среда представляет собой дискретную систему разномасштабных неодно-родностей;
- данная система имеет иерархическое строение. При этом отношение соседних характерных размеров находится в диапазоне 2-5 [1];
- с учетом того, что в литосфере происходит постоянное накопление упругой энергии, ее диссипация
определяет структурное строение отдельных элементов системы неоднородностей;
- способ диссипации энергии определяется величиной поступающей упругой энергии. В случае малых величин прироста упругой энергии она диссипирует на младших масштабах неоднородностей путем релаксации, микротресков. При значительных притоках энергии в сбросе упругой энергии принимают участие старшие масштабы, что влечет макросбросы и связанные с ними землетрясения;
- развиваясь на дискретном, иерархически самоподобном носителе (блоках среды), сейсмический процесс также имеет самоподобные свойства. Автомодельный вид имеют распределения землетрясений в энергетической, пространственной и временной областях [2]. На это указывает обобщенный закон повторяемости землетрясений [3]:
= + + + (1) где N - число землетрясений с энергией Е, произошедших в области с линейным размером Ь в течение времени Т; Ь, й - показатели самоподобия в энергетической и пространственной областях.
Исходя из данной концепции, иерархическая дискретность геофизической среды является естественной формой ее адаптации к условиям открытой динамической системы. Предлагаемая концепция определяет новый взгляд и на моделирование сейсмического режима исходя из фрактальный свойств геометрической системы разломов, блоков и плит. В этом случае фрактальная структура задана исходно. Путем моделирования геологической системы и затем моделирования ее однородного развития возможно получение
пространственно временной структуры напряженного состояния и сейсмического режима [4, 5]. Вопрос, который при этом остается за рамками модели, заключается в выяснении причин возникновения иерархических фрактальных структур неоднородностей земной коры.
Существует другой подход, при котором рассматривается динамика отдельных плит. В этом случае пространственная фрактальная структура не задана исходно. Организация сейсмичности является следствием некоторой эволюционной динамической системы. Такой подход развивается на основе концепции самоорганизации системы в устойчивое динамическое состояние со степенным характером распределения характеристик, являющееся аттрактором данной системы.
Простейшая физическая модель подобной динамической системы была предложена R. Burridge, L. Knopoff [6]. Она представляет собой двумерную пружинно-блоковую модель с постоянной упругой внешней силой и конечной силой трения покоя. В случае малости сил, действующих по направлению, перпендикулярному к направлению внешней силы, малости внутренних деформаций блоков и анизатропии упругих свойств системы, она сводится к к клеточной модели Олами-Федера-Кристиансена (OFC) [7].
Большое внимание исследователей к OFC-модели определяется тем, что при достаточно общих условиях она демонстрирует эволюцию к предельному устойчивому степенному распределению характеристик, свойственному системам в критическом состоянии. При этом можно выделить ряд принципиальных свойств, приближающих ее к поведению моделируемого процесса динамики геофизической среды:
- критическое состояние в модели OFC достигается при отсутствии управляющего параметра, т.е. самоорганизуется;
- в отличие от ранее полученных моделей с самоорганизованным критическим состоянием (напр. sand pile model) данная модель является диссипативной. При сбросе упругой энергии отдельным блоком осуществляется передача только части сброшенной энергии соседним блокам. Такое условие является принципиальным для достижения критического состояния;
вместе с тем необходимым условием возникновения самоорганизованного критического состояния (СОК) является существование взаимодействия между соседними блоками при сбросе упругой энергии одним из них. При отсутствии передачи энергии соседним блокам при сбросе состояние СОК также не наблюдается;
еще одним условием достижения критического состояния являются условия на границе двумерной системы блоков. При задании периодических условий на границе и отсутствии характерного пространственного размера в системе отсутствует сходимость к критическому состоянию. В модели обычно используются открытые или свободные граничные условия. В этом случае граничные блоки представляют собой неоднородные образования и структурная организация в системе происходит от граничных областей.
На основе многочисленных модельных экспериментов для системы с открытыми граничными условиями получены оценки распределения повторяемости сбросов по их размеру (под размером понимается количество блоков, участвующих в отдельном акте сброса). Полученный степенной показатель т«0,8 зависимости
функции плотности распределения Г 9] бли-
зок к значению 2у. где у~0.4-0.5 - эмпирическое значение степенного показателя функции плотности распределения землетрясений по энергии Л'С- ^ ;А (удвоенное значение берется исходя из того, что размер сброса пропорционален сброшенному значению упругой силы ~ А/7 ~ ).
Связь степенных показателей в модели OFC Рассмотрим выполнение еще одной статистической закономерности сейсмического процесса, выявленного в результате обработки имеющегося материала регистрации землетрясений в различных регионах. Связь между частотой землетрясений и их магни-тудой [10], сейсмическим моментом и магнитудой [11] и сейсмическим моментом и линейным размером очага [12] позволили получить соотношение между частотой землетрясений N и линейным размером очаговой зоны г в виде [13]
N = ¡3-г~1Ь , (2)
где Ь - степенной показатель в соотношении Гуттен-берга-Рихтера; Р - константа. Соотношение (2) определяет величину фрактальной размерности в пространстве как
£> = 26. (3)
С целью получения сравнительных данных с имеющимися геофизическими характеристиками сейсмического процесса, на основе модели ОБО проведена генерация пространственного распределения сбросов с заданием координат на поверхности событий различного масштаба.
Методика генерации заключается в получении для заданных параметров модели (Ртт. - максимально возможной величины параметра системы, эквивалентного силе или энергии; АЕ - величине приращения параметра за один шаг итерации; а - параметра диссипации в системе; Ь - линейного размера системы) критического состояния системы после заданного количества итераций N и дальнейшего накопления в течение определенного числа итераций АЛ' координат сбросов (х,у). Для определения координаты сброса масштаба 5, где 5 - число связанных сбросов при одиночном шаге приращения энергии использован принцип центра масс, в соответствии с которым координаты 1-го сброса масштаба 5 определяются как
xi = s , Уг = ЪУг,
где xy, yy - соответствую-
1 / м
щие координаты единичных элементов, входящих в 1-й сброс старшего масштаба. Полученный таким образом каталог представлен на рис. 1.
Для характеристики структурной неоднородности распределения и группируемости событий принято использовать фрактальную размерность полученного множества.
9,110-1 2,7 8,0 2,3-10-1 7,0-10-1 lg I -6,910-2Г
d^A
Рис. 1. Распределение сбросов на сетке Ь=200 в модели ОЕС при достижении критического состояния при длительности накопления (=300 итераций
Для измерения фрактальной размерности множества координат сбросов можно использовать два типа оценки: клеточную и корреляционную размерности. Клеточная размерность оценивается исходя непосредственно из определения фрактальной размерности как показателя самоподобия. Разобьем рассматриваемую пространственную область на клетки размером ^ и Тогда отношение числа клеток, содержащих хотя бы
одно событие для этих разбиений: N^N/ = ^ .
Следовательно, число непустых клеток Л' / ''.
В случае, если размерность й=г, где г - размерность евклидова пространства, в котором задано множество, то события распределены равномерно, если й<г, то это означает неравномерность распределения событий, т.е. с уменьшением размера клетки плотность событий р = п/1г возрастает (п - среднее число событий в клетке размером I). Исходя из этого, клеточная размерность оценивается как d0 = -limlgN/lgl. Корре-
/-> о
ляционная размерность определяется через корреляционный интеграл [14]: d?= 1пп ЬГХл'й/ •
/-> о
Для однородных фрактальных множеств имеет место равенство d0=d2. В общем случае d0 > дооценка й2 предпочтительна с точки зрения требуемых объемов данных, так как при длине выборки Ь анализируется I: пар событий.
На рис. 2 представлена зависимость 1оя,С2 .
Область линейного скейлинга определяется единичными размерами расстояний между ячейками. Значение корреляционной размерности по данным вычислений равно Л=1,4±0,1. Полученное значение й2<2 может служить характеристикой неоднородности распределения сбросов, их группируемости в пространстве.
-7,240-1
-1,5
-2,3
-3,1 -
1ёС20,')
Рис. 2. Оценка фрактальной размерности Аи1,4 распределения сбросов на сетке по наклону корреляционного интеграла С2
Рассмотрим еще одну важную характеристику поведения данной модели, являющейся простейшей моделью сейсмичности, - связь энергии сброса (масштаба сброса) и его размеров. В качестве размера сброса г будем рассматривать больший из линейный размеров сброса: /; =тах^гг,Дуг . В случае однородного распределения по масштабам естественно предполагается Е = /•.',„ахЛ' ~ Етахг2, где - площадь
сброса. На рис. 3. представлена зависимость Е(г) в двойном логарифмическом масштабе. Величина показателя степени у=-1,8±0,1.
1
-1,5
-2,5
7= / -1,8
у < Ч Г ff,
щ Г*
1,1
1,2
1,3
1,4
1,5
lg г
Рис. 3. Распределение числа сбросов по размеру сброса
Полученные в результате модельных расчетов величины показателей т, с/2, у позволяют сравнить (3) с
соотношением по модельным расчетам: — = 1,8 + 0,3 .
г
2
3
1
Данное значение с учетом погрешности совпадает с соотношением (3) с/2 =2-т. что указывает на близость модельных представлений к реальному сейсмическому процессу. При этом справедливо соотношение для модельных коэффициентов:
А?2 + т ■ у = 0 ± 0,3 . (4)
Погрешность оценки наклонов кривых в двойном логарифмическом масштабе в данном случае не позволяет более точно получить значение <Л2/т и требуются дальнейшие исследования для выяснения его связи с параметром у.
Возможность ускорения расчета эволюции системы к критическому состоянию
В соответствии с концепцией СОК, за конечное время т система достигает стационарного критического состояния, для которого характерен степенной характер функции плотности распределения масштабов сбросов я. Система находится в критическом состоянии, если для вероятности сброса справедливо масштабное соотношение [15]:
р^УГЫ^
(5)
где ¥(х) - масштабная функция, являющаяся константой до некоторого значения размера пространственного масштаба Ьтах, определяемого размером сетки Ь\ ¡5, 8- критические показатели распределения. Данное соотношение определяет степенной характер и зависимости Р($):
= s-ß/*F\ S} = ^f\ S
,1?) и°
где Т7 '''' - новая масштабная функция.
Таким образом, критический показатель распределения г], определяющий Р(й), можно представить как г) = Р/З. Следует отметить, что (4) формально совпадает с полученной связью степенных коэффициентов в модели состояния СОК.
Соотношение (6) позволяет построить самоафин-ное преобразование функции распределения вида
(7)
где А - показатель Гельдера.
Автомодельный характер полученной зависимости распределения сбросов по масштабу дает возможность построить алгоритм ускорения расчета критического состояния на больших сетках.
Процесс расчета модели заключается в исходном задании в каждом узле решетки некоторого случайного состояния и повторения итерационной процедуры, учитывающей сброс состояния отдельного узла в нулевое значение при достижении критического значения с передачей части значения параметра (силы или энергии) соседним узлам. Итерационный процесс продолжается до возникновения устойчивой во времени статистики возникновения сбросов различных масштабов - плотности распределения размеров сброса. Под размером сброса, произошедшим на к-м шаге итерационной процеду-
(6)
ры, принято понимать число соседних связанных ячеек, достигших критического состояния. Такое стационарное состояние достигается очень медленно: требуется около 108 - 109 сбросов в системе для достижения критического состояния на двумерной сетке размером Ь ~ 102. Причем необходимое число итераций возрастает по мере уменьшения диссипативного параметра а<1/4 [16]. Данная проблема становится серьезной при рассмотрении моделей с малым значением параметра а. В настоящее время нет общего мнения о поведении системы при значениях а, близких к 0. Так, критическое минимальное значение с^, при которых исчезает предельный степенной характер распределения, по разным данным варьируется от 0,05 до 0,16 [17, 18]. Большие затраты машинного времени для достижения состояния самоорганизованной критичности на больших сетках приводят также к противоречивым выводам о поведении системы в зависимости от величины Ь.
Пусть Р(я) - функция распределения размеров сброса при некотором фиксированном размере сетки
( \
Lo. PiyPi,LoyC-S-iF
KLo
Здесь константа С есть вероятность сброса при малых значениях s при заданных параметрах модели Ь0 и а (в реальных вычислениях обычно при s=smjn=1).
Увеличим размер сетки: a>1. Рассмот-
рим распределение размеров сбросов: s1 = a3 s . Для новых аргументов (s1,L1):
pO P<i, О <W F
f \ s
= С= С • а'РрЦ, где С = С^С .
1Л )
Сравнивая с (7), с учетом Ь = а3, для показателя Гельдера получаем значение Л= -г].
Таким образом, в случае, если получена СОК система на сетке размером /,,,. то путем ее растяжения в а
раз и увеличения величин сбросов в а8 раз, получаемая система по статистическим параметрам на старших масштабах сбросов также будет удовлетворять условию критической системы.
Проведенное большое количество модельных опытов [8, 19] показывает, что в двумерной ОРС-модели во всем диапазоне возможных значений а величина /7« 1,8, величина у подчиняется условию у < 2, оставаясь близкой к данному значению и характеризуя фрактальную размерность распределения сбросов на сетке.
Пусть Ь задано как Ь=2п, что определяется прежде всего удобством компьютерного представления. Зададим начальный размер сетки Ь0=2к (1<к<п). Оценка времени счета критического состояния составит [20]:
Tdir ~ А '
^тах | jïx+<f>„¡¿(г^сЩг ä
АЕ
а
(8)
где р да -1, % да 1,8, <р да -1,5 - степенные показатели, полученные в результате модельных расчетов; Т -среднее время счета одной итерации.
Совершим самоафинное преобразование с Ь=2 и Л= -2, перейдя к сетке размером Ь1=2Ь0. При таком преобразовании критическое состояние системы и соответственно степенной характер распределения будут нарушены за счет исчезновения сбросов амплитудой л/г и малого искажения степенного распределения при 2 (/7 «1,8).
Так как модель ОБО представляет самоорганизующуюся систему, то через конечное число итераций система вернется в критическое состояние. Оценим время восстановления полученной системы:
Tu «А-
Е„
АЕ
hiJ^e^-^T.t,.
а
Здесь учтено, что требуется конечное число итераций для настройки системы. С учетом того, что среднее время установления определяется достижением сброса величины 5=16±4, данная величина может быть представлена в виде
\ Ьтт да = 4 при Ь < 4 [ Ь при Ь> 4 Общее время расчета на сетке размером Ь1 соста-
hb
T,
Li
'11 L0
1
Е
/ Л\ max
i(-4a
2 х+Ч>
+ t,
n-k
i=1
(10)
Эффективность данного метода [20] можно охарактеризовать отношением времени полного расчета критического состояния на решетке размером /, и (8):
Тт ¿¡-к v-\-/p 2
1
п—к j
2
(11)
В зависимости от выбора к<п каждое из слагаемых является определяющим в параметре эффективности. На рис. 4 представлена зависимость значения эффективности от выбора начального размера сетки
Ьо(к)=2к.
Как следует из приведенной зависимости при выборе начального размера сетки с малым числом элементов эффективность определяется вторым слагаемым в (10), в противном случае основной вклад вносит первое слагаемое, и эффективность алгоритма становится существенно ниже.
Данный алгоритм по сравнению с прямым методом слабо зависит от значения диссипативного параметра а. В отличие от прямого метода полный расчет критического состояния с существенной зависимостью от а происходит только на малой сетке (первое слагаемое в (10)), а в дальнейшем используется самоафинное преобразование, не связанное с параметром а.
АЕ
ч-АСТ^С,,^. (9)
Здесь не учитывается время, необходимое на «растягивание» решетки в 2 раза в силу его малости по сравнению со слагаемыми в (9).
Проделывая аналогичную последовательность действий, состоящих из увеличения решетки в 2 раза и проведения корректировочного счета, для вычисления критического состояния на решетке размером Ь требуется время:
Рис. 4. Зависимость значения эффективности от выбора начального размера сетки Ь0(к)=2к. Кривая 1 соответствует полному значению эффективности, 2 - вкладу первого слагаемого (11), 3 - вкладу второго слагаемого (11). Параметры модели: а=ОД5,Ь=210
Рассмотренный метод предназначен прежде всего для ускорения расчета самоорганизованного критического состояния при больших значениях Ь. На рис. 5 представлена зависимость эффективности работы
метода от размера выбранной решетки Ь=2п при двух
к
различных значениях размера начальной сетки Ь0=2к.
Рис. 5. Зависимость эффективности метода от размера
выбранной решетки Ь=2п. Кривая 1 соответствует начальной сетке Ь0=2; кривая 2 - начальной сетке Ь0=25
Как следует из полученных зависимостей, эффективность данного метода существенно повышается с увеличением размера системы, оправдывая смысл его использования для поставленной задачи.
e
о
В качестве иллюстрации работы данного алгоритма на рис. 6 представлены расчетные зависимости ппр числа сбросов масштаба 5 во временном окне Т=3000 для последовательности увеличивающихся сеток.
Рис. 6. Распределение числа сбросов масштаба 5 во временном окне Т=3000 итераций в системе Ь=512, «=0,18. Кривые 1-6 - зависимости п(5) последовательности увеличивающихся сеток, начиная с Ь0= 16; 7- зависимость, полученная прямым методом после Л?=3-10о итераций; 8 - линия с наклоном 8 =1,8, характерным для систем в критическом состоянии. Расчеты получены во временном окне Т=1000 итераций при использовании оптимизирующего алгоритма с вычислением на каждом шаге N=1000 итераций
Хорошее совпадение полученных значений показателя степени зависимости п(5) при использовании рассмотренного алгоритма со значением при прямом расчете может служить подтверждением обоснован-
ности рассмотренной методики расчета критического состояния в модели OFC.
Литература
1. Садовский М.А., Болохвитинов Л.Г., Писаренко В.Ф. Деформирование геофизический среды и сейсмический процесс. М., 1987.
2. Садовский М.А., Писаренко В.Ф. Сейсмический процесс в блоковой среде. М., 1991.
3. Кейлис-Борок В.И., Кособоков В.Г., Мажкенов С.А. // Вычислительная сейсмология. 1989. № 22. С. 28-40.
4. Soloviev A.A., Vorobieva I.A., Panza G.F. // Pure and Appl. Geophys., 2000. Vol. 157. Р. 1-2, 97-110.
5. Бабешко В.А., Бабешко О.М. // Экологический вестн. науч. центров Черноморского экономического сотрудничества. 2005. № 2. С. 16-22.
6. Burridge R., Knopoff L. // Bull. Seism Soc. Am. 1967. Vol. 57. P. 341-371.
7. Olami Z., Feder H.J.S., Christensen K. // Phys. Rev. 1992. Lett. 68. Р. 1244-1247.
8. Lise S., PaczuskiM. // Phys. Rev. 2001. E 63. Р. 36111.
9. DrosselB. // Phys. Rev. Lett. 2002. Vol. 23. Р. 89.
10. Gutenberg B., Richter C.F. // Bull. Seism. Soc. Am. 1952. Vol. 32. P. 163-191.
11. Aki K. // J. Geophys. Res. 1967. Vol. 72. P. 1217-1231.
12. Kanamori H., Anderson D.L. // Bull. Seism. Soc. Am. 1975. Vol. 65. P. 1073-1095.
13. Turcotte D.L. Fractals and Chaos in Geology and Geophysics. Cambridge, 1992.
14. Черепанцев А.С. Черепанцев С.Ф. // Изв. ТРТУ. 2006. № 12. С. 155-161.
15. Family F., Vicsek T. // J. Phys. 1985. A 18. L-75-L81.
16. Grassberger P. // Phys. Rev. 1994. E 49. Р. 2436.
17. Christensen K., Olami Z. // Phys. Rev. 1992. A 46. Р. 1829.
18. Corral A. et al. // Phys. Rev. 1995. Lett. 74. Р. 118.
19. Pepke S.L., Carlson J.M. // Phys. Rev. 1994. E 50. Р. 236.
20. Черепанцев А.С. // Изв. вузов. Сев.-Кавк. регион. Ес-теств. науки. 2008. № 1. С. 72 - 77.
Таганрогский технологический институт Южного федерального университета
27 июня 2007 г.
