Profile of a Faceted Macrostep Caused by Anomalous Surface Tension

The height profile of amacrostep on a vicinal surface near equilibrium is studied numerically using a restricted solid-on-solidmodel with a point-contact-type step-step attraction (p-RSOSmodel).We calculate the surface tension of vicinal surfaces around the (001) surface inclined towards the ⟨111⟩ direction using the density-matrix-renormalization group method. We also calculate the height profiles of vicinal surfaces using the Monte Carlo method and study the connection between the height profile of the macrostep near equilibrium and the discontinuous surface tension. We find that the height profile of a macrostep on a vicinal surface near equilibrium can be classified depending on the zone in the faceting diagram where the system exists. We also find finite size effects both for the height profile and for the inhibition of the macrostep motion in the relaxation process to the equilibrium state.


Introduction
For 4H-SiC crystals grown by a solution method, macrostep formation is considered to be a major reason for grown crystal degradation [1]. The profile of a macrostep of 4H-SiC crystal formed on a vicinal surface is shown in Figure 1. To date, the formation of faceted macrosteps near equilibrium has been considered to be due to the effects of anomalous surface tension. 4H-SiC can be used as an example for understanding the importance of controlling the dynamics of macrosteps on a vicinal surface near equilibrium. To accomplish control of the dynamics of macrosteps, we have to clarify the fundamentals of the phenomenon step by step.
The instabilities with respect to macrostep formation near equilibrium have been considered phenomenologically to be the result of anomalous surface tension [2,3]. However, to date, such anomalous surface tension has not been calculated using a microscopic model. Since a crystal surface is a lowdimensional substance, it is difficult to obtain reliable results theoretically due to the large contribution from thermal fluctuations [4].
In previous papers [5][6][7][8][9][10], we calculated the surface tension numerically using a restricted solid-on-solid (RSOS) model with a point-contact-type step-step attraction (p-RSOS model, Figure 2). Here, "restricted" means that the height difference between nearest neighbor sites is restricted to {0, ±1}. We consider that the origin of the point-contacttype step-step attraction is the orbital overlap of the dangling bonds at the meeting point of the neighboring steps. The energy gained by forming the bonding state is regarded as the attractive energy between steps.
To obtain reliable results, we used the density-matrixrenormalization group (DMRG) method [11][12][13][14][15][16], and we showed numerically that the slope dependence of the surface tension is discontinuous (note that, in mathematical terms, "disconnected" is a more accurate description than "discontinuous") at low temperatures at equilibrium. This discontinuous surface tension affects the morphology of the crystal [5][6][7][8][9][10]. The equilibrium crystal shape (ECS), which is the crystal droplet shape with the least surface free energy, can be obtained from a polar graph of the surface tension using the Wulff theorem [17][18][19][20] or obtained using the Landau-Andreev method [21][22][23][24][25]. The discontinuous surface tension causes a jump on the surface slope at the facet edge on the ECS. This jump of the surface slope is referred to as a first-order shape transition, after the fashion of a phase transition [26].
The discontinuous surface tension also leads to the formation of a faceted macrostep on a vicinal surface, as two surfaces coexist [8] at equilibrium. Since the side surface  of a faceted macrostep is smooth, the kink density of the surface becomes small. Hence, under a small driving force for crystal growth or sublimation (dissolution), the velocity of the faceted macrostep is significantly smaller than the velocity of a single step [9,10]. In particular, at sufficiently low temperature, the macrostep moves intermittently. This inhomogeneous motion of macrosteps degrades the grown crystal.
In the present work, we study how the height profile of the macrostep is classified with the profile linked to the discontinuous surface tension. In a previous paper [5], we obtained the faceting diagram for the p-RSOS model using the DMRG method. The faceting diagram is shown in Figure 3, illustrating the various zones that the model can be divided into (see [5] for details). The macrostep is formed in the step-faceting zone and in the step droplet zone. We consider three cases in detail: a typical example from the step-faceting zone ( int / = −0.9 and / = 0.6), a typical example from step droplet zone I ( int / = −0.9 and / = 0.63), and a typical example from Gruber-Mullins-Pokrovsky-Talapov (GMPT) zone I [27,28], where and int are explained in the following section.
This paper is organized as follows. In Section 2, we explain the p-RSOS model. The DMRG calculated surface tensions Advances in Condensed Matter Physics Step-faceting Step  [5]. This figure is taken from [5]. are shown in Section 3. In Section 4, height profiles at equilibrium calculated by the Monte Carlo method for the step-faceting zone, step droplet zone, and GMPT universal zone are studied. In Section 5, we study the finite size effects both on the height profile and on the inhibition of the macrostep motion in the relaxation process to the equilibrium state. Discussions are given in Section 6 and conclusions are given in Section 7.

