[
Abstract
A new approximate scheme, DSUB, is described for the coupled cluster method. We then apply it to two wellstudied (spin1/2 Heisenberg antiferromagnet) spinlattice models, namely: the and the models on the square lattice in two dimensions. Results are obtained in each case for the groundstate energy, the sublattice magnetization and the quantum critical point. They are in good agreement with those from such alternative methods as spinwave theory, series expansions, quantum Monte Carlo methods and those from the CCM using the LSUB scheme.
Key words: Coupled cluster method; quantum antiferromagnet
PACS: 75.10.Jm, 75.30.Gw, 75.40.s, 75.50.Ee
The DSUB Approximation Scheme for the Coupled Cluster Method]The DSUB Approximation Scheme for the Coupled Cluster Method and Applications to Quantum Magnets R.F. Bishop et al.]R.F. Bishop, P.H.Y. Li, J. Schulenburg
[
[

School of Physics and Astronomy, Schuster Building, The University of Manchester, Manchester, M13 9PL, UK

Universitätsrechenzentrum, Universität Magdeburg, P.O. Box 4120, 39016 Magdeburg, Germany
1 Introduction
The coupled cluster method (CCM) is a universal microscopic technique of quantum manybody theory [1, 2, 3, 4, 5, 6, 7, 8, 9]. It has been applied successfully to many physical systems including:
A special characteristic of the CCM is that it deals with infinite systems from the outset, and hence one never needs to take the limit explicitly in the number, , of interacting particles or the number of sites in the lattice. However, approximations in the inherent cluster expansions for the correlation operators are required. Several efficient and systematic approximation schemes for the CCM have been specifically developed by us previously for use with lattice systems [16, 28, 29, 30]. Up till now the most favoured and most successful CCM approximation schemes for lattices have been the socalled LSUB and SUB– schemes that we describe more fully below in section 3. Although the LSUB scheme, in particular, has been shown to be highly successful in practice for a wide variety of both frustrated and unfrustrated spinlattice systems, a disadvantage of this scheme is that the number of spin configurations generally rises very rapidly with the truncation index . This motivates us to develop alternative schemes which satisfy one or both of the following two criteria:

that we are able to calculate more levels of approximation within the scheme, and hence have available more data points for the necessary extrapolations for calculated physical quantities to the exact limit where all configurations are retained; and

that one can capture the physically most important configurations at relatively low orders, so that the quantities of interest converge more rapidly with the truncation index.
A main aim of our work is thus to provide users of the CCM with more choices of approximation scheme. In this paper we outline the formal aspects, of and present preliminary results for some benchmark models for, a new CCM approximation scheme, called the DSUB scheme, for use with systems described on a spatial lattice. The rest of the paper is organized as follows. In section 2 we introduce the CCM formalism in general. In section 3 we then discuss the specific application of the CCM to spinlattice systems. We consider in section 4 both the existing truncation schemes and introduce the alternative new DSUB scheme. In order to evaluate the accuracy of the new approximation scheme we apply it to two very wellstudied antiferromagnetic spinlattice models [18, 20], namely: the and the models on the square lattice in two dimensions (2D). Both models have a quantum phase transition. They have also been very successfully investigated by the CCM within the LSUB scheme for the groundstate (gs) energy and the gs order parameter which is, in our case, the sublattice magnetization. All techniques applied to lattice spin systems need to be extrapolated in terms of some appropriate parameters. For exact diagonalization and quantum Monte Carlo methods, this is the lattice size . As noted above, one good aspect of the CCM is that we may work in the limit of infinite lattice size () from the outset. By contrast, the extrapolation for the CCM is in terms of some truncation index , where in the limit we retain all configurations. The CCM extrapolations that have been used up till now in the trunction index, e.g., for the LSUB approximation, are heuristic schemes, but we have considerable prior experience [20, 29, 31, 32, 33] in using and refining them, as described in section 5. The new DSUB approximation scheme is then applied to the squarelattice spin antiferromagnetic model in section 6 and the corresondingly model in section 7, respectively. The results for both models are compared critically with those from the corresponding LSUB scheme and also with the results from other methods. Finally, our conclusions are given in section 8 where we reiterate a brief summary of the results.
2 The CCM Formalism
This section briefly describes the CCM formalism (and see e.g., Refs. [8, 9] for further details). A first step in any CCM application is to choose a normzalized model (or reference) state that can act as a cyclic vector with respect to a complete set of mutually commuting multiconfigurational creation operators . The index here is a setindex that labels the manyparticle configuration created in the state . The exact ket and bra gs energy eigenstates and , of the manybody system are then parametrised in the CCM form as: Ketstate Brastate
(1) 
(2) 
(3) 
where we have defined . The requirements on the multiconfigurational creation operators are that any manyparticle state can be written exactly and uniquely as a linear combination of the states , which hence fulfill the completeness relation
(4) 
together with the conditions,
(5) 
(6) 
Approximations are necessary in practice to restrict the label set I to some finite (e.g., LSUB) or infinite (e.g., SUB) subset, as described more fully below. The correlation operator is a linkedcluster operator and is decomposed in terms of a complete set of creation operators . It creates excitations on the model state by acting on it to produce correlated cluster states. Although the manifest Hermiticity, , is lost, the normalization conditions are preserved. The CCM Schrödinger equations (3) are thus writtern as
(7) 
and its equivalent similaritytransformed form becomes
(8) 
We note that while the parametrizations of equations (1) and (2) are not manifestly Hermitian conjugate, they do preserve the important HellmannFeynman theorem at all levels of approximations (viz., when the complete set of manyparticle configurations is truncated) [9]. Furthermore, the amplitudes () form canonically conjugate pairs in a timedependent version of the CCM, by contrast with the pairs () coming from a manifestly Hermitianconjugate representation for , which are not canonically conjugate to one another.
The static gs CCM correlation operators and contain the real number correlation coefficients and that need to be calculated. Clearly, once the coefficients are known, all other gs properties of the manybody system can be derived from them. For example, the gs expectation value of an arbitrary operator can be expressed as
(9) 
To find the gs correlation coefficients we simply insert the parametrization of equation (2) into the Schrödinger equations (8), and project onto the complete sets of states and , respectively,
(10) 
(11) 
Equation (11) may also easily be rewritten, by premultiplying the ketstate equation (8) with the state and using the commutation relation (6), in the form
(12) 
Equations (10)–(12) may be equivalently derived by requiring that the gs energy expectation value, , is minimized with respect to the entire set . In practice we thus need to solve equations (10) and (12) for the set . We note that equations (9) and (10) show that the gs energy at the stationary point has the simple form
(13) 
which also follows immediately from the ketstate equation (8) by projecting it onto the state . It is important to note that this (bi)variational formulation does not necessarily lead to an upper bound for when the summations over the index set for and in equation (2) are truncated, due to the lack of manifest Hermiticity when such approximations are made. Nevertheless, as we have pointed out above, one can prove [9] that the important HellmannFeynman theorem is preserved in all such approximations.
We note that equation (10) represents a coupled set of multinomial equations for the number correlation coefficients . The nested commutator expansion of the similaritytransformed Hamiltonian,
(14) 
and the fact that all of the individual components of in the decomposition of equation (2) commute with one another by construction [and see equation (6)], together imply that each element of in equation (2) is linked directly to the Hamiltonian in each of the terms in equation (14). Thus, each of the coupled equations (10) is of Goldstone linkedcluster type. In turn this guarantees that all extensive variables, such as the energy, scale linearly with particle number . Thus, at any level of approximation obtained by truncation in the summations on the index in the parametrizations of equation (2), we may (and, in practice, do) work from the outset in the limit of an infinite system.
Furthermore, each of the seemingly infiniteorder (in ) linkedcluster equations (10) will actually be of finite length when expanded using equation (14), since the otherwise infinite series of equation (14) will actually terminate at a finite order, provided only (as is usually the case, including those for the Hamiltonians considered in this paper) that each term in the Hamiltonian contains a finite number of singleparticle destruction operators defined with respect to the reference (or generalized vacuum) state . Hence, the CCM parametrization naturally leads to a workable scheme that can be implemented computationally in an efficient manner to evaluate the set of configuration coefficients by solving the coupled sets of equations (10) and (12), once we have devised practical and systematic truncation hierarchies for limiting the set of multiconfigurational setindices to some suitable finite or infinite subset.
3 CCM for SpinLattice Systems
We will discuss various practical CCM truncation schemes that fulfill the criteria of being systematically improvable in some suitable truncation index , and that can be extrapolated accurately in practice to the exact, , limit for calculated physical quantities. Before doing so, however, we first describe how the general CCM formalism described in section 2 is implemented for spinlattice problems in practice. As is the case for any application of the CCM to a general quantum manybody system, the first step is to choose a suitable reference state in which the the state of the spin (viz., in practice, its projection onto a specific quantization axis in spin space) on every lattice site is characterized. Clearly, the choice of will depend on both the system being studied and, more importantly, which of its possible phases is being considered. We describe examples of such choices later for the particular models that we utilize here as test cases for our new truncation scheme.
Whatever the choice of we note firstly that is very convenient, in order to set up as universal a methodology as possible, to treat the spins on every lattice site in an arbitrarily given model state as being equivalent. A suitably simple way of so doing is to introduce a different local quantization axis and a correspondingly different set of spin coordinates on each lattice site , so that all spins, whatever their original orientation in in a global spincoordinate system, align along the same direction (which, to be specific, we henceforth choose as the negative direction) in these local spincoordinate frames. This can always be done in practice by defining a suitable rotation in spin space of the global spin coordinates at each lattice site . Such rotations are canonical transformations that leave unchanged the fundamental spin commutation relations,
(15) 
(16) 
among the usual SU(2) spin operators () on lattice site . Each spin has a total spin quantum number, , where is the SU(2) Casimir operator. For the models considered here, , .
In the local spin frames defined above the configuration indices simply become a set of lattice site indices, ). The corresponding generalized multiconfigurational creation operators thus become simple products of single spinraising operators, , and, for example, the ketstate CCM correlation operator is expressed as
(17) 
and is similarly defined in terms of the spinlowering operators . Since the operator acts on the state , in which all spins point along the negative axis in the local spincoordinate frames, every lattice site in equation (17) can be repeated up to no more than 2 times in each term where it is allowed, since a spin has only () possible projections along the quantization axis.
In practice the allowed configurations are often further constrained by symmetries in the problem and by conservation laws. An example of the latter is provided by the model considered below in section 6, for which we can easily show that the total component of spin, , in the original global spin coordinates, is a good quantum number since in this case. Finally, for the quasiclassical magnetically ordered states that we calculate here for both models in sections 6 and 7, the order parameter is the sublattice magnetization, , which is given within the local spin coordinates defined above as
(18) 
After the local spin axes have been chosen, the model state thus has all spins pointing downwards (i.e., in the negative direction, where is the quantization axis),
(19) 
where in the usual notation for single spin states.
The similaritytransformed Hamiltonian , and all of the corresponding matrix elements in equations (9)–(13) and equation (18), for example, may then be evaluated in the local spin coordinate frames by using the nested commutator expansion of equation (14), the commutator relations of equation (15), and the simple universal relations
(20) 
(21) 
that hold at all lattice sites in the local spin frames.
4 Approximation schemes
The CCM formalism is exact if all manybody configurations are included in the and operators. In practice, it is necessary to use approximation schemes to truncate the correlation operators.
4.1 Common previous truncation schemes
The main approximation scheme used to date for continuous systems is the socalled SUB scheme described below. For systems defined on a regular periodic spatial lattice, we have a further set of approximation schemes which are based on the discrete nature of the lattice, such as the SUB– and LSUB schemes described below. The various schemes and their definitions for spinlattice systems are:

the SUB scheme, in which only the correlations involving or fewer spinraising operators for are retained, but with no further restrictions on the spatial separations of the spins involved in the configurations;

the SUB– scheme wich includes only the subset of all spinflip configurations in the SUB scheme that are defined over all lattice animals of size , where a lattice animal is defined as a set of contiguous lattice sites, each of which is nearestneighbour to at least one other in the set; and

the LSUB scheme which includes all possible multispinflip configurations defined over all lattice animals of size . The LSUB scheme is thus equivalent to the SUB– scheme for for particles of spin quantum number . For example, for spin1/2 systems, for which no more than one spinraising operator, , can be applied at each site , LSUB SUB–.
4.2 The new DSUB scheme
Our new DSUB scheme is now defined to include in the correlation operator all possible configurations of spins involving spinraising operators where the maximum length or distance of any two spins apart is defined by , where is a vector joining sites on the lattice and the index labels lattice vectors in order of size. Hence DSUB includes only nearestneighbours, etc.
Table 1 shows how progresses in terms of and (which are the sides of the lattice in the and directions) for the case of a 2D square lattice with sides parallel to the and axes.
DSUB  

DSUB  1  0  1 
DSUB  1  1  
DSUB  2  0  2 
DSUB  1  2  
DSUB  2  2  
DSUB  3  0  3 
DSUB  1  3 
Similar tables can be constructed for an arbitrary regular lattice in any number of dimensions. Table 1 shows, for example, that the DSUB5 approximation on a 2D square lattice involves all clusters of spins (and their associated spinraising operators) for which the maximum distance between any two spins is . Clearly the DSUB scheme orders the multispin configurations in terms, roughly, of their compactness, whereas the LSUB scheme orders them according to the overall size of the lattice animals (or polyominos), defined as the number of contiguous lattice sites involved.
5 Extrapolation schemes
Any of the above truncated approximations clearly becomes exact when all possible multispin cluster configurations are retained, i.e., in the limit as and/or . We have considerable experience, for example, with the appropriate extrapolations for the LSUB scheme [20, 29, 31, 32, 33], that shows that the gs energy behaves in the large limit as a power series in , whereas the order parameter behaves as a power series in (at least for relatively unfrustrated systems). Initial experience with the new DSUB scheme shows, perhaps not unsurprisingly, that in the corresponding large limit the gs energy and order parameter behave as power series in and , respectively, as we show in more detail (below) for the two examples of the spin1/2 and models on the 2D square lattice. It is clear on physical grounds that the index should provide a better extrapolation variable than the index itself for the DSUB scheme, and so it turns out in practice. For the present, where we are interested primarily in a preliminary investigation of the power and accuracy of the DSUB scheme, we limit ourselves to retaining only the leading terms in the powerseries expansions,
(22) 
(23) 
Further subleading terms in each of the power series can easily be retained later should it prove useful to increase the accuracy of the extrapolations.
5.1 Three fundamental rules for the selection and extrapolation of the CCM raw data
We list below three fundamental rules as guidelines for the selection and extrapolation of the CCM raw data, using any approximation scheme:

RULE 1: In order to fit well to any fitting formula that contains unknown parameters, one should have at least () data points. This rule takes precedence over all other rules, and is vital to obtain a robust and stable fit.

RULE 2: Avoid using the lowest data points (e.g., LSUB2, SUB22, DSUB1, etc.) wherever possible, since these points are rather far from the large limit, unless it is necessary to do so to avoid breaking RULE 1, e.g., when only data points are available.

RULE 3: If RULE 2 has been broken (e.g., by including LSUB2 or SUB– data points), then do some other careful consistency checks on the robustness and accuracy of the results.
6 The spin antiferromagnetic model on the square lattice
In this section, we shall consider the spin1/2 model on the infinite square lattice. The Hamiltonian of the model, in global spin coordinates, is written as
(24) 
where the sum on runs over all nearestneighbour pairs of sites on the lattice and counts each pair only once. Since the square lattice is bipartite, we consider to be even, so that each sublattice contains spins, and we consider only the case where . The Néel state is the ground state (GS) in the trivial Ising limit , and a phase transition occurs at (or near to) . Indeed, the classical GS demonstrates perfect Néel order in the direction for , and a similar perfectly ordered  planar Néel phase for . For the classical GS is a ferromagnet.
The case is equivalent to the isotropic Heisenberg model, whereas is equivalent to the isotropic version of the model considered in section 7 below. The component of total spin, , is a good quantum number as it commutes with the Hamiltonian of equation (24). Thus one may readily check that . Our interest here is in those values of for which the GS is an antiferromagnet.
The CCM treatment of any spin system is started by choosing an appropriate model state (for a particular regime), so that a linear combinations of products of spinraising operators can be applied to this state and all possible spin configurations are determined. There is never a unique choice of model state . Our choice should clearly be guided by any physical insight that we can bring to bear on the system or, more specifically, to that particular phase of it that is under consideration. In the absence of any other insight into the quantum manybody system it is common to be guided by the behaviour of the corresponding classical system (i.e., equivalently, the system when the spin quantum number ). The model under consideration provides just such an illustrative example. Thus, for the classical Hamiltonian of equation (24) on the 2D square lattice (and, indeed, on any bipartite lattice) is minimized by a perfectly antiferromagnetically Néelordered state in the spin direction. However, the classical gs energy is minimized by a Néelordered state with spins pointing along any direction in the spin  plane (say, along the spin direction) for . Either of these states could be used as a CCM model state and both are likely to be of value in different regimes of appropriate to the particular quantum phases that mimic the corresponding classical phases. For present illustrative purposes we restrict ourselves to the aligned Néel state as our choice for , written schematically as
where and in the usual notation. Such a state is, clearly, likely to be a good startingpoint for all , down to the expected phase transition at from a aligned Néel phase to an  planar Néel phase.
As indicated in section 3 it is now convenient to perform a rotation of the axes for the uppointing spins (i.e., those on the sublattice with spins in the positive direction) by about the spin axis, so that takes the form given by equation (19). Under this rotation, the spin operators on the original up sublattice are transformed as
(25) 
The Hamiltonian of equation (24) may thus be rewritten in these local spin coordinate axes as
(26) 
As in any application of the CCM to spinlattice systems, we now include in our approximations at any given order only those fundamental configurations that are distinct under the point and space group symmetries of both the lattice and the Hamiltonian. The number, , of such fundamental configurations at any level of approximation may be further restricted whenever additional conservation laws come into play. For example, in our present case, the Hamiltonian of equation (24) commutes with the total uniform magnetization, , in the global spin coordinates, where runs over all lattice sites. The GS is known to lie in the subspace, and hence we exclude configurations with an odd number of spins or with unequal numbers of spins on the two equivalent sublattices of the bipartite square lattice. We show in figure 1
the fundamental configurations that are accordingly allowed for the DSUB approximations for this spin1/2 model on the 2D square lattice, with . We see that, for example, at the DSUB4 level of approximation. We also see that for both the DSUB2 and DSUB3 approximations, since the conservation law does not permit any additional configurations of spins with a maximum distance , apart from those already included in the DSUB2 approximation.
6.1 Groundstate energy and sublattice magnetization
The DSUB approximations can readily be implemented for the present spin1/2 model on the 2D square lattice for all values with reasonably modest computing power. By comparison, the LSUB scheme can be implemented with comparable computing resources for all values . Numerical results for the gs energy per spin and sublattice magnetization are shown in table 2
Method  Max. No.  
of spins  
DSUB1=LSUB2  1  1  0.64833  0.421  0  2 
DSUB2  2  0.65311  0.410  0.258  4  
DSUB3  2  2  0.65311  0.410  0.258  4 
DSUB4  12  0.66258  0.385  0.392  6  
DSUB5  20  0.66307  0.382  0.479  8  
DSUB6  3  43  0.66511  0.375  0.506  8 
DSUB7  135  0.66589  0.371  0.629  12  
DSUB8  831  0.66704  0.363  0.614  14  
DSUB9  4  1225  0.66701  0.363  0.654  14 
DSUB10  6874  0.66774  0.357  0.637  16  
DSUB11  14084  0.66785  0.356    16  
LSUB8  1287  0.66817  0.352  8  
LSUB10  29605  0.66870  0.345  10  
Extrapolations  
Based on  
DSUB  0.67082  0.308  1.009  
DSUB  0.67122  0.311  1.059  
DSUB  0.66978  0.319  
DSUB  0.67177  0.315  1.025  
DSUB  0.66967  0.325  0.979  
LSUB [34]  0.67029  0.305  
LSUB [34, 35]  0.66966  0.310  
LSUB  0.66962  0.308  
SWT [36]  0.66999  0.307  
SE [38]  0.6693  0.307  
ED [37]  0.6700  0.3173  
QMC [39]  0.669437(5)  0.3070(3) 
NOTES:
The magnetization point of inflexion for DSUB11 is not available since we only calculated at in this approximation.
The magnetization points of inflexion for the odd DSUB levels are extrapolated using .
The magnetization points of inflexion for the whole series of DSUB data are extrapolated as indicated, but without .
at the isotropic point at various levels of approximation, and corresponding results for the same quantities are displayed graphically in figures 2 and 3 as functions of the anisotropy parameter .
We also show in table 2 for the isotropic Heisenberg Hamiltonian () the results for the gs energy and sublattice magnetization using the leading (linear) extrapolation schemes of equations (22) and (23) respectively of the DSUB data, employing various subsets of results. Comparison is also made with corresponding LSUB extrapolation schemes for the same model [34, 35]. The results are generally observed to be in excellent agreement with each other, even though the DSUB extrapolations have employed the simple leading (linear) fits of equations (22) and (23), whereas the corresponding LSUB results shown [34, 35] have been obtained from the potentially more accurate quadratic fits , , to the LSUB data, in which the nextorder (quadratic) corrections have also been included in the relevant expansion parameters, and , respectively. Excellent agreement of all the CCM extrapolations is also obtained with the results from the best of the alternative methods for this model, including thirdorder spinwave theory (SWT) [36], linkedcluster series expansion techniques [38], and the extrapolations to infinite lattice size () from the exact diagonalization (ED) of small lattices [37], and quantum Monte Carlo (QMC) calculations for larger lattices [40].
We note that it has been observed and well documented in the past (and see, e.g., Ref. [34]) that the CCM LSUB results for this model (and many others) for both the gs energy and the sublattice magnetization show a distinct period2 ‘‘staggering’’ effect with index , according to whether is even or odd. As a consequency the LSUB data for both and converge differently for the even and the odd sequences, similar to what is observed very frequently in perturbation theory in corresponding even and odd orders [40]. As a rule, therefore, the LSUB data are generally extrapolated separately for even and for odd values of , since the staggering makes extrapolations using both odd and even values together extremely difficult. We show in figure 4
our DSUB results for the gs energy per spin and the sublattice magnetization plotted against and , respectively. The higher values clearly cluster well in both cases on straight lines, thereby justifying a posteriori our heuristic extrapolation fits of equations (22) and (23). Just as in the LSUB case a slight ‘‘evenodd staggering’’ effect is observed in the DSUB data (perhaps more so for the energy than for the sublattice magnetization), although it is less pronounced than for the corresponding LSUB data [34].
6.2 Termination or critical points
Before discussing our DSUB results further for this model we note that the comparable LSUB solutions actually terminate at a critical value , which itself depends on the truncation index [30]. Such LSUB termination points are very common for many spinlattice systems and are very well documented and understood (and see, e.g., Ref. [30]). In all such cases a termination point always arises due to the solution of the CCM equations becoming complex at this point, beyond which there exist two branches of entirely unphysical complex conjugate solutions [30]. In the region where the solution reflecting the true physical solution is real there actually also exists another (unstable) real solution. However, only the (shown) upper branch of these two solutions reflects the true (stable) physical GS, whereas the lower branch does not. The physical branch is usually easily identified in practice as the one which becomes exact in some known (e.g., perturbative) limit. This physical branch then meets the corresponding unphysical branch at some termination point beyond which no real solutions exist. The LSUB termination points are themselves also reflections of the quantum phase transitions in the real system, and may be used to estimate the position of the phase boundary [30].
It is interesting and intriguing to note that when the DSUB approximations are applied to the model, they do not terminate as do the corresponding LSUB approximations. We have no real explanation for this rather striking difference in behaviour for two apparently similar schemes applied to the same model. However, it is still possible to use our DSUB data to extract an estimate for the physical phase transition point at which the aligned Néel phase terminates. As has been justified and utilized elsewhere [18], a point of inflexion at in the sublattice magnetization as a function of clearly indicates the onset of an instability in the system. Such inflexion points occur for all DSUB approximations, as indicated in table 2 and figure 3. The DSUB approximations are thus expected to be unphysical for , and we hence show the corresponding results for the gs energy per spin in figure 2 only for values . Heuristically, we find that the magnetization inflexion points scale linearly with to leading order, and the extrapolated results shown in table 2 have been performed with the leading (linear) fit, , commensurate with the corresponding linear fits in and for the gs energy per spin and sublattice magnetization of equations (22) and (23), respectively. All of the various extrapolations shown in table 2 for in the limit are in good agreement with one another, thereby again demonstrating the robust quality of the heuristic extrapolation scheme. Furthermore, they are also in excellent agreement with the expected phase transition point at between two quasiclassical Néelordered phases aligned along the spin axis (for ) and in some arbitrary direction in the spin plane (for ).
Although we do not do so here, the  planar Néel phase could itself also easily be investigated by another CCM DSUB series of calculations based on a model state with perfect Néel ordering in, say, the direction.
Summarizing our results so far, we observe that the DSUB scheme has at, least partially, fulfilled the expectations placed upon it for the present model. Accordingly, we now apply it to the second test model of the spin1/2 model on the 2D square lattice.
7 The spin1/2 model on the square lattice
The Hamiltonian of the model [18] in global spin coordinates, is written as
(27) 
where the sum on again runs over all nearestneighbour pairs of lattice sites and counts each pair only once. We again consider the case of spin1/2 particles on each site of an infinite square lattice.
For the classical model described by equation (27), it is clear that the GS is a Néel state in the direction for and a Néel state in the direction for . Hence, since we only consider the case , we choose as our CCM model state for the quantum model a Néel state aligned along the direction, written schematically as,
Clearly the case is readily obtained from the case by interchange of the  and axes.
Once again we now perform our usual rotation of the spin axes on each lattice site so that takes the form given by equation (19) in the rotated local spin coordinate frame. Thus, for the spins on the sublattice where they point in the negative direction in the global spin axes (i.e., the leftpointing spins) we perform a rotation of the spin axes by about the spin axis. Similarly, for the spins on the other sublattice where they point in the positive direction in the global spin axes (i.e., the rightpointing spins) we perform a rotation of the spin axes by about the spin axis. Under these rotations the spin operators are transformed as
(28a)  
(28b) 
The Hamiltonian of equation (27) may thus be rewritten in the local spin coordinate axes defined above as
(29) 
As before, we now have to evaluate the fundamental configurations that are retained in the CCM correlation operators and at each DSUB level of approximation. Although the point and space group symmeries of the square lattice (common to both the and models considered here) and the two Hamiltonians of equations (26) and (29) are identical, the numbers of fundamental configurations for a given DSUB level are now larger (except for the case ) for the model than for the model, since the uniform magnetization is no longer a good quantum number for the model, . Nevertheless, we note from the form of equation (29), in which the spinraising and spinlowering operators appear only in combinations that either raise or lower the number of spin flips by two (viz., the and combinations, respectively) or leave them unchanged (viz., the and combinations), it is only necessary for the GS to consider fundamental configurations that contain an even number of spins. Thus, the main difference for the model over the model is that we must now also consider fundamental configurations in which we drop the restriction for the former case of having an equal number of spins on the two equivalent sublattices of the bipartite square lattice that was appropriate for the latter case. We show in figure 5 the fundamental configurations that are allowed for the spin1/2 model on the square lattice for the DSUB approximation with , and we invite the reader to compare with the corresponding fundamental configurations for the spin1/2 model on the same square lattice shown in figure 1.
The corresponding numbers of fundamental configurations for the model are also shown in table 3 for the higher DSUB approximations with for which we present results below.
7.1 Groundstate energy and sublattice magnetization
We present results for the spin1/2 model on the square lattice in the CCM DSUB approximations for all values that can be easily computed with very modest computing power. Comparable computing power enables the corresponding LSUB scheme to be implemented for all . Numerical results for the gs energy per spin and sublattice magnetization are shown in table 3 at the isotropic point at at various levels of approximation, and corresponding results for the same gs quantities are shown graphically in figures 6 and 7 as functions of the anisotropy parameter .
Method  Max. No.  
of spins  
DSUB1=LSUB2  1  1  0.54031  0.475  2  
DSUB2  3  0.54425  0.467  4  
DSUB3  2  6  0.54544  0.464  4  
DSUB4  21  0.54724  0.458  0.253  6  
DSUB5  44  0.54747  0.456  0.205  8  
DSUB6  3  78  0.54774  0.455  0.181  8 
DSUB7  388  0.54811  0.453  0.135  12  
DSUB8  1948  0.54829  0.451  0.107  14  
DSUB9  4  3315  0.54833  0.451  0.099  14 
LSUB6  131  0.54833  0.451  0.073  6  
LSUB8  2793  0.54862  0.447  0.04  8  
Extrapolations  
Based on  
DSUB  0.54879  0.442  0.036  
DSUB  0.54923  0.437  0.011  
DSUB  0.54884  0.441  0.029  
LSUB [18]  0.54892  0.435  0.00  
SE [41]  0.54883  0.43548  0.0  
QMC [42]  0.548824(2)  0.437(2) 
We also show in table 3 for the isotropic Hamiltonian () the results for the gs energy and sublattice magnetization using the leading (linear) extrapolation schemes of equations (22) and (23) respectively of the DSUB data, employing various subsets of our results, as for the model considered previously. We also compare in table 3 the present results with the corresponding CCM LSUB results [18] for the same model. All of the CCM results are clearly in excellent agreement both with one another and with the results of best of the alternative methods available for this model, including the linkedcluster series expansion (SE) techniques [41] and a quantum Monte Carlo (QMC) method [42].
We again show in figure 8 our DSUB results for the present model for the gs energy per spin and the sublattice magnetization, plotted respectively against and .
As previously for the model, the higher values cluster well on straight lines in both cases, thereby justifying once more our heursitic choice of extrapolation fits indicated in equations (22) and (23). Figures 8(a) and 8(b) again show an ‘‘evenodd’’ staggering effect in the termination index for the DSUB data, which is perhaps slightly more pronounced than that for the model shown in figures 4(a) and 4(b). For this reason we have again shown separate extrapolations of our DSUB results in table 3 for the even data and the odd data, as well as results using all (higher) values of .
7.2 Termination or critical points
It is interesting to note that for the present model the CCM DSUB solutions (with our choice of model state as a Néel state in the direction) now do physically terminate for all values of the truncation index at a critical value , exactly as commonly occurs (as for the present model) for the LSUB calculations, as we explained above in section 6.2. Why such DSUB terminations occur for the model but not for the previous model is not obvious to us. The corresponding termination points, , at various DSUB and LSUB levels of approximation are shown in table 3. It has been shown previously [29] that scales well with for the LSUB data, and the LSUB result [18] shown in table 3 was obtained by a leading (linear) fit, . We find heuristically that the best large asymptotic behaviour of the DSUB data for is against as the scaling parameter. Accordingly, the DSUB values for in table 3 are obtained with the leading (linear) fit, . We see that both the LSUB and DSUB results for agree very well with the value that is known to be the correct value for the phase transition in the onedimensional spin1/2 chain from the known exact solution [43], and which is believed also to be the phase transition point for higher dimensions, including the present 2D square lattice, on symmetry grounds.
8 Conclusions
From the two nontrivial benchmark spinlattice problems that we have investigated here, it is clear that the new DSUB approximation scheme works well for calculating their gs properties and phase boundaries. We have utilized here only the simplest leadingorder extrapolation schemes in the pertinent scaling variables, and have shown that these may be chosen, for example, as for the gs energy and for the order parameter. Clearly, in general, the results can be further improved by keeping higherorder terms in these asymptotic expansions (i.e., by retaining higher powers in the polynomial scaling expansions) although more data points may then be needed, especially in cases where the ‘‘evenodd’’ staggering effect is pronounced, as for model presented here. For further use of the scheme for more complex lattice models (e.g., those exhibiting geometric or dynamic frustration) it will be necessary to revisit the validity of these expansions, but a great deal of previous experience in such cases for the LSUB scheme will provide good guidance.
On the basis of the test results presented here, the DSUB scheme clearly fulfills the first of our two main criteria for introducing it, viz., that the number of fundamental configurations, , increases less rapidly with truncation index than for the corresponding LSUB series of approximations. At the same time our second criterion of capturing the physically most important configurations at relatively low levels of approximation also seems to be fulfilled, according to our experience with the convergence of the DSUB sequences for observable quantities. At the very least we now have two schemes (LSUB and DSUB) available to us for future investigations, each of which has its own merits, and which thus allows us more freedom in applications of the CCM to other spinlattice models in future.
The one slight drawback in the scheme which mitigates against our goal of obtaining more DSUB data points, for the same computing power than for the LSUB scheme applied to the same system, and that hence can be used together to attain more accuracy in the extrapolations, is the slight ‘‘evenodd’’ staggering in the data that is observed in the DSUB results, albeit that it is somewhat reduced from the similar stagerring in the corresponding LSUB results. We have some ideas on how the DSUB scheme might itself be modified to reduce this staggering and we hope to report results of these further investigations in a future paper.
References
 1. Coester F., Nucl. Phys., 1958, 7, 421.
 2. Čížek J., J. Chem. Phys., 1966, 45, 4256.
 3. Paldus J., Čížek J., Shavitt I., Phys. Rev. A, 1972, 5, 50.
 4. Kümmel H., Lührmann K.H., Zabolitzky J.G., Phys. Rep., 1978, 36C, 1.
 5. Arponen J.S., Ann. Phys. (N.Y.), 1983, 151, 311.
 6. Arponen J.S., Bishop R.F., Pajanne E., Phys. Rev. A, 1987, 36, 2539.
 7. Bartlett R.J., J. Phys. Chem., 1989, 93, 1697.
 8. Bishop R.F., Theor. Chim. Acta, 1991, 80, 95.
 9. Bishop R.F., in Microscopic Quantum ManyBody Theories and Their Applications, eds. Navarro J., Polls A., Lecture Notes in Physics, 1998, 510, SpringerVerlag, Berlin, p.1.
 10. Bishop R.F., Lührmann K.H., Phys. Rev. B, 1978, 17, 3757.
 11. Bishop R.F., Lührmann K.H., Phys. Rev. B, 1982, 26, 5523.
 12. Day B.D., Phys. Rev. Lett., 1981, 47, 226.
 13. Day B.D., Zabolitzky J.G., Nucl. Phys., 1981, A336, 221.
 14. Bartlett R.J., Ann. Rev. Phys. Chem., 1981 32, 359.
 15. Roger M., Hetherington J.H., Phys. Rev. B, 1990, 41, 200.
 16. Bishop R.F., Parkinson J.B., Xian Y., Phys. Rev. B, 1991, 44, 9425.
 17. Bursill R., Gehring G.A., Farnell D.J.J., Parkinson J.B., Xiang T., Zeng C., J. Phys.: Condens. Matter, 1995, 7, 8605.
 18. Farnell D.J.J., Krüger S.E., Parkinson J.B., J. Phys.: Condens. Matter, 1997, 9, 7601.
 19. Bishop R.F., Farnell D.J.J., Parkinson J.B., Phys. Rev. B, 1998, 58, 6394.
 20. Bishop R.F., Farnell D.J.J., Krüger S.E., Parkinson J.B., Richter J., Zeng C., J. Phys.: Condens. Matter, 2000, 12, 6887.
 21. Farnell D.J.J., Gernoth K.A., Bishop R.F., Phys. Rev. B, 2001, 64, 172409.
 22. Farnell D.J.J., Bishop R.F., Gernoth K.A., J. Stat. Phys., 2002, 108, 401.
 23. Bishop R.F., Li P.H.Y., Darradi R., Richter J., J. Phys.: Condens. Matter, 2008, 20, 255251.
 24. Bishop R.F., Li P.H.Y., Darradi R., Richter J., Europhys. Lett., 2008 83, 47004.
 25. Bishop R.F., Li P.H.Y., Darradi R., Schulenburg J., Richter J., Phys. Rev. B, 2008, 78, 054412.
 26. Bishop R.F., Li P.H.Y., Darradi R., Richter J., Campbell C.E., J. Phys.: Condens. Matter, 2008, 20, 415213.
 27. Bishop R.F., Li P.H.Y., Farnell D.J.J., Campbell C.E., Phys. Rev. B, 2009, 79 (in press).
 28. Bishop R.F., Parkinson J.B., Xian Yang, Phys. Rev. B, 1991, 43, 13782.
 29. Bishop R.F., Hale R.G., Xian Y., Phys. Rev. Lett., 1994, 73, 3157.
 30. Farnell D.J.J., Bishop R.F., in Quantum Magnetism, eds. Schollwöck U., Richter J., Farnell D.J.J., Bishop R.F., Lecture Notes in Physics, 2004, 645, SpringerVerlag, Berlin, p.307.
 31. Zeng C., Farnell D.J.J., Bishop R.F., J. Stat. Phys., 1998, 90, 327.
 32. Krüger S.E., Richter J., Schulenburg J., Farnell D.J.J., Bishop R.F., Phys. Rev. B, 2000, 61, 14607.
 33. Schmalfuß D., Darradi R., Richter J., Schulenburg J., Ihle D., Phys. Rev. Lett., 2006, 97, 157201.
 34. Farnell D.J.J., Bishop R.F., Int. J. Mod. Phys. B, 2008, 22, 3369.
 35. Richter J., Darradi R., Zinke R., Bishop R.F., Int. J. Mod. Phys. B, 2007, 21, 2273.
 36. Hamer C.J., Weihong Z., Arndt P., Phys. Rev. B, 1992, 46, 6276.
 37. Weihong Z., Oitmaa J., Hamer C.J., Phys. Rev. B, 1991, 43, 8321.
 38. Richter J., Schulenburg J., Honecker A., in Quantum Magnetism, eds. Schollwöck U., Richter J., Farnell D.J.J., Bishop R.F., Lecture Notes in Physics, 2004, 645, SpringerVerlag, Berlin, p.85.
 39. Sandvik A.W., Phys. Rev. B, 1997, 56, 11678.
 40. Morse P.M., Feshbach H., Methods of Theoretical Physics, Part II, 1953, McGrawHill, New York.
 41. Hamer C.J., Oitmaa J., Weihong Z., Phys. Rev. B, 1991, 43, 10789.
 42. Sandvik A.W., Hamer C.J., Phys. Rev. B, 1999, 60, 6588.
 43. Lieb E., Schultz T., Mattis D., Ann. Phys. (N.Y.), 1961, 16, 407.