PDE-constrained optimization problems arise in a broad number of applications. The resulting large-scale saddle-point systems are challenging to solve and acquiring a full solution is often infeasible. We present a new framework to find a low-rank approximation to the solution by reformulating the system into a system of Sylvester-like matrix equations. These matrix equations are subsequently projected onto a small subspace via rational Krlyov-subspace iterations and we obtain a reduced problem by imposing a Galerkin condition on its residual. In our presentation we discuss implementation details and dependence on the problem parameters. Numerical experiments will illustrate the performance of the new strategy.