Revisiting Conserved Charges in Higher Curvature Gravitational Theories

Restricting the covariant gravitational phase spaces to the manifold of parametrized families of solutions, the mass, angular momenta, entropies, and electric charges can be calculated by a single and simple method. In this method, which has been called"solution phase space method,"conserved charges are unambiguous and regular. Moreover, assuming the generators of the charges to be exact symmetries, entropies and other conserved charges can be calculated on almost arbitrary surfaces, not necessarily horizons or asymptotics. Hence, the first law of thermodynamics would be a local identity relating the exact symmetries to which the mass, angular momentum, electric charge, and entropy are attributed. In this paper, we apply this powerful method to the $f(R)$ gravitational theories accompanied by the terms quadratic in the Riemann and Ricci tensors. Furthermore, conserved charges and the first law of thermodynamics for some of their black hole solutions are exemplified. The examples include warped AdS$_3$, charged static BTZ, and 3-dimensional $z=3$ Lifshitz black holes.


Introduction
Since the realization of general relativity and Noether's theorems in 1915, there have been numerous attempts to attribute local (and later quasi-local) conserved charges to the symmetries in the presence of gravity. Nowadays, a century after that, the literature on this subject is rich and well established, but still in progress: local conserved charges have not been consistently formulated, while there are different successful formulations for quasi-local conserved charges (see Refs. [1,2] as reviews). Among the different approaches, two main lines of formulation can be distinguished: one is the Hamiltonian formulation which is based on space+time decomposition, and the other one is the Lagrangian formulation which is based on spacetime covariance. Reviewing the timeline of the major progress in these two formulations (which for sure might miss some interesting contributions) can give us an overview, in addition to clarify the motivations of the analysis in this paper.
Precursor of the Hamiltonian formulation was the introduction of quasi-local charges by Komar in 1959 [3]. In the Komar's method, quasi-local mass and angular momentum for asymptotic flat solutions could be found by an integration over a codimension-2 surface at constant time asymptotics. Soon after, the Hamiltonian formulation of the gravitational theories was elaborated in a series of works in 1959-62 by Arnowitt-Deser-Misner [4][5][6], known as the ADM formulation. Hence, in addition to introducing a sophisticated formulation for gravitational dynamics, the calculation of the mass and angular momentum at the constant time asymptotics was put on a firm basis (reviewed e.g. in Ref. [7]). After that, a similar formulation for the null asymptotics was proposed by Bondi et al. in 1962 [8, 9]. The Hamiltonian formulation for the asymptotic flat spacetimes reached its mature presentation by Regge-Teitelboim in 1974 [10], emphasizing the role of the surface terms in the Hamiltonian. Nonetheless, there was a shortcoming of the formulation, when the asymptotic flatness was relaxed to include asymptotic (anti) de Sitter solutions, mainly because of the appearance of divergent conserved charges. Later progress in this line of formulation has been mainly in the direction of ameliorating this problem (see the review [11]). Transferring to the Hamilton-Jacobi formulation by Brown-York in 1992 [12], ensued by addition of a surface counterterm to the Lagrangian [13], has been one of them. Reformulation of conserved charges based on covariantly defined conserved currents, and subtracting the contributions from a reference solution, has been another method proposed by Abbott-Deser in 1982 [14], and this was completed and extended to higher curvature theories by Deser-Tekin [15,16]. This method is known as the ADT method in the literature. Last but not least a contribution has been presented by Kim et al. in 2013 [17], and it is known as the quasi-local method. It is based on the ADT off-shell conserved current, while two major changes are considered: (1) instead of considering the difference between the solution and a reference solution, a one-parameter integration from the reference solution to the solution under consideration is performed (advocated in Refs. [18,19]), (2) the surface of integration is relaxed to be in the interior of the geometry. Hence, it provides a powerful method for calculating conserved charges, specifically for black hole solutions.
The second line of formulating conserved charges is based on the Lagrangian, which is covariant from the beginning. It was initiated by Ashtekar et al. [20,21] and Crnkovic-Witten [22] in 1987, and was consistently formulated in a series of works by Wald et al. [18,[23][24][25]. In this formulation, which is called covariant phase space formulation (see Refs. [26][27][28] for reviews), a covariant phase space was built without a space+time decomposition. The phase space manifold was constructed from dynamical fields all over the spacetime, without recruiting their momentum conjugates. The symplectic form was read from the Lagrangian, which entailed a concrete formulation for conserved charges associated with diffeomorphisms and gauge transformations. Besides, in this progress, the entropy of the nonextermal [24,25] (and later, of the extremal [29]) black holes was introduced as a conserved charge calculated on black hole horizons. Later in 2002, Barnich-Brandt reformulated the formalism in the language of variational bicomplex, in addition to proposing a version directly from the equation of motion (e.o.m) instead of the Lagrangian [30]. The conceptually strong point (but pragmatically weak point) of the covariant phase space formulation is that the phase space manifold and its tangent space, which are crucial for explicit calculation of conserved charges, are determined by some fall-off conditions. The usual fall-off conditions, although restrict the manifold, do not usually determine it such that calculation of the charges could be performed explicitly. To put the formulation into its full power of calculability, in Ref. [31] the manifold and its tangent space are constructed explicitly and directly from the beginning. The elaboration of the phase space and its tangent space, accompanied by relaxing the calculation of the entropy over horizons, has made the formulation to be a universal tool in the context of conserved charge calculations. This method can be dubbed a solution phase space method, because the phase space is constructed by some family of parametrized solutions. Interestingly, the recent independent progresses, the one by Kim et al. [17] in the Hamiltonian formulation, and the solution phase space method [31] in the Lagrangian formulation, have brought about the two lines of formulations to converge.
In this paper, we apply the solution phase space method to the higher curvature theories. For clarity, we will focus on f (R) gravity accompanied by quadratic terms in the Riemann and Ricci tensors (see Lagrangian (3.1)), although the generalization is straightforward. One of the motivations for this work is examining the method for gravitational theories beyond the Einstein-Hilbert gravity. Another motivation is providing detailed materials needed to perform calculations for general enough higher curvature theories. The analysis can be considered in parallel with higher curvature analysis in other methods, specifically the ADT and quasi-local methods studied e.g. in Refs. [15,16,[32][33][34][35][36][37][38][39][40][41][42][43]. In the following sections, first a review on the solution phase space method is presented. Then it is applied to the higher curvature theories. Finally, some interesting examples are provided and compared with the results of other methods.

Solution phase space method
Solution phase space method (SPSM) is a method for calculating conserved charges in gravitational theories. The SPSM is based on a powerful but not yet fully appreciated covariant formulation of gravitational phase spaces, which we are going to review in the next subsection. It is applicable to the solutions which are parametrized by some parameters p j . Specifically, it is a convenient method for calculating mass, angular momenta, entropies, and electric charges associated with the black hole solutions, although it is not exclusive to them [26,31].
Before delving into the details, it can be helpful to have a look at the big picture and the bottom line: calculation of variation of a conserved charge needs three pieces of information as inputs: (1) the theory in d dimensional spacetime, (2) the solution and some perturbation around it for which charge is calculated, and (3) the symmetry to which the charge is attributed. At the end of the day, integrating a d − 2-form k η (δΦ,Φ) over any codimension-2 surface yields the variation of the charge. Concerning the three inputs mentioned above, k is unambiguously determined by the theory. Its argumentsΦ andδΦ denote some elaborated solutions and perturbations. η carries information as regards the symmetry. If the result would be integrable, then an integration overδΦ produces the finite charge.

Covariant phase space formulation
Covariant phase space formulation is an appropriate and well-established construction of the gravitational phase spaces [20][21][22][23][24][25]. To have a self-contained document, we will review the basics of the formulation here, which might have some overlaps with the reviews in Refs. [31,44]. At the outset, it would be useful recalling some relevant elementary properties of a phase space. Phase space F(M, Ω) is a manifold M equipped with a closed nondegenerate symplectic form Ω. In order to introduce a physical phase space, one usually begins with a given Lagrangian. For example, dynamics of a particle in one dimension can be described by the Lagrangian L = mq 2 2 − V (q). M would be built by the position and its momentum conjugate (q, p), equipped with the canonical symplectic 2-form Ω = δp ∧ δq. A simple way to derive this symplectic form is to vary the Lagrangian Then, by recognizing the first term as the e.o.m, and the second term as a total derivative, d dt in the latter has to be dropped, and its exterior derivative on the phase space should be taken: Equipped with the Ω, and for a given vector field v on the M, a charge variation δH v can be defined by For instance, in our simple example, choosing v = ∂ q , then δH v = ∂ q · Ω = δp. Hence, H v = p, which is the mathematical manifestation of "momentum is the charge attributed to the translation in space". Notice that in order for a charge to be conserved, one needs extra conditions. In the case of the simple example mentioned above, p would be conserved if only V (q) would be a constant.
Covariant phase space: Let us begin with a given Lagrangian in d dimensional spacetime, with some classical dynamical fields, collectively denoted by Φ(x µ ). The fields might include the metric g αβ , some Abelian gauge fields A a µ , some scalars φ I , etc. It is usual to build the phase space canonically, i.e. to build the M from a subset of field configurations Φ( x) and their momentum conjugates defined on some privileged time foliation of spacetime. In this construction, solutions to the equation of motion are some curves on M parametrized by the time. Interestingly, in the context of generally covariant gravitational theories, there is a more suitable construction which does not break general covariance by specifying a time foliation. In this construction, M is composed of dynamical field configurations all over the spacetime Φ(x µ ). On the other hand, the field conjugates would not be needed to construct the manifold. As a result, any solution to the equation of motion in the phase space would be a point on M, instead of a curve. The tangent space of the manifold is also constituted from a subset of perturbations δΦ(x µ ).
Symplectic structure: The manifold M which is constructed is a phase space. The symplectic 2-form is constructed from the Lagrangian d-form L. To this end, using the same method as in Eqs. (2.1) and (2.2), first the Lee-Wald (d−1)-form Θ is picked up from the surface term appearing in the variation of Lagrangian: In the equation above, E Φ denotes equations of motion for the fields Φ, on which the summation convention should be understood. δ and d are exterior derivatives on M and on spacetime respectively. Then the pre-symplectic form can be defined as [23][24][25]] where The Σ is some codimension-1 (Cauchy) surface and δ 1,2 Φ are some members of the tangent space. The ω LW , which is a 2-form over the phase space and a d−1-form over the spacetime, is called pre-symplectic current. By construction, Ω in Eq. (2.5) is antisymmetric in δ 1 Φ ↔ δ 2 Φ, and is a closed form δΩ = 0. If it would be a non-degenerate form, it could be used to construct a symplectic structure. In this case, one drops the prefix in "pre-symplectic," and calls it a symplectic form.
Conservation: Apparently, Ω LW in Eq. (2.5) depends on a non-covariantly chosen surface Σ. In order to make Ω LW independent of Σ, which in this context is called "conservation of symplectic form," one needs dω LW = 0. Moreover, the flow of ω LW passing throughout the boundaries ∂Σ should vanish. The former is achieved if Φ and δΦ satisfy the e.o.m and the linearized e.o.m respectively. So, it is standard to restrict the phase space to the solutions, as we will do in the rest of the paper. On the other hand, achievement of the latter needs extra conditions, usually some boundary conditions on the perturbations.
Ambiguities: There are two kinds of ambiguities present in the covariant phase space formulation; one an irrelevant and another one a relevant. The irrelevant one originates from the fact that the formulation is based on Lagrangian formulation, which is itself ambiguous up to a surface term L → L + dK. Nonetheless, although it results to Θ → Θ + δK, but ω remains intact because of δ 2 K = 0. Another ambiguity, which is the relevant one, is an ambiguity originating from the definition of Θ in Eq. (2.4); one can add an exact This ambiguity entails corresponding ambiguities in Ω defined above, through Conserved charges: Let us consider a vector field ξ = ξ µ ∂ µ defined over the spacetime, which generates the diffeomorphism x µ → x µ −ξ µ . In addition, we might have some scalars on the spacetime λ a , generating the gauge transformations A a µ → A a µ + ∂ µ λ a . We can denote the generator of the combination diffeomorphism+gauge transformations by ǫ = {ξ, λ a } such that δ ǫ Φ ≡ L ξ Φ+δ λ a A a . Being equipped with the symplectic form, and motivated by the definition of charge Eq. (2.3), we might be able to associate a conserved charge (interchangeably called Hamiltonian generator) to the ǫ. To this end, variation of the charge is defined as [23][24][25]45] where the last equation follows from Stokes' theorem. The δ [Φ] emphasizes that δ acts on dynamical fields, not ǫ. In the equation above, the integrand in the first integration has been replaced by an exact (d−1)-form dk ǫ . This is the fundamental theorem of the covariant phase space formalism, which can be proved using the on-shell conditions [24,25]. The (d−2)-form k ǫ can be shown to be explicitly (see e.g. Appendix A in Ref. [31] for the detailed derivation) in which Q ǫ is the Noether-Wald charge density [24,25], defined by the relation Hence, by Eq. (2.10), k ǫ can be found for a given theory straightforwardly. Putting it into Eq. (2.9), the charge variation δH ǫ (Φ) can be calculated for an arbitrary generator ǫ. Concerning the conservation, by dω(δΦ, δ ǫ Φ, Φ) = d 2 k ǫ (δΦ, Φ) = 0 there would not be any source or sink in Σ. But vanishing of the flux/leakage through ∂Σ needs to be investigated.
Integrability: δH ǫ (Φ), which is calculated by the last integral in Eq. (2.9), might corresponds to the variation of a finite conserved charge H ǫ . In order to investigate this finite conserved charge, integrability over the phase space is needed. This condition is basically (δ 1 δ 2 − δ 2 δ 1 )H ǫ (Φ) = 0, in which Φs are any field configuration in the presumed phase space F, and δ 1,2 Φ are any arbitrary chosen member of its tangent space. It follows that this condition can be explained as [18,23,45] Symplectic symmetries: As far as conserved charges are concerned, conservation of δH ǫ can be guaranteed if ǫ is chosen such that on-shell. It is because there would not be any flow out of the boundaries locally, and hence globally. The family of ǫ's with this property, which has been dubbed "symplectic symmetry generators" [46], can be divided into two sets [31]: 1. Non-exact symplectic symmetries: The ǫ for which δ ǫ Φ = 0 at least on one point of the phase space. They constitute a closed algebraic structure, and have been proposed to be responsible for generating the phase space of a solution at some given constant thermodynamical variables [45][46][47]. Hence, by studying them and the phase space generated via exponentiating them, one might hope to understand the microstates of the system.

