Electromechanical Coupling Multiaxial Experimental and Micro-Constitutive Model Study of Pb(Mg1/3Nb2/3)O3- 0.32PbTiO3 Ferroelectric Single Crystal

Compared to polycrystalline ferroelectric ceramics such as Pb(Zr1/2Ti1/2)-O3 (PZT), domain en‐ gineered relaxor ferroelectric single crystals Pb(Zn1/3Nb2/3)O3-xPbTiO3 (PZN-xPT) and Pb(Mg1/3Nb2/3)O3-xPbTiO3 (PMN-xPT) show greatly enhanced electromechanical properties: the piezoelectric coefficient d33 and electrically induced strain of <001> oriented single crystals direction can respectively reach 2500pC/N and 1.7%, or even greater [1],[2]. PZN-xPT and PMN-xPT with outstanding properties are usually near the morphotropic phase boundary (MPB), separating rhombohedra (R) and tetragonal (T) phases [1], [4], [5]. For example, the MPB is located within a range of x=0.275-0.33 for PMN-xPT [6], [7]. Recent investigations have also shown that the MPB of PMN-xPT and PZN-xPT is actually in a multi-phase state at room temperature [8], [9]. e.g., the coexistence of rhombohedra (R) and monoclinic (M) phases at x=0.33 or rhombohedra and tetragonal phases at x=0.32. Other phases can also exist near MPB under stress and/or electric field loading. Experiment studies [10], [11], [12] and first-princi‐ ples calculations [4], [13]reveal that two different homogeneous polarization rotation path‐ ways are present between rhombohedra and tetragonal phases under an electric field. It is found that three kinds of intermediate monoclinic phases MA, MB and MC are associated with the two pathways [9], [10]. Orthorhombic phase has also been observed between the rhombo‐ hedra and tetragonal phases when PMN-PT is loaded with strong electric field along the <110> direction [9],[14]. Although the mechanism underlying the high performance of these crystals is not completely clear, the existence of the intermediate phases and the associated phase tran‐ sitions are believed to be one of the main reasons [1], [3], [4], [8], [10].

Note that PZN-PT and PMN-PT single crystals usually experience electric and/or mechanical loading during their in-service life.It has been shown that externally applied loading has significant effect on the properties of these crystals [1], [5]- [28].A number of studies have focused on the loading induced behavior of <001> and <110> oriented anisotropic PZN-xPT and PMN-xPT crystals [1], [5], [9], [15]- [19].It is has been shown loading in the form of electric field [1]- [21], and stress [21]- [26], can lead to polarization switching and phase transitions, which changes the crystal phase and domain structure of these single crystals, and hence dramatically alter their electromechanical properties.A mature level understanding of their responses to electrical, mechanical and temperature loading condition is thus essential to fulfill the applications of these crystals.
Rhombohedral phase of PZN-xPT(0<x<0.1) and PMN-xPT(0<x<0.35) exhibit excellent electromechanical properties along {001} and {110} orientations compared to the {111}(cubic cell reference) spontaneous polarization direction.A number of studies have focused on the loading induced behavior of {001}, {011} and {111} oriented PZN-xPT and PMN-xPT single crystals.It has been shown that loadings in the form of electric field and stress can lead to polarization rotation and phase transition, which change the crystal phase and domain structure of these single crystals, and hence dramatically alter their electromechanical properties.When an electric field is applied along the {001} direction of PZN-PT, polarization rotation occurs from {111} towards {001} via either M A or M B , depending the composition (e.g., the rotation is R-M A -T for PMN-4.5PT).For {110}-oriented PMN-PT, polarization rotates from {111} to {110} via M B when electric field is along {110} and ends up with an orthorhombic phase when the electric field exceeds a critical value.These single crystals can be polarized into a single domain state with {111} oriented electric field.Most available studies on the effect of bias stress on the crystal behavior focus only on uni-polar electric field loading.However, it has been shown that the uni-polar and bi-polar responses of these crystals can be very different.
A number of studies have focused on experimental, but the constitutive model of ferroelectric single is absent.Huber [29] and Bhattacharya [30], [31] et.al established a model based on the micromechanical method, but the simulation is not well compare to experimental result.
At present, a mature model to explain the stress-strain behavior of ferroelectric single is absent all along.In this study, electric field induced "butterfly" curves and polarization loops for a set of compressive bias stress of {001}, {011} and {111} poled PMN-0.32PTsingle crystals will be explored by systematical experiment study.The effects of the compressive bias stress on the material properties along these three crystallographic directions of PMN-0.32PTsingle crystals will be quantified.The underlying mechanisms for the observed feature will be explained in terms of phase transformation or domain switching,depending on the crystallographic direction.The stress-strain curves along <001> crystallographic direction of ferroelectric single crystals BaTiO 3 will be calculated in the first principle method to validate polarization rotation model.Finally, Based on the experimental phase transformation mechanism of ferroelectric single crystal, a constitutive model of ferroelectric single crystal is proposed based on micromechanical method.This constitutive model is facility and high computational efficiency.

