Fermionic functional renormalization group for first-order phase transitions: a mean-field model

First-order phase transitions in many-fermion systems are not detected in the susceptibility analysis of common renormalization-group (RG) approaches. Here we introduce a counterterm technique within the functional renormalization-group (fRG) formalism which allows access to all stable and metastable configurations. It becomes possible to study symmetry-broken states which occur through first-order transitions as well as hysteresis phenomena. For continuous transitions, the standard results are reproduced. As an example, we study discrete-symmetry breaking in a mean-field model for a commensurate charge-density wave. An additional benefit of the approach is that away from the critical temperature for the breaking of discrete symmetries large interactions can be avoided at all RG scales.


Introduction
Renormalization group (RG) methods are powerful systematic tools to investigate the emergence of the low-energy physics of many-particle systems. In the Wilsonian RG approach [1,2,3] for a many-fermion system, the quasiparticle excitations are integrated out successively with an energy scale decreasing toward the Fermi energy. This generates an energy-scale dependent effective interaction from which the leading interaction processes and correlations can be read off. Although most implementations involve a number of approximations, the flow of the vertex functions of the effective theory can be derived from an exact RG equation for a generating functional [4,5]. This powerful equation sparks the idea that a much more controlled approach to the low-energy phase than obtained previously is indeed possible.
In many applications for correlated electron systems one is primarily interested in changes of the ground state toward symmetry-breaking states. If there is such a symmetry change of the ground state, a certain bilinear B ij in the fermion operators c j with some quantum numbers i and j takes a nonzero anomalous (i.e. violating the orginal symmetry of the action) expectation value at a low RG scale. For a second order transition, this occurs continuously as RG scale or temperature are lowered. Then the tendency toward symmetry breaking is signaled in the RG flow by a divergence of a certain component of the effective interaction at a nonzero energy scale Λ c . In correspondence with this, the zero-frequency susceptibility for the coupling to the bilinear operator B ij diverges in power-law form at the same critical scale Λ c . If we view the expectation value B ij as a bosonic field φ (e.g. by a Hubbard-Stratonovich transformation), the inverse susceptibility appears in the quadratic part of the action for the φ-field. The symmetry breaking occurs at the temperature or energy scale where this quadratic part changes its sign, i.e. where the effective potential for φ reverses its curvature at φ = 0 from positive to negative, such that the potential minimum occurs at a nonzero value of φ. When the energy scale or temperature is lowered further, this minimum moves further away from zero and finally saturates unless collective fluctuations move it back to zero and restore the symmetric state (see, e.g., [6]).
As indicated in this description, the passage from the symmetric state at high energy scales or temperatures to the symmetry-broken state at low energies is most easily described in a bosonic language, and such methods are being developed ( [7,6,8]). However, in many interesting examples in many-fermion systems it is not clear from the outset which type of boson or which fermion bilinear acquires a nonzero expectation value. Hence, a purely fermionic description is a less biased approach, as it allows the system to develop its intrinsic correlations without any prejudice of grouping certain fermions in pairs and treating these in a saddle-point approach. Recently, some of us have presented a purely fermionic functional RG method based on [9] which allows one to continuously flow into the symmetry-broken state, requiring only an infinitesimal symmetry-breaking in the initial conditions. This approach has been shown to work correctly for models where mean-field solutions are exact, like the reduced Bardeen-Cooper-Schrieffer (BCS) model [10] or a similar charge-density-wave model with discrete-symmetry breaking [11,12]. Extensions for more general models will follow.
Yet, these new fRG extensions face two problems. First, there is still an almostdivergence of the interactions at the critical scale. Large values of the effective interaction are problematic for the extension to more general models, as they might cause a breakdown of the perturbative approach. A true divergence is avoided by the small explicit symmetry breaking put into the initial conditions. But in order to stay close to the symmetric model which one typically wants to study, some components of the interaction in the symmetry-breaking channel still reach large values whose feedback on other channels needs to be understood. While this problem may be solved either by a compromise with an adequate size of the initial symmetry breaking or by understanding the potentially rich and interesting physics of the effect of these large interactions on the flow, a second problem is more fundamental. The currently available fRG methods only test the stability of the symmetric phase locally around the origin in the space of anomalous expectation values. If there is no divergent susceptibility or no curvature change in the effective potential of a bosonic field, the fRG flow will remain finite and the anomalous expectation values will vanish continuously when the small initial symmetry-breaking is sent to zero. This in turn would be interpreted as stability of the symmetric phase and absence of spontaneous symmetry breaking. On the other hand, simple mean-field treatments of many models indicate first-order transitions into the symmetry-broken state. This occurs for example in the two-dimensional t-t ′ Hubbard model near half filling for the (commensurate) antiferromagnetic spin-density wave transition [13,14] and in a model with a particular forward scattering for the transition into a state with broken fourfold symmetry of the Fermi surface [15]. The fRG flow for these parameters in the almost symmetric state finds an enhanced but not a diverging susceptibility, hence missing the symmetry change of the ground state. The fRG method simply does not see the thermodynamic-potential minimum for nonzero order parameter if it is separated from the symmetric minimum by a finite energy barrier, as in the case of certain first-order transitions. What is needed is a method that, intuitively described, starts the flow further out in the symmetry-broken range and which then flows into the correct minimum corresponding to a spontaneously broken state of a model with a fully symmetric action.
In this work, we present a counterterm extension of the fRG method which is able to find symmetry-broken states not captured by the treatment of [11,12]. These would occur via first-order transitions upon varying temperature, chemical potential, or infrared cutoff. The procedure invokes a flow differing slightly from Wilson's momentumshell RG. The propagator of the fermionic excitations is increased continuously from zero to its full strength with the same rate for all excitation energies. This flow, also dubbed interaction flow, has been tested previously for the 2D Hubbard model without symmetry breaking and gives analogous results to those obtained with momentum-shell or temperature-flow schemes [16]. It can be regarded as an infinitely smooth cutoff, i.e. the opposite extreme to a sharp cutoff, for which all modes above a running energy scale Λ are included with full weight and all modes below Λ do not contribute.
Counterterms have been employed in the context of correlated fermions before, notably in calculations for the symmetric state of interacting electron systems. The early work [17] by Feldman and Trubowitz considers systems with spherical free Fermi surfaces. Unphysical infrared divergencies of perturbation theory Feynman graphs are removed by adjusting the chemical potential on the propagator lines to the interacting Fermi surface, countering this by an appropriate addition to the self-energy. The addition is determined by renormalization-group methods. This procedure can be viewed as a shift from the chemical potential into the interaction. Feldman, Salmhofer, and Trubowitz [18] employ it for more general Fermi surfaces which necessitate a counterterm that is a function in momentum space. They identify and find bounds for the most singular contributions to this couterterm function under slightly more restrictive conditions [19]. The authors prove regularity properties of the counterterm function in [20], which they use in [21] to study the invertibility of the mapping induced by the counterterm which transforms the interacting to the corresponding free Fermi surface. The analysis of [21] is restricted to space-inversion-symmetric dispersions.
Counterterms have also been used in the study of symmetry breaking. Feldman and Trubowitz [22] is devoted to flows with symmetry-breaking in the Cooper channel for a spherically symmetric electron-phonon system. Counterterms are again used to remove unphysical divergencies as in [17]. The BCS self-consistency equation is recovered as an approximation from the flow. Neumayr and Metzner [23] study the interplay of d-wave superconductivity and Fermi surface deformations by a counterterm technique. The effective interactions causing the symmetry breaking are constructed order by order in perturbation theory. The counterterm is determined by setting the self-energy on the Fermi surface equal to zero, giving rise to a self-consistency equation.
It is shown in [11,24] that such self-consistency equations are also naturally contained in the fRG method applied here and in [11,12,24]. This method is exact for the models considered, as is the self-consistency treatment. However, it should be noted that for general models, self-consistent calculations just as fRG calculations are not exact anymore and their results may differ. It may depend on the specific case whether self-consistent calculations or the broader and more fundamental fRG method give the more useful results.
The role of counterterms in the renormalization group literature has focussed on the removal of divergencies. In this work, the counterterm serves the purpose of selecting a certain symmetry-broken configuration. It also opens a gap in the spectrum of the system, which removes the divergencies on the Fermi surface treated in the literature. This paper is structured as follows. In section 2, we introduce the model, presenting and explaining the mean-field solution in section 3. We outline the renormalizationgroup scheme and introduce our conterterm extension in section 4, where we furthermore provide a numerical RG analysis of the model from section 2. We discuss our work and provide an outlook in section 5.

