SUMMARYThe main idea is a dual mathematical formulation and computational implementation of the iterative stochastic perturbation-based finite element method for both linear and nonlinear problems in solid mechanics. A general-order Taylor expansion with random coefficients serves here for the iterative determination of the basic probabilistic characteristics, where linearization procedure widely applicable in stochastic perturbation technique is replaced with the iterative one. The expected values and, in turn, the variances are derived first, and then, they are substituted into the equations for higher central probabilistic moments and additional probabilistic characteristics. The additional formulas for up to the fourth-order probabilistic characteristics are derived thanks to the 10th-order Taylor expansion. Computational implementation of this idea in the stochastic finite element method is provided by using the direct differentiation method and, independently, the response function method with polynomial basis. Numerical experiments include the simple tension of the elastic bar, nonlinear elasto-plastic analysis of the aluminum 2D truss, and solution to the homogenization problem of periodic fiber-reinforced composite with random elastic properties. The expected values, coefficients of variation, skewness, and kurtosis of the structural response determined via this iterative scheme are contrasted with these estimated by the Monte Carlo simulation as well as with the results of the semi-analytical probabilistic technique following the response function method itself. Although the entire methodology is illustrated here by using the Gaussian variables where all odd-order terms simply vanish, it can be extended towards non-Gaussian processes as well and completed with all the perturbation orders.