## 3 Results and Discussion  
To study the effects of film thickness on the domain switching pathways in the BFO thin film, we performed three cases of different film thicknesses, where \(h\) = 32, 16 and 8 nm, using the above phase-field model. The BFO films are assumed to be epitaxial (001) oriented on a rigid substrate. Firstly, the system was evolved, without any external electric field, from a paraelectric phase to an energy minimizing domain configuration in a ferroelectric rhombohedral phase. Then, the electric field was uniformly applied as 500 kV/cm to the pre-evolved film in the direction perpendicular to the film surface for all three cases. The switching behavior was characterized by switching percentage, which was defined as the ratio of the number of grids that took switch action in a period of time via one specific switching path to the total grid number in the system. These percentages were calculated every 100 steps. Figure 1: Schematic illustration of the ferroelectric BiFeO\({}_{3}\) (BFO) thin film model. In all three cases, it was observed that domains switch to electric-field-favorable orientations mainly through a 71 degree switching path and that 180 degree switching is barely observed due to the relatively high energy bias, as also reported in our previous work [14]. Figure 2a presents the typical time evolutions of the switching percentage for three paths in the 32 nm BFO film. It could be seen that, in the switching process, 180 degree switching, which consumed the most energy, remained at zero, and the percentage of the 109 degree switching path was even lower than 1%. The 71 degree path was the dominant one for polarization switching in the rhombohedral BFO films. To study the effect of thickness, switching percentages for this primary switching path in films with different thicknesses were compared (Figure 2b). For all samples, with an application of the constant electric field, polarization switching percentages burst into a maximum value in the beginning and drop to zero as time goes on. Films with different thickness have different decreasing tendencies, indicating different switching times. Here, we plotted the detailed volume fractions for upward variants R\({}^{+}\) (R\({}_{1}^{+}\), R\({}_{2}^{+}\), R\({}_{3}^{+}\), R\({}_{4}^{+}\)) and downward variants R\({}^{-}\) (R\({}_{1}^{-}\), R\({}_{2}^{-}\), R\({}_{3}^{-}\), R\({}_{4}^{-}\)) in these samples in Figure 3. As the system evolved under an applied downward electric field, energy-unfavorable variants (in this case, the upward variants R\({}^{+}\)) disappeared over time. We used the vanishing point of R\(+\) as the end of electric-field-induced polarization switching. In the 32 nm film, the volume fraction of R\(+\) dropped to zero around Step 1000. In the 16 nm film, this drop occurred at Step 3000, and, in the 8 nm film, at Step 7000. With the decrease of film thicknesses, the time that was needed for the system to finish the switching process became longer. Figure 3: Time evolution of volume fractions for downward (R- with solid line) and upward (R+ with dot line) variants. Films with different thicknesses (32, 16, and 8 nm) are indicated with different colors. Figure 2: Analysis of domain switching behavior for different film thicknesses using an identical electric field: \(E=500\) kV/cm. (**a**) Switching percentage of different switching paths as a function of time in the \(h=32\) nm film. (**b**) Percentage for 71 degree switching over time, with different film thicknesses indicated with different colors. The difference in the performance of domain switching between films with distinct thicknesses can be attributed to the influence of the coercive field. To study the relationship between film thickness and coercive field, ferroelectric hysteresis loops for films with different thicknesses--\(h=32,\,16,\,8,\,4\) nm--were calculated. Selected corresponding domain configurations were presented as well.