Exact symplectic symmetries:
The ones for which δ ǫ Φ = 0 all over the phase space. For clarity, let us denote such generators by η = {ζ, λ a }. They are considered as generators which by their conserved charges, the set of solutions in different thermodynamical variables can be labeled [31]. The phase space constructed by such field configurations has been called "solution phase space," [31] which will be described in a moment.
It has been conjectured that the phase space associated with the geometries without propagating degrees of freedom are direct product of these two families of phase spaces: the statistical phase space ⊗ the solution phase space [31,48].

Solution phase space method; conserved charges and the first law(s)
In the covariant phase space formulation reviewed above, it is standard to identify the phase space manifold by some asymptotic behaviors, usually through requiring some fall-off conditions. Fall-off conditions, although they delimit the phase space manifold, usually do not determine it completely. "Solution phase space method" is an alternative method for determining the phase space manifold. Restricting the covariant phase space formulation to some explicitly identified manifolds, empowers the calculability of this formulation. The specification is in three aspects, which will be described immediately: 1. Identifying the phase space manifold explicitly, 2. Crystallizing the tangent space of the specified manifolds, 3. Concentrating on the exact sympltectic symmetries of the proposed phase space.
Consider a family of (black hole) solutions to a generally covariant gravitational theory. Usually, such a family is identified by some isometries and some parameters p j . The field configuration of such a family can be denoted collectively byΦ(x µ , p j ). The parameters are some arbitrary (but maybe in some restricted domain of) real numbers appearing in the field configuration of the mentioned solutions. The p j can be reparametrized, but they cannot be removed by coordinate transformations. The manifoldM can be chosen to be composed of the members of the family, up to unphysical coordinate/gauge transformations. As an example, the set of all Schwarzschild black holes parametrized by one free parameter p 1 = m ≥ 0, construct a manifoldM. The symplectic 2-formΩ would be simply the Lee-Wald symplectic form (2.5), which is confined toM. Therefore, theF = (M,Ω) would constitute a phase space, the "solution phase space". Hence, any point of the manifold can be identified byΦ(x µ , p j ). The tangent space ofM is spanned (up to infinitesimal pure gauge transformations) by "parametric variations," which can be found simply by [49] δΦ = ∂Φ ∂p j δp j . In SPSM, the diff+gauge transformations, for which charges are calculated, should be restricted to the symplectic symmetries. However, our main focus would be on the set of exact symmetries. Denoting the generator of the exact symmetries by η = {ζ, λ a } such that δ ηΦ = 0, the ζ would be a Killing vector of all points of the phase spaceM. Besides, its action on the gauge fields has to be canceled by the action of λ a 's, i.e. L ζ A a µ + ∂ µ λ a = 0. As it was advertised in Sect. 2.2, conservation ofδH η is guaranteed. This is because of the relation ω LW (δΦ, δ ηΦ ,Φ) = 0, which itself is a result of linearity of ω LW in δ ηΦ = 0. Hence, there would not be any local and, therefore, any global flow of ω LW out of the boundaries ∂Σ.
Along with guaranteeing the conservation, focusing on the exact symmetries provides us some other nice features: • Independence ofδH η from the choice of ∂Σ: The relation ω LW (δΦ, δ ηΦ ,Φ) = 0 yields an interesting result:δH η would be independent of the chosen ∂Σ. It is because of vanishing of ω LW all over the Σ, and hence, vanishing of ω LW in the region enclosed between two different integrating surfaces ∂Σ 1 and ∂Σ 2 . Then, by the Stokes theorem, and taking the result of Eq. (2.9) into account, the claim is proved. Explaining this result in another way, although the integration in calculatingδH η is over a codimension-2 surface ∂Σ, but the result would be independent of all coordinates, including the two coordinates which are not integrated over.
• Discarding the ambiguity Y: This is because of Summarizing the last two paragraphs, the charges associated with exact symmetries are conserved, unambiguous, and independent of the chosen closed surfaces of integration ∂Σ. So far, the SPSM has provided all materials needed to calculate conserved charge variations. The final tasks would be checking integrability overM, and (if integrable) performing the integration. The integrability can be assessed by replacing Φ, δΦ, and ǫ in the integrability condition Eq. (2.12) byΦ,δΦ, and η, respectively. If integrable, then the integration over arbitrary path onM connecting a reference field configurationΦ(p j ) to the solution under considerationΦ(p 0 j ) yields the final result is the reference point (i.e. constant of integration) for the H η defined on some specific reference field configurationΦ(x µ ;p j ). It is worth mentioning that in order to perform the final tasks mentioned above, there is a shortcut: by the conservation+independence from ∂Σ, theδH η would only be a function of p j and δp j , not any coordinate of the spacetime. Hence, one can simply check whether it is a total derivative or not. Then, in the case of being a total derivative, the integration can be done by an appropriate choice of a reference field. For instance, if theM is parametrized by two parameters {p 1 , p 2 }, and one has foundδH η = p 1 δp 2 + p 2 δp 1 , then it is a total derivativê δH η =δ(p 1 p 2 ). So, the integrated charge would be simply H η = p 1 p 2 + const., where the constant would be fixed by the choice of a reference field configuration.
Before moving on to the next sections, which will provide us explicit examples, it can be useful recalling some remarks in the context of SPSM.
• Not all exact symmetries have integrable conserved charges. Hence, integrability puts constraint on the choice of exact symmetries to which mass, angular momenta etc. are attributed.
• Electric charge associated to the gauge field A a , denoted by Q a , is the conserved charge of the global gauge transformation η Qa = {0, 1 a } in which 1 a means λ a = 1 and λ b = 0 for b = a.
• Similar to the electric charge, mass and angular momenta are conserved charges which are attributed to the geometry as a whole. For the stationary solutions with some number of axial U (1) isometries (labeled by i), one can choose the coordinates such that the corresponding Killing vectors would be ∂ t and ∂ ϕ i , respectively. Then, up to a conventional normalization, the exact symmetries to which the mass M and angular momenta J i are attributed would be η . Ω i ∞ and Φ a ∞ are asymptotic angular velocities and electric potentials, which are usually adopted to be zero.
• In contrast with the charges mentioned above, to each one of the horizons in a geometry, one can associate an entropy. Hence, there might be more than one entropy in a single geometry, e.g. entropy of inner, outer or cosmological horizons. Entropies are considered to be conserved charges for the exact symmetries η H = 2π H are surface gravity, Killing vector, and electric potential of the horizon, respectively. Notice that the ζ H should be accompanied by the rigid gauge transformations λ a = −Φ a H , and be normalized by the surface gravity, in order to have an integrable charge.
• It is worth emphasizing again that assuming the generators of mass, angular momenta, electric charges and entropies to be exact symmetries, these charges can be calculated by integrations over almost arbitrary ∂Σ, and not necessarily the horizons or asymptotics. In this respect, the entropies are on equal footing with other conserved charges.
First law(s) of thermodynamics: To each one of the horizons denoted by "H", an entropy S H , temperature T H = κ H 2π , and some chemical potentials Ω i H , Φ a H etc. can be attributed [50][51][52]. The first law of thermodynamics for the chosen horizon relates δS H to the variations of other conserved charges attributed to the whole geometry. In the SPSM, derivation of the first law(s) is very simple, and originates from a local identity; η H is a linear combination of the generators of mass, angular momenta, and electric charges. From this identity, the first law follows by the linearity of the generic charge variations δH ǫ in terms of the generator ǫ (see Eq. (2.9)). Mathematically [24,25,31], where δS H ≡ δH η H , δM ≡ δH η M , δJ i ≡ δH η J i and δQ a ≡ δH η Qa . Notice that the δ in the proof is a generic perturbation which satisfies linearized e.o.m. So, it is not restricted to the parametric variations. Moreover, integration over horizons or asymptotics does not play any role in this proof.
3 Applying the method to higher curvature theories SPSM has reproduced successfully conserved charges and first law(s) for the standard (black hole) solutions to Einstein-Hilbert gravitational theories. Explicit examples can be found in Refs. [26,31,44]. The goal of this section is utilizing the SPSM for the gravitational theories with higher curvature terms. Explicitly, the Lagrangian which we will focus on, has the metric g αβ , some gauge fields A a µ , and some scalar fields φ I , in arbitrary dimension d: The R µ ναβ , R µν , and R are Riemann tensor, Ricci tensor, and Ricci scalar, respectively. F a = dA a are the field strengths. The coefficients a(φ), b(φ), c ab (φ), and d IJ (φ) are some functions of φ I . Notice that the f (R, φ) covers the Einstein-Hilbert gravity with a cosmological constant. Besides, the Gauss-Bonnet theory of gravity is also covered by the Lagrangian (3.1), hence the simplest Lanczos-Lovelock theories are also included [53][54][55][56]. Generalization to higher Lanczos-Lovelock theories is straightforward, and we will not consider in this paper. The Lagrangian d-form is the Hodge dual of (3.1), L = ⋆L, The ǫ µ 1 µ 2 ···µ d is the Levi-Civita symbol, i.e. ǫ 012···d−1 = +1 and its sign changes with the odd permutations of indices. We will use the conventions h µν ≡ δg µν = g µα g νβ δg αβ , δF µν ≡ g µα g νβ (δdA) αβ = g µα g νβ (dδA) αβ .  (3.6) where the notation f ′ ≡ ∂f ∂R is used. We need to find Θ LW , Q ǫ , and most importantly, the k ǫ for this theory. Their derivation and final results are standard practices in the literature. Hence we only report the final results here. Detailed analysis are similar to the simple Einstein-Hilbert Lagrangian, which can be found e.g. in Appendix A of Ref. [44].
By variation of Lagrangian δL and imposing the e.o.m, the surface d−1-form Θ LW can be read through Eq. (2.4) to be Θ LW = ⋆Θ, i.e.
in which 8) where h ≡ h α α . Having the Θ in our hand, for a generic ǫ = {ξ, λ a }, the Noether-Wald d−2-form Q ǫ can be read through Eq. (2.11) and imposing the e.o.m Eq. (3.4), as in which After varying the Q ǫ with respect to all dynamical fields and utilizing the textbook relations 11) one can find k ǫ by (2.10), to calculate variations of the conserved charges. The result, the final applicable tensor for calculation of conserved charges turns out to be where, using the notations f ′ ≡ ∂f ∂R , Having k ǫ , and equipped with the parametric variationsδΦ, calculation of the conserved charges associated with the exact symmetries η = {ζ, λ a } of the (black hole) solutionŝ Φ(x µ ; p j ) to the Lagrangian (3.1) can be performed.