Model
To illustrate the utility of the fermionic functional renormalization group (fRG) scheme introduced in section 4, we consider a tight-binding model on a d-dimensional hypercubic lattice with an interaction restricted to momentum transfers of Q = (π, π, . . .). The Hamiltonian reads In this work, ξ( k) = ǫ( k) − µ includes a chemical potential µ, which can be adjusted to produce a system exhibiting a first-order phase transition. We assume in the following that the dispersion fulfills the nesting relation ǫ( k) = −ǫ( k − Q). The kinetic energy ξ( k), interaction, and external-field parts related to ∆ ext are also introduced and briefly discussed in Chapter 2 of [12]. For V 0 > 0 the interaction term can lead to chargedensity-wave ordering with wavevector Q. This will be discussed below employing a mean-field approach. In fact, due to the restricted form of the interaction of the model, the mean-field treatment and the RG approach described below become exact in the thermodynamic limit, where the number of lattice sites N → ∞ (for a diagrammatic argument, see [24]). The same two-sublattice charge modulation is induced by the last term, which breaks the translational symmetry of the Hamiltonian by coupling an alternating charge field to the fermions. In the fRG treatment, Σ i will be used as the initial condition for the off-diagonal self-energy. The counterterm ∆ c will be included in the bare propagator to prevent this from affecting the physics of the system. We always set Σ i = ∆ c which guarantees a cancellation at the end of the fRG flow. ∆ ext allows us to study the effect of an external field and hysteresis phenomena.