Experimental methodology
At the test room temperature, PMN-0.32PTsingle crystals used in this study are of morphotropic composition, and in the rhombohedral phase, very close to MPB.The pseudocubic {001}, {011} and {111} directions of these crystals are determined by x-ray diffraction (XRD).Pellet-like specimens of dimensions 5 × 5 × 3mm 3 are then cut from these crystals, with the normal of the 5 × 5mm 2 major specimen surfaces along the pseu- do-cubic {001}, {011} or {111} direction.All specimens are electroded with silver on the 5 × 5mm major surfaces and poled along the {001}, {011} and {111} orientations (i.e., the specimen thickness direction) under a field of 1.5kV/mm.Note that there are eight possible dipole orientations along the body diagonal directions of unpoled PMN-0.32PTsingle crystals (i.e., the <111> direction).When an electric poling field is applied to the crystals along the {001} direction, a multi-domain structure can be produced, comprising four degenerate states.For {011} direction poled single crystals, the number of degenerate states is two.The single crystals can be poled into a single domain state when they are poled along the {111} direction.The states of {001}, {011} and {111} poled PMN-0.32PTsingle crystals are sketched in Fig. 1, where solid arrows refer to the domain states induced by poling (also labeled as type 1 in Figs.1b and 1c), dotted arrows represent possible domains switched from type 1 domains upon loading or unloading, and vice versa (also labeled as type 2, ).Since the focus of this study is to explore the effect of bias stress on the electromechanical properties of PMN-0.32PTsingle crystals along different crystallographic directions, experimental setup is adapted from Ref. [31] (Fig. 2) to allow simultaneously imposing uniaxial stress and electric field to the specimen along the thickness direction.Mechanical load is applied by a servo-hydraulic materials test system (MTS) and electric field is applied to the specimen using a high voltage power amplifier.Once the specimen is placed in the fixture, a compressive bias stress with magnitude of at least 0.4MPa is maintained throughout the test to ensure electrical contact.Stress controlled loading instead of displacement controlled loading is adopted during the test, so that the specimen is not clamped but is free to move longitudinally when electric field E 3 and mechanical stress σ 33 loading are applied, where subscript 3 refers to the thickness direction of the samples, corresponding to the {001}, {011} and {111} direction of {001}, {011} and {111} oriented PMN-0.32PTsingle crystals, respectively.During test, polarization P 3 (or electric displacement) is measured using a modified Sawyer-Tower bridge, and the deformation (i.e., the strain) is monitored by two pairs of strain gauges (in total four strain gauges used) mounted on the four 5 × 3mm 2 surfaces: one pair placed on two opposite 5 × 3mm 2 surfaces is applied to measure the longitudinal normal strain ε 33 along the {001}, {011} and {111} direction and another pair to measure the transverse normal strain ε 11 in the direction perpendicular to the {001}, {011} and {111} direction respectively.Output of the strain gauges during deformation is recorded by a computer through a multiple-channel analog-to-digital (AD) converter.
The first set of tests is performed for electric field loading of triangular wave form of magnitude 0.5kV/mm and frequency 0.02Hz, free of stress loading.The low frequency is chosen to mimic quasistatic electric loading, which is of particular interest in this study [32].Unless stated otherwise, this loading frequency for the electric field is used throughout the following test.The second set of tests consists of mechanical loading upon short circuited samples.The samples are compressed to -40MPa and unloaded to -0.4MPa at loading and unloading rate of 5MPa/min, followed by an electric field which is sufficiently large to remove the residual stress and strain to re-polarize the samples.In the third set of tests, triangular wave form electric field is applied to the samples which are simultaneously subject to co-axial constant compressive stress preload.The magnitude of the preload is varied from test to test and is in the range between 0 and -40MPa.Note that there is a time-dependent effect of the depolarization and strain responses under constant compressive stress.To minimize this effect, each electric field loading starts after a holding time of 150 seconds for a new stress preloading.It is found that three cycles of electric field loading and unloading are sufficient to produce stabilized response for each constant prestress, and the results for the last cycle are reported in the following.

