  We describe the numerical scheme for the discretization and solution of 2D
elliptic equations with strongly varying piecewise constant coefficients
arising in the stochastic homogenization of multiscale composite materials. An
efficient stiffness matrix generation scheme based on assembling the local
Kronecker product matrices is introduced. The resulting large linear systems of
equations are solved by the preconditioned CG iteration with a convergence rate
that is independent of the grid size and the variation in jumping coefficients
(contrast). Using this solver we numerically investigate the convergence of the
Representative Volume Element (RVE) method in stochastic homogenization that
extracts the effective behavior of the random coefficient field. Our numerical
experiments confirm the asymptotic convergence rate of systematic error and
standard deviation in the size of RVE rigorously established in [6]. The
asymptotic behavior of covariances of the homogenized matrix in the form of a
quartic tensor is also studied numerically. Our approach allows laptop
computation of sufficiently large number of stochastic realizations even for
large sizes of the RVE.