p-RSOS Model
The Hamiltonian of the (001) surface can be written as where N is the total number of lattice points, surf is the surface energy per unit cell on the planar (001) surface, is the microscopic ledge energy, ( , ) is the Kronecker delta, and int is the microscopic step-step interaction energy. The summation with respect to ( , ) is taken over all sites on the square lattice. The RSOS condition is required implicitly.
When int is negative, the step-step interaction becomes attractive (sticky steps).

Anomalous Surface Tension: DMRG Calculation
The surface tension is the surface free energy per normal area. In order to evaluate the surface free energy for a vicinal surface, we add the terms concerning the Andreev field [23,24]: = ( , ). The Hamiltonian for the grand canonical ensemble with respect to the number of steps is [29,30] The Andreev field behaves like a chemical potential with respect to a single step. From statistical mechanics, the grand partition function Z is calculated as follows: . The summation with respect to {ℎ( , )} is taken over all possible values of ℎ( , ). The Andreev free energỹ( ) [23,24] is the thermodynamic grand potential calculated from the grand partition function Z, as follows [29,30] where N is the number of lattice points on the square lattice.
In low-dimensional cases, more precision is required than that provided by a mean field calculation of the partition function [4]. Hence, to obtain reliable results, we adopt the DMRG method [11][12][13][14][15][16], which was developed for onedimensional quantum spin systems. We use the transfer matrix version of the DMRG method, which is known as the product wave-function renormalization group (PWFRG) method [14][15][16].
Examples of the Andreev free energy calculated by the DMRG method are shown in Figure 4 as a function of ( / , / ) for = . It should be noted that the profile of the Andreev free energỹ( , ) is similar to the profile of ECS = ( , ), wherẽ( , ) = ( , ), ( , ) = − ( , ), and represents the Lagrange multiplier related to the crystal volume [23][24][25].
Using the inverse Legendre transform with respect tõ ( ), we obtain the surface free energy per unit -area, (p). The surface tension is calculated from the surface free energy, as follows: The polar graphs of the surface tension for several typical temperatures are shown in Figure 5. In the figure, the angle is the tilt angle from the (001) surface towards the ⟨111⟩ direction. The surface gradient p is related to by |p| = ± tan .
A characteristic of the surface tension is a discontinuity near the (001) surface and the (111) surface. There are two transition temperatures ,1 and ,2 [8]. For < ,1 , the surface tension is discontinuous around the (111) surface ( Figure 5(b)). For < ,2 ( Figure 5(a)), the surface tension is discontinuous around the (001) surface. In the case of int / = −0.9, ,2 / = 0.613 ± 0.02 and ,1 / = 0.709 ± 0.02 [5]. As shown in Figure 3, we call the temperature range < ,2 the step-faceting zone, the range ,2 ≤ < ,1 the step droplet zone, and the range ,1 ≤ the GMPT universal zone. Furthermore, the step droplet zone and the GMPT zone are divided into two zones by the line of the roughening transition temperature of the (001) surface.
The ECS calculated as the Andreev free energy (shown by red lines) agrees well with the ECS obtained by the Wulff theorem [17][18][19][20]. In the step-faceting zone ( Figure 5(a)), the (001) surfaces directly contact the (111) surface without a curved area. In the step droplet zone ( Figure 5(b)), the {001} surfaces contact the curved areas without a jump of the slopes (point P in Figure 4), whereas the {111} surfaces contact the curved areas with a jump of the slopes (point Q in Figure 4). In the GMPT zone ( Figure 5(c)), both the (001) surfaces and (111) surfaces contact the curved areas without a jump of the slopes.