Crystallographic dependence of electric behavior and piezoelectric properties
The measured electric field induced polarization hysteresis loops and butterfly curves for {001}, {011} and {111} oriented poled PMN-0.32PTsingle crystals without stress loading are shown in Figs.3a and 3b.The remnant polarizations P r (defined as the polarization value at zero electric field), the coercive electric fields E c (defined as the electric field value at zero polarization) and the piezoelectric coefficients d 33 (defined as are 1828, 1049 and 200pC/N, respectively.It is noticed that there are eight possible polariza-tion orientations along the pseudo-cubic {111} for un-poled rhombohedral PMN-0.32PTsingle crystals.Upon poling, the dipoles switch as close as possible to the applied electric field direction: For {001} poled crystals, there are four equivalent polar vectors along the {111} orientation, with an inclined angle of -54.7º from the poling field (Fig. 1a); For {011} poled crystals, there are two equivalent polar vectors along the {111} direction (labeled as type 1 in Fig. 1b); For {111} poled crystals, there is one polar vector along {111} (type 1 in Fig. 1c).According to the domain configurations in Fig. 1, the remnant polarizations P r 001 and P r 011 are approximately related to P r 111 by P r 001 = P r 111 / 3 and P r 011 = 2P r 111 / =0.227C/m 2 and P r 011 =0.324C/m 2 ).Therefore, the measured results are consistent with the domain configurations shown in Fig. 1.
It is also seen from Figs. 3a and 3b that the coercive field is lowest for the {111} oriented crystals, and becomes successively higher for {001} and {110} orientations.This is same to Ref.25 except for {110} orientation.This trend of coercive field is due to two reasons: One is due to reorientation driving force being proportional to the component of electric field aligned with the rhombohedral direction; The other one is due to the domain switching process.In {111}oriented PMN-0.32PTsingle crystals, there are two types of domains (shown as type 1 and type 2 in Fig. 1c).When the electric field is decreased from 0.5kV/mm to -0.05kV/mm, the strain first decreases linearly (see Fig. 3b).When the electric field is decreased further, the type 1 domain switches to type 2 domains, leading to abrupt displacement change.When the electric field exceeds the coercive field -0.216kV/mm, the type 2 domains switch back to type 1 domain, recovering the deformation.In type 2 domain state, three equivalent polar vectors with an angle of 71º from the {111} direction can coexist and are separated by domain walls across which the normal components of electric displacement and displacement jump are zero.Ideally, this type of domain walls has no associated local stress or electric field.So the existing of type 2 domain state and the largest component of electric field along polarization direction induces the lowest coercive field in {111} orientation poled crystal.In the {011} orientation crystals, there are also two types of domains (i.e., type 1 and type 2 domain in Fig. 1b), the domain switching process is similar to {111}-poled crystals.In the type 2 domain state, however, the four possible polar vectors are perpendicular to the applied electric field (Fig. 1b).It is thus difficult to switch type 2 domain to type 1 domain only by applying electric field.Both the type 2 domain state and the smallest component of E contribute to the largest coercive field in {011} oriented crystals.In {001} oriented crystals, there will be no associated local stress or electric field, similar to type 2 domain state of {111} poled crystals.This feature again renders domain switching easy.On the other hand, the domain structure of {001} poled crystals is stable [4] and the component of electric field is smaller, giving the coercive field higher than that of {111} poled crystals but lower than that of {011} oriented crystals (see, Figs. 3a and 3b).Some researchers attribute the high piezoelectric coefficients along {001} and {011} oriented ferroelectric single crystals to the engineered domain state [4,11].It has also reported that, however, the piezoelectric coefficient along the {001} direction of single crystals with monodomain structure is comparable to that of crystals with multi-domain structure [17], implying the origins of the high piezoelectric constants of PMN-0.32PTsingle crystals may not be due to the engineered domain state.Instead, it could be due to the effect of crystal lattice properties.To further explore this issue, we follow Ref. [17] to calculate the piezoelectric co-  =1828 pC/N).The small discrepancy between the predictions and measurements is believed to be due to the fact that the predictions are based on ideal mono-domain single crystals while the material parameters used in Eq. ( 1) are actually from less ideal mono-domain crystals (in fact, they are more or less multi-domain crystals).Nevertheless, we can conclude that the dominant contribution to the large {001} and {011} piezoelectric response should be the crystal anisotropy other than engineered domain state.

