Self-Organized Criticality in an Anisotropic Earthquake Model
Li Bin-Quan, Wang Sheng-Jun
School of Physics and Information Technology, Shaanxi Normal University, Xi’an 710119, China

 

† Corresponding author. E-mail: wangshjun@snnu.edu.cn

Abstract
Abstract

We have made an extensive numerical study of a modified model proposed by Olami, Feder, and Christensen to describe earthquake behavior. Two situations were considered in this paper. One situation is that the energy of the unstable site is redistributed to its nearest neighbors randomly not averagely and keeps itself to zero. The other situation is that the energy of the unstable site is redistributed to its nearest neighbors randomly and keeps some energy for itself instead of reset to zero. Different boundary conditions were considered as well. By analyzing the distribution of earthquake sizes, we found that self-organized criticality can be excited only in the conservative case or the approximate conservative case in the above situations. Some evidence indicated that the critical exponent of both above situations and the original OFC model tend to the same result in the conservative case. The only difference is that the avalanche size in the original model is bigger. This result may be closer to the real world, after all, every crust plate size is different.

PACS: ;05.65.+b;;05.70.Ln;
1 Introduction

Self-organized criticality (SOC) is a key concept as a possible explanation for the widespread occurrence in many nature systems that long range correlations in space and time.[15] Earthquakes are probably the most relevant paradigm of SOC that can be observed by humans on earth. The relevance of SOC to earthquakes was first pointed out by Bak, Tang, and Wiesenfeld,[1] Sornette and Sornette.[6] According to this theory, plate tectonics provides energy input at a slow time scale into a spatially extended, dissipative system that can exhibit breakdown events via a chain reaction process of propagating instabilities in space and time. The empirical Gutenberg-Richter (GR) law[7] arises from the system of driven plates building up to a critical state with avalanches of all sizes. According to the GR law the distribution of earthquake events is scale free over many orders of magnitude in energy.

Then Olami, Feder and Christensen (OFC) introduced a nonconservative model on a lattice that displayed SOC.[8] The OFC model of earthquakes has played an important role in the context of SOC since 1992. However, the presence of criticality in the nonconservative version of the OFC model has been controversial since its introduction and it is still debated.[913]

OFC models on different topologies have been investigated in the papers, such as, the annealed random neighbor (ARN) graph model,[1417] the OFC model on a quenched random (QR) graph[18] and the effects of small-world and scale-free topologies on the criticality of the non-conservative OFC model.[1924]

Most work in this area is usually focused on their topological properties[2529] and homogeneous lattices with or without periodic boundary conditions. But the real systems modeled by these objects are not homogeneous. In a geological fault, for example, the local friction between the moving plates, which influences both the rate of motion and the redistribution on the neighbors in the OFC model, cannot be expected to be a constant value but should fluctuate according to local variations. Similarly, the local elasticity of the sheets, which determines how the energy is transferred from one point to another, is also expected to be variable. Therefore, a first step is to simply see how the introduction of quenched disorder in the simple coupled-map representation for these systems will affect their dynamical behavior. Some work has already been done along these lines.

A new earthquake model based on a random network was studied, on which the toppling mechanism of the system is that the force of the unstable site is redistributed to their nearest neighbors randomly.[3032] It is shown that when the system is conservative, the probability distribution displays power-law behavior. However, it displays no scaling behavior when the system is nonconservative. It is like the results in the ARN OFC model. But it is quite different to the model on quenched random graph[18] and the model on square lattice, both of which display criticality even the system is dissipated. It seems that the toppling mechanism of the system has affected the critical behavior of the system. They also compare the critical behavior of the model with different number of nearest neighbors. It is shown that different spatial topology does not alter the critical behavior of the system.

Mousseau studied the influence of quenched disorder on a coupled map model of earthquakes.[33] In his work, disorder is introduced in the redistribution fraction αi which now varies from site to site as αi = α+δi, where δi is a random number taken from a linear distribution [–δ, δ]. He said that the question of the role of disorder in dynamical systems is fundamental because most biological, neurological, or geological dynamical systems evolve in the presence of one or another type of disorder.

