higher terms truncated), the free energy is given by: where Px, Py, and Pz are the components of the polarization vector in the x, y, and z directions respectively, and the coefficients, must be consistent with the crystal symmetry. To investigate domain formation and other phenomena in ferroelectrics, these equations are often used in the context of a phase field model. Typically, this involves adding a gradient term, an electrostatic term and an elastic term to the free energy. The equations are then discretized onto a grid using the finite difference method or finite element method and solved subject to the constraints of Gauss's law and Linear elasticity. In all known ferroelectrics, and . These coefficients may be obtained experimentally or from ab-initio simulations. For ferroelectrics with a first order phase transition, , whereas for a second order phase transition. The spontaneous polarization, Ps of a ferroelectric for a cubic to tetragonal phase transition may be obtained by considering the 1D expression of the free energy which is: This free energy has the shape of a double well potential with two free energy minima at , the spontaneous polarization. We find the derivative of the free energy, and set it equal to zero in order to solve for : Since the Ps = 0 solution of this equation rather corresponds to a free energy maxima in the ferroelectric phase, the desired solutions for Ps correspond to setting the remaining factor to zero: whose solution is: and eliminating solutions which take the square root of a negative number (for either the first or second order phase transitions) gives: If , the solution for the spontaneous polarization reduces to: The hysteresis loop (Px versus Ex) may be obtained from the free energy expansion by including the term -Ex Px corresponding to the energy due to an external electric field Ex interacting with the polarization Px, as follows: We find the stable polarization values of Px under the influence of the external field, now denoted as Pe, again by setting the derivative of the energy with respect to Px to zero: Plotting Ex (on the X axis) as a function of Pe (but on the Y axis) gives an 'S' shaped curve which is multi-valued in Pe for some values of Ex. The central part of the 'S' corresponds to a free energy local maximum (since ). Elimination of this region, and connection