In this presentation, we present an efficient and robust technique for surface meshing based on harmonic maps. The overall procedure is implemented with a finite element discretization in the open-source mesh generator gmsh. We show how to ensure a one-to-one mapping for the discrete harmonic map and introduce a cubic representation of the geometry based on curved PN triangles. We show that with this approach we are able to recover high quality meshes from both low input STL triangulations and complex surfaces defined by many CAD patches.