The article is aimed to study the inverse problem of determining the right-hand side of a subdiffusion equation with Riemann-Liouville fractional derivative whose elliptic part has the most general form and is defined in an arbitrary multi-dimensional domain (with sufficiently smooth boundary). Namely, let \(\Omega \in \mathbb{R}^n\) be an arbitrary bouned domain with sufficiently smooth boundary \(\partial \Omega\), and let \(A(x,D)= \sum \limits_{|\alpha| < m} a_{\alpha}(x) D^{\alpha}\) be a positive formally self-adjoint elliptic differential operator of order \(m=2l\) with sufficiently smooth coefficients \(a_{\alpha}(x)\), \(\alpha = (\alpha_1, ..., \alpha_n)\) is a multiindex and \(D=(D_1, ..., D_n)\), \(D_i=\frac{\partial }{\partial x_i}\), \(0<p\leq 1\). Consider the subdiffusion equation \[ \partial^{p}_t u(x,t) + A(x, D) u(x,t)=f(x), \ \ x \in \Omega, \ \ 0<t \leq T.\] The direct initial-boundary value problem is to find the temperature distribution \(u(x,t)\) under the initial condition \[ \lim \limits_{t \to 0} \partial^{p-1}_t u(x,t) = \phi(x), \ \ x \in \overline{\omega}\] and the boundary conditions \[B_j u(x,t)= \sum \limits_{|\alpha| \leq m_j }b_{\alpha}(x) D^{\alpha}u(x,t)=0,\] \(0 \leq m_j \leq m-1\). The following inverse problem is considered: along with the solution \(u(x,t)\) to the initial-boundary value problem, find the heat source density \(f(x)\) such that the temperature distribution at time \(T\) achives a prescribed level, \(u(x,T)=\Phi(x)\) . The Fourier method is used to prove theorems on the existence and uniqueness of the classical solution of the initial-boundary value problem and on the unique reconstruction of the unknown right-hand side of the equation. The concept of the generalized solution is introduced and a theorem on its existence is proved. The stability of classical and generalized solutions is proved. Requirements for the initial function and for the additional condition are established under which the classical Fourier method can be applied to the inverse problem under consideration. The results obtained are also of interest for the classical diffusion equation.