Analysis of pattern formation in systems with competing range interactions

We analyzed pattern formation and identified various morphologies in a system of particles interacting through a non-monotonic potential with a competing range interaction characterized by a repulsive core (r < rc) and an attractive tail (r > rc), using molecular-dynamics simulations. Depending on parameters, the interaction potential models the inter-particle interaction in various physical systems ranging from atoms, molecules and colloids to vortices in low κ type-II superconductors and in recently discovered ‘type-1.5’ superconductors. We constructed a ‘morphology diagram’ in the plane ‘critical radius rc–density n’ and proposed a new approach to characterizing the different types of patterns. Namely, we elaborated a set of quantitative criteria in order to identify the different pattern types, using the radial distribution function (RDF), the local density function and the occupation factor.


Introduction
In a variety of systems pattern formation is governed by competing interactions [1]. Examples of such systems are: Langmuir monolayers [2,3], colloids and gels [4][5][6][7][8], ferrofluids [9][10][11][12], magnetic garnet thin films, type-I superconductors [13], the pasta phase in neutron stars [14], etc. In general, there is a strong correlation between pattern formation and inter-particle interaction. Thus, attraction favors aggregation while repulsion favors low local densities. This competition between repulsive and attractive interactions leads to very rich morphologies and phases, such as stripes, clusters, bubbles, etc [15]. Those different patterns were observed in many diverse systems. In ferrofluid systems, rich phases due to the competition between dipolar forces and short-range forces opposing density variations were found experimentally [9] and theoretically [10]. In neutron stars, the competition between the short-range nuclear attraction and the long-range Coulomb repulsion leads to complex pasta phases [14].
Pattern formation was extensively studied in colloidal systems where the colloid-colloid interaction is characterized by the competition of the hardcore excluded volume interaction, on the one hand, and the polarization of the particles, on the other. The interaction potential can be further controlled by introducing other contributions to the interaction. As a result, rich configurations, such as clusters, repulsive or attractive glassy states and gels, were found numerically [16][17][18][19][20][21], analytically [22,23] and experimentally [5,7] in colloidal systems with short-range attractive interaction. The properties of isolated clusters formed by the short-range attractive and long-range repulsive interaction were studied recently [24]. By controlling the interaction, the growth of the cluster was shown to change from nearly spherical to onedimensional (1D) patterns. The 1D growth of the clusters facilitated the collective packing into columnar or lamellar phases [24]. The columnar and lamellar phases in three dimensions were analyzed using molecular dynamics (MD) simulations [25].
In superconductors, the vortex-vortex interaction is usually considered to be either repulsive (in type-II superconductors, where the Ginzburg-Landau parameter κ, i.e. the ratio of the magnetic field penetration depth λ to the coherence length ξ , is >1/ √ 2) or attractive 3 (in type-I superconductors, where κ < 1/ √ 2 and vortices are unstable), while vortices do not interact with each other at the so-called 'dual point' when κ = 1/ √ 2 [26]. However, a somewhat deeper analysis of the inter-vortex interaction in type-II superconductors near the dual point revealed an attractive tail [27,28]. This repulsive-attractive inter-vortex interaction was used for the interpretation of unusual patterns in the intermediate state in low-κ superconductors (e.g. Nb): islands of the Meissner phase surrounded by the vortex phase or vice versa, i.e. vortex clusters surrounded by the Meissner phase [27,28].
The recent discovery of 'type-1.5' superconductors [29] induced a new wave of interest in systems with non-monotonic interactions (see, e.g., [30,31]), due to the fact that the observed vortex patterns in the new superconductors revealed a clear signature of repulsive-attractive inter-vortex interaction. In particular, some typical characteristics of the observed vortex patterns (e.g. vortex lattice with voids and the nearest-neighbor distribution) were explained using a simple model involving a non-monotonic inter-vortex interaction based on a more general approach to multi-order-parameter condensates [32].
In this paper, we consider a model competing range interaction potential which is repulsive for the short range and attractive for the longer range. In principle, this type of interaction potential could be used as a model for various systems with non-monotonic inter-particle interaction, e.g. atoms or molecules (i.e. the Lennard-Jones potential), or vortices in two-band superconductors, depending on specific parameters of the potential. We investigate the pattern formation in 2D systems for various interaction potential profiles; for example, we distinguish 'soft-core' and 'hard-core' interactions and analyze the transitions between different types of patterns. Based on this analysis, we construct a 'morphology diagram' for various interaction parameters and particle densities. We propose a new approach to characterizing the different morphologies: instead of qualitative characterization of the patterns (e.g. clusters or labyrinths), we introduce a number of quantitative criteria to distinguish these. In particular, the obtained patterns are analyzed in terms of the radial distribution function (RDF) and additional quantities characterizing, e.g., the local density of particles in clusters.
The paper is organized as follows. In section 2, we describe the model. The pattern formation for different interaction parameters is discussed in section 3. In section 4, we analyze different types of patterns using the RDF and discuss criteria for the identification of different patterns. The conclusions are given in section 5.

