Exact Expressions of Spin-Spin Correlation Functions of the Two-Dimensional Rectangular Ising Model on a Finite Lattice

We employ the spinor analysis method to evaluate exact expressions of spin-spin correlation functions of the two-dimensional rectangular Ising model on a finite lattice, special process enables us to actually carry out the calculation process. We first present some exact expressions of correlation functions of the model with periodic-periodic boundary conditions on a finite lattice. The corresponding forms in the thermodynamic limit are presented, which show the short-range order. Then, we present the exact expression of the correlation function of the two farthest pair of spins in a column of the model with periodic-free boundary conditions on a finite lattice. Again, the corresponding form in the thermodynamic limit is discussed, from which the long-range order clearly emerges as the temperature decreases.


Introduction
Since the exact solution of the partition function in the absence of a magnetic field of the two-dimensional rectangular Ising model with periodic-periodic boundary conditions is obtained in the thermodynamic limit [1] and in finite-size systems [2], many authors have contributed to the knowledge of various aspects of this model, such as different boundary conditions, the arrangement modes of the spin lattice, surfaces, or mathematical methods, etc. [3][4][5][6][7].
Besides the partition function of the model, the calculations of spin-spin correlation functions are an important subject in the research of the two-dimensional Ising model. Some expressions of correlation functions in the thermodynamic limit have been obtained [3][4][5]8,9], and the case in a finite lattice has been studied [10][11][12].
The determination of exact expressions of the partition function and spin-spin correlation functions of the model on a finite lattice is not only a theoretical subject; the results obtained can also be used in the research of finite-size scaling, finite-size corrections, and boundary effects [7,[13][14][15][16].
In this paper, we present some exact expressions of spin-spin correlation functions of the two-dimensional rectangular Ising model on a finite lattice by employing the spinor analysis method [2].
In Section 2, for the model with L rows and N columns and periodic-periodic boundary conditions (Onsager's lattice), we calculate some exact expressions of the correlation function σ 1 , 1 σ 1 , 1+Q and compare the corresponding forms in the thermodynamic limit obtained here with known results presented in Reference [9]. The investigation in Section 2 shows the main steps, key points, and problems of the approach used in this paper. Since the whole process are complex, here we outline the approach.
(1) Although any spin-spin correlation functions σ l , n σ l , n can be expressed by matrices, and the matrices belong to spin representatives [8], here we only consider σ l , 1 σ l , 1+Q , i.e., the correlation functions of pairwise spins in one column, of which exact expressions can be obtained by the spinor analysis method.

Some Results Concerning the Partition Function Z
From (1), we see that to obtain σ 1, 1 σ 1, 1+Q , we need knowledge about the partition function Z given by (2), from which we summarize some results presented in Reference [2] as follows.
Z can thus be obtained by finding the trace of a matrix, where Γ 2N+1 is the matrix U defined by (15.68) in Reference [17], and the matrices 1 can be diagonalized at the same time.
On the other hand, both matrices U (↑) and U (↓) are spin representatives. By ω (↑) and ω (↓) we denote the corresponding rotation matrices of U (↑) and U (↓) , respectively; both ω (↑) and ω (↓) can be diagonalized: where A + indicates a complex conjugate to a matrix A, In the above formulas, γ n and ϕ n are as the abbreviations for γ and, further, considering σ 1 , q σ 1 , 1+q = lim φ q →0 ∂e φ q σ 1 , q σ 1 , 1+q ∂φ q , according to (1), σ 1, 1 σ 1, 1+Q can be written in the form: Also, by a standard method [1,17] we can further write Y Q in the form: Y Q can thus be obtained by finding the trace of a matrix. Here we are not going to write out the explicit expressions of V (↑) and V (↓) , but only point out that: can be diagonalized at the same time.
When N = 2M, θ n have properties: In (15), the values of m and m are exactly the same as those in (7). We can first evaluate the eigenvalues of the rotation matrices ζHζ + , and then obtain Y Q in terms of the spinor analysis method. Finally, we obtain σ 1, 1 σ 1, 1+Q according to (8) and (9).
From (10) to (13), we see that H is a Hermitian conjugate matrix. Hence, the eigenvalues of both matrices ζHζ + and H are the same, and we therefore can only evaluate the eigenvalues of the rotation matrix H.

Basic Properties of the Eigenvalues and Eigenvectors of the Matrix H
Any rotation matrix A has the following property: If τ is an eigenvalue of A, then τ −1 is also an eigenvalue of A [2]. Since ζHζ + is a rotation matrix and the eigenvalues of ζHζ + and H are the same, H should have the same property. In this sub-section we prove this conclusion.
The eigen equation of H reads: where A T means the transpose of a matrix A, we denote the eigenvector of H, and, introducing C ∆ q , C ∇ q (q = 1, 2, · · · , Q) in terms of: we can prove that the eigen Equation (16) is equivalent to: where ) , for (↓), where we have used (7) and (15). We see that all A ± (k) and B(k) are real functions and satisfy: According to the above expressions and the properties of A ± (k) and B(k), we can conclude that if τ and C ∆ q (τ) , C ∇ q (τ) (q = 1, 2, · · · , Q) satisfy (17), then τ = τ −1 and also satisfy (17), where C 0 (τ) is an arbitrary function. From this discussion we not only prove the conclusion "If τ is an eigenvalue of the matrix H, then τ −1 is also an eigenvalue", but also obtain the relation (18) between C ∆ q (τ), C ∇ q (τ) and C ∆ q (τ −1 ) , C ∇ q (τ −1 ) . The conclusion and the relation (18) are useful to determine the forms of approximate eigenvalues and the expressions of the normalized eigenvectors of H, as well as to calculate the determinant of the matrix consisting of the eigenvectors in the actual calculation process.

Approximate Method for Solving the Eigen Equation (16)
It is very difficult to find the exact eigenvalues of H by solving the eigen Equation (16). On the other hand, the operator in (8) allows us to ignore all terms whose orders are higher than According to this key property, we can obtain the exact expressions of σ 1, 1 σ 1, 1+Q by only finding approximate eigenvalues of H.
Concretely, as the first step, the term 2 (10) is a diagonal matrix, whose eigenvalues and eigenvectors are summarized in the following formulas: , for the eigenvalue e Lγ n ; 0 0 , for the eigenvalue e −Lγ n , , for the eigenvalue e Lγ n ; 0 1 , for the eigenvalue e −Lγ n .
which are as the zeroth order approximation of the eigen Equation (16).
The term Q ∑ q=1 sinh2φ q H (q) with the factors sinh2φ q ≈ 2φ q (q = 1, 2, · · · , Q) in (10) can be regarded as a perturbation term. Then, by using RSPT, we can obtain the approximate eigenvalues of H.
However, although what eigenvalues we need are only corrected to the φ 1 q (= φ q ) order (q = 1, 2, · · · , Q), we must calculate the perturbation terms up to the Q-th order, not only for the first-order approximation, because all of the terms with the factor Q ∏ q=1 sinh2φ q appear in the Q-th order eigenvalues and are needed, which only include the φ 1 q order for every φ q . However, if we calculate the eigenvalues up to the Q-th order by using RSPT, then not only the actual calculation process is very complex, but there are also many unwanted terms with factors φ k q (k ≥ 2), for example, the term with the factor sinh 3 2φ 1 Q ∏ q=4 sinh2φ q , in the Q-th order eigenvalues.
To take out those terms with factors φ k q (k ≥ 2), we change "finding the eigenvalue up to the Q-th order" to "finding the eigenvalue through Q times first-order approximation".
We follow this approach up to sinh2φ Q H (Q) and every time we only calculate thr first-order approximation, which leads to the eigenvalues and eigenvectors τ (Q) , Ψ (Q) of all terms being only of the φ 1 q (= φ q ) order. All of the terms with Q ∏ q=1 sinh2φ q remain, and at the same time those unwanted terms with φ k q (k ≥ 2) do not appear. On the one hand, the above approximate method allows us to actually carry out the calculation process to find the eigenvalues and eigenvectors of H. In particular, once we obtain τ (1) , we can obtain Y 1 by the spinor analysis method, as well as obtain σ 1, 1 σ 1, 2 in terms of (8) and (9).
On the other hand, since RSPT is irregular, when Q is very large, e.g., Q ≈ N 2 , the above approach no longer functions. Hence, by this approach we can only obtain the exact expressions of correlation functions when Q is a small number, for example, σ 1, 1 σ 1, 2 , σ 1, 1 σ 1, 3 , σ 1, 1 σ 1, 4 , · · · , etc., which belong to the short-range order, but we cannot obtain the exact expressions of correlation functions when Q is larger, for example, , which belong to the long-range order.

Recurrence Formulas of the Eigenvalues and Eigenvectors τ (Q) , Ψ (Q)
According to the discussions in the above sub-section, we first regard sinh2φ 1 H (1) as a perturbation term, and, by using RSPT, evaluate eigenvalues and eigenvectors τ (1) , Ψ (1) of the matrix H (0) + sinh2φ 1 H (1) up to the first-order approximation. However, according to (7), all eigenvalues e ±Lγ (↑) n are doubly-degenerate; and, except e ±Lγ (↓) 2M and e ±Lγ (↓) M , all the remaining eigenvalues e ±Lγ (↓) n are also doubly-degenerate. Hence, for doubly-degenerate eigenvalues of H (0) , we must use the degenerate perturbation theory; the results obtained up to φ 1 1 order are as follows.
In principle, by following the above approach we obtain the eigenvalues τ (Q) . Furthermore, considering that up to the first-order approximation for φ q , we have 1 + Csinh2φ q ≈ e Csinh2φ q , the eigenvalues τ (Q) of H can be denoted by the forms: where the value of m is exactly the same as that in (7).
Based on the above forms of the eigenvalues τ (Q) and using the spinor analysis method, we obtain: Finally, according to (8) and (9), we obtain: where Z and Y Q are given by (6) and (23), respectively. Although in Section 2.5 we presented a simplified approximate method, the actual calculation process of σ 1, 1 σ 1, 1+Q is still complex; here, we only present the expressions of σ 1, 1 σ 1, 2 and σ 1, 1 σ 1, 3 directly.
When Q = 1, substituting α M , and δ 2M given by (20) into (23), we obtain Y 1 , and, further, we have: Substituting Z given by (6) and the above expression into (24), we obtain the expressions of σ 1, 1 σ 1, 2 of the model on a finite lattice: Then, using τ (1) , Ψ (1) presented in Table 1, (20) and (21), and according to (22), we obtain τ (2) , which can be denoted by the forms: In the above expressions, the values of m and m are exactly the same as those in (7). We see that all terms with the factor sinh2φ 1 sinh2φ 2 remain in the above expressions.

, tanh Lγ
Substituting the above two expressions into (26), we obtain: where the function θ(πx) in terms of (14) is defined by: The result (31) is in accordance with that in Reference [8].
According to the similar discussions, for σ 1, 1 σ 1, 3 we first have: (33) We discuss the first term in (28) as an example to show how to calculate lim m . First, using (7) and (15), the first term in (28) can be written in the form: According to (7), when k < m, γ k ) = ∞, and, thus, the first term in the above expression vanishes; when k > m, γ k ) = 0, we therefore obtain: Using this method to deal with the remaining terms in W (↑) m , we finally obtain: As M → ∞ , the second term in the above expression vanishes, and, according to the definition of the Riemann integral, we have: x dy(sin θ(πx) sin θ(πy) sin(πx) sin(πy) − cos θ(πx) cos θ(πy) cos(πx) cos(πy)) +2 x 0 dy cos(πx) cos(πy) , where the function θ(πx) is introduced by (32), Generally speaking, for the function f (u, v) and the domain D of the integration shown in Figure 1, we have: Substituting (31) and (36) into (33), we obtain the form of σ 1, 1 σ 1, 3 in the thermodynamic limit: On the other hand, the expressions of σ 1, 1 σ 1, 1+Q in the thermodynamic limit have been obtained [3,5,9]. Thus, we here cite the formulas (B6) and (B7) in Reference [9] for comparison. According to those two formulas: where θ(ω) is the function δ * (ω) in Reference [9]. We see that (31) and (37) obtained here are exactly the same as (38).

