Coupled diffusion with elasticity models are utilized in various applications, e.g., the phase decomposition and coarsening of solder alloys and microphase separation in hydrogels. Extensions to the classical Cahn-Hilliard equation incorporating elasticity have already been discussed and formalized in a more general sense, often referred to as the Cahn-Larche system. The arising set of coupled equations has been studied numerically for the cases of small and large elastic deformations. However, the coupling with gradient elasticity lacks numerical investigations in literature, despite its relevance for problems where classical elasticity is known to give qualitatively different predictions, e.g., under bending loading conditions. In the present contribution, we consider a coupled mechano-diffusion problem that combines Mindlin’s strain gradient elasticity with the Cahn–Hilliard diffusion model. This augments the modelling capabilities by incorporating elasticity-related length scales reflecting microstructural effects at the continuum level. The higher-order boundary value problem is formulated in a variational form within a proper Sobolev space setting. Conforming Galerkin discretizations for numerical results are obtained via an isogeometric analysis approach. We first assess the numerical formulation through a one-dimensional benchmark problem, verifying the implementation and establishing convergence properties. We then analyse size effects in 2D phase separation and coarsening, with particular emphasis on bending-type boundary conditions under which classical and strain gradient elasticity predict qualitatively different mechanical and microstructural behaviour. The simulations reveal how gradient elasticity modifies the evolution of the microstructure during coarsening and clarify the relation between the various length scales involved.