Profile of Macrosteps
To see the morphological effect of the discontinuous surface tension on the height profile of macrosteps, we studied the vicinal surface of the Hamiltonian equation (1) with a fixed number of steps step using the Monte Carlo method.
In this method, at the initial time, the steps were positioned at equal distances. Then, a lattice site to be updated is chosen randomly. The surface structure is updated nonconservatively by the Monte Carlo method with the Metropolis algorithm. With the RSOS restriction taken into consideration, the structure is updated with probability 1 when Δ ≤ 0 and with probability exp[−Δ / ] when Δ > 0, where Δ = − for the energy of the present configuration and the energy of the updated configuration. The energy is calculated based on the Hamiltonian equation (1). In the direction along the steps, we required a periodic boundary condition. For the direction normal to the steps, we connected the lowest side of the structure to the uppermost side of the square by adding a height with a number of steps step . The obtained snapshots of the surfaces and the height profiles of the vicinal surfaces at the bottom line of the top-downview of the surface after 2 × 10 8 Monte Carlo steps per site (MCS/site) are shown in Figure 6.
For the vicinal surface in the step-faceting zone (Figure 6(a)), there is one macrostep with an orientation of the side surface of (111). The structure of a regular train of steps is thermodynamically unstable for the direction between ⟨001⟩ and ⟨111⟩, except for the directions very close to ⟨001⟩ (see Figure 5(a)). In this manner, the steps selfassemble to form a faceted macrostep on the vicinal surface at equilibrium. Furthermore, the motion of the macrostep is almost inhibited. This is because the kinks in the vicinal surface are significantly small relative to the vicinal surface of the original RSOS model. It should be noted that the surface stiffness 2 (p)/ for the (001) surface and for the (111) is infinite since both the (001) and (111) surfaces are smooth.
For the vicinal surface in GMPT zone I (Figure 6(c)), the steps appear to be dissociating and the height profile is almost linear. With a small driving force for growth or for sublimation, the mean surface height increases or decreases linearly. However, the growth or sublimation rate is smaller than that of the original RSOS model because there remain locally merged steps (step droplets).
The structure of the vicinal surface in step droplet zone I ( Figure 6(b)) is complicated. The surface slope 1 at the point Q in Figure 4 is the key slope. For a low step density ( < 1 ), the steps are dissociated. No macrosteps are formed on the vicinal surface. For a high step density ( ≥ 1 ), the structure of a regular train of steps is thermodynamically unstable due to the discontinuity in the surface tension ( Figure 5(b)). At that location, the vicinal surface has two surfaces, as shown at point Q in Figure 4. At point Q, the (111) surface and the surface with slope 1 have the same value of the Andreev free energy. Therefore, the (111) surface and the surface with slope 1 coexist on the vicinal surface [8]. We calculated the averaged height profiles ℎ(̃) for the vicinal surface using the Monte Carlo method (Figure 7). Here, ℎ(̃) is defined by where is the number of lattice points along the step direction,̃is the location of the lattice points along the step Advances in Condensed Matter Physics direction, and̃= √ 2 + 2 is the location of steps in the direction normal to the step edge.
Comparing Figure 1 on 4H-SiC with the side views in Figure 6, the height profile in Figure 1(b) looks similar to the height profile in the step droplet zone, that is, Figure 6(b) or Figure 7. As for the case of Figure 1(a), it looks similar to Figure 6(a). According to [1], however, the terrace surface in Figure 1(a) is also slightly tilted. This similarity suggests that the surface tension of 4H-SiC around the side surface is discontinuous, whereas the surface tension around the terrace surface is continuous.