Janosi and Kertesz studied the effect of randomness in the threshold for the OFC rule.[34] They found that this type of disorder destroys criticality and changes the distribution of avalanche size from power law to exponential.

Ceva looked at uncorrelated and correlated disorder in the redistribution parameter α.[35] He was interested, however, in the effects of concentration of defects and not their amplitude. He found that SOC is stable under small concentration of defects.

In this work, we investigate the critical behavior on a modified anisotropic OFC model. Two situations are considered in this paper. One situation is that the energy of the unstable site is redistributed to its nearest neighbors randomly not averagely and keeps itself to zero. The other situation is that the energy of the unstable site is redistributed to its nearest neighbors randomly and keeps some energy for itself instead of reset to zero. Different boundary conditions are considered as well. The rest of the paper is organized as follows. In Sec. 2, we review the original OFC model and we point out the main reasons that have induced us to study the modify OFC model. In Sec. 3, we investigate the modified model and make comparison with other models by analyzing the distribution of earthquake sizes. Finally, in Sec. 4 we provide brief discussion and conclusion.

2 Model

The OFC model is defined on a two-dimensional square lattice of L × L sites. Each site is associated a real continuous energy Ei. To mimic that the system is driven continuously and uniformly, the value of Ei increases at the same rate. In simulations, find the largest value of energy Emax in the system and increase the energy of all sites by the same amount EthEmax. Therefore, the sites with the largest energy reaches the threshold value (EiEth) and becomes unstable. As soon as a site becomes unstable, i.e., Ei ≥ Eth, the global driving is stopped and the system evolves according to the following local relaxation rule

where “nn” stands for the collection of nearest neighbors to node i. In general, there are 4 of nearest neighbors, and only the isotropic situation is considered. The parameter α ∈ [0,1/4] controls the level of conservation of the dynamics, where α=1/4 corresponds to the conservative case, while α < 1/4 implies the model is nonconservative.

The toppling of one site triggers an avalanche, that is, neighbors of this site may become unstable and toppling propagates in the network. The avalanche is over until all the sites are below Eth. Then the driving to all sites recovers. The number of toppling sites during an earthquake is defined as the earthquake size S. Open boundary condition is used in OFC model.[36]

Here we modify the OFC model only in the toppling rule: when there is an unstable site toppling, the energy of the site is redistributed to its nearest neighbors randomly not averagely as follows.

If any EiEth then redistribute the energy on Ei to its neighbors randomly according to the following rule

where “nn” stands for the collection of nearest neighbors to node i, q is the number of nearest neighbors of every site. The anisotropic situation is considered in this modified OFC model. The size of parameter βnn is different for nearest neighbors. The parameter β ∈ [0,1] controls the level of conservation of the dynamics that is equivalent to 4α in the OFC model, where β = 1 corresponds to the conservative case, while β < 1 implies the model is nonconservative. The parameter βnn at each site is chosen randomly from a uniform distribution between 0 and β, and the sum is equal to β. Details as follow,

The other model that differs from above model only in the toppling rule: when there is an unstable site, the energy of the site is redistributed to its nearest neighbors randomly and keeps some energy for itself instead of reset to zero, that is not Eqs. (2), (3) but Eqs. (4), (5), as follow

where “nn” stands for the collection of nearest neighbors to node i, q is the number of nearest neighbors of every site. β00 is defined as follow

Table 1

The different of the two models.

.

To completely define the model, we need to consider the boundary conditions. We care about the open and periodic boundary conditions in our models.[3639] The energy αE of an unstable site at boundary is lost in the case of open boundary conditions. Periodic boundary conditions mean the energy is not lost, it is transferred to the other side, like a ring. In Table 1, we list the different of those models to be clearly understood.

In a system of SOC, the distribution of earthquake sizes is a power law function. The power-law exponent τ is defined as

