We present, to the best of our knowledge, the first complete three-dimensional numerical solutions to a broad range of boundary value problems for a general theory of finite strain gradient elasticity. We have chosen for our work, Toupin's theory [Arch. Rat. Mech. Anal., 11(1), 385-414, 1962]-one of the more general formulations of strain gradient elasticity. Our framework has three crucial ingredients: The first is iso-geometric analysis [Hughes et al., Comp. Meth. App. Mech. Engrg., 194(39-41), 4135-4195, 2005], which we have adopted for its straightforward and robust representation of C1-continuity. The second is a weak treatment of the higher-order Dirichlet boundary conditions in the formulation, which control the development of strain gradients in the solution. The third ingredient is algorithmic (automatic) differentiation, which eliminates the need for linearization "by hand" of the rather complicated geometric and material nonlinearities in gradient elasticity at finite strains. We present a number of numerical solutions to demonstrate that the framework is applicable to arbitrary boundary value problems in three dimensions. We discuss size effects, the role of higher-order boundary conditions, and perhaps most importantly, the relevance of the framework to problems with elastic free energy density functions that are non-convex in strain space. ; Comment: Appeared in Computer Methods in Applied Mechanics and Engineering (CMAME), Volume 278, 15 August 2014, Pages 705-728, ISSN 0045-7825, http://dx.doi.org/10.1016/j.cma.2014.06.015