D ec 2 01 9 Precise Wigner-Weyl calculus for lattice models

IV. Weyl symbol derivation 12 A. Lattice Wigner-Weyl calculus through the series expansion 12 1. Formal Definition 12 2. Integral representation for QW 14 3. Star identity I: AW ⋆ BB = (AB)B 15 4. Star identity II: AW ⋆ BW = (AB)W 16 5. Trace definition 17 6. Trace of a product 18 B. Wigner-Weyl formalism on a doubled lattice. One-dimensional case. 19 1. Operators on auxiliary lattice O and their Weyl symbols 19 2. W -symbol and its properties 20 3. Moyal product 20 4. Trace and its properties 22 5. Groenewold equation 23 C. Wigner-Weyl formalism on an extended lattice. The multidimensional case. 24 1. The auxiliary extended lattice 24 2. Moyal product 25 D. A brief word on motivation 26

Topology plays an important role both in relativistic quantum field theory and in condensed matter physics. One of its most well-known manifestations is nowadays associated with the Quantum Hall Effect (QHE). The first topological representation of the Hall conductivity has been proposed for homogeneous systems in the presence of constant magnetic field. It is called the TKNN invariant [1]. The TKNN invariant is also used for description of homogeneous Hall insulators that possess intrinsic anomalous quantum Hall effect (AQHE) in absence of external fields. The TKNN invariant is given by an integral of Berry curvature over the occupied energy levels. The corresponding value is not changed when the dependence of energy on momentum varies smoothly, see [2][3][4][5] and reviews [6,7].
Momentum space topological invariants were widely studied within the context of condensed matter physics since then. First of all, in the absence of interactions the Hall conductivity for the intrinsic AQHE has been expressed through the Green's function [8,9] (see also Chapter 21.2.1 in [10]). For two-dimensional topological insulators the Hall conductivity is given by It appeared later that the simplest possible topological invariant in momentum space is composed of the two-point Green's function of a homogeneous fermionic model, and is responsible for the stability of the Fermi surface [10] N 1 = tr C dp l 2πi G(p 0 , p)∂ l G −1 (p 0 , p).
Here ∂ ρ ≡ ∂/∂p ρ , while C is a closed path, which encloses the Fermi surface in momentum space. The topological stability of Fermi points is protected by another invariant, which reveals correspondence with the above mentioned expression for the Hall conductivity [10,11] Here S is a surface surrounding the Fermi point. These and more involved constructions are widely used in condensed matter physics theory [12][13][14][15][16]. Topological invariants in momentum space are responsible for the gapless nature of fermions at the edges of topological insulators [17,18] and in the bulk of Weyl semi-metals [10,16]. The fermion zero modes are related to various topological defects in 3 He-superfluid, and are protected by the similar topological invariants [19]. In the high energy theory the topology of momentum space has been discussed, for example, in [8,10,[20][21][22][23][24][25][26][27]. By construction, the TKNN invariant, as well as (1)(2)(3), is defined for systems without interactions. It allows to deal with physically homogeneous configurations even if the gauge potential corresponding to the constant magnetic field is coordinate dependent. In the latter case the magnetic Brillouin zone is to be introduced in order to define the TKNN invariant. The important question is what will happen to the topological representation of Hall conductivity if the system becomes non-homogeneous and interacting. The inhomogeneity can be given by the electric field of impurities, elastic deformations, or by the variations of magnetic field. Up to now the representation of Hall conductivity in the general case as a topological invariant has not been constructed. However, several particular cases have been investigated.
Although TKNN and (1) were obtained for the case when there are no interactions, according to the common lore they both remain valid for the models with weak interactions, if the two-point Green's function is taken with the interaction corrections. This statement was proved rigourously in 2 + 1D relativistic Quantum Electrodynamics [28,29]. Recently [30] the AQHE has been considered for the models of rather general types. There the mentioned statement has been proven, at least, in one-loop approximation.
For the systems with homogeneous magnetic field the influence of interactions on the Hall conductivity has been discussed widely (see, for example [31][32][33][34] and references therein). The general output of these studies is that the weak interactions do not affect the value of Hall conductivity. The case of varying magnetic field has been considered recently in the series of papers written with the participation of the authors on the present paper. In [35] a new expression for the Hall conductivity has been proposed. It was given in a topologically invariant form, composed, however, of approximate Wigner transformation of the two-point Green's functions. The adaptation of Weyl symbol used in [35] is described in Appendix B of the present paper. Graphene in the presence of inhomogeneities invoked by elastic deformations was considered in [36]. In [37] the proof was presented that in the presence of interactions the Hall conductivity is still given by the expression proposed in [35] with the complete two-point Green's function. The mentioned version of the Wigner-Weyl calculus has been also discussed in [38]. This approximate formalism allows for calculation of the Hall conductivity in the case when the magnitude of magnetic field is much smaller than several thousands Tesla, while the wavelengths are much larger than 1 Angstrom. Under these conditions the average Hall conductivity of a two dimensional systems is given by σ H = N 2π , with Here T is temperature that is assumed to be small, |A| is the (infinite) area of the system, G C (x, p) is a type of Wigner transformation (see Eq. B2) of the two-point Green's functionĜ, being in its turn the inverse of the Dirac operatorQ. As it was mentioned above those quantities are defined using the specific version of the Wigner-Weyl calculus designed especially for the consideration of the QHE in weak magnetic fields. In general case, it gives only an approximate expression for the conductivity. The original Wigner-Weyl calculus stemmed from the works of H. Groenewold [39] and J. Moyal [40], who considered the alternative formulations of quantum mechanics in infinite continuous space. This calculus accumulated the ideas of H. Weyl [41] and E. Wigner [42]. According to this formalism it is possible to use the Wigner distribution instead of the wave function, Weyl symbol instead of the operators of the observable quantities, and Moyal product in lieu of operator one [43,44]. The original version of Wigner-Weyl calculus has been widely applied to quantum mechanics [45,46]. Certain modifications of Wigner-Weyl formalism were proposed within the context of both high energy field theory and condensed matter physics [47][48][49][50][51][52][53]. The Wigner distribution has been used in QCD [54,55]. It has been considered also in quantum kinetic theory [56,57], in and in the noncommutative field theories [58,59]. There were several attempts to apply various versions of the Wigner-Weyl calculus in a variety of physical problems including cosmology [60][61][62].
Formulation of the Wigner-Weyl calculus on a lattice (or on a torus) faced certain difficulties, and still this field of research is not settled. It started from the works of Schwinger [63], and continued throughout the century. We would like to mention the works of F.Buot [53,64,65], Wooters [66], and Leonhardt [67] devoted to the physical applications, and the works of Kasperowitz [68], and Ligabó [69] on the mathematical formulations of the discreet Wigner-Weyl calculus. The closely related area of research is the deformational quantization [70,71]. Several important general facts were proven in this field, including the existence of a trace operation [72], and perturbative construction of the Moyal product [73]. However, the mentioned results are either too abstract or are possessing certain caveats, impeding the real physical applications.
The approximate version of the lattice Wigner-Weyl calculus, based on G C from above is described in Appendices B and C of the present paper. The application of this calculus, as already mentioned, is limited to slowly varying weak fields. However, even in this form it appears to be powerful enough to describe response of various nondissipative currents to external field strength [74][75][76][77][78][79]. In many cases such a response is expressed through topological invariants, which do not changed when the given model is modified smoothly. In this way the absence of the equilibrium chiral magnetic effect [80] has been proven [78], the AQHE has been investigated [79], the chiral separation effect (CSE) has been reestablished [81] within the properly regularized models (using lattice field theory) [74,76]. This formalism has also been applied to the high density QCD [77], and to the field systems in the presence of both gravity and magnetic field (the so-called scale magnetic effect) [75,82].
In the present paper we consummate the development of the Wigner-Weyl calculus for lattice models, aiming primarily at the systems in external magnetic field. We present the precise version of this calculus, which allows us to obtain topological expression for the Hall conductivity for the models defined on rectangular lattices without any limitations on the magnitude of external magnetic field and on the rate of its variations as a function of coordinates.
In the next Section we formulate our main result (upon introducing all necessary notation) in the form of four Theorems, and prove them one by one in the consequent Sections III and IV.

II. MAIN RESULT
Before the formulation of our main result let us introduce notations and definitions used throughout the paper. We assume the relativistic units with c = = 1, if not explicitly stated otherwise.

A. The Hilbert space
In the one-dimensional models the physical lattice is denoted by O, and for simplicity of the forthcoming considerations we write it as The first Brillouin zone is Apart from that, we also consider an extended lattice O Its first Brillouin zone is M = (−π/ℓ, π/ℓ]. The additional lattice O ′ can be considered as a translation of the physical one The usual properties of the physical states are assumed The latter is a Kronecker symbol since our coordinate is discreet. Accordingly, we have the Fourier decomposition In the present paper we will also consider a D-dimensional rectangular lattice. It will still be denoted by O. Its points and the corresponding Brillouin zone are We assume that the theory to be dealt with is Euclidean equilibrium one. For brevity we assume the discretization of imaginary time together with the discretization of space coordinates. However, in the applications to condensed matter physics typically imaginary time is continuous. Therefore, speaking of certain properties of our constructions we will take off the discretization of imaginary time. This will be pointed out explicitly. If the absence of the discretization of time is not mentioned explicitly, such a discretization is implied. The immediate generalization of (9) and (10) to multi-dimensional case is assumed to hold. The auxiliary extended lattice O in D-dimensional case, and the corresponding extended momentum space M are given by: To avoid the clatter of notation, we will not use bold-face notation in what follows. The physical Hilbert space H is defined using the basis |x (x ∈ O) that obeys Eq. (9). All operators to be discussed below are assumed to be defined on this Hilbert space.

B. Lattice Weyl symbol
An abstract lattice Weyl symbol is defined in the following way.

Definition 1 By the Weyl symbol of a linear operatorÂ acting on the Hilbert space H we understand the map
such that there could be defined a ⋆-product between Weyl symbols of any two operators, as well as trace operation Tr mapping L 2 (M × O) to complex numbers, C, satisfying the following conditions 1. Star product identity

Weyl symbol of identity operator
By tr we understand the trace of the operator itself in the original Hilbert space.
An explicit construction of such Weyl symbol leads to is a Weyl symbol of operatorQ in the sense of the Definition above, while In this case Here tr stands for the trace over the inner symmetries (if any), which we will omit in the future. The ⋆-product is the original Moyal one This pseudo-differential operator acting on functions of discreet coordinate x is understood either as an integral operator, see (132), or as acting on analytical continuation, see the end of the Section. We will show that the integral representation (17) is equivalent to a series representation for the Weyl symbol, which can be formulated as Theorem 2 For operators, possessing a series representation of a function of two variableŝ its Weyl symbol can be defined as with Moreover, provided that The notation above is evident in one-dimensional case. In D dimensions, the indices are to treated as multi-indices, n = (n 1 , ..., n D ), etc., with (i∂ p ) n = i (i∂ pi ) ni and C n k = i C ni ki . A sum over a multi-index is understand as a sum over all ot its components, and a condition of the type of n < k where n and k are the multi-indices, implies that the sum is over all those values of components that obey n i < k i for any i.
Altogether, throughout the paper we will use three different symbols of operators: Q W given by (17), Q W of (22), and also the symbol of operator proposed by Buot [64], The relation between the former two is given by the theorem above, while the latter does not satisfy our conditions for a proper Weyl symbol. However, it will be useful for establishing the properties of Q W and Q W . Buot's symbol is investigated in details in Appendix A.
We work with the model defined on the lattice. Therefore, the application of the derivatives in x present in the star product is questionable at a first look. The problem is that the analytical continuation of function defined on discreet (countable) set of points is not unique, generally speaking, [83]. However, it appears that the application of the ⋆-product to the symbols of operators does not depend on the choice of the analytical continuation, provided that we know their values on the extended lattice M.
Let us consider this problem for simplicity for one-dimensional lattice. Generalization to the multidimensional case is evident. The symbols we deal with (Q W , Q W and Q B ) are all π/ℓ-periodic in p. Thus they can be represented as with certain functions q n (x) defined on the lattice x ∈ M = ℓZ. Constructing an analytical continuation q n (z) : q n (ℓZ) = q n (ℓZ) (28) we investigate the ⋆-product of two such objects to see that it does not depend on the choice of the analytical continuation in (28), as the star acts by shifting the arguments of q n and s m in such a way that those arguments always belong to O. So, the application of ⋆-product is actually independent of the particularities of the continuation requiring only the knowledge of the function's values on x ∈ O.

C. Global current and averaged conductivity
In Euclidian space-time the partition function of a physical system (defined by its Dirac operatorQ) is given by with the action where the integration measure and normalization are understood to be chosen appropriately for the model under consideration. Note, thatψ stands for the second independent integration Grassman variable, and has nothing to do with the Dirac spinor conjugated operator, used in the operator formalism.
Using (1) and Definition 1, the action can be written as where ρ W (x, p) is the Weyl symbol of an operatorρ with Grassman-valued matrix elements ρ(x, y) ≡ψ(x)ψ(y) = x|ρ|y Then, we can formulate and prove the following two theorems Theorem 3 In the system described by (31) the total electric current of the sample is given bȳ where

Moreover, this expression is a topological invariant: it is not changed under small variations of the lattice Dirac operatorQ.
Theorem 4 In the system described by (31) subject to constant external electromagnetic field F (E) lm , the total electric current averaged over the volume |V| of the sample is given by The Hall conductivity, in turn,σ is a topological invariant.
Above β = 1/T is the inverse temperature, T is assumed to be small. We also assume that the system is gapped, so that the integral in Eq. (56) does not contain divergencies.

III. WIGNER-WEYL QUANTUM FIELD THEORY
In this section we describe an arbitrary non-interacting fermionic model defined on the lattice. The basic quantities are represented in terms of the Weyl symbols of the corresponding operators. We assume here, that there exists a Weyl symbol obeying the conditions given in the above Definition 1. Using these properties (but not the particular form of the Weyl symbol) we will derive in this section the expression for the quantum Hall conductivity. It is expressed through the Weyl symbol of the fermion propagator, and appears to be a topological invariant as soon as its expression does not contain divergences (in practice this occurs if the system is gapped, and the Fermi energy lies in the gap). Thus in this section we will prove Theorem 3 and Theorem 4. The explicit construction of the Weyl symbol of operator (and the proof of Theorem 1 and Theorem 2) will be given in the next Section.

A. Lattice model in momentum space
Following [35,38,78,79] we start with the lattice tight-binding fermionic model with the partition function of the following form Here D x,y is a matrix that depends on the discrete lattice coordinates x, y ∈ O. ψ,ψ are multi-component Grassmannvalued fields defined on the lattice sites. The corresponding indices are omitted here and below for brevity. Normalization of the functional integration measure is assumed to be chosen appropriately for the problem at hand. In the simplest case of a uniform system such a partition function can be rewritten in momentum space as follows where integration is over the fields defined in momentum space M. |M| is its volume, D is the dimensionality of space-time (we will omit it whenever this does not lead to misunderstanding),ψ(p) and ψ(p) are the anti-commuting multi-component Grassmann variables, now defined in momentum space. Without loss of generality we assume that time is discretized, so that momentum space is compact, and its volume is finite. In condensed matter physics the imaginary time typically is not discretized, the corresponding partition function can be obtained easily as the limit of Eq. (38) when the (imaginary) time spacing tends to zero. The partition function of Eq. (39) allows to describe noninteracting fermionic systems corresponding to matrix Q(p) (that is the Fourier transform of the lattice tight-binding operator D x,y ). The meaning of Q(p) for the lattice models of electrons in crystals is the inverse propagator of Bloch electron.
Introduction of an external gauge field A(x) defined as a function of coordinates effectively leads to the Peierls substitution (see, for example, [38,78,79]): where the products of operators inside expression Q(p − A(i∂ p )) are symmetrized. We relate operatorQ = Q(p − A(i∂ p )) and its inverseĜ =Q −1 defined in Hilbert space H of functions (on M) with their matrix elements Q(p, q) and G(p, q) correspondingly, in the usual way Q(p, q) = p|Q|q , G(p, q) = p|Q −1 |q .
Here the basis elements of H are normalized as in (9). These operators obey the following equation: In the non-uniform case (either due to the external gauge potential A(x), or because of any other reason) Eq. (40) can be rewritten as while the Green function of Bloch electron is given by Here indices a, b enumerate the components of the fermionic fields. In the following we will omit those indices for brevity.

B. Electric current in the Wigner-Weyl formalism
Let us suppose, that we modified the external gauge field as A → A + δA. The original external field A can vary arbitrarily, but δA is assumed to be slowly varying, i.e. its variation at the distance of the lattice spacing can be neglected. In the linear response theory, the functional derivative of partition function with respect to this extra contribution to the gauge potential gives the electric current. In the tight-binding model we can write, Here the sum is over all lattice links [xy], while j [xy] is the current along the link [xy]. Since the latter is constant along a link and δA can be considered constant as well, we can rewrite the above expression as The former sum is over the physical points O, while the latter is over the extended lattice O, where there are 2 D more points than in O.
For the variation of the partition function we have The trace over the fermionic indices (and any inner symmetries) is implied here. The definition and the properties of the Weyl symbols allow us to rewrite the last expression as a trace of the Weyl symbols of corresponding operators Now we use that δA(x) (unlike A(x)) varies slowly at the distances of the order of the lattice spacing. This permits us to represent Substituting this to (45), recalling the definition of Tr of (19), and comparing with (43), we come to the expression for the electric current Local current (47) is not a topological invariant. Let us calculate the total current It is worth mentioning, that this definition of the total current differs somehow from the conventional definition (which is the integral over the surface of the current density across the given surface). The quantity called in the present paper for brevity the total current according to Eq. (48) is the sum over the whole lattice (including the imaginary time direction) of the current density. The averaging in time is of no effect, since we assume the system to be time independent. Unlike (47),J k is a topological invariant. Indeed, under small variations of the Weyl symbol of the lattice Dirac operator, Q W ≈ Q W + δQ W , the Green's function varies accordingly, G W ≈ G W + δG W , and then Eq. 16 guarantees that the symbol of theQĜ =1 becomes the simplest Groenewold equation G W ⋆ Q W = 1. Then, δG W = −G W ⋆ δQ W ⋆ G W , the two terms of (49) become where we integrated by parts and used that ∂ p l G W = −G W ⋆ ∂ p l Q W ⋆ G W . Now, simple cyclic transformation inside the trace shows that providing for the proof of Theorem 3. A note on the allowed variations of the Dirac operator are in order. When varying the current we are only allowed to introduce such modifications of the Dirac operator (and its Green's function) that do not break the mere existence of the integration over momenta and coordinates in (48). Physically it means, that any variations of the system that give birth to zero (delocalized, propagating) modes are forbidden. These modes would naturally cause a pole inĜ and render the integration divergent. Such modifications of the system correspond to the topological phase transitions, where the invariant does indeed change (or may change, in principle) its value.

