- Email: [email protected]
YMnO3 Ising bilayer
Superlattices and Microstructures 136 (2019) 106293
Contents lists available at ScienceDirect
Superlattices and Microstructures journal homepage: www.elsevier.com/locate/superlattices
Monte Carlo simulation of polarization plateaus and hysteresis behaviors of an antiferroelectric/ferroelectric BiFeO3/YMnO3 Ising bilayer Ze-yuan Wang a, Qi Li a, Feng Wang b, Lei Sun a, Ming Tian a, Wei Wang a, * a b
School of Science, Shenyang University of Technology, Shenyang, 110870, China School of Materials Science & Engineering, Shenyang University of Technology, Shenyang, 110870, China
A R T I C L E I N F O
A B S T R A C T
Keywords: BiFeO3/YMnO3 bilayer Polarization plateaus Electric hysteresis loops Crystal-field Exchange coupling Monte Carlo
The mixed-spin Ising model with the crystal-field and the external electric field is used to describe an antiferroelectric/ferroelectric BiFeO3/YMnO3 bilayer structure. Using Monte Carlo simulation, we discuss the polarization plateaus and the hysteresis behaviors of the bilayer structure under the influence of different physical parameters. The results illustrate that the crystal-field, the exchange coupling and the external electric field have important roles in the low-temperature dielectric properties and multiple electric hysteresis loops of the antiferroelectric/ferroelectric BiFeO3/YMnO3 bilayer structure.
1. Introduction In recent years, the multiferroic bilayer structure has turned into the focus of interest in magnetoelectric materials due to its remarkable electrical, magnetic and thermodynamic properties. More particularly, the BiFeO3/YMnO3 bilayer has been one of the significant room temperature multiferroic materials since it has a higher antiferroelectric and ferroelectric transition temperature [1–4]. The research of multiferroic bilayer is developing towards practicality and has very attractive applications in many fields, such as phase shifter [5], magnetic field sensors [6], pyroelectric detector [7], magnetoelectric read head [8] and memory devices [9]. In experiment, the ferroelectric bilayers has aroused great interests. The various ferroelectric bilayers have been manufactured by diverse methods at present, such as pulsed laser deposition (PLD) [10], Sol-Gel method [11] and metal organic chemical vapor deposition (MOCVD) method [12]. In addition, many researchers have studied their physical properties experimentally. Using the magnetron sputtering technique, Zapata et al. have studied the electric and magnetic properties of the multiferroic BiFeO3/YMnO3 thin films [13]. By means of rapid liquid phase sintering method, Wang et al. have investigated the enhanced ferromagnetic properties of the BiFeO3/YMnO3 structure and solid solutions of (1-x) BiFeO3-x/YMnO3 (x ¼ 0, 0.05, 0.1) were synthesized [14]. Saturated electric hysteresis loop has been obtained in the sample with 5% of YMnO3 fraction at room temperature, while at 10% of YMnO3 content the polarization loops present a leakage characteristic. It is noteworthy that great progress has been made in the research of other mo lecular ferroelectric bilayers such as BiFeO3/YMnO3 [15,16], BaTiO3/CoFe2O4 [17], BiFeO3/CoFe2O4 [18]. On the other hand, theoretical explanations may be used to better understand the results of experimental observations. Various models and methods have been proposed to investigate the dielectric properties of the bilayer systems. By the effective field theory
* Corresponding author. E-mail address: [email protected] (W. Wang). https://doi.org/10.1016/j.spmi.2019.106293 Received 24 July 2019; Received in revised form 1 October 2019; Accepted 11 October 2019 Available online 16 October 2019 0749-6036/© 2019 Elsevier Ltd. All rights reserved.
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 1. Schematic sketch of an antiferroelectric/ferroelectric BiFeO3/YMnO3 Ising bilayer structure. J1 and J2 denote the intralayer antiferroelectric and ferroelectric exchange couplings of different pseudo-spins within the top layer (↑↓, σ ¼ 5/2) and within the bottom layer (↑↑, S ¼ 2), respec tively. J3 represents the interlayer exchange coupling between the top and bottom layers.
(EFT), Ainane et al. have deduced the effects of the interfacial coupling and the transverse field on the hysteresis behaviors and susceptibility of the ferroelectric/antiferroelectric bilayer described by the transverse Ising model with long-range interactions [19, 20]. By use of the EFT, Benhouria et al. have investigated the hysteresis loops and dielectric properties of compositionally graded (Ba, Sr) TiO3 thin films with the transverse Ising model [21]. By using quantum Monte Carlo simulation, Feraoun et al. have examined the influence of the crystal-field and the exchange coupling on the hysteresis behaviors of the ferroelectric bilayer with an antiferroelectric interfacial coupling [22]. Furthermore, the dielectric properties of other low-dimensional ferroelectric nanomaterials have been also widely studied such as nanoparticle [23], nanowire [24–27], nanoribbon [28]. These studies all reveal the important role of the crystal-field, the exchange coupling and the electric field in the dielectric properties of the low-dimensional system. It is worth to mention that there exists a close relationship between electric and magnetic properties of nano-materials. At first, one can realize the electric phenomena from the concept of magnetism. For example, there is a great similarity between electric and magnetic quantities such as polarization and magnetization, similar definition of step effects and hysteresis behaviors both in the electric and magnetic field. But one important difference between them is that their external fields are distinct, one is the electric field E and the other is the magnetic field h. In the last decade, the mixed-spin (2, 5/2) Ising model is one of the most widely studied mixedspin Ising models. Various mixed-spin (2, 5/2) Ising models have been proposed to understand the magnetic and thermodynamic properties of many low-dimensional magnetic materials such as nanoislands [29], magnetic bilayer Kekulene structure [30] and the molecular-based magnetic systems [31–34]. The following structures are theoretically constructed of the mixed-spin (2, 5/2) into a bilayer lattice system [35–37]. As is well-known, hysteresis loop is an important characteristic of magnetic and electric properties for low-dimensional nano-materials. Recently, as one of the low dimensional systems, hysteresis behaviors of graphene systems with bilayer honeycomb structure have also attracted extensive attention [38–41]. In particular, at low temperature, the system may exhibit interesting magnetization plateaus behaviors. By the use of the EFT, Jiang et al. have studied magnetization plateaus of a ferrimagnetic nanoparticle in the presence of single-ion anisotropies [42]. It is found that the single-ion anisotropy is the key factor to affect the number of magnetization plateaus. Performing the MC simulation, Feraoun et al. have investigated the magnetization plateaus be haviors of the spin 3/2 Ising model on a nanographene layer with higher-order exchange interaction couplings [43]. It is illustrated that the exchange interactions, the crystal and the magnetic fields have essential effects on the magnetization plateaus of the system. The effects of the crystal-field, the biquadratic exchange interaction and the external magnetic field on the critical behavior, the magnetization plateaus, the blocking temperature and the hysteresis behavior of the spin-5/2 Blume-Emery-Griffiths model on a nano-graphene layer have been studied by means of the MC simulation [44]. Using the MC simulation, Masrour et al. have researched size effect on magnetic properties of a nano-graphene bilayer structure [45]. The zero-field-cooled and the field cooled magnetization plateaus are obtained for different values of external magnetic field and exchange interaction between the two blocs. Ertas¸ et al. have employed dynamic mean-field calculations to discuss the dynamic hysteresis behaviors and the thermal behaviors including the coercivity and remance of the two-dimensional mixed-spin (2, 5/2) ferrimagnetic Ising model in an oscillating magnetic field [46]. By using the MC simulation, Masrour et al. have studied the hysteresis cycle on the Bethe lattice for different values of exchange in teractions, crystal-field and sizes [47,48]. Analysis on magnetization plateaus indicates the diversified change of the spin states of the system, originating from the competition among various physical parameters. By means of the density functional theory (DFT) and full potential linear augmented plane wave (FLAPW) method, Masrour et al. have researched electronic and magnetic properties of MnAg layers [49,50]. The excellent magnetic properties and the hysteresis loop behaviors have been also studied in other various bilayer magnetic systems with different structures and spin values, such as spin-1 and spin-2 bethe lattice [51], spin-5/2 Blume-Capel model [52], spin-3/2 and spin-5/2 kekulene structure [53], spin-1 and spin-3/2 square lattice [54], spin-7/2 and spin-3 graphene armchair nanoribbons [55]. Interestingly, the recent MC study has reported the electric polarization and hysteresis cycles of an anti ferroelectric/ferroelectric BiFeO3/YMnO3 bilayer with the different exchange interactions [56]. Because the spin values of the metallic 2
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 2. The total polarization P versus Monte Carlo steps (MCS) for D1 ¼ 1, D2 ¼ 1.5, J1 ¼ 1.5, |J3| ¼ 0.5, T ¼ 0.05 and E ¼ 3.
atoms Fe and Mn are 5/2 and 2 respectively, the ferroelectric/antiferroelectric BiFeO3/YMnO3 bilayer can be well described by the mixed-spin (2, 5/2) Ising bilayer systems. The results suggest that the crystal-field, the exchange coupling and the external field are crucial for the magnetic properties and hysteresis behaviors of the magnetic system. In our previous work, the MC simulations have been focused on the dielectric properties of the ferroelectric/antiferroelectric BiFeO3/YMnO3 bilayer. The effects of the crystal-fields, the exchange couplings and the external electric field on the polarization, the dielectric susceptibility, the internal energy, the specific heat and the phase diagrams of the blocking temperature have been examined in detail [57]. Despite these above studies [56,57], few studies have been concentrated on the polarization plateaus of the ferroe lectric/antiferroelectric BiFeO3/YMnO3 bilayer with MC simulation. Whether can the BiFeO3/YMnO3 bilayer display the rich polar ization plateaus in the external electric field. Similar to magnetic plateau behaviors observed in the Ising systems? Whether being there multiple electric hysteresis behaviors? Maybe it is interesting and inspires us to explore these issues for such mixed-spin bilayer system from the perspective of the theory. Thus, our aim in this paper is to investigate the effects of the crystal-field, exchange coupling and temperature on the polarization plateaus and electric hysteresis behaviors of the BiFeO3/YMnO3 bilayer. In Sec. 2, the Ising model of the crystal-field and external electric field is proposed and the process of the MC simulation is introduced. In Sec. 3, the polarization plateaus, phase diagrams and electric hysteresis loops of the system are presented. Finally, the summary is given in Sec. 4. 2. Model and Monte Carlo simulation We propose a three-dimensional antiferroelectric/ferroelectric BiFeO3/YMnO3 bilayer honeycomb structure in Fig. 1. The bilayer system consists of an antiferroelectric top layer and a ferroelectric bottom layer, which is constructed on the N � L honeycomb lattices. We define N and L as the number of spins in each layer and the bilayer thickness (L ¼ 2), respectively. Yellow and purple balls both represent Fe atoms with spin-5/2 on the top layer and the black balls denote Mn atoms with spin-2 on the bottom layer. It is worthy to mention that such antiferroelectric/ferroelectric bilayer has been also depicted in experiment [13–16]. The Hamiltonian of the system can be described as: X X X H ¼ J1 σziA σ zjB J2 SzmC SznD J3 σ ziAðjBÞ SzmCðnDÞ
X�
σ
D1 iðjÞ
z iAðjBÞ
�2
X D2 mðnÞ
�2 SzmCðnDÞ
2 μE
X
σ
iðjÞ
z iAðjBÞ
X þ SzmCðnDÞ
!
(1)
mðnÞ
σziA ; σzjB ; SzmC ; SznD are the pseudo-spin components of different atoms in each layer along the z-axis direction and σziðjÞ ¼ � 5=2;� 3=2;� 1=2 and SzmðnÞ ¼ �2; �1; 0 denote the pseudo-spin variables of sublattices A (B) on the top layer of and sublattices C (D) on the bottom
layer, respectively., , demonstrate the summation over all pairs of nearest-neighbor pseudo-spins within the top and bottom layers, as well as between the top and bottom layers, respectively. The first three items in the Hamiltonian (1) are exchange coupling energy between pseudo-spins. The exchange coupling comes from the interaction of strong electrostatic forces between the electrons of nearest atoms. J1 and J2 stand for the intralayer exchange coupling constants between sublattices A (B) on the top layer, and sublattices C (D) on the bottom layer. J3 is the interlayer exchange coupling constant between subalttices A(B) and C(D) between the top and bottom layers. The fourth and fifth items in the Hamiltonian (1) belong to the crystal-field energy for two sublattices. In general, the crystal-field reflects the anisotropy, which is related to the direction of polarization. Polarization is carried out at different angles from the direction of polarization crystal axis. Different polarization curves indicate that the free energy for polarization is different because of different directions of the crystal-field. For such an Ising system with only z-axis pseudo-spin component, we consider the numerical value of the crystal-field to judge the influence of the anisotropy of sublattice. Therefore, D1 and D2 denote the crystal-fields for sublattices A(B) and C(D), respectively. The last two items in the Hamiltonian (1) are the electric 3
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 3. Electric field E dependence of polarization (a, c, e, g) P and (b, d, f, h) P1, P2, P3, P4 of the antiferroelectric/ferroelectric Ising bilayer structure for different values of D1 with D2 ¼ 1.5, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05.
energy of two sublattices. E and μ represent the applied external electric field and the dipole moment on site, respectively. For simplicity, we take μ ¼ 1.0. The thermal quantities are given as following: The polarizations per site of sublattices A, B, C and D on the top and bottom layers are calculated by:
4
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 4. Critical and saturation electric field dependence of the crystal-field D1 for the antiferroelectric/ferroelectric Ising bilayer structure with D2 ¼ 1.5, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05.
+ * 2 X Z P1 ¼ σ N i iA P2 ¼
+ * 2 X Z σ jB N j
(2)
+ * 2 X z P3 ¼ S N m mC + * 2 X Z P4 ¼ S N n nD and the total polarization of the system per site P is: 1 P ¼ ðP1 þ P2 þ P3 þ P4 Þ 4
(3)
The standard Monte Carlo simulation with the Metropolis algorithm [58] has been applied to simulate the bilayer anti ferroelectric/ferroelectric Ising system. We should mention that some additional simulations have been conducted to examine the impact of N on the simulation results, but it is difficult to find any clear difference between N ¼ 50 to 100. As a result, we chose N ¼ 50 for the simulation. We use the periodic boundary condition in xy plane and the free boundary along z direction. In order to improve the accuracy of the data, every data we obtained was averaged by ten sample realization. Fig. 2 shows the total polarization (P) versus Monte Carlo steps (MCS) for D1 ¼ 1, D2 ¼ 1.5, J1 ¼ 1.5, |J3| ¼ 0.5, T ¼ 0.05 and E ¼ 3. From the figure, it is found that the P become relatively stable after MCS exceed 2 � 105. Therefore, in order to equilibrate the system, we use the rest 3 � 105 MCS per site to average various thermal quantities after approximately discarding previous 2 � 105 MCS to bring the system in steady state. 3. Results and discussion In the following calculation, take value of the exchange coupling J2 ¼ 1 as a unit of energy and temperature. Assuming the pseudospins of sublattices A, C, D are along the positive direction of the longitudinal electric field E and those of sublattice B in the opposite direction. We shall present typical results of the polarization plateaus, the phase diagrams and the hysteresis loops for the 5
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 5. Electric field E dependence of polarization (a, c) P and (b, d), P1, P2, P3, P4 for the antiferroelectric/ferroelectric Ising bilayer structure for different values of D2 with D1 ¼ 1, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05.
antiferroelectric/ferroelectric BiFeO3/YMnO3 bilayer. 3.1. Polarization plateaus The effects of various physical parameters on the polarization plateaus of the bilayer system are shown in Figs. 3–10. The polarization P and sublattice polarizations P1, P2, P3, P4 as the function of the external electric field E is plotted in Fig. 3. The physical parameters are fixed as D2 ¼ 1.5, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05 for selecting five typical values of the crystal-field D1 (¼ 0.4, 2.2, 2.6, 3.5 and 5.5). In Fig. 3(a), (c), (e), (g), the P curves show obvious step effects. It is found that six polarization plateaus with P ¼ 1.0, 1.25, 1.5, 1.75, 2.0 and 2.25 exist with the E increasing from 0. The existence of six polarization plateaus corresponds to six spin configurations of the system, namely, the spins of four sublattices (σ ziA , σzjB SzmC , SznD ) take the configurations
including (5/2, 5/2, 2, 2), (5/2, 3/2, 2, 2), (5/2, 1/2, 2, 2), (5/2, 1/2, 2, 2), (5/2, 3/2, 2, 2) and (5/2, 5/2, 2, 2) respectively. Therefore, according to Eq. (3), one can get the total polarization of the system by these above spin configurations. For example, P ¼ 1.25 comes from P ¼ 14 ðP1 þ P2 þ P3 þ P4 Þ ¼ 14 ð5 =2 3=2 þ 2 þ 2Þ ¼ 1:25.In addition, the most right dash lines arrow the values of the E that are called saturation electric field ES, where the P is saturated whatever the E increases. For example, when D1 ¼ 0.4, 2.2, 2.6, 3.5 and 5.5, ES ¼ 6.1, 9.7, 10.5, 12.3 and 16.3. The other dot lines point to the values of the E are the critical electric field EC, corresponding to every polarization plateau. For instance, when D1 ¼ 0.4, the polarization plateaus appear at EC1 ¼ 4.5, EC2 ¼ 4.9, EC3 ¼ 5.3 and EC4 ¼ 5.7 as the E increases. In Fig. 3(b), (d), (e) and (g), one can notice that the polarizations P3 and P4 always take P3 ¼ P4 ¼ 2 when D1 changes from 0.4 to 5.5, but the P1 and P2 are obviously sensitive to the variation of D1 and h. For small values of |D1| (¼0.4, 2.2) in Fig. 3(b), the P1 and P2 take the values P1 ¼ 5/2 and P2 ¼ �1/2, �3/2, �5/2, respectively and the width of the P2 increases with |D1| increasing, which is responsible for the variation in the width of the polarization plateaus in Fig. 3(a). This is 6
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 6. Critical and saturation electric field dependence of the crystal-field D2 on the antiferroelectric/ferroelectric Ising bilayer structure with D1 ¼ 1, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05.
because the spins of sublattice B in the top layer are difficult to flip from low spin state to a high one in order to overcome the stronger crystal-field D1 with |D1| increasing. Further increasing |D1| ¼ 2.6, an interesting phenomenon is observed in Fig. 3(d) that the P1 and P2 starts from 3/2 and 3/2 when h ¼ 0, suggesting that the spins of sublattices A and B stay at lower spin states for the stronger crystal-field D1. With h increasing, the P1 and P2 show P1 ¼ 5/2, 3/2 and P2 ¼ �1/2, �3/2, 5/2 and approach to saturated polarization, respectively. As observed in Fig. 3(f) and (h), for larger values of |D1| (¼3.5, 5.5), the P1 and P2 may take small values P1 ¼ 1/2 and P2 ¼ 1/2, respectively, and their widths are increasing with |D1| increasing. It demonstrates that the spins of sublattices A and B in the top layer would be fond of staying at the low spin states for a long time for the stronger crystal-field D1. Similar polarization plateaus behaviors can be comparable with magnetization plateaus in many magnetic systems such as nanoparticle [42], graphene monolayers [43,44] and graphene nanoribbon [59]. Fig. 4 shows the critical electric field (EC1, EC2, EC3 and EC4) and the saturation electric field (ES) as a function of D1 for D2 ¼ 1.5, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05. From this figure, we can find that the phase diagram is divided into six regions by five values of the EC1, EC2, EC3 and EC4 and one value of the ES. Every region illustrates the different spin configurations of the system. Note that two inflection points appear for D1 ¼ 2.6 and 3.4, which reflects the change of spin stated caused by the competition between crystalfield D1 and external electric field E. For example, when E < EC1, the spin states of sublattice A and B in the top layer are σziA ¼ 5=2 and σzjB ¼ 5=2 for |D1| < 2.6, σ ziA ¼ 3=2 and σzjB ¼ 3=2 for 2.6 � |D1| < 3.4, σ ziA ¼ 1=2 and σzjB ¼ 1=2 for |D1| �3.4, while the spins of sublattices C and D SzmCðnDÞ in the bottom layer always stay SzmCðnDÞ ¼ 2 state. More detailed spin configurations under the effects of the
crystal-field D1 and external the electric field E are presented in Table 1. In Figs. 5 and 6, the effect of the crystal-field D2 on the polarization plateaus behaviors of the bilayer structure is examined. The polarization P and P1, P2, P3, P4 as the function of the E is plotted in Fig. 5(a)–(d). The physical parameters are fixed as D1 ¼ 1, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05 for selecting two typical values of the D2 ¼ 0.5, 4.0. It is observed that the P curves have two separate numbers of the polarization plateaus. One is six polarization plateaus with P ¼ 1.0, 1.25, 1.5, 1.75, 2.0 and 2.25 for D2 ¼ 0.5 in Fig. 5(a) and the other is seven polarization plateaus with P ¼ 0.50, 1.00, 1.25, 1.50, 1.75, 2.00 and 2.25 when D2 ¼ 4.0 in Fig. 5(c). The more one plateau can be found from a comparison between Fig. 5(b) and (d) that for stronger crystal-field D2 (¼ 4.0), the P3 and P4 can take P3 ¼ P4 ¼1 under weak external electric field besides P3 ¼ P4 ¼2, which should be responsible for the increase in the number of polarization plateaus and is different from the behavior of P3 and P4 in Fig. 5(b). Fig. 6 shows the critical electric field EC1, EC2, EC3, EC4, EC5 and the saturation electric field ES as function of D2 in the system for D1 ¼ 1, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05. Comparing with Figs. 6 and 4, one can notice that more spin configurations may appear for stronger D2 (� 3.4) because the phase diagram is divided into seven regions from bottom to top by the EC1, EC2, EC3, EC4, EC5 and ES curves, corresponding to seven different spin configurations. This is mainly due to the fact that the spins of sublattice C and D can take another state SzmCðnDÞ ¼ 1 besides SzmCðnDÞ ¼ 2, while for weaker crystal-field D2 (¼ 0.5), they only take the larger state SzmCðnDÞ ¼ 2. 7
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 7. Electric field E dependence of polarization (a, c, e) P and (b, d, f) P1, P2, P3, P4 of the antiferroelectric/ferroelectric Ising bilayer structure for different values of J1 with D1 ¼ 1, D2 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05.
Corresponding spin configurations under the effects of the crystal-field D2 and external electric field E are presented in Table 2. Next, we will devote our attention to the influence of the exchange coupling J1 on the polarization plateaus behaviors of the bilayer structure. Fig. 7(a)–(f) shows the polarization P and P1, P2, P3, P4 as the function of the E for different values of J1 with fixed D1 ¼ 1, D2 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05. From Fig. 7(a), (c) and (e), one can notice that the number of the plateau of the system change when the value of J1 varies from 0.2 to 2.0. For instance, there are four plateaus with P ¼ 1.5, 1.75, 2.0 and 2.25 for J1 ¼ 0.2, five plateaus with P ¼ 1.25, 1.5, 1.75, 2.0 and 2.25 for J1 ¼ 0.6 and six plateaus with P ¼ 1.0, 1.25, 1.5, 1.75, 2.0 and 2.25 for J1 ¼ -0.8, 2.0, respectively. The result can be explained from Fig. 7(b), (d) and (f). It is found that the P3 and P4 always maintain P3 ¼P4 ¼ 2. The variations of the P1 and P2 should be responsible for the plateau behaviors in Fig. 7(a), (c) and (e). When E increases from 0 to 12, for the smaller value of |J1 |(¼ 0.2), the P1 and P2 take values P1 ¼ 3/2, 5/2 and P2 ¼ 1/2, 3/2, 5/2, respectively. Further increasing |J1 |, we can see that the P1 will always take value P1 ¼ 5/2, while the P2 can take negative values P2 ¼ 3/2, 1/2 for J1 ¼ 0.6, and P2 ¼ 5/2, 3/2, 1/2 for J1 ¼ 0.8, 2.0, respectively. This is because the spins of sublattices A and B more favor for being aligned in antiparallel under the stronger antiferromagnetic coupling J1, whereas the external electric field tends to force them to align towards its direction. As a result, the more plateaus may appear in the P2 when increasing the values of J1. It is worth noting that the plateaus move right overall with the width of the first plateau getting long when J1 changes from 0.8 to 2.0 in Fig. 7(e). The phenomenon can be explained physically from a point of energy as follows: with |J1| increasing, the exchange coupling energy between sublattices become more, so a larger external electric field is required to overcome the stronger exchange coupling energy to cause the appearance of next plateaus. Fig. 8 shows the critical electric field EC1, EC2, EC3, EC4 and the saturation electric field ES as function of J1 for D1 ¼ 1, D2 ¼ 1.5, | J3| ¼ 0.5 and T ¼ 0.05. The phase diagram is divided into four regions by three dash lines, corresponding to four different numbers of 8
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 8. Critical and saturation electric field dependence of the antiferroelectric coupling J1 on the top layer for the antiferroelectric/ferroelectric Ising bilayer structure with D1 ¼ 1, D2 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05.
Fig. 9. Electric field E dependence of polarization (a) P and (b) P1, P2, P3, P4 of the antiferroelectric/ferroelectric Ising bilayer structure for different values of |J3| with D1 ¼ 1, D2 ¼ 1.5, J1 ¼ 1.5 and T ¼ 0.05.
the plateaus. Detailed spin configurations under the effects of the exchange coupling J1 and external electric field E are presented in Table 3. We also investigate the effect of the interlayer exchange coupling |J3| on the polarization plateaus behaviors of the bilayer structure, as showed in Figs. 9 and 10. The polarization P and P1, P2, P3, P4 as the function of the E is plotted in Fig. 9(a) and (b). The physical parameters are fixed as D1 ¼ 1, D2 ¼ 1.5, J1 ¼ 1.5 and T ¼ 0.05 for selecting two typical values of the |J3| ¼ 0.4 and 2.0. In Fig. 9(a), The P curves remain six polarization plateaus with P ¼ 1.0, 1.25, 1.5, 1.75, 2.0 and 2.25 and shift right when |J3| increases from 0.4 to 2.0. This behavior is similar to those observed in Fig. 7(e). We can remark that the number of the plateaus is insensitive to the change of |J3|. From Fig. 9(b), we can observe that the plateaus of the total system come mainly from the contribution of the P2. Fig. 10 shows the critical electric field EC1, EC2, EC3, EC4 and the saturation electric field ES as function of |J3| in the system for D1 ¼ 1, D2 ¼ 1.5, J1 ¼ 1.5 and T ¼ 0.05. It is clearly seen that the EC1, EC2, EC3, EC4 and ES all increase linearly as |J3| increases. We can present the corresponding spin configurations under the effects of the exchange coupling |J3| and external electric field E in Table 4.
9
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 10. Critical and saturation electric field dependence of the interlayer coupling |J3| for the antiferroelectric/ferroelectric Ising bilayer structure with D1 ¼ 1, D2 ¼ 1.5, J1 ¼ 1.5 and T ¼ 0.05.
3.2. Electric hysteresis loops In this section, we investigate the effects of various physical parameters on the electric hysteresis loops of the total system. Fig. 11 shows the electric hysteresis loops for different values of the crystal-field D1 with fixed D2 ¼ 4.0, J1 ¼ 1.8, |J3| ¼ 0.3 and T ¼ 0.3. It is interesting that various multiple hysteresis loops may be discovered for different values of D1 such as sextuple loops (D1 ¼ 0.5), tenfold loops (D1 ¼ 2.2), quadruple loops (D1 ¼ 2.6) and double loops (D1 ¼ 2.7, 3.2), as shown in Fig. 11(a)–(f). The multiple hysteresis loops behavior reflects the coexistence of multiple spin states of four sublattices, due to the competition among the crystalfield, the exchange coupling and the electric field. Physically, on one hand, the external electric field tends to make the spins of all sublattices to align towards its direction, whereas the spins tend to be an antiparallel alignment because of the antiferroelectric ex change coupling between sublattices A and B in the top layer. On the other hand, the increasing of |D1| and |E| promote the change of the spin states of four sublattices. Therefore, together these factors determine the existence of the multiple hysteresis loops. What is interesting is that as D1 increasing, the system may undergo the first-order phase transition in the zero external electric field (E ¼ 0), which should be responsible for the presence of the double loops. In addition, the loops with even number can be found because the crystal-field is so strong that the spins become disordered at zero electric field, where a first-order phase transition occurs. Similar behaviors have been observed in theoretical studies of the ferroelectric and ferrielectric nanowires [27]. Furthermore, we can also notice that the system may present a paraelectric state when D1 ¼ 2.9, as observed in Fig. 11(e), suggesting that the P is completely controlled by the external electric field. The effect of the crystal-field D2 on the hysteresis loops for the total system is illustrated in Fig. 12. The values of physical parameter are fixed D1 ¼ 2.2, J1 ¼ 1.8, |J3| ¼ 0.3 and T ¼ 0.3, respectively. When we change the values of D2 increase from 0.6 to 6.0, the system can display the septuple loops (D2 ¼ 0.6, 2.8), tenfold loops (D2 ¼ 4.8) and sextuple loops (D2 ¼ 6.0), respectively. It is worth mentioning that the loops with even number appear usually for stronger D2. This is because the spins of sublattices C and D take integer spin values and favor for the lower spin states SzmCðnDÞ ¼ 0, which can explain why the polarization can take zero for E ¼ 0. Fig. 13 shows the hysteresis loops of the total system for various values of J1 with D1 ¼ 2.2, D2 ¼ 4, |J3| ¼ 0.3 and T ¼ 0.3. When J1 changes from 0.4 to 1.8, the double, quadruple and tenfold loops appear, with the area becoming greater. The increase of the number of the loops is attributed to the competition between the antiferroelectric exchange coupling and the external electric field. It is worth noting that the area of the hysteresis loop increases with the exchange coupling increasing. This hysteresis behavior can be comparable with that under the influence of exchange coupling [26,27,56]. However, it should be mentioned that we have not observed the exchange bias phenomenon similar to that in Ref. [56]. This is because the external magnetic field is not considered for the present system, which is different from that of Ref. [56]. Fig. 14 shows the hysteresis loops of the total system for various values of |J3| with D1 ¼ 0.8, D2 ¼ 0.3, J1 ¼ 2.0 and T ¼ 0.3. It is 10
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Table 1 Spin configurations and the total polarization with the effects of D1 and E. -2.6 < D1
3.4 < D1 � 2.6
σziA ¼
5 z σ ¼ 2 jB
σziA ¼
5 z σ ¼ 2 jB
3 2
5 4
EC1 � E < EC2
σziA ¼
5 z σ ¼ 2 jB
1 3 SzmC ¼ 2 SznD ¼ 2 P ¼ 2 2
EC2 � E < EC3
σziA
5 z ¼ σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
7 P¼ 4
EC3 � E < EC4
σziA
5 z ¼ σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼2
EC4 � E < ES
5 z σ ¼ 2 jB 3 z ¼ σ ¼ 2 jB
σziA ¼ σziA
P¼1
SzmC ¼ 2 SznD ¼ 2
P¼
E < EC1
5 9 SzmC ¼ 2 SznD ¼ 2 P ¼ ES � E 2 4 3 SzmC ¼ 2 SznD ¼ 2 P ¼ 1 E < EC1 2
σziA ¼
3 z σ ¼ 2 jB
1 2
5 4
EC1 � E < EC2
σziA ¼
5 z σ ¼ 2 jB
1 3 SzmC ¼ 2 SznD ¼ 2 P ¼ 2 2
EC2 � E < EC3
σziA ¼
5 z σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
P¼
EC3 � E < EC4
σziA ¼
5 z σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼2
5 z σ ¼ 2 jB 1 z ¼ σ ¼ 2 jB
σziA ¼
D1 � 3.4
5 SzmC ¼ 2 SznD ¼ 2 2
σziA
SzmC ¼ 2 SznD ¼ 2
P¼
7 4
EC4 � E < ES
5 9 SzmC ¼ 2 SznD ¼ 2 P ¼ ES � E 2 4 1 z z SmC ¼ 2 SnD ¼ 2 P ¼ 1 E < EC1 2
σziA ¼
1 z σ ¼ 2 jB
1 2
5 4
EC1 � E < EC2
σziA ¼
3 z σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
P¼
3 2
EC2 � E < EC3
σziA ¼
3 z σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼
7 4
EC3 � E < EC4
σziA
5 z ¼ σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼2
EC4 � E < ES
σziA
5 z ¼ σ ¼ 2 jB
5 SzmC ¼ 2 SznD ¼ 2 2
9 P¼ 4
SzmC ¼ 2 SznD ¼ 2 P ¼
ES � E
seen that the twelvefold loops appear for smaller value of |J3| (¼0.1) in Fig. 14(a). This is because the weak exchange coupling |J3| make all spins of sublattices flip more easily so as to cause more loops. On increasing |J3| to 0.4 and 0.6, tenfold loops and triple loops appear. The appearance of triple loops can be explained as follows: on one hand, the spins are constrained by strong |J3| so as to flip difficultly in the small external electric field. On the other hand, only the external field increases large enough to make spins align along its direction and the existence of the ferrielectric exchange coupling between sublattices exists, which causes two small outside loops appear in the large external field. Until |J3| ¼ 0.7, single loop comes into our sight. It is found that the area of the hysteresis loop increases obviously. Similar hysteresis behaviors have been observed in the ferroelectric/antiferroelectric bilayer superlattice [19] and ferroelectric bilayer superlattice with strong interlayer exchange coupling [22]. Finally, the effect of temperature on the hysteresis loops of the total system is illustrated in Fig. 15 with fixed values D1 ¼ 2.2, D2 ¼ 0.6, J1 ¼ 1.8 and |J3| ¼ 0.3. From this figure, we can find that the multiple hysteresis loops behaviors can be weakened by the increasing of temperature. Namely, the number of the loops changes from quintuple to triple loops and single loop when T increases from 0.0001 to 4, meanwhile the area of the hysteresis loop also decreases distinctly. Finally, the loop vanishes and turns into the paraelectric state when T � 6.1. The results are consistent with previous studies of ferroelectric or ferrielectric materials [20,21,26–28, 60], antiferroelectric/ferroelectric BiFeO3/YMnO3 bilayer [56] and magnetic nanostructures [61–71]. In experiment, similar the electric hysteresis behavior under the influence of the temperature has been also observed in the single crystal of stoichiometric lithium niobite [72]. In addition, we also find similar hysteresis behavior under the influence of temperature in the MC study of dynamic hysteresis behavior in an external oscillating field for the simple ferromagnetic Ising system. It is found that the system changes from a paramagnetic state to a single-loop hysteresis behavior with temperature decreasing [73]. In order to better understand the effect of temperature on the coercivity Ehc and whether it satisfies the N� eel-Arrhenius theory [74], we plot the T and T1/2 dependences of the coercivity for the system with D1 ¼ 2.2, D2 ¼ 0.6, J1 ¼ 1.8 and |J3| ¼ 0.3, as shown in Fig. 16. The red balls represent the simulation results of the coercivity in Fig. 16(a) and (b), and the blue dotted line represents the fitted curve according to Eq. (4) in N�eel-Arrhenius theory as follows: � Ehc ¼ Ehc0 1
�
T TB
�1=2 �
(4)
11
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Table 2 Spin configurations and the total polarization with the effects of D2 and E. -3.4 < D2
5.2 < D2 � 3.4
5.6 < D2 � 5.2
σziA ¼
5 z σ ¼ 2 jB
σziA ¼
5 z σ ¼ 2 jB
3 2
SzmC ¼ 2 SznD ¼ 2
P¼
5 4
EC1 � E < EC2
σziA ¼
5 z σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
P¼
3 2
EC2 � E < EC3
σziA
5 z ¼ σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
7 P¼ 4
EC3 � E < EC4
σziA
5 z ¼ σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼2
EC4 � E < ES
σziA
P¼1
E < EC1
5 9 SzmC ¼ 2 SznD ¼ 2 P ¼ ES � E 2 4 5 1 SzmC ¼ 1 SznD ¼ 1 P ¼ E < EC1 2 2
σziA ¼
5 z σ ¼ 2 jB
5 SzmC ¼ 2 SznD ¼ 2 2
P¼1
σziA ¼
5 z σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼
5 4
EC2 � E < EC3
σziA ¼
5 z σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
P¼
3 2
EC3 � E < EC4
σziA ¼
5 z 1 7 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4
EC4 � E < EC5
σziA ¼
5 z 3 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 2 jB 2
EC5 � E < ES
5 z σ ¼ 2 jB 5 z 5 ¼ σ ¼ 2 jB 2
σziA
EC1 � E < EC2
5 9 SzmC ¼ 2 SznD ¼ 2 P ¼ ES � E 2 4 1 z z SmC ¼ 1 SnD ¼ 1 P ¼ E < EC1 2
σziA ¼
σziA ¼
5 z σ ¼ 2 jB
3 SzmC ¼ 1 SznD ¼ 1 2
P¼
3 4
EC1 � E < EC2
σziA ¼
5 z σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼
5 4
EC2 � E < EC3
5 z ¼ σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
3 P¼ 2
EC3 � E < EC4
σziA
6.0 < D2 � 5.6
5 z σ ¼ 2 jB 5 z ¼ σ ¼ 2 jB
σziA ¼
5 SzmC ¼ 2 SznD ¼ 2 2
σziA ¼
5 z 1 7 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4
σziA ¼
5 z 3 σ ¼ SzmC ¼ 2 SznD ¼ 2 2 jB 2
5 z σ ¼ 2 jB 5 z 5 ¼ σ ¼ 2 jB 2
σziA ¼ σziA
EC4 � E < EC5
P¼2
EC5 � E < ES
5 9 SzmC ¼ 2 SznD ¼ 2 P ¼ ES � E 2 4 1 SzmC ¼ 1 SznD ¼ 1 P ¼ E < EC1 2
σziA ¼
5 z σ ¼ 2 jB
3 SzmC ¼ 1 SznD ¼ 1 2
P¼
σziA ¼
5 z σ ¼ 2 jB
1 SzmC ¼ 1 SznD ¼ 1 2
P¼1
σziA ¼
5 z σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
P¼
σziA ¼
5 z 1 7 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4
EC4 � E < EC5
σziA ¼
5 z 3 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 2 jB 2
EC5 � E < ES
σziA
3 4
3 2
EC1 � E < EC2 EC2 � E < EC3 EC3 � E < EC4
5 z 5 9 ¼ σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4
ES � E
where Ehc0 is the coercivity near zero temperature and TB is the blocking temperature corresponding to the paraelectric state of hysteresis loop. Their values can be obtained from Fig. 15(a) and (d), namely TB ¼ 6.1 for the paraelectric state for and Ehc0 ¼ 2.63 at T ¼ 0.0001. From Fig. 16(a), we can see that the coercivity of the system depends strongly on temperature, and decreases with the increase of temperature. This behavior can be explained by the fact that the bilayer system has higher thermal energy at a higher temperature so that a smaller electric field is required to reverse the polarization. Similar behavior has been also observed in the experimental results of the CoxFe3-xO4 nanoparticles, which can be described by the N�eel-Arrhenius model [75]. From Fig. 16(b), it is clearly seen that the Ehc decreases almost linearly with T1/2 and the fitted curve also changes linearly, which satisfies the rule of Ehc ∝ T1=2 in N�eel-Arrhenius theory [74]. Similar results have been found in experimental studies of highly crystalline superparamagnetic 12
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Table 3 Spin configurations and the total polarization with the effects of J1 and E. -0.2 < J1
0.4 < J1 � 0.2
0.8 < J1 � 0.4
σziA ¼
5 z σ ¼ 2 jB
1 7 SzmC ¼ 2 SznD ¼ 2 P ¼ 2 4
σziA ¼
5 z σ ¼ 2 jB
3 2
P ¼ 2 EC1 � E < ES
5 z 5 9 SzmC ¼ 2 SznD ¼ 2 P ¼ σ ¼ 2 jB 2 4 3 z 1 3 z z ¼ σ ¼ SmC ¼ 2 SnD ¼ 2 P ¼ 2 jB 2 2
σziA ¼ σziA
5 z 1 σ ¼ SzmC ¼ 2 SznD ¼ 2 2 jB 2
P¼
σziA ¼
5 z 3 σ ¼ SzmC ¼ 2 SznD ¼ 2 2 jB 2
P ¼ 2 EC2 � E < ES
σziA σziA
7 4
ES � E E < EC1
σziA ¼
EC1 � E < EC2
5 z 5 9 ¼ σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4 5 z 3 5 SzmC ¼ 2 SznD ¼ 2 P ¼ ¼ σ ¼ 2 jB 2 4 1 SzmC ¼ 2 SznD ¼ 2 2
ES � E E < EC1
σziA
5 z ¼ σ ¼ 2 jB
σziA
5 z 1 7 ¼ σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4
σziA ¼ 2.0 � J1 � 0.8
SzmC ¼ 2 SznD ¼ 2
E < EC1
EC1 � E < EC2 EC2 � E < EC3
5 z 3 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 EC3 � E < ES 2 jB 2
5 z σ ¼ 2 jB 5 z 5 ¼ σ ¼ 2 jB 2
σziA ¼ σziA
3 P¼ 2
5 9 SzmC ¼ 2 SznD ¼ 2 P ¼ 2 4 SzmC ¼ 2 SznD ¼ 2 P ¼ 1
ES � E E < EC1
σziA ¼
5 z σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼
5 4
EC1 � E < EC2
σziA ¼
5 z σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
P¼
3 2
EC2 � E < EC3
σziA ¼
5 z 1 7 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4
σziA ¼
5 z 3 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 EC4 � E < ES 2 jB 2
σziA ¼
5 z 5 9 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4
EC3 � E < EC4
ES � E
Table 4 Spin configurations and the total polarization with the effects of |J3| and E. 0 � |J3| � 2.0
5 z σ ¼ 2 jB
5 SzmC ¼ 2 SznD ¼ 2 2
P¼1
5 z ¼ σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
5 P¼ 4
σziA ¼
5 z σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
P¼
σziA ¼
σziA ¼ σziA
E < EC1 EC1 � E < EC2
3 2
EC2 � E < EC3
5 z 1 σ ¼ SzmC ¼ 2 SznD ¼ 2 2 jB 2
P¼
7 4
σziA ¼
5 z 3 σ ¼ SzmC ¼ 2 SznD ¼ 2 2 jB 2
P¼2
σziA ¼
5 z 5 σ ¼ SzmC ¼ 2 SznD ¼ 2 2 jB 2
P¼
9 4
EC3 � E < EC4 EC4 � E < ES ES � E
iron oxide nanoparticles [76]. 4. Conclusion In this paper, we investigate the polarization plateaus and the hysteresis behaviors of the antiferroelectric/ferroelectric BiFeO3/ YMnO3 bilayer structure by using Monte Carlo simulation. The effects of the crystal-field, the exchange coupling and the external electric field on the polarization plateaus and the hysteresis behaviors of the system have been clarified. It is found that the multiple plateaus and the hysteresis behaviors reflect multiple spin configurations of the system, originating from the competition among the crystal-field, the exchange coupling and the external electric field. We hope that these results may be helpful for the experimental application of the antiferroelectric/ferroelectric bilayer.
13
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 11. The hysteresis loops of the antiferroelectric/ferroelectric Ising bilayer structure for different values of D1 with D2 ¼ 4, J1 ¼ 1.8, | J3| ¼ 0.3 and T ¼ 0.3.
Fig. 12. The hysteresis loops of the antiferroelectric/ferroelectric Ising bilayer structure for different values of D2 with D1 ¼ 2.2, J1 ¼ 1.8, | J3| ¼ 0.3 and T ¼ 0.3.
Declaration of competing interest We declare that we do not have any commercial or associative interest that represents a conflict of interest in connection with the work submitted. 14
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 13. The hysteresis loops of the antiferroelectric/ferroelectric Ising bilayer structure for different values of J1 with D1 ¼ 2.2, D2 ¼ 4, | J3| ¼ 0.3 and T ¼ 0.3.
Fig. 14. The hysteresis loops of the antiferroelectric/ferroelectric Ising bilayer structure for different values of |J3| with D1 ¼ 2.2, D2 ¼ 4, J1 ¼ 1.8 and T ¼ 0.3.
15
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 15. The hysteresis loops of the antiferroelectric/ferroelectric Ising bilayer structure for different values of T with D1 ¼ 2.2, D2 ¼ 0.6, J1 ¼ 1.8 and |J3| ¼ 0.3.
Fig. 16. The T and T1/2 dependences of the coercivity for the antiferroelectric/ferroelectric Ising bilayer structure with D1 ¼ 2.2, D2 ¼ 0.6, J1 ¼ 1.8 and |J3| ¼ 0.3.
Acknowledgements This project was supported by the Natural Science Foundation of Liaoning Province, China (Grant No. 20170540672), “Xingliao Elite” Program Project of Liaoning Province (No. XLYC1807021), Youth Project of Liaoning Education Department (No. LQGD2017032).
16
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51] [52] [53] [54] [55] [56] [57] [58] [59] [60] [61] [62] [63] [64] [65] [66] [67] [68] [69] [70] [71] [72] [73] [74] [75] [76]
S.M. Selbach, T. Tybell, M. Einarsrud, T. Grande, Adv. Mater. 20 (19) (2008) 3692. S.R. Das, R.N.P. Choudhary, P. Bhattacharya, R.S. Katiyar, P. Dutta, A. Manivannan, M.S. Seehra, J. Appl. Phys. 101 (2007), 034104. G. Lescano, F.M. Figueiredo, F.M.B. Marques, J. Schmidt, J. Eur. Ceram. Soc. 21 (2001) 2037. Z.J. Huang, Y. Cao, Y.Y. Sun, Y.Y. Xue, C.W. Chu, Phys. Rev. B 56 (1997) 2623. M. Ouaddari, S. Delprat, F. Vidal, M. Chaker, K. Wu, IEEE Trans. Microw. Theory Tech. 53 (2005) 1390. M. Liu, X. Li, H. Imrane, Y. Chen, T. Goodrich, Z. Cai, K.S. Ziemer, J.Y. Huang, N.X. Sun, Appl. Phys. Lett. 90 (2007) 152501. S. Corkovic, Q. Zhang, J. Appl. Phys. 105 (2009), 061610. Y. Zhang, Z. Li, C. Deng, J. Ma, Y. Lin, C.W. Nan, Appl. Phys. Lett. 92 (2008) 152510. R. Thomas, J.F. Scott, D.N. Bose, R.S. Katiyar, J. Phys. Condens. Matter 22 (2010) 423201. L. Goux, M. Gervais, F. Gervais, C. Champeaux, A. Catherinot, Int. J. Inorg. Mater. 3 (2001) 839. T. Atsuki, N. Soyama, T. Yonezawa, K. Ogi, Jpn. J. Appl. Phys. 34 (1995) 5096. T. Kijima, H. Matsunaga, Jpn. J. Appl. Phys. 38 (1999) 2281. J. Zapata, J. Narv� aez, W. Lopera, M.E. G� omez, G.A. Mendoza, P. Prieto, Transactions on Magnetics 44 (11) (2008) 2895. B. Wang, Y. Yang, W.J. Qian, L.l. Liu, Y.p. Wang, Symposium on Piezoelectricity 47 (1) (2012) 298. P.X. Nie, Y.P. Wang, Y. Yang, G.L. Yuan, W. Li, X.T. Ren, Energy harvesting, Systems 2 (3–4) (2015) 157. S.N. Tripathy, K.K. Mishra, S. Sen, B.G. Mishra, Dhiren K. Pradhan, R. Palai, D.K. Pradhan, J. Appl. Phys. 114 (2013) 144104. D.X. Zhou, G. Jian, Y.X. Hu, Y.N. Zheng, S.p. Gong, H. Liu, Mater. Chem. Phys. 127 (2011) 316. M. Alam, S. Talukdar, K. Mandal, Mater. Lett. 210 (2018) 80. A. Ainane, I. Essaoudi, M. Saber, J. Magn. Magn. Mater. 315 (2007) 132. Y. Benhouria, I. Essaoudi, A. Ainane, A. Oubelkacem, M. Saber, F. Dujardin, Phys. Scr. 86 (2012), 045704. Y. Benhouria, I. Essaoudi, A. Ainane, R. Ahuja, F. Dujardin, Chin. J. Phys. 54 (4) (2016) 533. A. Feraoun, A. Zaim, M. Kerouad, Solid State Commun. 248 (2016) 88. Y. Benhouria, I. Essaoudi, A. Ainane, R. Ahuja, F. Dujardin, J. Supercond. Nov. Magnetism 27 (9) (2014) 2153. Y. Benhouria, I. Essaoudi, A. Ainane, R. Ahuja, F. Dujardin, Superlattice Microstruct. 73 (2014) 121. Y. Benhouria, I. Essaoudi, A. Ainane, R. Ahuja, F. Dujardin, J. Supercond. Nov. Magnetism 27 (9) (2014) 2153. Y. Benhouria, I. Essaoudi, A. Ainane, R. Ahuja, F. Dujardin, Physica A 506 (2018) 499. A. Feraoun, A. Zaim, M. Kerouad, J. Phys. Chem. Solids 96 (2016) 75. S. Zriouel, B. Taychour, F. El Yahyaoui, L.B. Drissi, Physica 91 (2017) 113. A. Jabar, R. Masrour, A. Benyoussef, M. Hamedoun, J. Supercond. Nov. Magnetism 30 (7) (2017) 1807. A. Jabar, R. Masrour, Phys. B Condens. Matter 539 (15) (2018) 21. J.L. Bi, W. Wang, Q. Li, Superlattice Microstruct. 107 (2017) 104. W. Wang, W. Jiang, D. Lv, F. Zhang, J. Phys. D Appl. Phys. 45 (2012) 475002. W. Wang, W. Jiang, D. Lv, Phys. Status Solidi B 249 (1) (2012) 190. W. Wang, D. Lv, F. Zhang, J.L. Bi, J.N. Chen, J. Magn. Magn. Mater. 385 (2015) 16. R. Masrour, A. Jabar, Physica A 491 (2018) 926. M. Ertas¸, M. Keskin, Phys. Lett. A 376 (2012) 2455. M. Ertas¸, M. Keskin, B. Deviren, Physica A 391 (2012) 1038. R. Masrour, L. Bahmad, A. Benyoussef, J. Magn. Magn. Mater. 324 (23) (2012) 3991. Y. Benhouria, I. Bouziani, I. Essaoudi, A. Ainane, R. Ahuja, J. Magn. Magn. Mater. 460 (2018) 223. R. Masrour, A. Jabar, J. Comput. Electron. 16 (1) (2017) 12. R. Masrour, L. Bahmad, A. Benyoussef, M. Hamedoun, E.K. Hlil, J. Supercond. Nov. Magnetism 26 (2013) 679. W. Jiang, L. Liu, X. Li, Q. Deng, H. Guan, F. Zhang, A. Guo, Physica B 407 (2012) 3933. A. Feraoun, M. Kerouad, Solid State Commun. 280 (2018) 56. A. Feraoun, M. Kerouad, Solid State Commun. 277 (2018) 25. A. Feraoun, M. Kerouad, Computational Condensed Matter 15 (2018) 7. M. Ertas¸, Superlattice Microstruct. 85 (2015) 734. A. Jabar, R. Masrour, A. Benyoussef, M. Hamedoun, J. Magn. Magn. Mater. 397 (2016) 287. R. Masrour, A. Jabar, A. Benyoussef, M. Hamedoun, J. Magn. Magn. Mater. 393 (2015) 151. R. Masrour, E.K. Hlil, M. Hamedoun, A. Benyoussef, O. Mounkachi, H. Moussaoui, Physica A 395 (2014) 128. R. Masrour, E.K. Hlil, M. Hamedoun, A. Benyoussef, O. Mounkachi, L. Bahmad, J. Magn. Magn. Mater. 326 (2013) 166. R. Masrour, A. Jabar, A. Benyoussef, M. Hamedoun, J. Magn. Magn. Mater. 401 (2016) 700. A. Jabar, R. Masrour, Solid State Commun. 268 (2017) 38. A. Jabar, R. Masrour, Chem. Phys. Lett. 700 (2018) 130. A. Jabar, R. Masrour, A. Benyoussef, M. Hamedoun, Chem. Phys. Lett. 670 (2017) 16. R. Masrour, A. Jabar, J. Magn. Magn. Mater. 426 (2017) 225. A. Jabar, R. Masrour, A. Benyoussef, M. Hamedoun, J. Supercond. Nov. Magnetism 29 (3) (2016) 733. Z.Y. Wang, W. Wang, Q. Li, M. Tian, Z. Gao, Y. Liu, Physica 110 (2019) 127. N. Metropolis, A.W. Rosenbluth, M.N. Rosenbluth, A.H. Teller, E. Teller, J. Chem. Phys. 21 (1953) 1087. W. Wang, Q. Li, D. Lv, R. Liu, Z. Peng, S. Yang, Carbon 120 (2017) 313. Y. Wu, D. Yao, Z. Li, J. Appl. Phys. 91 (2002) 1482. D. Lv, F. wang, R.J. Liu, Q. Xue, S.X. Li, J. Alloy. Comp. 701 (2016) 935. D. Lv, Y. Yang, W. Jiang, F. Wang, Z.Y. Gao, M. Tian, Physica A 514 (2019) 319. W. Wang, S.Q. Yang, Y. Yang, Z. Peng, B.C. Li, M. Yang, Physica 109 (2019) 30. Y. Yang, W. Wang, H. Ma, Q. Li, Z.Y. Gao, T. Huang, Physica 108 (2019) 358. H.J. Wu, W. Wang, B.C. Li, M. Yang, S.Q. Yang, F. Wang, Physica 112 (2019) 86. D. Lv, W. Jiang, Y. Ma, Z.Y. Gao, F. Wang, Physica 106 (2019) 101. D. Lv, W. Wang, J.P. Liu, D.Q. Guo, S.X. Li, J. Magn. Magn. Mater. 465 (2018) 348. W. Wang, D.D. Chen, D. Lv, J.P. Liu, Q. Li, Z. Peng, J. Phys. Chem. Solids 108 (2017) 39. S.Q. Yang, W. Wang, F. Wang, B.C. Li, J.H. Xu, J. Phys. Chem. Solids 135 (2019) 109110. W. Wang, Y. Liu, Z.Y. Gao, X.R. Zhao, Y. Yang, S. Yang, Physica E101 (2018) 110. W. Wang, J.L. Bi, R.J. Liu, X. Chen, J.P. Liu, Superlattice Microstruct. 98 (2016) 433. H. Bo, Q. Meng, H. Hu, H. Zhao, Z. Zhang, Q. Zhang, C. Zhang, Appl. Phys. A 124 (2018) 691. K. Bikas, M. Acharyya Chakrabarti, Rev. Mod. Phys. 71 (1999) 847. L. N�eel, Ann. Geophys. 5 (1949) 99. T.E. Torres, E. Lima, A. Mayoral, A. Ibarra, C. Marquina, M.R. Ibarra, G.F. Goya, J. Appl. Phys. 118 (2015) 183902. M. Tadi�c, V. Kusigerski, D. Markovi�c, M. Panjan, I. Milo�sevi�c, V. Spasojevi�c, J. Alloy. Comp. 525 (2012) 28.
17
Contents lists available at ScienceDirect
Superlattices and Microstructures journal homepage: www.elsevier.com/locate/superlattices
Monte Carlo simulation of polarization plateaus and hysteresis behaviors of an antiferroelectric/ferroelectric BiFeO3/YMnO3 Ising bilayer Ze-yuan Wang a, Qi Li a, Feng Wang b, Lei Sun a, Ming Tian a, Wei Wang a, * a b
School of Science, Shenyang University of Technology, Shenyang, 110870, China School of Materials Science & Engineering, Shenyang University of Technology, Shenyang, 110870, China
A R T I C L E I N F O
A B S T R A C T
Keywords: BiFeO3/YMnO3 bilayer Polarization plateaus Electric hysteresis loops Crystal-field Exchange coupling Monte Carlo
The mixed-spin Ising model with the crystal-field and the external electric field is used to describe an antiferroelectric/ferroelectric BiFeO3/YMnO3 bilayer structure. Using Monte Carlo simulation, we discuss the polarization plateaus and the hysteresis behaviors of the bilayer structure under the influence of different physical parameters. The results illustrate that the crystal-field, the exchange coupling and the external electric field have important roles in the low-temperature dielectric properties and multiple electric hysteresis loops of the antiferroelectric/ferroelectric BiFeO3/YMnO3 bilayer structure.
1. Introduction In recent years, the multiferroic bilayer structure has turned into the focus of interest in magnetoelectric materials due to its remarkable electrical, magnetic and thermodynamic properties. More particularly, the BiFeO3/YMnO3 bilayer has been one of the significant room temperature multiferroic materials since it has a higher antiferroelectric and ferroelectric transition temperature [1–4]. The research of multiferroic bilayer is developing towards practicality and has very attractive applications in many fields, such as phase shifter [5], magnetic field sensors [6], pyroelectric detector [7], magnetoelectric read head [8] and memory devices [9]. In experiment, the ferroelectric bilayers has aroused great interests. The various ferroelectric bilayers have been manufactured by diverse methods at present, such as pulsed laser deposition (PLD) [10], Sol-Gel method [11] and metal organic chemical vapor deposition (MOCVD) method [12]. In addition, many researchers have studied their physical properties experimentally. Using the magnetron sputtering technique, Zapata et al. have studied the electric and magnetic properties of the multiferroic BiFeO3/YMnO3 thin films [13]. By means of rapid liquid phase sintering method, Wang et al. have investigated the enhanced ferromagnetic properties of the BiFeO3/YMnO3 structure and solid solutions of (1-x) BiFeO3-x/YMnO3 (x ¼ 0, 0.05, 0.1) were synthesized [14]. Saturated electric hysteresis loop has been obtained in the sample with 5% of YMnO3 fraction at room temperature, while at 10% of YMnO3 content the polarization loops present a leakage characteristic. It is noteworthy that great progress has been made in the research of other mo lecular ferroelectric bilayers such as BiFeO3/YMnO3 [15,16], BaTiO3/CoFe2O4 [17], BiFeO3/CoFe2O4 [18]. On the other hand, theoretical explanations may be used to better understand the results of experimental observations. Various models and methods have been proposed to investigate the dielectric properties of the bilayer systems. By the effective field theory
* Corresponding author. E-mail address: [email protected] (W. Wang). https://doi.org/10.1016/j.spmi.2019.106293 Received 24 July 2019; Received in revised form 1 October 2019; Accepted 11 October 2019 Available online 16 October 2019 0749-6036/© 2019 Elsevier Ltd. All rights reserved.
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 1. Schematic sketch of an antiferroelectric/ferroelectric BiFeO3/YMnO3 Ising bilayer structure. J1 and J2 denote the intralayer antiferroelectric and ferroelectric exchange couplings of different pseudo-spins within the top layer (↑↓, σ ¼ 5/2) and within the bottom layer (↑↑, S ¼ 2), respec tively. J3 represents the interlayer exchange coupling between the top and bottom layers.
(EFT), Ainane et al. have deduced the effects of the interfacial coupling and the transverse field on the hysteresis behaviors and susceptibility of the ferroelectric/antiferroelectric bilayer described by the transverse Ising model with long-range interactions [19, 20]. By use of the EFT, Benhouria et al. have investigated the hysteresis loops and dielectric properties of compositionally graded (Ba, Sr) TiO3 thin films with the transverse Ising model [21]. By using quantum Monte Carlo simulation, Feraoun et al. have examined the influence of the crystal-field and the exchange coupling on the hysteresis behaviors of the ferroelectric bilayer with an antiferroelectric interfacial coupling [22]. Furthermore, the dielectric properties of other low-dimensional ferroelectric nanomaterials have been also widely studied such as nanoparticle [23], nanowire [24–27], nanoribbon [28]. These studies all reveal the important role of the crystal-field, the exchange coupling and the electric field in the dielectric properties of the low-dimensional system. It is worth to mention that there exists a close relationship between electric and magnetic properties of nano-materials. At first, one can realize the electric phenomena from the concept of magnetism. For example, there is a great similarity between electric and magnetic quantities such as polarization and magnetization, similar definition of step effects and hysteresis behaviors both in the electric and magnetic field. But one important difference between them is that their external fields are distinct, one is the electric field E and the other is the magnetic field h. In the last decade, the mixed-spin (2, 5/2) Ising model is one of the most widely studied mixedspin Ising models. Various mixed-spin (2, 5/2) Ising models have been proposed to understand the magnetic and thermodynamic properties of many low-dimensional magnetic materials such as nanoislands [29], magnetic bilayer Kekulene structure [30] and the molecular-based magnetic systems [31–34]. The following structures are theoretically constructed of the mixed-spin (2, 5/2) into a bilayer lattice system [35–37]. As is well-known, hysteresis loop is an important characteristic of magnetic and electric properties for low-dimensional nano-materials. Recently, as one of the low dimensional systems, hysteresis behaviors of graphene systems with bilayer honeycomb structure have also attracted extensive attention [38–41]. In particular, at low temperature, the system may exhibit interesting magnetization plateaus behaviors. By the use of the EFT, Jiang et al. have studied magnetization plateaus of a ferrimagnetic nanoparticle in the presence of single-ion anisotropies [42]. It is found that the single-ion anisotropy is the key factor to affect the number of magnetization plateaus. Performing the MC simulation, Feraoun et al. have investigated the magnetization plateaus be haviors of the spin 3/2 Ising model on a nanographene layer with higher-order exchange interaction couplings [43]. It is illustrated that the exchange interactions, the crystal and the magnetic fields have essential effects on the magnetization plateaus of the system. The effects of the crystal-field, the biquadratic exchange interaction and the external magnetic field on the critical behavior, the magnetization plateaus, the blocking temperature and the hysteresis behavior of the spin-5/2 Blume-Emery-Griffiths model on a nano-graphene layer have been studied by means of the MC simulation [44]. Using the MC simulation, Masrour et al. have researched size effect on magnetic properties of a nano-graphene bilayer structure [45]. The zero-field-cooled and the field cooled magnetization plateaus are obtained for different values of external magnetic field and exchange interaction between the two blocs. Ertas¸ et al. have employed dynamic mean-field calculations to discuss the dynamic hysteresis behaviors and the thermal behaviors including the coercivity and remance of the two-dimensional mixed-spin (2, 5/2) ferrimagnetic Ising model in an oscillating magnetic field [46]. By using the MC simulation, Masrour et al. have studied the hysteresis cycle on the Bethe lattice for different values of exchange in teractions, crystal-field and sizes [47,48]. Analysis on magnetization plateaus indicates the diversified change of the spin states of the system, originating from the competition among various physical parameters. By means of the density functional theory (DFT) and full potential linear augmented plane wave (FLAPW) method, Masrour et al. have researched electronic and magnetic properties of MnAg layers [49,50]. The excellent magnetic properties and the hysteresis loop behaviors have been also studied in other various bilayer magnetic systems with different structures and spin values, such as spin-1 and spin-2 bethe lattice [51], spin-5/2 Blume-Capel model [52], spin-3/2 and spin-5/2 kekulene structure [53], spin-1 and spin-3/2 square lattice [54], spin-7/2 and spin-3 graphene armchair nanoribbons [55]. Interestingly, the recent MC study has reported the electric polarization and hysteresis cycles of an anti ferroelectric/ferroelectric BiFeO3/YMnO3 bilayer with the different exchange interactions [56]. Because the spin values of the metallic 2
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 2. The total polarization P versus Monte Carlo steps (MCS) for D1 ¼ 1, D2 ¼ 1.5, J1 ¼ 1.5, |J3| ¼ 0.5, T ¼ 0.05 and E ¼ 3.
atoms Fe and Mn are 5/2 and 2 respectively, the ferroelectric/antiferroelectric BiFeO3/YMnO3 bilayer can be well described by the mixed-spin (2, 5/2) Ising bilayer systems. The results suggest that the crystal-field, the exchange coupling and the external field are crucial for the magnetic properties and hysteresis behaviors of the magnetic system. In our previous work, the MC simulations have been focused on the dielectric properties of the ferroelectric/antiferroelectric BiFeO3/YMnO3 bilayer. The effects of the crystal-fields, the exchange couplings and the external electric field on the polarization, the dielectric susceptibility, the internal energy, the specific heat and the phase diagrams of the blocking temperature have been examined in detail [57]. Despite these above studies [56,57], few studies have been concentrated on the polarization plateaus of the ferroe lectric/antiferroelectric BiFeO3/YMnO3 bilayer with MC simulation. Whether can the BiFeO3/YMnO3 bilayer display the rich polar ization plateaus in the external electric field. Similar to magnetic plateau behaviors observed in the Ising systems? Whether being there multiple electric hysteresis behaviors? Maybe it is interesting and inspires us to explore these issues for such mixed-spin bilayer system from the perspective of the theory. Thus, our aim in this paper is to investigate the effects of the crystal-field, exchange coupling and temperature on the polarization plateaus and electric hysteresis behaviors of the BiFeO3/YMnO3 bilayer. In Sec. 2, the Ising model of the crystal-field and external electric field is proposed and the process of the MC simulation is introduced. In Sec. 3, the polarization plateaus, phase diagrams and electric hysteresis loops of the system are presented. Finally, the summary is given in Sec. 4. 2. Model and Monte Carlo simulation We propose a three-dimensional antiferroelectric/ferroelectric BiFeO3/YMnO3 bilayer honeycomb structure in Fig. 1. The bilayer system consists of an antiferroelectric top layer and a ferroelectric bottom layer, which is constructed on the N � L honeycomb lattices. We define N and L as the number of spins in each layer and the bilayer thickness (L ¼ 2), respectively. Yellow and purple balls both represent Fe atoms with spin-5/2 on the top layer and the black balls denote Mn atoms with spin-2 on the bottom layer. It is worthy to mention that such antiferroelectric/ferroelectric bilayer has been also depicted in experiment [13–16]. The Hamiltonian of the system can be described as: X X X H ¼ J1 σziA σ zjB J2 SzmC SznD J3 σ ziAðjBÞ SzmCðnDÞ
X�
σ
D1 iðjÞ
z iAðjBÞ
�2
X D2 mðnÞ
�2 SzmCðnDÞ
2 μE
X
σ
iðjÞ
z iAðjBÞ
X þ SzmCðnDÞ
!
(1)
mðnÞ
σziA ; σzjB ; SzmC ; SznD are the pseudo-spin components of different atoms in each layer along the z-axis direction and σziðjÞ ¼ � 5=2;� 3=2;� 1=2 and SzmðnÞ ¼ �2; �1; 0 denote the pseudo-spin variables of sublattices A (B) on the top layer of and sublattices C (D) on the bottom
layer, respectively.
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 3. Electric field E dependence of polarization (a, c, e, g) P and (b, d, f, h) P1, P2, P3, P4 of the antiferroelectric/ferroelectric Ising bilayer structure for different values of D1 with D2 ¼ 1.5, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05.
energy of two sublattices. E and μ represent the applied external electric field and the dipole moment on site, respectively. For simplicity, we take μ ¼ 1.0. The thermal quantities are given as following: The polarizations per site of sublattices A, B, C and D on the top and bottom layers are calculated by:
4
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 4. Critical and saturation electric field dependence of the crystal-field D1 for the antiferroelectric/ferroelectric Ising bilayer structure with D2 ¼ 1.5, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05.
+ * 2 X Z P1 ¼ σ N i iA P2 ¼
+ * 2 X Z σ jB N j
(2)
+ * 2 X z P3 ¼ S N m mC + * 2 X Z P4 ¼ S N n nD and the total polarization of the system per site P is: 1 P ¼ ðP1 þ P2 þ P3 þ P4 Þ 4
(3)
The standard Monte Carlo simulation with the Metropolis algorithm [58] has been applied to simulate the bilayer anti ferroelectric/ferroelectric Ising system. We should mention that some additional simulations have been conducted to examine the impact of N on the simulation results, but it is difficult to find any clear difference between N ¼ 50 to 100. As a result, we chose N ¼ 50 for the simulation. We use the periodic boundary condition in xy plane and the free boundary along z direction. In order to improve the accuracy of the data, every data we obtained was averaged by ten sample realization. Fig. 2 shows the total polarization (P) versus Monte Carlo steps (MCS) for D1 ¼ 1, D2 ¼ 1.5, J1 ¼ 1.5, |J3| ¼ 0.5, T ¼ 0.05 and E ¼ 3. From the figure, it is found that the P become relatively stable after MCS exceed 2 � 105. Therefore, in order to equilibrate the system, we use the rest 3 � 105 MCS per site to average various thermal quantities after approximately discarding previous 2 � 105 MCS to bring the system in steady state. 3. Results and discussion In the following calculation, take value of the exchange coupling J2 ¼ 1 as a unit of energy and temperature. Assuming the pseudospins of sublattices A, C, D are along the positive direction of the longitudinal electric field E and those of sublattice B in the opposite direction. We shall present typical results of the polarization plateaus, the phase diagrams and the hysteresis loops for the 5
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 5. Electric field E dependence of polarization (a, c) P and (b, d), P1, P2, P3, P4 for the antiferroelectric/ferroelectric Ising bilayer structure for different values of D2 with D1 ¼ 1, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05.
antiferroelectric/ferroelectric BiFeO3/YMnO3 bilayer. 3.1. Polarization plateaus The effects of various physical parameters on the polarization plateaus of the bilayer system are shown in Figs. 3–10. The polarization P and sublattice polarizations P1, P2, P3, P4 as the function of the external electric field E is plotted in Fig. 3. The physical parameters are fixed as D2 ¼ 1.5, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05 for selecting five typical values of the crystal-field D1 (¼ 0.4, 2.2, 2.6, 3.5 and 5.5). In Fig. 3(a), (c), (e), (g), the P curves show obvious step effects. It is found that six polarization plateaus with P ¼ 1.0, 1.25, 1.5, 1.75, 2.0 and 2.25 exist with the E increasing from 0. The existence of six polarization plateaus corresponds to six spin configurations of the system, namely, the spins of four sublattices (σ ziA , σzjB SzmC , SznD ) take the configurations
including (5/2, 5/2, 2, 2), (5/2, 3/2, 2, 2), (5/2, 1/2, 2, 2), (5/2, 1/2, 2, 2), (5/2, 3/2, 2, 2) and (5/2, 5/2, 2, 2) respectively. Therefore, according to Eq. (3), one can get the total polarization of the system by these above spin configurations. For example, P ¼ 1.25 comes from P ¼ 14 ðP1 þ P2 þ P3 þ P4 Þ ¼ 14 ð5 =2 3=2 þ 2 þ 2Þ ¼ 1:25.In addition, the most right dash lines arrow the values of the E that are called saturation electric field ES, where the P is saturated whatever the E increases. For example, when D1 ¼ 0.4, 2.2, 2.6, 3.5 and 5.5, ES ¼ 6.1, 9.7, 10.5, 12.3 and 16.3. The other dot lines point to the values of the E are the critical electric field EC, corresponding to every polarization plateau. For instance, when D1 ¼ 0.4, the polarization plateaus appear at EC1 ¼ 4.5, EC2 ¼ 4.9, EC3 ¼ 5.3 and EC4 ¼ 5.7 as the E increases. In Fig. 3(b), (d), (e) and (g), one can notice that the polarizations P3 and P4 always take P3 ¼ P4 ¼ 2 when D1 changes from 0.4 to 5.5, but the P1 and P2 are obviously sensitive to the variation of D1 and h. For small values of |D1| (¼0.4, 2.2) in Fig. 3(b), the P1 and P2 take the values P1 ¼ 5/2 and P2 ¼ �1/2, �3/2, �5/2, respectively and the width of the P2 increases with |D1| increasing, which is responsible for the variation in the width of the polarization plateaus in Fig. 3(a). This is 6
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 6. Critical and saturation electric field dependence of the crystal-field D2 on the antiferroelectric/ferroelectric Ising bilayer structure with D1 ¼ 1, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05.
because the spins of sublattice B in the top layer are difficult to flip from low spin state to a high one in order to overcome the stronger crystal-field D1 with |D1| increasing. Further increasing |D1| ¼ 2.6, an interesting phenomenon is observed in Fig. 3(d) that the P1 and P2 starts from 3/2 and 3/2 when h ¼ 0, suggesting that the spins of sublattices A and B stay at lower spin states for the stronger crystal-field D1. With h increasing, the P1 and P2 show P1 ¼ 5/2, 3/2 and P2 ¼ �1/2, �3/2, 5/2 and approach to saturated polarization, respectively. As observed in Fig. 3(f) and (h), for larger values of |D1| (¼3.5, 5.5), the P1 and P2 may take small values P1 ¼ 1/2 and P2 ¼ 1/2, respectively, and their widths are increasing with |D1| increasing. It demonstrates that the spins of sublattices A and B in the top layer would be fond of staying at the low spin states for a long time for the stronger crystal-field D1. Similar polarization plateaus behaviors can be comparable with magnetization plateaus in many magnetic systems such as nanoparticle [42], graphene monolayers [43,44] and graphene nanoribbon [59]. Fig. 4 shows the critical electric field (EC1, EC2, EC3 and EC4) and the saturation electric field (ES) as a function of D1 for D2 ¼ 1.5, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05. From this figure, we can find that the phase diagram is divided into six regions by five values of the EC1, EC2, EC3 and EC4 and one value of the ES. Every region illustrates the different spin configurations of the system. Note that two inflection points appear for D1 ¼ 2.6 and 3.4, which reflects the change of spin stated caused by the competition between crystalfield D1 and external electric field E. For example, when E < EC1, the spin states of sublattice A and B in the top layer are σziA ¼ 5=2 and σzjB ¼ 5=2 for |D1| < 2.6, σ ziA ¼ 3=2 and σzjB ¼ 3=2 for 2.6 � |D1| < 3.4, σ ziA ¼ 1=2 and σzjB ¼ 1=2 for |D1| �3.4, while the spins of sublattices C and D SzmCðnDÞ in the bottom layer always stay SzmCðnDÞ ¼ 2 state. More detailed spin configurations under the effects of the
crystal-field D1 and external the electric field E are presented in Table 1. In Figs. 5 and 6, the effect of the crystal-field D2 on the polarization plateaus behaviors of the bilayer structure is examined. The polarization P and P1, P2, P3, P4 as the function of the E is plotted in Fig. 5(a)–(d). The physical parameters are fixed as D1 ¼ 1, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05 for selecting two typical values of the D2 ¼ 0.5, 4.0. It is observed that the P curves have two separate numbers of the polarization plateaus. One is six polarization plateaus with P ¼ 1.0, 1.25, 1.5, 1.75, 2.0 and 2.25 for D2 ¼ 0.5 in Fig. 5(a) and the other is seven polarization plateaus with P ¼ 0.50, 1.00, 1.25, 1.50, 1.75, 2.00 and 2.25 when D2 ¼ 4.0 in Fig. 5(c). The more one plateau can be found from a comparison between Fig. 5(b) and (d) that for stronger crystal-field D2 (¼ 4.0), the P3 and P4 can take P3 ¼ P4 ¼1 under weak external electric field besides P3 ¼ P4 ¼2, which should be responsible for the increase in the number of polarization plateaus and is different from the behavior of P3 and P4 in Fig. 5(b). Fig. 6 shows the critical electric field EC1, EC2, EC3, EC4, EC5 and the saturation electric field ES as function of D2 in the system for D1 ¼ 1, J1 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05. Comparing with Figs. 6 and 4, one can notice that more spin configurations may appear for stronger D2 (� 3.4) because the phase diagram is divided into seven regions from bottom to top by the EC1, EC2, EC3, EC4, EC5 and ES curves, corresponding to seven different spin configurations. This is mainly due to the fact that the spins of sublattice C and D can take another state SzmCðnDÞ ¼ 1 besides SzmCðnDÞ ¼ 2, while for weaker crystal-field D2 (¼ 0.5), they only take the larger state SzmCðnDÞ ¼ 2. 7
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 7. Electric field E dependence of polarization (a, c, e) P and (b, d, f) P1, P2, P3, P4 of the antiferroelectric/ferroelectric Ising bilayer structure for different values of J1 with D1 ¼ 1, D2 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05.
Corresponding spin configurations under the effects of the crystal-field D2 and external electric field E are presented in Table 2. Next, we will devote our attention to the influence of the exchange coupling J1 on the polarization plateaus behaviors of the bilayer structure. Fig. 7(a)–(f) shows the polarization P and P1, P2, P3, P4 as the function of the E for different values of J1 with fixed D1 ¼ 1, D2 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05. From Fig. 7(a), (c) and (e), one can notice that the number of the plateau of the system change when the value of J1 varies from 0.2 to 2.0. For instance, there are four plateaus with P ¼ 1.5, 1.75, 2.0 and 2.25 for J1 ¼ 0.2, five plateaus with P ¼ 1.25, 1.5, 1.75, 2.0 and 2.25 for J1 ¼ 0.6 and six plateaus with P ¼ 1.0, 1.25, 1.5, 1.75, 2.0 and 2.25 for J1 ¼ -0.8, 2.0, respectively. The result can be explained from Fig. 7(b), (d) and (f). It is found that the P3 and P4 always maintain P3 ¼P4 ¼ 2. The variations of the P1 and P2 should be responsible for the plateau behaviors in Fig. 7(a), (c) and (e). When E increases from 0 to 12, for the smaller value of |J1 |(¼ 0.2), the P1 and P2 take values P1 ¼ 3/2, 5/2 and P2 ¼ 1/2, 3/2, 5/2, respectively. Further increasing |J1 |, we can see that the P1 will always take value P1 ¼ 5/2, while the P2 can take negative values P2 ¼ 3/2, 1/2 for J1 ¼ 0.6, and P2 ¼ 5/2, 3/2, 1/2 for J1 ¼ 0.8, 2.0, respectively. This is because the spins of sublattices A and B more favor for being aligned in antiparallel under the stronger antiferromagnetic coupling J1, whereas the external electric field tends to force them to align towards its direction. As a result, the more plateaus may appear in the P2 when increasing the values of J1. It is worth noting that the plateaus move right overall with the width of the first plateau getting long when J1 changes from 0.8 to 2.0 in Fig. 7(e). The phenomenon can be explained physically from a point of energy as follows: with |J1| increasing, the exchange coupling energy between sublattices become more, so a larger external electric field is required to overcome the stronger exchange coupling energy to cause the appearance of next plateaus. Fig. 8 shows the critical electric field EC1, EC2, EC3, EC4 and the saturation electric field ES as function of J1 for D1 ¼ 1, D2 ¼ 1.5, | J3| ¼ 0.5 and T ¼ 0.05. The phase diagram is divided into four regions by three dash lines, corresponding to four different numbers of 8
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 8. Critical and saturation electric field dependence of the antiferroelectric coupling J1 on the top layer for the antiferroelectric/ferroelectric Ising bilayer structure with D1 ¼ 1, D2 ¼ 1.5, |J3| ¼ 0.5 and T ¼ 0.05.
Fig. 9. Electric field E dependence of polarization (a) P and (b) P1, P2, P3, P4 of the antiferroelectric/ferroelectric Ising bilayer structure for different values of |J3| with D1 ¼ 1, D2 ¼ 1.5, J1 ¼ 1.5 and T ¼ 0.05.
the plateaus. Detailed spin configurations under the effects of the exchange coupling J1 and external electric field E are presented in Table 3. We also investigate the effect of the interlayer exchange coupling |J3| on the polarization plateaus behaviors of the bilayer structure, as showed in Figs. 9 and 10. The polarization P and P1, P2, P3, P4 as the function of the E is plotted in Fig. 9(a) and (b). The physical parameters are fixed as D1 ¼ 1, D2 ¼ 1.5, J1 ¼ 1.5 and T ¼ 0.05 for selecting two typical values of the |J3| ¼ 0.4 and 2.0. In Fig. 9(a), The P curves remain six polarization plateaus with P ¼ 1.0, 1.25, 1.5, 1.75, 2.0 and 2.25 and shift right when |J3| increases from 0.4 to 2.0. This behavior is similar to those observed in Fig. 7(e). We can remark that the number of the plateaus is insensitive to the change of |J3|. From Fig. 9(b), we can observe that the plateaus of the total system come mainly from the contribution of the P2. Fig. 10 shows the critical electric field EC1, EC2, EC3, EC4 and the saturation electric field ES as function of |J3| in the system for D1 ¼ 1, D2 ¼ 1.5, J1 ¼ 1.5 and T ¼ 0.05. It is clearly seen that the EC1, EC2, EC3, EC4 and ES all increase linearly as |J3| increases. We can present the corresponding spin configurations under the effects of the exchange coupling |J3| and external electric field E in Table 4.
9
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 10. Critical and saturation electric field dependence of the interlayer coupling |J3| for the antiferroelectric/ferroelectric Ising bilayer structure with D1 ¼ 1, D2 ¼ 1.5, J1 ¼ 1.5 and T ¼ 0.05.
3.2. Electric hysteresis loops In this section, we investigate the effects of various physical parameters on the electric hysteresis loops of the total system. Fig. 11 shows the electric hysteresis loops for different values of the crystal-field D1 with fixed D2 ¼ 4.0, J1 ¼ 1.8, |J3| ¼ 0.3 and T ¼ 0.3. It is interesting that various multiple hysteresis loops may be discovered for different values of D1 such as sextuple loops (D1 ¼ 0.5), tenfold loops (D1 ¼ 2.2), quadruple loops (D1 ¼ 2.6) and double loops (D1 ¼ 2.7, 3.2), as shown in Fig. 11(a)–(f). The multiple hysteresis loops behavior reflects the coexistence of multiple spin states of four sublattices, due to the competition among the crystalfield, the exchange coupling and the electric field. Physically, on one hand, the external electric field tends to make the spins of all sublattices to align towards its direction, whereas the spins tend to be an antiparallel alignment because of the antiferroelectric ex change coupling between sublattices A and B in the top layer. On the other hand, the increasing of |D1| and |E| promote the change of the spin states of four sublattices. Therefore, together these factors determine the existence of the multiple hysteresis loops. What is interesting is that as D1 increasing, the system may undergo the first-order phase transition in the zero external electric field (E ¼ 0), which should be responsible for the presence of the double loops. In addition, the loops with even number can be found because the crystal-field is so strong that the spins become disordered at zero electric field, where a first-order phase transition occurs. Similar behaviors have been observed in theoretical studies of the ferroelectric and ferrielectric nanowires [27]. Furthermore, we can also notice that the system may present a paraelectric state when D1 ¼ 2.9, as observed in Fig. 11(e), suggesting that the P is completely controlled by the external electric field. The effect of the crystal-field D2 on the hysteresis loops for the total system is illustrated in Fig. 12. The values of physical parameter are fixed D1 ¼ 2.2, J1 ¼ 1.8, |J3| ¼ 0.3 and T ¼ 0.3, respectively. When we change the values of D2 increase from 0.6 to 6.0, the system can display the septuple loops (D2 ¼ 0.6, 2.8), tenfold loops (D2 ¼ 4.8) and sextuple loops (D2 ¼ 6.0), respectively. It is worth mentioning that the loops with even number appear usually for stronger D2. This is because the spins of sublattices C and D take integer spin values and favor for the lower spin states SzmCðnDÞ ¼ 0, which can explain why the polarization can take zero for E ¼ 0. Fig. 13 shows the hysteresis loops of the total system for various values of J1 with D1 ¼ 2.2, D2 ¼ 4, |J3| ¼ 0.3 and T ¼ 0.3. When J1 changes from 0.4 to 1.8, the double, quadruple and tenfold loops appear, with the area becoming greater. The increase of the number of the loops is attributed to the competition between the antiferroelectric exchange coupling and the external electric field. It is worth noting that the area of the hysteresis loop increases with the exchange coupling increasing. This hysteresis behavior can be comparable with that under the influence of exchange coupling [26,27,56]. However, it should be mentioned that we have not observed the exchange bias phenomenon similar to that in Ref. [56]. This is because the external magnetic field is not considered for the present system, which is different from that of Ref. [56]. Fig. 14 shows the hysteresis loops of the total system for various values of |J3| with D1 ¼ 0.8, D2 ¼ 0.3, J1 ¼ 2.0 and T ¼ 0.3. It is 10
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Table 1 Spin configurations and the total polarization with the effects of D1 and E. -2.6 < D1
3.4 < D1 � 2.6
σziA ¼
5 z σ ¼ 2 jB
σziA ¼
5 z σ ¼ 2 jB
3 2
5 4
EC1 � E < EC2
σziA ¼
5 z σ ¼ 2 jB
1 3 SzmC ¼ 2 SznD ¼ 2 P ¼ 2 2
EC2 � E < EC3
σziA
5 z ¼ σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
7 P¼ 4
EC3 � E < EC4
σziA
5 z ¼ σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼2
EC4 � E < ES
5 z σ ¼ 2 jB 3 z ¼ σ ¼ 2 jB
σziA ¼ σziA
P¼1
SzmC ¼ 2 SznD ¼ 2
P¼
E < EC1
5 9 SzmC ¼ 2 SznD ¼ 2 P ¼ ES � E 2 4 3 SzmC ¼ 2 SznD ¼ 2 P ¼ 1 E < EC1 2
σziA ¼
3 z σ ¼ 2 jB
1 2
5 4
EC1 � E < EC2
σziA ¼
5 z σ ¼ 2 jB
1 3 SzmC ¼ 2 SznD ¼ 2 P ¼ 2 2
EC2 � E < EC3
σziA ¼
5 z σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
P¼
EC3 � E < EC4
σziA ¼
5 z σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼2
5 z σ ¼ 2 jB 1 z ¼ σ ¼ 2 jB
σziA ¼
D1 � 3.4
5 SzmC ¼ 2 SznD ¼ 2 2
σziA
SzmC ¼ 2 SznD ¼ 2
P¼
7 4
EC4 � E < ES
5 9 SzmC ¼ 2 SznD ¼ 2 P ¼ ES � E 2 4 1 z z SmC ¼ 2 SnD ¼ 2 P ¼ 1 E < EC1 2
σziA ¼
1 z σ ¼ 2 jB
1 2
5 4
EC1 � E < EC2
σziA ¼
3 z σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
P¼
3 2
EC2 � E < EC3
σziA ¼
3 z σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼
7 4
EC3 � E < EC4
σziA
5 z ¼ σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼2
EC4 � E < ES
σziA
5 z ¼ σ ¼ 2 jB
5 SzmC ¼ 2 SznD ¼ 2 2
9 P¼ 4
SzmC ¼ 2 SznD ¼ 2 P ¼
ES � E
seen that the twelvefold loops appear for smaller value of |J3| (¼0.1) in Fig. 14(a). This is because the weak exchange coupling |J3| make all spins of sublattices flip more easily so as to cause more loops. On increasing |J3| to 0.4 and 0.6, tenfold loops and triple loops appear. The appearance of triple loops can be explained as follows: on one hand, the spins are constrained by strong |J3| so as to flip difficultly in the small external electric field. On the other hand, only the external field increases large enough to make spins align along its direction and the existence of the ferrielectric exchange coupling between sublattices exists, which causes two small outside loops appear in the large external field. Until |J3| ¼ 0.7, single loop comes into our sight. It is found that the area of the hysteresis loop increases obviously. Similar hysteresis behaviors have been observed in the ferroelectric/antiferroelectric bilayer superlattice [19] and ferroelectric bilayer superlattice with strong interlayer exchange coupling [22]. Finally, the effect of temperature on the hysteresis loops of the total system is illustrated in Fig. 15 with fixed values D1 ¼ 2.2, D2 ¼ 0.6, J1 ¼ 1.8 and |J3| ¼ 0.3. From this figure, we can find that the multiple hysteresis loops behaviors can be weakened by the increasing of temperature. Namely, the number of the loops changes from quintuple to triple loops and single loop when T increases from 0.0001 to 4, meanwhile the area of the hysteresis loop also decreases distinctly. Finally, the loop vanishes and turns into the paraelectric state when T � 6.1. The results are consistent with previous studies of ferroelectric or ferrielectric materials [20,21,26–28, 60], antiferroelectric/ferroelectric BiFeO3/YMnO3 bilayer [56] and magnetic nanostructures [61–71]. In experiment, similar the electric hysteresis behavior under the influence of the temperature has been also observed in the single crystal of stoichiometric lithium niobite [72]. In addition, we also find similar hysteresis behavior under the influence of temperature in the MC study of dynamic hysteresis behavior in an external oscillating field for the simple ferromagnetic Ising system. It is found that the system changes from a paramagnetic state to a single-loop hysteresis behavior with temperature decreasing [73]. In order to better understand the effect of temperature on the coercivity Ehc and whether it satisfies the N� eel-Arrhenius theory [74], we plot the T and T1/2 dependences of the coercivity for the system with D1 ¼ 2.2, D2 ¼ 0.6, J1 ¼ 1.8 and |J3| ¼ 0.3, as shown in Fig. 16. The red balls represent the simulation results of the coercivity in Fig. 16(a) and (b), and the blue dotted line represents the fitted curve according to Eq. (4) in N�eel-Arrhenius theory as follows: � Ehc ¼ Ehc0 1
�
T TB
�1=2 �
(4)
11
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Table 2 Spin configurations and the total polarization with the effects of D2 and E. -3.4 < D2
5.2 < D2 � 3.4
5.6 < D2 � 5.2
σziA ¼
5 z σ ¼ 2 jB
σziA ¼
5 z σ ¼ 2 jB
3 2
SzmC ¼ 2 SznD ¼ 2
P¼
5 4
EC1 � E < EC2
σziA ¼
5 z σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
P¼
3 2
EC2 � E < EC3
σziA
5 z ¼ σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
7 P¼ 4
EC3 � E < EC4
σziA
5 z ¼ σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼2
EC4 � E < ES
σziA
P¼1
E < EC1
5 9 SzmC ¼ 2 SznD ¼ 2 P ¼ ES � E 2 4 5 1 SzmC ¼ 1 SznD ¼ 1 P ¼ E < EC1 2 2
σziA ¼
5 z σ ¼ 2 jB
5 SzmC ¼ 2 SznD ¼ 2 2
P¼1
σziA ¼
5 z σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼
5 4
EC2 � E < EC3
σziA ¼
5 z σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
P¼
3 2
EC3 � E < EC4
σziA ¼
5 z 1 7 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4
EC4 � E < EC5
σziA ¼
5 z 3 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 2 jB 2
EC5 � E < ES
5 z σ ¼ 2 jB 5 z 5 ¼ σ ¼ 2 jB 2
σziA
EC1 � E < EC2
5 9 SzmC ¼ 2 SznD ¼ 2 P ¼ ES � E 2 4 1 z z SmC ¼ 1 SnD ¼ 1 P ¼ E < EC1 2
σziA ¼
σziA ¼
5 z σ ¼ 2 jB
3 SzmC ¼ 1 SznD ¼ 1 2
P¼
3 4
EC1 � E < EC2
σziA ¼
5 z σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼
5 4
EC2 � E < EC3
5 z ¼ σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
3 P¼ 2
EC3 � E < EC4
σziA
6.0 < D2 � 5.6
5 z σ ¼ 2 jB 5 z ¼ σ ¼ 2 jB
σziA ¼
5 SzmC ¼ 2 SznD ¼ 2 2
σziA ¼
5 z 1 7 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4
σziA ¼
5 z 3 σ ¼ SzmC ¼ 2 SznD ¼ 2 2 jB 2
5 z σ ¼ 2 jB 5 z 5 ¼ σ ¼ 2 jB 2
σziA ¼ σziA
EC4 � E < EC5
P¼2
EC5 � E < ES
5 9 SzmC ¼ 2 SznD ¼ 2 P ¼ ES � E 2 4 1 SzmC ¼ 1 SznD ¼ 1 P ¼ E < EC1 2
σziA ¼
5 z σ ¼ 2 jB
3 SzmC ¼ 1 SznD ¼ 1 2
P¼
σziA ¼
5 z σ ¼ 2 jB
1 SzmC ¼ 1 SznD ¼ 1 2
P¼1
σziA ¼
5 z σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
P¼
σziA ¼
5 z 1 7 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4
EC4 � E < EC5
σziA ¼
5 z 3 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 2 jB 2
EC5 � E < ES
σziA
3 4
3 2
EC1 � E < EC2 EC2 � E < EC3 EC3 � E < EC4
5 z 5 9 ¼ σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4
ES � E
where Ehc0 is the coercivity near zero temperature and TB is the blocking temperature corresponding to the paraelectric state of hysteresis loop. Their values can be obtained from Fig. 15(a) and (d), namely TB ¼ 6.1 for the paraelectric state for and Ehc0 ¼ 2.63 at T ¼ 0.0001. From Fig. 16(a), we can see that the coercivity of the system depends strongly on temperature, and decreases with the increase of temperature. This behavior can be explained by the fact that the bilayer system has higher thermal energy at a higher temperature so that a smaller electric field is required to reverse the polarization. Similar behavior has been also observed in the experimental results of the CoxFe3-xO4 nanoparticles, which can be described by the N�eel-Arrhenius model [75]. From Fig. 16(b), it is clearly seen that the Ehc decreases almost linearly with T1/2 and the fitted curve also changes linearly, which satisfies the rule of Ehc ∝ T1=2 in N�eel-Arrhenius theory [74]. Similar results have been found in experimental studies of highly crystalline superparamagnetic 12
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Table 3 Spin configurations and the total polarization with the effects of J1 and E. -0.2 < J1
0.4 < J1 � 0.2
0.8 < J1 � 0.4
σziA ¼
5 z σ ¼ 2 jB
1 7 SzmC ¼ 2 SznD ¼ 2 P ¼ 2 4
σziA ¼
5 z σ ¼ 2 jB
3 2
P ¼ 2 EC1 � E < ES
5 z 5 9 SzmC ¼ 2 SznD ¼ 2 P ¼ σ ¼ 2 jB 2 4 3 z 1 3 z z ¼ σ ¼ SmC ¼ 2 SnD ¼ 2 P ¼ 2 jB 2 2
σziA ¼ σziA
5 z 1 σ ¼ SzmC ¼ 2 SznD ¼ 2 2 jB 2
P¼
σziA ¼
5 z 3 σ ¼ SzmC ¼ 2 SznD ¼ 2 2 jB 2
P ¼ 2 EC2 � E < ES
σziA σziA
7 4
ES � E E < EC1
σziA ¼
EC1 � E < EC2
5 z 5 9 ¼ σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4 5 z 3 5 SzmC ¼ 2 SznD ¼ 2 P ¼ ¼ σ ¼ 2 jB 2 4 1 SzmC ¼ 2 SznD ¼ 2 2
ES � E E < EC1
σziA
5 z ¼ σ ¼ 2 jB
σziA
5 z 1 7 ¼ σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4
σziA ¼ 2.0 � J1 � 0.8
SzmC ¼ 2 SznD ¼ 2
E < EC1
EC1 � E < EC2 EC2 � E < EC3
5 z 3 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 EC3 � E < ES 2 jB 2
5 z σ ¼ 2 jB 5 z 5 ¼ σ ¼ 2 jB 2
σziA ¼ σziA
3 P¼ 2
5 9 SzmC ¼ 2 SznD ¼ 2 P ¼ 2 4 SzmC ¼ 2 SznD ¼ 2 P ¼ 1
ES � E E < EC1
σziA ¼
5 z σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
P¼
5 4
EC1 � E < EC2
σziA ¼
5 z σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
P¼
3 2
EC2 � E < EC3
σziA ¼
5 z 1 7 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4
σziA ¼
5 z 3 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 EC4 � E < ES 2 jB 2
σziA ¼
5 z 5 9 σ ¼ SzmC ¼ 2 SznD ¼ 2 P ¼ 2 jB 2 4
EC3 � E < EC4
ES � E
Table 4 Spin configurations and the total polarization with the effects of |J3| and E. 0 � |J3| � 2.0
5 z σ ¼ 2 jB
5 SzmC ¼ 2 SznD ¼ 2 2
P¼1
5 z ¼ σ ¼ 2 jB
3 SzmC ¼ 2 SznD ¼ 2 2
5 P¼ 4
σziA ¼
5 z σ ¼ 2 jB
1 SzmC ¼ 2 SznD ¼ 2 2
P¼
σziA ¼
σziA ¼ σziA
E < EC1 EC1 � E < EC2
3 2
EC2 � E < EC3
5 z 1 σ ¼ SzmC ¼ 2 SznD ¼ 2 2 jB 2
P¼
7 4
σziA ¼
5 z 3 σ ¼ SzmC ¼ 2 SznD ¼ 2 2 jB 2
P¼2
σziA ¼
5 z 5 σ ¼ SzmC ¼ 2 SznD ¼ 2 2 jB 2
P¼
9 4
EC3 � E < EC4 EC4 � E < ES ES � E
iron oxide nanoparticles [76]. 4. Conclusion In this paper, we investigate the polarization plateaus and the hysteresis behaviors of the antiferroelectric/ferroelectric BiFeO3/ YMnO3 bilayer structure by using Monte Carlo simulation. The effects of the crystal-field, the exchange coupling and the external electric field on the polarization plateaus and the hysteresis behaviors of the system have been clarified. It is found that the multiple plateaus and the hysteresis behaviors reflect multiple spin configurations of the system, originating from the competition among the crystal-field, the exchange coupling and the external electric field. We hope that these results may be helpful for the experimental application of the antiferroelectric/ferroelectric bilayer.
13
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 11. The hysteresis loops of the antiferroelectric/ferroelectric Ising bilayer structure for different values of D1 with D2 ¼ 4, J1 ¼ 1.8, | J3| ¼ 0.3 and T ¼ 0.3.
Fig. 12. The hysteresis loops of the antiferroelectric/ferroelectric Ising bilayer structure for different values of D2 with D1 ¼ 2.2, J1 ¼ 1.8, | J3| ¼ 0.3 and T ¼ 0.3.
Declaration of competing interest We declare that we do not have any commercial or associative interest that represents a conflict of interest in connection with the work submitted. 14
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 13. The hysteresis loops of the antiferroelectric/ferroelectric Ising bilayer structure for different values of J1 with D1 ¼ 2.2, D2 ¼ 4, | J3| ¼ 0.3 and T ¼ 0.3.
Fig. 14. The hysteresis loops of the antiferroelectric/ferroelectric Ising bilayer structure for different values of |J3| with D1 ¼ 2.2, D2 ¼ 4, J1 ¼ 1.8 and T ¼ 0.3.
15
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
Fig. 15. The hysteresis loops of the antiferroelectric/ferroelectric Ising bilayer structure for different values of T with D1 ¼ 2.2, D2 ¼ 0.6, J1 ¼ 1.8 and |J3| ¼ 0.3.
Fig. 16. The T and T1/2 dependences of the coercivity for the antiferroelectric/ferroelectric Ising bilayer structure with D1 ¼ 2.2, D2 ¼ 0.6, J1 ¼ 1.8 and |J3| ¼ 0.3.
Acknowledgements This project was supported by the Natural Science Foundation of Liaoning Province, China (Grant No. 20170540672), “Xingliao Elite” Program Project of Liaoning Province (No. XLYC1807021), Youth Project of Liaoning Education Department (No. LQGD2017032).
16
Superlattices and Microstructures 136 (2019) 106293
Z.-y. Wang et al.
References [1] [2] [3] [4] [5] [6] [7] [8] [9] [10] [11] [12] [13] [14] [15] [16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28] [29] [30] [31] [32] [33] [34] [35] [36] [37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50] [51] [52] [53] [54] [55] [56] [57] [58] [59] [60] [61] [62] [63] [64] [65] [66] [67] [68] [69] [70] [71] [72] [73] [74] [75] [76]
S.M. Selbach, T. Tybell, M. Einarsrud, T. Grande, Adv. Mater. 20 (19) (2008) 3692. S.R. Das, R.N.P. Choudhary, P. Bhattacharya, R.S. Katiyar, P. Dutta, A. Manivannan, M.S. Seehra, J. Appl. Phys. 101 (2007), 034104. G. Lescano, F.M. Figueiredo, F.M.B. Marques, J. Schmidt, J. Eur. Ceram. Soc. 21 (2001) 2037. Z.J. Huang, Y. Cao, Y.Y. Sun, Y.Y. Xue, C.W. Chu, Phys. Rev. B 56 (1997) 2623. M. Ouaddari, S. Delprat, F. Vidal, M. Chaker, K. Wu, IEEE Trans. Microw. Theory Tech. 53 (2005) 1390. M. Liu, X. Li, H. Imrane, Y. Chen, T. Goodrich, Z. Cai, K.S. Ziemer, J.Y. Huang, N.X. Sun, Appl. Phys. Lett. 90 (2007) 152501. S. Corkovic, Q. Zhang, J. Appl. Phys. 105 (2009), 061610. Y. Zhang, Z. Li, C. Deng, J. Ma, Y. Lin, C.W. Nan, Appl. Phys. Lett. 92 (2008) 152510. R. Thomas, J.F. Scott, D.N. Bose, R.S. Katiyar, J. Phys. Condens. Matter 22 (2010) 423201. L. Goux, M. Gervais, F. Gervais, C. Champeaux, A. Catherinot, Int. J. Inorg. Mater. 3 (2001) 839. T. Atsuki, N. Soyama, T. Yonezawa, K. Ogi, Jpn. J. Appl. Phys. 34 (1995) 5096. T. Kijima, H. Matsunaga, Jpn. J. Appl. Phys. 38 (1999) 2281. J. Zapata, J. Narv� aez, W. Lopera, M.E. G� omez, G.A. Mendoza, P. Prieto, Transactions on Magnetics 44 (11) (2008) 2895. B. Wang, Y. Yang, W.J. Qian, L.l. Liu, Y.p. Wang, Symposium on Piezoelectricity 47 (1) (2012) 298. P.X. Nie, Y.P. Wang, Y. Yang, G.L. Yuan, W. Li, X.T. Ren, Energy harvesting, Systems 2 (3–4) (2015) 157. S.N. Tripathy, K.K. Mishra, S. Sen, B.G. Mishra, Dhiren K. Pradhan, R. Palai, D.K. Pradhan, J. Appl. Phys. 114 (2013) 144104. D.X. Zhou, G. Jian, Y.X. Hu, Y.N. Zheng, S.p. Gong, H. Liu, Mater. Chem. Phys. 127 (2011) 316. M. Alam, S. Talukdar, K. Mandal, Mater. Lett. 210 (2018) 80. A. Ainane, I. Essaoudi, M. Saber, J. Magn. Magn. Mater. 315 (2007) 132. Y. Benhouria, I. Essaoudi, A. Ainane, A. Oubelkacem, M. Saber, F. Dujardin, Phys. Scr. 86 (2012), 045704. Y. Benhouria, I. Essaoudi, A. Ainane, R. Ahuja, F. Dujardin, Chin. J. Phys. 54 (4) (2016) 533. A. Feraoun, A. Zaim, M. Kerouad, Solid State Commun. 248 (2016) 88. Y. Benhouria, I. Essaoudi, A. Ainane, R. Ahuja, F. Dujardin, J. Supercond. Nov. Magnetism 27 (9) (2014) 2153. Y. Benhouria, I. Essaoudi, A. Ainane, R. Ahuja, F. Dujardin, Superlattice Microstruct. 73 (2014) 121. Y. Benhouria, I. Essaoudi, A. Ainane, R. Ahuja, F. Dujardin, J. Supercond. Nov. Magnetism 27 (9) (2014) 2153. Y. Benhouria, I. Essaoudi, A. Ainane, R. Ahuja, F. Dujardin, Physica A 506 (2018) 499. A. Feraoun, A. Zaim, M. Kerouad, J. Phys. Chem. Solids 96 (2016) 75. S. Zriouel, B. Taychour, F. El Yahyaoui, L.B. Drissi, Physica 91 (2017) 113. A. Jabar, R. Masrour, A. Benyoussef, M. Hamedoun, J. Supercond. Nov. Magnetism 30 (7) (2017) 1807. A. Jabar, R. Masrour, Phys. B Condens. Matter 539 (15) (2018) 21. J.L. Bi, W. Wang, Q. Li, Superlattice Microstruct. 107 (2017) 104. W. Wang, W. Jiang, D. Lv, F. Zhang, J. Phys. D Appl. Phys. 45 (2012) 475002. W. Wang, W. Jiang, D. Lv, Phys. Status Solidi B 249 (1) (2012) 190. W. Wang, D. Lv, F. Zhang, J.L. Bi, J.N. Chen, J. Magn. Magn. Mater. 385 (2015) 16. R. Masrour, A. Jabar, Physica A 491 (2018) 926. M. Ertas¸, M. Keskin, Phys. Lett. A 376 (2012) 2455. M. Ertas¸, M. Keskin, B. Deviren, Physica A 391 (2012) 1038. R. Masrour, L. Bahmad, A. Benyoussef, J. Magn. Magn. Mater. 324 (23) (2012) 3991. Y. Benhouria, I. Bouziani, I. Essaoudi, A. Ainane, R. Ahuja, J. Magn. Magn. Mater. 460 (2018) 223. R. Masrour, A. Jabar, J. Comput. Electron. 16 (1) (2017) 12. R. Masrour, L. Bahmad, A. Benyoussef, M. Hamedoun, E.K. Hlil, J. Supercond. Nov. Magnetism 26 (2013) 679. W. Jiang, L. Liu, X. Li, Q. Deng, H. Guan, F. Zhang, A. Guo, Physica B 407 (2012) 3933. A. Feraoun, M. Kerouad, Solid State Commun. 280 (2018) 56. A. Feraoun, M. Kerouad, Solid State Commun. 277 (2018) 25. A. Feraoun, M. Kerouad, Computational Condensed Matter 15 (2018) 7. M. Ertas¸, Superlattice Microstruct. 85 (2015) 734. A. Jabar, R. Masrour, A. Benyoussef, M. Hamedoun, J. Magn. Magn. Mater. 397 (2016) 287. R. Masrour, A. Jabar, A. Benyoussef, M. Hamedoun, J. Magn. Magn. Mater. 393 (2015) 151. R. Masrour, E.K. Hlil, M. Hamedoun, A. Benyoussef, O. Mounkachi, H. Moussaoui, Physica A 395 (2014) 128. R. Masrour, E.K. Hlil, M. Hamedoun, A. Benyoussef, O. Mounkachi, L. Bahmad, J. Magn. Magn. Mater. 326 (2013) 166. R. Masrour, A. Jabar, A. Benyoussef, M. Hamedoun, J. Magn. Magn. Mater. 401 (2016) 700. A. Jabar, R. Masrour, Solid State Commun. 268 (2017) 38. A. Jabar, R. Masrour, Chem. Phys. Lett. 700 (2018) 130. A. Jabar, R. Masrour, A. Benyoussef, M. Hamedoun, Chem. Phys. Lett. 670 (2017) 16. R. Masrour, A. Jabar, J. Magn. Magn. Mater. 426 (2017) 225. A. Jabar, R. Masrour, A. Benyoussef, M. Hamedoun, J. Supercond. Nov. Magnetism 29 (3) (2016) 733. Z.Y. Wang, W. Wang, Q. Li, M. Tian, Z. Gao, Y. Liu, Physica 110 (2019) 127. N. Metropolis, A.W. Rosenbluth, M.N. Rosenbluth, A.H. Teller, E. Teller, J. Chem. Phys. 21 (1953) 1087. W. Wang, Q. Li, D. Lv, R. Liu, Z. Peng, S. Yang, Carbon 120 (2017) 313. Y. Wu, D. Yao, Z. Li, J. Appl. Phys. 91 (2002) 1482. D. Lv, F. wang, R.J. Liu, Q. Xue, S.X. Li, J. Alloy. Comp. 701 (2016) 935. D. Lv, Y. Yang, W. Jiang, F. Wang, Z.Y. Gao, M. Tian, Physica A 514 (2019) 319. W. Wang, S.Q. Yang, Y. Yang, Z. Peng, B.C. Li, M. Yang, Physica 109 (2019) 30. Y. Yang, W. Wang, H. Ma, Q. Li, Z.Y. Gao, T. Huang, Physica 108 (2019) 358. H.J. Wu, W. Wang, B.C. Li, M. Yang, S.Q. Yang, F. Wang, Physica 112 (2019) 86. D. Lv, W. Jiang, Y. Ma, Z.Y. Gao, F. Wang, Physica 106 (2019) 101. D. Lv, W. Wang, J.P. Liu, D.Q. Guo, S.X. Li, J. Magn. Magn. Mater. 465 (2018) 348. W. Wang, D.D. Chen, D. Lv, J.P. Liu, Q. Li, Z. Peng, J. Phys. Chem. Solids 108 (2017) 39. S.Q. Yang, W. Wang, F. Wang, B.C. Li, J.H. Xu, J. Phys. Chem. Solids 135 (2019) 109110. W. Wang, Y. Liu, Z.Y. Gao, X.R. Zhao, Y. Yang, S. Yang, Physica E101 (2018) 110. W. Wang, J.L. Bi, R.J. Liu, X. Chen, J.P. Liu, Superlattice Microstruct. 98 (2016) 433. H. Bo, Q. Meng, H. Hu, H. Zhao, Z. Zhang, Q. Zhang, C. Zhang, Appl. Phys. A 124 (2018) 691. K. Bikas, M. Acharyya Chakrabarti, Rev. Mod. Phys. 71 (1999) 847. L. N�eel, Ann. Geophys. 5 (1949) 99. T.E. Torres, E. Lima, A. Mayoral, A. Ibarra, C. Marquina, M.R. Ibarra, G.F. Goya, J. Appl. Phys. 118 (2015) 183902. M. Tadi�c, V. Kusigerski, D. Markovi�c, M. Panjan, I. Milo�sevi�c, V. Spasojevi�c, J. Alloy. Comp. 525 (2012) 28.
17