The model
The model inter-particle interaction potential is chosen in the following form: Here, K 0 is the zeroth-order modified Bessel function, r i j = |r i − r j | is the inter-vortex distance, V 0 and λ are the units of energy and length, respectively. These units can be specified, depending on the system. For example, in the case of a type-II superconductor, the appropriate units of length λ and energy V 0 , correspondingly, are the magnetic field penetration depth λ and V 0 = 2 0 /8π 2 λ 2 , where 0 = hc/2e (see, e.g., [33]). In the dimensionless form, the interaction potential (1) reads as where the dimensionless length is defined as r i j = r i j /λ. Further on, we will omit the primes and use the dimensionless form of the potential (2). The interaction force is then given by where K 1 is the first-order modified Bessel function,r i j = (r i − r j )/r i j , and a and b are two positive coefficients. The first term of equation (3) is repulsive, while the second term describes an attractive interaction force. Indeed, for r → ∞, The interaction force (3) has a repulsive (attractive) tail when b < 1 (b > 1), while for r → 0, K 1 (x) → 1/x and thus F i j → (a/b − 1)/r . Therefore, for the short range the interaction force (3) is repulsive (attractive) when a > b (a < b), and we only consider the case a > b since an attractive interaction for short distances would result in a collapse of our system of point particles. When a > b and b < 1, the interaction is always repulsive, and particles form a Wigner crystal structure. The most interesting case is realized when a > b and b > 1, and the interaction has a repulsive core and attractive tail 1 . In this case, there exists a critical distance r c , where the inter-particle interaction energy (1) reaches a minimum (and the interaction force (3) changes sign). By setting the force equal to zero, the coefficient a is given by The pattern formation is determined by the coefficients b, r c and the particle density n.
We study pattern formation in a system of interacting point-like particles using the Langevin equations. The dynamics of a single particle i obeys the following overdamped equation of motion [16][17][18]33]: Here, η is the viscosity, which is set to unity, F i j is the inter-particle interaction force defined by equation (3) and F p i , given by is the pinning force (although often neglected, in the general case we introduced random pinning sites, see section 3). Here, N p is the number of pinning sites, f p is the maximum pinning force of each potential well, r p is the range of the pinning potential, is the Heaviside step function, and r The thermal stochastic term F T i in equation (5) obeys the following conditions: and We consider a 2D square simulation region L x × L y in the x y-plane and apply periodic boundary conditions in the xand y-directions. For the interaction force given by equation (3), which exponentially decays for large distances, we use a cut-off for r > 8. Note that in this way the interaction in our system is finite-range. The length of the square cell L x = L y = L has been varied from 60 to 180 to examine the finite-size effects, and we set L = 120 to optimize the calculations, without appreciable influence on the results. To obtain stable particle patterns, we performed simulated annealing simulations of interacting particles. For this purpose, particles were initially randomly distributed inside the simulation region at some suitable non-zero temperature. Then the temperature was gradually reduced to zero, and the simulation was continued until the total force acting on any single particle became much smaller than typical forces in the system.

Pattern formation
In this section, we study pattern formation and identify different morphologies depending on parameters of the interaction r c , b and the particle density n. In figure 1, we illustrate the change of the interaction force profile due to the increase of r c and b. The interaction force is very sensitive to r c , since its increase directly leads to an increase in the repulsive part in equation (3) (see the main panel of figure 1). On the other hand, an increase of b leads to an increase of both the repulsive and attractive parts (see the inset of figure 1). However, the repulsive interaction increases much faster for r < r c than the attractive interaction for r > r c .
Thus, increasing b facilitates stabilization of the minimum inter-particle distance at r c . In other words, an increase of b results in hardening of the core in the interaction force (3), i.e. the interaction changes continuously from a 'soft-core' to a 'hard-core' regime.