C. Calculation of Hall conductance
Let us start from Eq. (47) for the electric current. We represent the electromagnetic potential as a sum of two contributions: where A (E) is responsible for the electric field, while A (M) -for magnetic one. The former is assumed to be weak, and we expand (47) up to the terms linear in A (E) and its derivatives.
Following (46) the symbol of Dirac operator acquires the form The Groenewold equation for G W then can be solved iteratively. We will keep in this solution the terms linear in A (E) and in its first derivative. The zeroth order term (that does not contain Further expanding the stars in the above expression, which contains the derivatives in x acting on A (E) , we have where Upon substitution of (51) and (53) into (47) we obtain where the last product is the ordinary one. Notice that Q Assuming that the external field is constant across the system, F (E) lm = const (equivalent to neglecting higher derivatives of A (E) , which is already done), one can calculate the total current averaged over the volume |V| (and also averaged in time) Here β = 1/T is the inverse temperature, T is assumed to be small. In the above expression we restored the ⋆-product in the last factor using once again (15). In what follows we will omit the superscript (0) for brevity. The averaged Hall conductivity is given now bȳ while the other transport coefficients are given by similar expressions W i[mk] , i = 1, 2, m, k = 1, 2, 3. All of them are invariant under such variations ofQ that do not break the applicability of the W, see discussion in the previous section. Indeed, underQ →Q + δQ, we can write which after the integration by parts is identically zero, finalizing the proof of Theorem 4. For the two-dimensional system (D = 2 + 1) we come to the following representation of the average Hall current (i.e. the Hall current integrated over the whole area of the sample divided by this area |A|) in the presence of electric field along the x 2 axis (we take into account that the corresponding Euclidean field strength has the components F Here x = (R 1 , R 2 , τ ), |A| is the area of the system, T is temperature (that is assumed to be small). It is implied that the area |A| is much larger than the area of the elementary crystal lattice cell, so that we still can deal with continuous values of momenta. For the same reason the sum over the Matsubara frequencies ω = 2π(N + 1/2)T (N ∈ Z) can be substituted by an integral at T → 0. Recall, that from the very beginning we discretized imaginary time, so that the whole lattice is composed of the (D − 1)-dimensional crystal lattice and the discretized imaginary time (as it occurs in the lattice regularization of the field-theoretical models). In condensed matter physics time typically remains continuous, the Matsubara frequencies belong to the interval (−∞, +∞) (at T → 0). In this case Eq. (61) reads: Here x = (τ, x), x is the point in space, τ is imaginary time that varies between 0 and 1/T → ∞, V (2) is the area of the two -dimensional lattice cell. The consideration of the three-dimensional systems with D = 3 + 1 is completely similar. It gives the current density integrated over the whole volume divided by this volume and averaged with respect to time) in the presence of external electric field E i : Here |V| is the overall volume, |V (3) | is the volume of the three-dimensional lattice cell. Notice, that in the above expressions we deal with an equilibrium system. Therefore, function G W (p, x) entering the above expressions does not depend on imaginary time τ .

A. Lattice Wigner-Weyl calculus through the series expansion
The previous definitions of the Wigner-Weyl formalism (used earlier for the description of the quantum Hall effect), [78,79], given for completeness in Appendix B and Appendix C, are sufficient for the description of the systems in the presence of slowly varying fields. This assumes, in particular the requirement for the real crystals that the magnetic field is much smaller than about 10000 Tesla. The definition of Appendix A, originating from [64], in turn, although is valid for any magnitudes of external fields, does not satisfy the requirements needed to express the Hall conductivity through the topological invariant composed of the Weyl symbol ofQ.
In order to calculate Hall conductivity using the Wigner-Weyl formalism for the systems of general type (in particular, for those that are subject to strong magnetic fields) we need the precise version of the formalism, which satisfies the above given Definition 1 unlike the Buot's version (see also Eq. A23). We will modify the versions of Wigher-Weyl calculus of Appendices A-C. We will see, that the conditions of Definition 1 (and thus the precise Groenewold equation) hold true in this case, while the Weyl symbols of simple operators remain non-degenerate.

Formal Definition
Following Appendix C we define implicitly the symbol Q W of an operatorQ, but unlike the previous works [78,79] we apply this definition both to the Dirac operators and to their Green's functions.
Namely, if an operatorQ is given as an operator-valued function of operators p, i∂ p : its symbol (which will be shown to be a Weyl symbol in the sense of Definition 1) is given implicitly by relation which must hold true for arbitrary functions f (P, K) and h(P, K) defined on momentum space, P, K ∈ M. The derivatives ∂ P and ∂ K inside the arguments of Q W act only outside of this function, i.e. ∂ K acts on f (P, K) while ∂ P acts on h(P, K). At the same time the derivatives without arrows act as usual operators, i.e. not only right to the functionQ, but inside it as well. Recall that the argument p of function Q(x, p) belongs to M, while for Q (i∂ P + i∂ K , P/2 + K/2) the arguments P, K belong to M. This is the reason why the integrals in Eq. (64) are over M.
In the following we restrict ourselves to D-dimensional rectangular lattices only, bet generalization to triclinic ones are immediate. For the practical calculation of Q W we should representQ using the Taylor epxansion for the function Q(x, p) in vicinity of x = 0 Here n = (n 1 , ..., n D ) is a multi-index, while (i∂ p ) n = i (i∂ pi ) ni . The representation of Eq. For the further use, we might also write the operator aŝ Here and below the sums are understood for every component of the multi-indices. It is worth mentioning that for the given operatorQ its representation in form of the series of Eq. (66) is not unique. There may exist different sets of functions {Q n (p)|n = 0, 1, ...}) that correspond to the same operatorQ defined on the Hilbert space for the given lattice O. (That would not be so if we consider operators in continuous theories.) In order to use (64) we represent The original operator-valued function Q is periodic in p, the same can be assumed of Q n (p) and q lm . Substituting the above expressions into the RHS of (64) we are able to perform there the integration by parts in It is solved by x m q ml (p + q)y l , here m, n are the multi-indices, x n = i (x i ) ni . Alternatively the solution can be written as The later equality is true since Therefore, where The given set of functions {Q n (p)|n = 0, 1, ...} defines uniquely the operatorQ. However, as it was mentioned above, the given operatorQ does not define uniquely the set of functions {Q n (p)|n = 0, 1, ...}. Correspondingly, the above expressions Eq. (70) and (71) define a set of different symbols of an operator. Further we will formulate the restrictions on this set to be used for the construction of an appropriate unique Weyl symbol.