In simulation, we will be interested in the distribution of avalanche sizes P(S).

3 Simulations and Results
3.1 The One Case of the Modified OFC Model

Now we study the one case of the modified OFC model (Case 1). The energy of the unstable site is redistributed to its nearest neighbors randomly not averagely and keeps itself to zero. In comparison to the original OFC model we plot the distribution of avalanche sizes in Fig. 1. The statistics are collected in the critical state for 109 non-zero avalanches for each system size.

Fig. 1 (Color online) Avalanche size distribution with open boundary conditions for (a) different the value of the parameter β with the system size N = 352. Different curves correspond to β = 0.40, 0.60, 0.80, 0.90, 0.95, and 1.0, from left to right. For comparison purposes, the original OFC model of α = 0.25 is shown in the direction of the arrow. (b) Distribution of earthquake size for the different system size with β = 1.0.

In Fig. 1(a), we show that the earthquakes size distribution P(S) for different values of the dissipation parameter β, in the network with N = 352 and with open boundary conditions. Different curves correspond to β=0.40, 0.60, 0.80, 0.90, 0.95, and 1.0, from left to right. We can see that the model transits from non-SOC to SOC behavior with the increase of the parameter β. Power-law fit is shown as red solid line, the slope of the straight line is τ = 1.286 58 for β = 1.0. It is like the critical exponent of the original OFC model. For comparison purposes, the original OFC model of α = 0.25 is shown in the direction of the arrow.

In OFC model, the SOC states exhibit that the distribution is a power law function with an exponential cutoff. The largest avalanche size in OFC model is about 7000. In the modified OFC model, the distribution of avalanche size depends on the dissipation parameter β. The largest avalanche size in modified OFC model is about 1000, which is close to the system size N.

It is not like the result of the original OFC model. The original OFC model exhibits SOC behavior for a wide range of α values and the exponent τ depends on α. However, we find that this modified OFC model exhibits power law distribution when the value of β tends to 1. It is similar with the model on RN[14,16] and small-world networks.[19,2324]

In Fig. 1(b), we show that the earthquakes size distribution P(S) for different size of the system N, different curves correspond to N = 152, 252, 352, and 502, from left to right. Size effect is present in Fig. 1(b). It is like the result of the original OFC model. The scaling of the cutoff in the energy distribution as a function of the system size for β = 1.0.

Now we plot the distribution of earthquake size for different values of β with periodic boundary conditions. As shown in Fig. 2(a), we noticed that the model transits from non-SOC to SOC behavior with the increase of the dissipation parameter β. Power-law fit is shown as red solid line, the slope of the straight line is τ = 1.397 15 for N = 352 and β = 0.98. The largest avalanche size is about 103 in the modified OFC model with β < 1.0. Although the system size is N = 352, the largest size of avalanche is very large with periodic boundary conditions and β = 1.0. The largest size of avalanche is 108 much larger than 103. There is not much difference in the results under different boundary conditions, except in β = 1.0. The result is only a slight difference in the critical exponent. As shown in Fig. 2(b), the simulation results of avalanche size distribution in systems with dissipation parameter β = 0.95 and network size N = 152, 252, 352, 502, and 1002, respectively. Although system sizes range from 152 to 1002, the change of the largest size of avalanche is very small.

Fig. 2 (Color online) (a) Simulation result for the probability density of having an earthquake of energy E as a function of E for a dissipation parameter β with N = 352 and with periodic boundary conditions. Different curves correspond to β = 0.40, 0.60, 0.80, 0.90, 0.95, and 0.98, from left to right. The fitted curve is shown as red solid line. (b) Distribution of earthquake size for the different system size with β = 0.95. Different curves correspond to N=152, 252, 352, 502, and 1002.
3.2 The Other Case of the Modified OFC Model

Next, we study the other case of the modified OFC model (Case 2). The energy of the unstable site is redistributed to its nearest neighbors randomly and keeps some energy for itself instead of reset to zero.