Mean-field treatment
The low-temperature state of the model without explicit symmetry breaking, i.e. for (∆ ext + Σ i − ∆ c ) = 0, depends on the chemical potential µ. For µ = 0 and 0.5 particles (per spin orientation) per lattice site, the Fermi surface is perfectly nested and mean-field theory finds a second-order transition toward the charge-density wave state at a critical temperature T c . In the ordered state, ∆ = V 0 N k c † k c k+ Q becomes nonzero. Similarly, the RG without self-energy corrections would find a runaway flow at a nonzero critical scale Λ c . The runaway flow and the second-order mean-field transition are removed by a sufficiently large µ of the order of T c (µ = 0). For such µ, a local stability analysis of the symmetric state via an expansion around ∆ = 0 does not detect any instability. However, below a transition temperature T t , a mean-field search reveals global minima of the thermodynamic potential as function of ∆ which do not develop out of the local minimum at ∆ = 0 when the temperature is lowered. As the temperature is lowered, the system jumps from ∆ = 0 to the nonzero value, undergoing a first-order phase transition. The resulting phase diagram is shown in Fig. 1  potentials as functions of ∆ for µ = 0 and for first-order cases with µ = 0 are shown in Fig. 3.
To understand how such a simple model can have a first-order transition, we consider the self-consistency equation where f (x) = (exp(x/T ) + 1) −1 is the usual Fermi distribution and the eigenvalues are where R = (ε k ) 2 + ∆ 2 . (2) is obtained from the differential of the grand canonical potential Ω by demanding ∂ ∆ Ω = 0. To be specific, we base the following discussion on a two-dimensional square lattice with nearest neighbor hopping t (which is used as the unit of energy for the remainder of the paper) at T = 0. The self-consistency equation is plotted in Fig. 2. The RHS has a kink at ∆ = |µ| because one of the branches of the eigenvalues changes its sign on the whole umklapp surface for ∆ = |µ|, cutting out this part of the integration area of the RHS. It can be shown that at the kink the slope diverges as (1/ |µ| − ∆) log(|µ| − ∆) for ∆ ր |µ|. At T = 0, the RHS does not depend on µ for ∆ > |µ| since µ enters the RHS only via the Fermi functions.
is always positive (negative) for ∆ > |µ|, the thermodynamically stable solution is always half filled in this case. The dependence of the filling on ∆ should not be neglected, as is sometimes done in the literature [25]. In fact, mean-field theory already predicts the breakdown of the antiferromagnetic gap away from half filling in the homogeneous solution, in agreement with the quantum Monte-Carlo calculations of [25]. For a detailed discussion, see [14]. Since (2)  intersection points are extrema (or saddle-points) of Ω. If all intersection points are extrema of Ω(∆), the outermost one must be a minimum to guarantee the existence of a stable solution. By construction, the integral of the bisector minus the RHS of (2) is the energy gain through the opening of a finite gap. Thus, comparing the two shaded areas of Fig. 2 determines which of the two minima given by the intersections has lower energy. Plots of the thermodynamic potential are shown in Fig. 3. For small |µ| > 0, a local minimum develops at ∆ = 0, which becomes the global and finally the only minimum if µ is increased further. Numerical study of (2) and Ω for various T and µ leads to the phase diagram of Fig. 1

