Barium zirconate titanate [Ba(Zr x Ti 1−x)O 3 ] solid solutions are promising lead-free ferroelectric materials that have received substantial interest. Thermodynamic analysis based on phenomenological Landau-Devonshire theory is a powerful method for theoretical investigation of ferroelectric materials, but cannot be applied to Ba(Zr x Ti 1−x)O 3 because there is no thermodynamic potential. In this paper, a thermodynamic potential for Ba(Zr x Ti 1−x)O 3 (0 ≤ x ≤ 0.3) solid solutions is constructed, and then a thermodynamic analysis carried out. The results accurately reproduce known phase structures and their transition temperatures, with good agreement with experimentally measured polarization, dielectric, and piezoelectric constants. It is found that Ba(Zr x Ti 1−x) O 3 solid solutions at room temperature have three phase boundaries, including a tetragonal-orthorhombic phase boundary at x = 0.013, an orthorhombic-rhombohedral phase boundary at x = 0.0798, and a rhombohedral-paraelectric phase boundary at x = 0.2135. The results also indicate that the chemical composition-induced ferroelectric-paraelectric phase boundary has superior electromechanical properties, suggesting a new way to enhance electromechanical coupling in Ba(Zr x Ti 1−x)O 3 solid solutions.