A consistent treatment of the coupling of surface energy and elasticity within the multi-phasefield framework is presented. The model accurately reproduces stress distribution in a number of analytically tractable, yet non-trivial, cases including different types of spherical heterogeneities and a thin plate suspending in a gas environment. It is then used to study the stress distribution inside elastic bodies with non-spherical geometries, such as a solid ellipsoid and a sintered structure. In these latter cases, it is shown that the interplay between deformation and spatially variable surface curvature leads to heterogeneous stress distribution across the specimen.