Finite Size Effect
Since the faceted step moves very slowly at lower temperatures in the step-faceting zone, it takes an extremely long time to reach global equilibrium. The time evolution of the averaged surface height ⟨ℎ( )⟩ at temperatures / = 0.2 and / = 0.4 is shown in Figure 8, where ⟨ℎ( )⟩ = (1/ ) ∑ {̃} ℎ(̃) and is the number of lattice points normal to the mean running direction of the steps. In Figure 9, we show snapshots of the vicinal surface after 2 × 10 8 MCS/site. Since attachments and detachments on the crystal surface are balanced, local For the temperature / = 0.2, the surface structure changed during the initial few thousands of MCS/site but did not change for the rest of the period 2 × 10 8 MCS/site (Figure 8(a)). The surface froze due to the faceted macrosteps. The steps moved more frequently for / = 0.4 than for / = 0.2. From Figure 8(b), we can see that a step that detaches from the edge of the side surface of a macrostep due to thermal fluctuations can be readily captured by a neighboring macrostep. In this manner, the structure of the surface changes intermittently and takes far more time than 2 × 10 8 MCS/site to reach global equilibrium.
Two questions arise at this point: whether the system always reaches global equilibrium for / = 0.6 and why the edge of the macrostep looks smeared in Figure 6 when it should be sharp according to the surface tension of Figure 5(a). In the following, we present the Monte Carlo results, which demonstrate that local equilibrium is always achieved, while global equilibrium is not reached for a large system size, and that the edge of the macrostep converges to a sharp edge in the large size limit.
In order to study the size dependence of the height profile of macrosteps, we calculated the averaged height profile ℎ(̃) of vicinal surfaces. We show ℎ(̃) in Figure 10. As can be seen from the figures, even a time of 2 × 10 8 MCS/site is not sufficient to reach global equilibrium for sizes of 160 √ 2 × 160 √ 2 and 240 √ 2 × 240 √ 2.
However, the edges of the macrosteps become sharper as the size increases. If we express the width of the curved area by Δ(̃), Δ(̃) can be estimated by the following inequality relationship: 2 V ( ,̃)| =1 < Δ 2 (̃) < 2 (̃), where 2 V ( ,̃) is the fluctuation width of a single step on the vicinal surface with step density [31,32] and 2 (̃) is the fluctuation width of a single step on the (001) planar surface [33]. Then, we obtain the inequality expression: wherẽ( ) is the step stiffness defined bỹ( ) = ( ) + 2 ( )/ 2 , ( ) is the step tension, and is the in-plane tilt angle relative to the [010] direction. Hence, Δ(̃)/̃→ 0 in the thermodynamic limit̃→ ∞. This result is consistent with a discontinuous surface tension ( Figure 5). From the thermodynamics of two-surface coexistence [8], the surface slope jumps at the cross sectioned line of the two surfaces on the vicinal surface. Hence, in the thermodynamic limit, the curved area near the edge of the macrosteps converges to zero. Thus, we can conclude that the curved areas at the edge of a macrostep appear to be due to the finite size effect at equilibrium.
In the step droplet zone, global equilibrium is achieved more quickly. Though the faceted (111) side surface is difficult to move, the tilted surface (surface slope = 1 ) which contacts the (111) surface has many kinks. The thermal fluctuations on the tilted surface promote the achievement of global equilibrium.

