There are a few different ways of answering this one. For brevity, I'm going to be a bit hand-wavey. There is actually still some research going on with this.
Certain spacetimes will always have a conserved energy. These are the spacetimes that have what is called a global timelike (or, if you're wanting to be super careful and pedantic, perhaps null) Killing vector. Math-types will define this as a vector whose lowered form satisfies the Killing equation: $\nabla_{a}\xi_{b} + \nabla_{b} \xi_{a} = 0$. Physicists will just say that $\xi^{a}$ is a vector that generates time (or null) translations of the spacetime, and that Killing's equation just tells us that these translations are symmetries of the spacetime's geometry. If this is true, it is pretty easy to show that all geodesics will have a conserved quantity associated with the time component of their translation, which we can interpret as the gravitational potential energy of the observer (though there are some new relativistic effects--for instance, in the case of objects orbiting a star, you see a coupling between the mass of the star and the orbiting objects angular momentum that does not show up classically). The fact that you can define a conserved energy here is strongly associated with the fact that you can assign a conserved energy in any Hamiltonian system in which the time does not explicitly appear in the Hamiltonian--> time translation being a symmetry of the Hamiltonian means that there is a conserved energy associated with that symmetry. If time translation is a symmetry of the spacetime, you get a conserved energy in exactly the same way.
Secondly, you can have a surface in the spacetime (but not necessarily the whole spacetime) that has a conserved killing tangent vector. Then, the argument from above still follows, but that energy is a charge living on that surface. Since integrals over a surface can be converted to integrals over a bulk by Gauss's theorem, we can, in analogy with Gauss's Law, interpret these energies as the energy of the mass and energy inside the surface. If the surface is conformal spacelike infinity of an asymptotically flat spacetime, this is the ADM Energy. If it is conformal null infinity of an asymptotically flat spacetime, it is the Bondi energy. You can associate similar charges with Isolated Horizons, as well, as they have null Killing vectors associated with them, and this is the basis of the quasi-local energies worked out by York and Brown amongst others.
What you can't have is a tensor quantity that is globally defined that one can easily associate with 'energy density' of the gravitational field, or define one of these energies for a general spacetime. The reason for this is that one needs a time with which to associate a conserved quantity conjugate to time. But if there is no unique way of specifying time, and especially no way to specify time in such a way that it generates some sort of symmetry, then there is no way to move forward with this procedure. For this reason, a great many general spacetimes have quite pathological features. Only a very small proprotion of known exact solutions to Einstein's Equation are believed to have much to do with physics.
The Lagrangian formalism of physics is the way to start here. In this formulation, we define a function that maps all of the possible paths a particle takes to the reals, and call this the Lagrangian. Then, the [classical] path traveled by a particle is the path for which the Lagrangian has zero derivative with respect to small changes in each of the paths.
It turns out, due to a result known as Noether's theorem, that if the Lagrangian remains unchanged due to a symmetry, then the motion of the particles will necessarily have a conserved quantity.
Energy is a conserved quantity associated with a time translation symmetry in the Lagrangian of a system. So, if your Lagrangian is unchanged after substituting $t^{\prime} = t + c$ for $t$, then Noether's theorem tells us that the Lagrangian will have a conserved quantity. This quantity is the energy. If you know something about Lagrangians, you can explicitly calculate it. There are numerous googlable resources on all of these words, with links to how these calculations happen. I will answer further questions in edits.
Best Answer
We represent energy density and momentum and flux (which source the gravitational field) using a tensorial quantity $T^{\mu\nu}$. Conservation of energy is a ill defined notion because we cannot yet calculate all the the things which contribute to $T^{\mu\nu}$ - for example the gravitons, which we believe has to do with LHS (the gravitational/curvature part), of Einstein equation can interact with itself (higher order perturbations) which will contribute to $T^{\mu\nu}$.
Hence, we cannot directly define covariant conservation of Energy momentum tensor $\nabla_\mu T^{\mu\nu} = 0$. To our rescue, we have Bianchi identities which are mathematical way of implying conservation of curvature(LHS) tensor. From this identity we can deduce covariant conservation of energy momentum tensor. However, we can always define local conservation $\partial_\mu T^{\mu\nu} = 0$ independent of such identities.
Vectors and differential forms are defined locally hence $T^{\mu\nu}$ is a local quantity. In general, we do not have a globally conserved energy-momentum. On flat spacetime, we could add the local current vectors or energy densities to define global conservation of quantities. On the other hand, tensors (on curved spaces) which are defined in different tangent spaces cannot be added. Hence we define local covariant conservation.