{'fraction_non_alphanumeric': 0.04554973821989529, 'fraction_numerical': 0.018417000307976595, 'mean_word_length': 4.325733967525013, 'pattern_counts': {'":': 0, '<': 8, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
New perspectives on the Ising model

The Ising model, in presence of an external magnetic field, is isomorphic to a model of localized interacting particles satisfying the Fermi statistics. By using this isomorphism, we construct a general solution of the Ising model which holds for any dimensionality of the system. The Hamiltonian of the model is solved in terms of a complete finite set of eigenoperators and eigenvalues. The Green's function and the correlation functions of the fermionic model are exactly known and are expressed in terms of a finite small number of parameters that have to be self-consistently determined. By using the equation of the motion method, we derive a set of equations which connect different spin correlation functions. The scheme that emerges is that it is possible to describe the Ising model from a unified point of view where all the properties are connected to a small number of local parameters, and where the critical behavior is controlled by the energy scales fixed by the eigenvalues of the Hamiltonian. By using algebra and symmetry considerations, we calculate the self-consistent parameters for the one-dimensional case. All the properties of the system are calculated and obviously agree with the exact results reported in the literature.


I. INTRODUCTION
It is really very hard to say something new on the Ising model. The model, originally proposed by Lenz 1 in 1920, was exactly solved for the case of an infinite chain by Ising 2 in 1925. Since then, thousand and thousand of articles and several books have been published on the subject. The reason is that the model is very simple, but still can be considered as the prototype for systems subject to second order phase transitions and can be effectively used for studying critical phenomena. Moreover, the model and its generalizations and modifications can also be used for studying a large variety of physical systems. We do not attempt to summarize the enormous work done in these 80 years; it would go well beyond the purpose of this article. An excellent historical presentation of the Ising model can be found in Ref. 3, although it is old and obviously not updated. With no pretension of being exhaustive and complete, we here summarize the principal approaches used in these 80 years.
A basic tool is the transfer matrix method 4,5,6,7 . By means of this approach, Onsager 8 in 1944 succeed to give an exact solution of the model for a cubic twodimensional lattice in absence of external magnetic field. The theory of spinors and Lie algebra was used to simplify the Onsager solution 9, 10 . Among the exact results for the two-dimensional case, the calculation of the magnetization 11 and the writing of the spin correlation function in the form of a Toeplitz determinant 12 have to be mentioned. Other simplifications of the Onsager solution have been obtained by means of the Jordan-Wigner transformation and fermionization methods 13,14,15 . Different approaches are based on combinatorial methods 16,17,18 and pfaffian methods 19,20,21,22 . More recent approaches have seen the Ising Hamiltonian expressed as a Gaussian Grasmannian action 23,24 . Along this line, use of operatorial symmetries that simplify the algebra of the transfer matrix has led to the calculation of the partition function for a large class of lattices 25,26 .
Many approximation methods have been used with the goal of obtaining an expression for the partition function valid over a large temperature range: mean field theory, Bethe approximation 27 , cluster variational methods 28 , Monte Carlo simulations, series expansions. The spin correlation functions have been studied at the critical temperature 29 and in the asymptotic region 30,31 . To study critical phenomena and critical indices, tools like series expansions 32,33,34 , scaling 35,36,37 , renormalization group theory 38,39 have been used.
In spite of the tremendous work done, many problems remain unsolved. The exact partition function in a finite magnetic field is still unknown for dimensions higher than one. Very few exact results have been obtained for the three-dimensional model. There is no exact solution for the two-layer Ising model either. Most of all, a general approach which works in all dimensions and under general boundary conditions, although in some approximation, is needed.
In a recent work 40 , we have shown that there is a large class of models which are exactly solvable in terms of a finite number of parameters that have to be selfconsistently calculated. The purpose of the present paper is to apply the method proposed in Ref. 40 to the Ising model and to show that an exact solution of the model does exist for any dimension. In Section 2, we introduce the Ising model for a d -dimensional cubic lattice and show that the model is isomorphic to a system of localized spinless interacting particles, satisfying the Fermi statistics. In Section 3, the Hamiltonian of the latter model is solved, that is, a complete finite set of eigenoperators and the relative eigenvalues are determined. Then, as shown in Section 4, the exact form of the retarded Green's function (GF) and of the correlation function (CF) can be obtained. In Section 5, we derive a set of equations for determining the charge/spin correlation functions. As the composite operators do not satisfy a canonical al-gebra, the GF, the CF and the charge/spin correlation functions depend on a set of internal parameters not calculable by the dynamics. For the one-dimensional case, by means of the composite operator method 41,42,43 , we calculate these internal parameters (Section 6) and the charge/spin correlation functions (Section 7). Although obvious, it is worth noticing that all the results reproduce the exact solution known in the literature.
What are the advantages of the present method and what is new in the context of the Ising model? We present a new scheme of calculations for treating the model. The scheme is general and can be applied to any dimension. In the framework of this scheme we show that the model is always solvable for all dimensions. The energy spectra of the system are known. In the one dimensional case we show that the energy scales determined by the spectra rule the behavior at the critical temperature. It is reasonable to expect that this is true also for higher dimensions. General relations among different spin correlation functions have been obtained. These are exact relations and might be used to check the consistency of some approximate treatments or numerical calculations. In order to get quantitative results for the cases of two and three dimensions we have to determine a finite small number of parameters. All the properties of the Ising model, the magnetization, the thermodynamical quantities, the spin correlation functions, depend on these parameters that have to be self-consistently determined. By using algebra and symmetry considerations we calculate these parameters for the case d = 1. Extension of the calculations to higher dimensions is under investigation.

II. THE ISING MODEL
The Ising model, in presence of an uniform external magnetic field h, is described by the following Hamiltonian S(i) are spin variables, residing on a d-dimensional Bravais lattice of N sites spanned by the vectors R i = i . The variables S(i) takes only two values: up or down, or more simply S(i) = ±1. For a hypercubic lattice of lattice constant a with nearest neighbor interactions, the exchange matrix J ij is given by where d is the dimensionality of the system and k runs over the vectors in the first Brillouin zone. The exchange constant J can be positive or negative, and accordingly the coupling will be ferromagnetic or antiferromagnetic. According to (2), the Hamiltonian (1) can be rewritten as It is worth to recall that the Ising Hamiltonian (1) is invariant under the transformation Let us consider a system of N e interacting spinless fermions residing on the same lattice and let c(i) and c † (i) be the related annihilation and creation operators. These operators are Heisenberg fields [i = (i, t)] satisfying canonical anticommutation relations As a consequence of the algebra (6), each site can be occupied at most by a single particle. The occupation number of the site i, ν(i) = c † (i)c(i), takes only the values 0 and 1. By taking into account two-body interactions, the Hamiltonian for such a system reads as where µ is the chemical potential and V (i, j) is the potential. This model Hamiltonian can be connected to the Ising model by defining It is clear that By substituting (8) into (7) and by considering only a nearest-neighbor potential we can rewrite the Hamiltonian (7) in the following form where we defined Hamiltonian (10) is just the Ising Hamiltonian (3) as we have the equivalence The relation between the partition functions is Then, the thermal average of any operator A assumes the same value in both models According to this, we can choose to study either one or the other model and get both solutions at once. We decide to put attention to the model Hamiltonian (7), which for a nearest-neighbor potential reads as where The spin-inversion symmetry (5) of the Ising model (3) corresponds to the particle-hole symmetry exhibited by the Hamiltonian (15). In particular, we have that the chemical potential as a function of ν = ν(i) scales as

III. COMPOSITE OPERATORS AND EQUATIONS OF MOTION
It is immediate to see that the charge density operator ν(i) satisfies the equation of motion Then, standard methods based on the use of equations of motion and Green's function (GF) formalism are not immediately applicable. Indeed, it is easy to check that the causal propagator T [ν(i)ν(j)] [T is the chronological operator] and the correlation function ν(i)ν(j) assume the form where Γ(k) is the zero frequency function 41 which cannot be calculated by means of the dynamics 53 .
Then, in order to solve the Hamiltonian (15) let us consider the composite operator This field satisfies the equation of motion By taking higher-order time derivatives we generate a hierarchy of composite operators. However, we observe that for p ≥ 1 the number operator ν(i) = c † (i)c(i) satisfies the following algebra Therefore, the hierarchy of composite operators (20) must close for a certain value of p and we should be able to derive a finite closed set of eigenoperators of the Hamiltonian. To this purpose, on the basis of (22) the following fundamental property of the field [ν α (i)] p can be established where the coefficients A m are given in Appendix A for the cases d = 1, 2, 3. We now define the composite operator After (23), this field is an eigenoperator of the Hamiltonian (15) i where the (2d + 1) × (2d + 1) matrix ε (d) , the energy matrix, is defined in Appendix B. It is easy to see that the eigenvalues E (d) n of the energy matrix are given by E (d) n = −µ + (n − 1)V n = 1, 2, · · · , 2d + 1 (27) The Hamiltonian (15) has been solved since we know a complete set of eigenoperators and eigenvalues, and we can proceed to the calculations of observable quantities. This will be done in the next Sections by using the Green's function formalism .

IV. RETARDED AND CORRELATION FUNCTIONS
We define now the thermal retarded Green's function where · · · denotes the quantum-statistical average over the grand canonical ensemble. By introducing the Fourier transform (29) and by means of the Heisenberg equation (26) we obtain the equation where I (d) (k) is the Fourier transform of the normalization matrix, defined as The solution of Eq. (30) is The spectral density matrices σ [ where Ω (d) is the (2d+1)×(2d+1) matrix whose columns are the eigenvectors of the matrix ε (d) . The explicit expressions of Ω (d) are given in Appendix B. The spectral density matrices σ (d,n) (k) satisfy the sum rule where M (d,p) (k) are the spectral moments defined as  (34). The explicit expressions of I (d) (k) and σ (d,n) (k) are given in Appendices C and D, respectively, for the cases d = 1, 2, 3. The correlation function can be immediately calculated from (32) by using the spectral theorem and one obtains with Equations (32) and (38) are an exact solution of the model Hamiltonian (15). One is able to obtain an exact solution as the composite operators ψ p (i) = c(i)[ν α (i)] p−1 constitute a closed set of eigenoperators of the Hamiltonian. However, as stressed in Ref. 41, the knowledge of the GF is not fully achieved yet. The algebra of the field ψ (d) (i) is not canonical: as a consequence, the normalization matrix I (d) (k) in the equation ( 30) contains some unknown static correlation functions, correlators (see Appendix C for explicit calculations), that have to be self-consistently calculated. According to the scheme of calculations proposed by the composite operator method 41,42,43 (COM), one way of calculating these unknown correlators is by specifying the representation where the GF are realized. The knowledge of the Hamiltonian and of the operatorial algebra is not sufficient to completely determine the GF. The GF refer to a specific representation (i.e., to a specific choice of the Hilbert space) and this information must be supplied to the equations of motion that alone are not sufficient to completely determine the GF. Usually, the use of composite operators leads to an enlargement of the Hilbert space by the inclusion of some unphysical states. Since the GF depend on the unknown correlators, it is clear that the value of these parameters and the representation are intimately related. The procedure is the following. We set up some requirements on the representation and determine the correlators in order that these conditions be satisfied. From the algebra it is possible to derive several relations among the operators. We will call algebra constraints (AC) all possible relations among the operators dictated by the algebra. This set of relations valid at microscopic level must be satisfied also at macroscopic level, when expectations values are considered. Use of these considerations leads to some self-consistent equations which will be used to fix the unknown correlators appearing in the normalization matrix. An immediate set of rules is given by the equation where the l.h.s. is fixed by the AC and the boundary conditions compatible with the phase under investigation, while in the r.h.s. the correlation function C(k, ω) is computed by means of equation of motion [cfr. Eq. (38)]. Another important set of AC can be derived by observing that there exist some operators, O, which project out of the Hamiltonian a reduced part Another relation is the requirement of time translational invariance which leads to the condition that the spectral moments, defined by Eq. (35 ), must satisfy the following relation It can be shown that if (43) is violated, then states with a negative norm appear in the Hilbert space. Of course the above rules are not exhaustive and more conditions might be needed. According to the calculations given in appendices C and D, the GF and the correlation functions depend on the following parameters: external parameters (µ, T, V ), internal parameters (C 1,2d ), and (κ (1) , κ (1) , · · · κ (2d) ), defined as The parameters C (d)α µ,ν are determined by means of their own definitions (44), where the r.h.s. is calculated by means of ( 37)- (38). This equation gives From the results given in the Appendices C and D, we see that the spectral density matrices have the form where the matrices Λ 0 and Λ 1 do not depend on momentum k. Putting (47) into (46) we obtain Calculations given in the Appendices C and D show that the matrices Λ (d,n) 1 are linear combinations of the matrix elements C (d)α 1,p . Then, Eq. (48) gives a system of homogeneous linear equations. The determinant of this system is only function of the external parameters µ, T, V . This function will vanish only if there is a particular relation among these parameters. Since these parameters are independent variables the only solution is that all the matrix elements must vanish The matrices Λ (d,n) 1 are zero and the correlation function C (d) (k, ω) does not depend on momentum, as we expected. In the coordinate space the CF takes the expression The correlation function depends on 2d internal parameters: κ (1) , · · · , κ (2d) . In order to determine these parameters, we use the Pauli condition (40) which gives the self-consistent equations (50). New correlation functions appear and the set of self-consistent equations (51) is not sufficient to determine all unknown parameters. One needs more conditions. In the case of one-dimensional systems these extra conditions can be obtained by using the property (42).

V. CHARGE CORRELATIONS FUNCTIONS
In Sections 3 and 4, we have solved the problem of the Ising model in terms of a set of local parameters, defined by (45) and (52). In this Section, we want to show how we can calculate non-local correlation functions. Let us define the causal Green's function (for simplicity in this Section we drop the superindex (d)) the retarded and advanced functions where ψ(i) is the composite field defined in (25) and we used the fact the field operator ν(j) does not depend on time. The Fourier transforms of these quantities read as where Q = C, R, A. By means of the equation of motion (26) we have The most general solution of Eq. (57) is while the function g Q(n) (i, l, j) must be determined. P denotes the principal value. By recalling the retarded and advanced nature of F R,A (i, l, j), it is immediate to see that Therefore The solution of (58) is where the matrices d ψψ † (n) (i, l, j) and d ψ † ψ(n) (i, l, j) have to be determined. From the definitions (53)-(55) we can derive the following exact relations (65) A relation between the two correlation functions D ψψ † (i, l, j) and D ψ † ψ (i, l, j) can be established by means of trace properties. Indeed, it is straightforward to derive a KMS-like relation By recalling the definitions (55), this last equations can be written as (65) can be expressed in terms of the correlation functions as Analogous expression holds for the commutator. By means of (68) and by recalling that [see Eqs. (37)- (38)] we find that equations (65) have the following form the solution of (70) and (71,) is: By putting (72) and (73) into (60) and (64) we have From the study of the fermionic sector we have where σ (n) are the spectral functions given in Appendix D. Then, at equal time (74) becomes The system (78) gives a system of linear equations for the quantities D ψψ † (i, l, j). Since the inhomogeneous terms in this system are proportional to δ il , it is clear that D ψψ † (i, l, j) ∝ δ il . Then, we will take i = l and we write From (78) we have the system of equations From its own definition (79) and by using the recurrence relation ( 23), the matrix D(i, j) has the following structure.
(i) One dimension D 1,6 D 1,7 · · · D 5,7 D 6,7 D 1,7 D 2,7 · · · D 6,7 D 7,7 (87) With the definitions Then, we only need to calculate the matrix elements D 1,p (i, j) (p = 1, 2, · · · 2d + 1). The matrix J(i, j) can be obtained from the normalization matrix I(i, j) = {ψ(i, t), ψ † (j, t)} , calculated in Appendix C, by means of the following substitution Then, the matrices τ (n) (i, j) have the same expressions of the spectral matrices σ (n) when the following substitution is made. It can be seen that for j = i and j = i α the system (80) is exactly equivalent to the system (51). Then, it is enough to consider the case j = i, i α . In this case, the system (80) becomes with T n given by (39). The system (91) gives a set of exact relations among the correlation functions. We might think to solve this system by induction method, since some of the first correlation functions can be expressed in terms of the basic parameters κ (p) and λ (p) . However, when we do this, we immediately see that the number of equations is not sufficient to determine all the correlation functions and we need more equations. Once again, this can be done for the one dimensional system, as we shall see in the next Sections.

VI. SELF-CONSISTENT EQUATIONS FOR ONE-DIMENSIONAL SYSTEMS
Until now the analysis has been carried on in complete generality for a cubic lattice of d dimensions. We now consider one-dimensional systems, and in particular we will study an infinite chain in the homogeneous phase. For simplicity of notation we shall drop the superindex (d). By means of ( C6) and (D1)-(D2) the set of equations (51) gives the linear system where, because of translational invariance, we put It is immediate to see that for µ = V , the solution of the first equation in (92) for T > 0 is This is in agreement with the particle-hole symmetry enjoyed by the model [see (17)]. Recalling (8) and (11), this situation corresponds to the zero magnetization of the Ising model in absence of external magnetic field. Coming back to general value of µ, it is clear that Eqs.
(92) are not sufficient to specify completely the 4 parameters ν, κ (2) , λ (1) , λ (2) and we need another equation. A fourth equation can be easily obtained by means of the algebra. We observe that This relation leads to where By means of the requirement (42) the correlation function C 1,k = c(i)c † (i)[ν α (i)] k−1 can be written as where and · · · 0 denotes the thermal average with respect to H 0 . Let us define the retarded GF By means of the equation of motion we have Recalling the relation between retarded and correlation functions, from (102) we obtain By putting this result into (98) we have we obtain from (104) the relations Now, we observe 45 that H 0 describes a system where the original lattice is divided in two disconnected sublattices (the chains to the left and to the right of the site i). Then, in H 0 -representation, the correlation function which relates sites belonging to different sublattices can be decoupled: for j and m belonging to different sublattices. By using this property, invariance of H 0 under axis reflection and (106) we can write By putting (109) into (107), we obtain the following selfconsistent equation among the correlation functions By means of (37) and (38) and the results given in Appendices C and D, Eq. (110) takes the expression This equation together with equations (92) gives a system of 4 self-consistent equations for the 4 parameters ν, κ (2) , λ (1) , λ (2) as functions of µ, T, V . By solving the set of linear equations (92) with respect to κ (2) , λ (1) , λ (2) as functions of ν, we have To calculate the parameter ν let us put (112) into (111) and solve with respect to ν. We have two roots. One solution corresponds to an unstable state with negative compressibility and must be disregarded. By picking up the right root, and by using the relation we find As shown in Appendix E, the solutions (112)-(116) exactly correspond to the well-known solution of the 1D Ising model, obtained by means of the transfer matrix method. We could manipulate the expression (116) and the ones for κ (2) , λ (1) and λ (2) , obtained by substituting (116) into (112)-(114), in order to reproduce the expressions of the Ising model, given in Appendix E. However, we prefer to maintain the present expressions as the following discussion will be more transparent. In Section 3, we have seen that in the present model, in the onedimensional case, there are three energy scales At zero temperature, the three functions T 1 , T 2 and T 3 are not analytical functions at the points µ = 0, µ = V, and µ = 2V , respectively, and we expect that the parameters ν, κ (2) , λ (1) , λ (2) exhibit some discontinuous behavior at these points. As shown in Fig. 1, in the limit T → 0 the particle density ν has a discontinuity at µ = V for the case of negative V (i.e. J > 0, ferromagnetic coupling) and two discontinuities at µ = 0 and µ = 2V for the case of positive V (i.e. J < 0, antiferromagnetic coupling). Here and in the following, we take |V | = 1: all energies are measured in units of |V |. In particular, the particle density increases by increasing µ from zero to one. At zero temperature, in the ferromagnetic case ν is zero for µ < − |V | and equal to one for µ > − |V |; in the antiferromagnetic case ν is zero for µ < 0, jumps to 1/2 and exhibits a plateau, centered at µ = V , in the region 0 < µ < 2V , jumps to the value 1 for µ > 2V . The parameter κ (2) has a behavior similar to ν. In Fig. 2 we give the parameter λ (1) as a function of µ. For the ferromagnetic case the behavior is similar to that of ν . Instead, in the antiferromagnetic case λ (1) , at T = 0, exhibits only one discontinuity point at µ = 2V , where jumps from zero to one. The parameter λ (2) has a behavior similar to λ (1) . The different behavior exhibited by the pairs (ν, κ (2) ) and (λ (1) , λ (2) ) for V < 0 is naturally due to the antiferromagnetic correlations, when we recall that the pair (λ (1) , λ (2) ) describes a correlation between two first neighboring sites, while κ (2) describes correlations between two second neighboring sites. Of course in the point of discontinuity the two limits T → 0 and µ → µ c are not interchangeable. As we shall see in the next Section, the 4 local parameters ν, κ (2) , λ (1) , λ (2) are really basic since all the properties of the model are described in terms of them. It is worthwhile to note that some simple relations can be established among the parameters The Ising model in one dimension can be described in terms of only two parameters: ν and λ (1) .

VII. CHARGE CORRELATION FUNCTIONS FOR ONE-DIMENSIONAL SYSTEMS
The system of equations (91) establishes some relations among the non-local charge correlation functions As already discussed, the number of equations is not sufficient to determine completely the charge correlation functions, and one needs more equations to close the system. In the one-dimensional case a fourth equation can be easily obtained by algebraic considerations.
Recalling that c † (i)e −βH = c † (i)e −βH0 we can easily derive the following result Now, for j > i + a (because of invariance under axis reflection we could choose j < i − a as well), using the property (108) 122) Therefore, from (120) where we used Eq. (106). Equations (123) and (124) give Putting (127) into (125) we get the fourth self-consistent equation It is straightforward to verify that (128) is identically satisfied for j = i, while for j = i α coincides with equation (119). By adding (128) to (91), we have the following system of equations for the non-local spin correlation functions K (p) (i, j) and Λ (p) (i, j) where we put and we used K (0) (i, j) = ν(j) = ν. To calculate for general |i − j| we proceed by induction.
Let us put i − j = ma and concentrate the attention on the spin CF Λ (0) (m) = ν(m)ν(0) . We start by observing that where the two parameters ν, λ (1) have been calculated in the fermionic sector. By taking m = 2 we can calculate from the system (129)-(132) that The results (135)-(136) and use of the relation (118) show that Λ (0) (m) for m = 0, 1, 2 can be cast in the form where the parameter p is defined as By using the expressions of the basic parameters given in Section 5, it is possible to check that |p| < 1. Then, we can introduce the Fourier transform and reexpress (137) as Then, the CF Λ (0) (m + 1) can be calculated as In Fig. 3 we give p as a function of µ for V = −1 and V = 1 at various temperatures. We see that for negative V (ferromagnetic case), p is positive and various between zero and 1. For positive V (antiferromagnetic case), p is negative and various between −1 and zero. In particular, for negative V , p tends to 1 at µ = V in the limit T → 0. Instead, for positive V , p tends to −1 at µ = V in the limit T → 0. This is seen in Fig. 4 where p is plotted versus T at µ = V = −1 and at µ = V = 1.
Let us now discuss the correlation functions. Λ (0) (m) is plotted against m for µ = V = −1 [ Fig. 5 (top)] and for µ = V = 1 [Fig. 5 (bottom)] at various temperatures. We see that when zero temperature is approached a longrange order of ferromagnetic and antiferromagnetic type is established, respectively. Also, we can see from ( 141)-(147) that for m → ∞ and T = 0 (i.e. not at the critical temperature) the spin correlation functions assume the ergodic value. At the critical temperature T = 0 we have breakdown of the ergodicity.

VIII. CONCLUSIONS
The Ising model in presence of an external magnetic field is isomorphic to a model of localized spinless interacting particles, satisfying the Fermi statistics. The latter model belongs to a class of models always solvable, as shown in Ref. 40. On this basis, we have constructed a general solution of the Ising model which holds for any dimensionality of the system. The Hamiltonian of the model has been solved in terms of a complete finite set of eigenoperators and eigenvalues. The Green's function and the correlation functions of the fermionic model are exactly known and are expressed in terms of a finite small number of parameters that have to be selfconsistently determined. By using the equation of the motion method, we have derived a set of equations which connect different spin correlation functions. The scheme that emerges is that it is possible to describe the Ising model from a unified point of view where all the properties are connected to a small number of local parameters, and where the critical behavior is controlled by the energy scales fixed by the eigenvalues of the Hamiltonian. The latter considerations have been proved from a quantitative point of view in the one-dimensional case, where the equations which determine the self-consistent parameters and the spin correlation function have been solved. For d = 1 all the properties of the system have been calculated and obviously agree with the exact results reported in the literature. Extension of the calculations to higher dimensions is under investigation.
After the paper was completed and submitted for publication, the author learned that approaches to the spin-1/2 Ising model based on a fermionization of the model have been previously reported in Refs. 46 and 47. The author wishes to thank the referee and Prof. L. De Cesare for putting these papers to his attention. In particular, in Ref. 46,Tyablikov and Fedyanin showed that the chain of equations for the double-time GF closes and the number of equations is determined only by the co-ordination number, independently by the dimensionality. This conclusion agrees with the results given in Sections 3 and 4. In order to close the set of equations for the fermionic correlation functions in the one-dimensional case, the authors of Ref. 46 assumed ergodicity and solved the system, obtaining the exact solution of the 1D Ising model for an infinite chain. It should be remarked that ergodicity breaks down for finite systems and at the critical points. Kalashnikov and Fradkin in Ref. 47 used the spectral density method 48 to derive a system of equations for the correlation functions; also in this case the approach is valid for any dimension. However, the number of equations is less than the number of correlation functions.

APPENDIX A: ALGEBRAIC RELATIONS
As mentioned in Section 3, the number density operator ν(i) = c † (i)c(i) satisfies the algebra From this algebra an important relation can be derived. for the operator where i m are the first neighbors of the site i. We shall discuss separately the cases of different dimensions.

One dimension
We start from the equation After subtracting the terms m = 0 and m = p, we can use the algebraic relation (A1) to obtain From (A4), by putting p = 2 we obtain By substituting (A6) into (A4) we have the recurrence rule We note that the coefficients A where the operators Z m are defined as n n = 4(−1 + 3 · 2 p−1 − 3 p + 4 p−1 ) (A14) By solving the system (A12) with respect to variables Z m , we can obtain the recurrence rule where the coefficients A (p) m are defined as We note that for all p

APPENDIX C: THE NORMALIZATION MATRIX
We recall the definition of the normalization matrix It is straightforward to see that use of the hermiticity condition (43 ) leads to the fact that we have to calculate only the matrix elements I 1,m (k) (m = 1, 2, · · · 2d + 1). The calculations of these is very easy when one observes the following anticommutating rule (C2) By taking the expectation value of (C2) we obtain in momentum space with the definitions where 1,1 I 1,n+1 (k) (m = 2, · · · 5) (C8)