We demonstrate a fast numerical method of theoretical studies of skyrmion lattice or spiral order in magnetic materials with Dzyaloshinsky-Moriya interaction. The method is based on the Fourier expansion of the magnetisation combined with a minimisation of the free energy functional of the magnetic material in Fourier space, yielding the optimal configuration of the system for any given set of parameters. We employ a Lagrange multiplier technique in order to satisfy micromagnetic constraints. We apply this method to a system that exhibits, depending on the parameter choice, ferromagnetic, skyrmion lattice, or spiral (helical) order. Known critical fields corresponding to the helical-skyrmion as well as the skyrmion-ferromagnet phase transitions are reproduced with high precision. Using this numerical mehtod we predict new types of excited (metastable) states of the skyrmion lattice, which may be stabilized by coupling the skyrmion lattice with a superconducting vortex lattice. The method can be readily adapted to other micromagnetic systems.