Integral representation for QW
Our next purpose is derivation of an integral representation for Q W . First of all, for x ∈ O we can represent it as follows: Here we defined the D-dimensional periodic delta function as and used the following relations Recall, that p + k|p − k = δ [π/ℓ] (2k), see (9). Notice that where the sum is over all sets U , U ⊂ {1, 2, . . . D}. We plug it into the integrand of (72) and denoting, for brevity, p ± = p ± k, we represent each term of the sum over U using (66) as: Here the set of all indices is represented as a disjoint union of two subsets, {1, 2, ..., D} = U ⊔U ′ . The above expression allows us to obtain: i∂ p+,j nj and by n U we denote vector with n i = ℓ i for i ∈ U and n i = 0 for i ∈ U ′ . We come to the final result: where the sum is, once again, over all possible subsets of {1, ..., D}, while Q B is defined in (26) as It is discussed in Appendix A.
In the similar way we can represent Q W as follows One can see, that the given operatorQ defines uniquely the values of Q W (x, p) for x ∈ O. However, above it was argued that the given operator does not define uniquely the function Q W (x, p) for any real-valued x. One can see, that this ambiguity is related to the ambiguity of the analytical continuation of Q W (x, p) from its values at x ∈ O to the real-valued x.

Star identity I: AW ⋆ BB = (AB)B
Let us now consider the Moyal product of two different symbols of operatorsÂ andB, namely We substitute to this expression the definition of B B of Eq. (26) and obtain: Using (64) Thus we come to identity valid for In a similar way the following identity can be proved as well: Namely, let us calculate We substitute to this expression the above definition of B B , Eq. (26), and obtain: where A(p, q) = p|Â|q . We denote now A T (p, q) = A(q, p) = p|Â T |q , this gives The integration here is over M and the integration by parts can be performed (unlike Eq. (C10)). Without loss of generality we assume thatB(i∂ p , p) is written in the symmetric form: in each term the product of operators ∂ p and p is symmetrized. The corresponding representation can be obtained starting from Eqs. (66), (70) applying the commutation relation between operators ∂ p and p. ThenB T (i∂ p , p) =B(−i∂ p , p). We come to r = 1 2 D M dP e iP xB (i∂ P − i∂ p /2, p − P/2)A T (p − P/2, p + P/2)

