Elementary Cellular Automata along with delay sensitivity can model communal riot dynamics
aa r X i v : . [ n li n . C G ] J a n Elementary Cellular Automata along with delay sensitivity can model communal riotdynamics
Souvik Roy, ∗ Abhik Mukherjee, † and Sukanta Das ‡ (Dated: January 28, 2020)This work explores the potential of elementary cellular automata to model the dynamics of riot.Here, to model such dynamics, we introduce probabilistic loss of information and delay perturbationin the updating scheme of automata to capture sociological parameters - presence of anti-riot popu-lation and organizational presence of communal forces in the rioting society respectively. Moreover,delay has also been incorporated in the model to capture the non-local interaction of neighbours.Finally, the model is verified by a recent event of riot that occurred in Baduria of West Bengal,India. I. INTRODUCTION
Riots and their dynamics have been a popular topicto sociologists and historians [1–9]. In a parallel journey,computer scientists and mathematicians have found theirinterest in the study of riots, with a target to mathemat-ically model their dynamics [10–19]. The most popularapproach of developing such models is to adopt epidemi-ological framework [10–13, 17]. For example, 1960’s LosAngeles (1965) - Detroit (1967) - Washington D.C. (1968)riots [10], 2005 French riot [11], 2011 London riots [12]etc. have been modelled using this approach. Recently,non-local interactions along with neighbourhood depen-dency of elements of the system have been introducedin the models of riots [11, 13, 16, 20]. It is argued thatdue to the globalization and the advent of communica-tion technology, long range, that is, non-local commu-nications among elements are necessary to better modelthe dynamics of riots.In this scenario, we undertake this research to showthat the elementary cellular automata (ECAs) which relyonly on local neighbourhood dependency can efficientlymodel the dynamics of riots, if the neighbourhood de-pendency is delay sensitive . In particular, to model riotdynamics by ECAs, we introduce ‘probabilistic loss of in-formation perturbation’ and ‘delay perturbation’ in theupdating scheme of the automata. We observe that dueto this updating scheme, the ECAs show a new kind ofdynamical behaviour, which suggests us that some ECAscan be better model of riot. Finally, to validate our claim,we take into consideration a recent riot that happened inBaduria of West Bengal, India. Since media reports donot always reflect the ground realities of riots, we orga-nized an extensive field study in Baduria to get insight ∗ [email protected]; Department of Information Technology, In-dian Institute of Engineering Science and Technology, Shibpur,Howrah, West Bengal, India 711103. † [email protected]; Computer Science and Technology, In-dian Institute of Engineering Science and Technology, Shibpur,Howrah, West Bengal, India 711103. ‡ [email protected]; Department of Information Technology,Indian Institute of Engineering Science and Technology, Shibpur,Howrah, West Bengal, India 711103. about rioting dynamics.Here, in the proposed ECA based model, probabilisticloss of information perturbation rate is related to socio-logical factor such as the presence of anti-riot populationin the rioting society. Similarly, the presence of commu-nal elements in society, which plays a role to regeneraterioting spontaneity in the rioting society, indicates thephysical implication of delay in the system. However, theinherent property of CA is local interaction which contra-dicts the recent trends of considering non-locality in theage of globalization [11, 13, 16, 20]. The delay passivelyinduces a non-locality in the environment. To illustrate,the updated state information of a cell at time t reachesto its neighbouring cell at t + n time step where n depictsdelay for the cell and its neighbouring cell. This implies,non-local information from distance n reaches to the cor-responding neighbouring cell. The presence of communalorganization in society, which physically indicates delayin the system, induces this non-locality to regenerate therioting spontaneity. II. DELAY SENSITIVE CELLULARAUTOMATA
Here, we work with simple one-dimensional three-neighbouring two-state cellular automata, which arecommonly known as elementary cellular automata (ECA)[21]. The next state of each CA cell is determined as S t +1 i = f ( S ti − , S ti , S ti +1 ) where f is the next state function, S ti − , S ti , S ti +1 are the present states of left, self and rightneighbour of the i -th CA cell at time t , respectively. Thelocal transition function f : { , } → { , } can be ex-pressed as eight arguments of f . The decimal counterpartof eight next state is referred as ‘rule’. Each rule is asso-ciated with a ‘decimal code’ w , where w = f (0 , , + f (0 , , + · · · + f (1 , , , for the naming purpose.There are 2 (256) ECA rules, out of which 88 are mini-mal representative rules and the rest are their equivalent[22]. Classically, all the cells of a CA are updated simul-taneously. In the last decade, researchers have exploreddynamics of CA under asynchronous updating schemes[23–30].Classically, in ECA, delay and probabilistic loss of in-formation during information sharing among the neigh-bouring cells is not considered. In traditional cellularautomata, if a cell updates its state at time t , then thatstate information is available to neighbouring cell at time t + 1. To define the delay involved in sharing of infor-mation for two neighbouring cells i and j ( i = j ), weintroduce a non-negative integer function D( i , j ). In theproposed system, D( i , j ) = D( j , i ) > i and j . To illustrate, D( i , j ) = n inthe system implies, if cell i updates its state at time t ,then the updated state information is available to cell j at time t + n . In the proposed system, the delays arenon-uniform in space; i.e. D( i , j ) may be different fromD( i ′ , j ′ ), where i and j ; i ′ and j ′ are neighbouring cells,however, the delays are uniform in time. Practically, thedelay perturbation parameter d ∈ N assigns the maxi-mum possible delay in the proposed CA system. Everypair of neighbouring cells are randomly initialized withdelay between 1 to d following a uniform distribution.For the loss of information, one can consider that thedelay is ∞ (infinity). Here, ι (0 ≤ ι ≤
1) indicates theprobabilistic loss of information perturbation rate.Now, for introducing probabilistic loss of informationand delay in the system, each cell has to maintain stateinformation of neighbours to get a view of neighbour’sstate. In the proposed system, each cell has a view aboutthe states of its neighbours which may change from timeto time depending on the arrival of state informationabout neighbours. However, the cells act depending onthe current state information about neighbours at thattime. In this context, the state set is distinguished intotwo parts - the actualstate (self) of a cell and a vec-tor of neighbour’s viewstate. Now, the state set can bewritten as S ′ = S × S . Therefore, for a cell c , configu-ration at time t is distinguished into two parts - a tc and v tc where a tc ∈ S is the actualstate and v tc ∈ S is thevector of viewstate of left and right neighbours. Notethat, the actualstate set S is sufficient to represent tra-ditional CA. Here, in the proposed CA system, the localtransition function is also sub-divided into two parts - inthe first state update step, a cell changes its actualstate depending on the actualstate of the self and viewstates ofneighbours; and, in the second information sharing step,the cell shares its updated actualstate to its neighbouringcells. Now, the local transition function can be written as f ′ = f u ◦ f s , where, f u is the state update function, and f s is the information-sharing function. Here, the operator‘ ◦ ’ indicates that the functions are applied sequentiallyto represent the actual update.To illustrate, Fig 1 depicts a simple 3-cell ECA, where D ( i -1, i ) = D ( i +1, i -1) = 1 and D ( i , i +1) = 2. In Fig. 1,each cell has a view about the states of neighbours, i.e.left and right one for each cell. For every time step, thefirst step (dotted line) shows the state update functionand the second step (straight line) shows the informa-tion sharing function. Here, the information about statechange of cell i (resp. cell i +1) at first time step reachesto cell i +1 (resp. cell i ) at third time step due to de- i−1 i i+1tt+1t+2 FIG. 1. Example of delay and probabilistic loss of informationperturbation updating scheme. The applied rule is ECA 50. lay perturbation. In Fig 1, the information about statechange of cell i at first time step does not reach to cell i -1 at second time step due to probabilistic loss of infor-mation perturbation.To sum up, the proposed CA system depends on thefollowing two parameters : (i) The delay perturbationparameter d ∈ N indicates the maximum delay limit ofthe system; (ii) The probabilistic loss of information per-turbation rate ι (0 ≤ ι ≤
1) indicates the probabilisticloss of information during information sharing.
III. MODELING OF RIOTSA. Dynamic behaviour
To model the riot dynamics, we investigates the gen-erative behaviour of the proposed ECA system. Duringthis study, we start with the smallest possible seed pat-tern as initial configuration where a single cell is in state1, i.e. h· · · · · · i . From modelling and theoreti-cal research point of view, investigating dynamics of seedpatterns starting with single cell in state 1 is well estab-lish research approach [31–33]. From the physical impli-cation point of view, the initial seed with a single cell instate 1 represents the triggering event of riots.Here, we study the qualitative behaviour of the systemstarting from a single seed where we need to look at theevolution of the configuration, i.e. space-time diagrams,by inspection over a few time steps. Though this is nota formal method, but this approach can provide a goodcomparison. Note that, this generative behaviour studydoes not include 29 odd ECAs, out of 88 minimum rep-resentative ECAs, which have local transition function f (0 , ,
0) = 1. For odd ECAs, an empty background con-figuration, i.e. · · · · · · , evolves to · · · · · · which isunable to produce the generative behaviour of the sys-tem. Therefore, for even 59 rules, ECA depicts followingbehaviours - (i)
Evolution to zero : after one time step,the seed cell in state 1 has vanished; (ii)
Constant evolu-tion : the initial seed remains unchanged during the evo-lution of the system; (iii)
Left evolution : the seed shiftsor grows in the left side; (iv)
Right evolution : the seed
ECA 14 ECA 60 ECA 30 ECA 18FIG. 2. The samples of space time diagrams for the proposed updating schemes - (top) d = 1, ι = 0 .
0; (middle) d = 1, ι = 0 . d = 2, ι = 0 . shifts or grows in the right; (v) Growth behaviour : theseed cell develops into a pattern for both left and rightside. Fig. 2 depicts left evolution for ECA 14, right evo-lution for ECA 60 and growth behaviour for ECA 30.Table I shows the classification of ECAs depending onthe generative behaviour.
Evolution to Zero:
Constant evolution:
Left evolution:
Right evolution:
24 28 56 60 152 156 184
Growth behaviour:
18 22 26 30 50 54 58 90146 150 154 178TABLE I. Classification of ECA rules depending on the gen-erative behaviour.
Here, the target of this simple classification is to iden-tify the ECAs which develop into a pattern for both leftand right side, i.e. growth behaviour. Here, we makea sensible simple assumption that the riot propagationaffects every neighbour, i.e. both left and right for ECA.Therefore, 12 ECAs, out of 88, with growth behaviour are our target for modelling riot dynamics. In this context,note that, Redeker et. al. [31] classifies the behaviourof traditional synchronous CA starting from a singleseed into - ‘evolution to zero’,‘finite growth’,‘periodicpatterns’,‘Sierpi`nski patterns’ and ‘complex behaviour’which have no clear equivalence to Wolfram’s classes [34].‘Evolution to zero’ class only shows similarity with thisstudy. Here, the stable structure gets quickly destroyedin the presence of delay and probabilistic loss of informa-tion perturbation. As an evidence, in Fig. 2, fractal-likeSierpi`nski patterns [31] are destroyed for ECA 18 underthe proposed system. Therefore, now, the target of thisstudy is to identify candidate ECAs from 12 growth be-haviour ECAs to model the riot dynamics.
B. Candidate ECAs for modelling riots
Let us assume that the riot dynamics has two phases- spreading phase and diminishing phase. So, we choose4 candidate ECAs, out of 12, which show phase transi-tion under the proposed updating scheme, to model theriot dynamics (as an example, see ECA 18 in Fig 2).For these four ECAs 18, 26, 50 and 146, out of 88 min-imal representative rules, there exists a critical value of
FIG. 3. The samples of space time diagrams depicting phasetransition - (left) d = 1, ι = 0 .
3; (middle) d = 1, ι = 0 . d = 1, ι = 0 . probabilistic loss of information perturbation rate whichdistinguishes the behaviour of the system in two differ-ent phases - passive phase (i.e. the system converges to ahomogeneous fixed point of all 0’s) and active phase (i.e.the system oscillates around a non-zero density). As anexample for ECA 50, Fig. 3(left) depicts the active phasefor probabilistic loss of information perturbation rate 0 . ι = 0 . ι = 0 .
5. Asphysical implication with respect to rioting dynamics, thediminishing phase is not observed without the presenceof certain percentage of anti-riot population, i.e. prob-abilistic loss of information perturbation rate. However,presence of certain percentage of anti-riot population inthe society leads to passive phase in the diminishing riot-ing dynamics. According to [24],this phase transition be-longs to the directed percolation universality class. Notethat, in the literature of rioting dynamics research, theidea of critical threshold was also discussed in [13, 35] forunderstanding level of social tension to start a riot andsufficiently large number of protests to start a revolutionrespectively.
ECA ι cd =1 ι cd =2 ι cd =3 ι cd =4 . . . . . . . . . . . . . . . . Now, to understand the quantitative behaviour of these candidate rules, we let the system evolve through 2000time steps and average the density parameter value for100 time steps. Note that, for a configuration x ∈ S L ,the density can be defined as d ( x ) = x / | x | , where x is the number of 1’s in the actualstate for the configura-tion and | x | is the size of the configuration. Fig. 4 showsthe plot of the profile of density parameter starting froma single ‘1’ seed as a function of the probabilistic loss ofinformation perturbation rate with a fixed d parameterfor ECA rules which depicts phase transition behaviour.Table II depicts the critical value of probabilistic loss ofinformation perturbation rate for phase transition asso-ciated with these ECA rules where ι cd = k indicates thecritical value with d parameter value k . Note that, thecritical value for phase transition increases when the up-dating scheme is also associated with delay perturbation,see Table II for evidence. Moreover, the critical valueof probabilistic loss of information perturbation rate forphase transition proportionally increases with increasingvalue of delay. Table II justifies that the diminishingphase of riot needs more percentage of anti-riot popu-lation in the presence of sociological factor delay. Notethat, this phase transition result is not observed for onlydelay perturbation updating scheme. To sum up, ECAs18,26,50,146 are the final candidate rules for modelling ri-oting dynamics. Therefore, now, the target is to identifythe best candidate rule among those ECAs for valida-tion of Baduria riot dynamics. In this scenario, the nextsection depicts the case study comprising Baduria riot’sdataset. IV. BADURIA RIOT AND THE PROPOSEDSYSTEMA. Baduria riot dataset
Attracting nationwide media attention, Baduria riot isthe most well-exposed among recent Bengal’s riot events[36, 37]. The triggering event of Baduria riot took placeafter a social media religious post by a 17-year old studentin a village Baduria of West Bengal on 2 nd July, 2017.This social media post was seen as objectionable andwent viral in Baduria. Starting with this, violent clasheswere triggered between the two communities of Baduria,Basirhat, Swarupnagar, Deganga, Taki and Haroa sub-division.Here, we base our analysis here on reported incidentsin media reports during Baduria rioting time. The au-thenticity of media report data is cross verified with fieldstudy. We extract riot like events, as examples ‘attackon police vehicles’, ‘serious injuries’, ‘group clashed’ etc.,from 20 media reports [38–57] to build the data set. Inthe literature, the traditional methodology for quantify-ing the rioting activity is to study the daily crime reportsof police data for analysing the riot dynamics [11, 12].However, this methodology suffers due to lack of data onrioting event geographically located in third-world coun-
FIG. 4. The plot shows the profile of density parameter as a function of the probabilistic loss of information perturbation ratewith a fixed d parameter for ECA rules.
BaduriaDeganga Basirhat BasirhatTakiSwarupnagarHaroa Day1: 02/07/17 Day2: 03/07/17Day3: 04/07/17 Day4: 05/07/17 Day5: 06/07/17Day6: 07/07/17 Day7: 08/07/17 Day8: 09/07/17
FIG. 5. Graphical riot propagation dynamics of Baduria riot. try. Note that, arrest records, though available, do notindicate communal riot as explicit reason for the arrest.Here, we adopt two simple methodology for quantify-ing the rioting activity: Firstly, we define as a single eventany rioting-like act, as listed in the media reports afterambiguity checking, depending on its intensity. Thus,‘rail blockades at 3 places’ counts as 3 events, ‘threepolice vehicles have been torched’ indicates 3 riot likeevents. We thus get a dataset composed of number ofriot like events for every day from July 2 to July 9, 2017.Secondly, we quantify rioting activity from area (by sq.km.) affected by riot on a day to day basis from mediareports. Fig. 5 reflects the spatial propagation of riot dy-namics. It is not possible to quantify some important riotevents, like (number of) group clashes, road blockades,market/shop/school closure, from media reports. There- fore, we quantify those rioting events by area affected inriots day-wise. Figure 6 shows number of attack eventson (a) police (vehicles, stations, persons); (b) religiousplace, home, rail line and serious injuries; (c) affectedarea (in sq. km.) over the of time course.Now, we calculate the summation of percentage of in-tensity per rioting day, out of total intensity, for the ri-oting event datasets of Fig. 6(a),(b),(c). Hereafter, thepercentage of intensity per rioting day of this summa-rized intensity indicates the normalized overall intensity of Baduria riot which is reflected in Figure 6(d). Here,we work with this normalized overall intensity to under-stand rioting dynamics, which shows simple growth (up)and shrink (down) dynamics. Note that, this simple upand down dynamics, without any rebound, was also ob-served for 2005 French riots [11] and US ethnic riots [10].With a contradiction, the dynamics with up and suddendown (i.e. rise for four days and suddenly down on fifthday) was found in 2011 London riots [12, 13].On a different internal dynamics perspective, the prop-agation dynamics of Baduria riot depends on the localcommunication of violent events [58] and local rumourpropagation [59], [60]. However, during the riot event aparallel journey of religious harmony is also reported inmedia [39, 40] which finally converted the dynamics ofthe event as an early diminishing riot event. The fielddata also reflects this evidence [61], [62].In this context, anti-riot population of society does notparticipate in this rumour propagation. Moreover, theyplay an important role in the early diminishing dynamicsof the riots. Here, the term ‘anti-riot’ population is com-posed of following population - firstly, the secular pop-ulation of society [39, 40]; secondly, not ‘purely’ secularpopulation, however, due to the economical dependencythey play an anti-riot role during the riot time [63].
B. Verification
Now, this finding is mapped with the best candidaterule among ECA 18, 26, 50, 146 for modelling the Baduriarioting dynamics. To compare the CA dynamics andBaduria riot dynamics, here, we let the system evolvestarting from a single state ‘1’ seed and average the den-sity parameter value for every 100 time steps which de- (a) (b) (c) (d)FIG. 6. The plot shows quantified rioting activity for every day from July 2 to July 9, 2017 of ( a ) attack on police (vehi-cles,stations,persons); ( b ) attack on religious place, home, rail line and serious injuries; ( c ) affected area in riot (by sq. km.);( d ) normalized overall intensity. fines one single time unit ( ∴ ≈
15 minute).Therefore, the normalized density parameter and normal-ized overall intensity of Baduria riot are the parametersfor comparison in this study. Note that, normalized den-sity parameter is also calculated following the similarprocedure of calculating normalized overall intensity pa-rameter. Fig. 7 shows the profile of normalized densityparameter or normalized overall intensity of Baduria riotas a function of the time series. According to Fig. 7,ECAs 26 and 146 show similar dynamics with Baduriariot, however, ECAs 18 and 50 come up with ‘late’ con-vergence dynamics for critical probabilistic loss of infor-mation perturbation rate ( ι c ) where d = 1. We iden-tified the best fitting CA model using the standard rootmean square deviation = s P Tt =1 ( x ,t − x ,t ) T where x ,t and x ,t depict the time series normalized overall inten-sity of Baduria riots and density parameter value of CAmodel respectively. Here, ECA 26 ( d = 1 and ι c = . . . .
138 and 0 .
071 respec-tively. The best fitting CA with d = 1 depicts the evi-dence of no major presence of organized communal forcein the rioting society. Note that, the field data also re-flects this absence. Fig. 8 depicts the evidence associatedwith ‘late’ convergence with increasing value of delay per-turbation parameter for these ECA rules. Without pres-ence of organized communal force, the rioting society re-acts spontaneously and a simple up and down dynamicsis reflected. However, increasing rebound dynamics is ob-served for increasing value of delay perturbation param-eter, see Fig. 8 as evidence. As an insight, the rebounddynamics indicates that organized communal force playsrole for regenerating the spontaneity of rioting in the so-ciety.The proposed model only verifies rioting dynamics ofWest Bengal event. However, as a discussion, ECA 146depicts similar sudden-down dynamics of 2011 London ri-ots [12], for evidence see Fig. 8. Moreover, Berestycki et.al [13] have analysed sudden spike dynamics to represent FIG. 7. The plot compares normalized overall intensity ofBaduria riot and normalized density of ECA rules as a func-tion of time series. Here, d = 1 and ι = ι c . strong exogenous factor and slower but steady increaseto reflect endogenous factors in rioting dynamics. In thiscontext, ECA 18 and 50 shows similar sudden spike andsteady growth dynamics respectively, see Fig. 8 for evi-dence. However, proper understanding about this simi-lar signature behaviour of the proposed CA system andexogenous-endogenous factors is still an open questionfor us. Now, depending on the wide variety of results,the study strongly indicates that this model is relevantfor other rioting dynamics.To sumup, the study reflects the modelling aspects ofthe internal convergence rioting dynamics with respect tothe sociological factors - presence of anti-riot populationas well as organizational presence of communal forces.One may argue about the presence of other sociologicalfactors in the rioting dynamics, however, our target is topropose a simple model which can able to capture thecomplex rioting dynamics. To validate our argument, wequote from Burbeck et al. [10] which is the pioneeringwork on epidemiological mathematical riot modelling. “First efforts at model building in a new fieldnecessarily encounter the risk of oversimplifi-cation, yet if the models are not kept as simple FIG. 8. The plot shows the profile of density parameter as a function of time steps with changing delay parameter for ECArules. Here, the plot works with critical probabilistic loss of information perturbation rate ι c . as is practical they tend to become immune tofalsification.[10]” Moreover, this simple local interaction model can be ableto capture non-locality using delay perturbation. Here,the proposed CA system introduces non-uniformity withrespect to information sharing aspects, however, the sys-tem is associated with uniform rule. Several non-uniforminternal dynamics are also reflected in Baduria incident,i.e. here, three types of internal dynamics (first: loot-ing incident [64]; second: local resident vs refugees [65];third: local people vs police [66]) are reflected from thefield data. The current construction of CA model withuniform rule is unable to address this microscopic dynam- ics. Therefore, this research can be extended to explorethe dynamics of the CA model with non-uniform rule.In terms of dynamical systems, the model adopts timevarying system rather than linear time invariant systemin the journey towards full fledged non-linear dynamics.
Acknowledgements
This research is partially sup-ported by Impactful Policy Research in Social Science(IMPRESS) under Indian Council of Social Science Re-search, Govt. of India (P3271/220).The authors are grateful to Priyadarsini Sinha for herwork as a field data collecting team member. The authorsare also grateful to Prof. Subhasis Bandopadhyay forenlightening discussions on this topic. [1] M. I. Midlarsky, American Political Science Review ,996 (1978).[2] R. M. Govea and G. T. West, The Journal of ConflictResolution , 349 (1981).[3] J. Bohstedt and D. E. Williams, The Journal of Interdis-ciplinary History , 1 (1988).[4] A. Charlesworth, Rural History , 122 (1994).[5] D. Myers, Mobilization: An International Quarterly ,289 (2010).[6] D. J. Myers, American Sociological Review , 94 (1997).[7] S. Das, Social Scientist , 21 (1990).[8] S. Das, Communal Riots in Bengal, 1905-1947 (OxfordUniversity Press, 1993).[9] S. Das, Modern Asian Studies , 281 (2000).[10] S. L. Burbeck, W. J. Raine, and M. J. A. Stark, TheJournal of Mathematical Sociology , 1 (1978).[11] L. Bonnasse-Gahot, H. Berestycki, M.-A. Depuiset, M. B.Gordon, S. Roch, N. Rodriguez, and J.-P. Nadal, Scien-tific Reports (2018).[12] T. P. Davies, H. M. Fry, A. G. Wilson, and S. R. Bishop,Scientific Reports (2013).[13] H. Berestycki, J.-P. Nadal, and N. Rodiguez, Networks& Heterogeneous Media , 443 (2015).[14] B. L. Pitcher, R. L. Hamblin, and J. L. L. Miller, Amer-ican Sociological Review , 23 (1978).[15] D. Myers, American Journal of Sociology , 173(2000).[16] D. Braha, PLoS ONE (2012).[17] P. Baudains, S. D. Johnson, and A. M. Braithwaite,Applied Geography , 211 (2013). [18] M. Granovetter, American Journal of Sociology , 1420(1978).[19] M. J. A. Stark, W. J. Raine, S. L. Burbeck, and K. K.Davison, American Sociological Review , 865 (1974).[20] S. Gonzlez-Bailn, J. Borge-Holthoefer, A. Rivero, andY. Moreno, Scientific Reports (2011).[21] S. Wolfram, Theory and applications of cellular automata (World Scientific, Singapore, ISBN 9971-50-124-4 pbk,1986).[22] W. Li and N. Packard, Complex Systems , 281 (1990).[23] N. Fat`es, J. Cellular Automata , 387 (2014).[24] S. Roy, Physica A: Statistical Mechanics and its Appli-cations , 600 (2019).[25] B. Sethi, S. Roy, and S. Das, Complexity , 370 (2016).[26] B. Schnfisch and A. de Roos, Biosystems , 123 (1999).[27] H. J. Blok and B. Bergersen, Phys. Rev. E , 3876(1999).[28] S. M. Reia and O. Kinouchi, Phys. Rev. E , 042110(2015).[29] O. Bour´e, N. Fat`es, and V. Chevrier, Natural Computing , 553 (2012).[30] B. Sethi, S. Roy, and S. Das, Journal of Cellular Au-tomata , 527 (2018).[31] M. Redeker, A. Adamatzky, and G. J. Martinez, Interna-tional Journal of Modern Physics C , 1350010 (2013).[32] J. Gravner and D. Griffeath, Journal of StatisticalPhysics , 168 (2011).[33] J. Gravner and D. Griffeath, Theoretical Computer Sci-ence , 64 (2012).[34] S. Wolfram, Physica D: Nonlinear Phenomena , 1 (1984). [35] J. Lang and H. D. Sterck, Mathematical Social Sciences69