Microsoft Word - MATH990901.doc Science and Technology, 5 (2000) 77-84 © 2000 Sultan Qaboos University 77 Simple Monte Carlo Cellular Models for Surface Evolution H. H. Al- Barwani and A.Purnama Department of Mathematics and Statistics, College of Science, P.O.Box 36, Sultan Qaboos University, Al-Khod 123, Muscat, Sultanate of Oman. باستخدام طريقة مونت كارلو للتطور السطحينموذج خلوي مبسط حمدي البرواني وأنتون برناما تستعمل محاكاة الخلية المعروفة بمونت كارلو في بعض المجسمات البسيطة ذات السطح النشوئي حيث يحاكي النمو :خالصة إلى موقع معين في السطح على عدد السطحي بإضافة خاليا جديدة إلى السطح، هذا و يعتمد اندماج الخلية الجديدة التي تصل تم . الخاليا الموجودة حول هذا الموقع وعلى الخاليا الجديدة المحتمل التصاقها بالمواقع التي بها أقل عدد من الخاليا المطوقة مو في اتجاه تمثيل ذلك بتطبيقات معينه كمحاكاة انتشار مقدمة اللهب ، والترسيب األرضي والنمو البلوري المتباين حيث يكون الن .معين أكثر انتشارا من النمو في اتجاهات أخرى ABSTRACT : Monte Carlo cellular simulations are described for some simple surface evolution models. Surface growth is simulated by adding new cells to the surface. The bonding of a new cell arriving to a site on the surface depends on the number of cells present around that site; new cells are more likely to stick at sites with the fewest missing surrounding cells. Applications are given for simulating the propagation of a flame front and the formation of surface landforms; and the (anisotropic) growth of a crystal, where the surface may grow more rapidly in one direction than others. KEYWORDS: Cellular Automata, Crystal Growth, Geomorphology. ost computational methods for solving equations of motion of a surface can be subjected to problems such as convergence and stability of the algorithm. Simulations can also show unphysical features. For example, non-linear wave models simulations have shown how edges (folded surfaces) and other surface discontinuities can develop as a result of the erosion of an initially smooth surface (Al-Barwani, 1996). Without a detailed knowledge of non-linear wave theory, computational cellular model simulations can give the same results of surface growth (Smith, 1991a). Cellular Automata are based upon discrete dynamical systems where the state variables are chosen in a finite set and distributed on a discrete grid. The time evolution is run synchronously in all cells of the grid and each cell changes its state according to a local rule that only depends upon its neighbouring cell states. Models are simple to construct but generally display very complex behaviour. Cellular Automata supply useful models in trying to understand and describe nature (Manneville et al 1989; Toffoli and Margolus, 1987). Simple cellular growth models have been demonstrated successfully in simulating the most realistic crystal growth performance (Webb et al 1997) and erosion of landforms (Smith, 1991a). In cellular models, the challenge is to choose the local relationship for representing the growth mechanism. This rule specifies how each cell is going to be influenced by some nearby cells. Simple statistical rules are usually preferred but these deterministic rules may not necessarily give the shapes or growth rates observed experimentally; as from the same initial conditions, the simulations eventually lead to the same evolution. In contrast, with a probabilistic rule, the same current situations may lead to several different outcomes, each with a given probability. The effects of surface growth can be simulated by adding new cells (particles) to the surface. A new diffusing cell arriving on the surface is allowed to move around the surface until it finds a stable AL-BARWANI AND PURNAMA 78 location to stick and become part of the surface. The probability of a new cell staying at a stable site depends on the local properties of that particular site. For example, in a surface growing from a two- dimensional square cell, an edge is a likely stable site; since a site that had two adjacent edges could well have a higher probability of capturing a new cell Simple Monte Carlo cellular models for simulating the effects of surface growth are constructed here; and simulations are performed to animate the propagation of a flame front; the formation of surface landforms; and the (anisotropic) growth of a crystal. The surface growth mechanism is specified so that the bonding of a new cell arriving to a site on the surface depends on the number of cells present around that site; that is, new cells are more likely to stick at sites with the fewest missing surrounding cells. For the case of flame propagation, an unburnt cell is more likely to become burnt if it has many surrounding burnt cells. In crystal growth, the effects of directional surface growth dependence are included in the probability for adding a new cell to the surface. Simple Cellular Automata Models A cellular automaton model consists of a two-dimensional regular uniform grid with discrete variable states assigned to at each cell. For example, to simulate the flame front propagation, each cell can be represented as a three-state variable with values of 0, 1 and 2 (Figure 1). The model evolves in discrete timesteps with every new cell state being uniquely updated from the current states of its close neighbour cells states according to a definite set of rules, which are local (i.e. no action-at-a-distance is permitted) and uniform (i.e. same rules are applied to all cells). These rules are commonly deterministic, and from the same initial conditions, the simulations should invariably lead to the same final evolution. Figure 1. A two-dimensional grid used for simulating the propagation of a flame front. A cell is a neighbour of another one if it has a chance to directly affect its state, i.e. via the rule, in one timestep. So typically, the neighbourhood is taken to be the cell itself and all immediately adjacent cells. For a cell (i, j), where iNi ,,1 L= and jNj ,,1 L= , its neighbourhood Ui, j may be defined as ( ){ }c±= ji,U ji, , (1) where c is a constant. Symmetry of the neighbourhood is assumed with periodic boundary conditions so that for example, cell (1 - c, 1) is defined as equal to (Ni - c, 1). Figure 2 shows the commonly used neighbourhoods; namely, the square ( ) ( ) ( ) ( ){ }j1i1ji,j1i1ji,U ji, ,,,,, +−−+= , (2) SIMPLE MONTE CARLO CELLULAR MODELS FOR SURFACE EVOLUTION 79 and the hexagonal ( ) ( ) ( ) ( ) ( ) ( ){ }.,.,,,,.,,.,,,,. 1j50ij1i1j50i1-j50ij1-i1j50-iU ji, +++−+−+= (3) Figure 2. Cell neighbourhood in a two-dimensional grid. (a) Square and (b) Hexagonal. It is straightforward to construct simple Cellular Automata models, which will simulate the effects of surface evolution. To demonstrate the algorithm, consider as an example, a surface growing from one seed cell as illustrated in Figure 3. Suppose that initially only the seed cell is present; and the surface growth mechanism is postulated by adding a cell to all free edges of any existing cell on the surface. This rule is uniformly applied to all cells at each timestep, and then repeats the surface growth starting with the new (current) surface. As it can be seen from Figure 3 that for this simple Cellular Automata model, the surface grows uniformly in four directions. By a simple change in the surface growth rule, growth can be arranged, say parallel to the coordinate axes, as shown in Figure 4. That is, at even timesteps, the previous growth mechanism is prescribed on all free edges (Figure 3); and at odd timesteps, the growth is limited only to cell which has two edges adjoining the current surface. The above surface growth rules are deterministic, as at each timestep the new surface is updated cyclically. The algorithms will simulate the effects of isotropic surface growth from a single seed cell, but the growth may not be observed in nature; as the simulations always eventually lead to the same new surface. In contrast, if probabilistic growth rules are employed, the current surface may grow to several different new surfaces, each with a given probability. Simple Monte Carlo Cellular Models Assuming a simple power law relationship, the probability for adding a new cell to the surface is calculated as, for the hexagonal neighbourhood (3), nb p       = 6 , (4) where b is the number of missing neighbouring cells (calculated using equation (3)) and n is the assumed power law dependence which determines the strength of the influence of surrounding cells. If n = 0, all surface cells have the same probability of growth. If 1=n there is a 6:1 ratio of success between a cell which is completely surrounded by other cells and one which has only a single cell as its neighbour; whereas if n = 2, the ratio is 36:1. As n increases, a smoother surface evolution will develop but simulations take longer to achieve. Equation (4) also implies that the choice of an initial surface is another important factor; as the local geometry of the surface depends on n. If the surface is initially flat and isotropic (in the sense of equal probability), it will retain its flatness and hence, develop no surface structure. A Monte Carlo method is used to speed up the simulation, so that at each timestep, instead of updating all surface cells sequentially, a random surface cell is chosen. The growth rate of the Monte Carlo algorithm is about twice faster than that of the Cellular Automata, as demonstrated by Webb et al. AL-BARWANI AND PURNAMA 80 (1997, figure 1e). After 40 timesteps, the new surface almost completely covers the old surface for a growth rule that adds a new cell with a 50% probability if there is at least one free adjacent edge. It is also straightforward to construct simple Monte Carlo cellular models to simulate the effects of surface evolution. The procedure is, instead of updating all cells, to pick a cell randomly from the surface and check if its state has to be updated according to the surface growth rule prescribed. A random number generator is used to select cells from the surface. Figure 3. Cellular automaton animation of a single seed cell growing on all free edges of its surface, after timesteps (a) 1, (b) 2, (c) 3 and (d) 4. To update the state of a chosen surface cell, first another random number generator selects a number R between 0 and 1. If R ≤ p, then the chosen cell state is changed; but if pR > , the cell state remains unchanged. The algorithm then continues to select another surface cell at random. Simple Monte Carlo cellular models for simulating the effects of surface evolution are constructed here; and simulations are performed to animate the propagation of a flame front; the formation of surface landforms, and the (anisotropic) growth of a crystal. The surface growth mechanism is specified so that a new cell is more likely to stick at randomly chosen surface cells with the fewest missing surrounding cells. For the case of flame propagation, an unburnt cell is more likely to become burnt if it has many surrounding burnt cells. In crystal growth, the effects of directional growth dependence are included in the probability for adding a new cell to the surface. A FLAME FRONT PROPAGATION MODEL: Flame propagation is a highly complex phenomena. Physically, it depends on the fuel supplied and the media on which the flame is spreading. Basically, it can be illustrated as follows: a fire lit at a particular site on the surface will only spread if that site is SIMPLE MONTE CARLO CELLULAR MODELS FOR SURFACE EVOLUTION 81 surrounded by other sites which are already burnt. In this way, the effects of flame front propagation can be simulated using a simple Monte Carlo cellular model. For simplicity, the flame is assumed to propagate only from burnt to unburnt region. Once a cell is burnt it remains burnt, and only interface cells (Figure 1) that can be burnt. A simple algorithm used can be briefly described as follows. Start initially with a square burnt region. An interface unburnt cell is chosen at random and the number of its nearest neighbours is calculated using the hexagonal neighbourhood (3). The flame propagation mechanism is postulated so that, an interface unburnt cell is more likely to become burnt if it is surrounded by many neighbouring burnt cells. This probability is calculated using equation (4) with n = 2 Figure 4. Cellular automaton animation of a single seed cell growing in directions parallel to the coordinate axes, after timesteps (a) 1, (b) 2, (c) 3 and (d) 4. The simulated flame front propagation is shown in Figure 5. Simulations showed initially, some rough flames fronts emerge from the square burnt region. After 8000 timesteps (Figure 5b), the initial square shape of the burnt region still visible, but the corners start to disappear. By 20000 timesteps (Figure 5c), the burnt region has begun to take a circular form and corners have completely disappeared. Finally after 50000 timesteps (Figure 5d), the burnt regions new shape has assumed an approximately circular formation. A LAND DEPOSITION MODEL: Because of the complexity of phenomena in nature, it is a difficult task, if not impossible, to decide and single out one physical process which is dominant in reshaping surface landforms. Surface landscapes are generally formed through geological processes of erosion and deposition by wind and water, and volcanic eruptions. Basically, particles eroded from one location are transported along the surface until they are deposited at another location. So land erosion and deposition processes can be illustrated as a result of randomly removing and adding cells to the surface; which as observed, depends mainly on the geometry of surface. In this way, the evolution of surface landforms is another example that can be simulated using a simple Monte Carlo cellular model. For illustration, consider land deposition on a concave surface represented initially by a cosine curve. Again, the deposition of (adding) a new cell at a site on the surface is assumed to be dependentonly on the number of its neighbouring cells. The surface cells, which have more missing neighbouring cells, AL-BARWANI AND PURNAMA 82 Figure 5. Monte Carlo cellular model simulation of flame front propagation from an initial square burnt region, using the hexagonal neighbourhood with n = 2; after timesteps (a) 1000, (b) 8000, (c) 20000 and (d) 50000. Simulations are (a) – (d) in the order left to right and top to bottom. are more likely to have deposition occurring around them. The hexagonal neighbourhood (3) is employed and the probability of adding a new cell to the surface is calculated using equation (4) with n = 2. Figure 6. Monte Carlo cellular model simulation of land deposition on the surface represented Initially by a cosine curve, using the hexagonal neighbourhood with n = 2; after time steps (a) 10, (b) 5000, (c) 28000 and (d) 50000. (a) (b) (c) (d) (b)(a) (c) (d) SIMPLE MONTE CARLO CELLULAR MODELS FOR SURFACE EVOLUTION 83 The simulated land deposition on a cosine curve is shown in Figure 6. Simulations show that the deposition of new cells is dependent on the geometry of the surface. After 28000 timesteps (Figure 6c), the initial smooth surface represented by a cosine curve has completely disappeared. The more concave the surface the more material is being deposited. This in turn fills the shallow part until it reaches a stable elevation and the process reaches the steady state. After 50000 timesteps (Figure 6d), the surface evolution reaches its stable position. The simple Monte Carlo cellular models described above simulate only the effects of isotropic surface evolution. The simplicity of the models means that it is easy to incorporate different physical effects. Smith (1991a) pointed out that these effects could be incorporated simply by modifying the probability (equation (4)) of adding a new cell to the surface. Next, a simple model is presented to include a directional growth dependence, where the surface is growing more rapidly in one direction. AN ANISTROPIC CRYSTAL GROWTH MODEL: Physically, a crystal is an ordered lattice of atoms in a solid. For all types of atom, the particular ordering and the extent of the ordering have important effects on the macroscopic behaviour of the crystal. A crystal growth mechanism is considered here by adding new diffusing cells to the surface. Although the model of crystallization may seem to be simple, the actual surface simulated may be extremely complex. In anisotropic growth, a surface may grow more rapidly in one direction than in others. This growth mechanism can be prescribed simply by modifying the rule by weighting the added particles into specific orientation. Figure 7. Monte Carlo cellular model simulation of anisotropic crystal growth from an initial circular seed region, using the hexagonal neighbourhood with n = 2; after timesteps (a) 1000, (b) 30000, (c) 50000 and (d) 80000. Simulations are (a) – (d) in the order left to right and top to bottom. To demonstrate the anisotropic surface evolution, consider as an example a crystal growing from a spherical seed region. As before, the surface growth mechanism is prescribed so that new cells are more likely to stick at a site on the surface if the site has the fewest missing surrounding cells. A site is chosen randomly, and the number of its surrounding cells is calculated using the hexagonal neighbourhood (3). But instead of equation (4), the probability of adding a new cell to the surface is now calculated by nb p       γ= 6 , (5) where γ includes a directional growth dependence. For simplicity, take n = 2 and define AL-BARWANI AND PURNAMA 84 yx yx + + = 22 γ , (6) so that it has a maximum 1 along the directions of coordinate axes and a minimum 2/1 where | x | = | y |. The simulated effects of crystal growth from a spherical seed region is shown in Figure 7, where eventually growth occurs more rapidly along the coordinate axes. After 30000 timesteps (Figure 7b) the initial spherical shape has disappeared. By 50000 timesteps (Figure 7c) the crystal has grown to a new shape of a square rotated through 45° with respect to the coordinate axes. After 80000 timesteps, the final new crystal shape is shown in Figure 7c. Concluding Remarks As it has been demonstrated that simple Monte Carlo model simulations give the most realistic surface evolution performance. Cellular Automata models, although simple to construct, by its very nature do not create the surface roughening of an edge which is caused by different rates of cell growth. Monte Carlo algorithms have also been found much faster than that of Cellular Automata (Webb et al 1997). In addition to providing modeling tools, cellular models employ a scheme that avoids most of the problems associated with other surface tracking schemes, since it does not require the measurement of the derivatives and it can quickly detect topological changes. The models only used statistical approach. No particular dynamical laws are invoked to control the motion of the particles; rather the statistics and rules are presumed to contain all the necessary information. Other physical effects, if required, can also be added to the models as have been demonstrated by Smith (1991a,b) and Toffoli and Margolus (1987). References AL-BARWANI, H.H. 1996. Propagation of Fronts with Gradient and Curvature Dependent Velocities. Ph.D. Thesis, Department of Mathematical Sciences, Loughborough University, U.K. MANNEVILLE, P, BOCCARA, N, VICHNIAC, G.Y. and BIDAUX, R. (Editors). 1989. Cellular Automata and Modeling of Complex Physical Systems. Springer-Verlag, Berlin. SMITH, R. 1991a. The application of cellular automata to the erosion of landforms. Earth Surface Processes and Landforms, 16: 273-281. SMITH, R. 1991b. Modeling and Simulation of Particle-Surface Interactions, Diamond and Diamond-like Films and Coatings. Plemon Press, New York. TOFFOLI, T. and MARGOLUS, N. 1987. Cellular Automata Machines. MIT Press, Massachusetts. WEBB, R.P., SMITH, R., AL-BARWANI, H.H. and WILSON, I.H. 1997. Simple Cellular models for growth. Radiation Effects and Defects in Solids, 141: 211-222. Received 1 September 1999 Accepted 24 February 2000