Abstract. Recent experimental results on giant unilamellar vesicles (GUVs) show that mixed multiple lipid components on the surface of a membrane may decompose into coexisting phases with distinct compositions, with concomitant changes in the surface morphology. The driving forces for the evolution involves bending, line tension along the phase boundaries, inhomogeneous surface energy, and fluid forces. Here we are interested in exploring the emergent morphologies when the flow is present, and in particular when the surface tensions of the coexisting phases are different, which has not been considered previously. In this paper, we present a model capable of describing the nonlinear coupling among flow, membrane morphology, and the evolution of the surface phases. Using an energy variation approach, we derive a generalized surface tension and construct a constitutive equation connecting the surface tension and the phase variables. To investigate the nonlinear dynamics, we develop a numerical method that combines the immersed interface method to solve the flow equations, the level-set method to capture the interface motion, a non-stiff Eulerian algorithm to solve the phase field equations on the evolving surface, and a penalty term that enforces global inextensibility. Our numerical results suggest that the nonhomogeneous surface tension, together with the flow, introduces nontrivial vesicle dynamics including locomotion, wrinkling, and budding.