Crystallographic dependence of stress induced strain and polarization responses
Figure 4 shows the measured σ 33 − ε 33 and σ 33 − P 3 curves for {001}, {011} and {111} oriented short circuited samples.The results indicate obvious crystallographic anisotropy in stress inducing responses.From domain switching viewpoint, there should be no significant deformation in the thickness dimension of {001} poled crystal samples under stress loading, except for the elastic strain.However the longitudinal strain of {001} is contractive with a maximum magnitude of 0.3% under -40MPa.The contractive strain of {001} oriented crystals is about twice as much as those of {011} (about -0.18%) and {111} (about -0.17%) oriented crystals under a loading of -40MPa (Fig. 4a).This abnormal behavior of {001}-oriented crystals lies in that the mechanism underlying the stress induced response of {001}-oriented crystals is the R to O and T phase transition (see, Fig. 1a) rather than domain switching for {011} and {111} oriented crystals: the lattice distortion due to phase transition induces large deformation in the thickness dimension of {001} poled crystals.In {011} and {111}-poled crystals, on the contrary, the type 2 multi-domain state induced by compressive stress is the stable and preferred state (Figs.1b and 1c), and can form more easily by domain switching than phase transformation.Therefore, the stress induced strain and polarization curves in {011} and {111} orientation crystals are similar to those of ferroelectric polycrystals of which domain switching is also the dominant deformation mechanism (Figs.4a and 4b).During the unloading of {001} oriented crystals (corresponding to the returning to R phase of the unstable O phase), the σ 33 − ε 33 and σ 33 − P 3 curves show obvious nonlinear behavior.Note that the remnant strain and polarization at the end of unloading are attributed to the stable T phase which does not switches back to the R phase.For {011} and {111}-oriented crystals, however, there is only domain switching (i.e., switching between type 1 and type 2 domains) and no phase transformation occurs.During unloading, the σ 33 − ε 33 curves and σ 33 − P 3 curves of {011} and {111} show linear response since almost no domain switches back (Figs.4a and 4b).The stress cycles of {001} poled crystals can be explained by the polarization vector rotation mechanism sketched in Fig. 5 as follows.In general two most possible mechanisms responsible for the observed features of {001} poled crystals in Figs.4 are domain switching and polarization rotation associated phase transformation.Recall that the PMN-0.32PTsingle crystals considered here are in the R phase close to MPB (although in reality there could also be M or T phase in these crystals near MPB [8], [9] the R phase is nevertheless the dominant phase).It is noted that that, upon application of a field along the <001> poling axis of R phase domain engineered PMN-0.32PTsingle crystals, only four of the eight polarization orientations are possible, i.e., <111>, <1 ¯11>, <11 ¯1>, and <1 ¯1 ¯1>.Since the <001> components of these four polar vectors are completely equivalent, each domain wall cannot move under an external electric field along the <001> direction owing to the equivalent domain wall energies [32].In other words, no ferroelectric domain switching is possible.Equally, no ferroelastic domain switching is expected when a compressive stress is applied along the <001> axis.(For PMN-PT systems with orthorhombic 4O <001> domain engineered structure [33], 60º switching from <011> to <110> would be driven by an applied stress.But, this is not what we are considering).Such facts partly rationalize our hypothesis that the underlying mechanism associated with the observed behavior of the Rphased PMN-0.32PTsingle crystals considered in this paper is phase transformation.Nevertheless, it should be emphasized that since no in situ diffraction observation is conducted to properly identify the stress induced phase transformation in PMN-0.32PTsingle crystals, discussions in this paper on phase transition are solely inferred from the measured stress induced strain and polarization curves.).When the presence of an electric field along <110> direction, polarization rotation only occurs from R to O under compression along <001> direction.In the absence of electric field along <110>, however, almost all the polarization vectors may switch to T and O phases when the compression exceed 30MPa in magnitude (Fig. 5(c)).When the magnitude of the compressive loading is further increased, there is no more polarization rotation and the PMN-0.32PTsingle crystal shows linear response in Figs. 4. It is noted that O phase is usually unstable (its free energy balance between R and T phases depends on the electric and mechanical loading history [9], [14]).Upon unloading, the O phase switches back to R phase when the compression lower than -10MPa (Figs. 4, and 5(d)).On the other hand T phase is a stable phase.It will not switch back to R phase upon unloading, leading to remnant strain of about -0.07%for and remnant polarization of about -0.13C/m 2 at the end of the -40MPa stress cycle (Figs.4(a) and 4(b), Fig. 5(e)).After the stress cycle, a unipolar electric field of 0.5kV/mm in magnitude is applied to re-pole the single crystal to its initial state, i.e., the polarization vector completely switches back.As a result, the remnant strain and polarization diminish to zero at the end of the re-polarization process (Figs.5(a) and 5(e)).