Renormalization group
In this section, we outline an extension of the one-particle irreducible (1PI) fermionic functional renormalization group (fRG) scheme [4] with the goal of studying first-order phase transitions. We start from the flow equations as suggested by Katanin [9], studied for the case of the breaking of a continuous symmetry in [11] and for the case of the breaking of a discrete symmetry in [12]. We use the interaction-flow method [16], which permits the fruitful employment of counterterms. We show flows for first-order and second-order phase transitions, discussing the effect of the counterterm on possible approximations. We also study the effect of an external field on the flow, identifying a criterion for unphysical starting points.

Setup
The basic setup of a 1PI fRG calculation with Katanin's modification is outlined for two different models in [11,12]. To understand the problem solved in the following, we consider the order-parameter dependence of the thermodynamic potential in the plot for µ = 0 in Fig. 3. In [11,12], a small external field is included as the initial condition for the order parameter. This external field is relevant in Wilson's fixedpoint picture, but does not appreciably change the physics away from critical scales. It biases the fRG calculation toward the thermodynamic potential minimum at finite ∆ (see Fig. 3 µ = 0), which becomes the endpoint of the flow. This method is unable to deal with first-order phase transitions where the thermodynamically stable state is separated from the symmetric state by an energy barrier. A typical thermodynamic potential as function of the order parameter is depicted in the plot for µ = −0.2t in Fig. 3. Because the external field must be small, the flow always veers toward the symmetric thermodynamic potential minimum. Other minima are inaccessible, even if their thermodynamic potentials are smaller. Furthermore, the method cannot be used to study hysteresis effects appearing if systems exhibiting second-order phase transitions (see the graphs for µ = 0 and µ = −0.1t in Fig. 3) are placed in an external field. The external field destroys the axis symmetry of Ω(∆), but always biases the flow toward the global minimum of the thermodynamic potential, making it impossible to study the metastable configuration. Thus, the challenge encountered here is to set up an fRG scheme with a parameter that allows the selection of any stable or metastable configuration as the endpoint of the flow. The basic idea that will be exploited in the following is to include in the calculation a counterterm ∆ c of arbitrary strength that cancels an equally strong external field, but only at the end of the flow. In the fRG calculation, the external field is taken into account as the initial value Σ i of the order parameter Σ, while the counterterm ∆ c is added to the naked propagator Q 0 . Schematically, the matrix Green's function of the system reads where χ is the cutoff function and Q 0 is the inverse bare propagator without counterterm. For any χ ≡ 1, counterterm and self-energy do not cancel each other. However, the special case χ ≡ 1 is approached at the end of the RG flow. Thus, the vertex function of the symmetric model is obtained. We first consider the case of χ = χ Λ ( k) describing the sharp momentum-shell cutoff of conventional RG schemes (see [1]). Then, only momenta for which χ Λ ( k) = 1 are taken into account in the RG flow. For these momenta ∆ c and Σ cancel exactly at the start of the flow, causing Σ to remain stationary. This is due to the right-hand side of the 1PI flow equation for the order parameter being proportional to the anomalous part of the Green's function which vanishes if ∆ c and Σ cancel. This argument generalizes to the case of multiple order parameters, as can be seen by considering the adjugate matrix of the Green's function. Consequently, ∆ c and Σ also cancel at all stages of the flow. In the case of a second-order transition, the fRG will diverge at a nonzero scale since the spectrum remains unchanged. In the case of a first-order transition, the flow remains finite, signaling incorrectly the stability of the symmetric phase. We conclude that we have to employ a softer cutoff function.
Next, we consider the softest possible cutoff function, provided by the interaction flow scheme [16]. This scheme works by linearly turning on the interaction with the flow parameter (note the similarity to the treatment by Mahan [26]) while turning on the counterterm with the square root of the flow parameter. By rescaling the fields, this is found to be equivalent to employing a cutoff function χ which is constant at all momenta, zero at the start of the flow, and one at the end of the flow. With this, the counterterm and the external field do not cancel in (4) until the end of the flow, and all configurations become accessible, as will be illustrated in the following.
For a general cutoff function χ Λ with flow parameter Λ and a general model with two-particle bare interaction, the flow equations were obtained in [9]: where G Λ = Q/χ Λ − Σ Λ −1 and the traces run over all degrees of freedom of G Λ . The contribution of O(V 3 ) to (6) vanishes in the case of the Hamiltonian (1) when taking the thermodynamic limit due to the restricted momentum structure of the interaction. For the interaction flow, the cutoff function is where Λ runs from zero to one. Considering ε = ξ + µ as mapping the Brillouin zone onto R, denoting the effective gap ∆ f := −∆ c +χΣ, R := ε 2 + ∆ 2 f , and E ± := −µ±R, we obtain the Λ-dependent inverse free propagator and the Λ-dependent Matsubara Green's function We have omitted Λ-dependencies in the notation for brevity, as we continue to do in the following. (8) and (9) allow the calculation of the single-scale propagator S = −GQG.
Tracing over all degrees of freedom, analytically evaluating the Matsubara sum, and multiplying with the effective interaction and diagrammatic factors, we obtain the flow equatioṅ Similarly treating the right-most diagram of Fig. 3 from [12] or alternatively evaluating (6) for the special case (1) yields the interaction flow equatioṅ We furthermore need to calculate the thermodynamic potential using the functional renormalization group. From equation (48) of [4], using Ω i.a. = T γ 0 and equations (39) and (45) of [4], we determine the flow equatioṅ If the flow is started at the thermodynamic potential of the non-interacting system, the full thermodynamic potential is recovered. Since we will always do so, we have written Ω instead of Ω i.a. in (12). The right-hand side of (12) can again be drawn as a 1PI diagram. Evaluating the Matsubara sums, (12) readṡ (10), (11), and (13) constitute a closed system of integro-differential equations which can be numerically solved given appropriate initial conditions.

