We present an effective implementation of analytical calculations of the Lyα opacity distribution of the intergalactic medium (IGM) along multiple lines of sight (LOS) to distant quasars in a cosmological setting. The method assumes that the distribution of neutral hydrogen follows that of an underlying dark matter density field, and that the density distribution is a (local) lognormal distribution. It fully accounts for the expected correlations between LOS and the cosmic variance in the large-scale modes of the dark matter distribution. Strong correlations extending up to ~300kpc (proper) and more are found at redshifts z~2-3, in agreement with observations. These correlations are investigated using the cross-correlation coefficient and the cross-power spectrum of the flux distribution along different LOS, and by identifying coincident absorption features as fitted with a Voigt profile fitting routine. The cross-correlation coefficient between the LOS can be used to constrain the shape-parameter Γ of the power spectrum if the temperature and the temperature density relation of the IGM can be determined independently. We also propose a new technique to recover the 3D linear dark matter power spectrum by integrating over 1D flux cross-spectra, which is complementary to the usual `differentiation' of 1D auto-spectra. The cross-power spectrum suffers much less from errors uncorrelated in different LOS, like those introduced by continuum fitting. Investigations of the flux correlations in adjacent LOS should thus allow us to extend studies of the dark matter power spectrum with the Lyα forest to significantly larger scales than is possible with flux auto-power spectra. 30 pairs with separations of 1-2arcmin should be sufficient to determine the 1D cross-spectrum at scales of 60h-1Mpc with an accuracy of about 30 per cent (corresponding to a 15 per cent error of the rms density fluctuation amplitude) if the error is dominated by cosmic variance.