Crystallographic dependence of electric field induced behavior at constant bias compressive stress
The ε 33 − E 3 "butterfly" curves and P 3 − E 3 hysteresis loops of {001}, {011} and {111} poled crystals under different constant compressive bias stresses are shown in Fig. 6.Dependence of the electric coercive E c , remnant polarization P r , dielectric permittivity χ 33 , piezo- electric coefficient d 33 and aggregate strain Δε for {001}, {011} and {111} oriented crystals on the compressive bias stress is summarized in Fig. 7.The aggregate strain Δε is defined as the difference between the maximum and minimum strain for a complete responsive butterfly curves.Similar to d 33 , χ 33 are calculated by ΔP 3 / ΔE 3 where ΔP 3 is the polariza- tion difference between -0.05 and +0.05kV/mm.The calculated χ 33 within such a small field range is almost equal to the slope of P 3 − E 3 hysteresis loops when the electric field passes through zero.The calculated χ 33 includes both the reversible (intrinsic dielectric property) and irreversible (extrinsic domain switching and phase transformation related property) contributions of the material, which is generally higher than the permittivity measured by a dynamic method [32].
The influence of the preloaded compressive stress on the aggregate strain Δε, remnant polarization P r and piezoelectric coefficient d 33 seems to be similar for {001}, {011} and {111} oriented crystals: The remnant polarization P r decreases with increasing the magnitude of the compressive prestress; The aggregate strain Δε and piezoelectric coefficient d 33 first increase and then decrease with the magnitude of the prestress increasing.As suggested earlier, however, the underlying mechanisms for the electromechanical behavior of {001} {011}, and {111} oriented crystals are different.The change of the dielectric permittivity χ 33 {001} oriented crystals under compressive stress is similar to that of χ 33 induced by temperature: near the phase transformation temperature there is a peak in the χ 33 curves, As is shown in Fig. 7c, there is a small peak near -6MPa and a big peak near -20MPa in χ 33 curves for {001} oriented crystals.This may be due to R-M (near -6MPa) and M-R and M-O phase transformation (near -20MPa).On the other hand, there is no obvious peak in χ 33 curves for {011} and {111}-oriented crystals (Fig. 7c).The change of χ 33 curves induced by compressive stress and the aforementioned stress induced strain and polarization responses are consistent with the hypothesis that there is phase transformation for {001}-orientated crystals under compressive stress.
As suggested by Fu et al. [7], for {001} oriented crystal origins of large aggregate strain and high piezoelectric coefficient at a moderate compressive bias stress (i.e., around -20MPa for the single crystals considered here) may be attributed to the stress induced intermediate states between rhombohedral and tetragonal phases.Under a bias stress of about -20MPa, PMN-0.32PTsingle crystals, after a phase transformation, are in a state of monoclinic phase which has a larger c / a (c and a are lattice parameters) and a smaller polarization component along the field direction than those of rhombohedral phase [20], implying electric field induced greater aggregate strain and piezoelectric coefficient (Figs.7d and 7e).With the magnitude of the compressive bias stress increased further (i.e., σ 33 =-30 and -40MPa), the {001} oriented single crystals are in the state of a mixture of orthorhombic and tetragonal phases.
Although O and T phases also have large c to a ratio, the presence of a large compression has an opposite effect (i.e., preventing larger deformation induced by electric field) and results in small aggregate strain (Fig. 7e).
As is noted, the polarization rotation introduced by compression gives rise to θ in d 33 * (ϕ, θ) in the range between -54.7º and -90º (the angle between loading and polarization direction of T or O phase) under -20MPa, and results in larger d 33 * in accordance with Eq. ( 1) (see in Fig. 7d).When the magnitude of the compressive bias stress increases further, θ approaches 90º and d 33 * becomes smaller.Meanwhile, the remnant polarization P r and coercive field E c decrease monotonically with the applied compressive bias stress because of the decreased component polarization along the {001} direction under compression (Fig. 7b).
Under zero stress, the initial state of {011}-oriented crystals is of multi-domain with two equivalent polarization directions (Fig. 1b).For {111}-oriented crystals, the initial state is of mono-domain with polarization direction in the {111} direction (Fig. 1c).When the applied electric field decreases from 0.5kV/mm to -E c , the strain and polarization decrease linearly.
Upon approaching -E c , domain state switches to four polar domain state for {011} oriented crystals (type 2 in Fig. 1b) and three polar domain state for {111} oriented crystals (type 2 in Fig. 1c), giving the jumps in strain and polarization responses.When electric field exceeds -E c , domain complete second switching, from four polar domain state to two polar domain state for {011} oriented crystals and from three polar domain state to mono-domain state for {111} oriented crystals.When the electric field reaches -0.5kV/mm, the strain and polarization change linearly again (Fig. 7).
Note that compressive stress can induce domain switching in {011} and {111}-orientated crystals.When compressive stress is superimposed on the samples, the shapes of ε 33 − E 3 and P 3 − E 3 curves are different from those under the zero stress state (Fig. 6c-f).Due to compres- sion inducing depolarization, the remnant polarization decreases in {011} and {111} oriented crystals (Fig. 7b).A low compressive stress (e.g., -5MPa for {011}, -7MPa for {111} ) can lead to type 1 to type 2 domain switching.As a result, more domains take part in switching during the electric field loading cycle and larger aggregate strainΔε than under zero compressive stress is observed.This partly explains the observed features in Fig. 7e.

