For an arbitrary strong, spherically symmetric super-horizon curvature perturbation, we present analytical solutions of the Einstein equations in terms of asymptotic expansion over the ratio of the Hubble radius to the length-scale of the curvature perturbation under consideration. To obtain this solution we develop a recursive method of quasi-linearization which reduces the problem to a system of coupled ordinary differential equations for the $N$-th order terms in the asymptotic expansion with sources consisting of a non-linear combination of the lower order terms. We use this solution for setting initial conditions for subsequent numerical computations. For an arbitrary precision requirement predetermined by the intended accuracy and stability of the computer code, our analytical solution yields optimal truncated asymptotic expansion which can be used to find the upper limit on the moment of time when the initial conditions expressed in terms of such truncated expansion should be set. Examples of how these truncated (up to eighth order) solutions provide initial conditions with given accuracy for different radial profiles of curvature perturbations are presented.