In Fig. 3(a), we show that the earthquakes size distribution P(S) for different values of the dissipation parameter β, in the network with N = 352 and with open boundary conditions. Different curves correspond to β = 0.40, 0.60, 0.80, 0.90, 0.95, and 1.0, from left to right. It is shown that the distribution of avalanche size depends on the dissipation parameter β. We can see that there is criticality only in the conservative case. Power-law fit is shown as red solid line, the slope of the straight line is τ = 1.319 19 for β = 1.0. It is like the critical exponent of the original OFC model. For comparison purposes, the original OFC model of α = 0.25 is shown in the direction of the arrow.

Fig. 3 (Color online) Avalanche size distribution with open boundary conditions for (a) different the value of the parameter β with the system size N = 352. Different curves correspond to β = 0.40, 0.60, 0.80, 0.90, 0.95, and 1.0, from left to right. For comparison purposes, the original OFC model of α = 0.25 is shown in the direction of the arrow. (b) Distribution of earthquake size for the different system size with β = 1.0. The fitted curve is shown as red solid line.

The result of this case and the original OFC model tend to the same in the conservative case. The only difference is that the avalanche size in the original model is bigger. This result may be closer to the actual situation, after all, every crust plate size is different.

In Fig. 3(b), we show that the earthquakes size distribution P(S) for different size of the system N, different curves correspond to N = 152, 252, 352 and 502, from left to right. Size effect is strikingly present here in Fig. 3(b). As L or β increases, the behavior slowly converges to a power law distribution of earthquake sizes P(s) ∼ sτ with an exponent τ = 1.319 19.

Now we plot the distribution of earthquake size for different values of β with periodic boundary conditions. As shown in Fig. 4(a), we noticed that the model transits from non-SOC to SOC behavior with the increase of the dissipation parameter β. Power-law fit is shown as red solid line, the slope of the straights line is τ = 1.427 78 for N = 352 and β = 0.98. Similarly, the largest size of avalanche is 108 much larger than normal 103 for β = 1.0.

Fig. 4 (Color online) (a) Simulation result for the probability density of having an earthquake of energy E as a function of E for a dissipation parameter β with N = 352 and with periodic boundary conditions. Different curves correspond to β = 0.40, 0.60, 0.80, 0.90, 0.95, and 0.98, from left to right. The fitted curve is shown as red solid line. (b) Distribution of earthquake size for the different system size with β = 0.95. Different curves correspond to N = 152, 252, 352, and 502.

As shown in Fig. 4(b), the simulation results of avalanche size distribution in systems with dissipation parameter β = 0.95 and network size N = 152, 252, 352, and 502, respectively. Although system sizes range from 152 to 502, the change of the largest size of avalanche is very small. Some impact has produced on the distribution of avalanche size for the different system size N. We can see that size effect is exhibited but not particularly strong in Fig. 4(b). We find that, as L or β increases, the behavior slowly converges to a power law distribution of earthquake sizes P(s) ∼ sτ with an exponent τ = 1.427 78.

4 Conclusions

In summary, we have made an extensive numerical study of a modified anisotropic model proposed by Olami, Feder, and Christensen to describe earthquake behavior. The toppling rule is different with that of the OFC model. Two situations were considered in this paper. One case is that the energy of the unstable site is redistributed to its nearest neighbors randomly not averagely and keeps itself to zero. The other case is that the energy of the unstable site is redistributed to its nearest neighbors randomly and keeps some energy for itself instead of reset to zero. Different boundary conditions were considered as well. By analyzing the distribution of earthquake sizes, we found that both above cases can exhibit self-organized criticality only in the conservative case or the approximate conservative case. Some evidence indicated that the critical exponent of both above situations and the original OFC model tend to the same result in the conservative case. The only difference is that the avalanche size in the original model is bigger. It is different from the result of original OFC model. The original OFC model exhibits SOC behavior for a wide range of α values and the exponent τ depend on α. This result may be closer to the real world, after all, every crust plate size is different.