Some examples
To exemplify, in this section we will work out conserved charges and first law(s) of thermodynamics for some black hole solutions to the Lagrangian (3.1).
Example 1: z = 3 Lifshitz black hole in d= 3 Consider i.e. the new massive gravity (NMG) Lagrangian [58] in which Λ = − 13 2l 2 and m 2 = 1 2l 2 . We can have a family of black holesĝ αβ (x µ ; m) as solution to this theory in 3-dimensional case [59,60] for the z = 3. Let us analyze the thermodynamics of this family of black holes using SPSM. Putting Eq. (4.1) into the general result Eq. (3.12), k µν ǫ can be read for our specific theory. Then, choosing ∂Σ to be surfaces of constant (t, r) for simplicity, the conserved charge variations for an exact symmetry η can be simply read througĥ in which k tr η is the tr component of the k µν η . Inserting parametric variationsδg αβ = ∂ĝ αβ ∂m δm in it, conserved charges can be calculated, irrespective of the asymptotic Lifshitz behavior. Mass: We can choose the stationarity Killing −∂ t as the generator to which the mass is associated. The minus sign has been adopted to make the mass and entropy positive. Hence, by η M = {−∂ t , 0} the result of calculating Eq. (4.4) iŝ .
The reference point is chosen to be S H = 0 for the geometry identified by m = 0. Notice that the entropy is proportional to the area (here the length) of the horizon, but without the usual factor of 1 4 . The results above are in agreement with the results reported in Refs. [35,61]. First law: Having made the entropy free of being calculated on the horizons, the first law of thermodynamics would follow: The reference point for all of the charges above have been chosen to vanish on the geometry identified by m = j = 0. Our results would match exactly with the results reported in Refs. [64,65] if one replaces the parameter m → m − j l . The difference originates from considering the asymptotic angular velocity Ω ∞ in the definition of mass. Hence, the mass calculated here is different from the mass reported in [64,65] by a term, which is Ω ∞ J. First laws: For any generic perturbation which satisfies the linearized e.o.m, the first laws follow: which can be checked for the parametric variations in (4.12), (4.13), and (4.15).