Calculation methodology
PMN-PT and BaTiO 3 have the similar ABO 3 structure, so they have the similar ferroelectric properties.There is only Ti 4+ particle in B site of BaTiO 3 , however, there is not only Ti 4+ but also minim Mn 4+ and Ni 2+ particle in B site of PMN-PT.So the single cell of BaTiO 3 is convenient in calculation and it keeps the similar ferroelectric to PMN-PT.
In this paper, we calculate a single cell of BaTiO3, the single cell should be the smallest periodic reduplicate cell, it includes one Ti 4+ particle, three O 2-particles and one Ba 2+ particle (Fig. 8).It is suggested that the initialized state of BaTiO 3 ferroelectric single crystal is R phase after {001} oriented polarization.Refer to literature [34], the crystal lattice constant of R phase BaTiO 3 ferroelectric single crystal is 4.001A, the coordinate of particle in single cell is shown in table 1.The loading along {001} direction is carried out through application increasing strain by degrees.Stress and other parameters under each strain level is calculated by VASP.Calculation under each strain level include two steps.Firstly, the particle coordinate of single cell under each strain level is calcluated by first principle molecular dynamic method.Secondly, Stress and other parameters are obtained through relaxation that is based on the first result.Each increment of strain in this paper is 0.5%, untill 4%. 4 The first principle calculation of stress-strain 4.1 Calculation methodology PMN-PT and BaTiO 3 have the similar ABO 3 structure, so they have the similar ferroelectric properties.There is only Ti 4+ particle in B site of BaTiO 3 , however, there is not only Ti 4+ but also minim Mn 4+ and Ni 2+ particle in B site of PMN-PT.So the single cell of BaTiO 3 is convenient in calculation and it keeps the similar ferroelectric to PMN-PT.
In this paper, we calculate a single cell of BaTiO3, the single cell should be the smallest periodic reduplicate cell, it includes one Ti 4+ particle, three O 2-particles and one Ba 2+ particle(Fig.8).It i suggested that the initialized state of BaTiO 3 ferroelectric single crystal is R phase after {001} oriented polarization.Refer to literatureError!Reference source not found., the crystal lattice constant of R phase BaTiO 3 ferroelectric single crystal is 4.001  A , the coordinate of particle in single cell is shown in table 1.The loading along {001} direction is carried out through application increasing strain by degrees.Stress and other parameters under each strain level i calculated by VASP.Calculation under each strain level include two steps.Firstly, the particle coordinate of single cell under each strain level is calcluated by first principle molecula dynamic method.Secondly, Stress and other parameters are obtained through relaxation that i based on the first result.Each increment of strain in this paper is 0.5%, untill 4%.

Discussion
In order to validate that the method a mentioned in this paper is correct, we calculate the elasticity constant C 33 of T phase BaTiO 3 using the method mentioned before.Fig. 9 shows that the C 33 is 180GPa, which is approach to the experimental value 189± 8Gpa reported in refer [34].This means that the method used in this paper is trusty. -4  33 (%) Stress-strain curve of {001} orientated R phase BaTiO 3 calculated in first principle method is shown in Fig. 10(a).Fig. 10(a) shows that the result has simlar nonlinear behavior to the experimental result of PMN-0.32PT during the loading that sketched in Fig. 4(a), namely, there is obvious "a,b,c" steps during loading.We know that the nonlinear behavior of PMN-0.32PTshown in Fig. 4(a)should be polarization rotation (R→M→O and R→M→T).
The R→M→O is corresponding to the processing of polarization vector P R switching to P O , which is shown in Fig. 10(b).table 2 is calculated coordinate of particle in BaTiO 3 , the coordinate is correspond to point A in fig.10(a).Table 3 is coordinate of particle in O phase BaTiO 3 ferroelectric single crystal from refer [35].Data in table 2 are equal to those in table 3, which indicates that "A" point in fig.10(a) should be O phase.With increasing strain, the coordinates of particle are unchangeable after "A", this indicates BaTiO 3 ferroelectric single crystal is stable in O phase after "A".From the first principle calculation, we testify that the PMN-0.32PTferroelectric single crystal undergoes polarization rotation(R→M→O), this prove the polarization rotation model is reasonable.