Flows for first-order phase transitions
We first consider a system exhibiting a first-order phase transition at a transition temperature T t . Studying flows below T t (see Fig. 4), we notice two strong attractors. By its lower thermodynamic potential at the end of the flow, one of them can be identified with the stable, symmetry-broken configuration. The values for the order parameter, effective interaction, and thermodynamic potential difference reproduce the exact mean-field results. Note that the final values do not depend on the magnitude of the counterterm. The dependence of the results on the external field, known from [12,27,11], is eliminated. Considering Fig. 4(a), we see that there is a separatrix between the effective gap flows to the metastable configuration and the flows to the stable configuration. For more complicated problems than (1), the flow probably has to be stopped at Λ < 1 because of the rise of large effective interactions which are not correctly taken into account by our method. Therefore, flows which exhibit weak changes of the effective interaction away from Λ = 1 and weak changes of the order parameter close to Λ = 1 are most promising for determining order parameters in more complicated models. We see that the flows closest to the separatrix are not optimal flows since ∆(Λ) still shows a large slope at Λ = 1. This final slope is smallest for a counterterm of roughly twice the magnitude of the order parameter. Such a flow would yield an excellent approximation of the order parameter even if stopped at Λ = 0.75.
Considering the flows of the effective interaction in diagram 4(b), we do not find a separatrix as in (a). Instead, the flows corresponding to the ones close to the separatrix in (a) develop a shoulder which becomes a maximum for certain values of the counterterm (for an example of such a maximum, see diagram 5(b)). The flows exhibiting the strongest effective interactions correspond to the flows giving the worst approximations for the order parameter if terminated prematurely. This entails for applications to models where this method is not exact that carefully chosing the counterterm can significantly simplify the calculation and improve the accuracy of the results.
The flows of the thermodynamic potential shown in 4(c) again exhibit separatrix behaviour. This is in contrast to the flows above T t . For the thermodynamic potential, the flows closest to the separatrix yield the best approximations if terminated prematurely. However, it is apparent from 4(c) that the flow must be continued until the  external field is completely compensated by the counterterm to obtain reliable results for the thermodynamic potential. Flows above T t are illustrated in Fig. 5. Again, we clearly discern two attractors. In this case, however, the thermodynamic potential flows of the symmetry-broken-phase attractor cross the thermodynamic potential flows of the symmetric-phase attractor. The symmetric phase is therefore thermodynamically more stable. Apart from this, the flows behave similarly as those below T t . For inconveniently chosen counterterms, the maxima in the effective-interaction flows are clearly visible in 5(b). Nevertheless, the flows reproduce the exact mean-field results.