Discussions
It would be helpful if real numbers for surface tension and temperatures were known for 4H-SiC. However, in order to determine these, we need to know the values of and int , and to determine and int , we need real values for the kink energies and temperature dependence of the step stiffness of a single step for the interface of 4H-SiC. To date, these real values are not accurately known. Therefore, more experimental observations on and int are required.
Fortunately, for step-faceting the atomic composition of the substance is less relevant because it is a first-order phase transition. From thermodynamics, the free energy is expressed as = − , where is an internal energy and is entropy. At = 0, = , and depends on the atomic composition of the substance. In the phase transition ( > 0), however, entropy plays a dominant role. The entropy is hardly affected by the atomic species of the substance but depends on the symmetry, the interacting force range, the number of components, and the dimensionality (as seen in Trouton's rule for molar entropy of vaporization, e.g., the molar entropy of vaporization for nonpolar materials is almost the same value for various kinds of liquids).
In the present work, we have shown that the characteristics of the height profile of a macrostep can be classified by the zone in the faceting diagram where the vicinal surface exists. In the step-faceting zone, the terrace for a macrostep is (001), whereas, in the step droplet zone, the surface in contact with a macrostep has the slope 1 ̸ = 0. This classification is a direct result of the connectivity of the surface tension. Since the p-RSOS model is constructed on a slightly coarse grained vicinal surface [8] with a square mesh, this result is applicable to a wide variety of materials to qualitatively explain the effect of a discontinuous surface tension. From the classification scheme, we can find a method of dissociation of macrosteps for 4H-SiC at equilibrium. As mentioned at the end of Section 4, 4H-SiC in Figure 1 is thought to be in the step droplet zone. To dissociate the macrostep, the temperature should be increased so that the state of the system moves to the GMPT zone. When the stepstep attraction is strong (a large | int |), the area of the step droplet zone is large. Hence, we have to raise the temperature substantially to reach the GMPT zone. In that case, the introduction of a small amount of "surfactant" like Al in 4H-SiC (Figure 1(b)) may make it easier to dissociate the macrostep [34][35][36], more so than increasing the temperature.
Comparing the height profiles of macrosteps observed for 4H-SiC (Figure 1) with the height profiles simulated by the p-RSOS model in detail, we find that the shape of the profile at the upper edge of a macrostep differs slightly from that at the lower edge of a macrostep for a real system, whereas, for the p-RSOS model at equilibrium, the shapes of the height profiles at the upper and lower edges are the same. Since the height profile of 4H-SiC (Figure 1) is obtained during the crystal growth, some nonequilibrium effects may apply. In the case of vapor growth, surface diffusion and the Ehrlich-Schwoebel barrier effect [37][38][39][40] are known to cause step bunching or step meandering.
It is important to distinguish between the relaxation process from an initial configuration to the equilibrium state and the steady-state motion of a vicinal surface with a driving force for crystal growth or crystal dissolution. For the relaxation to the equilibrium state near equilibrium on the vicinal surface, the faceted step is not able to move due to the divergence of the surface stiffness [9,10,41]. We can see this from the time-dependent Ginzburg-Landau (TDGL) theory for the surface (interface) motion [42,43]. Applying the finite size effect, the results in Section 5 are consistent with the expectation from the TDGL equations. Macrosteps with a small number of steps are easier to move than macrosteps with a larger number of steps. In real observations near equilibrium, the same phenomenon will be encountered. The dynamics of the height profile of a macrostep for steady-state motion will be studied elsewhere.
Furthermore, in a real system there exists a step-step repulsion originating from the elastic interaction among the steps [44][45][46], which is expressed by the equation ∑ ̸ = ∑̃0/(̃(̃)−̃(̃)) 2 , where 0 is the coupling constant of the elastic repulsion,̃is the location of the th step on a vicinal surface normal to the mean running direction of steps, and̃is the location along the steps.
The elastic step-step repulsion causes three effects: (1) an increase of the (111) surface energy, (2) a weakening of the step-step attraction [5], and (3) the development of a regular array of macrosteps [47]. Regarding (1), the step-step repulsion increases the surface energy of the (111) surface, which increases the area of the (001) surface on the ECS while decreasing the area of the (111) surface. Hence, wider (001) terraces are preferred on the vicinal surface. For (2), if we consider an averaged field for interactions among the steps, the step-step repulsion reduces the step-step attraction overall. Then, the step-faceting zone shrinks and the zone boundary between the step-faceting zone and the step droplet zone shifts to a lower temperature [5]. For (3), due to the difference in the force range between the step-step attraction (point-contact type) and the step-step repulsion (long range), the balance between the accumulated interactions depends on the distance between steps. For short distances between steps, the attractive interaction dominates, whereas for long distances between steps, the repulsive interaction dominates. Hence, a regular array of -merged steps develops [47], where is selected as the least free energy of the system. We studied the height profile of macrosteps for 4H-SiC because faceted macrosteps for 4H-SiC can be explained by the anomaly of the surface (interface) tension [1]. No other explanation for the formation of the faceted macrostep has yet been found. For 6H-SiC or 3C-SiC the classification of the height profile of a macrostep in the present work is applicable Advances in Condensed Matter Physics 9 to explain the effect of the surface tension qualitatively. If faceted macrosteps in the mesoscopic scale for 6H-SiC or 3C-SiC are observed at equilibrium, the surface tension is suggested to be discontinuous. When the slope of a terrace 1 is zero (a smooth flat surface), the vicinal surface is in the step-faceting zone, whereas when the slope of a terrace is not zero, the vicinal surface is in the step droplet zone.
Finally, we note that, formally, the Andreev free energy equals the negative pressure − [5]. We can see this from the analogy between (3) and the grand potential of a gas-liquid system. This pressure represents the tendency of the surface area to become wider. On the mesoscopic scale, the surfaces on a crystal determine their locations so that the pressure of neighboring planar surfaces balances out. Since this pressure is a kind of force, it propagates rapidly around the crystal surface.

Conclusions
The ECS and surface tension of the vicinal surface for the p-RSOS model with int / = −0.9 were calculated ( Figures  4 and 5). To ensure the reliability of the calculations, the DMRG method was used. The microscopic structure of vicinal surfaces was studied using the Monte Carlo method.
We showed that the characteristics of the height profile for a macrostep can be classified by the zone in the faceting diagram where the vicinal surface exists. In the step-faceting zone, the terrace for a macrostep is (001), whereas in the step droplet zone, the surface which contacts a macrostep has the slope 1 ̸ = 0. This classification is a direct result of the connectivity of the surface tension.
Furthermore, in the step-faceting zone, the macrostep is so stiff that it takes an infinitely long time to reach global equilibrium in the thermodynamic (large system size) limit. In the Monte Carlo simulation, global equilibrium is reached because of the finite size. The relaxation time increases as the temperature decreases (Figures 8 and 9). Smearing at the edge of the macrostep is observed, which is also caused by the finite size effect ( Figure 10).
The observed height profile of a macrostep of 4H-SiC [1] has the characteristics of the height profile in the step droplet zone. This suggests the surface tension is discontinuous around the side surface of the macrostep.