Soft-core interaction
In this subsection, we analyze pattern formation in the absence of additional pinning. The influence of random pinning on the pattern formation will be discussed in the next subsection.
To study pattern formation in the soft-core regime, we set the coefficient b = 1.1. We define the critical density as n c = 8r −2 c / √ 3, which is the density of an ideal (hexagonal) Wigner crystal with the lattice constant a = r c . The density is defined as n = N /S, where S = 120 × 120 is the area of the simulation region and N is the number of particles (further on, in our analysis of various patterns, we will refer either to the number of particles N in the simulation cell or to the density which is n = N /14 400). For the case of n < n c which is considered here, the Wigner crystal is not stable due to the attractive interaction. In figure 2, we present patterns formed by N = 1500 particles when r c increases from 1.9 (a) to 2.9 (f). Note that the condition n < n c is always fulfilled for all the values of r c in this range. We found that for r c < 2.1 particles form clusters similar to the formation of the 'clump' phase found in colloids and 2D electrons [1,[16][17][18] (see figures 2(a) and (b)). The main difference to the patterns found in [16][17][18] is that the relatively softer core in our case is compressed due to the attractive interaction, and the clusters acquire a circular shape. In addition, the interaction between the clusters (decaying exponentially for long distances within the interaction range) becomes negligible for intercluster distance of the order of a few to several r 0 . Therefore, the clusters can be considered as non-interacting for low densities (although they still do not approach each other), contrary to the situation described in [16,18] where a super-lattice is formed due to the long-range cluster-cluster repulsion. As r c increases, the clusters expand. In particular, for r c > 2.1, the clusters start elongating, which is an indication of the instability with respect to the transition to stripes. For 2.1 < r c < 2.3, a mixed state with both stripes and clusters is observed (see figure 2(c)). Further increase of r c gradually destroys the clusters and leads to the formation of labyrinths (see figures 2(d)-(f)).
In order to reveal the influence of the density on the pattern formation, we gradually increase the number of particles from 1500 to 10 500 in our simulation cell. First, in figure 3, we present patterns formed by N = 5500 particles for varying r c . Note that n > n c for N = 5500. As compared to the lower density case shown in figure 2, in the case of clusters, the additional particles lead to the expansion of the clusters (see figures 3(a) and (b)). A mixture of clusters and stripes is formed when the additional particles form bridges connecting the clusters (see figure 3(c)), which will be discussed in detail below. For labyrinths (shown in figures 2(d)-(f) for N = 1500), the additional particles fill the empty regions (voids), resulting in the formation of the hexagonal lattice with varying local density (figure 3(d)) and, finally, a regular hexagonal lattice (see figures 3(e) and (f)). Note that the lattice with varying local density (figure 3(d)) is not stable in systems with pure repulsive interaction such as vortices in type-II superconductors.
We analyzed in detail the intermediate regime (i.e. corresponding to the transition from clusters to stripes, see figure 2(c)) where r c ≈ 2.3. The resulting patterns for varying density are shown in figure 4. Figure 4(a) displays a configuration at low density with N = 1500 when many individual clusters are formed. With increasing the density (see figure 4(b)), the clusters connect to each other and form stripes. A further increase of the number of particles (figure 4(c)) results in the formation of a complex mixture of interconnected stripes with voids. A very interesting and counter-intuitive evolution is observed when the number of particles increases from N = 5500 to N = 6500 (see figure 4(d)): in contrast to the gradual transition from a void-rich configuration to a lattice-rich configuration when N changes from N = 3500 (b) to N = 5500 (c), in the case of N = 6500 we observe 're-entrant' behavior, i.e. the void-rich configuration starts recovering, which is compensated for by an increase of the local density in the stripes. However, further increasing the density results in the expansion of the stripes to the empty regions, which is accompanied by a decrease of the local density in the stripes, as shown in figure 4(e). The distribution of particles becomes more uniform, with only a few small voids. Finally, for N = 9500 (figure 4(f)), we obtain a deformed hexagonal lattice characterized by a varying local density with only one small void.
Our calculations show that the obtained patterns are very sensitive to variations in r c . Thus, if r c slightly decreases (e.g. r c = 2.25), the number of clusters greatly increases as compared to the case r c = 2.3, for the same density of particles. For r c = 2.25, we also observe the transition from a void-rich configuration to a lattice-rich configuration. However, since the decrease of r c increases the attractive component of the inter-particle interaction this occurs at a much higher density. For even smaller values of r c , i.e. r c < 2.1, we do not observe a lattice pattern, even for an extremely large number of particles (up to N = 20 000).

Analysis of the stability of the patterns
To further study the stability of the patterns, here we analyze the influence of the cut-off procedure and the annealing time used in our simulations. Alternatively, we consider the effect of additional weak pinning and modify the inter-particle interaction by introducing a weak repulsive tail.
First, we consider the case when clusters are formed (i.e. r c < 2.3). We choose r c = 1.9 and increase the cut-off distance from 8 to 16λ. We also increase the annealing time from typically (10 000-20 000) t to 200 000 t, where t = 0.04 is the time step of the MD simulation. In figure 5, we plot the evolution of the energy with time, and we show some representative snapshots in the insets. We start from a random configuration of particles (see figure 5(b)). During the annealing, the energy first drops rapidly due to the clusterization of the particles. After a relatively short time, small clusters start to merge and form larger clusters (see figure 5(c)). The average inter-cluster distance increases when the total number of clusters decreases. The merging becomes more difficult and considerably slows down when the clusters are well separated. As a result, short-living metastable states are formed (see figure 5(d)). However, during the merging of the clusters, the lifetime of the metastable states increases very rapidly, and gradually the short-living metastable states evolve into long-living metastable states when the clusters are well separated (see figure 5(e)). Comparing with the relatively faster annealing case (i.e. 20 000 t), we see that the slower annealing process results in a decrease of the number of clusters inside the simulation cell, which increases the lifetime of the metastable states. However, even for the shorter annealing times (i.e. 20 000 t), we arrive at long-living metastable configurations which can have either slightly higher energy or even the same energy as those states obtained during the long-time annealing (note that the energy of the configuration does not change over the time interval from ≈15 000 t to 50 000 t). It is worth noting that the morphology is similar to that found in the 'fast' annealing process. Thermodynamically stable states can be obtained in the presence of weak random pinning. To verify this, we introduce weak random pinning sites in the system, with the maximum pinning force f p = 0.2 and radius r p = 0.3 and the density of random pinning sites being the same as the density of the particles. We found that such a weak pinning produces no appreciable impact on stripes, labyrinths or lattices. In the case of clusters, the bonding energy due to the inter-particle interaction within every individual cluster is strong, and the weak pinning is not sufficient to destroy the structure of the clusters. However, the inter-cluster interaction is weak, therefore, those clusters can be easily stabilized by weak pinning. Thus, comparing to the patterns obtained without pinning, we find that slightly smaller clusters are formed in the presence of pinning, and the merging of clusters is prevented even for two relatively close but small neighbor clusters (see figure 6). Note that in the absence of pinning, clusters have nearly circular shape while their shape is less regular in the case of weak pinning. In the latter case, some of the particles were pinned individually, i.e. outside the clusters, see figure 6(b).
An alternative way to stabilize clusters (and approach the thermodynamically stable configurations) consists in using a modified force with an additional weak repulsive tail, e.g. taken in the form: F i j = F i j cosh(r e − r ), where r e is a parameter chosen as r e = 10, and the cut-off of the interaction is set at 16λ. Then the resulting modified force coincides with that given by equation (3) for r < 8 but the weak attractive tail smoothly evolves to repulsive when r > 10. The modified interaction force leads to very similar patterns as the original one (given by equation (3)) because the weak repulsive tail is not efficient for the clusters to form superstructures (see, e.g., [16,18,21]), but is sufficient for stabilizing the configurations of particles. When we increase the inter-cluster repulsive interaction by increasing the particle density, the resulting patterns are similar to the previously calculated long-living metastable states. The reason for this behavior is that the radius of a cluster is insensitive to the particle numbers within the cluster. As a result, the inter-cluster repulsion increases very slowly as the density increases.

Morphology diagram
Based on the above analysis of the patterns, we constructed a 'morphology diagram' (as we explained above, the calculated patterns are long-living metastable states, i.e. in the absence of additional pinning or the long-range repulsive tail-rather than real thermodynamic phases; thus the diagram presents various morphologies) in the plane of r c and the number density n (see figure 7). For extremely low density (n < 0.03), particles form small clusters that are well separated, and the patterns are rather insensitive to variations in r c . However, morphologies become richer when density n increases (i.e. n > 0.1). Low values of r c (i.e. r c < 2.1) still favor cluster formation for a broad range of densities, although the particle density in clusters increases for large n. With increasing r c , clusters become unstable with respect to the formation of the 'bridges' between separate clusters, which is a precursor of the formation of stripes. The 'morphology diagram' in the plane 'critical radius r c -density n' (n) and representative patterns (a-m). For r c < 2.1 and extremely low density (n < 0.03), particles form small, well-separated clusters. For larger r c , the size of the clusters increases (see the change from (d) to (e) and from (b) to (c)). Increase of the density results in increasing size of the clusters (see the change from (d) to (b) and from (e) to (c)). Further increase of r c leads to elongation of the clusters and the formation of 'bridges' between them (a). Thus, the configurations change gradually from clusters to stripes (i). At larger r c or n, stripes interconnect and form patterns with voids (j). Thus the stripes are divided into two sub-groups, I and II. When r c becomes even larger, labyrinth-like configurations are formed for n < n c (see (g) and (f)). With further increase of r c or n, deformed hexagonal lattices with (k) or without voids (l) are formed. The deformation of a hexagonal lattice is reduced for even larger values of r c and n (m).
Stripes are formed in a rather narrow range of r c when r c > 2.1 (see figure 7). The stripe patterns can be divided into two sub-groups, I and II, i.e. void-rich (see figure 7(i)) and void-poor (latticerich) stripes ( figure 7(j)). For larger r c and n < n c , particles form labyrinth structures. However, when increasing the density, additional particles fill the empty regions and finally they form a deformed hexagonal lattice with varying density. Depending on r c and n < n c , deformed lattice is characterized by a varying local density or by the appearance of voids. Correspondingly, we distinguish two regions (1 and 2 in figure 7). Note that in the vicinity of the morphology boundaries, patterns are always mixtures of the two morphologies (e.g. clusters and stripes) except the transition from the stripes to the deformed hexagonal lattice. Near this morphology boundary, particles form either stripes with high local density or deformed hexagonal lattices with lower local density. The probability of finding the lattice configuration greatly increases when r c or n increases. Finally, for very high density, when the average inter-particle distance becomes smaller than r c the repulsive interaction prevails, resulting in the formation of an almost perfect hexagonal lattice.
The stability of clusters for varying parameter of the inter-particle interaction, r c , can be examined by applying an additional weak confinement potential and calculating the energy corresponding to the varying number of particles in the system. For this purpose, we use a parabolic potential, similar to F i p in equation (6), to confine all the particles and gradually decrease f p to zero during the annealing. This allows us to obtain one large cluster for various r c . Then we study the stability of these induced clusters, by analyzing the average energy per particle U for various r c and N c , where N c is the total number of particles composing the cluster. We found that for r c < 2.3, the clusters are stable: U is approximately proportional to N 2 c and decreases while N c increases (see figure 8(a)). The size of the cluster grows with the number of constituent particles: it goes to infinity for N c → ∞, which corresponds to phase separation. Note that this result could be understood from a comparison of the model potential (equation (2)) with the Lennard-Jones potential: for small r c , both the potentials become very similar, in particular, for r c ≈ 1. On the other hand, it is well known that the Lennard-Jones potential leads to phase separation.
For r c = 2.3, which is a typical value for stripe formation, there is a saddle point in the function U versus N c , which shows that the clusters have an optimal size, e.g. N c ≈ 280 (see figure 8(b)). Finally, for r c > 2.3, U changes sign when N c > 7, resulting in a positive energy, indicating that the clusters in this case are overheated by the parabolic potential and therefore are not stable. Labyrinths are more favorable in this case.
The patterns presented in figure 7 are found in many real systems. For example, clusters are found in such systems as colloids [4,6,19,20] and 2D electrons [16][17][18]. The obtained stripe patters are very similar to those in Langmuir monolayers [3]. Furthermore, several of the calculated morphologies were found in superconductors, e.g. the stripe phase in the mixed state of type-I superconductors, clusters of the Meissner phase or vortex clusters in the intermediate state in low-κ type-II superconductors [27,28], the labyrinth phase (i.e. vortex lattice with voids) or vortex clusters in type-1.5 superconductors [29]. In spite of the variety of physical systems and length scales ranging from nano-and micro-objects to cosmic objects, the common feature of all these systems is a competing attractive-repulsive inter-particle interaction, which allows analyzing them within the same approach.

Hard-core interaction
Let us now analyze the influence of the b parameter on the pattern formation. As mentioned above (see figure 1), an increase in parameter b changes the potential from soft-core to hardcore. Hardening the repulsive core in the interaction potential generally leads to a decreasing compressibility of the inner parts of the patterns where particles are closely packed. Clearly, in this case the patterns change as compared to the soft-core regime.
As shown in figure 9, for b = 4 and low density (N = 1500 per simulation region), all the patterns are clusters of different shapes. Thus for r c = 1.15 and r c = 1.3, the clusters are of a circular shape similar to that found in the soft-core regime, although for r c = 1.3 some clusters, composed of a small number of particles, have a polygonal shape. With increasing r c , polygonal-shaped clusters become more favorable, which allows us to identify them as 'short stripes'. For even larger r c , the clusters become much larger. They represent separate islands of the hexagonal lattice. For higher density (N = 6500), the variety of patterns is much richer. Although for r c = 1.15 circular-shaped clusters are still observed (see figure 10(a)), which grow in size for r c = 1.3 and slightly deviate in shape from circular ( figure 10(b)), for larger r c = 1.45, the deviations of the cluster shape from circular become more pronounced (see figure 10(c)). With further increasing r c , i.e. r c = 1.6, stripes start to form (see figure 10(d)) followed by lattices with voids for r c = 1.75 and r c = 1.9. In the hard-core regime, deformations of lattices occur via the appearance of voids instead of varying local density (as happens in the soft-core regime).
The obtained patterns presented in figures 9 and 10 resemble the nanoscale domain patterns formed by Cu and Pb atoms on a surface of Cu [34,35]. In particular, the lattice islands (figures 9(d)-(f)), the stripes (figure 10(d)) and void-rich configurations (figure 10(e)) are similar to the patterns evolving from islands to stripes and 'inverted' islands formed by Cu and Pb atoms arranged in a single atomic layer on a Cu(111) surface [34,35], with increasing atomic coverage.

Radial distribution function
Although the different patterns and morphologies studied above are qualitatively well distinguished, it is highly desirable to build a set of solid criteria that would allow us to identify them in a quantitative manner. For this purpose, we analyze here different morphologies by means of the RDF. The RDF, g(r ), describes the variation of the atomic (particle) density as a function of the distance from one particular atom (particle). If we define an average density as n = N /V (where V is the volume or surface area in the 2D case), then the local density at distance r is ng(r ). The knowledge of RDF is important since one can measure g(r ) experimentally using digital video microscopy [36]. Moreover, macroscopic thermodynamic quantities can be calculated using g(r ) in thermodynamics [37].
In our calculations, we define the RDF g i (r ) as follows: Here, the lower index indicates that the RDF centers at the position of the ith particle; N is the number of particles whose distance to the ith particles is between r and r + r . The average g(r ) is given by In figure 11(a), we plot the function g(r ) calculated for the low-density cluster configuration shown in figure 2(b). The function g(r ) has two well-pronounced peaks. The first peak corresponds to the average distance to the first coordination sphere (nearest neighbors), while the second one is located at a distance approximately twice the distance to the first peak, which shows short-range periodicity. For r larger than the size of the cluster and smaller than the inter-cluster distance, g(r ) < 1, since only a few particles are located inside this range. The decreasing tail in the RDF shows strong aggregation of the particles forming clusters. The minimum of g(r ) is very close to zero, indicating that most of the clusters have circular symmetry and are well separated. The position of the minimum of g(r ) (marked by the gray arrow in figure 11(a)) gives an estimate of the average diameter of the clusters.
For stripes (see figure 11(b) for the pattern shown in figure 2(d)), there are also two peaks indicating the short-range periodicity, similar to clusters. However, unlike the case of clusters, the minimum of g(r ) is non-zero since the stripes are generally continuous. For labyrinths (see figure 11(c) for the pattern shown in figure 2(f)), only short-range periodic ordering exists, and g(r ) becomes uniform for larger r .
In figure 12, we plot g(r ) for the patterns formed at high density (N = 5500). For clusters, increase in the density leads to an increasing variation of the cluster size. Note, however, that the sizes of individual clusters are rather insensitive to moderate variations of the number of particles in the clusters, due to the strong compressibility of the core. As a result, the function g(r ) for the high-density clusters misses the short-range periodicity and the second peak disappears (see figure 12(a)). Strikingly, this peak still exists for the void-rich patterns (see figure 12(b)). The decrease of the tail of g(r ) becomes very slow, which shows a minor aggregation of the particles.
For lattices (see figure 12(b)), the position of the first peak is r 1 ≈ a = 2/( √ 3n), where a is the distance between two neighboring particles in the ideal hexagonal lattice with the density n. The positions of the second and third peaks are at r 2 ≈ √ 3a and at r 2 ≈ 2a, respectively, corresponding to those in a hexagonal lattice. The fourth, fifth and even sixth peaks appear, which show that this structure is much more ordered. However, due to the variation of the local density, these peaks are inhomogeneously broadened: they are actually a combination of many neighboring peaks. These peaks become more pronounced for larger r c or n, which implies that the lattices become more regular.
It is interesting to compare the results for the function g(r ) for hard-core interaction with those for soft-core interaction. Although the clusters shown in figures 10(a) and (b) have shapes very similar to those of the clusters shown in figures 3(a) and (b), the analysis using the RDF shows that most of the clusters shown in figure 10(b) have hexagonal ordered cores. The peaks in g(r ) appear to be much better separated than, e.g., for the deformed hexagonal lattice formed in the soft-core case (see figure 12(c)). For the stripe and void-rich patterns, the RDF shows a much clearer hexagonal lattice ordering.
Therefore, the analysis using the RDF allowed us to establish quantitative criteria for the different patterns and reveal the differences in the structure of the patterns (which look similar, e.g. clusters) in the case of soft-and hard-core inter-particle interactions ( figure 13).

Local density
Using the RDF, let us define another useful quantity, the 'order parameter', to characterize various types of patterns. By definition, the local density is Here, N ξ ≈ π ξ 2 n is the average number of particles within the circle centered at the ith particle with radius R = ξ . In an ideal hexagonal lattice, one particle has six nearest neighbors, therefore N ξ = 7 and ξ = N ξ /π n = √ 7/π n. Thus, for any configuration characterized by a small local density fluctuation, the average local density I = I i = 6. From the definition given by equation (11), we can see that the presence of a large fraction of empty regions can result in a considerable increase of I since only a 'shell' of those regions of thickness ξ is considered. In figure 14(a), we plot the function I for different r c and N for the soft-core interaction with b = 1.1. We see that the function I can be used to characterize the different types of patterns. Thus we found that, for clusters, I is always large (I > 20) since aggregation is strong. For stripes, aggregation is weaker, and thus I is smaller (10 < I < 17). For labyrinths, the regions free of particles are relatively small. Therefore, I ranges between 6 and 9. Finally, Figure 14. The average local density I versus the particle density n for b = 1.1 (soft-core) (a). Curves A, B, C, D, E and F correspond to r c = 1.9, 2.1, 2.3, 2.5, 2.7 and 2.9, respectively. For clusters, I > 20, for stripes, 10 < I < 20, for labyrinths, 6 < I < 10; I = 6 for lattice and deformed lattice configurations. The local density I versus the particle density n for b = 4 (hardcore) (b). Curves A, B, C, D, E and F correspond to r c = 1.15, 1.3, 1.45, 1.6, 1.75 and 1.9, respectively. for a hexagonal lattice, I is always 6. Therefore, the function I serves as a measure of the aggregation of particles and allows us to effectively distinguish the different types of patterns.
The function I also provides a tool for analyzing the stability of patterns with increasing density. This is demonstrated in figure 14(b), where we plot I for the hard-core case when b = 4. For r c = 1.15, I > 20 for all the values of the density, and thus the clusters are stable. However, the situation changes for r c = 1.3: while clusters are formed for low densities up to N ≈ 7500, for larger density I decreases below the value I = 20, which means that clusters elongate and interconnect, giving rise to the onset of stripe formation. For r c = 1.45, the particles start to form stripes at even lower density (N = 3500). For higher density, the additional particles fill in the empty regions and finally they form a regular lattice. For r c > 1.6, the island-like clusters formed at low densities are very unstable with respect to an increase in density. The additional particles rapidly occupy the empty regions, and the patterns change from clusters to stripes and from stripes to lattice with increasing density. Therefore, the phenomenological description of pattern evolution with increasing density we revealed above in section 3.1 has been verified in a quantitative manner using the solid basis of the local density function approach.

Occupation factor
We demonstrated that the (average) RDF g(r ) (equation (10)) and the local density function I (equation (11)) allowed us to unambiguously characterize various types of patterns. At the same time, it is still hard to distinguish, using the above tools, between labyrinths (or a lattice with voids) and lattices, especially in the soft-core regime when the corresponding RDF do not differ much from each other, and both types of patterns are characterized by small values of I . However, these can be easily distinguished by employing another additional criterion, namely, the occupation factor which characterizes the fraction of the space occupied by particles to the total space.
Let us assume that all the patterns are formed by islands of hexagonal lattice with the average distance between two nearest-neighbor particles r = r 1 , where r 1 is the position of the first peak of the RDF. Then the ratio of the area occupied by the particles and the whole simulation area is A = (r 1 /a) 2 . As shown in figure 15(a), this ratio can be used to distinguish between labyrinths and lattices, since for lattices A ≈ 1.
We also found that circular clusters in the case of soft-core interaction are stable with respect to an increase in the density. In this case, the additional particles cannot efficiently increase the occupation factor due to the increase in the local density of the core. However, for stripes and labyrinths, the occupation factor A increases with the density (see figure 15(a)). These patterns are not stable for high density. Thus the large jump in A(N ) at N = 6500 in curve C in figure 15(a) shows the transition from stripes to lattice.
As one can expect, in the case of hard-core interaction, the occupation factor increases nearly linearly versus the density (see figure 15(b)). Therefore, the polygon shaped and islandlike clusters are not stable: with increasing density, they evolve to the lattice configuration (the plateau in A(N ) in figure 15(b)).
In addition, let us introduce another useful quantity which characterizes the degree of 'perfection' of a lattice. Let us define the particles with inter-particle distance shorter than ξ as neighboring particles. Then is the average displacement of the inter-distance between two neighbor particles (measured in units of a), which is independent of the density. The function ε is non-zero for a deformed Table 1. The set of criteria used to quantitatively identify various morphologies in terms of the RDF g(r ), the local density function I , the occupation factor A and the parameter ε (for details, see the text).

Lattice
Several peaks: r 1 , r 2 , . . . ≈6 ≈1 0 Deformed Several peaks: r 1 , r 2 , . . . ≈6 ≈1 >0 lattice hexagonal lattice and zero for the ideal hexagonal lattice. Thus, ε quantitatively measures the degree of perfection of a lattice. Note that ε is only used for lattices as an auxiliary tool to distinguish ordered lattices from less ordered ones (see table 1).

Conclusions
Using MD simulations, we analyzed the pattern formation and identified different morphologies in a system of particles interacting via a non-monotonic competing range potential, with a repulsive short-range part and an attractive long-range part. The form of the interacting potential is rather general: it describes, depending on the specific parameters, the inter-particle interaction in a variety of physical systems ranging from, e.g., atoms and molecules (Lennard-Jones potential) to colloids and vortices in superconductors. In particular, it can be used as a model of inter-vortex interaction in low κ type-II superconductors and in the recently discovered so-called 'type-1.5' superconductors.
Since the non-monotonic range potential involves interaction frustration, the obtained patterns are in general long-living metastable states rather than equilibrium configurations. We demonstrated that thermodynamically stable patterns can be obtained if we either use a modified potential with an additional repulsive tail or introduce weak pinning.
The obtained morphologies were summarized in a 'morphology diagram' in the plane 'critical radius r c -density n' (r c is the critical radius where the interaction force changes its sign). We also analyzed the influence of the hardness of the 'core', i.e. the strength of the repulsive core part of the interaction potential, on the pattern formation.
We developed a set of criteria in order to unambiguously identify the obtained types of patterns using the following approaches: (i) the RDF g(r ), (ii) the local density function I and (iii) the occupation factor A. In addition, we introduced a parameter which characterizes the degree of perfection of a lattice ε. Employing these approaches, we elaborated a set of criteria for the identification of the different morphologies summarized in table 1.