Domain switching in the vicinity of a crack tip is known as one of the major aspects of local nonlinear behavior of ferroelectrics, and it plays an important role in the fracture behavior. In the present paper, a fracture model based on a phase field continuum and a damage variable is presented to study the fracture behavior of ferroelectrics and its interaction with the domain structures. In this model the energy of fracture is regularized by the damage variable. When the damage variable equals one, it represents undamaged material. In this case, the energy reduces to the phase field potential with the spontaneous polarization being an order parameter, and the system of equations becomes the same as that of a conventional phase field continuum. When the damage variable becomes zero, it represents a crack region, and the potential becomes the energy density stored in the crack medium. The evolution of the damage variable is governed by a Ginzburg-Landau type equation. In this way, the fracture model can simulate the fracture behavior such as crack growth, kinking and formation, with no a priori assumption on fracture criteria and predefined crack paths. The model is implemented in a 2D Finite Element Method in combination with implicit time integration and non-linear Newton iteration. As example, the fracture model is used to simulate the fracture of an edge crack in a ferroelectric single crystal under mechanical mode-I loading. In the simulation crack propagation, kinking and formation are observed. In particular, the results show the interaction between the domain structure evolution and the crack propagation.