Flows for second-order phase transitions
The flows for second-order phase transitions behave similarly to the flows for first-order phase transitions, but only a single attractor appears, as can be seen in Fig. 6(a). counterterm of twice the non-approximated value of the gap.

External field
If the external field ∆ ext in (1) is zero, the half-filled system below T c exhibits two degenerate stable configurations (local minima of the thermodynamic potential) distinguished by the sign of the order parameter. If ∆ ext is non-zero, this degeneracy is lifted. The fRG scheme outlined above allows us to select the endpoint of the flow independently of the external field by appropriately setting the counterterm, in contrast to the scheme from [11,12]. The flows toward the metastable configuration are shown in Fig. 7 (broken lines). As the external field is increased toward a critical value, the order parameter for the metastable configuration vanishes and the corresponding local minimum of the thermodynamic potential disappears. Beyond this critical value, the flows for negative counterterms become divergent as the effective gap reaches values close to zero, exposing low-energy modes. Such a divergence can be used as an indicator for an adversely chosen starting point. The flows for positive counterterms (full lines in Fig. 7) always attain the stable solution. Calculating the flows as for Fig. 7, but for negative values of the external field ∆ ext , we obtain the hysteresis curve of Fig. 8.

Discussion and outlook
We have described a counter-term scheme within the framework of the functional renormalization group concept for interacting fermions. The scheme is able to detect ordered ground states which, in terms of the effective bosonic potential for the ordering field, occur as first-order transitions. Second-order transitions are captured as well. Therefore, the new scheme remedies one shortcoming of conventional functional RG techniques (e.g., [28,29,30]) for interacting fermions which are blind with respect to symmetry-broken states which are separated by an energy barrier from the symmetric state. A second shortcoming that is addressed is the need for a finite external field [11,12], which does not arise in our approach. As a third improvement, at least for discrete-symmetry breaking away from the critical temperature, our scheme can circumvent the development of large interactions in the flow. The validity of the new scheme has been demonstrated in a mean-field model, where our truncated version of the fRG equations and mean-field theory are known to be exact. More generally, the scheme can be used to target individual minima of the effective bosonic potential, for example to study metastable states in the presence of an external field. As an example, we have shown the calculation of a hysteresis curve. The fRG scheme can be applied to more realistic models such as the Hubbard model as well. Here a new complication arises, as the truncated form of the RG equations used is no longer exact. Hence large values of the coupling constants could spoil the validity. This is a general problem for perturbative functional renormalization group schemes which can for example be circumvented by partial bosonization at the expense of introducing a dynamical bosonic field [6,7,8,31]. The counter-term technique presented here allows to avoid this problem in certain cases. As demonstrated and described above, in the case of the breaking of a discrete symmetry, no large couplings occur unless one is close to a curvature change of the effective bosonic potential. Hence for these cases, an inclusion of non-mean-field-type interactions should be possible without much trouble. For broken continuous symmetries, Goldstone modes may again lead to large couplings. It will be interesting to see whether the Goldstone physics only enters the flow in a critical way toward the end of the flow, when the off-diagonal selfenergy has already converged to the final value. The results above indicate a rather quick convergence for most values of the counter-term. These issues will be the subject of future work.
We emphasize however that even in the case that the flow cannot be followed over the full range, the present scheme should be useful for the detection of symmetry-broken states which are separated by an energy barrier from the symmetric state. For example, if the conventional analysis signals a locally stable symmetric state, the counter-term technique can be used to scan non-zero values of the suspected order parameter. If there is a lower thermodynamic potential minimum further out, the self-energy shows a flow toward non-zero values. Further, there should be a local curvature change of the effective bosonic potential somewhere between 0 and the value related to the global minimum, which is detected as a runaway flow of the interactions in the corresponding channel. Hence, the method is useful at least for qualitative scans of the bosonic potential landscape. More generally, variations of the proposed scheme can be used for the perturbative analysis around states which differ from the free ground state of the model. This as well should lead to numerous new applications.