Star identity II: AW ⋆ BW = (AB)W
The above definition (68) of the symbol Q W of an operatorQ through the series in powers of x works not only for the lattice models, but for the continuous theories as well (it actually originated there). Then in Eq. (64) instead of M one integrates over infinite momentum space R D . Correspondingly, one again represents Q W as Now Q n (p) depends on p ∈ R D , and is not periodic. The derivation of Sect. IV A 2 can be easily adopted to a continuum theory to give for any values of x the following representation of the symbol of an operator: Notice, that unlike the lattice case in the continuum theory the representation through the series defines the Weyl symbol uniquely. The above representation allows to prove the following identity (in the way similar to the corresponding proof given in Appendix A for the Buot's symbol) This identity in continuum theory is valid for any continuous values of x. Looking at Eq. (92) one recognizes that when the symbol of operators is written as a series, Eq. 92 has purely algebraic nature, i.e. it contains the series in the derivatives (entering the star operation) and the series in powers of x for both A W and B W . Therefore, Eq. (92) should be valid for both continuous and lattice systems.
Let us now give a direct proof of this identity in one-dimensional case. We have from (78) At the same time in Appendix A, Eq. (A6) it is shown that A B ⋆B B = (AB) B . Therefore, we conclude, that A −1,B (x+ ℓ, p)⋆ B B (x, p) = 0 at x ∈ O for arbitrary operatorsÂ andB. In a similar way we obtain A B (x, p)⋆ B −1,B (x+ ℓ, p) = 0. Then, The condition, which we encountered this way, A −1,B (x + ℓ, p) ⋆ B B (x, p) = 0, is deeply non-trivial. We note that a more general form must hold true

Trace definition
Let us study the trace of an operator, given by its series (66) trQ ≡ M dp p|Q|p = n M dp p|(i∂ p ) n Q n (p)|p = x∈O n M dp p|(i∂ p ) n |x x|p Q n (p) We consider here one-dimensional case for simplicity. The trace over indices corresponding to any inner symmetries is implied, as before. The cross check in coordinate space gives, using (93) and (B2) 1 |M| x∈O M dp Q W (x, p) = 1 |M| x,z1,z2∈O M dp z 1 |Q|z 2 δ 2x,z1+z2 e ip(z1−z2) + 1 |M| x,z1,z2∈O M dp z 1 |Q −1 |z 2 δ 2(x+ℓ),z1+z2 e ip(z1−z2) = x,z1,z2∈O z 1 |Q|z 2 δ 2x,z1+z2 e ip(z1−z2) δ z1,z2 + x,z1,z2∈O We might define the trace operation as However, it does not satisfy (15), as we will show in the next subsection. An alternative definition, which does respect (15), is where we used (A13) and (A17). Thus, for the class of operators satisfyingQ the Weyl trace is given by Again, the consideration in multi-dimensional case is similar, and for arbitrary D we obtain for operators satisfyingQ Note that 2 D |M| = |M|.

Trace of a product
Finally, we will check Eq. (15), i.e. if the star between two symbols of operators can be removed from an expression standing inside the trace. In the one-dimensional case (for simplicity) we can use (93) for each operator's symbol of a product to write where we used one of the remarkable properties of the B-symbol, namely (A20), to get rid of the mixed terms of the type of A B (x, p)B B (x + ℓ, p). Applying next (A14) for each of the terms above, along with (A15) for the second one, we obtain On the other hand, using (102) forQ =ÂB (and thus Q W = A W ⋆ B B due to Eq. 95) we immediately establish We come to the identity Its validity in the multi-dimensional case can be proved in similar way. At the same time ifÂ −U =Â andB −U =B for any U ⊂ {1, 2, . . . D}, then This compleets the proof of Theorem 2.
B. Wigner-Weyl formalism on a doubled lattice. One-dimensional case.

Operators on auxiliary lattice O and their Weyl symbols
In Eq. (7) we introduced the extended auxiliary lattice O that contains our physical lattice O as a subset, To build a lattice Wigner transformation we now enlarge operators acting on O to act on O. To this end we first introduce additional position and momenta eigenstates, and continue physical operators to the auxiliary lattice by the following relations, Moreover, we demand that the inter-lattice matrix elements vanish The extended Fourier decomposition reads Then, the matrix elements in momentum space of such operators become where by definition Here |M| = |M|/2 = π/ℓ (for the lattice in 1D case). This definition, in particular, implies that Q(p, q) is periodic with the period π/ℓ in each of the arguments, since O is a lattice with the link length being equal to 2ℓ, but the original matrix element, p|Q |q , is still periodic with the period 2π/ℓ due to the additional exponential factor in (116). Here by tr we denoted the ordinary, physical, trace of an operator.

W -symbol and its properties
We now use the B-symbol of operators, but defined on the extended lattice O. Using Eq. (116) we can express it as follows Formally speaking, Q W (x, p) is defined for any x ∈ R. However, for the discrete values, x ∈ O we get a simpler expression. Indeed, due to the structure of the argument of the exponential factor in (116), the integrand of (120) is periodic with the period π/ℓ, and thus the integration can be reduced from M to M, This is, therefore, the new definition of the Weyl symbol of an operatorQ defined on the lattice O. The second sublattice O ′ is an auxiliary instrument that may actually be omitted in the following as well as the extended lattice O.
Further, using Eq. (79) in (121), we obtain One can see, that Q W = Q W for the operators of the class (106) considered in the previous Section, withQ 1 =Q. The inverse transformation reads where by definition δ [a] (q) = N ∈Z δ(q − aN ). It gives Q(P 1 , P 2 ) = 1 It is worth mentioning, that the definition of the Weyl symbol of operatorQ given by Eq. (121) obeys the following important property:

Moyal product
We know that the B-symbol does map the product of operators into the star product of their symbols. Thus, the same property is naturally expected from Eq. (121) based on the definition of the B-symbol, where we denoted f (q) = 1 + e −2iqℓ for brevity. Both f (q) and A, B are periodic with the period π/ℓ, while M = (−π/2ℓ, π/2ℓ]. It is instructive to prove the above equation without referring to the extended lattice O, i.e. working within O only. The RHS of Eq. (125) reads Transforming the integration area to a rhombus, see Fig. 1, and changing the variables we have We used that the Jacobian of the transformation is |J| = 1/2. Notice that we cannot reduce back the integration area from M = (−π/ℓ, π/ℓ] to M = (−π/(2ℓ), π/(2ℓ)] because f ( P±R 2 ) is periodic with the period 2π/ℓ as a function of P or R, while A, B are periodic with the period π/ℓ. We introduce now then Noting that g(p ± π/(2ℓ)) = −g(p), we see that the first and the last terms in the square parenthesis do not change their signs when either P or R is shifted by π/ℓ. Thus, these terms give rise to 4 times the integral over M = (−π/(2ℓ), π/(2ℓ)].
The behavior of the remaining terms in the square parenthesis is more complicated. The shift of squares of type 1 on Fig. 2 by π/ℓ in the corresponding argument (i.e. to the position denoted by 1 ′ ) changes the sign of the g-factor, while shifting of the squares of type 2 into position 2 ′ does not. Consequently, the net result of the integration of these two terms gives zero -eight positives squares (4 diagonal and 4 inner ones), and eight negatives ones (the outer off diagonal squares). Thus, we come to In 1D we trivially have since g(q) = e −2iqℓ . To finalize the procedure we study This finishes the proof of (125).

Trace and its properties
As before, we can introduce two possible trace operations. One is using the summation over the points of lattice O only, = M dp Q(p, p) = trQ.

Groenewold equation
Let us check now that introduction of the auxiliary lattice O does not spoil anything. We consider two operatorsQ,Ĝ that obey Eqs. (113) and (114), which are inverse to each other: The Fourier representation gives whose LHS using (117) can be rewritten as To check the consistency, we calculate in the last line we expanded the summation in x 2 to O = O ∪ O ′ using (113). And now M dq Q(p, q)G(q, k) = 2 π/ℓ x1,y2∈O x 1 |QĜ |y 2 e i(y2k−x1p) = 2 π/ℓ x1,y2∈O δ x1y2 e i(y2k−x1p) Thus we have We recall now that M = (−π/ℓ, π/ℓ] while Q(p, q) is periodic with the period π/ℓ, and the exponential factor changes its sign under the shift by π/ℓ. Therefore, the last term in the above expression vanishes. We arrive at which is actually a representation of δ [2π/ℓ] (k − p) as in (141), while Thus the Weyl transform given by (121) does result in C. Wigner-Weyl formalism on an extended lattice. The multidimensional case.

The auxiliary extended lattice
Let us extend our consideration to the multidimensional case. We are going to obtain an analogue of Eq. (121). Rectangular lattice and its first Brillouin zone are It can also be written as 2 M copies of O: here U is any subset of {1, 2, . . . M } including the empty one, the number of such subsets is 2 M . The empty subset corresponds to O, evidently. The relation between the sublattices is now where and we assume that ℓ ∅ = 0. Along the lines of the 1D-case, we consider Q-class of operators, which satisfies the following condition which guarantees that the matrix elements are non-zero only for the transitions within the same sublattice. The second condition becomes First, we note that in the D-dimensional case all states |p are periodic with the periods g (j) . The Fourier transform will read Using (154) we obtain Here we introduced similarly to (117) This object is periodic with the periods g (j) /2. We also represent

Moyal product
The multidimensional Weyl symbol shall be defined as It satisfies all the properties we require for the proper Weyl symbol in Sect. II B Let us prove, as an example of the D-dimensional considerations, that In other words, we shall prove that (162) The RHS of this expression reads (163) Then we introduce notation and transform the integration area to the form of a multi-dimensional rhombus As in the 1D case we write Then For each term of g( P±R 2 ) there exists such j that the shift P → P ± g (j) /2 (or similar with R) changes its sign. Thus such a term will be cancelled when integrated in the j-th direction of M (see Fig. 2).
The product g( P+R 2 )g( P−R 2 ) requires an additional analysis. Indeed, it will contain the exponents of the type now shifting P by g (j) will affect it by the following factor and similarly for the shifts of R. So, under D shifts we transform M to the smaller M, and only the terms with survive. Thus The R integration goes exactly as in (134), and we arrive at where each Weyl symbol is understood in terms of the operators acting in O only, as in the RHS of (160).
In the complete analogy with the one-dimensional case we can also demonstrate all other relevant properties. In particular, for the operators that obeyÂ −U =Â we have A W (x, p) = A W (x, p).

D. A brief word on motivation
To develop the double/extended lattice approach we used as a motivation the following consideration. Let us consider the following B-symbol in the extended lattice and apply it to an exponential operatorQ where b and l are certain constant vectors,x ≡ −i∂ p . Then using the Hausdorff formula we obtain [pb,xl] = ibl ⇒ e i(pb+xl) = e ipb e ixl e −ibl/2 (175) the B-symbol becomes (since e ixl = e +l∂p ) here δ g (·) is a delta function modulo the reciprocal lattice vectors. To be more specific where det {ℓ (j) } = det (...ℓ (j) ...) is the determinant of matrix composed of vectors ℓ (j) . We took into account that g (j) ℓ (j) = 2π. Given that q ∈ M, there are the following non-vanishing contributions of where g U = j∈U (−1) αj g (j) , α j = 0, 1 and U is any subset of {1, 2, . . . M }, i.e., for some particular choice of signs α j depending on l, the combination l/2 + g U /2 will also be inside M for any choice of U . Thus we have We took into account that for any lattice vector x ∈ O we have gx = 2πk, k ∈ Z. Now to have Q W (x, p) = e i(pb+lx) it is sufficient to require that b is an even combination of the lattice basis vectors Recalling that e ipb is simply the translation operator in the coordinate representation, we come to the requirement that our class of operators must describe only the even jumps. Hence, we separate the original lattice into the even and the odd ones, O = O ∪ O ′ , and associate the even one, O, with the physical crystal lattice.

V. DISCUSSION
It is clear that in all real experiments on the Quantum Hall Effect the actual magnetic field is non-homogeneous (though, possibly, only slightly). Still, the quantization of the Hall conductivity follows the same rigorous rules as derived for the magnetic fields strictly constant over the sample. Up to now this apparent contradiction was not studied, and we presented here the first rigorous demonstration that in the non-homogeneous systems the Hall conductivity is given by a topological invariant.
To this end, we proposed the precise Wigner-Weyl calculus for lattice models of condensed matter physics or quantum field theory. We restricted ourselves to consideration of the models on rectangular lattices, but generalization is immediate to general triclinic ones. All operators in such models are defined on Hilbert space described in Sect. II A. The constructed Weyl symbols of these operators obey the properties listed in Sect. II B. Our Weyl symbols differ from those proposed in the other papers, and allow us to derive the topological expression for the Hall conductivity presented in Sect. II C. This expression remains valid for arbitrary magnetic field unlike our previous considerations [35][36][37], which were limited to magnetic fields much smaller than several thousands Tesla and of wavelengths much larger than 1 Angstrom. Thus, our present construction can be also used for consideration of artificial crystal lattices [84][85][86], in which external magnetic field of order of several Tesla may be sufficient to observe the Hofstadter butterfly. The corresponding non-homogeneous systems can be investigated using the formalism proposed here.
The extension of the presented formalism may be relevant for the description of other non-dissipative transport phenomena (Chiral Separation Effect, Chiral Torsional Effect, Spin Hall Effect, etc), which were previously considered for homogeneous magnetic fields. With the aid of our precise Wigner-Weyl calculus the corresponding transport coefficients can be written in closed form, and their topological properties investigated.

It gives
Notice that the summation here goes over the more dense lattice O. Therefore, to restore the operator by its symbol we need to know the values of the latter not only on the physical points of O, but on the intermediate ones as well, O ′ . It appears, that with this definition of the symbols the Moyal product expression is precise for x ∈ O: The proof is as follows The factor 2 D in the third line results from the Jacobian. In addition, we took into account that the integration in the third line is over q, k ∈ M, that is the integration over P = q + k, K = q − k is over the region of rhomboid form, that is contained inside M. However, due to the periodicity the integral over this figure is equal to 2 D times the integral over P, K ∈ M. This results in the factor 1 2 D in the third line. It is similarly easy to check, that all of the properties (13)-(16) of Definition 1 for this version of the Wigner-Weyl technique are satisfied, except for the last one, see following subsections.
Certain properties of the symbol are exceptional. In particular, for any two operatorsÂ andB. Let us check (A7) directly, using the x-representation (B2) of the B-symbol valid for x ∈ O. The sum does not vanish if both Kronecker symbols are non zero, i.e. z 1 − z 2 +z 1 +z 2 = 2ℓ −z 1 +z 2 + z 1 + z 2 . This leads to equation which is never satisfied for z 2 ,z 1 ∈ O, leading indeed to (A7).

Star product without differentiation
Let us represent the star product of A W (x, p) and B W (x, p) for x ∈ O through the matrix elements ofÂ andB Next, using Eq. (B2) we obtain (for x ∈ O) One can see, that in order to define the star product of the Weyl symbols A B (x, p) and B B (x, p) for x ∈ O we do not need to know the values of these functions for all real values of x. It is enough to know the values of the Weyl symbols for x ∈ O, which is the extension of the physical lattice (it is doubled in one-dimensional models, and enlarged 2 D times in general case).

Trace and its properties
There are two tentative definitions of the trace distinct in the summation set: O in the former, and O in the latter. Both satisfy (14), but only "Tr O " also cater for (15) Tr For the future notice, it is worth mentioning that since the shift of the argument of the symbol only results in the change of the initial point of the sum in x ∈ O. We obtain immediately simply because A B (x + ℓ, p) ⋆ B B (x + ℓ, p) = (ÂB) B (x + ℓ, p). On the other hand, And a similar property holds for a ⋆-product, From the basic trace and ⋆ properties, (A5) and (A13), it can be trivially shown that The trace without a star, on the other hand, gives (with an arbitrary shift κ proportional to the distance between the adjacent points of the extended lattice O) For the sum to be non-zero, it must hold, in particular, that which has solutions in O only for κ given by an even integer times ℓ (in every component, if we are dealing with a D-dimensional problem). For κ = 0 we get then which together with (A19) proves that Tr O of (A12) does satisfy (A14)

B-symbol of unity and the other examples
On the other hand, it is a straightforward calculation to show that (16) is not satisfied for B-symbol. In 1D case (for simplicity) it reads, dq e 2ixq δ(q) + 1 2 δ(q − π/2ℓ) + 1 2 δ(q + π/2ℓ) = 1 2 (1 + cos(πx/ℓ)). (A22) The Groenewold equation in its turn (forQĜ =1 and x ∈ O) becomes Thus, the Groenewold equation acquires the form, which is not convenient for the derivation of the proper expression for the Hall conductivity.
In a similar fashion we can consider another example: an exponential operator, corresponding to a jump between the adjacent points in O,Q = e 2ipℓ , e 2ipℓ B (x, p) = π/2ℓ −π/2ℓ dq e 2iqx e 2i(p−q)ℓ δ [π/ℓ] (2q) dq e 2iqx e 2i(p−q)ℓ δ(q) + 1 2 δ(q − π/2ℓ) + 1 2 δ(q + π/2ℓ) = 1 2 e 2ipℓ + 1 2 e 2i(p−π/2+πx/ℓ) + 1 2 e 2i(p+π/2−πx/ℓ) = 1 2 e 2ipℓ (1 − cos(πx/ℓ)). (A24) One can see, that for x ∈ O the B-symbol of this operator vanishes at all. In a more general case of the homogeneous system withQ = f (p), where f (p) is a function defined in momentum space with periodic boundary conditions, we have Q B (p, x) = 1 2 (f (p) + f (p + π/2ℓ)cos(πx/ℓ)). This property of the B-symbol of an operator does not allow to obtain slowly varying function of x for Q B (x, p) in the case, when operatorQ(x,p) depends on the operatorx slowly. This does not allow to apply in practice the given version of Wigner-Weyl formalism to a wide class of lattice models. In particular, in the simplest tight-binding models the Dirac operator Q(p) is proportional to cos 2pℓ, and the corresponding B-symbol vanishes identically for x ∈ O.
The proof is given in [38]. We repeat it here briefly Here the bra-and ket-vectors in momentum space are defined modulo vectors of reciprocal lattice. In the second line we change variables with the Jacobian This results in the factor 2 D in the third line. Here D is the dimension of space-time.
The transition between the second and the third lines of Eq. (B4) is only approximate whenever the operators are not diagonal. If, however, the off-diagonal matrix elements are small, it becomes possible to substitute the region of the integration in q and k (that corresponds to P, K/2 ∈ M) by M ⊗ M.
This version of the Wigner-Weyl technique may be applied successfully to the lattice Dirac operator (the inverse Green function) and the Green function G(p, q) itself (to be considered as matrix elements of an operatorĜ: G(p + q/2, p − q/2) = p + q/2|Ĝ|p − q/2 ). Both are almost diagonal if the external electromagnetic field varies slowly, i.e. if its variation on the distance of the order of lattice spacing may be neglected. This occurs for the magnitudes of external magnetic field much smaller than thousands Tesla, and for the wavelengths much larger than 1 Angstrom. One has in this approximation G C (x, p) ⋆ Q C (x, p) = (ĜQ) C (x, p) = q∈M dqe ixq p + q/2|ĜQ|p − q/2 = q∈M dqe ixq δ(q) = 1 . (B5) Thus we have the Groenewold equation With this definition of Weyl symbols forQ = e 2ipℓ we have (in the one-dimensional model with the Brillouin zone (−π, π]): dqe iqx e i(p−q)ℓ δ [π/ℓ] (q) = e 2ipℓ In a more general case of the homogeneous system withQ = f (p), where f (p) is the function defined in momentum space with periodic boundary conditions, we have Q C (p, x) = f (p). An alternative version of Wigner-Weyl formalism has been proposed in [78,79], although still giving only approximate results for lattice models.
At the same time, for the operators defined by their functional symbols Q(x, p),Q = Q(i∂ p , p), an implicit definition was proposed based on the integral equation M dP dK f (P, K)Q W (−i ∂ K + i ∂ P , P/2 + K/2) h(P, K) = M dP dK f (P, K)Q i∂ P + i∂ K , P/2 + K/2 h(P, K), which has to hold for arbitrary functions f (P, K) and h(P, K) defined on momentum space P, K ∈ M. Notice, that ∂ ∂(P/2+K/2) = ∂ P + ∂ K and ∂ ∂(P/2−K/2) = ∂ P − ∂ K . The practical calculation of Q W is as follows. We should represent the integrand of the RHS in (C1) in the following way (for that we use the commutation relations to rearrange the order of operators): Q i∂ P + i∂ K , P/2 + K/2 = i,j (i∂ P ) i q ij (P/2 + K/2)(i∂ K ) j .
If the functional symbol was given by its Taylor series around zero, then q ab (p) = ∞ j=0 C a+j a+b+j C j a+j 2 −j (i∂ p ) j Q a+b+j (p) .
The operator-valued function Q is periodic in P and K. However, this does not necessarily apply to the functions q ij . Nevertheless, we have: Q i∂ P + i∂ K , P/2 + K/2 = i,j (i∂ P ) i q ij ([P/2 + K/2] mod g j )(i∂ K ) j (C4) where g j are vectors of the reciprocal lattice. Substituting this into the RHS of (C1) we are able to perform the integration by parts, obtaining M dP dK f (P, K)Q W (−i ∂ K + i ∂ P , P/2 + K/2) h(P, K) = M dP dK f (P, K) It is solved by alternatively written as The later equality is true since it holds that q a0 (p) = 1 2 a a b=0 q b,a−b (p). (C7) Let us consider the product assuming thatĜ andQ are inverse operators,ĜQ =QĜ =1. We substitute to this expression the definition of G C and obtain: r = M dP Q W (x + i ∂ p /2, p − i ∂ x /2)e iP x G(p + P/2, p − P/2).
In this expression the derivatives ∂ p and ∂ x inside the arguments of Q W act only outside of this function, i.e. on e iP x G(p + P/2, p − P/2) and do not act inside the function Q W , i.e. on p and x in its arguments. This gives r = M dP e iP x Q W (−i ∂ P + i ∂ p /2, p + P/2)G(p + P/2, p − P/2).
Now let us suppose, that the function G(p + P/2, p − P/2) is nonzero only in the small vicinity of P = 0 mod π/ℓ, i.e. G is almost diagonal operator. Then the integration by parts may be applied, although only approximately, which upon substituting Q W using (C6) gives r ≈ M dP e iP x   i,j (i∂ P ) i q ij (p + P/2)(i∂ p /2) j   G(p + P/2, p − P/2) = M dP e iP x Q(i∂ P + i∂ p /2, p + P/2)G(p + P/2, p − P/2) = M dP e iP x δ [π/ℓ] (P ) = 1 ∀x ∈ O, and we come to the Groenewold equation, although only approximate. Thus, in the same way, as for the case of Q C , the given definition allows to deal effectively with the systems, where external fields vary slowly. The advantage of the formalism discussed in the current section is that the Weyl symbol may be calculated relatively easily for a certain class of operatorsQ.