  We are interested in numerically approximating the solution ${\bf U}(t)$ of
the large dimensional semilinear matrix differential equation $\dot{\bf U}(t) =
{ \bf A}{\bf U}(t) + {\bf U}(t){ \bf B} + {\cal F}({\bf U},t)$, with
appropriate starting and boundary conditions, and $ t \in [0, T_f]$. In the
framework of the Proper Orthogonal Decomposition (POD) methodology and the
Discrete Empirical Interpolation Method (DEIM), we derive a novel
matrix-oriented reduction process leading to an effective, structure aware low
order approximation of the original problem. The reduction of the nonlinear
term is also performed by means of a fully matricial interpolation using left
and right projections onto two distinct reduction spaces, giving rise to a new
two-sided version of DEIM. By maintaining a matrix-oriented reduction, we are
able to employ first order exponential integrators at negligible costs.
  Numerical experiments on benchmark problems illustrate the effectiveness of
the new setting.