The viscoplastic model
From the experimental analysis, it is known that the <001>-oriented PMN-0.32PTsingle crystal undergo the R→M→T and R→M→O phase transformation under the compression.In order to establish a compact model and keep the essence of experiment, R→O phase transformation of PMN-0.32PT is considered in this study.Polarization rotation cause bulk deformation and shear deformation associated with the slip planes.It is assumed that the corresponding deformation associated with the slip planes is shear dominated, a feature similar to that of the multi-slip system of the crystal plasticity.This similarity renders possible to use crystal plasticity models to describe the transformation deformation of PMN-0.32PTsingle crystal.
The poled ferroelectric single is four domain state, the polarization vector is along <111>, <1 11>, <11 ¯1> and <1 ¯1 ¯1> respectively.In this study, it is assumed that the polarization vector of R phase can switch to the polarization vector of O phase, such as <111> vector switches to <001>, the vector also can switch back from O phase to R phase.The polarization rotation has an analogy to crystal plasticity slip, so we suggest there is eight slip systems in PMN-0.32PTferroelectric single crystal.According to the crystalgraphic theory, we suggest ferroelectric single crystal has possible 8 variants; transformation of variants can be characterized by the habit planes illustrated in Fig. 11.where n denotes the unit normal to the habit plane and s refers to the direction of transformation.
In general, a criterion (i.e., the phase transformation criterion) exists for the phase transformation of ferroelectric single crystal, and the material is assumed to undergo phase transformation when at least one of the 8 variants satisfies the phase transformation criterion.This is detailed as follows.Upon loading, the condition to produce R→O polarization rotation on a specified habit plane is that the driving force G of that plane reaches the critical value G 0O .The driving force is composed of the chemical driving force G chem and mechanical driving force G mech [36] A similar condition holds for reverse transformation from O→R polarization rotation with a critical value G 0A In Eqs. ( 2) and ( 3), the mechanical driving force can be expressed by Where τ a and E a are the resolved stress on the "a" transformation system, and γ * and P * denote the associated transformation strain and transformation polarization.Following the crystal theory of plasticity, the resolved stress τ a of the variant "a" is related to the stress tensor σ ij and the Schmid factor a ij a by 3 , 1 where the Schmid factor is defined as follow with m i and n i being the unit normal to the habit plane and shear direction to the variant "a", respectively.The resolved stress E a of the variant "a" is related to the electrical field tensor E i The chemical driving force in Eqs. ( 2) and ( 3) is assumed to be a linear function of the temperature where T and T 0 denote the temperature in the single crystal and the equilibrium tempera- ture respectively, and β is the stress-temperature coefficient.Note that the equilibrium tem- perature T 0 is defined as the average of the starting temperature of O phase transformation and that of R phase transformation, that is When applying the rate independent crystal theory based model, Eqs. ( 2) and ( 3), to simulate the behavior of ferroelectric single crystals, one of the most computationally consuming tasks is to determine the set of instantaneously active transformation systems among the 8 possible variants at crystal level.This determination is usually achieved by an iterative procedure and must be carried out at each loading step, requiring extensive computation.Note that in the crystal theory of plasticity, a similar problem exists whilst in a rate dependent viscoplastic version of crystal theory of plasticity, determination of the set of active transformation systems is not necessary.As a result, computation effort can be reduced significantly.Following this idea, a viscoplastic version of Eqs. ( 21) and ( 32) are proposed and employed in this study.In the viscoplastic crystal model for PMN-0.32PTferroelectric single crystals all transformation systems are assumed to be instantaneously active of varying extent, which is governed by a rate dependent viscoplastic law.In this paper, the phase transformation of variant "a" is assumed to comply with the following power law of viscoplasticity, where f ˙a is phase fraction transformation rate of variant "a", f ˙0 is reference phase fraction transformation rate, G a is the driving force of variant "a", G 0 is the critical driving force, ex- ponents k and m are material parameters dictating the rate effect, c depends upon the phase fraction, c 0 the reference value at initial state.
The phase fraction of variant "a" is calculated as, Summation over all possible variants provides the phase fraction for the whole ferroelectric single crystal, that is The incremental form of Eq. ( 12) can be written as The transformation strain tensor ε ij tr , which is associated with d f a , can be obtained as where α ij is Schmid factor.d ε ij tr is increment of phase train.γ * is the maximal strain transfor- mation during loading.d f a is increment of phase fraction.
Elastic strain and electrical field induced strain during loading can expressed as Increment of the total strain d ε ij can then expressed as the sum of an elastic component dε ij e , a electrical field induced component d ε ij E and a transformation induced component

Numerical results
To validate the model, the corresponding calculation result is compared to experimental stress-strain curve.The material parameters used in this study are as follow [37,38].Materi-  Table 3.The coordinate of particle in O phase BaTiO 3 ferroelectric single crystal cell from refer.

Conclusion and discussion
In this paper, Stress induced strain and polarization, and electric field induced "butterfly" curves and polarization loops for a set of compressive bias stress for {001}, {011} and {111} poled PMN-0.32PTsingle crystals are experimentally explored.Obtained results indicate that high piezoelectric responses of PMN-0.32PTsingle crystals are controlled by the anisotropy of the crystals and the multi-domain structure (i.e., engineered domain structure) has a relatively minor effect.Analysis shows that in all three directions the electric field induced aggregate strain Δε and piezoelectric constant d 33 increase with increasing the magnitude of the compressive bias stress.However, when the magnitude of the compressive bias stress is further increased the electric field induce Δε and d 33 decrease.As a result, an optimized compressive bias stress exists for the purpose of enhancing the electromechanical properties of {001} and {111} oriented PMN-0.32PTsingle crystals.These results have apparent importance in the design of actuators and sensors using PMN-0.32PTsingle crystals.It is found that the observed stress induced strain and polarization in {001}-oriented PMN-0.32PTcan be described by a polarization rotation mechanism, i.e., polarization rotates from rhombohedral (R) to orthohombic (O) and tetragonal (T) phases through the intermediate Monoclinic (M) phase during loading, and O to R transition during unloading.However, domain switching is believed to be the main mechanism dictating the electromechanical behavior of {011} and {111} oriented PMN-0.32PTsingle crystals.polarization rotation model is developed to explain the observed behaviors of PMN-0.32PT.The stress-strain curve along <001> crystallographic direction of ferroelectric single crystal BaTiO 3 is calculated with the first principle method.Obtained results show that the R→M→O polarization rotation (phase transformation) takes place in rhombohedral BaTiO 3 ferroelectric single crystals under compression, which is consistent with the polarization rotation model.Based on the polarization rotation model, a constitutive model of PMN-0.32PT is proposed based on micromechanical model.It is shown that the developed model can faithfully capture the key characteristic of the observed constitutive behavior of <001> oriented PMN-0.32PT.