Long Range-Order in the Model with Periodic-Free Boundary Conditions
For the model with L rows and N columns and periodic boundary condition in the horizontal direction and free boundary condition in the vertical direction, we consider σ l , n σ l , n , i.e., correlation functions of pairwise spins in one column, periodic boundary condition in the horizontal direction leads to: where is the partition function of the system in absence of a magnetic field, where σ L+1 , m = σ 1 , m , J (> 0) and J(> 0) are the interaction constants for the horizontal and vertical directions, respectively.

Some Results Concerning the Partition Function Z 0
We summarize some results concerning Z 0 given by (40), some of which are obtained in Reference [12]. However, the approximate values of some quantities presented here show improvement over those given by Reference [12].
By using the spinor analysis method, Z 0 is obtained in Reference [12]: where γ n−1 (n = 1, 2, · · · , N) are determined by: where (n = 2, 3, · · · , N) (n = 1, 2, · · · , N) are N roots of the N-th order algebraic equation in x: where is an n-th degree polynomial in x. If by x ≡ d + d −1 2 we introduce the quantity d, then g n (x) can be written in the form: The expression in (45) is not only simple but also convenient for investigating the properties of g n (x), especially if we assume x = cos ϕ, then g n (x) = sin(n + 1)ϕ sin ϕ . Substituting these forms of g n (x) into (43), for the N − 1 roots of the N-th order algebraic Equation (43) we obtain: Further, θ n−1 can be determined by solving a transcendental equation about θ; if N is finite, then the evaluation of θ n−1 is complex because of the so-called "finite size effect"; for the limit case σ l, 1 σ l, N , we can assume θ n−1 = ∞ ∑ k=0 θ (k) n−1 N k and obtain θ n−1 by the iterative method. Further, we obtain γ 1 , γ 2 , · · · , γ N−1 in terms of (42). Concretely, correcting to 1 N order, we have: where γ (0) n−1 and θ (0) n−1 (n = 2, 3, · · · , N) are introduced by: More important are the values of x 0 and γ 0 ; to present x 0 and γ 0 , we first introduce a temperature T c in terms of: When T ≥ T c , 0 < x 0 ≤ 1; however, when T < T c , K < K and 1 < x 0 < tanh2K tanh2K . For the limit case N → ∞ , we can obtain the approximate values of x 0 and γ 0 , whose low-order approximations are: We can thus make a comparison between Onsager's lattice and the model with periodic-free boundary conditions. For Onsager's lattice, when the system crosses its critical temperature, (7) changes sign; however, from (51) we see that, for the model with periodic-free boundary conditions, when T ≥ T c , γ 0 ≈ 2(K − K). Once the system crosses its critical temperature T c , γ 0 becomes exponentially smaller and then vanishes rapidly as N → ∞ . This property of γ 0 plays a key role for the correlation function σ 1, 1 σ 1, N .
3. An Exact Expression of σ 1, 1 σ 1, N on a Finite Lattice In Reference [12], all roots of the equation F ± (τ) = 0 are obtained by correcting to the e −LC 0 order of magnitude (C 0 is a positive constant). These approximate roots can lead to the exact expression of σ l, 1 σ l, N in the thermodynamic limit, since lim L→∞ e −LC 0 = 0, but cannot lead to the exact expression of σ l, 1 σ l, N on a finite lattice. Hence, the expression of σ l, 1 σ l, N presented in Reference [12] is only an approximate result.
Substituting the above result and (41) into (54), we obtain the exact expression of σ 1, 1 σ 1, N of the model on a finite lattice: Although the whole calculation process is complex, the final result (60) is simple.
As for γ 0 , when T ≥ T c , from (51) we see that (61) still holds for γ 0 ; however, when T < T c , from the last expression in (51) we see that maybe lim