Example 3: Schwarzschild-AdS black holes in d-dimensions
The family of black holes are solutions to the theories where . (4.21) The reference point has been chosen to be M = 0 for the geometry which is identified by m = 0. Angular momentum: By η J = {−∂ ϕ , 0}, angular momentum is calculated to bê Entropy: Surface gravity is a property of solutions, and it is independent of the theory. For the event horizon of the solutions (4.17), it is {ζ H , 0} in which ζ H = ∂ t , the entropy variation attributed to the event horizon, which is calculated on arbitrary surfaces of integration, would beδ Noticing the linearity of δH ǫ (2.9) in ǫ, this result can also be found by multiplication ofδM , which is calculated in Eq. (4.20), by the factor 2π κ H . Hence, using First law: It is simply which can be checked for the parametric variations by the results in Eqs. (4.20) and (4.24). We finish this example by mentioning two remarks: -One can add any factor of Gauss-Bonnet Lagrangian L GB ∝ R 2 −4R µν R µν +R µναβ R µναβ in d = 3, 4 to the Lagrangian (4.18), without affecting the e.o.m and conserved charges.
-The family of black holes in this example has the property R µν = R d g µν . The geometries with such a property are called Einstein geometries. The following theorem (see Ref. [66] and references therein) sheds light on the results of the calculations above. Theorem: Any theory which is described by a Lagrangian L = L(g αβ , R µν ), with an Einstein geometry g αβ as a solution, can be mapped to the Einstein-Hilbert theory with the solutionḡ This theorem clarifies the observation that the mass and entropy calculated above are the mass and entropy in the Einstein-Hilbert theory multiplied by the factor X .
as a solution to the theory described by for Λ = −1 l 2 . k µν ǫ for this theory can be read through Eq. (3.12) by putting f = R − 2Λ, a = b = d IJ = 0, and c ab = δ ab . Making the simplifying choice of taking ∂Σ to be the lines of constant (t, r), conserved charge variations for an exact symmetry η would bê δH η = 2π 0 −ĝ k tr η (δΦ,Φ) dϕ. (4.30) The dynamical fieldsΦ are the metricĝ αβ and gauge fieldÂ in Eq. (4.28), parametrized by p j = {m, q}. So, the parametric variations would bê δg αβ = ∂ĝ αβ ∂m δm + ∂ĝ αβ ∂q δq ,δA µ = ∂Â µ ∂m δm + ∂Â µ ∂q δq . (4.34) Entropies: For any horizon present in this geometry, one can associate an entropy. Surface gravities, electric potentials, and horizon Killing vectors for different horizons would be collectively Reference points for the charges above are chosen to vanish for the pure AdS 3 geometry, i.e. the geometry identified by m = q = 0. In comparison with the calculations done in the literature (see e.g. Ref. [68]), the charges above are finite, and one does not need to regularize any divergent result. Notice that by replacing log( r l ) → log( r r 0 ) for some r 0 , the solution would remain a solution, but M → M + q 2 16 G log( l r 0 ). First laws: For any one of the horizons, the first law would be At the end, it is worth mentioning that in the definition of η H in the examples above, we tacitly assumed κ H = 0. For κ H = 0 cases, which are called extremal black holes, one can find infinite number of exact symmetries in their near horizon regions, as generators of the entropy [29,49]. Using any one of these generators, the SPSM would reproduce the entropy for the extremal black holes too. An explicit example for such an analysis can be found in Ref. [31], where the near horizon of the extremal Kerr-Newman black hole is studied.

Conclusion
In this work, after reviewing the solution phase space method, we applied it to a family of higher curvature gravitational theories. The family which we focused on, contained f (R) gravity, quadratic Riemann and Ricci terms, an arbitrary number of Abelian gauge fields, and arbitrary scalar fields. After elaborating the k ǫ , which is pragmatically the most important differential form for the calculations, four families of black hole solutions were analyzed. Specifically, their conserved charges were calculated, confirming the results formerly calculated by the other methods. By the way, the results ameliorated the divergence appearing in the calculation of mass for the charged static BTZ black hole. The main advantages of the method are: (1) it works for any higher curvature theory in any dimension, (2) asymptotics and horizons are unimportant in the charge calculations, (3) conserved charges are automatically regular, (4) conserved charges are unambiguous, (5) all the charges, including the entropy and electric charge, are calculated by a single machinery, (6) the proof of the first law(s) is very simple.