Figure 1 .
Figure 1.(a) Multidomain rhombohedral crystal obtained by poling along the {001} orientation.The solid line arrows show possible directions of the polarization vector in a fully {001} poled rhombohedral crystal, the dashed line arrows show possible directions of the polarization vector in orthorhombic and tetragonal phases crystal.(b) Two domain rhombohedral crystal obtained by {011} poling.The solid line arrows indicate type 1 domains in a fully {011} poled crystal and dashed line arrows indicate possible type 2 domains under compressive stress.(c) Monodomain rhombohedral crystal poled along the {111} orientation.The solid and dashed line arrows have similar means to (b).The hollow arrows in Fig. a ~ Fig. c show the electric field and stress loading.The small rhombohedral distortion is neglected and all numbers and notations refer to a quasi-cubic unit cell.

Figure 2 .
Figure 2. The sketch of experimental set

Figure 3 .
Figure 3. Electric field induced polarization and strain responses for {001}, {011} and {111}-oriented crystals of PMN-0.32PT:(a) P 3 − E 3 curves and (b) ε 33 − E 3 curves.(c) Relationship between rhombohedral (solid lines) and cubic (dashed lines) coordinate systems.Subscripts r and c denote directions with respect to the rhombohedral and cubic coordinate systems, respectively.(d) Orientational dependence of piezoelectric coefficient d * 33 of PMN-0.32PT in the

Figure 5 .
Figure 5. Schematic drawing of the polarization rotation from R to T and O phases due to <001> direction compression: (a) the initial state, (b) polarization starts to rotate from R to O and T via two pathways M A and M C , (c) R phase completely transformed, (d) O phase starts to switch back to R phase upon unloading, (e) O phase completely switched back to R phase.The states (a) to (e) correspond to the stress levels marked in Figs.2(a) and 2(b).At the initial poled state (a), PMN-0.32PTpossesses four equivalent <111>-oriented polarizations vector close to <001> direction though only one polarization rotation processes is shown here.

Fig. 5
Fig.5illustrates the polarization rotation mechanism of PMN-0.32PTsingle crystals under stress loading and unloading cycle.The polarization vector states shown in Figs.5(a)-5(e) correspond to the stress levels marked on the curves in Figs.4 for the -40MPa stress cycle.In the initial poled state, PMN-0.32PTsingle crystals possess four equivalent <111> polarizations, with only one shown in Fig.5afor the sake of clarity.Upon loading, polarization vector starts to rotate from R to O and T phases through intermediate phases M A and M C , when the compressive stress exceeds about 15MPa in magnitude (Fig.5(b)).This gives a mixture of R and T phases and remarkable augment in the polarization and strain change around σ 33 =-20MPa (Figs. 4).When the presence of an electric field along <110> direction, polarization rotation only occurs from R to O under compression along <001> direction.In the absence of electric field along <110>, however, almost all the polarization vectors may switch to T and O phases when the compression exceed 30MPa in magnitude (Fig.5(c)).When the magnitude of the compressive loading is further increased, there is no more polarization rotation and the PMN-0.32PTsingle crystal shows linear response in Figs.4.It is noted that O phase is usually unstable (its free energy balance between R and T phases depends on the electric Figure.6 Electric field induced 33

Figure. 7
Figure.7 Effect of compressive bias stress on (a) the coercive field c E , (b) remnant polarization r P , (c

Figure. 8 Figure 8 .
Figure.8The cell calculation model of BaTiO 3 in R phase

Figure 9 .
Figure9.The stress-train curve of BaTiO 3 ferrorlctric single crystal obtained by the method amentioned in this paper, the {001} oriented BaTiO 3 ferroelectric single crystal is in T phase.

Figure 10 .
Figure 10.(a)The stress-strain curve of {001} oriented BaTiO 3 ferroelectric single crystal calculated in first principle calculation, the scatter point is calculation value and the broken line fit the scatter point.(b)The sketch map of polarization rotation during the calculation.The "a,b,c" letter in (a) is correspond to the "a,b,c" letter in (b), the letter indicates different polarization states during loading.

Figure 11 . 3 )Figure 12 .Table 1 .Table 2 .
Figure 11.The comparison between ferroelectric phase transformation and plastic slip of single crystal, n is direction normal to slip surface, s is direction along phase transformation.
33 =Δε 33 / ΔE 3 where ΔE 3 is 3, respective- ly.By taking the measured value for P r 111 (0.395C/m 2 ), P r 001 and P r 011 are predicted to be 0.228 and 0.322C/m 2 , respectively, which are very close to the measured ones (P r 001 With d 33 =200pC/N taken from our measured value and d 15 =4100, d 31 =-90 and d 22 =1340pC/N,