  The multi-configurational self-consistent field theory is considered the
standard starting point for almost all multireference approaches required for
strongly-correlated molecular problems. The limitation of the approach is
generally given by the number of strongly-correlated orbitals in the molecule,
as its cost will grow exponentially with this number. We present a new
multi-configurational self-consistent field approach, wherein linear
determinant coefficients of a multi-configurational wavefunction are optimized
via the stochastic full configuration interaction quantum Monte Carlo technique
at greatly reduced computational cost, with non-linear orbital rotation
parameters updated variationally based on this sampled wavefunction. This
extends this approach to strongly-correlated systems with far larger active
spaces than it is possible to treat by conventional means. By comparison with
this traditional approach, we demonstrate that the introduction of stochastic
noise in both the determinant amplitudes and the gradient and Hessian of the
orbital rotations does not preclude robust and reliable convergence of the
orbital optimization. It can even improve the ability to avoid convergence to
local minima in the orbital space, and therefore aid in finding variationally
lower-energy solutions. We consider the effect on the convergence of the
orbitals as the number of walkers and the sampling time within the active space
increases, as well as the effect on the final energy and error. The scope of
the new protocol is demonstrated with a study of the increasingly strongly
correlated electronic structure in a series of polycyclic aromatic
hydrocarbons, up to the large coronene molecule in a complete active space of
24 pi electrons in 24 orbitals, requiring only modest computational resources.