Reference
[1] Bak P. Tang C. Wiesenfeld K. Phys. Rev. Lett. 59 1987 381
[2] Haldeman C. Beggs J. M. Phys. Rev. Lett. 94 2005 058101
[3] Wang S. J. Zhou C. New J. Phys. 14 2012 023005
[4] Plenz D. Schuster H. G. Criticality in Neural Systems Wiley New York 2014
[5] Wang S. J. Ouyang G. Guang J. et al. Phys. Rev. Lett. 116 2016 018101
[6] Sornette A. Sornette D. Europhys. Lett. 9 1989 197
[7] Gutenberg B. Richter C. F. Ann. Geofis. 9 1956 1
[8] Olami Z. Feder H. J. S. Christensen K. Phys. Rev. Lett. 68 1992 1244
[9] Klein W. Rundle J. Phys. Rev. Lett. 71 1993 1288
[10] Christensen K. Phys. Rev. Lett. 71 1993 1289
[11] Carvalho J. X. Prado C. P. C. Phys. Rev. Lett. 84 2000 4006
[12] Carvalho J. X. Prado C. P. C. Phys. Rev. Lett. 87 2001 039802
[13] Christensen K. Hamon D. Jensen H.J. Lise S. Phys. Rev. Lett. 87 2001 039801
[14] Lise S. Jensen H. J. Phys. Rev. Lett. 76 1996 2326
[15] Chabanol M. L. Hakim V. Phys. Rev. E 56 1997 R2343
[16] Broker H. M. Grassberger P. Phys. Rev. E 56 1997 3944
[17] Kinouchi O. Pinho S. T. R. Prado C. P. C. Phys. Rev. E 58 1998 3997
[18] Lise S. Paczuski M. Phys. Rev. Lett. 88 2002 228301
[19] Caruso F. Latora V. Pluchino A. et al. Eur. Phys. J. B 50 2006 243
[20] Caruso F. Latora V. Rapisarda A. Complexity, Metastability and Nonextensivity World Scientific Singapore 2005 355
[21] Masuda N. Miwa H. Konno N. Phys. Rev. E 71 2005 036108
[22] Rozenfeld A. F. Cohen R. ben Avraham D. Havlin S. Phys. Rev. Lett. 89 2002 218701
[23] Lin Min Zhao Xiao-Wei Chen Tian-Lun Commun. Theor. Phys. 41 2004 557
[24] Lin Min Wang Gang Chen Tian-Lun Commun. Theor. Phys. 46 2006 1011
[25] Rattana P. Berthouze L. Kiss I. Z. Phys. Rev. E 90 2014 052806
[26] Dominguez R. Tiampo K. Serino C. A. Klein W. Phys. Rev. E 87 2013 022809
[27] Markovic D. Gros C. Phys. Rep. 536 2014 41
[28] De Arcangelis L. Godano C. Grasso J. R. Lippiello E. Phys. Rep. 628 2016 1
[29] Perkins A. A. Galeano J. Pastor J. M. Phys. Rev. E 94 2016 052304
[30] Zhang Duan-Ming Sun Fan et al. Commun. Theor. Phys. 45 2006 293
[31] Zhang Duan-Ming Sun Fan et al. Commun. Theor. Phys. 46 2006 261
[32] Sun Fan Zhang Duan-Ming Commun. Theor. Phys. 50 2008 417
[33] Mousseau N. Phys. Rev. Lett. 77 1996 968
[34] Jánosi I. M. Kertész J. Physica (Amsterdam) 200A 1993 179
[35] Ceva H. Phys. Rev. E 52 1995 154
[36] Lise S. Paczuski M. Phys. Rev. E 63 2001 036111
[37] Socolar J. E. S. Grinstein G. Jayaprakash C. Phys. Rev. E 47 1993 2366
[38] Grassberger P. Phys. Rev. E 49 1994 2436
[39] Middleton A. A. Tang C. Phys. Rev. Lett. 74 1995 742